Skip to content

Commit

Permalink
fix logic for timestep size when not substepping (erf-model#1794)
Browse files Browse the repository at this point in the history
  • Loading branch information
@asalmgren
asalmgren authored Sep 8, 2024
1 parent e0d3772 commit e7e76f3
Show file tree
Hide file tree
Showing 2 changed files with 53 additions and 36 deletions.
39 changes: 21 additions & 18 deletions Source/TimeIntegration/ERF_ComputeTimestep.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -173,27 +173,30 @@ ERF::estTimeStep (int level, long& dt_fast_ratio) const
}
}

if (fixed_dt > 0. && fixed_fast_dt > 0.) {
dt_fast_ratio = static_cast<long>( fixed_dt / fixed_fast_dt );
} else if (fixed_dt > 0.) {
dt_fast_ratio = static_cast<long>( std::ceil((fixed_dt/estdt_comp)) );
} else {
dt_fast_ratio = (estdt_lowM_inv > 0.0) ? static_cast<long>( std::ceil((estdt_lowM/estdt_comp)) ) : 1;
}
if (!l_no_substepping) {
if (fixed_dt > 0. && fixed_fast_dt > 0.) {
dt_fast_ratio = static_cast<long>( fixed_dt / fixed_fast_dt );
} else if (fixed_dt > 0.) {
dt_fast_ratio = static_cast<long>( std::ceil((fixed_dt/estdt_comp)) );
} else {
dt_fast_ratio = (estdt_lowM_inv > 0.0) ? static_cast<long>( std::ceil((estdt_lowM/estdt_comp)) ) : 1;
}

// Force time step ratio to be an even value
if (solverChoice.force_stage1_single_substep) {
if ( dt_fast_ratio%2 != 0) dt_fast_ratio += 1;
} else {
if ( dt_fast_ratio%6 != 0) {
Print() << "mri_dt_ratio = " << dt_fast_ratio
<< " not divisible by 6 for N/3 substeps in stage 1" << std::endl;
dt_fast_ratio = static_cast<int>(std::ceil(dt_fast_ratio/6.0) * 6);
// Force time step ratio to be an even value
if (solverChoice.force_stage1_single_substep) {
if ( dt_fast_ratio%2 != 0) dt_fast_ratio += 1;
} else {
if ( dt_fast_ratio%6 != 0) {
Print() << "mri_dt_ratio = " << dt_fast_ratio
<< " not divisible by 6 for N/3 substeps in stage 1" << std::endl;
dt_fast_ratio = static_cast<int>(std::ceil(dt_fast_ratio/6.0) * 6);
}
}
}

if (verbose && !l_no_substepping)
Print() << "smallest even ratio is: " << dt_fast_ratio << std::endl;
if (verbose) {
Print() << "smallest even ratio is: " << dt_fast_ratio << std::endl;
}
} // if substepping

if (fixed_dt > 0.0) {
return fixed_dt;
Expand Down
50 changes: 32 additions & 18 deletions Source/TimeIntegration/ERF_MRI.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -190,16 +190,14 @@ public:
using namespace amrex;

// *******************************************************************************
// version == 0: we only update the fast variables every fast timestep, then update
// the slow variables after the acoustic sub-stepping. This has
// two calls to slow_rhs so that we can update the slow variables
// with the velocity field after the acoustic substepping using
// the time-averaged velocity from the substepping
// version == 1: we don't do any acoustic subcyling so we only make one call per RK
// stage to slow_rhs
// !no_substepping: we only update the fast variables every fast timestep, then update
// the slow variables after the acoustic sub-stepping. This has
// two calls to slow_rhs so that we can update the slow variables
// with the velocity field after the acoustic substepping using
// the time-averaged velocity from the substepping
// no_substepping: we don't do any acoustic subcyling so we only make one call per RK
// stage to slow_rhs
// *******************************************************************************
int version = no_substepping ? 1 : 0;

timestep = time_step;

const int substep_ratio = get_slow_fast_timestep_ratio();
Expand All @@ -208,8 +206,6 @@ public:
AMREX_ALWAYS_ASSERT(substep_ratio > 1 && substep_ratio % 2 == 0);
}

const amrex::Real sub_timestep = timestep / substep_ratio;

// Assume before advance() that S_old is valid data at the current time ("time" argument)
// And that if data is a MultiFab, both S_old and S_new contain ghost cells for evaluating a stencil based RHS
// We need this from S_old. This is convenient for S_new to have so we can use it
Expand Down Expand Up @@ -245,6 +241,8 @@ public:
amrex::Real time_stage = time;
amrex::Real old_time_stage;

const amrex::Real sub_timestep = timestep / substep_ratio;

if (!incompressible) {
// RK3 for compressible integrator
for (int nrk = 0; nrk < 3; nrk++)
Expand All @@ -255,13 +253,29 @@ public:
old_time_stage = time_stage;

if (nrk == 0) {
if (force_stage1_single_substep)
{ nsubsteps = 1; dtau = timestep / 3.0; time_stage = time + timestep / 3.0;}
else
{ nsubsteps = substep_ratio/3; dtau = sub_timestep ; time_stage = time + timestep / 3.0;}
if (force_stage1_single_substep || no_substepping) {
nsubsteps = 1; dtau = timestep / 3.0;
} else {
nsubsteps = substep_ratio/3; dtau = sub_timestep ;
}
time_stage = time + timestep / 3.0;
}
if (nrk == 1) {
if (no_substepping) {
nsubsteps = 1; dtau = 0.5 * timestep;
} else {
nsubsteps = substep_ratio/2; dtau = sub_timestep;
}
time_stage = time + timestep / 2.0;
}
if (nrk == 2) {
if (no_substepping) {
nsubsteps = 1; dtau = timestep;
} else {
nsubsteps = substep_ratio; dtau = sub_timestep;
}
time_stage = time + timestep;
}
if (nrk == 1) { nsubsteps = substep_ratio/2; dtau = sub_timestep ; time_stage = time + timestep / 2.0;}
if (nrk == 2) { nsubsteps = substep_ratio; dtau = sub_timestep ; time_stage = time + timestep ;}

// step 1 starts with S_stage = S^n and we always start substepping at the old time
// step 2 starts with S_stage = S^* and we always start substepping at the old time
Expand All @@ -285,7 +299,7 @@ public:
// ****************************************************
// Acoustic substepping
// ****************************************************
if (version == 0)
if (!no_substepping)
{
// *******************************************************************************
// Update the fast variables
Expand Down

0 comments on commit e7e76f3

Please sign in to comment.