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