File: N:\mfix\model\write_res1.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: WRITE_RES1                                             C
4     !  Purpose: write out the time-dependent restart records               C
5     !                                                                      C
6     !  Author: P. Nicoletti                               Date: 13-DEC-91  C
7     !  Reviewer: P. Nicoletti, W. Rogers, M. Syamlal      Date: 24-JAN-92  C
8     !                                                                      C
9     !  Revision Number:                                                    C
10     !  Purpose:                                                            C
11     !  Author:                                            Date: dd-mmm-yy  C
12     !  Reviewer:                                          Date: dd-mmm-yy  C
13     !                                                                      C
14     !  Literature/Document References:                                     C
15     !                                                                      C
16     !  Variables referenced: TIME, NSTEP, EP_g, P_g, P_star, RO_g, ROP_g   C
17     !                        T_g, T_s, U_g, V_g, W_g, ROP_s, U_s    C
18     !                        V_s, W_s, IJKMAX2                             C
19     !  Variables modified: None                                            C
20     !                                                                      C
21     !  Local variables: LC, N, NEXT_REC                                    C
22     !                                                                      C
23     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
24     !
25           SUBROUTINE WRITE_RES1
26     !...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
27     !...Switches: -xf
28     !
29     !-----------------------------------------------
30     !   M o d u l e s
31     !-----------------------------------------------
32           USE cdist
33           USE compar           !//
34           USE energy
35           USE fldvar
36           USE funits
37           USE geometry
38           USE machine, only: flush_res
39           USE mpi_utility      !//d pnicol : for gather
40           USE output
41           USE param
42           USE param1
43           USE physprop
44           USE run
45           USE rxns
46           USE scalars
47           USE sendrecv         !//d pnicol : for gather
48     !//d pnicol  ... not needed    USE tmp_array
49           IMPLICIT NONE
50     !-----------------------------------------------
51     !   G l o b a l   P a r a m e t e r s
52     !-----------------------------------------------
53     !-----------------------------------------------
54     !   L o c a l   P a r a m e t e r s
55     !-----------------------------------------------
56     !-----------------------------------------------
57     !   L o c a l   V a r i a b l e s
58     !-----------------------------------------------
59     !
60     !
61     !
62     !
63     !//d pnicol : allocate arrays for gather/convert_to_io_dp
64           double precision, allocatable :: array1(:)
65           double precision, allocatable :: array2(:)
66     
67     
68     !             loop counter
69           INTEGER :: LC, NN
70     !
71     !             pointer to first time-dependent record in restart file
72           INTEGER :: NEXT_REC
73     !-----------------------------------------------
74     !
75     
76     !//d pnicol
77     !      if (myPE.eq.PE_IO .and. .not.distio) then
78              allocate (array1(ijkmax2))
79              allocate (array2(ijkmax3))
80     !      else
81     !         allocate (array1(1))
82     !         allocate (array2(1))
83     !      end if
84     
85     
86           if (myPE.eq.PE_IO .or. bDist_IO) then
87              READ (UNIT_RES, REC=3) NEXT_REC
88              WRITE (UNIT_RES, REC=NEXT_REC) TIME, DT, NSTEP
89              NEXT_REC = NEXT_REC + 1
90           end if
91     !
92     !\\SP Local Send Receive - need to be moved to source later!!
93     
94           if (.not. bDist_IO) then
95     
96               call send_recv(EP_g,2)
97               call send_recv(P_g,2)
98               call send_recv(P_star,2)
99               call send_recv(RO_g,2)
100               call send_recv(ROP_g,2)
101               call send_recv(X_g,2)
102               call send_recv(T_g,2)
103               call send_recv(U_g,2)
104               call send_recv(V_g,2)
105               call send_recv(W_g,2)
106               call send_recv(ROP_S,2)
107               call send_recv(T_S,2)
108               call send_recv(U_S,2)
109               call send_recv(V_S,2)
110               call send_recv(W_S,2)
111               call send_recv(THETA_M,2)
112               call send_recv(X_S,2)
113               if(NScalar > 0)call send_recv(Scalar,2)
114               if(K_Epsilon) THEN
115                   call send_recv(K_Turb_G,2)
116                   call send_recv(E_Turb_G,2)
117               endif
118               call send_recv(GAMA_RG,2)
119               call send_recv(T_RG,2)
120               call send_recv(GAMA_RS,2)
121               call send_recv(T_RS,2)
122               if(nRR > 0)call send_recv(ReactionRates,2)
123     
124           end if
125     
126     
127           call gatherWriteRes (EP_g,array2, array1, NEXT_REC)  !//d pnicol
128     !
129           call gatherWriteRes (P_g,array2, array1, NEXT_REC)  !//d pnicol
130     !
131           call gatherWriteRes (P_star,array2, array1, NEXT_REC)  !//d pnicol
132     !
133           call gatherWriteRes (RO_g,array2, array1, NEXT_REC)  !//d pnicol
134     !
135           call gatherWriteRes (ROP_g,array2, array1, NEXT_REC)  !//d pnicol
136     !
137           call gatherWriteRes (T_g,array2, array1, NEXT_REC)  !//d pnicol
138     !
139           DO NN = 1, NMAX(0)
140                 call gatherWriteRes (X_g(:,nn),array2, array1, NEXT_REC)  !//d pnicol
141           END DO
142     !
143           call gatherWriteRes (U_g,array2, array1, NEXT_REC)  !//d pnicol
144     !
145           call gatherWriteRes (V_g,array2, array1, NEXT_REC)  !//d pnicol
146     !
147           call gatherWriteRes (W_g,array2, array1, NEXT_REC)  !//d pnicol
148     !
149           DO LC = 1, MMAX
150     !
151             call gatherWriteRes (ROP_s(:,LC),array2, array1, NEXT_REC)  !//d pnicol
152     
153             IF(ANY(SOLVE_ROs)) &
154                 call gatherWriteRes (RO_S(:,LC),array2, array1, NEXT_REC)  !//d pnicol
155     !
156             call gatherWriteRes (T_s(:,LC),array2, array1, NEXT_REC)  !//d pnicol
157     !
158             call gatherWriteRes (U_s(:,LC),array2, array1, NEXT_REC)  !//d pnicol
159     !
160             call gatherWriteRes (V_s(:,LC),array2, array1, NEXT_REC)  !//d pnicol
161     !
162             call gatherWriteRes (W_s(:,LC),array2, array1, NEXT_REC)  !//d pnicol
163     !
164             call gatherWriteRes (THETA_M(:,LC),array2, array1, NEXT_REC)  !//d pnicol
165     !
166              DO NN = 1, NMAX(LC)
167                 call gatherWriteRes (X_s(:,LC,NN),array2, array1, NEXT_REC)  !//d pnicol
168              END DO
169           END DO
170     !
171     !     Version 1.3
172     
173           DO LC = 1, NScalar
174                 call gatherWriteRes (Scalar(:,LC),array2, array1, NEXT_REC)  !//d pnicol
175           END DO
176     !
177     !     Version 1.4 -- write radiation variables in write_res1
178           call gatherWriteRes (GAMA_RG,array2, array1, NEXT_REC)  !//d pnicol
179     
180           call gatherWriteRes (T_RG,array2, array1, NEXT_REC)  !//d pnicol
181     
182           DO LC = 1, MMAX
183             call gatherWriteRes (GAMA_RS(1,LC),array2, array1, NEXT_REC)  !//d pnicol
184     
185             call gatherWriteRes (T_RS(1,LC),array2, array1, NEXT_REC)  !//d pnicol
186           ENDDO
187     
188     !
189     !     Version 1.5
190           DO LC = 1, nRR
191                 call gatherWriteRes (ReactionRates(:,LC),array2, array1, NEXT_REC)  !//d pnicol
192           END DO
193     !
194     !     Version 1.6
195     
196           if (K_epsilon) then
197                 call gatherWriteRes (K_turb_G,array2, array1, NEXT_REC)  !//d pnicol
198               call gatherWriteRes (E_turb_G,array2, array1, NEXT_REC)
199           endif
200     !---------------------------------------------------------------------
201     
202           if ( (myPE.eq.PE_IO .and. .not.bDist_IO) .or. bDist_IO) then
203                CALL FLUSH_res (UNIT_RES)
204             end if
205     
206     !      call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
207     !
208           deallocate (array1)  !//d pnicol
209           deallocate (array2)  !//d pnicol
210     !     call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
211     
212           call write_res1_netcdf
213     !
214           RETURN
215           END SUBROUTINE WRITE_RES1
216     
217           subroutine gatherWriteRes(VAR, array2, array1, NEXT_REC)
218     
219           USE geometry
220           USE funits
221           USE cdist
222           USE compar           !//
223           USE mpi_utility      !//d pnicol : for gather
224           USE sendrecv         !//d pnicol : for gather
225     
226           USE cutcell
227           USE in_binary_512
228           USE param, only: dimension_3
229     
230           IMPLICIT NONE
231     
232           double precision, dimension(ijkmax2) :: array1
233           double precision, dimension(ijkmax3) :: array2
234           double precision, dimension(DIMENSION_3) :: VAR,TMP_VAR
235     
236           INTEGER :: NEXT_REC
237     
238     !     call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
239           if (.not.bDist_IO) then
240     
241     
242              IF(RE_INDEXING) THEN
243                 CALL UNSHIFT_DP_ARRAY(VAR,TMP_VAR)
244                 CALL gather (TMP_VAR,array2,root)
245              ELSE
246                 CALL gather (VAR,array2,root)
247              ENDIF
248     !         call gather (VAR,array2,root)  !//d pnicol
249     
250     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
251              if (myPE.eq.PE_IO) then
252                 call convert_to_io_dp(array2,array1,ijkmax2)
253                 CALL OUT_BIN_512 (UNIT_RES, array1, IJKMAX2, NEXT_REC)
254              end if
255     
256           else
257     
258              IF(RE_INDEXING) THEN
259                 CALL UNSHIFT_DP_ARRAY(VAR,TMP_VAR)
260                 CALL OUT_BIN_512 (UNIT_RES, TMP_VAR, size(TMP_VAR), NEXT_REC)
261              ELSE
262                 CALL OUT_BIN_512 (UNIT_RES, var, size(var), NEXT_REC)
263              ENDIF
264     !         CALL OUT_BIN_512 (UNIT_RES, var, size(var), NEXT_REC)
265     
266           end if
267     
268           End subroutine gatherWriteRes
269     
270     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
271     !                                                                              !
272     !                                      write_res1_netcdf                       !
273     !                                                                              !
274           subroutine write_res1_netcdf
275     
276             USE param
277             USE param1
278             USE fldvar
279             USE geometry
280             USE physprop
281             USE run
282             USE scalars
283             USE rxns
284             USE cdist
285             USE compar
286             USE mpi_utility
287             USE MFIX_netcdf
288             USE energy
289             USE in_binary_512
290     !       USE tmp_array
291     
292     
293             implicit none
294     
295             integer :: I , nn
296     
297             integer   :: ncid , xyz_dimid
298             integer   :: varid_time , t_dimid
299             integer   :: dimids(1) , varid_epg , varid_pg,dims_time(1)
300             integer   :: varid_pstar  , varid_ug , varid_vg , varid_wg
301             integer   :: varid_tg  , varid_ropg
302             integer   :: varid_rog , varid_gamaRG , varid_TRG
303             integer   :: varid_gamaRS(20) , varid_TRS(20)
304     
305             integer   :: varid_us(20) , varid_vs(20) , varid_ws(20)  !! MMAX
306             integer   :: varid_rops(20)  , varid_ts(20) !! mmax
307             integer   :: varid_thetam(20) !! mmax
308     
309             integer   :: varid_xg(20)  ! nmax(0)
310             integer   :: varid_xs(20,20)  ! mmax , MAX(nmax(1:mmax))
311     
312             integer   :: varid_scalar(20)  ! nscalar
313             integer   :: varid_rr(20)      ! nRR
314     
315             integer   :: varid_kturbg , varid_eturbg
316     
317     
318             character(LEN=80) :: fname, var_name
319     
320             double precision, allocatable :: arr1(:)
321             double precision, allocatable :: arr2(:)
322     
323     ! bWrite_netcdf(1)  : EP_g
324     ! bWrite_netcdf(2)  : P_g
325     ! bWrite_netcdf(3)  : P_star
326     ! bWrite_netcdf(4)  : U_g / V_g / W_g
327     ! bWrite_netcdf(5)  : U_s / V_s / W_s
328     ! bWrite_netcdf(6)  : ROP_s
329     ! bWrite_netcdf(7)  : T_g
330     ! bWrite_netcdf(8)  : T_s
331     ! bWrite_netcdf(9)  : X_g
332     ! bWrite_netcdf(10) : X_s
333     ! bWrite_netcdf(11) : Theta_m
334     ! bWrite_netCDF(12) : Scalar
335     ! bWrite_netCDF(13) : ReactionRates
336     ! bWrite_netCDF(14) : k_turb_g , e_turb_g
337     
338             if (.not. MFIX_usingNETCDF()) return
339             if (.not. bGlobalNetcdf) return
340     
341             return ! for now ... do not write restart files using netCDF
342     
343             if (myPE .eq. PE_IO) then
344                allocate (arr1(ijkmax2))
345                allocate (arr2(ijkmax3))
346             else
347                allocate (arr1(1))
348                allocate (arr2(1))
349             end if
350     
351     !       call mpi$barrier(MPI_COMM_WORLD,mpierr)
352             if (myPE .ne. PE_IO) goto 1234
353     
354             fname = trim(run_name) // "_RES1.nc"
355             call MFIX_check_netcdf( MFIX_nf90_create(fname, NF90_CLOBBER, ncid) )
356     
357             call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "xyz", ijkmax2, xyz_dimid) )
358             call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "t"  ,         1        ,   t_dimid) )
359     
360             ! The dimids array is used to pass the IDs of the dimensions of
361             ! the variables. Note that in fortran arrays are stored in
362             ! column-major format.
363             dims_time(1) = t_dimid
364             dimids =  (/ xyz_dimid /)
365     
366             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid ,"time"  , NF90_DOUBLE ,dims_time , varid_time ) )
367             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "EP_g"  , NF90_DOUBLE, dimids    , varid_epg  ) )
368             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "P_g"   , NF90_DOUBLE, dimids    , varid_pg   ) )
369             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "P_star", NF90_DOUBLE, dimids    , varid_pstar) )
370             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "U_g"   , NF90_DOUBLE, dimids    , varid_ug   ) )
371             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "V_g"   , NF90_DOUBLE, dimids    , varid_vg   ) )
372             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "W_g"   , NF90_DOUBLE, dimids    , varid_wg   ) )
373             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "T_g"   , NF90_DOUBLE, dimids    , varid_tg   ) )
374             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "ROP_g" , NF90_DOUBLE, dimids    , varid_ropg ) )
375             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "RO_g"  , NF90_DOUBLE, dimids    , varid_rog  ) )
376             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "gamaRG"  , NF90_DOUBLE, dimids    , varid_gamaRG  ) )
377             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "TRG"  , NF90_DOUBLE, dimids       , varid_TRG  ) )
378     
379             do i = 1,1   ! mmax
380                var_name = 'U_s_xxx'
381                write (var_name(5:7),'(i3.3)') I
382                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_us(I)) )
383     
384                var_name = 'V_s_xxx'
385                write (var_name(5:7),'(i3.3)') I
386                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_vs(I)) )
387     
388                var_name = 'W_s_xxx'
389                write (var_name(5:7),'(i3.3)') I
390                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_ws(I)) )
391     
392                var_name = 'ROP_s_xxx'
393                write (var_name(7:10),'(i3.3)') I
394                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_rops(I)) )
395     
396                var_name = 'T_s_xxx'
397                write (var_name(5:7),'(i3.3)') I
398                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_ts(I)) )
399     
400                var_name = 'Theta_m_xxx'
401                write (var_name(9:11),'(i3.3)') I
402                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_thetam(I)) )
403     
404                var_name = 'gamaRS_xxx'
405                write (var_name(8:10),'(i3.3)') I
406                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_gamaRS(I)) )
407     
408                var_name = 'TRS_xxx'
409                write (var_name(5:7),'(i3.3)') I
410                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_trs(I)) )
411     
412                DO NN = 1, NMAX(i)
413                   var_name = 'X_s_xxx_xxx'
414                   write (var_name(5:7) ,'(i3.3)') I
415                   write (var_name(9:11),'(i3.3)') nn
416                   call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_xs(I,nn)) )
417                END DO
418     
419     
420             end do
421     
422             do i = 1,nmax(0)
423                var_name = 'X_g_xxx'
424                write (var_name(5:7),'(i3.3)') I
425                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_xg(I)) )
426             end do
427     
428             do i = 1,nscalar
429                var_name = 'Scalar_xxx'
430                write (var_name(8:10),'(i3.3)') I
431                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_scalar(I)) )
432             end do
433     
434             do i = 1,nRR
435                var_name = 'RRates_xxx'
436                write (var_name(8:10),'(i3.3)') I
437                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_rr(I)) )
438             end do
439     
440     
441             if (k_Epsilon) then
442                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, 'k_turb_g', NF90_DOUBLE, dimids, varid_kturbg) )
443                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, 'e_turb_g', NF90_DOUBLE, dimids, varid_eturbg) )
444             end if
445     
446     
447             call MFIX_check_netcdf( MFIX_nf90_enddef(ncid) )
448     
449      1234   continue
450     !       call mpi$barrier(MPI_COMM_WORLD,mpierr)
451     
452             if (myPE .eq. PE_IO) then
453     !           call check_netcdf( nf90_put_var(ncid,varid_time,the_time) )   ! ???
454                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_time,time) )       ! ???
455             end if
456     
457             call gather(EP_g,arr2,root)
458     
459             if (myPE .eq. PE_IO) then
460                call convert_to_io_dp(arr2,arr1,ijkmax2)
461                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_epg, arr1) )
462             end if
463     
464             call gather(P_g,arr2,root)
465             if (myPE .eq. PE_IO) then
466                call convert_to_io_dp(arr2,arr1,ijkmax2)
467                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_pg, arr1) )
468             end if
469     
470             call gather(P_Star,arr2,root)
471             if (myPE .eq. PE_IO) then
472                call convert_to_io_dp(arr2,arr1,ijkmax2)
473                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_pstar, arr1) )
474             end if
475     
476             call gather(U_g,arr2,root)
477             if (myPE .eq. PE_IO) then
478                call convert_to_io_dp(arr2,arr1,ijkmax2)
479                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_ug, arr1) )
480             end if
481     
482             call gather(V_g,arr2,root)
483             if (myPE .eq. PE_IO) then
484                call convert_to_io_dp(arr2,arr1,ijkmax2)
485                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_vg, arr1) )
486             end if
487     
488             call gather(W_g,arr2,root)
489             if (myPE .eq. PE_IO) then
490                call convert_to_io_dp(arr2,arr1,ijkmax2)
491                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_wg, arr1) )
492             end if
493     
494             call gather(T_g,arr2,root)
495             if (myPE .eq. PE_IO) then
496                call convert_to_io_dp(arr2,arr1,ijkmax2)
497                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_tg, arr1) )
498             end if
499     
500             call gather(gama_rg,arr2,root)
501             if (myPE .eq. PE_IO) then
502                call convert_to_io_dp(arr2,arr1,ijkmax2)
503                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_gamarg, arr1) )
504             end if
505     
506             call gather(ro_g,arr2,root)
507             if (myPE .eq. PE_IO) then
508                call convert_to_io_dp(arr2,arr1,ijkmax2)
509                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_rog, arr1) )
510             end if
511     
512             call gather(rop_g,arr2,root)
513             if (myPE .eq. PE_IO) then
514                call convert_to_io_dp(arr2,arr1,ijkmax2)
515                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_ropg, arr1) )
516             end if
517     
518             call gather(t_rg,arr2,root)
519             if (myPE .eq. PE_IO) then
520                call convert_to_io_dp(arr2,arr1,ijkmax2)
521                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_trg, arr1) )
522             end if
523     
524             do i = 1,1   ! mmax
525     
526                call gather(U_s(:,i),arr2,root)
527                if (myPE .eq. PE_IO) then
528                   call convert_to_io_dp(arr2,arr1,ijkmax2)
529                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_us(i), arr1) )
530                end if
531     
532                call gather(V_s(:,i),arr2,root)
533                if (myPE .eq. PE_IO) then
534                   call convert_to_io_dp(arr2,arr1,ijkmax2)
535                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_vs(i), arr1) )
536                end if
537     
538     
539                call gather(W_s(:,i),arr2,root)
540                if (myPE .eq. PE_IO) then
541                   call convert_to_io_dp(arr2,arr1,ijkmax2)
542                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_ws(i), arr1) )
543                end if
544     
545                call gather(ROP_s(:,i),arr2,root)
546                if (myPE .eq. PE_IO) then
547                   call convert_to_io_dp(arr2,arr1,ijkmax2)
548                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_rops(i), arr1 ) )
549                end if
550     
551                call gather(T_s(:,i),arr2,root)
552                if (myPE .eq. PE_IO) then
553                   call convert_to_io_dp(arr2,arr1,ijkmax2)
554                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_ts(i), arr1 ) )
555                end if
556     
557                call gather(Theta_m(:,i),arr2,root)
558                if (myPE .eq. PE_IO) then
559                   call convert_to_io_dp(arr2,arr1,ijkmax2)
560                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_thetam(i), arr1 ) )
561                end if
562     
563                call gather(gama_rs(:,i),arr2,root)
564                if (myPE .eq. PE_IO) then
565                   call convert_to_io_dp(arr2,arr1,ijkmax2)
566                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_gamaRS(i), arr1 ) )
567                end if
568     
569                call gather(T_rs(:,i),arr2,root)
570                if (myPE .eq. PE_IO) then
571                   call convert_to_io_dp(arr2,arr1,ijkmax2)
572                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_TRS(i) , arr1 ) )
573                end if
574     
575                do nn = 1,nmax(i)
576                   call gather(X_s(:,i,NN),arr2,root)
577                   if (myPE .eq. PE_IO) then
578                      call convert_to_io_dp(arr2,arr1,ijkmax2)
579                      call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_xs(i,NN), arr1 ) )
580                   end if
581                end do
582     
583             end do
584     
585             do i = 1,nmax(0)
586                call gather(X_g(:,i),arr2,root)
587                if (myPE .eq. PE_IO) then
588                   call convert_to_io_dp(arr2,arr1,ijkmax2)
589                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_xg(i), arr1 ) )
590                end if
591             end do
592     
593             do i = 1,nscalar
594                call gather(Scalar(:,i),arr2,root)
595                if (myPE .eq. PE_IO) then
596                   call convert_to_io_dp(arr2,arr1,ijkmax2)
597                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_scalar(i), arr1 ) )
598                end if
599             end do
600     
601             do i = 1,nRR
602                call gather(ReactionRates(:,i),arr2,root)
603                if (myPE .eq. PE_IO) then
604                   call convert_to_io_dp(arr2,arr1,ijkmax2)
605                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_rr(i), arr1 ) )
606                end if
607             end do
608     
609             if (k_epsilon) then
610                call gather(k_turb_g,arr2,root)
611                if (myPE .eq. PE_IO) then
612                   call convert_to_io_dp(arr2,arr1,ijkmax2)
613                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_kturbg, arr1) )
614                end if
615     
616                call gather(e_turb_g,arr2,root)
617                if (myPE .eq. PE_IO) then
618                   call convert_to_io_dp(arr2,arr1,ijkmax2)
619                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_eturbg, arr1) )
620                end if
621             end if
622     
623      2222   continue
624     
625             ! Close the file. This frees up any internal netCDF resources
626             ! associated with the file, and flushes any buffers.
627             if (myPE .eq. PE_IO) then
628                call MFIX_check_netcdf( MFIX_nf90_close(ncid) )
629             end if
630     
631             deallocate (arr1)
632             deallocate (arr2)
633     
634             return
635     
636           end subroutine write_res1_netcdf
637     
638     
639     
640