diff --git a/src/hydro/Gizmo/hydro.h b/src/hydro/Gizmo/hydro.h index ccbeff6743de5a6097e9361137711eb97dbb040d..6273f5009a051e17d2a2954ac994c3f189110464 100644 --- a/src/hydro/Gizmo/hydro.h +++ b/src/hydro/Gizmo/hydro.h @@ -24,6 +24,7 @@ #include "approx_math.h" #include "equation_of_state.h" #include "hydro_gradients.h" +#include "hydro_properties.h" #include "hydro_space.h" #include "hydro_unphysical.h" #include "hydro_velocities.h" @@ -130,8 +131,9 @@ __attribute__((always_inline)) INLINE static void hydro_first_init_part( /* remember that we store the total thermal energy, not the specific thermal energy (as in Gadget) */ #if defined(EOS_ISOTHERMAL_GAS) - /* this overwrites the internal energy from the initial condition file */ - p->conserved.energy = mass * eos.isothermal_internal_energy; + /* this overwrites the internal energy from the initial condition file + * Note that we call the EoS function just to get the constant u here. */ + p->conserved.energy = mass * gas_internal_energy_from_entropy(0.f, 0.f); #else p->conserved.energy *= mass; #endif @@ -608,7 +610,9 @@ __attribute__((always_inline)) INLINE static void hydro_kick_extra( p->conserved.momentum[1] += p->conserved.flux.momentum[1]; p->conserved.momentum[2] += p->conserved.flux.momentum[2]; #if defined(EOS_ISOTHERMAL_GAS) - p->conserved.energy = p->conserved.mass * eos.isothermal_internal_energy; + /* We use the EoS equation in a sneaky way here just to get the constant u */ + p->conserved.energy = + p->conserved.mass * gas_internal_energy_from_entropy(0.f, 0.f); #else p->conserved.energy += p->conserved.flux.energy; #endif @@ -707,15 +711,11 @@ __attribute__((always_inline)) INLINE static void hydro_kick_extra( __attribute__((always_inline)) INLINE static float hydro_get_internal_energy( const struct part* restrict p) { - if (p->primitives.rho > 0.) { -#ifdef EOS_ISOTHERMAL_GAS - return eos.isothermal_internal_energy; -#else - return p->primitives.P / hydro_gamma_minus_one / p->primitives.rho; -#endif - } else { + if (p->primitives.rho > 0.) + return gas_internal_energy_from_pressure(p->primitives.rho, + p->primitives.P); + else return 0.; - } } /** @@ -727,7 +727,7 @@ __attribute__((always_inline)) INLINE static float hydro_get_entropy( const struct part* restrict p) { if (p->primitives.rho > 0.) { - return p->primitives.P / pow_gamma(p->primitives.rho); + return gas_entropy_from_pressure(p->primitives.rho, p->primitives.P); } else { return 0.; } @@ -741,16 +741,10 @@ __attribute__((always_inline)) INLINE static float hydro_get_entropy( __attribute__((always_inline)) INLINE static float hydro_get_soundspeed( const struct part* restrict p) { - if (p->primitives.rho > 0.) { -#ifdef EOS_ISOTHERMAL_GAS - return sqrtf(eos.isothermal_internal_energy * hydro_gamma * - hydro_gamma_minus_one); -#else - return sqrtf(hydro_gamma * p->primitives.P / p->primitives.rho); -#endif - } else { + if (p->primitives.rho > 0.) + return gas_soundspeed_from_pressure(p->primitives.rho, p->primitives.P); + else return 0.; - } } /** diff --git a/src/hydro/Gizmo/hydro_io.h b/src/hydro/Gizmo/hydro_io.h index 09d667d8ff22136c7472045d3b6188ab5ec557cc..ea0a083f37078859b236ca0799268066f6fa0385 100644 --- a/src/hydro/Gizmo/hydro_io.h +++ b/src/hydro/Gizmo/hydro_io.h @@ -18,7 +18,7 @@ ******************************************************************************/ #include "adiabatic_index.h" -#include "equation_of_state.h" +#include "hydro.h" #include "hydro_flux_limiters.h" #include "hydro_gradients.h" #include "hydro_slope_limiters.h" @@ -72,15 +72,7 @@ void hydro_read_particles(struct part* parts, struct io_props* list, * @return Internal energy of the particle */ float convert_u(struct engine* e, struct part* p) { - if (p->primitives.rho > 0.) { -#ifdef EOS_ISOTHERMAL_GAS - return eos.isothermal_internal_energy; -#else - return p->primitives.P / hydro_gamma_minus_one / p->primitives.rho; -#endif - } else { - return 0.; - } + return hydro_get_internal_energy(p); } /** @@ -91,11 +83,7 @@ float convert_u(struct engine* e, struct part* p) { * @return Entropic function of the particle */ float convert_A(struct engine* e, struct part* p) { - if (p->primitives.rho > 0.) { - return p->primitives.P / pow_gamma(p->primitives.rho); - } else { - return 0.; - } + return hydro_get_entropy(p); } /**