From e78c4d7a52a7b51aaa485fc4d73661237be8400a Mon Sep 17 00:00:00 2001
From: Matthieu Schaller <matthieu.schaller@durham.ac.uk>
Date: Sun, 24 Apr 2016 23:06:46 +0100
Subject: [PATCH] Read snapshot basename and interval between snapshots from
 yaml file.

---
 examples/UniformBox/uniformBox.yml |  6 ++++++
 examples/main.c                    | 22 ++++++++++++----------
 examples/parameter_example.yml     |  6 ++++++
 src/common_io.c                    | 23 +++++++++++++++--------
 src/common_io.h                    |  4 ++--
 src/engine.c                       | 15 +++++++++++++++
 src/engine.h                       |  4 ++++
 src/parallel_io.c                  | 27 ++++++++++++++-------------
 src/parallel_io.h                  |  6 +++---
 src/serial_io.c                    | 19 +++++++++----------
 src/serial_io.h                    |  5 +++--
 src/single_io.c                    | 19 +++++++++----------
 src/single_io.h                    |  3 ++-
 13 files changed, 100 insertions(+), 59 deletions(-)

diff --git a/examples/UniformBox/uniformBox.yml b/examples/UniformBox/uniformBox.yml
index 0474b0f820..dcde96be9e 100644
--- a/examples/UniformBox/uniformBox.yml
+++ b/examples/UniformBox/uniformBox.yml
@@ -21,6 +21,12 @@ TimeIntegration:
   dt_min:     1e-6  # The minimal time-step size of the simulation (in internal units).
   dt_max:     1e-2  # The maximal time-step size of the simulation (in internal units).
 
+# Parameters governing the snapshots
+Snapshots:
+  basename:   uniformBox # Common part of the name of output files
+  time_first: 0.         # Time of the frist output (in internal units)
+  delta_time: 0.01       # Time difference between consecutive outputs (in internal units)
+
 # Parameters for the hydrodynamics scheme
 SPH:
   resolution_eta:        1.2349   # Target smoothing length in units of the mean inter-particle separation (1.2349 == 48Ngbs with the cubic spline kernel).
diff --git a/examples/main.c b/examples/main.c
index 481a6bf43a..e8b4e08a06 100644
--- a/examples/main.c
+++ b/examples/main.c
@@ -422,24 +422,22 @@ int main(int argc, char *argv[]) {
     fflush(stdout);
   }
 
-  /* Now that everything is ready, no need for the parameters any more */
-  free(params);
-  params = NULL;
-
   int with_outputs = 1;
+  char baseName[200];
+  parser_get_param_string(params, "Snapshots:basename", baseName);
   if (with_outputs && !dry_run) {
     /* Write the state of the system before starting time integration. */
     if (myrank == 0) clocks_gettime(&tic);
 #if defined(WITH_MPI)
 #if defined(HAVE_PARALLEL_HDF5)
-    write_output_parallel(&e, &us, myrank, nr_nodes, MPI_COMM_WORLD,
+    write_output_parallel(&e, baseName, &us, myrank, nr_nodes, MPI_COMM_WORLD,
                           MPI_INFO_NULL);
 #else
-    write_output_serial(&e, &us, myrank, nr_nodes, MPI_COMM_WORLD,
+    write_output_serial(&e, baseName, &us, myrank, nr_nodes, MPI_COMM_WORLD,
                         MPI_INFO_NULL);
 #endif
 #else
-    write_output_single(&e, &us);
+    write_output_single(&e, baseName, &us);
 #endif
     if (myrank == 0 && verbose) {
       clocks_gettime(&toc);
@@ -449,6 +447,10 @@ int main(int argc, char *argv[]) {
     }
   }
 
+  /* Now that everything is ready, no need for the parameters any more */
+  free(params);
+  params = NULL;
+
 /* Init the runner history. */
 #ifdef HIST
   for (k = 0; k < runner_hist_N; k++) runner_hist_bins[k] = 0;
@@ -505,7 +507,7 @@ int main(int argc, char *argv[]) {
     /* Take a step. */
     engine_step(&e);
 
-    if (with_outputs && j % 100 == 0) {
+    if (with_outputs && j % 10 == 0) {
 
       if (myrank == 0) clocks_gettime(&tic);
 #if defined(WITH_MPI)
@@ -517,7 +519,7 @@ int main(int argc, char *argv[]) {
                           MPI_INFO_NULL);
 #endif
 #else
-      write_output_single(&e, &us);
+      write_output_single(&e, baseName, &us);
 #endif
       if (myrank == 0 && verbose) {
         clocks_gettime(&toc);
@@ -621,7 +623,7 @@ int main(int argc, char *argv[]) {
                         MPI_INFO_NULL);
 #endif
 #else
-    write_output_single(&e, &us);
+    write_output_single(&e, baseName, &us);
 #endif
     if (myrank == 0 && verbose) {
       clocks_gettime(&toc);
diff --git a/examples/parameter_example.yml b/examples/parameter_example.yml
index 228e748f3b..38b0a7c76f 100644
--- a/examples/parameter_example.yml
+++ b/examples/parameter_example.yml
@@ -21,6 +21,12 @@ TimeIntegration:
   dt_min:     1e-6  # The minimal time-step size of the simulation (in internal units).
   dt_max:     1e-2  # The maximal time-step size of the simulation (in internal units).
 
+# Parameters governing the snapshots
+Snapshots:
+  basename:   uniformBox # Common part of the name of output files
+  time_first: 0.         # Time of the frist output (in internal units)
+  delta_time: 0.01       # Time difference between consecutive outputs (in internal units)
+
 # Parameters for the hydrodynamics scheme
 SPH:
   resolution_eta:        1.2349   # Target smoothing length in units of the mean inter-particle separation (1.2349 == 48Ngbs with the cubic spline kernel).
diff --git a/src/common_io.c b/src/common_io.c
index 2a635723d5..143c1ca941 100644
--- a/src/common_io.c
+++ b/src/common_io.c
@@ -335,11 +335,15 @@ void writeCodeDescription(hid_t h_file) {
  *
  * @todo Use a proper XML library to avoid stupid copies.
  */
-FILE* prepareXMFfile() {
+FILE* prepareXMFfile(const char* baseName) {
   char buffer[1024];
 
-  FILE* xmfFile = fopen("output.xmf", "r");
-  FILE* tempFile = fopen("output_temp.xmf", "w");
+  char fileName[FILENAME_BUFFER_SIZE];
+  char tempFileName[FILENAME_BUFFER_SIZE];
+  snprintf(fileName, FILENAME_BUFFER_SIZE, "%s.xmf", baseName);
+  snprintf(tempFileName, FILENAME_BUFFER_SIZE, "%s_temp.xmf", baseName);
+  FILE* xmfFile = fopen(fileName, "r");
+  FILE* tempFile = fopen(tempFileName, "w");
 
   if (xmfFile == NULL) error("Unable to open current XMF file.");
 
@@ -355,8 +359,8 @@ FILE* prepareXMFfile() {
   fclose(xmfFile);
 
   /* We then copy the XMF file back up to the closing lines */
-  xmfFile = fopen("output.xmf", "w");
-  tempFile = fopen("output_temp.xmf", "r");
+  xmfFile = fopen(fileName, "w");
+  tempFile = fopen(tempFileName, "r");
 
   if (xmfFile == NULL) error("Unable to open current XMF file.");
 
@@ -369,7 +373,7 @@ FILE* prepareXMFfile() {
   }
   fprintf(xmfFile, "\n");
   fclose(tempFile);
-  remove("output_temp.xmf");
+  remove(tempFileName);
 
   return xmfFile;
 }
@@ -380,8 +384,11 @@ FILE* prepareXMFfile() {
  * @todo Exploit the XML nature of the XMF format to write a proper XML writer
  *and simplify all the XMF-related stuff.
  */
-void createXMFfile() {
-  FILE* xmfFile = fopen("output.xmf", "w");
+void createXMFfile(const char* baseName) {
+
+  char fileName[FILENAME_BUFFER_SIZE];
+  snprintf(fileName, FILENAME_BUFFER_SIZE, "%s.xmf", baseName);
+  FILE* xmfFile = fopen(fileName, "w");
 
   fprintf(xmfFile, "<?xml version=\"1.0\" ?> \n");
   fprintf(xmfFile, "<!DOCTYPE Xdmf SYSTEM \"Xdmf.dtd\" []> \n");
diff --git a/src/common_io.h b/src/common_io.h
index b7f3a1a317..70ed25993c 100644
--- a/src/common_io.h
+++ b/src/common_io.h
@@ -97,8 +97,8 @@ void writeAttribute_i(hid_t grp, char* name, int data);
 void writeAttribute_l(hid_t grp, char* name, long data);
 void writeAttribute_s(hid_t grp, char* name, const char* str);
 
-void createXMFfile();
-FILE* prepareXMFfile();
+void createXMFfile(const char* baseName);
+FILE* prepareXMFfile(const char* baseName);
 void writeXMFoutputheader(FILE* xmfFile, char* hdfFileName, float time);
 void writeXMFoutputfooter(FILE* xmfFile, int outputCount, float time);
 void writeXMFgroupheader(FILE* xmfFile, char* hdfFileName, size_t N,
diff --git a/src/engine.c b/src/engine.c
index 7af7c05661..ab250aaea4 100644
--- a/src/engine.c
+++ b/src/engine.c
@@ -2381,6 +2381,11 @@ void engine_init(struct engine *e, struct space *s,
   e->ti_old = 0;
   e->ti_current = 0;
   e->timeStep = 0.;
+  e->timeBase = 0.;
+  e->timeFirstSnapshot =
+      parser_get_param_double(params, "Snapshots:time_first");
+  e->deltaTimeSnapshot =
+      parser_get_param_double(params, "Snapshots:delta_time");
   e->dt_min = parser_get_param_double(params, "TimeIntegration:dt_min");
   e->dt_max = parser_get_param_double(params, "TimeIntegration:dt_max");
   e->file_stats = NULL;
@@ -2557,6 +2562,16 @@ void engine_init(struct engine *e, struct space *s,
     error("Maximal time-step size larger than the simulation run time t=%e",
           e->timeEnd - e->timeBegin);
 
+  /* Deal with outputs */
+  if (e->deltaTimeSnapshot < 0.)
+    error("Time between snapshots (%e) must be positive.",
+          e->deltaTimeSnapshot);
+
+  if (e->timeFirstSnapshot < e->timeBegin)
+    error(
+        "Time of first snapshot (%e) must be after the simulation start t=%e.",
+        e->timeFirstSnapshot, e->timeBegin);
+
 /* Construct types for MPI communications */
 #ifdef WITH_MPI
   part_create_mpi_types();
diff --git a/src/engine.h b/src/engine.h
index 82e76481f8..3acc53e012 100644
--- a/src/engine.h
+++ b/src/engine.h
@@ -131,6 +131,10 @@ struct engine {
   /* Time base */
   double timeBase;
 
+  /* Snapshot times */
+  double timeFirstSnapshot;
+  double deltaTimeSnapshot;
+
   /* File for statistics */
   FILE *file_stats;
 
diff --git a/src/parallel_io.c b/src/parallel_io.c
index d1c739b590..4579be8f04 100644
--- a/src/parallel_io.c
+++ b/src/parallel_io.c
@@ -509,8 +509,12 @@ void read_ic_parallel(char* fileName, double dim[3], struct part** parts,
  *its XMF descriptor
  *
  * @param e The engine containing all the system.
- * @param us The UnitSystem used for the conversion of units
- *in the output
+ * @param baseName The common part of the snapshot file name.
+ * @param us The UnitSystem used for the conversion of units in the output.
+ * @param mpi_rank The MPI rank of this node.
+ * @param mpi_size The number of MPI ranks.
+ * @param comm The MPI communicator.
+ * @param info The MPI information object
  *
  * Creates an HDF5 output file and writes the particles
  *contained
@@ -522,9 +526,9 @@ void read_ic_parallel(char* fileName, double dim[3], struct part** parts,
  * Calls #error() if an error occurs.
  *
  */
-void write_output_parallel(struct engine* e, struct UnitSystem* us,
-                           int mpi_rank, int mpi_size, MPI_Comm comm,
-                           MPI_Info info) {
+void write_output_parallel(struct engine* e, const char* baseName,
+                           struct UnitSystem* us, int mpi_rank, int mpi_size,
+                           MPI_Comm comm, MPI_Info info) {
   hid_t h_file = 0, h_grp = 0, h_grpsph = 0;
   const size_t Ngas = e->s->nr_parts;
   const size_t Ntot = e->s->nr_gparts;
@@ -536,22 +540,19 @@ void write_output_parallel(struct engine* e, struct UnitSystem* us,
   static int outputCount = 0;
   FILE* xmfFile = 0;
 
-  /* Number of particles of each type */
-  // const size_t Ndm = Ntot - Ngas;
-
-  /* MATTHIEU: Temporary fix to preserve master */
+  /* Number of unassociated gparts */
   const size_t Ndm = Ntot > 0 ? Ntot - Ngas : 0;
-  /* MATTHIEU: End temporary fix */
 
   /* File name */
   char fileName[FILENAME_BUFFER_SIZE];
-  snprintf(fileName, FILENAME_BUFFER_SIZE, "output_%03i.hdf5", outputCount);
+  snprintf(fileName, FILENAME_BUFFER_SIZE, "%s_%03i.hdf5", baseName,
+           outputCount);
 
   /* First time, we need to create the XMF file */
-  if (outputCount == 0 && mpi_rank == 0) createXMFfile();
+  if (outputCount == 0 && mpi_rank == 0) createXMFfile(baseName);
 
   /* Prepare the XMF file for the new entry */
-  if (mpi_rank == 0) xmfFile = prepareXMFfile();
+  if (mpi_rank == 0) xmfFile = prepareXMFfile(baseName);
 
   /* Open HDF5 file */
   hid_t plist_id = H5Pcreate(H5P_FILE_ACCESS);
diff --git a/src/parallel_io.h b/src/parallel_io.h
index f3691cb29b..970ad8c41d 100644
--- a/src/parallel_io.h
+++ b/src/parallel_io.h
@@ -36,9 +36,9 @@ void read_ic_parallel(char* fileName, double dim[3], struct part** parts,
                       int* periodic, int mpi_rank, int mpi_size, MPI_Comm comm,
                       MPI_Info info, int dry_run);
 
-void write_output_parallel(struct engine* e, struct UnitSystem* us,
-                           int mpi_rank, int mpi_size, MPI_Comm comm,
-                           MPI_Info info);
+void write_output_parallel(struct engine* e, const char* baseName,
+                           struct UnitSystem* us, int mpi_rank, int mpi_size,
+                           MPI_Comm comm, MPI_Info info);
 
 #endif
 
diff --git a/src/serial_io.c b/src/serial_io.c
index 10eab97f1b..f2c0dcce1c 100644
--- a/src/serial_io.c
+++ b/src/serial_io.c
@@ -590,6 +590,7 @@ void read_ic_serial(char* fileName, double dim[3], struct part** parts,
  * @brief Writes an HDF5 output file (GADGET-3 type) with its XMF descriptor
  *
  * @param e The engine containing all the system.
+ * @param baseName The common part of the snapshot file name.
  * @param us The UnitSystem used for the conversion of units in the output.
  * @param mpi_rank The MPI rank of this node.
  * @param mpi_size The number of MPI ranks.
@@ -604,8 +605,9 @@ void read_ic_serial(char* fileName, double dim[3], struct part** parts,
  * Calls #error() if an error occurs.
  *
  */
-void write_output_serial(struct engine* e, struct UnitSystem* us, int mpi_rank,
-                         int mpi_size, MPI_Comm comm, MPI_Info info) {
+void write_output_serial(struct engine* e, const char* baseName,
+                         struct UnitSystem* us, int mpi_rank, int mpi_size,
+                         MPI_Comm comm, MPI_Info info) {
   hid_t h_file = 0, h_grp = 0, h_grpsph = 0;
   const size_t Ngas = e->s->nr_parts;
   const size_t Ntot = e->s->nr_gparts;
@@ -617,16 +619,13 @@ void write_output_serial(struct engine* e, struct UnitSystem* us, int mpi_rank,
   static int outputCount = 0;
   FILE* xmfFile = 0;
 
-  /* Number of particles of each type */
-  // const size_t Ndm = Ntot - Ngas;
-
-  /* MATTHIEU: Temporary fix to preserve master */
+  /* Number of unassociated gparts */
   const size_t Ndm = Ntot > 0 ? Ntot - Ngas : 0;
-  /* MATTHIEU: End temporary fix */
 
   /* File name */
   char fileName[FILENAME_BUFFER_SIZE];
-  snprintf(fileName, FILENAME_BUFFER_SIZE, "output_%03i.hdf5", outputCount);
+  snprintf(fileName, FILENAME_BUFFER_SIZE, "%s_%03i.hdf5", baseName,
+           outputCount);
 
   /* Compute offset in the file and total number of particles */
   size_t N[NUM_PARTICLE_TYPES] = {Ngas, Ndm, 0};
@@ -646,10 +645,10 @@ void write_output_serial(struct engine* e, struct UnitSystem* us, int mpi_rank,
   if (mpi_rank == 0) {
 
     /* First time, we need to create the XMF file */
-    if (outputCount == 0) createXMFfile();
+    if (outputCount == 0) createXMFfile(baseName);
 
     /* Prepare the XMF file for the new entry */
-    xmfFile = prepareXMFfile();
+    xmfFile = prepareXMFfile(baseName);
 
     /* Write the part corresponding to this specific output */
     writeXMFoutputheader(xmfFile, fileName, e->time);
diff --git a/src/serial_io.h b/src/serial_io.h
index 74ab8326db..b7ed6eb62d 100644
--- a/src/serial_io.h
+++ b/src/serial_io.h
@@ -36,8 +36,9 @@ void read_ic_serial(char* fileName, double dim[3], struct part** parts,
                     int* periodic, int mpi_rank, int mpi_size, MPI_Comm comm,
                     MPI_Info info, int dry_run);
 
-void write_output_serial(struct engine* e, struct UnitSystem* us, int mpi_rank,
-                         int mpi_size, MPI_Comm comm, MPI_Info info);
+void write_output_serial(struct engine* e, const char* baseName,
+                         struct UnitSystem* us, int mpi_rank, int mpi_size,
+                         MPI_Comm comm, MPI_Info info);
 
 #endif
 
diff --git a/src/single_io.c b/src/single_io.c
index 1dc71087e1..d545fb086f 100644
--- a/src/single_io.c
+++ b/src/single_io.c
@@ -456,7 +456,8 @@ void read_ic_single(char* fileName, double dim[3], struct part** parts,
  * @brief Writes an HDF5 output file (GADGET-3 type) with its XMF descriptor
  *
  * @param e The engine containing all the system.
- * @param us The UnitSystem used for the conversion of units in the output
+ * @param baseName The common part of the snapshot file name.
+ * @param us The UnitSystem used for the conversion of units in the output.
  *
  * Creates an HDF5 output file and writes the particles contained
  * in the engine. If such a file already exists, it is erased and replaced
@@ -466,7 +467,8 @@ void read_ic_single(char* fileName, double dim[3], struct part** parts,
  * Calls #error() if an error occurs.
  *
  */
-void write_output_single(struct engine* e, struct UnitSystem* us) {
+void write_output_single(struct engine* e, const char* baseName,
+                         struct UnitSystem* us) {
 
   hid_t h_file = 0, h_grp = 0;
   const size_t Ngas = e->s->nr_parts;
@@ -478,25 +480,22 @@ void write_output_single(struct engine* e, struct UnitSystem* us) {
   struct gpart* dmparts = NULL;
   static int outputCount = 0;
 
-  /* Number of particles of each type */
-  // const size_t Ndm = Ntot - Ngas;
-
-  /* MATTHIEU: Temporary fix to preserve master */
+  /* Number of unassociated gparts */
   const size_t Ndm = Ntot > 0 ? Ntot - Ngas : 0;
-  /* MATTHIEU: End temporary fix */
 
   long long N_total[NUM_PARTICLE_TYPES] = {Ngas, Ndm, 0};
 
   /* File name */
   char fileName[FILENAME_BUFFER_SIZE];
-  snprintf(fileName, FILENAME_BUFFER_SIZE, "output_%03i.hdf5", outputCount);
+  snprintf(fileName, FILENAME_BUFFER_SIZE, "%s_%03i.hdf5", baseName,
+           outputCount);
 
   /* First time, we need to create the XMF file */
-  if (outputCount == 0) createXMFfile();
+  if (outputCount == 0) createXMFfile(baseName);
 
   /* Prepare the XMF file for the new entry */
   FILE* xmfFile = 0;
-  xmfFile = prepareXMFfile();
+  xmfFile = prepareXMFfile(baseName);
 
   /* Write the part corresponding to this specific output */
   writeXMFoutputheader(xmfFile, fileName, e->time);
diff --git a/src/single_io.h b/src/single_io.h
index 587ebe07b6..adfc5b4394 100644
--- a/src/single_io.h
+++ b/src/single_io.h
@@ -30,7 +30,8 @@ void read_ic_single(char* fileName, double dim[3], struct part** parts,
                     struct gpart** gparts, size_t* Ngas, size_t* Ndm,
                     int* periodic, int dry_run);
 
-void write_output_single(struct engine* e, struct UnitSystem* us);
+void write_output_single(struct engine* e, const char* baseName,
+                         struct UnitSystem* us);
 
 #endif
 
-- 
GitLab