parthenon-hpc-lab / athenapk

AthenaPK: a performance portable version of Athena++ built on Parthenon and Kokkos
BSD 3-Clause "New" or "Revised" License
55 stars 13 forks source link

Split of cooling time cfl #11

Open pgrete opened 1 year ago

pgrete commented 1 year ago

The adaptive and Townsend integrators do not require a limitation on the timestep a priori. Nevertheless, from a physical point of view it may be desired to limit the timestep anyway so that the thermal dynamics and other dynamics don't decouple. While this is currently supported/tied to the adaptive rk12, and rk45 cooling integrators, we should separate that option and properly document it.

BenWibking commented 1 year ago

I tried this for my cloud-crushing simulations and it was extraordinarily expensive.

I think it would be worth considering either fixed-point iteration to reduce the operator splitting error, e.g.: https://amrex-combustion.github.io/PeleC/Algorithms.html#standard-time-advance.

Is the current method first-order operator splitting, or Strang splitting? Strang should reduce the operator-splitting error a bit if that's not already done.

pgrete commented 1 year ago

I tried this for my cloud-crushing simulations and it was extraordinarily expensive.

What exactly? Limiting dt to e/edot = 0.1?

I think it would be worth considering either fixed-point iteration to reduce the operator splitting error, e.g.: https://amrex-combustion.github.io/PeleC/Algorithms.html#standard-time-advance.

Would that result in just calling the source function in every stage of the integrator?

Is the current method first-order operator splitting, or Strang splitting? Strang should reduce the operator-splitting error a bit if that's not already done.

Currently, it's first order, but the infrastructure to swap it to Strang splitting is already in place, if required.

BenWibking commented 1 year ago

I tried this for my cloud-crushing simulations and it was extraordinarily expensive.

What exactly? Limiting dt to e/edot = 0.1?

Yes. I was perhaps in an unfavorable parameter regime.

I think it would be worth considering either fixed-point iteration to reduce the operator splitting error, e.g.: https://amrex-combustion.github.io/PeleC/Algorithms.html#standard-time-advance.

Would that result in just calling the source function in every stage of the integrator?

One could do that, but the way it's implemented in PeleC is to save the source term from the last timestep and use that in each stage of the hydro RK2 update. Then at the end, compute the final stage hydro update and compute an implicit update with the hydro flux as a source term, repeating until it converges.

Is the current method first-order operator splitting, or Strang splitting? Strang should reduce the operator-splitting error a bit if that's not already done.

Currently, it's first order, but the infrastructure to swap it to Strang splitting is already in place, if required.

BenWibking commented 1 year ago

I realized the PeleC algorithm could cause the internal energy to become negative for simulations with strong cooling, so it might not be applicable unless we recompute the source term each stage.

forrestglines commented 1 year ago

cooling/cfl already limits the time step, i.e. cooling/cfl=0.1 will set dt>= 0.1*min(e/edot)

However, it depends on a working DeDt in the cooling module which assumes an evenly spaced cooling table and uses linear interpolation in log space. It's not very tied to any integration scheme.

There's also the issue when the cooling times are very short like in the ICM simulations as the gas is cooling down to the temperature floor. However, since the cooling goes to zero at the temperature floor, we should only be taking on the order of 1/cooling cfl more time steps as the temperature goes to the floor.