# Re: cooling, heating and temperature oscillation

From: Volker Springel <volker_at_MPA-Garching.MPG.DE>
Date: Fri, 30 Mar 2007 12:13:44 +0200

Hi Yves,

You need to solve the cooling part in the energy equation with an
implicit time-integration scheme to cure this problem.

For example, the new energy (or entropy) at the end of the timestep can
be computed as

u' = u + (du/dt)_visc * dt + (du/dt)'_cool * dt

Here (du/dt)'_cool is the cooling rate evaluated for the *new*
temperature u'. It's a bit cumbersome to solve this equation (will
require some iterative scheme), but then your oscillations should be gone.

Volker

Yves Revaz wrote:
>
> I encounter problems in reproducing the Temperature-Density distribution
> of gas particles in a LCDM model (see for example Fig. 3 of Katz et al
> 96 or Fig.11 of
> Springel & Hernquist 2002).
>
> The overall distribution is correct, however, I fail to reproduce the
> thinness
> of the horizontal branch, corresponding to high density regions
> ((rho/rhom)>1e4, T<1e5K).
> In my simulations, the temperature dispersion of the horizontal branch
> is high,
> with some particles having temperature up to 1e6 K !
>
> This problem comes from the competition between cooling (dA/dt)_rad and
> viscosity heating (dA/dt)_visc, where A is the entropy.
> For a particle with a density (rho/rhom)>1e4 and temperature > 1e4K,
> we have:
> |(dA/dt)_rad| >> |(dA/dt)_visc| => dA/dt)_tot << 0,
>
> the cooling dominates and the temperature quickly decreases. When the
> temperature of
> the particle goes below 1e4K, (dA/dt)_rad drops nearly to zero (cutoff
> in the cooling function),
> and the entropy variation is only due to the (dA/dt)_visc therm, which,
> in some cases is so high
> that the particle temperature instantaneously rises up to 1e6K !!!
> In summary, in the horizontal branch, instead of being more or less
> constant at 1e4K (equilibrium between
> viscosity heating and radiative cooling), the temperature of the
> particles oscillate between 1e4 and 1e5-1e6K.
>
> This behavior is the result of the cooling and heating time scale, much
> shorter than
> the time-step imposed by the currant condition. The cooling is limited
> by the condition
> that :
> dA/dt > -0.5 A.
>
> Imposing also
>
> dA/dt < A,
>
> in not sufficient to damp the temperature oscillation.
> There is probably a well known solution to this problem,
> but I haven't found it in the literature.
>
> Does anyone has a solution ?
>
>
>
> Yves
>
>
> MY PARAMETERS
> ---------------
>
> The simulation test contains 2*64^3 particles in 20 Mpc^3 h^-3
>
> I use the following parameters :
>
> ErrTolIntAccuracy 0.025
> CourantFac 0.15
> MaxSizeTimestep 0.03
> MinSizeTimestep 0
>
> ErrTolTheta 0.8
> TypeOfOpeningCriterion 0
> ErrTolForceAcc 0.005
>
> DesNumNgb 32
> MaxNumNgbDeviation 2
> ArtBulkViscConst 0.8
>
> MinGasHsmlFractional 0.25
> SofteningGas 8.
> SofteningHalo 8.
>