diff --git a/Makefile.am b/Makefile.am index c724ebb85360fa6543b66d27adf3bf1bdaffcf57..9623ffc9e203d7d0645aed89400c7183b5b6501a 100644 --- a/Makefile.am +++ b/Makefile.am @@ -80,12 +80,13 @@ endif SUBDIRS = tests # List required headers -include_HEADERS = logger_header.h logger_loader_io.h logger_particle.h logger_time.h logger_tools.h logger_reader.h \ - logger_logfile.h logger_index.h quick_sort.h logger_python_tools.h logger_interpolation.h +include_HEADERS = logger_header.h logger_loader_io.h logger_particle.h logger_time.h logger_tools.h +include_HEADERS += logger_reader.h logger_logfile.h logger_index.h quick_sort.h logger_python_tools.h +include_HEADERS += logger_interpolation.h logger_parameters.h # Common source files AM_SOURCES = logger_header.c logger_loader_io.c logger_particle.c logger_time.c logger_tools.c logger_reader.c -AM_SOURCES += logger_logfile.c logger_index.c quick_sort.c +AM_SOURCES += logger_logfile.c logger_index.c quick_sort.c logger_parameters.c AM_SOURCES += $(GRAVITY_SRC) $(STARS_SRC) $(HYDRO_SRC) $(CHEMISTRY_SRC) $(STAR_FORMATION_SRC) if HAVEPYTHON diff --git a/chemistry/GEAR/logger_chemistry.h b/chemistry/GEAR/logger_chemistry.h index f0fbdd6119fc9689e7dfd301e31f97378f621ea9..635e0d65d5cd4d104f4de2cc052b590cf9583639 100644 --- a/chemistry/GEAR/logger_chemistry.h +++ b/chemistry/GEAR/logger_chemistry.h @@ -25,6 +25,7 @@ #include "chemistry_logger.h" #include "logger_interpolation.h" #include "logger_loader_io.h" +#include "logger_parameters.h" #include "logger_python_tools.h" /* Index of the mask in the header mask array */ @@ -70,12 +71,14 @@ chemistry_logger_reader_link_derivatives_spart(struct header *head) {} * @param t Requested time. * @param field The field to reconstruct (follows the order of * #chemistry_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void chemistry_logger_interpolate_field_part( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t, const int field) { + void *restrict output, const double t, const int field, + const struct logger_parameters *params) { #ifdef SWIFT_DEBUG_CHECKS /* Check the times */ if (t_before > t || t_after < t) { @@ -112,12 +115,14 @@ chemistry_logger_interpolate_field_part( * @param t Requested time. * @param field The field to reconstruct (follows the order of * #chemistry_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void chemistry_logger_interpolate_field_spart( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t, const int field) { + void *restrict output, const double t, const int field, + const struct logger_parameters *params) { #ifdef SWIFT_DEBUG_CHECKS /* Check the times */ @@ -132,7 +137,8 @@ chemistry_logger_interpolate_field_spart( switch (field) { case chemistry_logger_field_spart_metal_mass_fractions: interpolate_linear_double_ND(t_before, before, t_after, after, output, t, - GEAR_CHEMISTRY_ELEMENT_COUNT); + GEAR_CHEMISTRY_ELEMENT_COUNT, + /* periodic= */ 0, params); break; default: diff --git a/chemistry/none/logger_chemistry.h b/chemistry/none/logger_chemistry.h index e776cb352d514e3bbf0a50fdaaed682941b0f18d..8a2a60d161ba77d0ba10ae07f0012aa2507aa294 100644 --- a/chemistry/none/logger_chemistry.h +++ b/chemistry/none/logger_chemistry.h @@ -25,6 +25,7 @@ #include "chemistry_logger.h" #include "logger_interpolation.h" #include "logger_loader_io.h" +#include "logger_parameters.h" #include "logger_python_tools.h" /* Index of the mask in the header mask array */ @@ -70,12 +71,14 @@ chemistry_logger_reader_link_derivatives_spart(struct header *head) {} * @param t Requested time. * @param field The field to reconstruct (follows the order of * #chemistry_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void chemistry_logger_interpolate_field_part( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t, const int field) {} + void *restrict output, const double t, const int field, + const struct logger_parameters *params) {} /** * @brief Interpolate a field of the #spart at the given time. @@ -92,12 +95,14 @@ chemistry_logger_interpolate_field_part( * @param t Requested time. * @param field The field to reconstruct (follows the order of * #chemistry_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void chemistry_logger_interpolate_field_spart( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t, const int field) {} + void *restrict output, const double t, const int field, + const struct logger_parameters *params) {} #ifdef HAVE_PYTHON /** diff --git a/gravity/MultiSoftening/logger_gravity.h b/gravity/MultiSoftening/logger_gravity.h index c5cfcdc07e476f8c04191fe139b446ed281fd107..07154fc24a3a55e33df8b3269a8a13493f120209 100644 --- a/gravity/MultiSoftening/logger_gravity.h +++ b/gravity/MultiSoftening/logger_gravity.h @@ -25,6 +25,7 @@ #include "gravity_logger.h" #include "logger_interpolation.h" #include "logger_loader_io.h" +#include "logger_parameters.h" #include "logger_python_tools.h" /* Index of the mask in the header mask array */ @@ -72,6 +73,7 @@ gravity_logger_reader_link_derivatives(struct header *head) { * @param t Requested time. * @param field The field to reconstruct (follows the order of * #gravity_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void gravity_logger_interpolate_field(const double t_before, @@ -79,7 +81,8 @@ gravity_logger_interpolate_field(const double t_before, const double t_after, const struct logger_field *restrict after, void *restrict output, const double t, - const int field) { + const int field, + const struct logger_parameters *params) { #ifdef SWIFT_DEBUG_CHECKS /* Check the times */ @@ -94,18 +97,21 @@ gravity_logger_interpolate_field(const double t_before, switch (field) { case gravity_logger_field_coordinates: interpolate_quintic_double_float_ND(t_before, before, t_after, after, - output, t, /* dimension= */ 3); + output, t, /* dimension= */ 3, + params->periodic, params); break; case gravity_logger_field_velocities: interpolate_cubic_float_ND(t_before, before, t_after, after, output, t, - /* dimension= */ 3); + /* dimension= */ 3, /* periodic= */ 0, params); break; case gravity_logger_field_accelerations: interpolate_linear_float_ND(t_before, before, t_after, after, output, t, - /* dimension= */ 3); + /* dimension= */ 3, /* periodic= */ 0, + params); break; case gravity_logger_field_masses: - interpolate_linear_float(t_before, before, t_after, after, output, t); + interpolate_linear_float(t_before, before, t_after, after, output, t, + /* periodic= */ 0, params); break; case gravity_logger_field_particle_ids: interpolate_ids(t_before, before, t_after, after, output, t); diff --git a/hydro/Gadget2/logger_hydro.h b/hydro/Gadget2/logger_hydro.h index affa61e8ec7bb6491f99ffdb134a0fa3bd0bf84e..c8f5efe2c4897c78ff81979aab1d7c6ac89ec2e0 100644 --- a/hydro/Gadget2/logger_hydro.h +++ b/hydro/Gadget2/logger_hydro.h @@ -25,6 +25,7 @@ #include "hydro_logger.h" #include "logger_interpolation.h" #include "logger_loader_io.h" +#include "logger_parameters.h" #include "logger_python_tools.h" /* Index of the mask in the header mask array */ @@ -72,6 +73,7 @@ hydro_logger_reader_link_derivatives(struct header *head) { * @param t Requested time. * @param field The field to reconstruct (follows the order of * #hydro_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void hydro_logger_interpolate_field(const double t_before, @@ -79,7 +81,8 @@ hydro_logger_interpolate_field(const double t_before, const double t_after, const struct logger_field *restrict after, void *restrict output, const double t, - const int field) { + const int field, + const struct logger_parameters *params) { #ifdef SWIFT_DEBUG_CHECKS /* Check the times */ @@ -95,23 +98,25 @@ hydro_logger_interpolate_field(const double t_before, /* Do the position */ case hydro_logger_field_coordinates: interpolate_quintic_double_float_ND(t_before, before, t_after, after, - output, t, /* dimesion= */ 3); + output, t, /* dimesion= */ 3, + params->periodic, params); break; /* Do the velocity */ case hydro_logger_field_velocities: interpolate_cubic_float_ND(t_before, before, t_after, after, output, t, - /* dimesion= */ 3); + /* dimesion= */ 3, /* periodic= */ 0, params); break; case hydro_logger_field_accelerations: interpolate_linear_float_ND(t_before, before, t_after, after, output, t, - /* dimesion= */ 3); + /* dimesion= */ 3, /* periodic= */ 0, params); break; /* Do the linear interpolation of float. */ case hydro_logger_field_smoothing_lengths: case hydro_logger_field_entropies: case hydro_logger_field_densities: case hydro_logger_field_masses: - interpolate_linear_float(t_before, before, t_after, after, output, t); + interpolate_linear_float(t_before, before, t_after, after, output, t, + /* periodic= */ 0, params); break; /* Check the ids */ case hydro_logger_field_particle_ids: diff --git a/hydro/SPHENIX/logger_hydro.h b/hydro/SPHENIX/logger_hydro.h index 7f4288f24e269ced02d4d066fde43a1e01155538..bb60efd67e58ab4c54594543d00eb6e14a55496e 100644 --- a/hydro/SPHENIX/logger_hydro.h +++ b/hydro/SPHENIX/logger_hydro.h @@ -26,6 +26,7 @@ #include "hydro_logger.h" #include "logger_interpolation.h" #include "logger_loader_io.h" +#include "logger_parameters.h" #include "logger_python_tools.h" /* Index of the mask in the header mask array */ @@ -73,6 +74,7 @@ hydro_logger_reader_link_derivatives(struct header *head) { * @param t Requested time. * @param field The field to reconstruct (follows the order of * #hydro_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void hydro_logger_interpolate_field(const double t_before, @@ -80,7 +82,8 @@ hydro_logger_interpolate_field(const double t_before, const double t_after, const struct logger_field *restrict after, void *restrict output, const double t, - const int field) { + const int field, + const struct logger_parameters *params) { #ifdef SWIFT_DEBUG_CHECKS /* Check the times */ @@ -96,23 +99,26 @@ hydro_logger_interpolate_field(const double t_before, /* Do the position */ case hydro_logger_field_coordinates: interpolate_quintic_double_float_ND(t_before, before, t_after, after, - output, t, /* dimension= */ 3); + output, t, /* dimension= */ 3, + params->periodic, params); break; /* Do the velocity */ case hydro_logger_field_velocities: interpolate_cubic_float_ND(t_before, before, t_after, after, output, t, - /* dimension= */ 3); + /* dimension= */ 3, /* periodic= */ 0, params); break; case hydro_logger_field_accelerations: interpolate_linear_float_ND(t_before, before, t_after, after, output, t, - /* dimension= */ 3); + /* dimension= */ 3, /* periodic= */ 0, + params); break; /* Do the linear interpolation of float. */ case hydro_logger_field_masses: case hydro_logger_field_smoothing_lengths: case hydro_logger_field_internal_energies: case hydro_logger_field_densities: - interpolate_linear_float(t_before, before, t_after, after, output, t); + interpolate_linear_float(t_before, before, t_after, after, output, t, + /* periodic= */ 0, params); break; /* Check the ids */ case hydro_logger_field_particle_ids: @@ -136,10 +142,12 @@ hydro_logger_interpolate_field(const double t_before, /* Div v */ x[n_linear] = interpolate_cubic_hermite_spline( t_before, bef[n_linear], bef[n_linear + 1], t_after, aft[n_linear], - aft[n_linear + 1], t); + aft[n_linear + 1], t, /* periodic= */ 0, params); /* d Div v / dt */ x[n_linear + 1] = wa * aft[n_linear + 1] + wb * bef[n_linear + 1]; + /* Use the linear interpolation */ + x[1] = wa * div_aft[1] + wb * div_bef[1]; break; } diff --git a/logger_index.c b/logger_index.c index b4bfe7adb5f62b70cc46cf73e07a09e97bdb13d3..606f7cd8c059b313ac71efab51a79c8bc0aaaa45 100644 --- a/logger_index.c +++ b/logger_index.c @@ -208,7 +208,7 @@ struct index_data *logger_index_get_removed_history(struct logger_index *index, */ int logger_index_contains_time_array(struct logger_index *index) { /* Only the first index file should have a time array */ - if (index->time != 0.) { + if (index->integer_time != 0) { error_python("Only the first index file can have a time array."); } diff --git a/logger_interpolation.h b/logger_interpolation.h index 7bc125e1c38fe742d19e7c1a34d9d295d5b793c1..775059aa3086d9ab288425e925331a7c82cb0d7d 100644 --- a/logger_interpolation.h +++ b/logger_interpolation.h @@ -20,6 +20,7 @@ #ifndef LOGGER_LOGGER_INTERPOLATION_H #define LOGGER_LOGGER_INTERPOLATION_H +#include "logger_parameters.h" #include "logger_tools.h" /** @@ -151,6 +152,8 @@ interpolate_cubic_hermite_spline(const double t0, const float v0, * @param t_before Time of field_before (< t). * @param t_after Time of field_after (> t). * @param t Requested time. + * @param periodic Should the periodic boundary be applied? + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void interpolate_quintic_double_float_ND(const double t_before, @@ -158,35 +161,72 @@ interpolate_quintic_double_float_ND(const double t_before, const double t_after, const struct logger_field *restrict after, void *restrict output, const double t, - const int dimension) { - /* Compute the interpolation scaling. */ - const double wa = (t - t_before) / (t_after - t_before); - const double wb = 1. - wa; + const int dimension, int periodic, + const struct logger_parameters *params) { - for (int i = 0; i < dimension; i++) { - double *x = (double *)output; - const double *x_bef = (double *)before->field; - const float *v_bef = (float *)before->first_deriv; - const float *a_bef = (float *)before->second_deriv; + /* Get the arrays */ + const double *x_bef = (double *)before->field; + const float *v_bef = (float *)before->first_deriv; + const float *a_bef = (float *)before->second_deriv; + + const double *x_aft = (double *)after->field; + const float *v_aft = (float *)after->first_deriv; + const float *a_aft = (float *)after->second_deriv; - const double *x_aft = (double *)after->field; - const float *v_aft = (float *)after->first_deriv; - const float *a_aft = (float *)after->second_deriv; + double *x = (double *)output; + + /* Loop over the dimensions */ + for (int i = 0; i < dimension; i++) { + /* Get the current variables */ + double x0 = x_bef[i]; + double x1 = x_aft[i]; + float v0 = v_bef ? v_bef[i] : -1; + float v1 = v_aft ? v_aft[i] : -1; + float a0 = a_bef ? a_bef[i] : -1; + float a1 = a_aft ? a_aft[i] : -1; + + /* Apply the periodic boundaries */ + if (periodic) { + /* Move towards the origin for more accuracy */ + const double half = 0.5 * params->box_size[i]; + if (x0 > x1 + half) + x0 -= params->box_size[i]; + else if (x1 > x0 + half) + x1 -= params->box_size[i]; + } /* Use quintic hermite spline. */ if (v_bef && v_aft && a_bef && a_aft) { - x[i] = interpolate_quintic_hermite_spline(t_before, x_bef[i], v_bef[i], - a_bef[i], t_after, x_aft[i], - v_aft[i], a_aft[i], t); + x[i] = interpolate_quintic_hermite_spline(t_before, x0, v0, a0, t_after, + x1, v1, a1, t); } /* Use cubic hermite spline. */ else if (v_bef && v_aft) { - x[i] = interpolate_cubic_hermite_spline(t_before, x_bef[i], v_bef[i], - t_after, x_aft[i], v_aft[i], t); + x[i] = interpolate_cubic_hermite_spline(t_before, x0, v0, t_after, x1, v1, + t); } /* Use linear interpolation. */ else { - x[i] = wa * x_aft[i] + wb * x_bef[i]; + const double wa = (t - t_before) / (t_after - t_before); + const double wb = 1. - wa; + x[i] = wa * x1 + wb * x0; + } + + /* Apply the periodic boundaries to the output */ + if (periodic) { + /* We suppose that only 1 wrapping is enough + otherwise something really bad happened */ + box_wrap(x[i], 0, params->box_size[i]); + +#ifdef SWIFT_DEBUG_CHECKS + /* Check if the periodic box is correctly applied */ + if (x[i] >= params->box_size[i] || x[i] < 0) { + error_python( + "A particle is outside the periodic box (%g):" + "before=%g, after=%g interpolation=%g", + params->box_size[i], x0, x1, x[i]); + } +#endif } } } @@ -203,11 +243,21 @@ interpolate_quintic_double_float_ND(const double t_before, * @param t_before Time of field_before (< t). * @param t_after Time of field_after (> t). * @param t Requested time. + * @param periodic Should the periodic boundary be applied? + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void interpolate_cubic_float_ND( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t, const int dimension) { + void *restrict output, const double t, const int dimension, int periodic, + const struct logger_parameters *params) { + +#ifdef SWIFT_DEBUG_CHECKS + /* The position is expected to be a double => error with periodic */ + if (periodic) { + error_python("Not implemented yet"); + } +#endif /* Compute the interpolation scaling. */ const float wa = (t - t_before) / (t_after - t_before); @@ -244,11 +294,21 @@ __attribute__((always_inline)) INLINE static void interpolate_cubic_float_ND( * @param t_before Time of field_before (< t). * @param t_after Time of field_after (> t). * @param t Requested time. + * @param periodic Should the periodic boundary be applied? + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void interpolate_linear_float_ND( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t, const int dimension) { + void *restrict output, const double t, const int dimension, int periodic, + const struct logger_parameters *params) { + +#ifdef SWIFT_DEBUG_CHECKS + /* The position is expected to be a double => error with periodic */ + if (periodic) { + error_python("Not implemented yet"); + } +#endif /* Compute the interpolation scaling. */ const float wa = (t - t_before) / (t_after - t_before); @@ -274,11 +334,22 @@ __attribute__((always_inline)) INLINE static void interpolate_linear_float_ND( * @param t_before Time of field_before (< t). * @param t_after Time of field_after (> t). * @param t Requested time. + * @param periodic Should the periodic boundary be applied? + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void interpolate_linear_double_ND( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t, const int dimension) { + void *restrict output, const double t, const int dimension, int periodic, + const struct logger_parameters *params) { + +#ifdef SWIFT_DEBUG_CHECKS + /* The velocity and acceleration are supposed to be provided => quintic => + * error here */ + if (periodic) { + error_python("Not implemented yet"); + } +#endif /* Compute the interpolation scaling. */ const double wa = (t - t_before) / (t_after - t_before); @@ -302,11 +373,21 @@ __attribute__((always_inline)) INLINE static void interpolate_linear_double_ND( * @param t_before Time of field_before (< t). * @param t_after Time of field_after (> t). * @param t Requested time. + * @param periodic Should the periodic boundary be applied? + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void interpolate_linear_float( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t) { + void *restrict output, const double t, int periodic, + const struct logger_parameters *params) { + +#ifdef SWIFT_DEBUG_CHECKS + /* The position is expected to be a double => error with periodic */ + if (periodic) { + error_python("Not implemented yet"); + } +#endif /* Compute the interpolation scaling. */ const float wa = (t - t_before) / (t_after - t_before); @@ -337,4 +418,5 @@ __attribute__((always_inline)) INLINE static void interpolate_ids( } *(long long *)output = *(long long *)after->field; } + #endif // LOGGER_LOGGER_INTERPOLATION_H diff --git a/logger_parameters.c b/logger_parameters.c new file mode 100644 index 0000000000000000000000000000000000000000..cb78d7cea44a750a2dc2cfd43faee51b6bd22905 --- /dev/null +++ b/logger_parameters.c @@ -0,0 +1,70 @@ +/******************************************************************************* + * This file is part of SWIFT. + * Copyright (c) 2019 Loic Hausammann (loic.hausammann@epfl.ch) + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU Lesser General Public License as published + * by the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. + * + ******************************************************************************/ + +/* Include corresponding header */ +#include "logger_parameters.h" + +/* Include the logger */ +#include "logger_reader.h" +#include "logger_tools.h" + +/* Include swift */ +#include "parser.h" + +/** + * @brief Initialize the parameters from the yaml file. + * + * @param reader The #logger_reader. + * @param basename The basename of the logger files. + * @param verbose The verbose level. + */ +void logger_parameters_init(struct logger_parameters *params, + const struct logger_reader *reader) { + + /* Generate filename */ + char filename[STRING_SIZE]; + strcpy(filename, reader->basename); + size_t len = strlen(filename); + /* Remove the rank number */ + filename[len - 5] = '\0'; + strcat(filename, ".yml"); + + /* Initialize the parser */ + struct swift_params swift_params; + parser_read_file(filename, &swift_params); + +#ifdef SWIFT_DEBUG_CHECKS + /* Print the parameters */ + if (reader->verbose > 0) parser_print_params(&swift_params); +#endif + + /* Read the number of dimension */ + params->dimension = parser_get_param_int(&swift_params, "Header:Dimension"); + + /* Read the box size */ + parser_get_param_double_array(&swift_params, "Header:BoxSize", 3, + params->box_size); + + /* Read the periodicity */ + params->periodic = parser_get_param_int(&swift_params, "Header:Periodic"); + + /* Read if we are running with cosmology */ + params->with_cosmology = + parser_get_param_int(&swift_params, "Policy:cosmological integration"); +} diff --git a/logger_parameters.h b/logger_parameters.h new file mode 100644 index 0000000000000000000000000000000000000000..8bf8bd2f55a7f85e570b48fb68f9e8ef1095e702 --- /dev/null +++ b/logger_parameters.h @@ -0,0 +1,51 @@ +/******************************************************************************* + * This file is part of SWIFT. + * Copyright (c) 2019 Loic Hausammann (loic.hausammann@epfl.ch) + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU Lesser General Public License as published + * by the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. + * + ******************************************************************************/ + +/** + * @file logger_parameters.h + * @brief This file contains the functions that deals with the yaml file. + */ + +#ifndef LOGGER_LOGGER_PARAMETERS_H +#define LOGGER_LOGGER_PARAMETERS_H + +struct logger_reader; + +/** + * @brief Structure containing the simulation's parameters. + */ +struct logger_parameters { + + /* Number of dimension */ + int dimension; + + /* Simulation volume */ + double box_size[3]; + + /* Is the volume periodic? */ + int periodic; + + /* Are we doing a cosmological simulation? */ + int with_cosmology; +}; + +void logger_parameters_init(struct logger_parameters *params, + const struct logger_reader *reader); + +#endif // LOGGER_LOGGER_PARAMETERS_H diff --git a/logger_particle.c b/logger_particle.c index 6a5084763ea90320f68a2a3ae8332bef00e0185d..e3b02d61e17313eb9e304e5a8d7dbe815045bde5 100644 --- a/logger_particle.c +++ b/logger_particle.c @@ -147,12 +147,14 @@ logger_particle_read_special_flag(const struct logger_reader *reader, * @param time_after Time of field_after (> t). * @param time Requested time. * @param field The field to reconstruct. + * @param params The simulation's #logger_parameters. */ void logger_particle_interpolate_field( const double time_before, const struct logger_field *restrict before, const double time_after, const struct logger_field *restrict after, void *restrict output, const double time, - const struct field_information *field, enum part_type type) { + const struct field_information *field, enum part_type type, + const struct logger_parameters *params) { /* Select the correct interpolation */ switch (type) { @@ -162,13 +164,14 @@ void logger_particle_interpolate_field( switch (field->module) { case field_module_default: hydro_logger_interpolate_field(time_before, before, time_after, after, - output, time, field->local_index); + output, time, field->local_index, + params); break; case field_module_chemistry: - chemistry_logger_interpolate_field_part(time_before, before, - time_after, after, output, - time, field->local_index); + chemistry_logger_interpolate_field_part( + time_before, before, time_after, after, output, time, + field->local_index, params); break; default: @@ -183,7 +186,8 @@ void logger_particle_interpolate_field( if (field->module != field_module_default) error("Module not implemented"); gravity_logger_interpolate_field(time_before, before, time_after, after, - output, time, field->local_index); + output, time, field->local_index, + params); break; /* Stars */ @@ -191,19 +195,20 @@ void logger_particle_interpolate_field( switch (field->module) { case field_module_default: stars_logger_interpolate_field(time_before, before, time_after, after, - output, time, field->local_index); + output, time, field->local_index, + params); break; case field_module_chemistry: - chemistry_logger_interpolate_field_spart(time_before, before, - time_after, after, output, - time, field->local_index); + chemistry_logger_interpolate_field_spart( + time_before, before, time_after, after, output, time, + field->local_index, params); break; case field_module_star_formation: - star_formation_logger_interpolate_field(time_before, before, - time_after, after, output, - time, field->local_index); + star_formation_logger_interpolate_field( + time_before, before, time_after, after, output, time, + field->local_index, params); break; default: diff --git a/logger_particle.h b/logger_particle.h index ce73e9e05a40550d3d624a94bd07d11acc9682d4..c20aabcf82f4d46a9af0d27082dd783d4723bb5d 100644 --- a/logger_particle.h +++ b/logger_particle.h @@ -27,6 +27,7 @@ #include "logger_gravity.h" #include "logger_header.h" #include "logger_hydro.h" +#include "logger_parameters.h" #include "logger_star_formation.h" #include "logger_stars.h" #include "logger_time.h" @@ -53,7 +54,8 @@ void logger_particle_interpolate_field( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, void *restrict output, const double t, - const struct field_information *field, enum part_type type); + const struct field_information *field, enum part_type type, + const struct logger_parameters *params); enum logger_special_flags logger_particle_read_special_flag( const struct logger_reader *reader, size_t offset, size_t *mask, diff --git a/logger_python_wrapper.c b/logger_python_wrapper.c index 6a69b8b95dfcf0bc3afd6b09f941f41941d9b72b..b696022712ef5428f8719c52c8f9b8887d80cb26 100644 --- a/logger_python_wrapper.c +++ b/logger_python_wrapper.c @@ -181,6 +181,30 @@ static PyObject *getTimeLimits(PyObject *self, PyObject *Py_UNUSED(ignored)) { return (PyObject *)out; } +/** + * @brief Read the box size. + * + * <b>returns</b> tuple containing the box size. + */ +static PyObject *get_box_size(PyObject *self, PyObject *Py_UNUSED(ignored)) { + if (!((PyObjectReader *)self)->ready) { + error_python( + "The logger is not ready yet." + "Did you forget to open it with \"with\"?"); + } + + /* initialize the reader. */ + struct logger_reader *reader = &((PyObjectReader *)self)->reader; + + /* Create the output */ + PyObject *out = PyTuple_New(reader->params.dimension); + for (int i = 0; i < reader->params.dimension; i++) { + PyTuple_SetItem(out, i, PyFloat_FromDouble(reader->params.box_size[i])); + } + + return (PyObject *)out; +} + #define find_field_in_module_internal(MODULE, PART) \ for (int local = 0; local < MODULE##_logger_field##PART##_count; local++) { \ const int global = MODULE##_logger_local_to_global##PART[local]; \ @@ -871,6 +895,12 @@ static PyMethodDef libloggerReaderMethods[] = { "-------\n\n" "times: tuple\n" " time min, time max\n"}, + {"get_box_size", get_box_size, METH_NOARGS, + "Gives the box size of the simulation.\n\n" + "Returns\n" + "-------\n\n" + "box_size: tuple\n" + " The box size.\n"}, {"get_data", (PyCFunction)pyGetData, METH_VARARGS | METH_KEYWORDS, "Read some fields from the logfile at a given time.\n\n" "Parameters\n" diff --git a/logger_reader.c b/logger_reader.c index 8875860878730d57199675ec3b61b5912fe5e083..b7713f7e158799b752a68b4dc39f71f8f231d655 100644 --- a/logger_reader.c +++ b/logger_reader.c @@ -20,6 +20,8 @@ /* Include corresponding header */ #include "logger_reader.h" +#include "logger_parameters.h" + /* Include standard library */ #include <sys/sysinfo.h> #include <unistd.h> @@ -50,6 +52,7 @@ void logger_reader_init(struct logger_reader *reader, const char *basename, /* Initialize the reader variables. */ reader->verbose = verbose; + logger_parameters_init(&reader->params, reader); /* Generate the logfile filename */ char logfile_name[STRING_SIZE]; @@ -539,7 +542,8 @@ int logger_reader_read_field(struct logger_reader *reader, double time, /* Interpolate the data. */ logger_particle_interpolate_field(time_before, &before, time_after, &after, - output, time, field_wanted, type); + output, time, field_wanted, type, + &reader->params); } return 0; } diff --git a/logger_reader.h b/logger_reader.h index d770d8abf41949081559c354302ed391a8969808..8189589d7d47480471a8383163e5c498577300f9 100644 --- a/logger_reader.h +++ b/logger_reader.h @@ -50,6 +50,7 @@ #include "logger_index.h" #include "logger_loader_io.h" #include "logger_logfile.h" +#include "logger_parameters.h" #include "logger_particle.h" /** @@ -101,6 +102,9 @@ struct logger_reader { size_t index; } time; + /* Information from the yaml file */ + struct logger_parameters params; + /* Level of verbosity. */ int verbose; }; diff --git a/star_formation/GEAR/logger_star_formation.h b/star_formation/GEAR/logger_star_formation.h index bcdd80e4dd8d16b39cb771cd31e87a663c1d324d..729d914cec8cb6ef70fa229d26e34b2d359c5785 100644 --- a/star_formation/GEAR/logger_star_formation.h +++ b/star_formation/GEAR/logger_star_formation.h @@ -24,6 +24,7 @@ /* local includes */ #include "logger_interpolation.h" #include "logger_loader_io.h" +#include "logger_parameters.h" #include "logger_python_tools.h" #include "star_formation_particle_logger.h" @@ -58,12 +59,14 @@ star_formation_logger_reader_link_derivatives(struct header *head) {} * @param t Requested time. * @param field The field to reconstruct (follows the order of * #star_formation_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void star_formation_logger_interpolate_field( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t, const int field) { + void *restrict output, const double t, const int field, + const struct logger_parameters *params) { #ifdef SWIFT_DEBUG_CHECKS /* Check the times */ @@ -78,7 +81,7 @@ star_formation_logger_interpolate_field( switch (field) { case star_formation_logger_field_all: interpolate_linear_float_ND(t_before, before, t_after, after, output, t, - /* Dimension */ 2); + /* Dimension */ 2, /* periodic= */ 0, params); /* Deal with the progenitor id */ long long *proj_id = diff --git a/star_formation/none/logger_star_formation.h b/star_formation/none/logger_star_formation.h index ad5721b94b6cf758b76abcb6df36186590f28920..2c1ed94887f4b9f0abc98cf607edfdd867b4d6c1 100644 --- a/star_formation/none/logger_star_formation.h +++ b/star_formation/none/logger_star_formation.h @@ -24,6 +24,7 @@ /* local includes */ #include "logger_interpolation.h" #include "logger_loader_io.h" +#include "logger_parameters.h" #include "logger_python_tools.h" #include "star_formation_particle_logger.h" @@ -58,12 +59,14 @@ star_formation_logger_reader_link_derivatives(struct header *head) {} * @param t Requested time. * @param field The field to reconstruct (follows the order of * #star_formation_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void star_formation_logger_interpolate_field( const double t_before, const struct logger_field *restrict before, const double t_after, const struct logger_field *restrict after, - void *restrict output, const double t, const int field) {} + void *restrict output, const double t, const int field, + const struct logger_parameters *params) {} #ifdef HAVE_PYTHON __attribute__((always_inline)) INLINE static void diff --git a/stars/Basic/logger_stars.h b/stars/Basic/logger_stars.h index 65d6b37dc58b069d1aa44e4be59d7ba34a7489af..a979d74f6406392e4c7e3acc652c2f996e6a3bd3 100644 --- a/stars/Basic/logger_stars.h +++ b/stars/Basic/logger_stars.h @@ -24,6 +24,7 @@ /* local includes */ #include "logger_interpolation.h" #include "logger_loader_io.h" +#include "logger_parameters.h" #include "logger_python_tools.h" #include "stars_logger.h" @@ -72,6 +73,7 @@ stars_logger_reader_link_derivatives(struct header *head) { * @param t Requested time. * @param field The field to reconstruct (follows the order of * #stars_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void stars_logger_interpolate_field(const double t_before, @@ -79,7 +81,8 @@ stars_logger_interpolate_field(const double t_before, const double t_after, const struct logger_field *restrict after, void *restrict output, const double t, - const int field) { + const int field, + const struct logger_parameters *params) { #ifdef SWIFT_DEBUG_CHECKS /* Check the times */ @@ -94,19 +97,22 @@ stars_logger_interpolate_field(const double t_before, switch (field) { case stars_logger_field_coordinates: interpolate_quintic_double_float_ND(t_before, before, t_after, after, - output, t, /* dimension= */ 3); + output, t, /* dimension= */ 3, + params->periodic, params); break; case stars_logger_field_velocities: interpolate_cubic_float_ND(t_before, before, t_after, after, output, t, - /* dimension= */ 3); + /* dimension= */ 3, /* periodic= */ 0, params); break; case stars_logger_field_accelerations: interpolate_linear_float_ND(t_before, before, t_after, after, output, t, - /* dimension= */ 3); + /* dimension= */ 3, /* periodic= */ 0, + params); break; case stars_logger_field_smoothing_lengths: case stars_logger_field_masses: - interpolate_linear_float(t_before, before, t_after, after, output, t); + interpolate_linear_float(t_before, before, t_after, after, output, t, + /* periodic= */ 0, params); break; case stars_logger_field_particle_ids: interpolate_ids(t_before, before, t_after, after, output, t); diff --git a/stars/GEAR/logger_stars.h b/stars/GEAR/logger_stars.h index 2e7f8ad19253fbe435b7a383f31d700c7d471525..70ad1938ac867d3493c78a15228880effde67846 100644 --- a/stars/GEAR/logger_stars.h +++ b/stars/GEAR/logger_stars.h @@ -24,6 +24,7 @@ /* local includes */ #include "logger_interpolation.h" #include "logger_loader_io.h" +#include "logger_parameters.h" #include "logger_python_tools.h" #include "stars_logger.h" @@ -72,6 +73,7 @@ stars_logger_reader_link_derivatives(struct header *head) { * @param t Requested time. * @param field The field to reconstruct (follows the order of * #stars_logger_fields). + * @param params The simulation's #logger_parameters. */ __attribute__((always_inline)) INLINE static void stars_logger_interpolate_field(const double t_before, @@ -79,7 +81,8 @@ stars_logger_interpolate_field(const double t_before, const double t_after, const struct logger_field *restrict after, void *restrict output, const double t, - const int field) { + const int field, + const struct logger_parameters *params) { #ifdef SWIFT_DEBUG_CHECKS /* Check the times */ @@ -94,20 +97,23 @@ stars_logger_interpolate_field(const double t_before, switch (field) { case stars_logger_field_coordinates: interpolate_quintic_double_float_ND(t_before, before, t_after, after, - output, t, /* dimension= */ 3); + output, t, /* dimension= */ 3, + params->periodic, params); break; case stars_logger_field_velocities: interpolate_cubic_float_ND(t_before, before, t_after, after, output, t, - /* dimension= */ 3); + /* dimension= */ 3, /* periodic= */ 0, params); break; case stars_logger_field_accelerations: interpolate_linear_float_ND(t_before, before, t_after, after, output, t, - /* dimension= */ 3); + /* dimension= */ 3, /* periodic= */ 0, + params); break; case stars_logger_field_smoothing_lengths: case stars_logger_field_masses: case stars_logger_field_birth_scale_factors: - interpolate_linear_float(t_before, before, t_after, after, output, t); + interpolate_linear_float(t_before, before, t_after, after, output, t, + /* periodic= */ 0, params); break; case stars_logger_field_particle_ids: interpolate_ids(t_before, before, t_after, after, output, t);