Skip to content
Snippets Groups Projects
Commit a9df8822 authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Clarify the physical nature of the minimal energy.

parent c5130bee
No related branches found
No related tags found
1 merge request!641Fixes to the cosmo hydro
...@@ -656,12 +656,13 @@ __attribute__((always_inline)) INLINE static void hydro_kick_extra( ...@@ -656,12 +656,13 @@ __attribute__((always_inline)) INLINE static void hydro_kick_extra(
xp->entropy_full += p->entropy_dt * dt_therm; xp->entropy_full += p->entropy_dt * dt_therm;
/* Apply the minimal energy limit */ /* Apply the minimal energy limit */
const float density = p->rho * cosmo->a3_inv; const float physical_density = p->rho * cosmo->a3_inv;
const float min_energy = hydro_props->minimal_internal_energy; const float min_physical_energy = hydro_props->minimal_internal_energy;
const float min_entropy = const float min_physical_entropy =
gas_entropy_from_internal_energy(density, min_energy); gas_entropy_from_internal_energy(physical_density, min_physical_energy);
if (xp->entropy_full < min_entropy) { const float min_comoving_entropy = min_physical_entropy; /* A' = A */
xp->entropy_full = min_entropy; if (xp->entropy_full < min_comoving_entropy) {
xp->entropy_full = min_comoving_entropy;
p->entropy_dt = 0.f; p->entropy_dt = 0.f;
} }
...@@ -693,20 +694,21 @@ __attribute__((always_inline)) INLINE static void hydro_convert_quantities( ...@@ -693,20 +694,21 @@ __attribute__((always_inline)) INLINE static void hydro_convert_quantities(
struct part *restrict p, struct xpart *restrict xp, struct part *restrict p, struct xpart *restrict xp,
const struct cosmology *cosmo, const struct hydro_props *hydro_props) { const struct cosmology *cosmo, const struct hydro_props *hydro_props) {
/* We read u in the entropy field. We now get (comoving) S from (physical) u /* We read u in the entropy field. We now get (comoving) A from (physical) u
* and (physical) rho. Note that comoving S == physical S */ * and (physical) rho. Note that comoving A (A') == physical A */
xp->entropy_full = xp->entropy_full =
gas_entropy_from_internal_energy(p->rho * cosmo->a3_inv, p->entropy); gas_entropy_from_internal_energy(p->rho * cosmo->a3_inv, p->entropy);
p->entropy = xp->entropy_full; p->entropy = xp->entropy_full;
/* Apply the minimal energy limit */ /* Apply the minimal energy limit */
const float density = p->rho * cosmo->a3_inv; const float physical_density = p->rho * cosmo->a3_inv;
const float min_energy = hydro_props->minimal_internal_energy; const float min_physical_energy = hydro_props->minimal_internal_energy;
const float min_entropy = const float min_physical_entropy =
gas_entropy_from_internal_energy(density, min_energy); gas_entropy_from_internal_energy(physical_density, min_physical_energy);
if (xp->entropy_full < min_entropy) { const float min_comoving_entropy = min_physical_entropy; /* A' = A */
xp->entropy_full = min_entropy; if (xp->entropy_full < min_comoving_entropy) {
p->entropy = min_entropy; xp->entropy_full = min_comoving_entropy;
p->entropy = min_comoving_entropy;
p->entropy_dt = 0.f; p->entropy_dt = 0.f;
} }
......
...@@ -662,10 +662,10 @@ __attribute__((always_inline)) INLINE static void hydro_kick_extra( ...@@ -662,10 +662,10 @@ __attribute__((always_inline)) INLINE static void hydro_kick_extra(
xp->u_full += p->u_dt * dt_therm; xp->u_full += p->u_dt * dt_therm;
/* Apply the minimal energy limit */ /* Apply the minimal energy limit */
const float min_energy = const float min_comoving_energy =
hydro_props->minimal_internal_energy / cosmo->a_factor_internal_energy; hydro_props->minimal_internal_energy / cosmo->a_factor_internal_energy;
if (xp->u_full < min_energy) { if (xp->u_full < min_comoving_energy) {
xp->u_full = min_energy; xp->u_full = min_comoving_energy;
p->u_dt = 0.f; p->u_dt = 0.f;
} }
...@@ -704,11 +704,11 @@ __attribute__((always_inline)) INLINE static void hydro_convert_quantities( ...@@ -704,11 +704,11 @@ __attribute__((always_inline)) INLINE static void hydro_convert_quantities(
xp->u_full = p->u; xp->u_full = p->u;
/* Apply the minimal energy limit */ /* Apply the minimal energy limit */
const float min_energy = const float min_comoving_energy =
hydro_props->minimal_internal_energy / cosmo->a_factor_internal_energy; hydro_props->minimal_internal_energy / cosmo->a_factor_internal_energy;
if (xp->u_full < min_energy) { if (xp->u_full < min_comoving_energy) {
xp->u_full = min_energy; xp->u_full = min_comoving_energy;
p->u = min_energy; p->u = min_comoving_energy;
p->u_dt = 0.f; p->u_dt = 0.f;
} }
......
...@@ -69,13 +69,13 @@ struct hydro_props { ...@@ -69,13 +69,13 @@ struct hydro_props {
/*! Minimal temperature allowed */ /*! Minimal temperature allowed */
float minimal_temperature; float minimal_temperature;
/*! Minimal internal energy per unit mass */ /*! Minimal physical internal energy per unit mass */
float minimal_internal_energy; float minimal_internal_energy;
/*! Initial temperature */ /*! Initial temperature */
float initial_temperature; float initial_temperature;
/*! Initial internal energy per unit mass */ /*! Initial physical internal energy per unit mass */
float initial_internal_energy; float initial_internal_energy;
/*! Primordial hydrogen mass fraction for initial energy conversion */ /*! Primordial hydrogen mass fraction for initial energy conversion */
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment