File: N:\mfix\model\check_data_30.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !  Module name: CHECK_DATA_30                                          !
3     !  Author: M. Syamlal                                 Date: 27-OCT-92  !
4     !                                                                      !
5     !  Purpose: Check whether the sum of reaction rates is zero and the    !
6     !  sum of mass fractions is 1.0 and EP_g >= EP_Star.                   !
7     !           and EP_g >= EP_star. Set miscellaneous constants           !
8     !                                                                      !
9     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
10           SUBROUTINE CHECK_DATA_30
11     
12     ! Global variables: (common to sub-functions)
13     !---------------------------------------------------------------------//
14           use run, only: ENERGY_EQ, SPECIES_EQ
15           use compar, only: ISTART2, IEND2
16           use compar, only: JSTART2, JEND2
17           use compar, only: KSTART2, KEND2
18           use physprop, only: SMAX, NMAX
19           use param1, only: ZERO, ONE
20           use rxns, only: USE_RRATES
21           use mms, only: USE_MMS
22           use param, only: DIMENSION_M
23     
24           use error_manager
25     
26           IMPLICIT NONE
27     
28     ! Local variables:
29     !---------------------------------------------------------------------//
30     ! Flag for error message headers.
31           INTEGER :: ERR_TAG
32     ! Flags to report stats on species equations.
33           LOGICAL :: REPORT_SPECIES(0:DIMENSION_M)
34     !......................................................................!
35     
36     
37     ! Check physical properties in inflow/outflow cells.
38           IF(.NOT.USE_MMS) CALL CHECK_FLOW_CELL_PROPS
39     ! Verify physical values for field variables.
40           CALL CHECK_PHYSICAL_BOUNDS
41     ! Generate species report if needed.
42           IF(ANY(REPORT_SPECIES)) CALL REPORT_SPECIES_STATS
43     ! Check interphase mass transfer.
44           IF(USE_RRATES) CALL CHECK_RXN_MASS_BALANCE
45     
46           RETURN
47     
48           CONTAINS
49     
50     !----------------------------------------------------------------------!
51     !  Module name: CHECK_FLOW_PROPS                                       !
52     !  Author: M. Syamlal                                 Date: 27-OCT-92  !
53     !                                                                      !
54     !  Purpose: Verify that inflow/outflow cells do not contain physical   !
55     !  properties for specified variables.                                 !
56     !                                                                      !
57     !----------------------------------------------------------------------!
58           SUBROUTINE CHECK_FLOW_CELL_PROPS
59     
60     ! Global variables:
61     !---------------------------------------------------------------------//
62           use visc_s, only: MU_S, LAMBDA_S
63           use visc_g, only: MU_GT, LAMBDA_GT
64           use physprop, only: K_G, K_S
65           use physprop, only: DIF_s, DIF_g
66     
67           use functions, only: FUNIJK
68           use functions, only: FLOW_AT
69           use compar, only: DEAD_CELL_AT
70     
71           IMPLICIT NONE
72     
73     ! Local variables:
74     !---------------------------------------------------------------------//
75     ! Loop counters
76           INTEGER :: I, J, K, IJK
77           INTEGER :: M, N
78     ! Integer error flag.
79           INTEGER :: IER
80     !......................................................................!
81     
82           CALL INIT_ERR_MSG("CHECK_FLOW_CELL_PROPS")
83     
84           IER = 0
85           ERR_TAG=2000
86     
87           DO K = KSTART2, KEND2
88           DO J = JSTART2, JEND2
89           DO I = ISTART2, IEND2
90     
91              IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
92     
93              IJK = FUNIJK(I,J,K)
94     
95              IF(FLOW_AT(IJK)) THEN
96     
97     ! Turbulent viscosity of fluid phase.
98                 IF(MU_gt(IJK) /= ZERO) CALL REPORT_ERROR                   &
99                    (IER, I, J, K, MU_GT(IJK), '/=', ZERO, 'MU_GT')
100     
101     ! Granular second coefficient of viscosity.
102                 IF(LAMBDA_gt(IJK) /= ZERO) CALL REPORT_ERROR               &
103                    (IER, I, J, K, LAMBDA_gt(IJK), '/=', ZERO, 'LAMBDA_GT')
104     ! Gas conductivity.
105                 IF(K_g(IJK) /= ZERO) CALL REPORT_ERROR                     &
106                    (IER, I, J, K, K_g(IJK),'/=',ZERO,'K_G')
107     ! Diffusivity of gas species N.
108                 DO N = 1, NMAX(0)
109                    IF(DIF_g(IJK, N) /= ZERO) CALL REPORT_ERROR             &
110                       (IER, I, J, K, DIF_g(IJK,N),'/=',ZERO,'DIF_G',N)
111                 ENDDO
112     
113                 DO M = 1, SMAX
114     ! Granular first coefficient of (shear) viscosity.
115                    IF(MU_s(IJK, M) /= ZERO) CALL REPORT_ERROR              &
116                       (IER, I, J, K, MU_s(IJK,M),'/=',ZERO,'MU_S',M)
117     ! Granular second coefficient of viscosity.
118                    IF(LAMBDA_s(IJK, M) /= ZERO) CALL REPORT_ERROR          &
119                       (IER, I, J, K, LAMBDA_s(IJK,M),'/=',ZERO,'LAMBDA_S',M)
120     ! Solids thermal conductivity.
121                    IF(K_s(IJK, M) /= ZERO) CALL REPORT_ERROR               &
122                       (IER, I, J, K, K_s(IJK,M),'/=',ZERO,'K_S', M)
123     ! Diffusivity of solids phase M, species N.
124                    DO N = 1, NMAX(M)
125                       IF(DIF_s(IJK,M,N) /= ZERO) CALL REPORT_ERROR         &
126                          (IER, I,J,K, DIF_s(IJK,M,N),'/=',ZERO,'DIF_S',M,N)
127                    ENDDO
128                 ENDDO
129              ENDIF
130     
131           ENDDO
132           ENDDO
133           ENDDO
134     
135           IF(IER /= 0) THEN
136              WRITE(ERR_MSG,"('End of Report.')")
137              CALL FLUSH_ERR_MSG(HEADER=.FALSE., ABORT=.TRUE.)
138     ! Close DMP logs when running in interactive mode.
139              CALL CLOSE_PE_LOG
140           ENDIF
141     
142           CALL FINL_ERR_MSG
143     
144           RETURN
145           END SUBROUTINE CHECK_FLOW_CELL_PROPS
146     
147     
148     
149     !----------------------------------------------------------------------!
150     !                                                                      !
151     !  Module name: CHECK_PHYSICAL_BOUNDS                                  !
152     !  Author: M. Syamlal                                 Date: 27-OCT-92  !
153     !                                                                      !
154     !  Purpose: Verify that fluid cells have physical values for the       !
155     !  specified variables.                                                !
156     !                                                                      !
157     !----------------------------------------------------------------------!
158           SUBROUTINE CHECK_PHYSICAL_BOUNDS
159     
160     ! Global variables:
161     !---------------------------------------------------------------------//
162           use toleranc, only: TMIN, TMAX, TOL_COM
163           use fldvar, only: X_G, X_S
164           use fldvar, only: T_G, T_S
165           use fldvar, only: ROP_s
166           use rxns, only: RoX_GC, RoX_SC
167           use rxns, only: R_GP, R_SP
168           use visc_s, only: MU_S
169           use physprop, only: K_G, K_S
170           use physprop, only: C_PG, C_PS
171           use physprop, only: DIF_s, DIF_g
172           use physprop, only: MU_G
173           use physprop, only: MW_MIX_G
174     
175           use mpi_utility, only: GLOBAL_ALL_SUM
176           use mpi_utility, only: GLOBAL_ALL_OR
177     
178           use functions, only: FUNIJK
179           use functions, only: WALL_AT
180           use compar, only: DEAD_CELL_AT
181     
182           IMPLICIT NONE
183     
184     ! Local variables:
185     !---------------------------------------------------------------------//
186           INTEGER :: I, J, K, IJK
187           INTEGER :: M, N
188     ! Integer error flag.
189           INTEGER :: IER
190     
191           CALL INIT_ERR_MSG("CHECK_PHYSICAL_BOUNDS")
192     
193           IER = 0
194           ERR_TAG = 3000
195           REPORT_SPECIES = .FALSE.
196     
197           DO K = KSTART2, KEND2
198           DO J = JSTART2, JEND2
199           DO I = ISTART2, IEND2
200              IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
201              IJK = FUNIJK(I,J,K)
202              IF (.NOT.WALL_AT(IJK)) THEN
203     
204     ! Gas viscosity must be non-negative.
205                 IF(MU_G(IJK) < ZERO) CALL REPORT_ERROR                     &
206                    (IER, I, J, K, MU_G(IJK), '<', ZERO, 'MU_G')
207     
208     ! Mixture molecular weight must be positive.
209                 IF(MW_MIX_G(IJK) <= ZERO) CALL REPORT_ERROR                &
210                    (IER, I, J, K,MW_MIX_G(IJK), '<=', ZERO, 'MW_MIX_G')
211     
212     ! Verify thermodynamic properties when solving energy equations:
213                 IF(ENERGY_EQ) THEN
214     ! Gas conductivity must be non-negative.
215                    IF(K_G(IJK) < ZERO) CALL REPORT_ERROR                   &
216                       (IER, I, J, K, K_G(IJK), '<', ZERO, 'K_G')
217     ! Gas phase specific heat must be positive.
218                    IF(C_PG(IJK) <= ZERO) CALL REPORT_ERROR                 &
219                       (IER, I, J, K, C_PG(IJK), '<=', ZERO, 'C_PG')
220     ! Verify that the gas phase temperature is within the bounds.
221                    IF(T_G(IJK) <= TMIN) CALL REPORT_ERROR                  &
222                       (IER, I, J, K, T_G(IJK), '<=', TMIN, 'T_G')
223                    IF(T_G(IJK) >= TMAX) CALL REPORT_ERROR                  &
224                       (IER, I, J, K, T_G(IJK), '>=', TMAX, 'T_G')
225     
226     ! Diffusivity of gas species N must be non-negative.
227                    DO N = 1, NMAX(0)
228                       IF(DIF_g(IJK, N) < ZERO) CALL REPORT_ERROR           &
229                          (IER, I,J,K, DIF_g(IJK,N),'<',ZERO, 'DIF_G',N)
230                       IF(X_g(IJK,N) < ZERO) CALL REPORT_ERROR              &
231                          (IER, I, J, K, X_g(IJK,N), '<', ZERO, 'X_G',N)
232                       IF(X_g(IJK,N) > ONE) CALL REPORT_ERROR               &
233                          (IER, I, J, K, X_g(IJK,N), '>', ONE, 'X_G',N)
234                    ENDDO
235                 ENDIF
236     
237     ! Verify that the gas phase mass fractons sum to one.
238                 IF(SPECIES_EQ(0)) THEN
239                    IF(ABS(ONE - sum(X_G(IJK,1:NMAX(0)))) > TOL_COM) &
240                       REPORT_SPECIES(0) = .TRUE.
241     
242     ! Verify that the rates of formation and consumption adhear to the
243     ! expected coding restraints. (non-negative)
244                    IF(USE_RRATES) THEN
245                       DO N = 1, NMAX(0)
246                          IF(R_GP(IJK,N) < ZERO) CALL REPORT_ERROR          &
247                            (IER,I,J,K, R_GP(IJK,N),'<',ZERO,'R_GP',N)
248                          IF(ROX_GC(IJK,N) < ZERO) CALL REPORT_ERROR        &
249                            (IER,I,J,K,RoX_GC(IJK,N),'<',ZERO,'RoX_GC',N)
250                       ENDDO
251                    ENDIF
252     
253                 ENDIF
254     
255                 DO M = 1, SMAX
256     
257     ! Solids viscosity should be non-negativel.
258                    IF(MU_S(IJK,M) < ZERO) CALL REPORT_ERROR                &
259                       (IER, I, J, K, MU_S(IJK,M), '<', ZERO, 'MU_S',M)
260     
261                    IF(ENERGY_EQ) THEN
262     
263     ! Thermal conductivity must be non-negative.
264                       IF(K_S(IJK,M) < ZERO) CALL REPORT_ERROR              &
265                          (IER, I, J, K, K_S(IJK,M), '<', ZERO, 'K_S',M)
266     
267     ! Solids specific heat must be positive.
268                       IF(C_PS(IJK,M) <= ZERO) CALL REPORT_ERROR            &
269                          (IER, I, J, K, C_PS(IJK,M), '<=', ZERO, 'C_PS',M)
270     
271     ! Verify that the solids temperature is within the required bounds.
272                       IF(T_S(IJK,M) <= TMIN) CALL REPORT_ERROR             &
273                          (IER, I, J, K, T_S(IJK,M), '<=', TMIN, 'T_S',M)
274                       IF(T_S(IJK,M) >= TMAX) CALL REPORT_ERROR             &
275                          (IER, I, J, K, T_S(IJK,M), '>=', TMAX, 'T_S',M)
276     
277                    ENDIF
278     
279                    DO N = 1, NMAX(M)
280                       IF(DIF_s(IJK, M, N) < ZERO) CALL REPORT_ERROR        &
281                          (IER, I, J, K, DIF_s(IJK,M,N),'<',ZERO,'DIF_S',M,N)
282                    ENDDO
283     
284     ! Sum of solids mass fractions should be one
285                    IF(SPECIES_EQ(M)) THEN
286                       IF(ROP_S(IJK,M) /= ZERO) THEN
287                          IF(ABS(ONE-sum(X_S(IJK,M,1:NMAX(M)))) > TOL_COM) &
288                             REPORT_SPECIES(M) = .TRUE.
289                       ENDIF
290     
291                       DO N = 1, NMAX(M)
292                          IF(X_s(IJK, M, N) < ZERO) CALL REPORT_ERROR       &
293                             (IER, I,J,K, X_s(IJK,M,N),'<',ZERO, 'X_S',M,N)
294                          IF(X_s(IJK, M, N) > ONE) CALL REPORT_ERROR        &
295                             (IER, I,J,K, X_s(IJK,M,N),'>',ONE, 'X_S',M,N)
296                       ENDDO
297                    ENDIF
298     
299     ! Verify that the rates of formation and consumption adhear to the
300     ! expected coding restraints. (non-negative)
301                    IF(SPECIES_EQ(M) .AND. USE_RRATES) THEN
302                       DO N = 1, NMAX(0)
303                          IF(R_SP(IJK,M,N) < ZERO) CALL REPORT_ERROR        &
304                            (IER,I,J,K,R_SP(IJK,M,N),'<',ZERO,'R_SP',M,N)
305                          IF(ROX_GC(IJK,N) < ZERO) CALL REPORT_ERROR        &
306                            (IER,I,J,K,RoX_SC(IJK,M,N),'<',ZERO,'RoX_SC',M,N)
307                       ENDDO
308                    ENDIF
309                 ENDDO
310     
311              ENDIF ! IF(.NOT.WALL_AT(IJK))
312           ENDDO ! DO I = ISTART2, IEND2
313           ENDDO ! DO J = JSTART2, JEND2
314           ENDDO ! DO K = KSTART2, KEND2
315     
316           CALL GLOBAL_ALL_SUM(IER)
317     
318           IF(IER /= 0) THEN
319              WRITE(ERR_MSG,"('End of Report.')")
320              CALL FLUSH_ERR_MSG(HEADER=.FALSE., ABORT=.TRUE.)
321              CALL CLOSE_PE_LOG
322           ENDIF
323     
324           CALL GLOBAL_ALL_OR(REPORT_SPECIES)
325     
326           CALL FINL_ERR_MSG
327     
328           RETURN
329           END SUBROUTINE CHECK_PHYSICAL_BOUNDS
330     
331     
332     
333     !----------------------------------------------------------------------!
334     !                                                                      !
335     !  Module name: REPORT_SPECIES_STATS                                   !
336     !  Author: M. Syamlal                                 Date: 27-OCT-92  !
337     !                                                                      !
338     !  Purpose: Collect information on the sums of species mass fractions. !
339     !  This routine is call as-needed.                                     !
340     !                                                                      !
341     !----------------------------------------------------------------------!
342           SUBROUTINE REPORT_SPECIES_STATS
343     
344     ! Global variables:
345     !---------------------------------------------------------------------//
346           use fldvar, only: X_G, X_S
347           use fldvar, only: ROP_S
348           use run, only: TIME
349           use mpi_utility, only: GLOBAL_ALL_SUM
350           use mpi_utility, only: GLOBAL_ALL_MAX
351           use mpi_utility, only: GLOBAL_ALL_MIN
352     
353           use functions, only: FUNIJK
354           use functions, only: WALL_AT
355     
356           use param1, only: UNDEFINED
357           use compar, only: DEAD_CELL_AT
358     
359           IMPLICIT NONE
360     
361     ! Local variables:
362     !---------------------------------------------------------------------//
363           INTEGER :: I, J, K, IJK, L
364           INTEGER :: M
365     
366           INTEGER :: COUNT(0:DIMENSION_M, 9)
367           INTEGER :: MIN_LOC(0:DIMENSION_M)
368           INTEGER :: MAX_LOC(0:DIMENSION_M)
369     
370           DOUBLE PRECISION :: MIN_VAL(0:DIMENSION_M)
371           DOUBLE PRECISION :: MAX_VAL(0:DIMENSION_M)
372     
373           DOUBLE PRECISION :: lSUM
374     !......................................................................!
375     
376           CALL INIT_ERR_MSG('REPORT_SPECIES_STATS')
377     
378           WRITE(ERR_MSG, 4000) TIME
379           CALL FLUSH_ERR_MSG(FOOTER=.FALSE.)
380     
381      4000 FORMAT('Time = ',G12.5,/'Warning: The sum of mass fractions ',   &
382           'is not equal to one.')
383     
384     ! Initialize the counters
385           COUNT = 0
386           MAX_VAL = -UNDEFINED
387           MIN_VAL =  UNDEFINED
388     
389     ! Collect stats on species across the domain.
390           DO K = KSTART2, KEND2
391           DO J = JSTART2, JEND2
392           DO I = ISTART2, IEND2
393              IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
394              IJK = FUNIJK(I,J,K)
395              IF (.NOT.WALL_AT(IJK)) THEN
396     
397                 IF(SPECIES_EQ(0)) THEN
398                    lSUM = sum(X_G(IJK,1:NMAX(0)))
399                    IF (lSUM < MIN_VAL(0)) THEN
400                       MIN_VAL(0) = lSUM
401                       MIN_LOC(0) = IJK
402                    ENDIF
403                    IF (lSUM > MAX_VAL(0)) THEN
404                       MAX_VAL(0) = lSUM
405                       MAX_LOC(0) = IJK
406                    ENDIF
407                    IF (lSUM < 0.9) THEN
408                       COUNT(0,1) = COUNT(0,1) + 1 ! < 0.9pP
409                    ELSE IF (lSUM < 0.99) THEN
410                       COUNT(0,2) = COUNT(0,2) + 1 ! 0.9    - 0.99
411                    ELSE IF (lSUM < 0.999) THEN
412                       COUNT(0,3) = COUNT(0,3) + 1 ! 0.99   - 0.999
413                    ELSE IF (lSUM < 0.9999) THEN
414                       COUNT(0,4) = COUNT(0,4) + 1 ! 0.999  - 0.9999
415                    ELSE IF (lSUM < 1.0001) THEN
416                       COUNT(0,5) = COUNT(0,5) + 1 ! 0.9999 - 1.0001
417                    ELSE IF (lSUM < 1.001) THEN
418                       COUNT(0,6) = COUNT(0,6) + 1 ! 1.0001 - 1.001
419                    ELSE IF (lSUM < 1.01) THEN
420                       COUNT(0,7) = COUNT(0,7) + 1 ! 1.001  - 1.01
421                    ELSE IF (lSUM < 1.1) THEN
422                       COUNT(0,8) = COUNT(0,8) + 1 ! 1.01   - 1.1
423                    ELSE
424                       COUNT(0,9) = COUNT(0,9) + 1 ! > 1.1
425                    ENDIF
426                 ENDIF
427     
428     
429                 DO M = 1, SMAX
430     !  Sum of solids mass fractions should be one
431                    IF(SPECIES_EQ(M)) THEN
432                       lSUM = sum(X_S(IJK,M,1:NMAX(M)) )
433                       IF(ROP_S(IJK,M) /= ZERO) THEN
434                          IF(lSUM < MIN_VAL(M)) THEN
435                             MIN_VAL(M) = lSUM
436                             MIN_LOC(M) = IJK
437                          ENDIF
438                          IF (lSUM > MAX_VAL(M)) THEN
439                             MAX_VAL(M) = lSUM
440                             MAX_LOC(M) = IJK
441                          ENDIF
442                          IF (lSUM < 0.9) THEN
443                             COUNT(M,1) = COUNT(M,1) + 1 ! < 0.9
444                          ELSE IF (lSUM < 0.99) THEN
445                             COUNT(M,2) = COUNT(M,2) + 1 ! 0.9    - 0.99
446                          ELSE IF (lSUM < 0.999) THEN
447                             COUNT(M,3) = COUNT(M,3) + 1 ! 0.99   - 0.999
448                          ELSE IF (lSUM < 0.9999) THEN
449                             COUNT(M,4) = COUNT(M,4) + 1 ! 0.999  - 0.9999
450                          ELSE IF (lSUM < 1.0001) THEN
451                             COUNT(M,5) = COUNT(M,5) + 1 ! 0.9999 - 1.0001
452                          ELSE IF (lSUM < 1.001) THEN
453                             COUNT(M,6) = COUNT(M,6) + 1 ! 1.0001 - 1.001
454                          ELSE IF (lSUM < 1.01) THEN
455                             COUNT(M,7) = COUNT(M,7) + 1 ! 1.001  - 1.01
456                          ELSE IF (lSUM < 1.1) THEN
457                             COUNT(M,8) = COUNT(M,8) + 1 ! 1.01   - 1.1
458                          ELSE
459                             COUNT(M,9) = COUNT(M,9) + 1 ! > 1.1
460                          ENDIF
461                       ENDIF
462                    ENDIF
463                 ENDDO
464     
465     
466              ENDIF ! IF(.NOT.WALL_AT(IJK))
467           ENDDO ! DO I = ISTART2, IEND2
468           ENDDO ! DO J = JSTART2, JEND2
469           ENDDO ! DO K = KSTART2, KEND2
470     
471     
472           CALL GLOBAL_ALL_SUM(COUNT)
473           CALL GLOBAL_ALL_MAX(MAX_VAL)
474           CALL GLOBAL_ALL_MIN(MIN_VAL)
475     
476     
477      4100 FORMAT(//'Statistics of sum of gas species mass fraction',/      &
478              1X,'Minimum sum of X_g=',G12.5,/1X,'Maximum sum of X_g=',     &
479              G12.5,2/,3x,'Sum of X_g',T20,'No of Cells  Distribution')
480     
481           IF(REPORT_SPECIES(0)) THEN
482              lSUM = SUM(COUNT(0,:))
483              WRITE(ERR_MSG,4100) MIN_VAL(0), MAX_VAL(0)
484              CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
485              WRITE(ERR_MSG,4999) (COUNT(0,L),DBLE(COUNT(0,L))/lSUM,L=1,9)
486              CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
487           ENDIF
488     
489     
490      4200 FORMAT(//'Statistics of sum of solids phase ',I2,' species ',    &
491              'mass fractions',/1X,'Minimum sum of X_s=',G12.5,/1X,         &
492              'Maximum sum of X_s=',G12.5,2/,3x,'Sum of ',A,T20,'No of ',   &
493              'Cells',2x,'Distribution')
494     
495           DO M=1,SMAX
496              IF(REPORT_SPECIES(M)) THEN
497                 lSUM = SUM(COUNT(M,:))
498                 WRITE(ERR_MSG,4200) M,  MIN_VAL(M), MAX_VAL(M), &
499                    trim(iVAR('X_s',M))
500                 CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
501                 WRITE(ERR_MSG,4999) (COUNT(M,L),DBLE(COUNT(M,L))/lSUM,L=1,9)
502                 CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
503              ENDIF
504           ENDDO
505     
506           WRITE(ERR_MSG,"(/'End of report.')")
507           CALL FLUSH_ERR_MSG(HEADER=.FALSE.)
508     
509      4999 FORMAT(                                    &
510              1X,'<0.9            ',T20,I9,T33,G12.5,/&
511              1X,' 0.9    - 0.99  ',T20,I9,T33,G12.5,/&
512              1X,' 0.99   - 0.999 ',T20,I9,T33,G12.5,/&
513              1X,' 0.999  - 0.9999',T20,I9,T33,G12.5,/&
514              1X,' 0.9999 - 1.0001',T20,I9,T33,G12.5,/&
515              1X,' 1.0001 - 1.001 ',T20,I9,T33,G12.5,/&
516              1X,' 1.001  - 1.01  ',T20,I9,T33,G12.5,/&
517              1X,' 1.01   - 1.1   ',T20,I9,T33,G12.5,/&
518              1X,'>1.1            ',T20,I9,T33,G12.5)
519     
520           CALL FINL_ERR_MSG
521     
522           RETURN
523           END SUBROUTINE REPORT_SPECIES_STATS
524     
525     
526     
527     !----------------------------------------------------------------------!
528     !                                                                      !
529     !  Module name: CHECK_RXN_MASS_BALANCE                                 !
530     !  Author: M. Syamlal                                 Date: 27-OCT-92  !
531     !                                                                      !
532     !  Purpose: Verify that the net interphase mass transfer rates sum to  !
533     !  zero. This check is not needed with updated rates implementation.   !
534     !                                                                      !
535     !----------------------------------------------------------------------!
536           SUBROUTINE CHECK_RXN_MASS_BALANCE
537     !
538     !-----------------------------------------------
539     !   M o d u l e s
540     !-----------------------------------------------
541           USE param
542           USE param1
543           USE toleranc
544           USE fldvar
545           USE rxns
546           USE visc_s
547           USE visc_g
548           USE geometry
549           USE run
550           USE constant
551           USE physprop
552           USE indices
553           USE funits
554           USE mpi_utility
555           USE discretelement
556           USE mms
557           USE functions
558           use compar, only: DEAD_CELL_AT
559     
560           IMPLICIT NONE
561     !-----------------------------------------------
562     !   G l o b a l   P a r a m e t e r s
563     !-----------------------------------------------
564     !-----------------------------------------------
565     !   L o c a l   P a r a m e t e r s
566     !-----------------------------------------------
567     !-----------------------------------------------
568     !   L o c a l   V a r i a b l e s
569     !-----------------------------------------------
570     !
571     ! Loop counters:
572           INTEGER :: I, J, K, IJK
573           INTEGER :: M, L, LM
574     ! Temp variable for summations.
575           DOUBLE PRECISION :: lSUM
576     
577           INTEGER :: COUNT(DIMENSION_M+2, 2)
578           INTEGER :: MAX_LOC(DIMENSION_M+2)
579           DOUBLE PRECISION :: MAX_VAL(DIMENSION_M+2)
580     !......................................................................!
581     
582           CALL INIT_ERR_MSG('CHECK_RXN_MASS_BALANCE')
583     
584     
585           COUNT = 0
586           MAX_LOC = 0
587           MAX_VAL = ZERO
588     
589     
590           DO K = KSTART2, KEND2
591           DO J = JSTART2, JEND2
592           DO I = ISTART2, IEND2
593              IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
594              IJK = FUNIJK(I,J,K)
595              IF (.NOT.WALL_AT(IJK)) THEN
596     
597     !---------------------------------------------------------------------//
598     
599     ! The rate of interphase mass transfer must sum to zero over all phases.
600                 lSUM = SUM_R_G(IJK)
601                 IF(SMAX > 0) lSUM = lSUM + sum(SUM_R_S(IJK,1:SMAX))
602     
603                 IF (ABS(lSUM) > SMALL_NUMBER) THEN
604                    IF(abs(lSUM) > abs(MAX_VAL(1))) THEN
605                       MAX_VAL(1) = lSUM
606                       MAX_LOC(1) = IJK
607                    ENDIF
608     ! Bin the fluid cells: Good/Bad
609                    IF (ABS(lSUM) > TOL_COM) THEN
610                       COUNT(1,2) = COUNT(1,2) + 1
611                    ELSE
612                       COUNT(1,1) = COUNT(1,1) + 1
613                    ENDIF
614                 ENDIF
615     
616     ! Verify that the net rate of production (SUM_R_x) matches the the total
617     ! amount of mass transferred from other phases.
618                 DO L = 0, SMAX
619                    IF (L == 0) THEN
620                       lSUM = SUM_R_G(IJK)
621                    ELSE
622                       lSUM = SUM_R_S(IJK,L)
623                    ENDIF
624                    DO M = 0, SMAX
625                       IF (M > L) THEN
626                          LM = L + 1 + (M - 1)*M/2
627                          lSUM = lSUM - R_PHASE(IJK,LM)
628                       ELSE IF (L > M) THEN
629                          LM = M + 1 + (L - 1)*L/2
630                          lSUM = lSUM + R_PHASE(IJK,LM)
631                       ENDIF
632                    ENDDO
633     
634                    IF(ABS(lSUM) > SMALL_NUMBER) THEN
635                       IF(abs(lSUM) > abs(MAX_VAL(L+2))) THEN
636                          MAX_VAL(L+2) = lSUM
637                          MAX_LOC(L+2) = IJK
638                       ENDIF
639     
640     ! Force an exit if an imbalance exceeds the tollerance.
641                       IF(ABS(lSUM) > TOL_COM) then
642                          COUNT(L+2,2) = COUNT(L+2,2) + 1
643                       ELSE
644     ! Count the number of cells that do not sum to one, but are below the
645     ! tollerance for force and exit.
646                          COUNT(L+2,1) = COUNT(L+2,1) + 1
647                       ENDIF
648                    ENDIF
649                 END DO
650     
651              ENDIF ! IF(.NOT.WALL_AT(IJK))
652           ENDDO ! DO I = ISTART2, IEND2
653           ENDDO ! DO J = JSTART2, JEND2
654           ENDDO ! DO K = KSTART2, KEND2
655     
656           CALL GLOBAL_ALL_SUM(COUNT)
657     
658           IF(sum(COUNT(:,2)) > 0) THEN
659     
660              CALL GLOBAL_ALL_MAX(MAX_VAL)
661     
662              WRITE(ERR_MSG,5000)
663              CALL FLUSH_ERR_MSG(FOOTER=.FALSE.)
664     
665              IF(COUNT(1,2) > 0) THEN
666                 WRITE(ERR_MSG, 5100) COUNT(1,1), COUNT(1,2), MAX_VAL(1)
667                 CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
668              ENDIF
669     
670              DO L=0,SMAX
671                 IF(COUNT(L+2,2) > 0) THEN
672                    WRITE(ERR_MSG, 5200) L, COUNT(L+2,1), COUNT(L+2,2),&
673                       MAX_VAL(L+2)
674                    CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
675                 ENDIF
676              ENDDO
677     
678              WRITE(ERR_MSG,"(/'End of report.')")
679              CALL FLUSH_ERR_MSG(HEADER=.FALSE., ABORT=.TRUE.)
680           ENDIF
681     
682           CALL FINL_ERR_MSG
683     
684           RETURN
685      5000 FORMAT(5X,'Time = ',G12.5,/&
686              'Message: One or more of the following errors detected:',/    &
687              '  1. Discrepancies in the reaction rates.',/                 &
688              '  2. The rate of production of phases (SUM_R_g or SUM_R_s)',/&
689              '     and the interphase mass transfer rates (R_Phase) are',/ &
690              '     inconsistent (in subroutine RRATES).',/4X,'I',T14,'J',  &
691              T24,'K',T34,'M',T45,'Value')
692     
693      5100 FORMAT(//'Sum of all the reaction rates is not zero!',/,         &
694              'Number of cells with discrepancy < error tolerance = ',I5,/, &
695              'Number of cells with discrepancy > error tolerance = ',I5,/, &
696              'Maximum discrepancy = ',G12.5)
697     
698      5200 FORMAT(//'Mesage: Production of phase ',I2,' not equal to ',     &
699              'total mass transfer',/'from other phases!',/                 &
700              'Number of cells with discrepancy < error tolerance = ',I9,/  &
701              'Number of cells with discrepancy > error tolerance = ',I9,/  &
702              'Maximum discrepancy = ',G12.4)
703     
704     
705           END SUBROUTINE CHECK_RXN_MASS_BALANCE
706     
707     
708     
709     !----------------------------------------------------------------------!
710     !  Subroutine: REPORT_ERROR                                            !
711     !                                                                      !
712     !  Purpose: Manage error messages for CHECK_DATA_30.                   !
713     !----------------------------------------------------------------------!
714           SUBROUTINE REPORT_ERROR(pIER, pI, pJ, pK, VAL,  RELATION, BND, &
715              VAR, LC1, LC2)
716     
717           INTEGER, INTENT(INOUT) :: pIER
718           INTEGER, INTENT(IN) :: pI, pJ, pK
719           DOUBLE PRECISION, INTENT(IN) :: BND, VAL
720           CHARACTER(LEN=*), INTENT(IN) :: RELATION
721           CHARACTER(LEN=*), INTENT(IN) :: VAR
722           INTEGER, INTENT(IN), OPTIONAL :: LC1, LC2
723           CHARACTER(LEN=32) :: VAR_FULL
724     
725           INTEGER :: lIER
726     
727           VAR_FULL=''
728           IF(PRESENT(LC2)) THEN
729              VAR_FULL = iVAR(VAR,LC1,LC2)
730           ELSEIF(PRESENT(LC1)) THEN
731              VAR_FULL = iVAR(VAR,LC1)
732           ELSE
733              VAR_FULL = VAR
734           ENDIF
735     
736           IF(pIER == 0) THEN
737              CALL OPEN_PE_LOG(lIER)
738              SELECT CASE(ERR_TAG)
739              CASE(2000); WRITE(ERR_MSG,2000)
740              CASE(3000); WRITE(ERR_MSG,3000)
741              END SELECT
742              CALL FLUSH_ERR_MSG(FOOTER=.FALSE.)
743              pIER = 1
744           ENDIF
745     
746      2000 FORMAT('Error 2000: Physical properties detected in flow cells.',&
747              2/3X,'I',6x,'J',6x,'K',5x,'Value',8x,A,2x,'Bound',5x,'Variable')
748     
749      3000 FORMAT('Error 3000: Unphysical field variables detected.',&
750              2/3X,'I',6x,'J',6x,'K',5x,'Value',8x,A,2x,'Bound',5x,'Variable')
751     
752           WRITE(ERR_MSG,9000) pI, pJ, pK, VAL, RELATION, BND, trim(VAR_FULL)
753           CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
754     
755      9000 FORMAT(3(I6,1X),G12.4,1X,A,G12.4,1X,A)
756     
757           RETURN
758           END SUBROUTINE REPORT_ERROR
759           END SUBROUTINE CHECK_DATA_30
760