Re: Should timesteps for hydrodynamics be assigned after the force is computed?

From: Leonard Romano <leonard.romano_at_tum.de>
Date: Wed, 11 Aug 2021 12:20:14 +0200

Hi Volker,

I don't quite agree that the Courant condition would not be violated by
the treatment in Gadget-4.
In particular, I would like to note that the signal velocity for the
Courant timestep is computed during the force loop. Since major changes
in the temperature due to cooling or feedback would typically happen --
even on small time steps -- within the function
calculate_non_standard_physics_end_of_step(), by the time new time steps
are assigned the signal velocity hasn't been updated yet and so the
Courant condition can be violated for one full time step of the particle
in question.
The simple fix that made my simulations run much more stable (and got
rid of spurious high velocity particles), was to move the function call
"find_hydro_timesteps()" to the function "hydro_kick()" right after the
force-loop as

    if(step_indicator == FIRST_HALF_STEP)
         find_hydro_timesteps();

This way the most recent value of "MaxSignalVel" for the Courant time
step criterion and the most recent hydro force for the kinematic time
step criterion are available, which can substantially reduce numerical
instability.
For a while I also believed that the tree based time step criterion
would be sufficient for detecting feedback bombs, but I came to the
conclusion that this is only partially true as energy injections due to
feedback typically happen instantaneous and so it is basically
impossible to predict when and where they would happen (though it does
help to predict incoming waves at the location of particles that are far
away). Close to the energy injections I think that a time step limiter
similar to the one proposed by Saitoh & Makino (2009) and Durier & Dalla
Vecchia (2012) might be necessary. However, of course this closely
depends on how one injects energy.

Cheers,
Leonard

On 10.08.21 19:49, Volker Springel wrote:
> Hi Leonard,
>
>> On 31. Jul 2021, at 13:22, Leonard Romano <leonard.romano_at_tum.de> wrote:
>>
>> Dear gadget-list members,
>>
>> Recently Julianne Goddard raised an issue with crashes in cosmological simulations with subgrid physics (Cooling, star formation and feedback) due to exceedingly small timesteps
>> (https://wwwmpa.mpa-garching.mpg.de/gadget/gadget-list/0926.html).
>> I am running similar simulations with a custom feedback model, and often times run into similar problems. Now I noticed that a major difference between previous versions of Gadget and Gadget-4 is, that in Gadget-4 the time steps for hydrodynamics are assigned before the forces are computed. Without sub-grid physics this is probably fine as changes in quantities like Entropy and velocities are small, but with sub-grid physics there can be sudden changes by orders of magnitude and so the time steps before the force calculation might be heavily overestimated possibly leading to dramatic violations of the Courant-condition.
> Normally the Courant condition doesn't depend on the pressure forces, only on the sound speed and the bulk velocity of the gas (or more accurately the wave speeds). So computing hydro forces before deciding on the timestep doesn't directly influence this.
>
> I note that the code uses the tree_based_timesteps() function to determine the first arrival of waves at a given location, in principle from everywhere in the simulation volume, i.e. if you set off a "feedback bomb" somewhere, then the timestepping should react to this appropriately before a blast wave from this arrives, thus that the Courant condition is always respected. This is however only guaranteed to work if the external energy injection happens on a globally synchronized timestep. If you apply your custom feedback on arbitrary timesteps, including ones where only a small number of particles are synchronized, then it can indeed happen that other neighbouring particles of the feedback-affected particle have already been assigned to longer, partially completed timesteps, and this can indeed cause problems. One way to avoid this would be to add strong feedback only on global synchronization points, or at least only at times when the neighboring particles of a particle th!
> at receives feedback are also active.
>
>> This is why suggest to move the assignment of time steps in between the first force calculation and the first kick, in order to protect against such rapid changes in these quantities.
> It would be possible to split up the hydro force computation and the kick such that one has also the pressure force available before setting the hydro step. (Note however that the originating particle of the high pressure force should then also be hot, so induce a short timestep via the Courant condition too). This hasn't been done in Gadget4 because of an intended analogy/similarity with Arepo (where the tree-based timestep assignment scheme has first been described), and the way one usually does mesh-based hydro.
>
> Best,
> Volker
>
>
>> Best,
>> Leonard
>>
>> --
>> ===================================================
>> Leonard Romano, B.Sc.(レオナルド・ロマノ)
>> Physics Department
>> Technical University of Munich (TUM), Germany
>> Theoretical Astrophysics Group
>> Department of Earth and Space Science
>> Graduate School of Science, Osaka University, Japan
>> he / him / his
>> ===================================================
>>
>>
>> -----------------------------------------------------------
>>
>> If you wish to unsubscribe from this mailing, send mail to
>> minimalist_at_MPA-Garching.MPG.de with a subject of: unsubscribe gadget-list
>> A web-archive of this mailing list is available here:
>> http://www.mpa-garching.mpg.de/gadget/gadget-list
>
>
>
> -----------------------------------------------------------
>
> If you wish to unsubscribe from this mailing, send mail to
> minimalist_at_MPA-Garching.MPG.de with a subject of: unsubscribe gadget-list
> A web-archive of this mailing list is available here:
> http://www.mpa-garching.mpg.de/gadget/gadget-list

-- 
===================================================
Leonard Romano, B.Sc.(レオナルド・ロマノ)
Physics Department
Technical University of Munich (TUM), Germany
Theoretical Astrophysics Group
Department of Earth and Space Science
Graduate School of Science, Osaka University, Japan
he / him / his
===================================================
Received on 2021-08-11 12:20:26

This archive was generated by hypermail 2.3.0 : 2022-09-01 14:03:43 CEST