File: RELATIVE:/../../../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     
228     !      call MPI_barrier(MPI_COMM_WORLD,mpierr)
229           deallocate( array1 )
230           deallocate( array2 )
231     !      call MPI_barrier(MPI_COMM_WORLD,mpierr)
232     
233           if (.not.bDist_IO .or. bStart_with_one_RES) then
234              call send_recv(rop_g)
235              call send_recv(ro_g)
236              call send_recv(rop_s)
237           end if
238     
239     
240           IF (DT_FAC == ONE) DT = DT_SAVE
241     !
242     
243     !     We may no longer need PATCH_AFTER_RESTART
244     !     CALL PATCH_AFTER_RESTART
245     
246           RETURN
247           END SUBROUTINE READ_RES1
248     
249           subroutine readScatterRes(VAR, array2, array1, init, NEXT_REC)
250             use param1, only: zero, undefined
251             use param, only: dimension_3
252             USE geometry
253             USE funits
254             USE compar
255             USE cdist
256             USE mpi_utility
257             USE sendrecv
258             USE in_binary_512
259             IMPLICIT NONE
260             double precision, dimension(ijkmax2) :: array1
261             double precision, dimension(ijkmax3) :: array2
262             double precision, dimension(DIMENSION_3) :: VAR
263             INTEGER :: init  ! define VAR initialization, 0: undefin, 1: zero
264             INTEGER :: NEXT_REC
265     
266     !// Reset global scratch arrays
267             if( init==0 ) then
268               array1(:) = Undefined
269               array2(:) = Undefined
270             else
271               array1(:) = zero
272               array2(:) = zero
273             endif
274     
275             if (.not.bDist_IO .or. bStart_with_one_RES) then
276              if (myPE == PE_IO) then
277                 CALL IN_BIN_512 (UNIT_RES, array1, IJKMAX2, NEXT_REC)
278                 CALL convert_from_io_dp(array1, array2, IJKMAX2)
279              end if
280     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
281              call scatter(VAR, array2, PE_IO)
282     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
283           else
284              CALL IN_BIN_512 (UNIT_RES, var, size(var) , NEXT_REC)
285             end if
286     
287           End subroutine readScatterRes
288     
289     
290           subroutine readScatterRes_netcdf(VAR, array2, array1, ncid , varid)
291           USE param, only: dimension_3
292           USE geometry
293           USE compar
294           USE cdist
295           USE mpi_utility
296           USE sendrecv
297           USE MFIX_netcdf
298           USE in_binary_512
299     
300           IMPLICIT NONE
301     
302           double precision, dimension(ijkmax2)     :: array1
303           double precision, dimension(ijkmax3)     :: array2
304           double precision, dimension(DIMENSION_3) :: VAR
305     
306           integer :: ncid , varid
307     
308     
309     
310           if (myPE .eq. PE_IO) then
311              call MFIX_check_netcdf( MFIX_nf90_get_var(ncid , varid   , array1   ) )
312              call convert_from_io_dp(array1,array2,ijkmax2)
313           end if
314     
315     !     call MPI_barrier(MPI_COMM_WORLD,mpierr)
316           call scatter(var,array2,PE_IO)
317     !     call MPI_barrier(MPI_COMM_WORLD,mpierr)
318     
319           End subroutine readScatterRes_netcdf
320     
321     
322     
323     
324     
325     !// Comments on the modifications for DMP version implementation
326     !// 001 Include header file and common declarations for parallelization
327     !// 020 New local variables for parallelization: array1, array2
328     !// 400 Added sendrecv module and send_recv calls for COMMunication
329     !// 400 Added mpi_utility module and other global reduction (bcast) calls
330     
331     
332     
333     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
334     !                                                                              !
335     !                                      read_res1_netcdf                       !
336     !                                                                              !
337           subroutine read_res1_netcdf
338     
339             USE param
340             USE param1
341             USE fldvar
342             USE geometry
343             USE physprop
344             USE run
345     !       USE funits
346             USE scalars
347     !       USE output
348             USE rxns
349             USE cdist
350             USE compar
351             USE mpi_utility
352             USE MFIX_netcdf
353     !       USE tmp_array
354             USE energy
355     
356     
357             implicit none
358     
359             integer :: I , n
360     
361             integer   :: ncid
362             integer   :: varid_time
363             integer   :: varid_epg , varid_pg
364             integer   :: varid_pstar  , varid_ug , varid_vg , varid_wg
365             integer   :: varid_tg
366             integer   :: varid_rog , varid_gamaRG , varid_TRG
367             integer   :: varid_gamaRS(20) , varid_TRS(20) , varid_ropg
368     
369             integer   :: varid_us(20) , varid_vs(20) , varid_ws(20)  !! MMAX
370             integer   :: varid_rops(20)  , varid_ts(20) !! mmax
371             integer   :: varid_thetam(20) !! mmax
372     
373             integer   :: varid_xg(20)  ! nmax(0)
374             integer   :: varid_xs(20,20)  ! mmax , MAX(nmax(1:mmax))
375     
376             integer   :: varid_scalar(20)  ! nscalar
377             integer   :: varid_rr(20)      ! nRR
378     
379             integer   :: varid_kturbg , varid_eturbg
380     
381     
382             character(LEN=80) :: fname, var_name
383     
384             integer nDim , nVars , nAttr , unID , formatNUM
385     
386             integer xyz_id , xyz_dim
387     
388             character(LEN=80) :: varname
389             integer vartype,nvdims,vdims(10),nvatts,rcode
390     
391             double precision, allocatable :: array1(:)
392             double precision, allocatable :: array2(:)
393     
394     
395     !        integer :: MFIX_nf90_create
396     !        integer :: MFIX_nf90_def_dim
397     !        integer :: MFIX_nf90_close
398     !        integer :: MFIX_nf90_open
399     !       integer :: MFIX_nf90_inquire
400     !       integer :: MFIX_nf90_inq_dimid
401     !       integer :: MFIX_nf90_inquire_dimension
402     !       integer :: MFIX_nf90_inq_varid
403     
404     
405     
406     
407     
408     
409     
410     
411     !
412     !
413     
414     
415     ! bWrite_netcdf(1)  : EP_g
416     ! bWrite_netcdf(2)  : P_g
417     ! bWrite_netcdf(3)  : P_star
418     ! bWrite_netcdf(4)  : U_g / V_g / W_g
419     ! bWrite_netcdf(5)  : U_s / V_s / W_s
420     ! bWrite_netcdf(6)  : ROP_s
421     ! bWrite_netcdf(7)  : T_g
422     ! bWrite_netcdf(8)  : T_s
423     ! bWrite_netcdf(9)  : X_g
424     ! bWrite_netcdf(10) : X_s
425     ! bWrite_netcdf(11) : Theta_m
426     ! bWrite_netCDF(12) : Scalar
427     ! bWrite_netCDF(13) : ReactionRates
428     ! bWrite_netCDF(14) : k_turb_g , e_turb_g
429     
430       !    return
431     
432             if (.not. MFIX_usingNETCDF()) return
433     
434           if (myPE .eq. PE_IO) then
435              allocate (array1(ijkmax2))
436              allocate (array2(ijkmax3))
437           else
438              allocate (array1(1))
439              allocate (array2(1))
440           end if
441     
442     
443     
444     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
445           if (myPE .eq. PE_IO) then
446              fname = trim(run_name) // "_RES1.nc"
447     
448              call MFIX_check_netcdf( MFIX_nf90_open(fname, NF90_NOWRITE, ncid) )
449     
450              call MFIX_check_netcdf( MFIX_nf90_inquire(ncid, nDim , nVars , nAttr , unID , formatNUM) )
451     
452     !         write (*,*) ' nDim      = ' , nDim
453     !         write (*,*) ' nVars     = ' , nVars
454     !         write (*,*) ' nAttr     = ' , nAttr
455     !         write (*,*) ' unID      = ' , unID
456     !         write (*,*) ' formatNum = ' , formatNum
457     
458              call MFIX_check_netcdf( MFIX_nf90_inq_dimid(ncid,"xyz",xyz_id) )
459              call MFIX_check_netcdf( MFIX_nf90_inquire_dimension(ncid,xyz_id,len=xyz_dim) )
460     
461     
462              do i = 1,nVars
463                 call MFIX_ncvinq(ncid,i,varname,vartype,nvdims,vdims,nvatts,rcode)
464              end do
465     
466              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid ,"time"  , varid_time  ) )
467              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "EP_g"  , varid_epg   ) )
468              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "P_g"   , varid_pg    ) )
469              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "P_star", varid_pstar ) )
470              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "U_g"   , varid_ug    ) )
471              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "V_g"   , varid_vg    ) )
472              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "W_g"   , varid_wg    ) )
473              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "T_g"   , varid_tg    ) )
474              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "RO_g"   , varid_rog  ) )
475              call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, "ROP_g"  , varid_ropg ) )
476              call MFIX_check_netcdf( MFIX_nf90_get_var(ncid,varid_time,time) )
477     
478           end if
479     
480     
481     !      if (myPE .eq. PE_IO) then
482     !        call check_netcdf( nf90_get_var(ncid , varid_epg   , array1   ) )
483     !         call convert_from_io_dp(array1,array2,ijkmax2)
484     !      end if
485     !      call scatter(EP_g,array2,PE_IO)
486     
487           call readScatterRes_netcdf(EP_g   , array2, array1, ncid , varid_epg)
488           call readScatterRes_netcdf(P_g    , array2, array1, ncid , varid_pg)
489           call readScatterRes_netcdf(P_star , array2, array1, ncid , varid_pstar)
490           call readScatterRes_netcdf(ro_g   , array2, array1, ncid , varid_rog)
491           call readScatterRes_netcdf(rop_g  , array2, array1, ncid , varid_ropg)
492           call readScatterRes_netcdf(u_g    , array2, array1, ncid , varid_ug)
493           call readScatterRes_netcdf(v_g    , array2, array1, ncid , varid_vg)
494           call readScatterRes_netcdf(w_g    , array2, array1, ncid , varid_wg)
495           call readScatterRes_netcdf(t_g    , array2, array1, ncid , varid_tg)
496     
497           do i = 1,1   ! mmax
498     
499              if (myPe .eq. PE_IO) then
500                 var_name = 'U_s_xxx'
501                 write (var_name(5:7),'(i3.3)') I
502                call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_us(I)) )
503              end if
504              call readScatterRes_netcdf(U_s(:,i) , array2, array1, ncid , varid_us(i))
505     
506              if (myPe .eq. PE_IO) then
507                 var_name = 'V_s_xxx'
508                 write (var_name(5:7),'(i3.3)') I
509                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_vs(I)) )
510              end if
511              call readScatterRes_netcdf(V_s(:,i) , array2, array1, ncid , varid_vs(i))
512     
513              if (myPe .eq. PE_IO) then
514                 var_name = 'W_s_xxx'
515                 write (var_name(5:7),'(i3.3)') I
516                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_ws(I)) )
517              end if
518              call readScatterRes_netcdf(W_s(:,i) , array2, array1, ncid , varid_ws(i))
519     
520              if (myPe .eq. PE_IO) then
521                 var_name = 'ROP_s_xxx'
522                 write (var_name(7:10),'(i3.3)') I
523                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_rops(I)) )
524              end if
525              call readScatterRes_netcdf(ROP_s(:,i) , array2, array1, ncid , varid_rops(i))
526     
527              if (myPe .eq. PE_IO) then
528                 var_name = 'T_s_xxx'
529                 write (var_name(5:7),'(i3.3)') I
530                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_ts(I)) )
531              end if
532              call readScatterRes_netcdf(T_s(:,i) , array2, array1, ncid , varid_ts(i))
533     
534              if (myPe .eq. PE_IO) then
535                 var_name = 'Theta_m_xxx'
536                 write (var_name(9:11),'(i3.3)') I
537                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_thetam(I)) )
538              end if
539              call readScatterRes_netcdf(theta_m(:,i) , array2, array1, ncid , varid_thetam(i))
540     
541              if (myPe .eq. PE_IO) then
542                 var_name = 'gamaRS_xxx'
543                 write (var_name(8:10),'(i3.3)') I
544                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_gamaRS(I)) )
545              end if
546              call readScatterRes_netcdf(gama_rs(:,i) , array2, array1, ncid , varid_gamaRS(i))
547     
548              if (myPe .eq. PE_IO) then
549                 var_name = 'TRS_xxx'
550                 write (var_name(5:7),'(i3.3)') I
551                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_TRS(I)) )
552              end if
553              call readScatterRes_netcdf(T_rs(:,i) , array2, array1, ncid , varid_TRS(i))
554     
555              DO N = 1, NMAX(i)
556                 if (myPe .eq. PE_IO) then
557                    var_name = 'X_s_xxx_xxx'
558                    write (var_name(5:7) ,'(i3.3)') I
559                    write (var_name(9:11),'(i3.3)') n
560                    call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_xs(I,n)) )
561                 end if
562                 call readScatterRes_netcdf(X_s(:,i,n) , array2, array1, ncid , varid_xs(i,n))
563              END DO
564     
565           end do
566     
567           do i = 1,nmax(0)
568              if (myPe .eq. PE_IO) then
569                 var_name = 'X_g_xxx'
570                 write (var_name(5:7),'(i3.3)') I
571                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_xg(I)) )
572              end if
573              call readScatterRes_netcdf(X_g(:,i) , array2, array1, ncid , varid_xg(i))
574           end do
575     
576           do i = 1,nscalar
577              if (myPe .eq. PE_IO) then
578                 var_name = 'Scalar_xxx'
579                 write (var_name(8:10),'(i3.3)') I
580                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_scalar(I)) )
581              end if
582              call readScatterRes_netcdf(scalar(:,i) , array2, array1, ncid , varid_scalar(i))
583           end do
584     
585           do i = 1,nRR
586              if (myPe .eq. PE_IO) then
587                 var_name = 'RRates_xxx'
588                 write (var_name(8:10),'(i3.3)') I
589                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, var_name, varid_rr(I)) )
590              end if
591              call readScatterRes_netcdf(ReactionRates(:,i) , array2, array1, ncid , varid_rr(i))
592           end do
593     
594           if (k_Epsilon) then
595              if (myPe .eq. PE_IO) then
596                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, 'k_turb_g', varid_kturbg) )
597              end if
598              call readScatterRes_netcdf(k_turb_g , array2, array1, ncid , varid_kturbg)
599     
600              if (myPe .eq. PE_IO) then
601                 call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, 'e_turb_g', varid_eturbg) )
602              end if
603              call readScatterRes_netcdf(e_turb_g , array2, array1, ncid , varid_eturbg)
604           end if
605     
606          if (myPe .eq. PE_IO) then
607             call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, 'gamaRG', varid_gamaRG) )
608          end if
609          call readScatterRes_netcdf(gama_rg , array2, array1, ncid , varid_gamaRG)
610     
611          if (myPe .eq. PE_IO) then
612             call MFIX_check_netcdf( MFIX_nf90_inq_varid(ncid, 'TRG', varid_TRG) )
613          end if
614          call readScatterRes_netcdf(T_rg , array2, array1, ncid , varid_TRG)
615     
616             ! Close the file. This frees up any internal netCDF resources
617             ! associated with the file, and flushes any buffers.
618     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
619             if (myPE .eq. PE_IO) then
620                call MFIX_check_netcdf( MFIX_nf90_close(ncid) )
621            end if
622     !        call MPI_barrier(MPI_COMM_WORLD,mpierr)
623     
624       !      stop
625     
626     
627             return
628     
629           end subroutine read_res1_netcdf
630     
631     
632     
633     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
634     !                                                                      C
635     !  Module name: PATCH_AFTER_RESTART                                    C
636     !  Purpose: Patch new fluid cells after a restart                      C
637     !           This could occur when restarting with a different          C
638     !           grid partition when the RESTART file was generated         C
639     !           prior to the Dec. 4th 2014 bug fix.                        C
640     !                                                                      C
641     !  Author: Jeff Dietiker                              Date: 14-APR-15  C
642     !                                                                      C
643     !  Revision Number:                                                    C
644     !  Purpose:                                                            C
645     !  Author:                                            Date: dd-mmm-yy  C
646     !  Reviewer:                                          Date: dd-mmm-yy  C
647     !                                                                      C
648     !  Literature/Document References:                                     C
649     !                                                                      C
650     !                                                                      C
651     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
652     
653           SUBROUTINE PATCH_AFTER_RESTART
654     !...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
655     !...Switches: -xf
656     !
657     !-----------------------------------------------
658     !   M o d u l e s
659     !-----------------------------------------------
660           USE param
661           USE param1
662           USE fldvar
663           USE geometry
664           USE physprop
665           USE run
666           USE rxns
667           USE scalars
668           USE funits
669           USE energy
670           USE compar
671           USE cdist
672           USE mpi_utility
673           USE sendrecv
674           USE cutcell
675           use functions
676     
677           IMPLICIT NONE
678     !-----------------------------------------------
679     !   G l o b a l   P a r a m e t e r s
680     !-----------------------------------------------
681     !-----------------------------------------------
682     !   L o c a l   P a r a m e t e r s
683     !-----------------------------------------------
684     !-----------------------------------------------
685     !   L o c a l   V a r i a b l e s
686     !-----------------------------------------------
687     !
688     !
689     !-----------------------------------------------
690     ! Local variables
691     !-----------------------------------------------
692     ! indices
693           INTEGER :: I,J,K, IJK, IJKNB
694           INTEGER :: M,N
695           INTEGER :: NB
696           INTEGER, DIMENSION(6) :: NBCELL
697           LOGICAL :: NB_FOUND
698     
699     
700     !-----------------------------------------------
701     
702           DO IJK = ijkstart3, ijkend3
703     
704              IF (FLUID_AT(IJK).AND.EP_G(IJK)==UNDEFINED) THEN
705     
706     ! Detects new fluid cells that used to be blocked cells with undefined
707     ! values. When a fluid cell has undefined void fraction, this means all
708     ! variables need to be patched. Typically, this will be a fairly small
709     ! cut cell that was flagged as BLOCKED cell with a different partition.
710     ! If a valid fluid cell is found next to this undefined cell, all field
711     ! variables will be copied over. If no valid fluid cell is found, the
712     ! code will continue and will likely stop during the check_data_30
713     ! (zero species mass fractions will yield a zero specific heat).
714                 I = I_OF(IJK)
715                 J = J_OF(IJK)
716                 K = K_OF(IJK)
717     
718                 NBCELL(1) = IM_OF(IJK)
719                 NBCELL(2) = JM_OF(IJK)
720                 NBCELL(3) = KM_OF(IJK)
721                 NBCELL(4) = IP_OF(IJK)
722                 NBCELL(5) = JP_OF(IJK)
723                 NBCELL(6) = KP_OF(IJK)
724     
725                 NB_FOUND = .FALSE.
726     
727                 DO NB = 1,6
728     
729                    IJKNB = NBCELL(NB)
730     
731                    IF(FLUID_AT(IJKNB).AND.EP_G(IJKNB)/=UNDEFINED) THEN
732                       NB_FOUND = .TRUE.
733                       WRITE (*, 1010) MyPE, I,J,K
734     
735                       EP_G(IJK)   = EP_G(IJKNB)
736                       P_G(IJK)    = P_G(IJKNB)
737                       P_STAR(IJK) = P_STAR(IJKNB)
738                       RO_G(IJK)   = RO_G(IJKNB)
739                       ROP_G(IJK)  = ROP_G(IJKNB)
740                       T_G(IJK)    = T_G(IJKNB)
741     
742                       T_s(IJK,1:MMAX) = T_s(IJKNB,1:MMAX)
743     
744                       U_s(IJK,1:MMAX) = U_s(IJKNB,1:MMAX)
745                       V_s(IJK,1:MMAX) = V_s(IJKNB,1:MMAX)
746                       W_s(IJK,1:MMAX) = W_s(IJKNB,1:MMAX)
747     
748                       X_g(IJK,1:NMAX(0)) = X_g(IJKNB,1:NMAX(0))
749     
750                       U_G(IJK) = U_G(IJKNB)
751                       V_G(IJK) = V_G(IJKNB)
752                       W_G(IJK) = W_G(IJKNB)
753     
754                       ROP_S(IJK,1:MMAX) = ROP_S(IJKNB,1:MMAX)
755     
756                       IF(ANY(SOLVE_ROs)) RO_S(IJK,1:MMAX) = RO_S(IJKNB,1:MMAX)
757     
758     
759                       THETA_M(IJK,1:MMAX) = THETA_M(IJKNB,1:MMAX)
760     
761                       DO M = 1,MMAX
762                          DO N = 1, NMAX(M)
763                             X_S(IJK,M,N)= X_S(IJKNB,M,N)
764                          ENDDO
765                       ENDDO
766     
767     
768                       DO N = 1, NScalar
769                          Scalar(IJK,N) = Scalar(IJKNB,N)
770                       END DO
771     
772                       GAMA_RG(IJK) = GAMA_RG(IJKNB)
773                       T_RG(IJK)    = T_RG(IJKNB)
774     
775                       GAMA_RS(IJK,1:MMAX) = GAMA_RS(IJKNB,1:MMAX)
776                       T_RS(IJK,1:MMAX)    = T_RS(IJKNB,1:MMAX)
777     
778     
779                       DO N = 1, nRR
780                          ReactionRates(IJK,N) = ReactionRates(IJKNB,N)
781                       END DO
782     
783                       IF (K_Epsilon) THEN
784                          K_Turb_G(IJK) = K_Turb_G(IJKNB)
785                          E_Turb_G(IJK) = E_Turb_G(IJKNB)
786                       ENDIF
787     
788                       EXIT ! Exit as soon as first valid neighbor cell is found
789                    ENDIF  ! NB is a fluid cell
790     
791                 ENDDO ! NB Loop
792     
793                 IF(.NOT.NB_FOUND) WRITE (*, 1020) MyPE, I,J,K   ! NO FLUID CELL AMONG NEIGBHORS
794     
795              ENDIF ! New fuid cell
796     
797           ENDDO ! IJK loop
798     
799     1010  FORMAT(1X,'PATCHING NEW FLUID CELL UPON RESTART: MyPE,I,J,K =' ,I6,I6,I6,I6)
800     1020  FORMAT(1X,'UNABLE TO PATCH NEW FLUID CELL UPON RESTART: MyPE,I,J,K =' ,I6,I6,I6,I6)
801           END SUBROUTINE PATCH_AFTER_RESTART
802