File: RELATIVE:/../../../mfix.git/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, N
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 N = 1, NMAX(0)
140                 call gatherWriteRes (X_g(:,n),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 N = 1, NMAX(LC)
167                 call gatherWriteRes (X_s(:,LC,N),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     
229           IMPLICIT NONE
230     
231           double precision, dimension(ijkmax2) :: array1
232           double precision, dimension(ijkmax3) :: array2
233           double precision, dimension(DIMENSION_3) :: VAR,TMP_VAR
234     
235           INTEGER :: NEXT_REC
236     
237     !     call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
238           if (.not.bDist_IO) then
239     
240     
241              IF(RE_INDEXING) THEN
242                 CALL UNSHIFT_DP_ARRAY(VAR,TMP_VAR)
243                 CALL gather (TMP_VAR,array2,root)
244              ELSE
245                 CALL gather (VAR,array2,root)
246              ENDIF
247     !         call gather (VAR,array2,root)  !//d pnicol
248     
249     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
250              if (myPE.eq.PE_IO) then
251                 call convert_to_io_dp(array2,array1,ijkmax2)
252                 CALL OUT_BIN_512 (UNIT_RES, array1, IJKMAX2, NEXT_REC)
253              end if
254     
255           else
256     
257              IF(RE_INDEXING) THEN
258                 CALL UNSHIFT_DP_ARRAY(VAR,TMP_VAR)
259                 CALL OUT_BIN_512 (UNIT_RES, TMP_VAR, size(TMP_VAR), NEXT_REC)
260              ELSE
261                 CALL OUT_BIN_512 (UNIT_RES, var, size(var), NEXT_REC)
262              ENDIF
263     !         CALL OUT_BIN_512 (UNIT_RES, var, size(var), NEXT_REC)
264     
265           end if
266     
267           End subroutine gatherWriteRes
268     
269     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
270     !                                                                              !
271     !                                      write_res1_netcdf                       !
272     !                                                                              !
273           subroutine write_res1_netcdf
274     
275             USE param
276             USE param1
277             USE fldvar
278             USE geometry
279             USE physprop
280             USE run
281             USE scalars
282             USE rxns
283             USE cdist
284             USE compar
285             USE mpi_utility
286             USE MFIX_netcdf
287             USE energy
288             USE in_binary_512
289     !       USE tmp_array
290     
291     
292             implicit none
293     
294             integer :: I , n
295     
296             integer   :: ncid , xyz_dimid
297             integer   :: varid_time , t_dimid
298             integer   :: dimids(1) , varid_epg , varid_pg,dims_time(1)
299             integer   :: varid_pstar  , varid_ug , varid_vg , varid_wg
300             integer   :: varid_tg  , varid_ropg
301             integer   :: varid_rog , varid_gamaRG , varid_TRG
302             integer   :: varid_gamaRS(20) , varid_TRS(20)
303     
304             integer   :: varid_us(20) , varid_vs(20) , varid_ws(20)  !! MMAX
305             integer   :: varid_rops(20)  , varid_ts(20) !! mmax
306             integer   :: varid_thetam(20) !! mmax
307     
308             integer   :: varid_xg(20)  ! nmax(0)
309             integer   :: varid_xs(20,20)  ! mmax , MAX(nmax(1:mmax))
310     
311             integer   :: varid_scalar(20)  ! nscalar
312             integer   :: varid_rr(20)      ! nRR
313     
314             integer   :: varid_kturbg , varid_eturbg
315     
316     
317             character(LEN=80) :: fname, var_name
318     
319             double precision, allocatable :: arr1(:)
320             double precision, allocatable :: arr2(:)
321     
322     ! bWrite_netcdf(1)  : EP_g
323     ! bWrite_netcdf(2)  : P_g
324     ! bWrite_netcdf(3)  : P_star
325     ! bWrite_netcdf(4)  : U_g / V_g / W_g
326     ! bWrite_netcdf(5)  : U_s / V_s / W_s
327     ! bWrite_netcdf(6)  : ROP_s
328     ! bWrite_netcdf(7)  : T_g
329     ! bWrite_netcdf(8)  : T_s
330     ! bWrite_netcdf(9)  : X_g
331     ! bWrite_netcdf(10) : X_s
332     ! bWrite_netcdf(11) : Theta_m
333     ! bWrite_netCDF(12) : Scalar
334     ! bWrite_netCDF(13) : ReactionRates
335     ! bWrite_netCDF(14) : k_turb_g , e_turb_g
336     
337             if (.not. MFIX_usingNETCDF()) return
338             if (.not. bGlobalNetcdf) return
339     
340             return ! for now ... do not write restart files using netCDF
341     
342             if (myPE .eq. PE_IO) then
343                allocate (arr1(ijkmax2))
344                allocate (arr2(ijkmax3))
345             else
346                allocate (arr1(1))
347                allocate (arr2(1))
348             end if
349     
350     !       call mpi$barrier(MPI_COMM_WORLD,mpierr)
351             if (myPE .ne. PE_IO) goto 1234
352     
353             fname = trim(run_name) // "_RES1.nc"
354             call MFIX_check_netcdf( MFIX_nf90_create(fname, NF90_CLOBBER, ncid) )
355     
356             call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "xyz", ijkmax2, xyz_dimid) )
357             call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "t"  ,         1        ,   t_dimid) )
358     
359             ! The dimids array is used to pass the IDs of the dimensions of
360             ! the variables. Note that in fortran arrays are stored in
361             ! column-major format.
362             dims_time(1) = t_dimid
363             dimids =  (/ xyz_dimid /)
364     
365             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid ,"time"  , NF90_DOUBLE ,dims_time , varid_time ) )
366             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "EP_g"  , NF90_DOUBLE, dimids    , varid_epg  ) )
367             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "P_g"   , NF90_DOUBLE, dimids    , varid_pg   ) )
368             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "P_star", NF90_DOUBLE, dimids    , varid_pstar) )
369             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "U_g"   , NF90_DOUBLE, dimids    , varid_ug   ) )
370             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "V_g"   , NF90_DOUBLE, dimids    , varid_vg   ) )
371             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "W_g"   , NF90_DOUBLE, dimids    , varid_wg   ) )
372             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "T_g"   , NF90_DOUBLE, dimids    , varid_tg   ) )
373             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "ROP_g" , NF90_DOUBLE, dimids    , varid_ropg ) )
374             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "RO_g"  , NF90_DOUBLE, dimids    , varid_rog  ) )
375             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "gamaRG"  , NF90_DOUBLE, dimids    , varid_gamaRG  ) )
376             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "TRG"  , NF90_DOUBLE, dimids       , varid_TRG  ) )
377     
378             do i = 1,1   ! mmax
379                var_name = 'U_s_xxx'
380                write (var_name(5:7),'(i3.3)') I
381                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_us(I)) )
382     
383                var_name = 'V_s_xxx'
384                write (var_name(5:7),'(i3.3)') I
385                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_vs(I)) )
386     
387                var_name = 'W_s_xxx'
388                write (var_name(5:7),'(i3.3)') I
389                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_ws(I)) )
390     
391                var_name = 'ROP_s_xxx'
392                write (var_name(7:10),'(i3.3)') I
393                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_rops(I)) )
394     
395                var_name = 'T_s_xxx'
396                write (var_name(5:7),'(i3.3)') I
397                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_ts(I)) )
398     
399                var_name = 'Theta_m_xxx'
400                write (var_name(9:11),'(i3.3)') I
401                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_thetam(I)) )
402     
403                var_name = 'gamaRS_xxx'
404                write (var_name(8:10),'(i3.3)') I
405                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_gamaRS(I)) )
406     
407                var_name = 'TRS_xxx'
408                write (var_name(5:7),'(i3.3)') I
409                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_trs(I)) )
410     
411                DO N = 1, NMAX(i)
412                   var_name = 'X_s_xxx_xxx'
413                   write (var_name(5:7) ,'(i3.3)') I
414                   write (var_name(9:11),'(i3.3)') n
415                   call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_xs(I,n)) )
416                END DO
417     
418     
419             end do
420     
421             do i = 1,nmax(0)
422                var_name = 'X_g_xxx'
423                write (var_name(5:7),'(i3.3)') I
424                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_xg(I)) )
425             end do
426     
427             do i = 1,nscalar
428                var_name = 'Scalar_xxx'
429                write (var_name(8:10),'(i3.3)') I
430                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_scalar(I)) )
431             end do
432     
433             do i = 1,nRR
434                var_name = 'RRates_xxx'
435                write (var_name(8:10),'(i3.3)') I
436                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_rr(I)) )
437             end do
438     
439     
440             if (k_Epsilon) then
441                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, 'k_turb_g', NF90_DOUBLE, dimids, varid_kturbg) )
442                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, 'e_turb_g', NF90_DOUBLE, dimids, varid_eturbg) )
443             end if
444     
445     
446             call MFIX_check_netcdf( MFIX_nf90_enddef(ncid) )
447     
448      1234   continue
449     !       call mpi$barrier(MPI_COMM_WORLD,mpierr)
450     
451             if (myPE .eq. PE_IO) then
452     !           call check_netcdf( nf90_put_var(ncid,varid_time,the_time) )   ! ???
453                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_time,time) )       ! ???
454             end if
455     
456             call gather(EP_g,arr2,root)
457     
458             if (myPE .eq. PE_IO) then
459                call convert_to_io_dp(arr2,arr1,ijkmax2)
460                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_epg, arr1) )
461             end if
462     
463             call gather(P_g,arr2,root)
464             if (myPE .eq. PE_IO) then
465                call convert_to_io_dp(arr2,arr1,ijkmax2)
466                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_pg, arr1) )
467             end if
468     
469             call gather(P_Star,arr2,root)
470             if (myPE .eq. PE_IO) then
471                call convert_to_io_dp(arr2,arr1,ijkmax2)
472                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_pstar, arr1) )
473             end if
474     
475             call gather(U_g,arr2,root)
476             if (myPE .eq. PE_IO) then
477                call convert_to_io_dp(arr2,arr1,ijkmax2)
478                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_ug, arr1) )
479             end if
480     
481             call gather(V_g,arr2,root)
482             if (myPE .eq. PE_IO) then
483                call convert_to_io_dp(arr2,arr1,ijkmax2)
484                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_vg, arr1) )
485             end if
486     
487             call gather(W_g,arr2,root)
488             if (myPE .eq. PE_IO) then
489                call convert_to_io_dp(arr2,arr1,ijkmax2)
490                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_wg, arr1) )
491             end if
492     
493             call gather(T_g,arr2,root)
494             if (myPE .eq. PE_IO) then
495                call convert_to_io_dp(arr2,arr1,ijkmax2)
496                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_tg, arr1) )
497             end if
498     
499             call gather(gama_rg,arr2,root)
500             if (myPE .eq. PE_IO) then
501                call convert_to_io_dp(arr2,arr1,ijkmax2)
502                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_gamarg, arr1) )
503             end if
504     
505             call gather(ro_g,arr2,root)
506             if (myPE .eq. PE_IO) then
507                call convert_to_io_dp(arr2,arr1,ijkmax2)
508                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_rog, arr1) )
509             end if
510     
511             call gather(rop_g,arr2,root)
512             if (myPE .eq. PE_IO) then
513                call convert_to_io_dp(arr2,arr1,ijkmax2)
514                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_ropg, arr1) )
515             end if
516     
517             call gather(t_rg,arr2,root)
518             if (myPE .eq. PE_IO) then
519                call convert_to_io_dp(arr2,arr1,ijkmax2)
520                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_trg, arr1) )
521             end if
522     
523             do i = 1,1   ! mmax
524     
525                call gather(U_s(:,i),arr2,root)
526                if (myPE .eq. PE_IO) then
527                   call convert_to_io_dp(arr2,arr1,ijkmax2)
528                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_us(i), arr1) )
529                end if
530     
531                call gather(V_s(:,i),arr2,root)
532                if (myPE .eq. PE_IO) then
533                   call convert_to_io_dp(arr2,arr1,ijkmax2)
534                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_vs(i), arr1) )
535                end if
536     
537     
538                call gather(W_s(:,i),arr2,root)
539                if (myPE .eq. PE_IO) then
540                   call convert_to_io_dp(arr2,arr1,ijkmax2)
541                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_ws(i), arr1) )
542                end if
543     
544                call gather(ROP_s(:,i),arr2,root)
545                if (myPE .eq. PE_IO) then
546                   call convert_to_io_dp(arr2,arr1,ijkmax2)
547                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_rops(i), arr1 ) )
548                end if
549     
550                call gather(T_s(:,i),arr2,root)
551                if (myPE .eq. PE_IO) then
552                   call convert_to_io_dp(arr2,arr1,ijkmax2)
553                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_ts(i), arr1 ) )
554                end if
555     
556                call gather(Theta_m(:,i),arr2,root)
557                if (myPE .eq. PE_IO) then
558                   call convert_to_io_dp(arr2,arr1,ijkmax2)
559                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_thetam(i), arr1 ) )
560                end if
561     
562                call gather(gama_rs(:,i),arr2,root)
563                if (myPE .eq. PE_IO) then
564                   call convert_to_io_dp(arr2,arr1,ijkmax2)
565                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_gamaRS(i), arr1 ) )
566                end if
567     
568                call gather(T_rs(:,i),arr2,root)
569                if (myPE .eq. PE_IO) then
570                   call convert_to_io_dp(arr2,arr1,ijkmax2)
571                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_TRS(i) , arr1 ) )
572                end if
573     
574                do N = 1,nmax(i)
575                   call gather(X_s(:,i,N),arr2,root)
576                   if (myPE .eq. PE_IO) then
577                      call convert_to_io_dp(arr2,arr1,ijkmax2)
578                      call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_xs(i,N), arr1 ) )
579                   end if
580                end do
581     
582             end do
583     
584             do i = 1,nmax(0)
585                call gather(X_g(:,i),arr2,root)
586                if (myPE .eq. PE_IO) then
587                   call convert_to_io_dp(arr2,arr1,ijkmax2)
588                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_xg(i), arr1 ) )
589                end if
590             end do
591     
592             do i = 1,nscalar
593                call gather(Scalar(:,i),arr2,root)
594                if (myPE .eq. PE_IO) then
595                   call convert_to_io_dp(arr2,arr1,ijkmax2)
596                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_scalar(i), arr1 ) )
597                end if
598             end do
599     
600             do i = 1,nRR
601                call gather(ReactionRates(:,i),arr2,root)
602                if (myPE .eq. PE_IO) then
603                   call convert_to_io_dp(arr2,arr1,ijkmax2)
604                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_rr(i), arr1 ) )
605                end if
606             end do
607     
608             if (k_epsilon) then
609                call gather(k_turb_g,arr2,root)
610                if (myPE .eq. PE_IO) then
611                   call convert_to_io_dp(arr2,arr1,ijkmax2)
612                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_kturbg, arr1) )
613                end if
614     
615                call gather(e_turb_g,arr2,root)
616                if (myPE .eq. PE_IO) then
617                   call convert_to_io_dp(arr2,arr1,ijkmax2)
618                   call MFIX_check_netcdf( MFIX_nf90_put_var(ncid, varid_eturbg, arr1) )
619                end if
620             end if
621     
622      2222   continue
623     
624             ! Close the file. This frees up any internal netCDF resources
625             ! associated with the file, and flushes any buffers.
626             if (myPE .eq. PE_IO) then
627                call MFIX_check_netcdf( MFIX_nf90_close(ncid) )
628             end if
629     
630             deallocate (arr1)
631             deallocate (arr2)
632     
633             return
634     
635           end subroutine write_res1_netcdf
636     
637     
638     
639