Commit 33a58b65 authored by Matthieu Schaller's avatar Matthieu Schaller

Merge branch 'time_base_in_snap_header' into 'master'

Time base in snap header

Closes #704

See merge request !1161
parents 3e7143c1 279470f5
......@@ -33,6 +33,15 @@ The ``RunName`` field contains the name of the simulation that was specified as
the ``run_name`` in the :ref:`Parameters_meta_data` section of the YAML
parameter file.
The ``TimeBase_dloga`` field contains the change in logarithm of the
scale-factor corresponding to a time-step of length 1 on the integer
time-line. This is the smallest time-step size that the code can use. This field
is zero in non-cosmological runs. Similarly, the field ``TimeBase_dt`` contains
the smallest time-step size (in internal units) that the code can take. This
would be the increase in time a particle in the time-bin one would have. Note
that in cosmological runs this quantity evolves with redhsift as the (logarithm
of the) scale-factor is used on the integer time-line.
Meta-data about the code and run
--------------------------------
......
......@@ -853,6 +853,54 @@ double cosmology_get_delta_time_from_scale_factors(const struct cosmology *c,
return t2 - t1;
}
/**
* @brief Compute the time corresponding to the timebase interval at the current
* redshift
*
* This function is slow as it performs the actual integral.
*
* @param c The comology model.
* @param ti_current The current point on the time-line.
*
* @return The time corresponding to c->time_base in internal time units.
*/
double cosmology_get_timebase(struct cosmology *c,
const integertime_t ti_current) {
#ifdef HAVE_LIBGSL
/* We are going to average over a number of time-bins
* as in some cases the smallest bin is too small for double accuracy */
const int num_bins = 24;
/* Range in log(a) */
const double log_a_start = c->log_a_begin + (ti_current)*c->time_base;
const double log_a_end =
c->log_a_begin + (ti_current + (1LL << num_bins)) * c->time_base;
/* Range in (a) */
const double a_start = exp(log_a_start);
const double a_end = exp(log_a_end);
/* Initalise the GSL workspace and function */
gsl_integration_workspace *space =
gsl_integration_workspace_alloc(GSL_workspace_size);
gsl_function F = {&time_integrand, c};
/* Perform the integral */
double result, abserr;
gsl_integration_qag(&F, a_start, a_end, 0, 1.0e-10, GSL_workspace_size,
GSL_INTEG_GAUSS61, space, &result, &abserr);
result /= (1LL << num_bins);
return result;
#else
error("Code not compiled with GSL. Can't compute cosmology integrals.");
#endif
}
/**
* @brief Compute scale factor from time since big bang (in internal units).
*
......
......@@ -213,6 +213,9 @@ double cosmology_get_delta_time_from_scale_factors(const struct cosmology *c,
const double a_start,
const double a_end);
double cosmology_get_timebase(struct cosmology *c,
const integertime_t ti_current);
double cosmology_get_scale_factor(const struct cosmology *cosmo, double t);
double cosmology_get_time_since_big_bang(const struct cosmology *c, double a);
......
......@@ -394,6 +394,16 @@ void write_output_distributed(struct engine* e,
io_write_attribute_s(h_grp, "Code", "SWIFT");
io_write_attribute_s(h_grp, "RunName", e->run_name);
/* Write out the time-base */
if (with_cosmology) {
io_write_attribute_d(h_grp, "TimeBase_dloga", e->time_base);
const double delta_t = cosmology_get_timebase(e->cosmology, e->ti_current);
io_write_attribute_d(h_grp, "TimeBase_dt", delta_t);
} else {
io_write_attribute_d(h_grp, "TimeBase_dloga", 0);
io_write_attribute_d(h_grp, "TimeBase_dt", e->time_base);
}
/* Store the time at which the snapshot was written */
time_t tm = time(NULL);
struct tm* timeinfo = localtime(&tm);
......
......@@ -1175,6 +1175,16 @@ void prepare_file(struct engine* e, const char* fileName,
io_write_attribute_s(h_grp, "Code", "SWIFT");
io_write_attribute_s(h_grp, "RunName", e->run_name);
/* Write out the time-base */
if (with_cosmology) {
io_write_attribute_d(h_grp, "TimeBase_dloga", e->time_base);
const double delta_t = cosmology_get_timebase(e->cosmology, e->ti_current);
io_write_attribute_d(h_grp, "TimeBase_dt", delta_t);
} else {
io_write_attribute_d(h_grp, "TimeBase_dloga", 0);
io_write_attribute_d(h_grp, "TimeBase_dt", e->time_base);
}
/* Store the time at which the snapshot was written */
time_t tm = time(NULL);
struct tm* timeinfo = localtime(&tm);
......
......@@ -1058,6 +1058,17 @@ void write_output_serial(struct engine* e,
io_write_attribute_s(h_grp, "Code", "SWIFT");
io_write_attribute_s(h_grp, "RunName", e->run_name);
/* Write out the time-base */
if (with_cosmology) {
io_write_attribute_d(h_grp, "TimeBase_dloga", e->time_base);
const double delta_t =
cosmology_get_timebase(e->cosmology, e->ti_current);
io_write_attribute_d(h_grp, "TimeBase_dt", delta_t);
} else {
io_write_attribute_d(h_grp, "TimeBase_dloga", 0);
io_write_attribute_d(h_grp, "TimeBase_dt", e->time_base);
}
/* Store the time at which the snapshot was written */
time_t tm = time(NULL);
struct tm* timeinfo = localtime(&tm);
......
......@@ -893,6 +893,16 @@ void write_output_single(struct engine* e,
io_write_attribute_s(h_grp, "Code", "SWIFT");
io_write_attribute_s(h_grp, "RunName", e->run_name);
/* Write out the time-base */
if (with_cosmology) {
io_write_attribute_d(h_grp, "TimeBase_dloga", e->time_base);
const double delta_t = cosmology_get_timebase(e->cosmology, e->ti_current);
io_write_attribute_d(h_grp, "TimeBase_dt", delta_t);
} else {
io_write_attribute_d(h_grp, "TimeBase_dloga", 0);
io_write_attribute_d(h_grp, "TimeBase_dt", e->time_base);
}
/* Store the time at which the snapshot was written */
time_t tm = time(NULL);
struct tm* timeinfo = localtime(&tm);
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment