diff --git a/src/common_io.c b/src/common_io.c index 1c6e9434d0091ead5cb6072a72a57fd1688d484f..8066a305fdcddfaea41d5673b46d479937e2b264 100644 --- a/src/common_io.c +++ b/src/common_io.c @@ -968,7 +968,8 @@ void io_write_output_field_parameter(const char* filename) { } -void time_array_read_file(struct time_array *times, const char* filename) { +void time_array_read_file(struct time_array *times, const char* filename, + struct cosmology *cosmo) { /* initialize times */ times->size = 0; @@ -976,10 +977,35 @@ void time_array_read_file(struct time_array *times, const char* filename) { FILE* file = fopen(filename, "r"); if (file == NULL) error("Error opening file '%s'", filename); - /* Read file */ + /* Read header */ ssize_t read; size_t len = 0; char *line = NULL; + if((read = getline(&line, &len, file)) == -1) + error("Unable to read header in file '%s'", filename); + + /* Remove end of line character */ + line[strcspn(line, "\n")] = 0; + + /* Find type of data in file */ + int type; + if (!strcmp(line, "# Redshift")) + type = TIME_ARRAY_REDSHIFT; + else if (!strcmp(line, "# Age")) + type = TIME_ARRAY_AGE; + else if (!strcmp(line, "# Scale Factor")) + type = TIME_ARRAY_SCALE_FACTOR; + else + error("Unable to interpret the header (%s) in file '%s'", + line, filename); + + if (!cosmo && + (type == TIME_ARRAY_SCALE_FACTOR || type == TIME_ARRAY_REDSHIFT)) + error("Unable to compute a redshift or a scale factor without cosmology. " + "Please change the header in '%s'", filename); + + + /* Read file */ while ((read = getline(&line, &len, file)) != -1) { /* Check data size */ @@ -988,19 +1014,29 @@ void time_array_read_file(struct time_array *times, const char* filename) { "Please decrease the number of time required in '%s'.", filename); + double *time = ×->times[times->size]; /* Write data to time array */ - if (sscanf(line, "%lf", ×->times[times->size]) != 1) { + if (sscanf(line, "%lf", time) != 1) { error( "Tried parsing double but found '%s' with illegal double " "characters in file '%s'.", line, filename); } + /* Transform input into correct time (e.g. ages or scale factor) */ + if (type == TIME_ARRAY_REDSHIFT) + *time = cosmology_get_a_from_z(cosmo, *time); + + if (cosmo && type == TIME_ARRAY_AGE) + error("Not implemented"); + /* Update size */ times->size += 1; } time_array_print(times); + + fclose(file); } void time_array_print(const struct time_array *times) { diff --git a/src/common_io.h b/src/common_io.h index 039a794e866b86f7a4f143c4ab0f021f415f3f42..cb6308c46d37c3da5f59a4fdcf5d8ead258d1e3e 100644 --- a/src/common_io.h +++ b/src/common_io.h @@ -33,6 +33,7 @@ #define TIME_ARRAY_MAX_SIZE 8192 /* Avoid cyclic inclusion problems */ +struct cosmology; struct part; struct gpart; struct spart; @@ -63,6 +64,12 @@ enum IO_DATA_TYPE { */ enum IO_STF_OUTPUT_FORMAT { STEPS = 0, TIME }; +enum TIME_ARRAY_TYPE { + TIME_ARRAY_AGE, + TIME_ARRAY_REDSHIFT, + TIME_ARRAY_SCALE_FACTOR, +}; + struct time_array { double times[TIME_ARRAY_MAX_SIZE]; size_t size; @@ -119,7 +126,7 @@ void io_check_output_fields(const struct swift_params* params, void io_write_output_field_parameter(const char* filename); -void time_array_read_file(struct time_array *times, const char* filename); +void time_array_read_file(struct time_array *times, const char* filename, struct cosmology *cosmo); void time_array_print(const struct time_array *times); #endif /* SWIFT_COMMON_IO_H */ diff --git a/src/engine.c b/src/engine.c index cba6888ca831a9111c00784ce57f94bce14f9c5a..2b1b4804f2889c83c4abc4151ae1828e5f6d1080 100644 --- a/src/engine.c +++ b/src/engine.c @@ -6890,6 +6890,9 @@ void engine_struct_restore(struct engine *e, FILE *stream) { void engine_read_time_files(struct engine *e, const struct swift_params *params) { char filename[PARSER_MAX_LINE_SIZE]; + struct cosmology *cosmo = NULL; + if (e->policy & engine_policy_cosmology) + cosmo = e->cosmology; /* Read snapshot time array */ e->time_array_snapshots = (struct time_array*) malloc(sizeof(struct time_array)); @@ -6900,7 +6903,7 @@ void engine_read_time_files(struct engine *e, const struct swift_params *params) if (strcmp(filename, "")) { message("Reading snaplist file."); - time_array_read_file(e->time_array_snapshots, filename); + time_array_read_file(e->time_array_snapshots, filename, cosmo); } } @@ -6912,9 +6915,9 @@ void engine_read_next_snapshot_time(struct engine *e) { /* Find upper-bound on last output */ double time_end; if (is_cosmo) - time_end = e->cosmology->a_end * e->delta_time_snapshot; + time_end = e->cosmology->a_end; else - time_end = e->time_end + e->delta_time_snapshot; + time_end = e->time_end; /* Find next snasphot above current time */ double time;