diff --git a/src/cooling/Compton/cooling.h b/src/cooling/Compton/cooling.h index f440cd03455c07d2eeb64c37189aed36efe78e09..d6ec24f02338a44ac4144968fba9b080b979cf17 100644 --- a/src/cooling/Compton/cooling.h +++ b/src/cooling/Compton/cooling.h @@ -273,6 +273,49 @@ __attribute__((always_inline)) INLINE static void cooling_first_init_part( xp->cooling_data.radiated_energy = 0.f; } +/** + * @brief Compute the temperature of a #part based on the cooling function. + * + * @param phys_const #phys_const data structure. + * @param hydro_props The properties of the hydro scheme. + * @param us The internal system of units. + * @param cosmo #cosmology data structure. + * @param cooling #cooling_function_data struct. + * @param p #part data. + * @param xp Pointer to the #xpart data. + */ +INLINE static float cooling_get_temperature( + const struct phys_const* restrict phys_const, + const struct hydro_props* restrict hydro_props, + const struct unit_system* restrict us, + const struct cosmology* restrict cosmo, + const struct cooling_function_data* restrict cooling, + const struct part* restrict p, const struct xpart* restrict xp) { + + /* Physical constants */ + const double m_H = phys_const->const_proton_mass; + const double k_B = phys_const->const_boltzmann_k; + + /* Gas properties */ + const double T_transition = hydro_props->hydrogen_ionization_temperature; + const double mu_neutral = hydro_props->mu_neutral; + const double mu_ionised = hydro_props->mu_ionised; + + /* Particle temperature */ + const double u = hydro_get_physical_internal_energy(p, xp, cosmo); + + /* Temperature over mean molecular weight */ + const double T_over_mu = hydro_gamma_minus_one * u * m_H / k_B; + + /* Are we above or below the HII -> HI transition? */ + if (T_over_mu > (T_transition + 1.) / mu_ionised) + return T_over_mu * mu_ionised; + else if (T_over_mu < (T_transition - 1.) / mu_neutral) + return T_over_mu * mu_neutral; + else + return T_transition; +} + /** * @brief Returns the total radiated energy by this particle. * diff --git a/src/cooling/Compton/cooling_io.h b/src/cooling/Compton/cooling_io.h index d020587c920f781450a5183954bc6c429e461512..8fa3944ff78e7592da3978ee9465451c96e1d533 100644 --- a/src/cooling/Compton/cooling_io.h +++ b/src/cooling/Compton/cooling_io.h @@ -23,6 +23,7 @@ #include "../config.h" /* Local includes */ +#include "cooling.h" #include "io_properties.h" #ifdef HAVE_HDF5 @@ -41,11 +42,18 @@ __attribute__((always_inline)) INLINE static void cooling_write_flavour( } #endif +INLINE static void convert_part_T(const struct engine* e, const struct part* p, + const struct xpart* xp, float* ret) { + + ret[0] = cooling_get_temperature(e->physical_constants, e->hydro_properties, + e->internal_units, e->cosmology, + e->cooling_func, p, xp); +} + /** * @brief Specifies which particle fields to write to a dataset * - * Nothing to write for this scheme. - * + * @param parts The particle array. * @param xparts The extended particle array. * @param list The list of i/o properties to write. * @param cooling The #cooling_function_data @@ -53,10 +61,14 @@ __attribute__((always_inline)) INLINE static void cooling_write_flavour( * @return Returns the number of fields to write. */ __attribute__((always_inline)) INLINE static int cooling_write_particles( - const struct xpart* xparts, struct io_props* list, + const struct part* parts, const struct xpart* xparts, struct io_props* list, const struct cooling_function_data* cooling) { - return 0; + list[0] = io_make_output_field_convert_part("Temperature", FLOAT, 1, + UNIT_CONV_TEMPERATURE, parts, + xparts, convert_part_T); + + return 1; } #endif /* SWIFT_COOLING_IO_COMPTON_H */ diff --git a/src/cooling/const_lambda/cooling_io.h b/src/cooling/const_lambda/cooling_io.h index abd586b61a0355663a961174fc51712e5a5dd291..9437f0f94db41725d6715cf349843bf079137305 100644 --- a/src/cooling/const_lambda/cooling_io.h +++ b/src/cooling/const_lambda/cooling_io.h @@ -61,8 +61,6 @@ INLINE static void convert_part_T(const struct engine* e, const struct part* p, /** * @brief Specifies which particle fields to write to a dataset * - * Nothing to write for this scheme. - * * @param parts The particle array. * @param xparts The extended particle array. * @param list The list of i/o properties to write.