diff --git a/configure.ac b/configure.ac index 9d0797d8b2b2d0be5a6927836ef4fbbf8c5f1b2d..bfd611d5683af93d1ba0eeb1af28c61119d7e668 100644 --- a/configure.ac +++ b/configure.ac @@ -533,7 +533,7 @@ if test "x$with_grackle" != "xno"; then [grackle], [initialize_chemistry_data], [AC_DEFINE([HAVE_GRACKLE],1,[The GRACKLE library appears to be present.]) - AC_DEFINE([CONFIG_BFLOAT_4],1,[Use double in grackle]) + AC_DEFINE([CONFIG_BFLOAT_8],1,[Use double in grackle]) ], [AC_MSG_ERROR(Cannot find grackle library!)], [$GRACKLE_LIBS $GRACKLE_INCS $FCLIBS] diff --git a/src/cooling/grackle/cooling.h b/src/cooling/grackle/cooling.h index ea4bf86ea3e7262544f42f7ba5c2911818c32000..ba67d35e3ea5acb7d08960b196152b5a4e6006e3 100644 --- a/src/cooling/grackle/cooling.h +++ b/src/cooling/grackle/cooling.h @@ -267,6 +267,180 @@ __attribute__((always_inline)) INLINE static void cooling_print_backend( #endif } +/** + * @brief allocate required field + * + * @param data the #grackle_field_data + */ +__attribute__((always_inline)) INLINE static void cooling_malloc_data(grackle_field_data *data) { + +#if COOLING_GRACKLE_MODE >= 1 + /* primordial chemistry >= 1 */ + data->HI_density = malloc(sizeof(gr_float)); + data->HII_density = malloc(sizeof(gr_float)); + data->HeI_density = malloc(sizeof(gr_float)); + data->HeII_density = malloc(sizeof(gr_float)); + data->HeIII_density = malloc(sizeof(gr_float)); + data->e_density = malloc(sizeof(gr_float)); +#endif // MODE >= 1 + +#if COOLING_GRACKLE_MODE >= 2 + /* primordial chemistry >= 2 */ + data->HM_density = malloc(sizeof(gr_float)); + data->H2I_density = malloc(sizeof(gr_float)); + data->H2II_density = malloc(sizeof(gr_float)); +#endif // MODE 2 + +#if COOLING_GRACKLE_MODE >= 3 + /* primordial chemistry >= 3 */ + data->DI_density = malloc(sizeof(gr_float)); + data->DII_density = malloc(sizeof(gr_float)); + data->HDI_density = malloc(sizeof(gr_float)); +#endif // MODE >= 3 + + /* metal cooling = 1 */ + data->metal_density = malloc(sizeof(gr_float)); + + /* /\* volumetric heating rate *\/ */ + /* data->volumetric_heating_rate = NULL; */ + + /* /\* specific heating rate *\/ */ + /* data->specific_heating_rate = NULL; */ + +} + +/** + * @brief free the allocated memory + * + * @param data the #grackle_field_data + */ + +__attribute__((always_inline)) INLINE static void cooling_free_data(grackle_field_data *data) { + +#if COOLING_GRACKLE_MODE >= 1 + /* primordial chemistry >= 1 */ + free(data->HI_density); + free(data->HII_density); + free(data->HeI_density); + free(data->HeII_density); + free(data->HeIII_density); + free(data->e_density); +#endif // MODE >= 1 + +#if COOLING_GRACKLE_MODE >= 2 + /* primordial chemistry >= 2 */ + free(data->HM_density); + free(data->H2I_density); + free(data->H2II_density); +#endif // MODE 2 + +#if COOLING_GRACKLE_MODE >= 3 + /* primordial chemistry >= 3 */ + free(data->DI_density); + free(data->DII_density); + free(data->HDI_density); +#endif // MODE >= 3 + + /* metal cooling = 1 */ + free(data->metal_density); + + /* /\* volumetric heating rate *\/ */ + /* data->volumetric_heating_rate = NULL; */ + + /* /\* specific heating rate *\/ */ + /* data->specific_heating_rate = NULL; */ + +} + +/** + * @brief copy xp to data + * + * requires the particle to have been transformed into density + * + * @param data the #grackle_field_data + * @param xp the #xpart + */ +__attribute__((always_inline)) INLINE static void cooling_copy_to_data(grackle_field_data *data, const struct xpart *xp) { + +#if COOLING_GRACKLE_MODE >= 1 + /* primordial chemistry >= 1 */ + data->HI_density[0] = xp->cooling_data.HI_frac; + data->HII_density[0] = xp->cooling_data.HII_frac; + data->HeI_density[0] = xp->cooling_data.HeI_frac; + data->HeII_density[0] = xp->cooling_data.HeII_frac; + data->HeIII_density[0] = xp->cooling_data.HeIII_frac; + data->e_density[0] = xp->cooling_data.e_frac; +#endif // MODE >= 1 + +#if COOLING_GRACKLE_MODE >= 2 + /* primordial chemistry >= 2 */ + data->HM_density[0] = xp->cooling_data.HM_frac; + data->H2I_density[0] = xp->cooling_data.H2I_frac; + data->H2II_density[0] = xp->cooling_data.H2II_frac; +#endif // MODE 2 + +#if COOLING_GRACKLE_MODE >= 3 + /* primordial chemistry >= 3 */ + data->DI_density[0] = xp->cooling_data.DI_frac; + data->DII_density[0] = xp->cooling_data.DII_frac; + data->HDI_density[0] = xp->cooling_data.HDI_frac; +#endif // MODE >= 3 + + /* metal cooling = 1 */ + data->metal_density[0] = xp->cooling_data.metal_frac; + + /* volumetric heating rate */ + data->volumetric_heating_rate = NULL; + + /* specific heating rate */ + data->specific_heating_rate = NULL; +} + + +/** + * @brief copy data to xp + * + * @param data the #grackle_field_data + * @param xp the #xpart + */ +__attribute__((always_inline)) INLINE static void cooling_copy_to_particle(const grackle_field_data *data, struct xpart *xp) { + +#if COOLING_GRACKLE_MODE >= 1 + /* primordial chemistry >= 1 */ + xp->cooling_data.HI_frac = data->HI_density[0]; + xp->cooling_data.HII_frac = data->HII_density[0]; + xp->cooling_data.HeI_frac = data->HeI_density[0]; + xp->cooling_data.HeII_frac = data->HeII_density[0]; + xp->cooling_data.HeIII_frac = data->HeIII_density[0]; + xp->cooling_data.e_frac = data->e_density[0]; +#endif // MODE >= 1 + +#if COOLING_GRACKLE_MODE >= 2 + /* primordial chemistry >= 2 */ + xp->cooling_data.HM_frac = data->HM_density[0]; + xp->cooling_data.H2I_frac = data->H2I_density[0]; + xp->cooling_data.H2II_frac = data->H2II_density[0]; +#endif // MODE 2 + +#if COOLING_GRACKLE_MODE >= 3 + /* primordial chemistry >= 3 */ + xp->cooling_data.DI_frac = data->DI_density[0]; + xp->cooling_data.DII_frac = data->DII_density[0]; + xp->cooling_data.HDI_frac = data->HDI_density[0]; +#endif // MODE >= 3 + + /* metal cooling = 1 */ + xp->cooling_data.metal_frac = data->metal_density[0]; + + /* /\* volumetric heating rate *\/ */ + /* data->volumetric_heating_rate = NULL; */ + + /* /\* specific heating rate *\/ */ + /* data->specific_heating_rate = NULL; */ + +} + + /** * @brief Compute the cooling rate and update the particle chemistry data * @@ -327,48 +501,26 @@ __attribute__((always_inline)) INLINE static double cooling_rate( /* transform gas fraction to densities */ cooling_compute_density(xp, *data.density); -#if COOLING_GRACKLE_MODE >= 1 - /* primordial chemistry >= 1 */ - data.HI_density = &xp->cooling_data.HI_frac; - data.HII_density = &xp->cooling_data.HII_frac; - data.HeI_density = &xp->cooling_data.HeI_frac; - data.HeII_density = &xp->cooling_data.HeII_frac; - data.HeIII_density = &xp->cooling_data.HeIII_frac; - data.e_density = &xp->cooling_data.e_frac; -#endif // MODE >= 1 - -#if COOLING_GRACKLE_MODE >= 2 - /* primordial chemistry >= 2 */ - data.HM_density = &xp->cooling_data.HM_frac; - data.H2I_density = &xp->cooling_data.H2I_frac; - data.H2II_density = &xp->cooling_data.H2II_frac; -#endif // MODE 2 - -#if COOLING_GRACKLE_MODE >= 3 - /* primordial chemistry >= 3 */ - data.DI_density = &xp->cooling_data.DI_frac; - data.DII_density = &xp->cooling_data.DII_frac; - data.HDI_density = &xp->cooling_data.HDI_frac; -#endif // MODE >= 3 - + /* allocate grackle data */ + cooling_malloc_data(&data); - /* metal cooling = 1 */ - data.metal_density = &xp->cooling_data.metal_frac; - - /* volumetric heating rate */ - data.volumetric_heating_rate = NULL; - - /* specific heating rate */ - data.specific_heating_rate = NULL; + /* copy data from particle to grackle data */ + cooling_copy_to_data(&data, xp); /* solve chemistry with table */ if (solve_chemistry(&units, &data, dt) == 0) { error("Error in solve_chemistry."); } + /* copy from grackle data to particle */ + cooling_copy_to_particle(&data, xp); + /* transform densities to gas fraction */ cooling_compute_fraction(xp, *data.density); + /* free allocated memory */ + cooling_free_data(&data); + /* compute rate */ return (energy - energy_before) / dt; } diff --git a/src/cooling/grackle/cooling_struct.h b/src/cooling/grackle/cooling_struct.h index 362c5a04aea75f1af740baca0cca606c6912c932..3d6feeffc83865510e569d14511e854db726c777 100644 --- a/src/cooling/grackle/cooling_struct.h +++ b/src/cooling/grackle/cooling_struct.h @@ -61,24 +61,24 @@ struct cooling_xpart_data { /* here all fractions are mass fraction */ #if COOLING_GRACKLE_MODE >= 1 /* primordial chemistry >= 1 */ - gr_float HI_frac; - gr_float HII_frac; - gr_float HeI_frac; - gr_float HeII_frac; - gr_float HeIII_frac; - gr_float e_frac; + float HI_frac; + float HII_frac; + float HeI_frac; + float HeII_frac; + float HeIII_frac; + float e_frac; #if COOLING_GRACKLE_MODE >= 2 /* primordial chemistry >= 2 */ - gr_float HM_frac; - gr_float H2I_frac; - gr_float H2II_frac; + float HM_frac; + float H2I_frac; + float H2II_frac; #if COOLING_GRACKLE_MODE >= 3 /* primordial chemistry >= 3 */ - gr_float DI_frac; - gr_float DII_frac; - gr_float HDI_frac; + float DI_frac; + float DII_frac; + float HDI_frac; #endif // MODE >= 3 #endif // MODE >= 2 @@ -86,7 +86,7 @@ struct cooling_xpart_data { #endif // MODE >= 1 /* metal cooling = 1 */ - gr_float metal_frac; + float metal_frac; }; #endif /* SWIFT_COOLING_STRUCT_NONE_H */