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

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Subroutine: READ_RES0                                               C
4     !  Purpose: read the initial restart records (namelist data)           C
5     !                                                                      C
6     !  Author: P. Nicoletti                               Date: 13-DEC-91  C
7     !  Reviewer: P. Nicoletti, W. Rogers, M. Syamlal      Date:            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     !                                                                      C
15     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
16     
17           SUBROUTINE READ_RES0
18     
19     !-----------------------------------------------
20     ! Modules
21     !-----------------------------------------------
22           USE param
23           USE param1
24           USE geometry
25           USE physprop
26           USE run
27           USE ic
28           USE bc
29           USE is
30           USE constant
31           USE funits
32           USE output
33           USE scales
34           USE ur_facs
35           USE toleranc
36           USE leqsol
37           USE scalars
38           USE rxns
39           USE compar
40           USE mpi_utility
41           USE fldvar
42           USE stiff_chem
43           USE in_binary_512
44           USE in_binary_512i
45     
46           IMPLICIT NONE
47     !-----------------------------------------------
48     ! Local variables
49     !-----------------------------------------------
50     ! loop counters
51           INTEGER    LC, L , N, M
52     ! Pointer to the next record
53           INTEGER    NEXT_RECA
54     ! file version id
55           CHARACTER(LEN=512) :: VERSION
56     ! version number
57           REAL       VERSION_NUMBER
58     ! Temporary arrays
59           DOUBLE PRECISION IC_Tmp(DIMENSION_IC), BC_Tmp(DIMENSION_BC)
60     !
61           INTEGER    DIM_IC , DIM_BC , DIM_C , DIM_IS
62     ! Read an array dimension
63           INTEGER :: DIM_tmp
64     !
65           logical :: doingPost
66           integer :: n_spx_res
67     
68     ! Place holder for deprecated keywords:
69           LOGICAL :: CALL_ISAT
70     
71     ! declare global scratch arrays:
72     ! declare integer Global SCRatch array
73           INTEGER, ALLOCATABLE, DIMENSION(:) :: IGTEMP1, iGTEMP2
74     ! declare real*4 Global SCRatch array
75           REAL, ALLOCATABLE, DIMENSION(:) :: rGTEMP
76     ! declare real*8 Global SCRatch array
77           DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: dGTEMP
78     ! declaring following arrays to pack scalar variables when
79     ! BCASTing to reduce the number of BCAST calls:
80     ! packing array for integers
81           INTEGER, ALLOCATABLE, DIMENSION(:) :: INTPACK
82     ! packing array for doubles
83           DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: DBLPACK
84     !-----------------------------------------------
85     
86           doingPost = .false.
87     
88     !  1) Check to ensure that this subroutine was updated.
89     !  2) Initialize missing constants from earlier versions.
90     !  3) Add new read statements at the end of the file.
91     
92     ! only PE_IO reads the restart file
93         if (myPE == PE_IO ) then
94           READ (UNIT_RES, REC=1) VERSION
95           READ (VERSION(6:512), *) VERSION_NUMBER
96     
97           IF (VERSION_NUMBER > 1.8) THEN
98              WRITE (*, *) ' Update Subroutine read_res0'
99              CALL SLUMBER
100     !         STOP
101              call exitMPI(myPE)  ! Abort all PEs, not only the current one
102           ENDIF
103         endif
104     
105     ! Initialize required constants missing from earlier versions
106           P_REF = ZERO
107           P_SCALE = ONE
108           DIM_IC = 5
109           DIM_BC = 5
110           DIM_C = 5
111           DIM_IS = 5
112           C_E = 1.0
113           C_F = 0.0
114           PHI = 0.0
115           PHI_W = 0.0
116     
117     ! only PE_IO reads the restart
118         if (myPE == PE_IO ) then
119           READ (UNIT_RES, REC=2) RUN_NAME, ID_MONTH, ID_DAY, ID_YEAR, ID_HOUR, &
120              ID_MINUTE, ID_SECOND
121           READ (UNIT_RES, REC=3) NEXT_RECA
122     
123           IF (VERSION == 'RES = 01.00') THEN
124              READ (UNIT_RES, REC=4) IMIN1, JMIN1, KMIN1, IMAX, JMAX, KMAX, IMAX1, &
125                 JMAX1, KMAX1, IMAX2, JMAX2, KMAX2, IJMAX2, IJKMAX2, MMAX, DT, &
126                 XLENGTH, YLENGTH, ZLENGTH
127           ELSEIF (VERSION=='RES = 01.01' .OR. VERSION=='RES = 01.02') THEN
128              READ (UNIT_RES, REC=4) IMIN1, JMIN1, KMIN1, IMAX, JMAX, KMAX, IMAX1, &
129                 JMAX1, KMAX1, IMAX2, JMAX2, KMAX2, IJMAX2, IJKMAX2, MMAX, DIM_IC, &
130                 DIM_BC, DT, XLENGTH, YLENGTH, ZLENGTH
131           ELSEIF (VERSION == 'RES = 01.03') THEN
132              READ (UNIT_RES, REC=4) IMIN1, JMIN1, KMIN1, IMAX, JMAX, KMAX, IMAX1, &
133                 JMAX1, KMAX1, IMAX2, JMAX2, KMAX2, IJMAX2, IJKMAX2, MMAX, DIM_IC, &
134                 DIM_BC, DT, XMIN, XLENGTH, YLENGTH, ZLENGTH
135           ELSEIF (VERSION == 'RES = 01.04') THEN
136              READ (UNIT_RES, REC=4) IMIN1, JMIN1, KMIN1, IMAX, JMAX, KMAX, IMAX1, &
137                 JMAX1, KMAX1, IMAX2, JMAX2, KMAX2, IJMAX2, IJKMAX2, MMAX, DIM_IC, &
138                 DIM_BC, DIM_C, DT, XMIN, XLENGTH, YLENGTH, ZLENGTH
139           ELSEIF (VERSION == 'RES = 01.05') THEN
140              READ (UNIT_RES, REC=4) IMIN1, JMIN1, KMIN1, IMAX, JMAX, KMAX, IMAX1, &
141                 JMAX1, KMAX1, IMAX2, JMAX2, KMAX2, IJMAX2, IJKMAX2, MMAX, DIM_IC, &
142                 DIM_BC, DIM_C, DIM_IS, DT, XMIN, XLENGTH, YLENGTH, ZLENGTH
143           ELSE
144              READ (UNIT_RES, REC=4) IMIN1, JMIN1, KMIN1, IMAX, JMAX, KMAX, IMAX1, &
145                 JMAX1, KMAX1, IMAX2, JMAX2, KMAX2, IJMAX2, IJKMAX2, MMAX, DIM_IC, &
146                 DIM_BC, DIM_C, DIM_IS, DT, XMIN, XLENGTH, YLENGTH, ZLENGTH, C_E, &
147                 C_F, PHI, PHI_W
148           ENDIF
149         endif
150     
151         Allocate( INTPACK(30))  ! ALLOCate packing array
152         Allocate( DBLPACK(20))  ! ALLOCate packing array
153         INTPACK = 0
154         DBLPACK = 0.d0
155     
156     ! Root PE (PE_IO) : Pack and broadcast ; others unpack
157     ! WARNING: This implementation assumes VERSION > 1.05, filtering needed for
158     !          earlier versions
159         if (myPE == PE_IO ) then
160            call bcast(VERSION_NUMBER, PE_IO)  ! BCAST0r
161            INTPACK(1) = ID_MONTH
162            INTPACK(2) = ID_DAY
163            INTPACK(3) = ID_YEAR
164            INTPACK(4) = ID_HOUR
165            INTPACK(5) = ID_MINUTE
166            INTPACK(6) = ID_SECOND
167            INTPACK(7) = IMIN1
168            INTPACK(8) = JMIN1
169            INTPACK(9) = KMIN1
170            INTPACK(10) = IMAX
171            INTPACK(11) = JMAX
172            INTPACK(12) = KMAX
173            INTPACK(13) = IMAX1
174            INTPACK(14) = JMAX1
175            INTPACK(15) = KMAX1
176            INTPACK(16) = IMAX2
177            INTPACK(17) = JMAX2
178            INTPACK(18) = KMAX2
179            INTPACK(19) = IJMAX2
180            INTPACK(20) = IJKMAX2
181            INTPACK(21) = MMAX
182            INTPACK(22) = DIM_IC
183            INTPACK(23) = DIM_BC
184            INTPACK(24) = DIM_C
185            INTPACK(25) = DIM_IS
186     ! In spite of the overhead, using MPI_Pack/MPI_Unpack may be worth using here
187            call bcast(INTPACK(1:25),PE_IO)  ! BCAST1i
188            DBLPACK(1) = DT
189            DBLPACK(2) = XMIN
190            DBLPACK(3) = XLENGTH
191            DBLPACK(4) = YLENGTH
192            DBLPACK(5) = ZLENGTH
193            DBLPACK(6) = C_E
194            DBLPACK(7) = C_F
195            DBLPACK(8) = PHI
196            DBLPACK(9) = PHI_W
197            call bcast(DBLPACK,PE_IO) ! BCAST1d
198            call bcast(VERSION,PE_IO) ! BCAST0c
199            call bcast(RUN_NAME,PE_IO) ! BCAST0c
200         else
201            call bcast(VERSION_NUMBER, PE_IO)  ! BCAST0r
202            call bcast(INTPACK(1:25),PE_IO)    ! BCAST1i (receive)
203            ID_MONTH = INTPACK(1)
204            ID_DAY = INTPACK(2)
205            ID_YEAR = INTPACK(3)
206            ID_HOUR = INTPACK(4)
207            ID_MINUTE = INTPACK(5)
208            ID_SECOND = INTPACK(6)
209            IMIN1 = INTPACK(7)
210            JMIN1 = INTPACK(8)
211            KMIN1 = INTPACK(9)
212            IMAX = INTPACK(10)
213            JMAX = INTPACK(11)
214            KMAX = INTPACK(12)
215            IMAX1 = INTPACK(13)
216            JMAX1 = INTPACK(14)
217            KMAX1 = INTPACK(15)
218            IMAX2 = INTPACK(16)
219            JMAX2 = INTPACK(17)
220            KMAX2 = INTPACK(18)
221            IJMAX2 = INTPACK(19)
222            IJKMAX2 = INTPACK(20)
223            MMAX = INTPACK(21)
224            DIM_IC = INTPACK(22)
225            DIM_BC = INTPACK(23)
226            DIM_C = INTPACK(24)
227            DIM_IS = INTPACK(25)
228            call bcast(DBLPACK,PE_IO)  ! BCAST1d (recv)
229            DT = DBLPACK(1)
230            XMIN = DBLPACK(2)
231            XLENGTH = DBLPACK(3)
232            YLENGTH = DBLPACK(4)
233            ZLENGTH = DBLPACK(5)
234            C_E = DBLPACK(6)
235            C_F = DBLPACK(7)
236            PHI = DBLPACK(8)
237            PHI_W = DBLPACK(9)
238            call bcast(VERSION,PE_IO)  ! BCAST0c (recv)
239            call bcast(RUN_NAME,PE_IO) ! BCAST0c (recv)
240         endif
241     
242         DeAllocate (INTPACK)
243     
244     
245     
246     ! CHECK DIMENSIONS
247     
248           IF (.NOT. ( (DIM_IC <= DIMENSION_IC)  .AND. &
249                       (DIM_BC <= DIMENSION_BC)  .AND. &
250                       (DIM_C <= DIMENSION_C)  .AND. &
251                       (DIM_IS <= DIMENSION_IS) )) GOTO 900
252     
253           IF (MMAX + 1 > 0) THEN
254              NMAX(:MMAX) = 1
255           ENDIF
256     
257           NEXT_RECA = 5
258     
259     
260           IF (VERSION_NUMBER >= 1.04) THEN
261     
262              if (myPE == PE_IO) then
263                 CALL IN_BIN_512 (UNIT_RES, C, DIM_C, NEXT_RECA)
264     ! work around for -O3 compiler bug
265                  NEXT_RECA = 1 + NEXT_RECA
266                  NEXT_RECA = NEXT_RECA - 1
267              endif
268              call bcast(C,PE_IO)   ! BCAST1d user defined constants,C
269     
270              if (myPE == PE_IO) then
271                 DO LC = 1, DIM_C
272                    READ (UNIT_RES, REC=NEXT_RECA) C_NAME(LC)
273                    NEXT_RECA = NEXT_RECA + 1
274                 ENDDO
275              endif
276              call bcast(C_NAME,PE_IO)  ! BCAST1c user defined constant names
277     
278              IF (myPE == PE_IO) THEN
279                 IF (VERSION_NUMBER < 1.12) THEN
280                    CALL IN_BIN_512I (UNIT_RES, NMAX, MMAX + 1, NEXT_RECA)
281                 ELSE
282                    READ (UNIT_RES, REC=NEXT_RECA) (NMAX(L),L=0, MMAX)
283                    NEXT_RECA = NEXT_RECA + 1
284                 ENDIF
285              ENDIF
286              call bcast(NMAX,PE_IO)   !BCAST1i total # of gas OR solid species
287           ENDIF
288     
289     ! The following occurs when mfix.dat has not been read and the
290     ! dimensions are from the .RES file.  Note that check on NMAX for
291     ! solids is not done.
292           IF ( .NOT. ( (IMAX2 <= DIMENSION_I)    .AND. &
293                        (JMAX2 <= DIMENSION_J)    .AND. &
294                        (KMAX2 <= DIMENSION_K)    .AND. &
295                        (IJKMAX2 <= DIMENSION_3)  .AND. &
296                        (MMAX <= DIMENSION_M)     .AND. &
297                        (NMAX(0) <= DIMENSION_N_G) ) ) then
298     
299              IF(IMAX2 == 1)NO_I=.TRUE.
300              IF(JMAX2 == 1)NO_J=.TRUE.
301              IF(KMAX2 == 1)NO_K=.TRUE.
302     
303     ! This is a "fix" for post_mfix ... assumes only 1 processor
304              write (*,*) ' '
305              write (*,*) ' ********************************************************'
306              write (*,*) ' '
307              write (*,*) ' read_res0 : code valid for running on 1 processor only'
308              write (*,*) ' '
309              write (*,*) ' ********************************************************'
310              write (*,*) ' '
311              IMAX3 = imax2
312              JMAX3 = jmax2
313              KMAX3 = kmax2
314              kend3 = kmax2
315              kstart3 = 1
316              jend3 = jmax2
317              jstart3 = 1
318              iend3 = imax2
319              istart3 = 1
320              ijkmax3 = imax3*jmax3*kmax3
321              allocate (ijksize3_all(0:1))
322              ijksize3_all(:) = ijkmax3
323     
324              nScalar = 0                ! since NScalar
325              nRR     = 0                ! and NRR not read
326              doingPost = .true.         ! until later
327     
328              call set_parameters
329              call allocate_arrays_geometry
330              call allocate_arrays_increments
331              call allocate_arrays       ! do for mfix/post_mfix
332              deallocate(ijksize3_all)   ! post_mfix "fix"
333     
334           ENDIF
335     
336           if (myPE == PE_IO) then
337              CALL IN_BIN_512 (UNIT_RES, DX, IMAX2, NEXT_RECA)
338              CALL IN_BIN_512 (UNIT_RES, DY, JMAX2, NEXT_RECA)
339              CALL IN_BIN_512 (UNIT_RES, DZ, KMAX2, NEXT_RECA)
340           endif
341     
342     
343     ! the RES file for version <= 1.4 write out : dx(1) to dx(imax2)
344     ! but reads starting at dx(0). Need to shift the arrays for
345     ! proper pos-processing etc.
346           if (version_number < 1.41) then
347              do L = imax2,1,-1
348                 dx(L) = dx(L-1)
349              end do
350              do L = jmax2,1,-1
351                 dy(L) = dy(L-1)
352              end do
353              do L = kmax2,1,-1
354                 dz(L) = dz(L-1)
355              end do
356           endif
357           call bcast(dx, PE_IO)   !//PAR_I/O BCAST1d
358           call bcast(dy, PE_IO)   !//PAR_I/O BCAST1d
359           call bcast(dz, PE_IO)   !//PAR_I/O BCAST1d
360     
361     
362           IF (myPE == PE_IO) THEN
363              READ (UNIT_RES, REC=NEXT_RECA) RUN_NAME, &
364                    DESCRIPTION, UNITS, RUN_TYPE, &
365                    COORDINATES
366              NEXT_RECA = NEXT_RECA + 1
367           ENDIF
368           call bcast(RUN_NAME,PE_IO)    ! BCAST0c
369           call bcast(DESCRIPTION,PE_IO) ! BCAST0c
370           call bcast(UNITS,PE_IO)       ! BCAST0c
371           call bcast(RUN_TYPE,PE_IO)    ! BCAST0c
372           call bcast(COORDINATES,PE_IO) ! BCAST0c
373     
374           IF (myPE == PE_IO) THEN
375              IF (VERSION=='RES = 01.00' .OR. VERSION=='RES = 01.01') THEN
376                 READ (UNIT_RES, REC=NEXT_RECA) (D_P0(L),L=1,&
377                    MMAX), (RO_S0(L),L=1,MMAX), EP_STAR, &
378                    MU_G0, MW_AVG
379              ELSEIF (VERSION == 'RES = 01.02') THEN
380                 READ (UNIT_RES, REC=NEXT_RECA) (D_P0(L),L=1,&
381                    MMAX), (RO_S0(L),L=1,MMAX), EP_STAR, &
382                    RO_G0, MU_G0, MW_AVG
383              ELSEIF (VERSION == 'RES = 01.03') THEN
384                 READ (UNIT_RES, REC=NEXT_RECA) (D_P0(L),L=1,&
385                     MMAX), (RO_S0(L),L=1,MMAX), EP_STAR, &
386                     RO_G0, MU_G0, MW_AVG
387              ELSEIF (VERSION_NUMBER >= 1.04) THEN
388                 READ (UNIT_RES, REC=NEXT_RECA) (D_P0(L),L=1,&
389                     MMAX), (RO_S0(L),L=1,MMAX), EP_STAR, &
390                     RO_G0, MU_G0, MW_AVG
391              ENDIF
392              NEXT_RECA = NEXT_RECA + 1
393           ENDIF
394           call bcast(D_P0, PE_IO)      ! BCAST1d
395           call bcast(RO_S0, PE_IO)     ! BCAST1d
396           call bcast(EP_STAR, PE_IO)   ! BCAST0d
397           call bcast(RO_G0, PE_IO)     ! BCAST0d
398           call bcast(MU_G0, PE_IO)     ! BCAST0d
399           call bcast(MW_AVG, PE_IO)    ! BCAST0d
400     
401           IF (VERSION_NUMBER >= 1.04) THEN
402              IF (myPE == PE_IO) THEN
403                 CALL IN_BIN_512 (UNIT_RES, MW_G, NMAX(0), NEXT_RECA)
404                 DO LC = 1, MMAX
405                    READ (UNIT_RES, REC=NEXT_RECA) (MW_S(LC,N),&
406                       N=1,NMAX(LC))
407                    NEXT_RECA = NEXT_RECA + 1
408                 ENDDO
409              ENDIF
410              call bcast(MW_G, PE_IO)     ! BCAST1d
411              call bcast(MW_S, PE_IO)     ! BCAST2d
412           ENDIF
413     
414           if (myPE == PE_IO) then
415              CALL IN_BIN_512 (UNIT_RES, IC_X_W, DIM_IC, NEXT_RECA)
416              CALL IN_BIN_512 (UNIT_RES, IC_X_E, DIM_IC, NEXT_RECA)
417              CALL IN_BIN_512 (UNIT_RES, IC_Y_S, DIM_IC, NEXT_RECA)
418              CALL IN_BIN_512 (UNIT_RES, IC_Y_N, DIM_IC, NEXT_RECA)
419           endif
420           call bcast(IC_X_W, PE_IO)     ! BCAST1d
421           call bcast(IC_X_E, PE_IO)     ! BCAST1d
422           call bcast(IC_Y_S, PE_IO)     ! BCAST1d
423           call bcast(IC_Y_N, PE_IO)     ! BCAST1d
424     
425           if (myPE == PE_IO)then
426              CALL IN_BIN_512 (UNIT_RES, IC_Z_B, DIM_IC, NEXT_RECA)
427              CALL IN_BIN_512 (UNIT_RES, IC_Z_T, DIM_IC, NEXT_RECA)
428           endif
429           call bcast(IC_Z_B, PE_IO)     ! BCAST1d
430           call bcast(IC_Z_T, PE_IO)     ! BCAST1d
431     
432           if (myPE == PE_IO)then
433              CALL IN_BIN_512I (UNIT_RES, IC_I_W, DIM_IC, NEXT_RECA)
434              CALL IN_BIN_512I (UNIT_RES, IC_I_E, DIM_IC, NEXT_RECA)
435              CALL IN_BIN_512I (UNIT_RES, IC_J_S, DIM_IC, NEXT_RECA)
436              CALL IN_BIN_512I (UNIT_RES, IC_J_N, DIM_IC, NEXT_RECA)
437           endif
438           call bcast(IC_I_W, PE_IO)     ! BCAST1i
439           call bcast(IC_I_E, PE_IO)     ! BCAST1i
440           call bcast(IC_J_S, PE_IO)     ! BCAST1i
441           call bcast(IC_J_N, PE_IO)     ! BCAST1i
442     
443           if (myPE == PE_IO)then
444              CALL IN_BIN_512I (UNIT_RES, IC_K_B, DIM_IC, NEXT_RECA)
445              CALL IN_BIN_512I (UNIT_RES, IC_K_T, DIM_IC, NEXT_RECA)
446           endif
447           call bcast(IC_K_B, PE_IO)     ! BCAST1i
448           call bcast(IC_K_T, PE_IO)     ! BCAST1i
449     
450           if (myPE == PE_IO) then
451              CALL IN_BIN_512 (UNIT_RES, IC_EP_G, DIM_IC, NEXT_RECA)
452              CALL IN_BIN_512 (UNIT_RES, IC_P_G, DIM_IC, NEXT_RECA)
453              CALL IN_BIN_512 (UNIT_RES, IC_T_G, DIM_IC, NEXT_RECA)
454           endif
455           call bcast(IC_EP_G, PE_IO)    ! BCAST1d
456           call bcast(IC_P_G, PE_IO)     ! BCAST1d
457           call bcast(IC_T_G, PE_IO)     ! BCAST1d
458     
459           IF (VERSION_NUMBER < 1.15) THEN
460              IF (myPE == PE_IO) THEN
461                 CALL IN_BIN_512 (UNIT_RES, IC_T_S(1,1), DIM_IC, NEXT_RECA)
462                 IF (MMAX >= 2) THEN
463                    CALL IN_BIN_512 (UNIT_RES, IC_T_S(1,2), DIM_IC, NEXT_RECA)
464                 ELSE
465                    CALL IN_BIN_512 (UNIT_RES, IC_TMP, DIM_IC, NEXT_RECA)
466                 ENDIF
467              ENDIF
468                 call bcast(IC_T_S, PE_IO)    ! BCAST2d
469           ENDIF
470     
471           IF (VERSION_NUMBER >= 1.04) THEN
472              IF (myPE == PE_IO) THEN
473                 DO N = 1, NMAX(0)
474                    CALL IN_BIN_512 (UNIT_RES, IC_X_G(1,N), DIM_IC, &
475                       NEXT_RECA)
476                 ENDDO
477              ENDIF
478              call bcast(IC_X_G, PE_IO)    ! BCAST2d
479           ENDIF
480     
481           if (myPE == PE_IO)then
482              CALL IN_BIN_512 (UNIT_RES, IC_U_G, DIM_IC, NEXT_RECA)
483              CALL IN_BIN_512 (UNIT_RES, IC_V_G, DIM_IC, NEXT_RECA)
484              CALL IN_BIN_512 (UNIT_RES, IC_W_G, DIM_IC, NEXT_RECA)
485           ENDIF
486           call bcast(IC_U_G, PE_IO)    ! BCAST1d
487           call bcast(IC_V_G, PE_IO)    ! BCAST1d
488           call bcast(IC_W_G, PE_IO)    ! BCAST1d
489     
490           IF (myPE == PE_IO) THEN
491              DO LC = 1, MMAX
492                 CALL IN_BIN_512 (UNIT_RES, IC_ROP_S(1,LC), &
493                      DIM_IC, NEXT_RECA)
494                 CALL IN_BIN_512 (UNIT_RES, IC_U_S(1,LC), &
495                      DIM_IC, NEXT_RECA)
496                 CALL IN_BIN_512 (UNIT_RES, IC_V_S(1,LC), &
497                      DIM_IC, NEXT_RECA)
498                 CALL IN_BIN_512 (UNIT_RES, IC_W_S(1,LC), &
499                      DIM_IC, NEXT_RECA)
500                 IF (VERSION_NUMBER >= 1.15) THEN
501                    CALL IN_BIN_512(UNIT_RES, IC_T_S(1,LC), DIM_IC, &
502                      NEXT_RECA)
503                 ENDIF
504     
505                 IF (VERSION_NUMBER >= 1.04) THEN
506                    DO N = 1, NMAX(LC)
507                      CALL IN_BIN_512 (UNIT_RES, IC_X_S(1,LC,N), &
508                            DIM_IC, NEXT_RECA)
509                    ENDDO
510                 ENDIF
511              ENDDO
512           ENDIF
513           call bcast(IC_ROP_S, PE_IO)  ! BCAST2d
514           call bcast(IC_U_S, PE_IO)    ! BCAST2d
515           call bcast(IC_V_S, PE_IO)    ! BCAST2d
516           call bcast(IC_W_S, PE_IO)    ! BCAST2d
517     
518           if (VERSION_NUMBER >= 1.15) call bcast(IC_T_S, PE_IO)  ! BCAST2d
519           if (VERSION_NUMBER >= 1.04) call bcast(IC_X_S, PE_IO)  ! BCAST3d
520     
521           if (myPE == PE_IO) then
522              CALL IN_BIN_512 (UNIT_RES, BC_X_W, DIM_BC, NEXT_RECA)
523              CALL IN_BIN_512 (UNIT_RES, BC_X_E, DIM_BC, NEXT_RECA)
524              CALL IN_BIN_512 (UNIT_RES, BC_Y_S, DIM_BC, NEXT_RECA)
525              CALL IN_BIN_512 (UNIT_RES, BC_Y_N, DIM_BC, NEXT_RECA)
526           endif
527           call bcast(BC_X_W, PE_IO)    ! BCAST1d
528           call bcast(BC_X_E, PE_IO)    ! BCAST1d
529           call bcast(BC_Y_S, PE_IO)    ! BCAST1d
530           call bcast(BC_Y_N, PE_IO)    ! BCAST1d
531     
532           if (myPE == PE_IO) then
533              CALL IN_BIN_512 (UNIT_RES, BC_Z_B, DIM_BC, NEXT_RECA)
534              CALL IN_BIN_512 (UNIT_RES, BC_Z_T, DIM_BC, NEXT_RECA)
535           endif
536           call bcast(BC_Z_B, PE_IO)    ! BCAST1d
537           call bcast(BC_Z_T, PE_IO)    ! BCAST1d
538     
539           IF (myPE == PE_IO) THEN
540              CALL IN_BIN_512I (UNIT_RES, BC_I_W, DIM_BC, NEXT_RECA)
541               CALL IN_BIN_512I (UNIT_RES, BC_I_E, DIM_BC, NEXT_RECA)
542               CALL IN_BIN_512I (UNIT_RES, BC_J_S, DIM_BC, NEXT_RECA)
543               CALL IN_BIN_512I (UNIT_RES, BC_J_N, DIM_BC, NEXT_RECA)
544           ENDIF
545           call bcast(BC_I_W, PE_IO)    ! BCAST1i
546           call bcast(BC_I_E, PE_IO)    ! BCAST1i
547           call bcast(BC_J_S, PE_IO)    ! BCAST1i
548           call bcast(BC_J_N, PE_IO)    ! BCAST1i
549     
550           IF (myPE == PE_IO) THEN
551                CALL IN_BIN_512I (UNIT_RES, BC_K_B, DIM_BC, NEXT_RECA)
552                CALL IN_BIN_512I (UNIT_RES, BC_K_T, DIM_BC, NEXT_RECA)
553           ENDIF
554           call bcast(BC_K_B, PE_IO)    ! BCAST1i
555           call bcast(BC_K_T, PE_IO)    ! BCAST1i
556     
557           IF (myPE == PE_IO) THEN
558              CALL IN_BIN_512 (UNIT_RES, BC_EP_G, DIM_BC, NEXT_RECA)
559              CALL IN_BIN_512 (UNIT_RES, BC_P_G, DIM_BC, NEXT_RECA)
560              CALL IN_BIN_512 (UNIT_RES, BC_T_G, DIM_BC, NEXT_RECA)
561           ENDIF
562           call bcast(BC_EP_G, PE_IO)   ! BCAST1d
563           call bcast(BC_P_G, PE_IO)    ! BCAST1d
564           call bcast(BC_T_G, PE_IO)    ! BCAST1d
565     
566           IF (VERSION_NUMBER < 1.15) THEN
567              IF (myPE == PE_IO) THEN
568                 CALL IN_BIN_512 (UNIT_RES, BC_T_S(1,1), DIM_BC, NEXT_RECA)
569                 IF (MMAX >= 2) THEN
570                    CALL IN_BIN_512 (UNIT_RES, BC_T_S(1,2), DIM_BC, NEXT_RECA)
571                 ELSE
572     ! dummy read, no need to broadcast
573                    CALL IN_BIN_512 (UNIT_RES, BC_TMP, DIM_BC, NEXT_RECA)
574                 ENDIF
575              ENDIF
576              call bcast(BC_T_S, PE_IO)   ! BCAST2d
577           ENDIF
578     
579           IF (VERSION_NUMBER >= 1.04) THEN
580              IF (myPE == PE_IO) THEN
581                 DO N = 1, NMAX(0)
582                    CALL IN_BIN_512 (UNIT_RES, BC_X_G(1,N), &
583                         DIM_BC, NEXT_RECA)
584                 ENDDO
585              ENDIF
586              call bcast(BC_X_G, PE_IO)   ! BCAST2d
587           ENDIF
588     
589           IF (myPE == PE_IO) THEN
590              CALL IN_BIN_512 (UNIT_RES, BC_U_G, DIM_BC, NEXT_RECA)
591              CALL IN_BIN_512 (UNIT_RES, BC_V_G, DIM_BC, NEXT_RECA)
592              CALL IN_BIN_512 (UNIT_RES, BC_W_G, DIM_BC, NEXT_RECA)
593              CALL IN_BIN_512 (UNIT_RES, BC_RO_G, DIM_BC, NEXT_RECA)
594              CALL IN_BIN_512 (UNIT_RES, BC_ROP_G, DIM_BC, NEXT_RECA)
595              CALL IN_BIN_512 (UNIT_RES, BC_VOLFLOW_G, DIM_BC, NEXT_RECA)
596              CALL IN_BIN_512 (UNIT_RES, BC_MASSFLOW_G, DIM_BC, NEXT_RECA)
597           ENDIF
598           call bcast(BC_U_G, PE_IO)       ! BCAST1d
599           call bcast(BC_V_G, PE_IO)       ! BCAST1d
600           call bcast(BC_W_G, PE_IO)       ! BCAST1d
601           call bcast(BC_RO_G, PE_IO)      ! BCAST1d
602           call bcast(BC_ROP_G, PE_IO)     ! BCAST1d
603           call bcast(BC_VOLFLOW_G, PE_IO) ! BCAST1d
604           call bcast(BC_MASSFLOW_G, PE_IO)! BCAST1d
605     
606           IF (myPE == PE_IO) THEN
607              DO LC = 1, MMAX
608                 CALL IN_BIN_512 (UNIT_RES, BC_ROP_S(1,LC), DIM_BC, &
609                    NEXT_RECA)
610                 CALL IN_BIN_512 (UNIT_RES, BC_U_S(1,LC), DIM_BC, &
611                    NEXT_RECA)
612                 CALL IN_BIN_512 (UNIT_RES, BC_V_S(1,LC), DIM_BC, &
613                    NEXT_RECA)
614     
615     ! Note : previous versions did not write out BC_W_s, X_S
616                 IF (VERSION_NUMBER >= 1.04) THEN
617                    CALL IN_BIN_512 (UNIT_RES, BC_W_S(1,LC), DIM_BC, &
618                       NEXT_RECA)
619     ! Note : previous versions did not write out BC_T_s
620                    IF (VERSION_NUMBER >= 1.15) &
621                       CALL IN_BIN_512 (UNIT_RES, BC_T_S(1,LC), DIM_BC, &
622                          NEXT_RECA)
623                    DO N = 1, NMAX(LC)
624                       CALL IN_BIN_512 (UNIT_RES, BC_X_S(1,LC,N), &
625                          DIM_BC, NEXT_RECA)
626                    ENDDO
627                 ENDIF
628                 CALL IN_BIN_512 (UNIT_RES, BC_VOLFLOW_S(1,LC), DIM_BC,&
629                    NEXT_RECA)
630                 CALL IN_BIN_512 (UNIT_RES, BC_MASSFLOW_S(1,LC), DIM_BC,&
631                     NEXT_RECA)
632              ENDDO
633           ENDIF
634           call bcast(BC_ROP_S, PE_IO) ! BCAST2d
635           call bcast(BC_U_S, PE_IO)   ! BCAST2d
636           call bcast(BC_V_S, PE_IO)   ! BCAST2d
637           if (VERSION_NUMBER >= 1.04) then
638              call bcast(BC_W_S, PE_IO)   ! BCAST2d
639              if (VERSION_NUMBER >= 1.15)  call bcast(BC_T_S, PE_IO)   ! BCAST2d
640              call bcast(BC_X_S, PE_IO)   ! BCAST2d
641           endif
642           call bcast(BC_VOLFLOW_S, PE_IO)   ! BCAST2d
643           call bcast(BC_MASSFLOW_S, PE_IO)  ! BCAST2d
644     
645           IF (myPE == PE_IO) THEN
646              IF (VERSION == 'RES = 01.00') THEN
647                 L = 10
648              ELSE
649                 L = DIM_BC
650              ENDIF
651              DO LC = 1, L
652                 READ (UNIT_RES, REC=NEXT_RECA) BC_TYPE(LC)
653                 NEXT_RECA = NEXT_RECA + 1
654              ENDDO
655           ENDIF
656           call bcast(BC_TYPE, PE_IO)   ! BCAST1c
657     
658           if (myPE == PE_IO) then
659              Allocate(IGTEMP1(IJKMAX2))   ! ALLOCate INT Global scratch
660              Allocate(iGTEMP2(IJKMAX3))   ! ALLOCate INT Global scratch
661              CALL IN_BIN_512I (UNIT_RES, IGTEMP1, IJKMAX2, NEXT_RECA)
662              call convert_from_io_i(IGTEMP1,iGTEMP2,ijkmax2)
663           else
664              Allocate(IGTEMP1(1))   ! ALLOCate INT Global scratch
665              Allocate(iGTEMP2(1))   ! ALLOCate INT Global scratch
666           endif
667           call scatter(flag,iGTEMP2,PE_IO)        !//PAR_I/O SCATTER1d
668           DeAllocate (IGTEMP1, iGTEMP2)
669     
670     
671     ! ------------------------------------------------------------------------
672           IF (VERSION_NUMBER >= 1.04) THEN
673              if (myPE == PE_IO) then
674                 CALL IN_BIN_512 (UNIT_RES, IS_X_W, DIM_IS, NEXT_RECA)
675                 CALL IN_BIN_512 (UNIT_RES, IS_X_E, DIM_IS, NEXT_RECA)
676                 CALL IN_BIN_512 (UNIT_RES, IS_Y_S, DIM_IS, NEXT_RECA)
677                 CALL IN_BIN_512 (UNIT_RES, IS_Y_N, DIM_IS, NEXT_RECA)
678                 CALL IN_BIN_512 (UNIT_RES, IS_Z_B, DIM_IS, NEXT_RECA)
679                 CALL IN_BIN_512 (UNIT_RES, IS_Z_T, DIM_IS, NEXT_RECA)
680                 CALL IN_BIN_512I (UNIT_RES, IS_I_W, DIM_IS, NEXT_RECA)
681                 CALL IN_BIN_512I (UNIT_RES, IS_I_E, DIM_IS, NEXT_RECA)
682                 CALL IN_BIN_512I (UNIT_RES, IS_J_S, DIM_IS, NEXT_RECA)
683                 CALL IN_BIN_512I (UNIT_RES, IS_J_N, DIM_IS, NEXT_RECA)
684                 CALL IN_BIN_512I (UNIT_RES, IS_K_B, DIM_IS, NEXT_RECA)
685                 CALL IN_BIN_512I (UNIT_RES, IS_K_T, DIM_IS, NEXT_RECA)
686                 CALL IN_BIN_512 (UNIT_RES, IS_PC(1,1), DIM_IS, NEXT_RECA)
687                 CALL IN_BIN_512 (UNIT_RES, IS_PC(1,2), DIM_IS, NEXT_RECA)
688              endif
689              call bcast(IS_X_W, PE_IO)   ! BCAST1d
690              call bcast(IS_X_E, PE_IO)   ! BCAST1d
691              call bcast(IS_Y_S, PE_IO)   ! BCAST1d
692              call bcast(IS_Y_N, PE_IO)   ! BCAST1d
693              call bcast(IS_Z_B, PE_IO)   ! BCAST1d
694              call bcast(IS_Z_T, PE_IO)   ! BCAST1d
695              call bcast(IS_I_W, PE_IO)   ! BCAST1i
696              call bcast(IS_I_E, PE_IO)   ! BCAST1i
697              call bcast(IS_J_S, PE_IO)   ! BCAST1i
698              call bcast(IS_J_N, PE_IO)   ! BCAST1i
699              call bcast(IS_K_B, PE_IO)   ! BCAST1i
700              call bcast(IS_K_T, PE_IO)   ! BCAST1i
701              call bcast(IS_PC, PE_IO)    ! BCAST1i
702     
703              IF (VERSION_NUMBER >= 1.07) THEN
704                 if (myPE == PE_IO) then
705                    DO LC = 1, MMAX
706                      CALL IN_BIN_512 (UNIT_RES, IS_VEL_S(1,LC), &
707                            DIM_IS, NEXT_RECA)
708                    ENDDO
709                 endif
710                 call bcast(IS_VEL_S, PE_IO)   ! BCAST2d
711              ENDIF
712     
713              if (myPE == PE_IO) then
714                 DO LC = 1, DIM_IS
715                    READ (UNIT_RES, REC=NEXT_RECA) IS_TYPE(LC)
716                    NEXT_RECA = NEXT_RECA + 1
717                  ENDDO
718              endif
719              call bcast(IS_TYPE, PE_IO)   ! BCAST1c
720           ENDIF
721     
722     ! ------------------------------------------------------------------------
723     ! Additions from new versions of .RES file
724           IF (VERSION_NUMBER >= 1.08) THEN
725              IF (myPE == PE_IO) THEN
726                 READ (UNIT_RES, REC=NEXT_RECA) CYCLIC_X, &
727                    CYCLIC_Y, CYCLIC_Z, CYCLIC_X_PD, &
728                    CYCLIC_Y_PD, CYCLIC_Z_PD, DELP_X, DELP_Y, &
729                    DELP_Z, U_G0, U_S0, V_G0, V_S0, W_G0, W_S0
730                 NEXT_RECA = NEXT_RECA + 1
731              ENDIF
732              call bcast(CYCLIC_X,PE_IO)     ! BCAST0l
733              call bcast(CYCLIC_Y,PE_IO)     ! BCAST0l
734              call bcast(CYCLIC_Z,PE_IO)     ! BCAST0l
735              call bcast(CYCLIC_X_PD,PE_IO)  ! BCAST0l
736              call bcast(CYCLIC_Y_PD,PE_IO)  ! BCAST0l
737              call bcast(CYCLIC_Z_PD,PE_IO)  ! BCAST0l
738              call bcast(DELP_X,PE_IO)       ! BCAST1d
739              call bcast(DELP_Y,PE_IO)       ! BCAST1d
740              call bcast(DELP_Z,PE_IO)       ! BCAST1d
741              call bcast(U_G0,PE_IO)         ! BCAST1d
742              call bcast(U_S0,PE_IO)         ! BCAST1d
743              call bcast(V_G0,PE_IO)         ! BCAST1d
744              call bcast(V_S0,PE_IO)         ! BCAST1d
745              call bcast(W_G0,PE_IO)         ! BCAST1d
746              call bcast(W_S0,PE_IO)         ! BCAST1d
747           ENDIF
748     
749     ! ------------------------------------------------------------------------
750           IF (VERSION_NUMBER >= 1.09) THEN
751              IF (myPE == PE_IO) THEN
752                 READ (UNIT_RES, REC=NEXT_RECA) TIME, TSTOP, &
753                     ENERGY_EQ, RES_DT, OUT_DT, NLOG, &
754                     L_SCALE0, NO_I, NO_J, NO_K, CALL_USR
755                  NEXT_RECA = NEXT_RECA + 1
756              ENDIF
757              call bcast(TIME,PE_IO)      ! BCAST0d
758              call bcast(TSTOP,PE_IO)     ! BCAST0d
759              call bcast(ENERGY_EQ,PE_IO) ! BCAST0l
760              call bcast(RES_DT,PE_IO)    ! BCAST0d
761              call bcast(OUT_DT,PE_IO)    ! BCAST0d
762              call bcast(NLOG,PE_IO)      ! BCAST0i
763              call bcast(L_SCALE0,PE_IO)  ! BCAST0d
764              call bcast(NO_I,PE_IO)      ! BCAST0l
765              call bcast(NO_J,PE_IO)      ! BCAST0l
766              call bcast(NO_K,PE_IO)      ! BCAST0l
767              call bcast(CALL_USR,PE_IO)  ! BCAST0l
768     
769              IF (myPE == PE_IO) THEN
770                 if (version_number >= 1.50) then
771                    read (unit_res,rec=next_reca) n_spx_res
772                    next_reca = next_reca + 1
773                    if (n_spx_res > n_spx) then
774                       write (*,*) ' n_spx too small '
775                       write (*,*) ' n_spx = ' , n_spx
776                       write (*,*) ' n_spx must equal ' , n_spx_res
777                       call exitMPI(myPE)  ! Abort all PEs, not only the current one
778                    endif
779                 else
780                    n_spx_res = 9
781                 endif
782     
783                 DO LC = 1, N_SPX_RES
784                    READ (UNIT_RES, REC=NEXT_RECA) SPX_DT(LC)
785                    NEXT_RECA = NEXT_RECA + 1
786                 ENDDO
787                 DO LC = 0, MMAX
788                    READ (UNIT_RES, REC=NEXT_RECA) SPECIES_EQ(LC)
789                    NEXT_RECA = NEXT_RECA + 1
790                 ENDDO
791              ENDIF
792              call bcast(SPX_DT,PE_IO)     ! BCAST1d
793              call bcast(SPECIES_EQ,PE_IO) ! BCAST1l (recv)
794     
795              IF (myPE == PE_IO) THEN
796                 CALL IN_BIN_512 (UNIT_RES, USR_DT, DIMENSION_USR, NEXT_RECA)
797                 CALL IN_BIN_512 (UNIT_RES, USR_X_W, DIMENSION_USR, NEXT_RECA)
798                 CALL IN_BIN_512 (UNIT_RES, USR_X_E, DIMENSION_USR, NEXT_RECA)
799                 CALL IN_BIN_512 (UNIT_RES, USR_Y_S, DIMENSION_USR, NEXT_RECA)
800                 CALL IN_BIN_512 (UNIT_RES, USR_Y_N, DIMENSION_USR, NEXT_RECA)
801                 CALL IN_BIN_512 (UNIT_RES, USR_Z_B, DIMENSION_USR, NEXT_RECA)
802                 CALL IN_BIN_512 (UNIT_RES, USR_Z_T, DIMENSION_USR, NEXT_RECA)
803              ENDIF
804              call bcast(USR_DT,PE_IO)  ! BCAST1d
805              call bcast(USR_X_W,PE_IO) ! BCAST1d
806              call bcast(USR_X_E,PE_IO) ! BCAST1d
807              call bcast(USR_Y_S,PE_IO) ! BCAST1d
808              call bcast(USR_Y_N,PE_IO) ! BCAST1d
809              call bcast(USR_Z_B,PE_IO) ! BCAST1d
810              call bcast(USR_Z_T,PE_IO) ! BCAST1d
811     
812              IF (myPE == PE_IO) THEN
813                 DO LC = 1, DIMENSION_USR
814                    READ (UNIT_RES, REC=NEXT_RECA) USR_FORMAT(LC), &
815                       USR_EXT(LC), USR_TYPE(LC), USR_VAR(LC)
816                    NEXT_RECA = NEXT_RECA + 1
817                 ENDDO
818              ENDIF
819              call bcast(USR_FORMAT,PE_IO) ! BCAST1c
820              call bcast(USR_EXT,PE_IO)    ! BCAST1c
821              call bcast(USR_TYPE,PE_IO)   ! BCAST1c
822              call bcast(USR_VAR,PE_IO)    ! BCAST1c
823     
824              if (myPE == PE_IO) then
825                 CALL IN_BIN_512 (UNIT_RES, IC_P_STAR, DIM_IC, NEXT_RECA)
826                 CALL IN_BIN_512 (UNIT_RES, IC_L_SCALE, DIM_IC, NEXT_RECA)
827              endif
828              call bcast(IC_P_STAR,PE_IO)  ! BCAST1d
829              call bcast(IC_L_SCALE,PE_IO) ! BCAST1d
830     
831              IF (myPE == PE_IO) THEN
832                 DO LC = 1, DIM_IC
833                    READ (UNIT_RES, REC=NEXT_RECA) IC_TYPE(LC)
834                    NEXT_RECA = NEXT_RECA + 1
835                 ENDDO
836              ENDIF
837              call bcast(IC_TYPE,PE_IO) ! BCAST1c
838     
839              IF (myPE == PE_IO) THEN
840                 CALL IN_BIN_512 (UNIT_RES, BC_DT_0, DIM_BC, NEXT_RECA)
841                 CALL IN_BIN_512 (UNIT_RES, BC_JET_G0, DIM_BC, NEXT_RECA)
842                 CALL IN_BIN_512 (UNIT_RES, BC_DT_H, DIM_BC , NEXT_RECA)
843                 CALL IN_BIN_512 (UNIT_RES, BC_JET_GH, DIM_BC, NEXT_RECA)
844                 CALL IN_BIN_512 (UNIT_RES, BC_DT_L, DIM_BC , NEXT_RECA)
845                 CALL IN_BIN_512 (UNIT_RES, BC_JET_GL, DIM_BC, NEXT_RECA)
846              ENDIF
847              call bcast(BC_DT_0,PE_IO) !//PAR_I/O BCAST1d
848              call bcast(BC_JET_G0,PE_IO) !//PAR_I/O BCAST1d
849              call bcast(BC_DT_H,PE_IO) !//PAR_I/O BCAST1d
850              call bcast(BC_JET_GH,PE_IO) !//PAR_I/O BCAST1d
851              call bcast(BC_DT_L,PE_IO) !//PAR_I/O BCAST1d
852              call bcast(BC_JET_GL,PE_IO) !//PAR_I/O BCAST1d
853           ENDIF   ! endif (version_number >=1.09)
854     
855     ! ------------------------------------------------------------------------
856           IF (VERSION_NUMBER >= 1.10) THEN
857              IF (myPE == PE_IO) THEN
858                 READ (UNIT_RES, REC=NEXT_RECA) MU_GMAX
859                 NEXT_RECA = NEXT_RECA + 1
860              ENDIF
861              call bcast(MU_GMAX,PE_IO) ! BCAST0d
862           ENDIF
863     
864     ! ------------------------------------------------------------------------
865          IF (VERSION_NUMBER >= 1.11) THEN
866             IF (myPE == PE_IO) THEN
867                READ (UNIT_RES, REC=NEXT_RECA) V_EX, MODEL_B
868                NEXT_RECA = NEXT_RECA + 1
869             ENDIF
870             call bcast(V_EX,PE_IO)    ! BCAST0d
871             call bcast(MODEL_B,PE_IO) ! BCAST0l
872          ENDIF
873     
874     ! ------------------------------------------------------------------------
875          IF (VERSION_NUMBER >= 1.12) THEN
876             IF (myPE == PE_IO) THEN
877                READ (UNIT_RES, REC=NEXT_RECA) P_REF, &
878                   P_SCALE, UR_FAC, TOL_RESID, DT_MAX, &
879                   DT_MIN, DT_FAC, CLOSE_PACKED, GRAVITY, &
880                   MU_S0
881                 NEXT_RECA = NEXT_RECA + 1
882             ENDIF
883             call bcast(P_REF,PE_IO)        ! BCAST0d
884             call bcast(P_SCALE,PE_IO)      ! BCAST0d
885             call bcast(UR_FAC,PE_IO)       ! BCAST0d
886             call bcast(TOL_RESID,PE_IO)    ! BCAST0d
887             call bcast(DT_MAX,PE_IO)       ! BCAST0d
888             call bcast(DT_MIN,PE_IO)       ! BCAST0d
889             call bcast(DT_FAC,PE_IO)       ! BCAST0d
890             call bcast(CLOSE_PACKED,PE_IO) ! BCAST0l
891             call bcast(GRAVITY,PE_IO)      ! BCAST0d
892             call bcast(MU_S0,PE_IO)        ! BCAST0d
893     
894             IF (myPE == PE_IO) THEN
895                READ (UNIT_RES, REC=NEXT_RECA) LEQ_IT, LEQ_METHOD
896                NEXT_RECA = NEXT_RECA + 1
897             ENDIF
898             call bcast(LEQ_IT,PE_IO)     ! BCAST1i
899             call bcast(LEQ_METHOD,PE_IO) ! BCAST1i
900     
901             IF (myPE == PE_IO) THEN
902                CALL IN_BIN_512 (UNIT_RES, BC_HW_G, DIM_BC, NEXT_RECA)
903                CALL IN_BIN_512 (UNIT_RES, BC_UW_G, DIM_BC, NEXT_RECA)
904                CALL IN_BIN_512 (UNIT_RES, BC_VW_G, DIM_BC, NEXT_RECA)
905                 CALL IN_BIN_512 (UNIT_RES, BC_WW_G, DIM_BC, NEXT_RECA)
906             ENDIF
907             call bcast(BC_HW_G,PE_IO) ! BCAST1d
908             call bcast(BC_UW_G,PE_IO) ! BCAST1d
909             call bcast(BC_VW_G,PE_IO) ! BCAST1d
910             call bcast(BC_WW_G,PE_IO) ! BCAST1d
911     
912             IF (myPE == PE_IO) THEN
913                DO LC = 1, MMAX
914                   CALL IN_BIN_512 (UNIT_RES, BC_HW_S(1,LC), DIM_BC, NEXT_RECA)
915                   CALL IN_BIN_512 (UNIT_RES, BC_UW_S(1,LC), DIM_BC, NEXT_RECA)
916                   CALL IN_BIN_512 (UNIT_RES, BC_VW_S(1,LC), DIM_BC, NEXT_RECA)
917                   CALL IN_BIN_512 (UNIT_RES, BC_WW_S(1,LC), DIM_BC, NEXT_RECA)
918                 ENDDO
919             ENDIF
920             call bcast(BC_HW_S,PE_IO) ! BCAST2d
921             call bcast(BC_UW_S,PE_IO) ! BCAST2d
922             call bcast(BC_VW_S,PE_IO) ! BCAST2d
923             call bcast(BC_WW_S,PE_IO) ! BCAST2d
924          ENDIF   ! endif (version_number >=1.12)
925     
926          LC = 0
927          IF (MMAX + 1 > 0) THEN
928             MOMENTUM_X_EQ(:MMAX) = .TRUE.
929             MOMENTUM_Y_EQ(:MMAX) = .TRUE.
930             MOMENTUM_Z_EQ(:MMAX) = .TRUE.
931             LC = MMAX + 1
932          ENDIF
933          TOL_DIVERGE = 1.E+4
934     
935     ! ------------------------------------------------------------------------
936          IF (VERSION_NUMBER >= 1.13) THEN
937             IF (myPE == PE_IO) THEN
938                READ (UNIT_RES, REC=NEXT_RECA) MOMENTUM_X_EQ, &
939                   MOMENTUM_Y_EQ, MOMENTUM_Z_EQ, TOL_DIVERGE, &
940                   DISCRETIZE, FULL_LOG
941                NEXT_RECA = NEXT_RECA + 1
942             ENDIF
943             call bcast(MOMENTUM_X_EQ,PE_IO) ! BCAST1l
944             call bcast(MOMENTUM_Y_EQ,PE_IO) ! BCAST1l
945             call bcast(MOMENTUM_Z_EQ,PE_IO) ! BCAST1l
946             call bcast(TOL_DIVERGE,PE_IO)   ! BCAST0d
947             call bcast(DISCRETIZE,PE_IO)    ! BCAST1i
948             call bcast(FULL_LOG,PE_IO)      ! BCAST0l
949          ENDIF
950     
951     ! ------------------------------------------------------------------------
952          IF (VERSION_NUMBER >= 1.14) THEN
953             IF (myPE == PE_IO) THEN
954                READ (UNIT_RES, REC=NEXT_RECA) DETECT_STALL
955                NEXT_RECA = NEXT_RECA + 1
956             ENDIF
957             call bcast(DETECT_STALL,PE_IO) ! BCAST0l
958          ENDIF
959     
960     ! ------------------------------------------------------------------------
961          IF (VERSION_NUMBER >= 1.15) THEN
962             IF (myPE == PE_IO) THEN
963                READ (UNIT_RES, REC=NEXT_RECA) K_G0, K_S0(1), &
964                   C_PG0, C_PS0(1), TOL_RESID_T, TOL_RESID_X
965                NEXT_RECA = NEXT_RECA + 1
966             CALL IN_BIN_512 (UNIT_RES, IC_GAMA_RG, DIM_IC, NEXT_RECA)
967             CALL IN_BIN_512 (UNIT_RES, IC_T_RG, DIM_IC, NEXT_RECA)
968             ENDIF
969             call bcast(K_G0,PE_IO)        ! BCAST0d
970             call bcast(K_S0(1),PE_IO)        ! BCAST0d
971             call bcast(C_PG0,PE_IO)       ! BCAST0d
972             call bcast(C_PS0(1),PE_IO)    ! BCAST0d
973             call bcast(TOL_RESID_T,PE_IO) ! BCAST0d
974             call bcast(TOL_RESID_X,PE_IO) ! BCAST0d
975             call bcast(IC_GAMA_RG,PE_IO)  ! BCAST1d
976             call bcast(IC_T_RG,PE_IO)     ! BCAST1d
977     
978             IF (myPE == PE_IO) THEN
979                DO LC = 1, MMAX
980                   CALL IN_BIN_512 (UNIT_RES, IC_GAMA_RS(1,LC), DIM_IC, &
981                      NEXT_RECA)
982                   CALL IN_BIN_512 (UNIT_RES, IC_T_RS(1,LC), DIM_IC, &
983                      NEXT_RECA)
984                ENDDO
985             ENDIF
986             call bcast(IC_GAMA_RS,PE_IO) ! BCAST2d
987             call bcast(IC_T_RS,PE_IO)    ! BCAST2d
988          ENDIF
989     
990     ! ------------------------------------------------------------------------
991          IF (VERSION_NUMBER >= 1.2) THEN
992             IF (myPE == PE_IO) THEN
993                READ (UNIT_RES, REC=NEXT_RECA) NORM_G, NORM_S
994                NEXT_RECA = NEXT_RECA + 1
995             ENDIF
996             call bcast(NORM_G,PE_IO) ! BCAST0d
997             call bcast(NORM_S,PE_IO) ! BCAST0d
998          ENDIF
999     
1000     ! ------------------------------------------------------------------------
1001          IF (VERSION_NUMBER >= 1.3) THEN
1002             IF (myPE == PE_IO) THEN
1003                READ (UNIT_RES, REC=NEXT_RECA) NScalar, TOL_RESID_Scalar, DIM_tmp
1004                NEXT_RECA = NEXT_RECA + 1
1005                CALL IN_BIN_512I (UNIT_RES, Phase4Scalar, DIM_tmp, NEXT_RECA)
1006     
1007     ! post mfix fix ...
1008                if (doingPost .and. nscalar.gt.0) then
1009                   DIMENSION_Scalar = NScalar
1010                   Allocate(  Scalar (DIMENSION_3,  DIMENSION_Scalar) )
1011                   Allocate(  Scalaro (DIMENSION_3, DIMENSION_Scalar) )
1012                endif
1013     
1014             ENDIF
1015             call bcast(NScalar,PE_IO)          ! BCAST0d
1016             call bcast(TOL_RESID_Scalar,PE_IO) ! BCAST0d
1017             call bcast(Phase4Scalar,PE_IO)     ! BCAST0d
1018          ELSE
1019             NScalar = 0
1020          ENDIF
1021     
1022     ! ------------------------------------------------------------------------
1023     ! Version 1.4 -- read radiation variables in read_res1
1024           IF (VERSION_NUMBER >= 1.499) THEN
1025              IF (myPE == PE_IO) THEN
1026                 READ (UNIT_RES, REC=NEXT_RECA) nRR
1027                 NEXT_RECA = NEXT_RECA + 1
1028                 if (doingPost .and. nRR.gt.0) then
1029                    Allocate( ReactionRates(DIMENSION_3, nRR) )
1030                 endif
1031              ENDIF
1032              call bcast(nRR,PE_IO) ! BCAST0d
1033           ELSE
1034              nRR = 0
1035           ENDIF
1036     
1037     ! ------------------------------------------------------------------------
1038     ! Version 1.6 -- read K_Epsilon and dqmom
1039           IF (VERSION_NUMBER >= 1.599) THEN
1040              IF (myPE == PE_IO) THEN
1041                 READ (UNIT_RES, REC=NEXT_RECA) K_Epsilon, Call_DQMOM
1042                 NEXT_RECA = NEXT_RECA + 1
1043                 if (doingPost .and. K_epsilon) then
1044                    Allocate( K_Turb_G(DIMENSION_3) )
1045                    Allocate( E_Turb_G(DIMENSION_3) )
1046                 end if
1047              ENDIF
1048              call bcast(K_Epsilon,PE_IO) !//PAR_I/O BCAST0d
1049              call bcast(Call_DQMOM,PE_IO) !//PAR_I/O BCAST0d
1050           ELSE
1051              K_Epsilon = .FALSE.
1052              Call_DQMOM =.FALSE.
1053           ENDIF
1054     
1055     ! ------------------------------------------------------------------------
1056     ! Version 1.7 -- Stiff Chemistry
1057           IF (VERSION_NUMBER >= 1.699) THEN
1058              IF (myPE == PE_IO) THEN
1059                 READ (UNIT_RES, REC=NEXT_RECA) STIFF_CHEMISTRY, CALL_ISAT
1060                 NEXT_RECA = NEXT_RECA + 1
1061              ENDIF
1062              call bcast(STIFF_CHEMISTRY,PE_IO)
1063              !call bcast(CALL_ISAT,PE_IO)
1064           ELSE
1065              STIFF_CHEMISTRY = .FALSE.
1066           ENDIF
1067     
1068     ! ------------------------------------------------------------------------
1069     ! Version 1.8 -- Variable solid density and  each solids species
1070           IF (VERSION_NUMBER >= 1.799) THEN
1071              IF (myPE == PE_IO) THEN
1072                 READ (UNIT_RES, REC=NEXT_RECA) (SOLVE_ROs(LC),LC=1,MMAX)
1073                 NEXT_RECA = NEXT_RECA + 1
1074                 DO LC = 1, MMAX
1075                    READ (UNIT_RES, REC=NEXT_RECA) INERT_SPECIES(LC), &
1076                       BASE_ROs(LC), (X_s0(LC,N),N=1,NMAX(LC))
1077                    NEXT_RECA = NEXT_RECA + 1
1078                 ENDDO
1079              ENDIF
1080              call bcast(SOLVE_ROs,PE_IO)
1081              call bcast(INERT_SPECIES,PE_IO)
1082              call bcast(BASE_ROs,PE_IO)
1083              call bcast(X_s0,PE_IO)
1084     
1085              ANY_SOLVE_ROs = ANY(SOLVE_ROs)
1086              IF (doingPost .AND. (.NOT.ANY_SOLVE_ROs)) THEN
1087                 DO LC = 1, MMAX
1088                    RO_S(:,LC) = RO_S0(LC)
1089                 END DO
1090              ENDIF
1091           ELSE
1092              SOLVE_ROs = .FALSE.
1093              ANY_SOLVE_ROs = .FALSE.
1094              DO LC = 1, MMAX
1095                 RO_S(:,LC) = RO_S0(LC)
1096              END DO
1097           ENDIF
1098     
1099     ! Add new read statements above this line.  Remember to update NEXT_RECA.
1100     ! Remember to update the version number check near beginning of this subroutine.
1101     !------------------------------------------------------------------------------
1102     
1103           READ (UNIT_RES, REC=3) NEXT_RECA
1104     
1105     ! Since the value of UNDEFINED was changed ...
1106           IF (RO_G0 >= 1E30) RO_G0 = UNDEFINED
1107           IF (MU_G0 >= 1E30) MU_G0 = UNDEFINED
1108           IF (MW_AVG >= 1E30) MW_AVG = UNDEFINED
1109           IF (C_E >= 1E30) C_E = UNDEFINED
1110     
1111           RETURN
1112     
1113     ! HERE IF DIMENSION ERROR
1114     
1115       900 CONTINUE
1116           WRITE (*, *) ' '
1117           WRITE (*, *) ' **************************************'
1118           WRITE (*, "('(PE ',I6,'): From: READ_RES0')") myPE
1119           WRITE (*, *) ' DIMENSION ERROR ---'
1120           WRITE (*, *) ' '
1121           WRITE (*, *) ' DIMENSION_IC = ', DIMENSION_IC, ' DIM_IC       = ', DIM_IC
1122           WRITE (*, *) ' DIMENSION_BC = ', DIMENSION_BC, ' DIM_BC       = ', DIM_BC
1123           WRITE (*, *) ' DIMENSION_IS = ', DIMENSION_IS, ' DIM_IS       = ', DIM_IS
1124           WRITE (*, *) ' DIMENSION_C  = ', DIMENSION_C, ' DIM_C        = ', DIM_C
1125           WRITE (*, *) ' '
1126     
1127           call exitMPI(myPE)
1128     
1129           END SUBROUTINE READ_RES0
1130     
1131