check_data_30.f

Go to the documentation of this file.

!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
!  Module name: CHECK_DATA_30                                          !
!  Author: M. Syamlal                                 Date: 27-OCT-92  !
!                                                                      !
!  Purpose: Check whether the sum of reaction rates is zero and the    !
!  sum of mass fractions is 1.0 and EP_g >= EP_Star.                   !
!           and EP_g >= EP_star. Set miscellaneous constants           !
!                                                                      !
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
      SUBROUTINE check_data_30

! Global variables: (common to sub-functions)
!---------------------------------------------------------------------//
      use run, only: energy_eq, species_eq
      use compar, only: istart2, iend2
      use compar, only: jstart2, jend2
      use compar, only: kstart2, kend2
      use physprop, only: smax, nmax
      use param1, only: zero, one
      use rxns, only: use_rrates
      use mms, only: use_mms
      use param, only: dimension_m

      use error_manager

      IMPLICIT NONE

! Local variables:
!---------------------------------------------------------------------//
! Flag for error message headers.
      INTEGER :: ERR_TAG
! Flags to report stats on species equations.
      LOGICAL :: REPORT_SPECIES(0:dimension_m)
!......................................................................!


! Check physical properties in inflow/outflow cells.
      IF(.NOT.use_mms) CALL check_flow_cell_props
! Verify physical values for field variables.
      CALL check_physical_bounds
! Generate species report if needed.
      IF(any(report_species)) CALL report_species_stats
! Check interphase mass transfer.
      IF(use_rrates) CALL check_rxn_mass_balance

      RETURN

      CONTAINS

!----------------------------------------------------------------------!
!  Module name: CHECK_FLOW_PROPS                                       !
!  Author: M. Syamlal                                 Date: 27-OCT-92  !
!                                                                      !
!  Purpose: Verify that inflow/outflow cells do not contain physical   !
!  properties for specified variables.                                 !
!                                                                      !
!----------------------------------------------------------------------!
      SUBROUTINE check_flow_cell_props

! Global variables:
!---------------------------------------------------------------------//
      use visc_s, only: mu_s, lambda_s
      use visc_g, only: mu_gt, lambda_gt
      use physprop, only: k_g, k_s
      use physprop, only: dif_s, dif_g

      use functions, only: funijk
      use functions, only: flow_at
      use compar, only: dead_cell_at

      IMPLICIT NONE

! Local variables:
!---------------------------------------------------------------------//
! Loop counters
      INTEGER :: I, J, K, IJK
      INTEGER :: M, N
! Integer error flag.
      INTEGER :: IER
!......................................................................!

      CALL init_err_msg("CHECK_FLOW_CELL_PROPS")

      ier = 0
      err_tag=2000

      DO k = kstart2, kend2
      DO j = jstart2, jend2
      DO i = istart2, iend2

         IF (dead_cell_at(i,j,k)) cycle  ! skip dead cells

         ijk = funijk(i,j,k)

         IF(flow_at(ijk)) THEN

! Turbulent viscosity of fluid phase.
            IF(mu_gt(ijk) /= zero) CALL report_error                   &
               (ier, i, j, k, mu_gt(ijk), '/=', zero, 'MU_GT')

! Granular second coefficient of viscosity.
            IF(lambda_gt(ijk) /= zero) CALL report_error               &
               (ier, i, j, k, lambda_gt(ijk), '/=', zero, 'LAMBDA_GT')
! Gas conductivity.
            IF(k_g(ijk) /= zero) CALL report_error                     &
               (ier, i, j, k, k_g(ijk),'/=',zero,'K_G')
! Diffusivity of gas species N.
            DO n = 1, nmax(0)
               IF(dif_g(ijk, n) /= zero) CALL report_error             &
                  (ier, i, j, k, dif_g(ijk,n),'/=',zero,'DIF_G',n)
            ENDDO

            DO m = 1, smax
! Granular first coefficient of (shear) viscosity.
               IF(mu_s(ijk, m) /= zero) CALL report_error              &
                  (ier, i, j, k, mu_s(ijk,m),'/=',zero,'MU_S',m)
! Granular second coefficient of viscosity.
               IF(lambda_s(ijk, m) /= zero) CALL report_error          &
                  (ier, i, j, k, lambda_s(ijk,m),'/=',zero,'LAMBDA_S',m)
! Solids thermal conductivity.
               IF(k_s(ijk, m) /= zero) CALL report_error               &
                  (ier, i, j, k, k_s(ijk,m),'/=',zero,'K_S', m)
! Diffusivity of solids phase M, species N.
               DO n = 1, nmax(m)
                  IF(dif_s(ijk,m,n) /= zero) CALL report_error         &
                     (ier, i,j,k, dif_s(ijk,m,n),'/=',zero,'DIF_S',m,n)
               ENDDO
            ENDDO
         ENDIF

      ENDDO
      ENDDO
      ENDDO

      IF(ier /= 0) THEN
         WRITE(err_msg,"('End of Report.')")
         CALL flush_err_msg(header=.false., abort=.true.)
! Close DMP logs when running in interactive mode.
         CALL close_pe_log
      ENDIF

      CALL finl_err_msg

      RETURN
      END SUBROUTINE check_flow_cell_props



!----------------------------------------------------------------------!
!                                                                      !
!  Module name: CHECK_PHYSICAL_BOUNDS                                  !
!  Author: M. Syamlal                                 Date: 27-OCT-92  !
!                                                                      !
!  Purpose: Verify that fluid cells have physical values for the       !
!  specified variables.                                                !
!                                                                      !
!----------------------------------------------------------------------!
      SUBROUTINE check_physical_bounds

! Global variables:
!---------------------------------------------------------------------//
      use toleranc, only: tmin, tmax, tol_com
      use fldvar, only: x_g, x_s
      use fldvar, only: t_g, t_s
      use fldvar, only: rop_s
      use rxns, only: rox_gc, rox_sc
      use rxns, only: r_gp, r_sp
      use visc_s, only: mu_s
      use physprop, only: k_g, k_s
      use physprop, only: c_pg, c_ps
      use physprop, only: dif_s, dif_g
      use physprop, only: mu_g
      use physprop, only: mw_mix_g

      use mpi_utility, only: global_all_sum
      use mpi_utility, only: global_all_or

      use functions, only: funijk
      use functions, only: wall_at
      use compar, only: dead_cell_at

      IMPLICIT NONE

! Local variables:
!---------------------------------------------------------------------//
      INTEGER :: I, J, K, IJK
      INTEGER :: M, N
! Integer error flag.
      INTEGER :: IER

      CALL init_err_msg("CHECK_PHYSICAL_BOUNDS")

      ier = 0
      err_tag = 3000
      report_species = .false.

      DO k = kstart2, kend2
      DO j = jstart2, jend2
      DO i = istart2, iend2
         IF (dead_cell_at(i,j,k)) cycle  ! skip dead cells
         ijk = funijk(i,j,k)
         IF (.NOT.wall_at(ijk)) THEN

! Gas viscosity must be non-negative.
            IF(mu_g(ijk) < zero) CALL report_error                     &
               (ier, i, j, k, mu_g(ijk), '<', zero, 'MU_G')

! Mixture molecular weight must be positive.
            IF(mw_mix_g(ijk) <= zero) CALL report_error                &
               (ier, i, j, k,mw_mix_g(ijk), '<=', zero, 'MW_MIX_G')

! Verify thermodynamic properties when solving energy equations:
            IF(energy_eq) THEN
! Gas conductivity must be non-negative.
               IF(k_g(ijk) < zero) CALL report_error                   &
                  (ier, i, j, k, k_g(ijk), '<', zero, 'K_G')
! Gas phase specific heat must be positive.
               IF(c_pg(ijk) <= zero) CALL report_error                 &
                  (ier, i, j, k, c_pg(ijk), '<=', zero, 'C_PG')
! Verify that the gas phase temperature is within the bounds.
               IF(t_g(ijk) <= tmin) CALL report_error                  &
                  (ier, i, j, k, t_g(ijk), '<=', tmin, 'T_G')
               IF(t_g(ijk) >= tmax) CALL report_error                  &
                  (ier, i, j, k, t_g(ijk), '>=', tmax, 'T_G')

! Diffusivity of gas species N must be non-negative.
               DO n = 1, nmax(0)
                  IF(dif_g(ijk, n) < zero) CALL report_error           &
                     (ier, i,j,k, dif_g(ijk,n),'<',zero, 'DIF_G',n)
                  IF(x_g(ijk,n) < zero) CALL report_error              &
                     (ier, i, j, k, x_g(ijk,n), '<', zero, 'X_G',n)
                  IF(x_g(ijk,n) > one) CALL report_error               &
                     (ier, i, j, k, x_g(ijk,n), '>', one, 'X_G',n)
               ENDDO
            ENDIF

! Verify that the gas phase mass fractons sum to one.
            IF(species_eq(0)) THEN
               IF(abs(one - sum(x_g(ijk,1:nmax(0)))) > tol_com) &
                  report_species(0) = .true.

! Verify that the rates of formation and consumption adhear to the
! expected coding restraints. (non-negative)
               IF(use_rrates) THEN
                  DO n = 1, nmax(0)
                     IF(r_gp(ijk,n) < zero) CALL report_error          &
                       (ier,i,j,k, r_gp(ijk,n),'<',zero,'R_GP',n)
                     IF(rox_gc(ijk,n) < zero) CALL report_error        &
                       (ier,i,j,k,rox_gc(ijk,n),'<',zero,'RoX_GC',n)
                  ENDDO
               ENDIF

            ENDIF

            DO m = 1, smax

! Solids viscosity should be non-negativel.
               IF(mu_s(ijk,m) < zero) CALL report_error                &
                  (ier, i, j, k, mu_s(ijk,m), '<', zero, 'MU_S',m)

               IF(energy_eq) THEN

! Thermal conductivity must be non-negative.
                  IF(k_s(ijk,m) < zero) CALL report_error              &
                     (ier, i, j, k, k_s(ijk,m), '<', zero, 'K_S',m)

! Solids specific heat must be positive.
                  IF(c_ps(ijk,m) <= zero) CALL report_error            &
                     (ier, i, j, k, c_ps(ijk,m), '<=', zero, 'C_PS',m)

! Verify that the solids temperature is within the required bounds.
                  IF(t_s(ijk,m) <= tmin) CALL report_error             &
                     (ier, i, j, k, t_s(ijk,m), '<=', tmin, 'T_S',m)
                  IF(t_s(ijk,m) >= tmax) CALL report_error             &
                     (ier, i, j, k, t_s(ijk,m), '>=', tmax, 'T_S',m)

               ENDIF

               DO n = 1, nmax(m)
                  IF(dif_s(ijk, m, n) < zero) CALL report_error        &
                     (ier, i, j, k, dif_s(ijk,m,n),'<',zero,'DIF_S',m,n)
               ENDDO

! Sum of solids mass fractions should be one
               IF(species_eq(m)) THEN
                  IF(rop_s(ijk,m) /= zero) THEN
                     IF(abs(one-sum(x_s(ijk,m,1:nmax(m)))) > tol_com) &
                        report_species(m) = .true.
                  ENDIF

                  DO n = 1, nmax(m)
                     IF(x_s(ijk, m, n) < zero) CALL report_error       &
                        (ier, i,j,k, x_s(ijk,m,n),'<',zero, 'X_S',m,n)
                     IF(x_s(ijk, m, n) > one) CALL report_error        &
                        (ier, i,j,k, x_s(ijk,m,n),'>',one, 'X_S',m,n)
                  ENDDO
               ENDIF

! Verify that the rates of formation and consumption adhear to the
! expected coding restraints. (non-negative)
               IF(species_eq(m) .AND. use_rrates) THEN
                  DO n = 1, nmax(0)
                     IF(r_sp(ijk,m,n) < zero) CALL report_error        &
                       (ier,i,j,k,r_sp(ijk,m,n),'<',zero,'R_SP',m,n)
                     IF(rox_gc(ijk,n) < zero) CALL report_error        &
                       (ier,i,j,k,rox_sc(ijk,m,n),'<',zero,'RoX_SC',m,n)
                  ENDDO
               ENDIF
            ENDDO

         ENDIF ! IF(.NOT.WALL_AT(IJK))
      ENDDO ! DO I = ISTART2, IEND2
      ENDDO ! DO J = JSTART2, JEND2
      ENDDO ! DO K = KSTART2, KEND2

      CALL global_all_sum(ier)

      IF(ier /= 0) THEN
         WRITE(err_msg,"('End of Report.')")
         CALL flush_err_msg(header=.false., abort=.true.)
         CALL close_pe_log
      ENDIF

      CALL global_all_or(report_species)

      CALL finl_err_msg

      RETURN
      END SUBROUTINE check_physical_bounds



!----------------------------------------------------------------------!
!                                                                      !
!  Module name: REPORT_SPECIES_STATS                                   !
!  Author: M. Syamlal                                 Date: 27-OCT-92  !
!                                                                      !
!  Purpose: Collect information on the sums of species mass fractions. !
!  This routine is call as-needed.                                     !
!                                                                      !
!----------------------------------------------------------------------!
      SUBROUTINE report_species_stats

! Global variables:
!---------------------------------------------------------------------//
      use fldvar, only: x_g, x_s
      use fldvar, only: rop_s
      use run, only: time
      use mpi_utility, only: global_all_sum
      use mpi_utility, only: global_all_max
      use mpi_utility, only: global_all_min

      use functions, only: funijk
      use functions, only: wall_at

      use param1, only: undefined
      use compar, only: dead_cell_at

      IMPLICIT NONE

! Local variables:
!---------------------------------------------------------------------//
      INTEGER :: I, J, K, IJK, L
      INTEGER :: M

      INTEGER :: COUNT(0:dimension_m, 9)
      INTEGER :: MIN_LOC(0:dimension_m)
      INTEGER :: MAX_LOC(0:dimension_m)

      DOUBLE PRECISION :: MIN_VAL(0:dimension_m)
      DOUBLE PRECISION :: MAX_VAL(0:dimension_m)

      DOUBLE PRECISION :: lSUM
!......................................................................!

      CALL init_err_msg('REPORT_SPECIES_STATS')

      WRITE(err_msg, 4000) time
      CALL flush_err_msg(footer=.false.)

 4000 FORMAT('Time = ',g12.5,/'Warning: The sum of mass fractions ',   &
      'is not equal to one.')

! Initialize the counters
      count = 0
      max_val = -undefined
      min_val =  undefined

! Collect stats on species across the domain.
      DO k = kstart2, kend2
      DO j = jstart2, jend2
      DO i = istart2, iend2
         IF (dead_cell_at(i,j,k)) cycle  ! skip dead cells
         ijk = funijk(i,j,k)
         IF (.NOT.wall_at(ijk)) THEN

            IF(species_eq(0)) THEN
               lsum = sum(x_g(ijk,1:nmax(0)))
               IF (lsum < min_val(0)) THEN
                  min_val(0) = lsum
                  min_loc(0) = ijk
               ENDIF
               IF (lsum > max_val(0)) THEN
                  max_val(0) = lsum
                  max_loc(0) = ijk
               ENDIF
               IF (lsum < 0.9) THEN
                  count(0,1) = count(0,1) + 1 ! < 0.9pP
               ELSE IF (lsum < 0.99) THEN
                  count(0,2) = count(0,2) + 1 ! 0.9    - 0.99
               ELSE IF (lsum < 0.999) THEN
                  count(0,3) = count(0,3) + 1 ! 0.99   - 0.999
               ELSE IF (lsum < 0.9999) THEN
                  count(0,4) = count(0,4) + 1 ! 0.999  - 0.9999
               ELSE IF (lsum < 1.0001) THEN
                  count(0,5) = count(0,5) + 1 ! 0.9999 - 1.0001
               ELSE IF (lsum < 1.001) THEN
                  count(0,6) = count(0,6) + 1 ! 1.0001 - 1.001
               ELSE IF (lsum < 1.01) THEN
                  count(0,7) = count(0,7) + 1 ! 1.001  - 1.01
               ELSE IF (lsum < 1.1) THEN
                  count(0,8) = count(0,8) + 1 ! 1.01   - 1.1
               ELSE
                  count(0,9) = count(0,9) + 1 ! > 1.1
               ENDIF
            ENDIF


            DO m = 1, smax
!  Sum of solids mass fractions should be one
               IF(species_eq(m)) THEN
                  lsum = sum(x_s(ijk,m,1:nmax(m)) )
                  IF(rop_s(ijk,m) /= zero) THEN
                     IF(lsum < min_val(m)) THEN
                        min_val(m) = lsum
                        min_loc(m) = ijk
                     ENDIF
                     IF (lsum > max_val(m)) THEN
                        max_val(m) = lsum
                        max_loc(m) = ijk
                     ENDIF
                     IF (lsum < 0.9) THEN
                        count(m,1) = count(m,1) + 1 ! < 0.9
                     ELSE IF (lsum < 0.99) THEN
                        count(m,2) = count(m,2) + 1 ! 0.9    - 0.99
                     ELSE IF (lsum < 0.999) THEN
                        count(m,3) = count(m,3) + 1 ! 0.99   - 0.999
                     ELSE IF (lsum < 0.9999) THEN
                        count(m,4) = count(m,4) + 1 ! 0.999  - 0.9999
                     ELSE IF (lsum < 1.0001) THEN
                        count(m,5) = count(m,5) + 1 ! 0.9999 - 1.0001
                     ELSE IF (lsum < 1.001) THEN
                        count(m,6) = count(m,6) + 1 ! 1.0001 - 1.001
                     ELSE IF (lsum < 1.01) THEN
                        count(m,7) = count(m,7) + 1 ! 1.001  - 1.01
                     ELSE IF (lsum < 1.1) THEN
                        count(m,8) = count(m,8) + 1 ! 1.01   - 1.1
                     ELSE
                        count(m,9) = count(m,9) + 1 ! > 1.1
                     ENDIF
                  ENDIF
               ENDIF
            ENDDO


         ENDIF ! IF(.NOT.WALL_AT(IJK))
      ENDDO ! DO I = ISTART2, IEND2
      ENDDO ! DO J = JSTART2, JEND2
      ENDDO ! DO K = KSTART2, KEND2


      CALL global_all_sum(count)
      CALL global_all_max(max_val)
      CALL global_all_min(min_val)


 4100 FORMAT(//'Statistics of sum of gas species mass fraction',/      &
         1x,'Minimum sum of X_g=',g12.5,/1x,'Maximum sum of X_g=',     &
         g12.5,2/,3x,'Sum of X_g',t20,'No of Cells  Distribution')

      IF(report_species(0)) THEN
         lsum = sum(count(0,:))
         WRITE(err_msg,4100) min_val(0), max_val(0)
         CALL flush_err_msg(header=.false., footer=.false.)
         WRITE(err_msg,4999) (count(0,l),dble(count(0,l))/lsum,l=1,9)
         CALL flush_err_msg(header=.false., footer=.false.)
      ENDIF


 4200 FORMAT(//'Statistics of sum of solids phase ',i2,' species ',    &
         'mass fractions',/1x,'Minimum sum of X_s=',g12.5,/1x,         &
         'Maximum sum of X_s=',g12.5,2/,3x,'Sum of ',a,t20,'No of ',   &
         'Cells',2x,'Distribution')

      DO m=1,smax
         IF(report_species(m)) THEN
            lsum = sum(count(m,:))
            WRITE(err_msg,4200) m,  min_val(m), max_val(m), &
               trim(ivar('X_s',m))
            CALL flush_err_msg(header=.false., footer=.false.)
            WRITE(err_msg,4999) (count(m,l),dble(count(m,l))/lsum,l=1,9)
            CALL flush_err_msg(header=.false., footer=.false.)
         ENDIF
      ENDDO

      WRITE(err_msg,"(/'End of report.')")
      CALL flush_err_msg(header=.false.)

 4999 FORMAT(                                    &
         1x,'<0.9            ',t20,i9,t33,g12.5,/&
         1x,' 0.9    - 0.99  ',t20,i9,t33,g12.5,/&
         1x,' 0.99   - 0.999 ',t20,i9,t33,g12.5,/&
         1x,' 0.999  - 0.9999',t20,i9,t33,g12.5,/&
         1x,' 0.9999 - 1.0001',t20,i9,t33,g12.5,/&
         1x,' 1.0001 - 1.001 ',t20,i9,t33,g12.5,/&
         1x,' 1.001  - 1.01  ',t20,i9,t33,g12.5,/&
         1x,' 1.01   - 1.1   ',t20,i9,t33,g12.5,/&
         1x,'>1.1            ',t20,i9,t33,g12.5)

      CALL finl_err_msg

      RETURN
      END SUBROUTINE report_species_stats



!----------------------------------------------------------------------!
!                                                                      !
!  Module name: CHECK_RXN_MASS_BALANCE                                 !
!  Author: M. Syamlal                                 Date: 27-OCT-92  !
!                                                                      !
!  Purpose: Verify that the net interphase mass transfer rates sum to  !
!  zero. This check is not needed with updated rates implementation.   !
!                                                                      !
!----------------------------------------------------------------------!
      SUBROUTINE check_rxn_mass_balance
!
!-----------------------------------------------
!   M o d u l e s
!-----------------------------------------------
      USE param
      USE param1
      USE toleranc
      USE fldvar
      USE rxns
      USE visc_s
      USE visc_g
      USE geometry
      USE run
      USE constant
      USE physprop
      USE indices
      USE funits
      USE mpi_utility
      USE discretelement
      USE mms
      USE functions
      use compar, only: dead_cell_at

      IMPLICIT NONE
!-----------------------------------------------
!   G l o b a l   P a r a m e t e r s
!-----------------------------------------------
!-----------------------------------------------
!   L o c a l   P a r a m e t e r s
!-----------------------------------------------
!-----------------------------------------------
!   L o c a l   V a r i a b l e s
!-----------------------------------------------
!
! Loop counters:
      INTEGER :: I, J, K, IJK
      INTEGER :: M, L, LM
! Temp variable for summations.
      DOUBLE PRECISION :: lSUM

      INTEGER :: COUNT(dimension_m+2, 2)
      INTEGER :: MAX_LOC(dimension_m+2)
      DOUBLE PRECISION :: MAX_VAL(dimension_m+2)
!......................................................................!

      CALL init_err_msg('CHECK_RXN_MASS_BALANCE')


      count = 0
      max_loc = 0
      max_val = zero


      DO k = kstart2, kend2
      DO j = jstart2, jend2
      DO i = istart2, iend2
         IF (dead_cell_at(i,j,k)) cycle  ! skip dead cells
         ijk = funijk(i,j,k)
         IF (.NOT.wall_at(ijk)) THEN

!---------------------------------------------------------------------//

! The rate of interphase mass transfer must sum to zero over all phases.
            lsum = sum_r_g(ijk)
            IF(smax > 0) lsum = lsum + sum(sum_r_s(ijk,1:smax))

            IF (abs(lsum) > small_number) THEN
               IF(abs(lsum) > abs(max_val(1))) THEN
                  max_val(1) = lsum
                  max_loc(1) = ijk
               ENDIF
! Bin the fluid cells: Good/Bad
               IF (abs(lsum) > tol_com) THEN
                  count(1,2) = count(1,2) + 1
               ELSE
                  count(1,1) = count(1,1) + 1
               ENDIF
            ENDIF

! Verify that the net rate of production (SUM_R_x) matches the the total
! amount of mass transferred from other phases.
            DO l = 0, smax
               IF (l == 0) THEN
                  lsum = sum_r_g(ijk)
               ELSE
                  lsum = sum_r_s(ijk,l)
               ENDIF
               DO m = 0, smax
                  IF (m > l) THEN
                     lm = l + 1 + (m - 1)*m/2
                     lsum = lsum - r_phase(ijk,lm)
                  ELSE IF (l > m) THEN
                     lm = m + 1 + (l - 1)*l/2
                     lsum = lsum + r_phase(ijk,lm)
                  ENDIF
               ENDDO

               IF(abs(lsum) > small_number) THEN
                  IF(abs(lsum) > abs(max_val(l+2))) THEN
                     max_val(l+2) = lsum
                     max_loc(l+2) = ijk
                  ENDIF

! Force an exit if an imbalance exceeds the tollerance.
                  IF(abs(lsum) > tol_com) then
                     count(l+2,2) = count(l+2,2) + 1
                  ELSE
! Count the number of cells that do not sum to one, but are below the
! tollerance for force and exit.
                     count(l+2,1) = count(l+2,1) + 1
                  ENDIF
               ENDIF
            END DO

         ENDIF ! IF(.NOT.WALL_AT(IJK))
      ENDDO ! DO I = ISTART2, IEND2
      ENDDO ! DO J = JSTART2, JEND2
      ENDDO ! DO K = KSTART2, KEND2

      CALL global_all_sum(count)

      IF(sum(count(:,2)) > 0) THEN

         CALL global_all_max(max_val)

         WRITE(err_msg,5000)
         CALL flush_err_msg(footer=.false.)

         IF(count(1,2) > 0) THEN
            WRITE(err_msg, 5100) count(1,1), count(1,2), max_val(1)
            CALL flush_err_msg(header=.false., footer=.false.)
         ENDIF

         DO l=0,smax
            IF(count(l+2,2) > 0) THEN
               WRITE(err_msg, 5200) l, count(l+2,1), count(l+2,2),&
                  max_val(l+2)
               CALL flush_err_msg(header=.false., footer=.false.)
            ENDIF
         ENDDO

         WRITE(err_msg,"(/'End of report.')")
         CALL flush_err_msg(header=.false., abort=.true.)
      ENDIF

      CALL finl_err_msg

      RETURN
 5000 FORMAT(5x,'Time = ',g12.5,/&
         'Message: One or more of the following errors detected:',/    &
         '  1. Discrepancies in the reaction rates.',/                 &
         '  2. The rate of production of phases (SUM_R_g or SUM_R_s)',/&
         '     and the interphase mass transfer rates (R_Phase) are',/ &
         '     inconsistent (in subroutine RRATES).',/4x,'I',t14,'J',  &
         t24,'K',t34,'M',t45,'Value')

 5100 FORMAT(//'Sum of all the reaction rates is not zero!',/,         &
         'Number of cells with discrepancy < error tolerance = ',i5,/, &
         'Number of cells with discrepancy > error tolerance = ',i5,/, &
         'Maximum discrepancy = ',g12.5)

 5200 FORMAT(//'Mesage: Production of phase ',i2,' not equal to ',     &
         'total mass transfer',/'from other phases!',/                 &
         'Number of cells with discrepancy < error tolerance = ',i9,/  &
         'Number of cells with discrepancy > error tolerance = ',i9,/  &
         'Maximum discrepancy = ',g12.4)


      END SUBROUTINE check_rxn_mass_balance



!----------------------------------------------------------------------!
!  Subroutine: REPORT_ERROR                                            !
!                                                                      !
!  Purpose: Manage error messages for CHECK_DATA_30.                   !
!----------------------------------------------------------------------!
      SUBROUTINE report_error(pIER, pI, pJ, pK, VAL,  RELATION, BND, &
         var, lc1, lc2)

      INTEGER, INTENT(INOUT) :: pIER
      INTEGER, INTENT(IN) :: pI, pJ, pK
      DOUBLE PRECISION, INTENT(IN) :: BND, VAL
      CHARACTER(LEN=*), INTENT(IN) :: RELATION
      CHARACTER(LEN=*), INTENT(IN) :: VAR
      INTEGER, INTENT(IN), OPTIONAL :: LC1, LC2
      CHARACTER(LEN=32) :: VAR_FULL

      INTEGER :: lIER

      var_full=''
      IF(PRESENT(lc2)) THEN
         var_full = ivar(var,lc1,lc2)
      ELSEIF(PRESENT(lc1)) THEN
         var_full = ivar(var,lc1)
      ELSE
         var_full = var
      ENDIF

      IF(pier == 0) THEN
         CALL open_pe_log(lier)
         SELECT CASE(err_tag)
         CASE(2000); WRITE(err_msg,2000)
         CASE(3000); WRITE(err_msg,3000)
         END SELECT
         CALL flush_err_msg(footer=.false.)
         pier = 1
      ENDIF

 2000 FORMAT('Error 2000: Physical properties detected in flow cells.',&
         2/3x,'I',6x,'J',6x,'K',5x,'Value',8x,a,2x,'Bound',5x,'Variable')

 3000 FORMAT('Error 3000: Unphysical field variables detected.',&
         2/3x,'I',6x,'J',6x,'K',5x,'Value',8x,a,2x,'Bound',5x,'Variable')

      WRITE(err_msg,9000) pi, pj, pk, val, relation, bnd, trim(var_full)
      CALL flush_err_msg(header=.false., footer=.false.)

 9000 FORMAT(3(i6,1x),g12.4,1x,a,g12.4,1x,a)

      RETURN
      END SUBROUTINE report_error
      END SUBROUTINE check_data_30