diff --git a/Source/chem.f90 b/Source/chem.f90 index caad873b03a..8998f5d56c8 100644 --- a/Source/chem.f90 +++ b/Source/chem.f90 @@ -807,6 +807,7 @@ SUBROUTINE CVODE_SERIAL(CC,ZZ_0, TMP_IN, TMP_UNMIX, PR_IN, ZETA0, TAU_MIX, CELL_ TMP_OUT, CHEM_TIME, WRITE_SUBSTEPS, CVODE_CALL_OPTION) USE PHYSICAL_FUNCTIONS, ONLY : MOLAR_CONC_TO_MASS_FRAC, CALC_EQUIV_RATIO, GET_ENTHALPY, GET_MOLECULAR_WEIGHT USE COMP_FUNCTIONS, ONLY: GET_FILE_NUMBER +USE CHEMCONS, ONLY: CVODE_WARNING_CELLS,CVODE_ERR_CODE_MIN,CVODE_ERR_CODE_MAX USE GLOBAL_CONSTANTS USE FCVODE_MOD ! FORTRAN INTERFACE TO CVODE USE FSUNDIALS_CONTEXT_MOD ! FORTRAN INTERFACE TO SUNCONTEXT @@ -1015,15 +1016,18 @@ SUBROUTINE CVODE_SERIAL(CC,ZZ_0, TMP_IN, TMP_UNMIX, PR_IN, ZETA0, TAU_MIX, CELL_ ENDIF IF (IERR_C .NE. CV_SUCCESS) THEN - IF (IERR_C == CV_TOO_MUCH_WORK) THEN - WRITE(LU_ERR,'(A, 2E18.8, I8, A)')" WARN: CVODE took all internal substeps. CUR_CFD_TIME, DT, MAXTRY=", CUR_CFD_TIME, & - (TEND-TCUR), MAXTRY, ". If the warning persists, reduce the timestep." - ELSE - WRITE(LU_ERR,'(A, I4, A, 2E18.8, A)')" WARN: CVODE didn't finish ODE solution with message code:", IERR_C, & - " and CUR_CFD_TIME, DT=", CUR_CFD_TIME, (TEND-TCUR), ". If the warning persists, reduce the timestep." + IF(IERR_C>=CVODE_ERR_CODE_MIN .AND. IERR_C<=CVODE_ERR_CODE_MAX) THEN + CVODE_WARNING_CELLS(IERR_C) = CVODE_WARNING_CELLS(IERR_C) + 1 ENDIF - IF (DEBUG) THEN + IF (IERR_C == CV_TOO_MUCH_WORK) THEN + WRITE(LU_ERR,'(A, 2E18.8, I8, A)')" WARN: CVODE took all internal substeps. CUR_CFD_TIME, DT, MAXTRY=", CUR_CFD_TIME, & + (TEND-TCUR), MAXTRY, ". If the warning persists, reduce the timestep." + ELSE + WRITE(LU_ERR,'(A, I4, A, 2E18.8, A)')" WARN: CVODE didn't finish ODE solution with message code:", IERR_C, & + " and CUR_CFD_TIME, DT=", CUR_CFD_TIME, (TEND-TCUR), ". If the warning persists, reduce the timestep." + ENDIF + CALL MOLAR_CONC_TO_MASS_FRAC(CC(1:N_TRACKED_SPECIES), ZZ(1:N_TRACKED_SPECIES)) CALL CALC_EQUIV_RATIO(ZZ(1:N_TRACKED_SPECIES), EQUIV) DO NS = 1, N_TRACKED_SPECIES diff --git a/Source/cons.f90 b/Source/cons.f90 index 3421ee34133..3d8df244ef7 100644 --- a/Source/cons.f90 +++ b/Source/cons.f90 @@ -912,6 +912,10 @@ MODULE CHEMCONS INTEGER :: MAX_CVODE_SUBSTEPS=100000 INTEGER :: CVODE_MAX_TRY=4 INTEGER :: CVODE_ORDER=0 +INTEGER :: CVODE_ERR_CODE_MIN=-100 +INTEGER :: CVODE_ERR_CODE_MAX=100 +INTEGER :: CVODE_WARNING_CELLS(-100:100)! Index of the array is error code and value is Cell count +CHARACTER(LEN=100) :: CVODE_WARN_MESSAGES(-100:100) ! FOR WRITING CVODE SUBSTEPS LOGICAL :: WRITE_CVODE_SUBSTEPS = .FALSE. @@ -932,4 +936,12 @@ MODULE CHEMCONS INTEGER :: N_IGNITION_ZONES = 0 TYPE(IGNITION_ZONE_TYPE), DIMENSION(MAX_IGNITION_ZONES) :: IGNITION_ZONES !< Coordinates of ignition zones +CONTAINS + SUBROUTINE INIT_CVODE_WARN_MESSAGES() + CVODE_WARN_MESSAGES = 'CVODE didn''t finish ODE solution with this code.' + CVODE_WARN_MESSAGES(-1) = 'CVODE took all internal substeps.' + CVODE_WARN_MESSAGES(-3) = 'Minimum step size was reached.' + CVODE_WARN_MESSAGES(-4) = 'Convergence test failure.' + END SUBROUTINE INIT_CVODE_WARN_MESSAGES + END MODULE CHEMCONS diff --git a/Source/fire.f90 b/Source/fire.f90 index aca4c8f6d38..44a03dc6c7e 100644 --- a/Source/fire.f90 +++ b/Source/fire.f90 @@ -41,6 +41,7 @@ SUBROUTINE COMBUSTION_LOAD_BALANCED(T,DT) USE SOOT_ROUTINES, ONLY: SOOT_SURFACE_OXIDATION USE COMP_FUNCTIONS, ONLY: CURRENT_TIME +USE CHEMCONS, ONLY: CVODE_WARNING_CELLS,INIT_CVODE_WARN_MESSAGES REAL(EB), INTENT(IN) :: T,DT INTEGER :: NM,ICC,JCC REAL(EB) :: TNOW @@ -55,12 +56,14 @@ SUBROUTINE COMBUSTION_LOAD_BALANCED(T,DT) COMBUSTION_INIT = .TRUE. NVAR_TO_SEND = N_TRACKED_SPECIES +7 ! NS, TEMP, RHO, PRES, MU, DELTA, VOL, IGN_ZN NVAR_TO_RECEIVE = N_TRACKED_SPECIES +4 ! NS, Q_OUT, MIX_TIME_OUT, CHI_R_OUT, CHEM_SUBIT_TMP_OUT + CALL INIT_CVODE_WARN_MESSAGES() ENDIF DO NM=LOWER_MESH_INDEX,UPPER_MESH_INDEX CALL POINT_TO_MESH(NM) Q = 0._EB CHI_R = 0._EB + CVODE_WARNING_CELLS = 0 IF (CC_IBM) THEN DO ICC=1,MESHES(NM)%N_CUTCELL_MESH DO JCC=1,CUT_CELL(ICC)%NCELL @@ -109,6 +112,7 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT) REAL(EB) :: MIX_TIME_OUT, Q_OUT, CHI_R_OUT, MYTEMP, MYRHO, MYMU, DELTA, VOL, TNOW2 INTEGER :: IGN_ZN + Q_EXISTS = .FALSE. !------ @@ -329,6 +333,13 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT) ENDDO ENDIF +! Output accumulated CVODE warnings (instead of cell by cell) in verbose mode +#ifdef WITH_SUNDIALS +IF (COMBUSTION_ODE_SOLVER == CVODE_SOLVER .AND. VERBOSE) THEN + CALL DUMP_CVODE_WARNING_SUMMARY(NCHEM_ACTIVE_CELLS_AND_CC) +ENDIF +#endif + END SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED @@ -1269,6 +1280,60 @@ SUBROUTINE CALC_AFT_REAC_AND_PROD(ZZ,ZZ_REAC,ZZ_PROD) END SUBROUTINE CALC_AFT_REAC_AND_PROD +!> \brief Dump CVODE warning summary +SUBROUTINE DUMP_CVODE_WARNING_SUMMARY(NCHEM_ACTIVE_CELLS_AND_CC) + USE CHEMCONS, ONLY: CVODE_WARNING_CELLS,CVODE_ERR_CODE_MIN,CVODE_ERR_CODE_MAX,CVODE_WARN_MESSAGES + + INTEGER, INTENT(IN) :: NCHEM_ACTIVE_CELLS_AND_CC + INTEGER :: NCHEM_ACTIVE_CELLS_AND_CC_GLOBAL,IERR, CVODE_ERR_CODE, TOT_WARN_CELLS + REAL(EB) :: PCNT + CHARACTER(LEN=16) :: CELL_STR, PCNT_STR + CHARACTER(LEN=40) :: CELL_INFO + + IF (N_MPI_PROCESSES > 1) THEN + CALL MPI_ALLREDUCE(MPI_IN_PLACE, CVODE_WARNING_CELLS, SIZE(CVODE_WARNING_CELLS), & + MPI_INTEGER, MPI_SUM, MPI_COMM_WORLD, IERR) + + CALL MPI_ALLREDUCE(NCHEM_ACTIVE_CELLS_AND_CC, NCHEM_ACTIVE_CELLS_AND_CC_GLOBAL, & + 1, MPI_INTEGER, MPI_SUM, MPI_COMM_WORLD, IERR) + ELSE + NCHEM_ACTIVE_CELLS_AND_CC_GLOBAL = NCHEM_ACTIVE_CELLS_AND_CC + END IF + TOT_WARN_CELLS = SUM(CVODE_WARNING_CELLS) + + IF (TOT_WARN_CELLS > 0 .AND. NCHEM_ACTIVE_CELLS_AND_CC_GLOBAL > 0 .AND. MY_RANK == 0) THEN + + WRITE(LU_ERR, '(A)') '--------------- CVODE Warning Summary ----------------' + WRITE(LU_ERR, '(A)') ' Code Cells (% Chem) Message' + + DO CVODE_ERR_CODE = CVODE_ERR_CODE_MIN, CVODE_ERR_CODE_MAX + IF (CVODE_WARNING_CELLS(CVODE_ERR_CODE) > 0) THEN + + PCNT = REAL(CVODE_WARNING_CELLS(CVODE_ERR_CODE)) / & + REAL(NCHEM_ACTIVE_CELLS_AND_CC_GLOBAL) * 100.0_EB + + ! Convert to strings and left-align + WRITE(CELL_STR, '(I0)') CVODE_WARNING_CELLS(CVODE_ERR_CODE) + WRITE(PCNT_STR, '(F6.2)') PCNT + CELL_STR = ADJUSTL(CELL_STR) + PCNT_STR = ADJUSTL(PCNT_STR) + + ! Combine count + percentage into one left-aligned string + WRITE(CELL_INFO, '(A, " (", A, " %)")') TRIM(CELL_STR), TRIM(PCNT_STR) + CELL_INFO = ADJUSTL(CELL_INFO) + + ! Print one line of summary + WRITE(LU_ERR, '(1X, I4, 5X, A, 5X, A)') CVODE_ERR_CODE, & + TRIM(CELL_INFO), TRIM(CVODE_WARN_MESSAGES(CVODE_ERR_CODE)) + END IF + END DO + + WRITE(LU_ERR, '(A)') '------------------------------------------------------' + + END IF + +END SUBROUTINE DUMP_CVODE_WARNING_SUMMARY + #endif SUBROUTINE CHECK_AUTO_IGNITION(EXTINCT,TMP_IN,AIT,IIC,JJC,KKC,REAC_INDEX)