diff --git a/integration/VODE/actual_integrator.H b/integration/VODE/actual_integrator.H index 22b341dfe..05e36496d 100644 --- a/integration/VODE/actual_integrator.H +++ b/integration/VODE/actual_integrator.H @@ -22,7 +22,8 @@ void actual_integrator (BurnT& state, amrex::Real dt, bool is_retry=false) auto istate = dvode(state, vode_state); - integrator_cleanup(vode_state, state, istate, state_save, dt); + integrator_cleanup(vode_state, state, istate, state_save, dt, + vode_final_state_species_failure_tolerance_factor); } diff --git a/integration/VODE/actual_integrator_sdc.H b/integration/VODE/actual_integrator_sdc.H index 8a1a362b9..c52683eca 100644 --- a/integration/VODE/actual_integrator_sdc.H +++ b/integration/VODE/actual_integrator_sdc.H @@ -27,7 +27,8 @@ void actual_integrator (BurnT& state, amrex::Real dt, bool is_retry=false) auto istate = dvode(state, vode_state); state.error_code = istate; - integrator_cleanup(vode_state, state, istate, state_save, dt); + integrator_cleanup(vode_state, state, istate, state_save, dt, + vode_final_state_species_failure_tolerance_factor); } diff --git a/integration/VODE/vode_type.H b/integration/VODE/vode_type.H index a1df1f9bf..b8c73344b 100644 --- a/integration/VODE/vode_type.H +++ b/integration/VODE/vode_type.H @@ -22,6 +22,10 @@ constexpr amrex::Real HMIN = 0.0_rt; constexpr amrex::Real vode_increase_change_factor = 4.0_rt; constexpr amrex::Real vode_decrease_change_factor = 0.25_rt; +// The interpolation back to tout is not monotonic, so the final state +// allows slightly more negativity than internal integration nodes. +constexpr amrex::Real vode_final_state_species_failure_tolerance_factor = 1.5_rt; + // For the backward differentiation formula (BDF) integration // the maximum order should be no greater than 5. constexpr int VODE_MAXORD = 5; diff --git a/integration/integrator_setup_sdc.H b/integration/integrator_setup_sdc.H index 2eac92bef..65067572b 100644 --- a/integration/integrator_setup_sdc.H +++ b/integration/integrator_setup_sdc.H @@ -135,7 +135,8 @@ state_backup_t integrator_backup (const BurnT& state) { template AMREX_GPU_HOST_DEVICE AMREX_INLINE void integrator_cleanup (IntegratorT& int_state, BurnT& state, - int istate, const state_backup_t& state_save, amrex::Real dt) + int istate, const state_backup_t& state_save, amrex::Real dt, + const amrex::Real final_state_species_failure_tolerance_factor = 1.0_rt) { // Copy the integration data back to the burn state. @@ -185,8 +186,11 @@ void integrator_cleanup (IntegratorT& int_state, BurnT& state, state.success = false; } + const amrex::Real final_state_species_failure_tolerance = + final_state_species_failure_tolerance_factor * state.rho * integrator_rp::species_failure_tolerance; + for (int n = 0; n < NumSpec; ++n) { - if (state.y[SFS+n] < -state.rho * integrator_rp::species_failure_tolerance) { + if (state.y[SFS+n] < -final_state_species_failure_tolerance) { state.success = false; } diff --git a/integration/integrator_setup_strang.H b/integration/integrator_setup_strang.H index 60f874ab4..6793630e1 100644 --- a/integration/integrator_setup_strang.H +++ b/integration/integrator_setup_strang.H @@ -119,7 +119,8 @@ state_backup_t integrator_backup (const BurnT& state) { template AMREX_GPU_HOST_DEVICE AMREX_INLINE void integrator_cleanup (IntegratorT& int_state, BurnT& state, - int istate, const state_backup_t& state_save, amrex::Real dt) + int istate, const state_backup_t& state_save, amrex::Real dt, + const amrex::Real final_state_species_failure_tolerance_factor = 1.0_rt) { // Copy the integration data back to the burn state. @@ -156,8 +157,11 @@ void integrator_cleanup (IntegratorT& int_state, BurnT& state, state.success = false; } + const amrex::Real final_state_species_failure_tolerance = + final_state_species_failure_tolerance_factor * integrator_rp::species_failure_tolerance; + for (int n = 1; n <= NumSpec; ++n) { - if (int_state.y(n) < -integrator_rp::species_failure_tolerance) { + if (int_state.y(n) < -final_state_species_failure_tolerance) { state.success = false; }