File: /nfs/home/0/users/jenkins/mfix.git/model/read_res1.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: READ_RES1                                              C
4     !  Purpose: read in the time-dependent restart records                 C
5     !                                                                      C
6     !  Author: P. Nicoletti                               Date: 03-JAN-92  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: IJKMAX2, MMAX, DT                             C
17     !  Variables modified: TIME, NSTEP, EP_g, P_g, P_star, RO_g            C
18     !                      ROP_g, T_g, T_s,  U_g, V_g, W_g, ROP_s    C
19     !                      U_s, V_s, W_s                                   C
20     !                                                                      C
21     !  Local variables: TIME_READ, LC, NEXT_REC                            C
22     !                                                                      C
23     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
24     !
25           SUBROUTINE READ_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 rxns
39           USE scalars
40           USE funits
41           USE energy
42           USE compar
43           USE cdist
44           USE mpi_utility
45           USE sendrecv
46           USE in_binary_512
47     
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     !             loop counter
60           INTEGER LC
61     !
62     !                      Local species index
63           INTEGER          N
64     !
65     !             pointer to the next record
66           INTEGER NEXT_REC
67     !
68     !                file version id
69           CHARACTER(LEN=512) :: VERSION
70     !
71     !                version number
72           REAL       VERSION_NUMBER
73     !
74     !                      Dummy array
75           DOUBLE PRECISION DT_SAVE
76     
77     !//PAR_I/O declare global scratch arrays
78           double precision, allocatable :: array1(:)
79           double precision, allocatable :: array2(:)
80     !-----------------------------------------------
81     !
82           if (myPE .eq. PE_IO .or. .not.bStart_with_one_res) then
83              allocate (array1(ijkmax2))
84              allocate (array2(ijkmax3))
85           else
86              allocate (array1(1))
87              allocate (array2(1))
88           end if
89     
90     !      call MPI_barrier(MPI_COMM_WORLD,mpierr)
91     !
92     !     Use DT from data file if DT_FAC is set to 1.0
93           IF (DT_FAC == ONE) DT_SAVE = DT
94     !
95     
96     !//PAR_I/O only PE_IO reads the restart file
97           if (myPE == PE_IO .or. (bDist_IO .and. .not.bStart_with_one_RES)) then
98              READ (UNIT_RES, REC=1) VERSION
99              READ (VERSION(6:512), *) VERSION_NUMBER
100     
101              READ (UNIT_RES, REC=3) NEXT_REC
102              IF (VERSION_NUMBER >= 1.12) THEN
103                 READ (UNIT_RES, REC=NEXT_REC) TIME, DT, NSTEP
104              ELSE
105                 READ (UNIT_RES, REC=NEXT_REC) TIME, NSTEP
106              ENDIF
107              NEXT_REC = NEXT_REC + 1
108           end if
109     
110     
111     
112           if (.not.bDist_IO  .or. bStart_with_one_RES) then
113     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
114     
115              call bcast(VERSION, PE_IO)        !//PAR_I/O BCAST0c
116              call bcast(VERSION_NUMBER, PE_IO) !//PAR_I/O BCAST0r
117              call bcast(TIME, PE_IO)           !//PAR_I/O BCAST0d
118              call bcast(NSTEP, PE_IO)          !//PAR_I/O BCAST0i
119              if (VERSION_NUMBER >= 1.12) call bcast(DT, PE_IO)   !//PAR_I/O BCAST0d
120             end if
121     !      call MPI_barrier(MPI_COMM_WORLD,mpierr)
122     
123     !AE TIME 091501 Store the timestep counter level at the begin of RESTART run
124             NSTEPRST = NSTEP
125     
126     ! for now ... do not do RES1 file in netCDF format
127     !       call read_res1_netcdf
128     !        goto 999
129     !
130     
131     
132     !
133           call readScatterRes(EP_G,array2, array1, 0, NEXT_REC)
134     
135           call readScatterRes(P_G,array2, array1, 0, NEXT_REC)
136     
137           call readScatterRes(P_STAR,array2, array1, 1, NEXT_REC)
138     
139           call readScatterRes(RO_G,array2, array1, 0, NEXT_REC)
140     
141           call readScatterRes(ROP_G,array2, array1, 0, NEXT_REC)
142     
143           call readScatterRes(T_G,array2, array1, 1, NEXT_REC)
144     !
145     
146           IF (VERSION_NUMBER < 1.15) THEN
147              call readScatterRes (T_s(:,1),array2, array1, 1, NEXT_REC)
148     
149              IF (MMAX >= 2) THEN
150                 call readScatterRes (T_s(:,2),array2, array1, 1, NEXT_REC)
151              ELSE
152                 if (myPE == PE_IO) &
153                        CALL IN_BIN_512 (UNIT_RES, array1, IJKMAX2, NEXT_REC)
154              ENDIF
155           ENDIF
156     
157           IF (VERSION_NUMBER >= 1.05) THEN
158              DO N = 1, NMAX(0)
159                 call readScatterRes (X_g(:,n),array2, array1, 1, NEXT_REC)
160              END DO
161           ENDIF
162     
163           call readScatterRes(U_G, array2, array1, 0, NEXT_REC)
164           call readScatterRes(V_G, array2, array1, 0, NEXT_REC)
165           call readScatterRes(W_G, array2, array1, 0, NEXT_REC)
166     !
167           DO LC = 1, MMAX
168              call readScatterRes(ROP_S(:,LC), array2, array1, 0, NEXT_REC)
169     
170              IF(ANY(SOLVE_ROs)) &
171                 CALL readScatterRes(RO_S(:,LC), array2, array1, 0, NEXT_REC)
172     
173              IF (VERSION_NUMBER >= 1.15) THEN
174                 call readScatterRes(T_S(:,LC), array2, array1, 1, NEXT_REC)
175              END IF
176              call readScatterRes(U_S(:,LC), array2, array1, 0, NEXT_REC)
177              call readScatterRes(V_S(:,LC), array2, array1, 0, NEXT_REC)
178              call readScatterRes(W_S(:,LC), array2, array1, 0, NEXT_REC)
179     
180              IF (VERSION_NUMBER >= 1.2) then
181                 call readScatterRes(THETA_M(:,LC), array2, array1, 1, NEXT_REC)
182              end if
183              IF (VERSION_NUMBER >= 1.05) THEN
184                 DO N = 1, NMAX(LC)
185                    call readScatterRes(X_S(:,LC,N), array2, array1, 1, NEXT_REC)
186                 END DO
187              ENDIF
188           END DO
189     
190           IF (VERSION_NUMBER >= 1.3) THEN
191             DO N = 1, NScalar
192               call readScatterRes(Scalar(:,N), array2, array1, 1, NEXT_REC)
193             END DO
194           ENDIF
195     
196           IF (VERSION_NUMBER >= 1.4) THEN
197             call readScatterRes(GAMA_RG, array2, array1, 1, NEXT_REC)
198     
199             call readScatterRes(T_RG, array2, array1, 0, NEXT_REC)
200     
201             DO LC = 1, MMAX
202               call readScatterRes(GAMA_RS(1,LC), array2, array1, 1, NEXT_REC)
203     
204               call readScatterRes(T_RS(1,LC), array2, array1, 0, NEXT_REC)
205     
206             ENDDO
207           ELSE
208             GAMA_RG(:)   = ZERO
209             T_RG (:)     = ZERO
210             GAMA_RS(:,:) = ZERO
211             T_RS(:,:)    = ZERO
212           ENDIF
213     
214     
215           IF (VERSION_NUMBER >= 1.5) THEN
216             DO N = 1, nRR
217               call readScatterRes(ReactionRates(:,N), array2, array1, 1, NEXT_REC)
218             END DO
219           ENDIF
220     
221           IF (VERSION_NUMBER >= 1.6 .AND. K_Epsilon) THEN
222               call readScatterRes(K_Turb_G, array2, array1, 1, NEXT_REC)
223               call readScatterRes(E_Turb_G, array2, array1, 1, NEXT_REC)
224           ENDIF
225     !------------------------------------------------------------------------
226     !
227      999  continue
228     
229     !      call MPI_barrier(MPI_COMM_WORLD,mpierr)
230           deallocate( array1 )
231           deallocate( array2 )
232     !      call MPI_barrier(MPI_COMM_WORLD,mpierr)
233     
234           if (.not.bDist_IO .or. bStart_with_one_RES) then
235              call send_recv(rop_g)
236              call send_recv(ro_g)
237              call send_recv(rop_s)
238           end if
239     
240     
241           IF (DT_FAC == ONE) DT = DT_SAVE
242     !
243     
244           RETURN
245           END SUBROUTINE READ_RES1
246     
247           subroutine readScatterRes(VAR, array2, array1, init, NEXT_REC)
248             USE geometry
249             USE funits
250             USE compar
251               USE cdist
252             USE mpi_utility
253             USE sendrecv
254             USE in_binary_512
255             IMPLICIT NONE
256             double precision, dimension(ijkmax2) :: array1
257             double precision, dimension(ijkmax3) :: array2
258             double precision, dimension(DIMENSION_3) :: VAR
259             INTEGER :: init  ! define VAR initialization, 0: undefin, 1: zero
260             INTEGER :: NEXT_REC
261     
262     !// Reset global scratch arrays
263             if( init==0 ) then
264               array1(:) = Undefined
265               array2(:) = Undefined
266             else
267               array1(:) = zero
268               array2(:) = zero
269             endif
270     
271             if (.not.bDist_IO .or. bStart_with_one_RES) then
272              if (myPE == PE_IO) then
273                 CALL IN_BIN_512 (UNIT_RES, array1, IJKMAX2, NEXT_REC)
274                 CALL convert_from_io_dp(array1, array2, IJKMAX2)
275              end if
276     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
277              call scatter(VAR, array2, PE_IO)
278     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
279           else
280              CALL IN_BIN_512 (UNIT_RES, var, size(var) , NEXT_REC)
281             end if
282     
283           End subroutine readScatterRes
284     
285     
286           subroutine readScatterRes_netcdf(VAR, array2, array1, ncid , varid)
287     
288           USE geometry
289           USE compar
290           USE cdist
291           USE mpi_utility
292           USE sendrecv
293           USE MFIX_netcdf
294           USE in_binary_512
295     
296           IMPLICIT NONE
297     
298           double precision, dimension(ijkmax2)     :: array1
299           double precision, dimension(ijkmax3)     :: array2
300           double precision, dimension(DIMENSION_3) :: VAR
301     
302           integer :: ncid , varid
303     
304     
305     
306           if (myPE .eq. PE_IO) then
307              call MFIX_check_netcdf( MFIX_nf90_get_var(ncid , varid   , array1   ) )
308              call convert_from_io_dp(array1,array2,ijkmax2)
309           end if
310     
311     !     call MPI_barrier(MPI_COMM_WORLD,mpierr)
312           call scatter(var,array2,PE_IO)
313     !     call MPI_barrier(MPI_COMM_WORLD,mpierr)
314     
315           End subroutine readScatterRes_netcdf
316     
317     
318     
319     
320     
321     !// Comments on the modifications for DMP version implementation
322     !// 001 Include header file and common declarations for parallelization
323     !// 020 New local variables for parallelization: array1, array2
324     !// 400 Added sendrecv module and send_recv calls for COMMunication
325     !// 400 Added mpi_utility module and other global reduction (bcast) calls
326     
327     
328     
329     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
330     !                                                                              !
331     !                                      read_res1_netcdf                       !
332     !                                                                              !
333           subroutine read_res1_netcdf
334     
335             USE param
336             USE param1
337             USE fldvar
338             USE geometry
339             USE physprop
340             USE run
341     !       USE funits
342             USE scalars
343     !       USE output
344             USE rxns
345             USE cdist
346             USE compar
347             USE mpi_utility
348             USE MFIX_netcdf
349     !       USE tmp_array
350             USE energy
351     
352     
353             implicit none
354     
355             integer :: I , n
356     
357             integer   :: ncid
358             integer   :: varid_time
359             integer   :: varid_epg , varid_pg
360             integer   :: varid_pstar  , varid_ug , varid_vg , varid_wg
361             integer   :: varid_tg
362             integer   :: varid_rog , varid_gamaRG , varid_TRG
363             integer   :: varid_gamaRS(20) , varid_TRS(20) , varid_ropg
364     
365             integer   :: varid_us(20) , varid_vs(20) , varid_ws(20)  !! MMAX
366             integer   :: varid_rops(20)  , varid_ts(20) !! mmax
367             integer   :: varid_thetam(20) !! mmax
368     
369             integer   :: varid_xg(20)  ! nmax(0)
370             integer   :: varid_xs(20,20)  ! mmax , MAX(nmax(1:mmax))
371     
372             integer   :: varid_scalar(20)  ! nscalar
373             integer   :: varid_rr(20)      ! nRR
374     
375             integer   :: varid_kturbg , varid_eturbg
376     
377     
378             character(LEN=80) :: fname, var_name
379     
380             integer nDim , nVars , nAttr , unID , formatNUM
381     
382             integer xyz_id , xyz_dim
383     
384             character(LEN=80) :: varname
385             integer vartype,nvdims,vdims(10),nvatts,rcode
386     
387             double precision, allocatable :: array1(:)
388             double precision, allocatable :: array2(:)
389     
390     
391     !        integer :: MFIX_nf90_create
392     !        integer :: MFIX_nf90_def_dim
393     !        integer :: MFIX_nf90_close
394     !        integer :: MFIX_nf90_open
395     !       integer :: MFIX_nf90_inquire
396     !       integer :: MFIX_nf90_inq_dimid
397     !       integer :: MFIX_nf90_inquire_dimension
398     !       integer :: MFIX_nf90_inq_varid
399     
400     
401     
402     
403     
404     
405     
406     
407     !
408     !
409     
410     
411     ! bWrite_netcdf(1)  : EP_g
412     ! bWrite_netcdf(2)  : P_g
413     ! bWrite_netcdf(3)  : P_star
414     ! bWrite_netcdf(4)  : U_g / V_g / W_g
415     ! bWrite_netcdf(5)  : U_s / V_s / W_s
416     ! bWrite_netcdf(6)  : ROP_s
417     ! bWrite_netcdf(7)  : T_g
418     ! bWrite_netcdf(8)  : T_s
419     ! bWrite_netcdf(9)  : X_g
420     ! bWrite_netcdf(10) : X_s
421     ! bWrite_netcdf(11) : Theta_m
422     ! bWrite_netCDF(12) : Scalar
423     ! bWrite_netCDF(13) : ReactionRates
424     ! bWrite_netCDF(14) : k_turb_g , e_turb_g
425     
426       !    return
427     
428             if (.not. MFIX_usingNETCDF()) return
429     
430           if (myPE .eq. PE_IO) then
431              allocate (array1(ijkmax2))
432              allocate (array2(ijkmax3))
433           else
434              allocate (array1(1))
435              allocate (array2(1))
436           end if
437     
438     
439     
440     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
441           if (myPE .eq. PE_IO) then
442              fname = trim(run_name) // "_RES1.nc"
443     
444              call MFIX_check_netcdf( MFIX_nf90_open(fname, NF90_NOWRITE, ncid) )
445     
446              call MFIX_check_netcdf( MFIX_nf90_inquire(ncid, nDim , nVars , nAttr , unID , formatNUM) )
447     
448     !         write (*,*) ' nDim      = ' , nDim
449     !         write (*,*) ' nVars     = ' , nVars
450     !         write (*,*) ' nAttr     = ' , nAttr
451     !         write (*,*) ' unID      = ' , unID
452     !         write (*,*) ' formatNum = ' , formatNum
453     
454              call MFIX_check_netcdf( MFIX_nf90_inq_dimid(ncid,"xyz",xyz_id) )
455              call MFIX_check_netcdf( MFIX_nf90_inquire_dimension(ncid,xyz_id,len=xyz_dim) )
456     
457     
458              do i = 1,nVars
459                 call MFIX_ncvinq(ncid,i,varname,vartype,nvdims,vdims,nvatts,rcode)
460              end do
461     
462              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid ,"time"  , varid_time  ) )
463              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "EP_g"  , varid_epg   ) )
464              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "P_g"   , varid_pg    ) )
465              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "P_star", varid_pstar ) )
466              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "U_g"   , varid_ug    ) )
467              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "V_g"   , varid_vg    ) )
468              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "W_g"   , varid_wg    ) )
469              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "T_g"   , varid_tg    ) )
470              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "RO_g"   , varid_rog  ) )
471              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "ROP_g"  , varid_ropg ) )
472              call MFIX_check_netcdf( MFIX_nf90_get_var(ncid,varid_time,time) )
473     
474           end if
475     
476     
477     !      if (myPE .eq. PE_IO) then
478     !        call check_netcdf( nf90_get_var(ncid , varid_epg   , array1   ) )
479     !         call convert_from_io_dp(array1,array2,ijkmax2)
480     !      end if
481     !      call scatter(EP_g,array2,PE_IO)
482     
483           call readScatterRes_netcdf(EP_g   , array2, array1, ncid , varid_epg)
484           call readScatterRes_netcdf(P_g    , array2, array1, ncid , varid_pg)
485           call readScatterRes_netcdf(P_star , array2, array1, ncid , varid_pstar)
486           call readScatterRes_netcdf(ro_g   , array2, array1, ncid , varid_rog)
487           call readScatterRes_netcdf(rop_g  , array2, array1, ncid , varid_ropg)
488           call readScatterRes_netcdf(u_g    , array2, array1, ncid , varid_ug)
489           call readScatterRes_netcdf(v_g    , array2, array1, ncid , varid_vg)
490           call readScatterRes_netcdf(w_g    , array2, array1, ncid , varid_wg)
491           call readScatterRes_netcdf(t_g    , array2, array1, ncid , varid_tg)
492     
493           do i = 1,1   ! mmax
494     
495              if (myPe .eq. PE_IO) then
496                 var_name = 'U_s_xxx'
497                 write (var_name(5:7),'(i3.3)') I
498                call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_us(I)) )
499              end if
500              call readScatterRes_netcdf(U_s(:,i) , array2, array1, ncid , varid_us(i))
501     
502              if (myPe .eq. PE_IO) then
503                 var_name = 'V_s_xxx'
504                 write (var_name(5:7),'(i3.3)') I
505                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_vs(I)) )
506              end if
507              call readScatterRes_netcdf(V_s(:,i) , array2, array1, ncid , varid_vs(i))
508     
509              if (myPe .eq. PE_IO) then
510                 var_name = 'W_s_xxx'
511                 write (var_name(5:7),'(i3.3)') I
512                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_ws(I)) )
513              end if
514              call readScatterRes_netcdf(W_s(:,i) , array2, array1, ncid , varid_ws(i))
515     
516              if (myPe .eq. PE_IO) then
517                 var_name = 'ROP_s_xxx'
518                 write (var_name(7:10),'(i3.3)') I
519                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_rops(I)) )
520              end if
521              call readScatterRes_netcdf(ROP_s(:,i) , array2, array1, ncid , varid_rops(i))
522     
523              if (myPe .eq. PE_IO) then
524                 var_name = 'T_s_xxx'
525                 write (var_name(5:7),'(i3.3)') I
526                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_ts(I)) )
527              end if
528              call readScatterRes_netcdf(T_s(:,i) , array2, array1, ncid , varid_ts(i))
529     
530              if (myPe .eq. PE_IO) then
531                 var_name = 'Theta_m_xxx'
532                 write (var_name(9:11),'(i3.3)') I
533                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_thetam(I)) )
534              end if
535              call readScatterRes_netcdf(theta_m(:,i) , array2, array1, ncid , varid_thetam(i))
536     
537              if (myPe .eq. PE_IO) then
538                 var_name = 'gamaRS_xxx'
539                 write (var_name(8:10),'(i3.3)') I
540                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_gamaRS(I)) )
541              end if
542              call readScatterRes_netcdf(gama_rs(:,i) , array2, array1, ncid , varid_gamaRS(i))
543     
544              if (myPe .eq. PE_IO) then
545                 var_name = 'TRS_xxx'
546                 write (var_name(5:7),'(i3.3)') I
547                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_TRS(I)) )
548              end if
549              call readScatterRes_netcdf(T_rs(:,i) , array2, array1, ncid , varid_TRS(i))
550     
551              DO N = 1, NMAX(i)
552                 if (myPe .eq. PE_IO) then
553                    var_name = 'X_s_xxx_xxx'
554                    write (var_name(5:7) ,'(i3.3)') I
555                    write (var_name(9:11),'(i3.3)') n
556                    call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_xs(I,n)) )
557                 end if
558                 call readScatterRes_netcdf(X_s(:,i,n) , array2, array1, ncid , varid_xs(i,n))
559              END DO
560     
561           end do
562     
563           do i = 1,nmax(0)
564              if (myPe .eq. PE_IO) then
565                 var_name = 'X_g_xxx'
566                 write (var_name(5:7),'(i3.3)') I
567                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_xg(I)) )
568              end if
569              call readScatterRes_netcdf(X_g(:,i) , array2, array1, ncid , varid_xg(i))
570           end do
571     
572           do i = 1,nscalar
573              if (myPe .eq. PE_IO) then
574                 var_name = 'Scalar_xxx'
575                 write (var_name(8:10),'(i3.3)') I
576                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_scalar(I)) )
577              end if
578              call readScatterRes_netcdf(scalar(:,i) , array2, array1, ncid , varid_scalar(i))
579           end do
580     
581           do i = 1,nRR
582              if (myPe .eq. PE_IO) then
583                 var_name = 'RRates_xxx'
584                 write (var_name(8:10),'(i3.3)') I
585                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_rr(I)) )
586              end if
587              call readScatterRes_netcdf(ReactionRates(:,i) , array2, array1, ncid , varid_rr(i))
588           end do
589     
590           if (k_Epsilon) then
591              if (myPe .eq. PE_IO) then
592                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, 'k_turb_g', varid_kturbg) )
593              end if
594              call readScatterRes_netcdf(k_turb_g , array2, array1, ncid , varid_kturbg)
595     
596              if (myPe .eq. PE_IO) then
597                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, 'e_turb_g', varid_eturbg) )
598              end if
599              call readScatterRes_netcdf(e_turb_g , array2, array1, ncid , varid_eturbg)
600           end if
601     
602          if (myPe .eq. PE_IO) then
603             call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, 'gamaRG', varid_gamaRG) )
604          end if
605          call readScatterRes_netcdf(gama_rg , array2, array1, ncid , varid_gamaRG)
606     
607          if (myPe .eq. PE_IO) then
608             call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, 'TRG', varid_TRG) )
609          end if
610          call readScatterRes_netcdf(T_rg , array2, array1, ncid , varid_TRG)
611     
612             ! Close the file. This frees up any internal netCDF resources
613             ! associated with the file, and flushes any buffers.
614      1234   continue
615     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
616             if (myPE .eq. PE_IO) then
617                call MFIX_check_netcdf( MFIX_nf90_close(ncid) )
618            end if
619     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
620     
621       !      stop
622     
623     
624             return
625     
626           end subroutine read_res1_netcdf
627