add hotstart from checkpoint

7 files changed, 100 insertions, 16 deletions

diff --git a/input/options.input b/input/options.input
index 2ec212f..d6db646 100644
--- a/input/options.input
+++ b/input/options.input
@@ -3,8 +3,8 @@
 n_agent         : 1000						// Number of agents
 fuzz            : 0.01						// Fuzzy equivalence constant, used for geometry
 seed            : 108						// Random seed
-n_iteration     : 1000						// Number of iteration steps
 subgrid_dx      : 5.0						// Subgrid cell size
+start_time      : 0.0                                           // Simulation start time (s)
 
 surfacemesh	: ./example/ats_vis_surface_mesh.h5		// Surface mesh file
 surfacedata	: ./example/ats_vis_surface_data.h5		// Surface data file
@@ -16,4 +16,5 @@ subsurfacedata  : ./example/ats_vis_data.h5			// Subsurface data file
 subsurfacetime  : ./input/subtimes.txt			 	// Time levels for subsurface
 subsurfacegroup : ./input/subgroups.txt				// Group names for subsurface
 
-restartfile : none
\ No newline at end of file
+restartfile : ./input/checkpoint.input				// Restart file
+
diff --git a/main.cpp b/main.cpp
index 76277cf..96ac25f 100644
--- a/main.cpp
+++ b/main.cpp
@@ -89,7 +89,7 @@ int main(int argc, char* argv[])
                 step_surface ++;
                 t_surface = parser.surface_time[step_surface]; // get new time step
                 pmanager.environment.update_surface(); // update environment
-                t_combined = std::min(t_surface, t_subsurface);
+                t_combined = std::min(t_surface, t_subsurface); // update the time horizon
 
                 if (step_surface == nstep_surface - 1) {
                     fin_surface = true;
diff --git a/parser.hpp b/parser.hpp
index f0142a0..8d2e3dc 100644
--- a/parser.hpp
+++ b/parser.hpp
@@ -81,20 +81,22 @@ public:
     std::vector<double> subsurface_time;
 
     int n_agent;
-    int n_iteration;
     int seed;
 
-    double start_time = 0.0;
-    double end_time   = 1.0;
+    double start_time = -9999.0;
 
     double fuzz;
     double subgrid_dx;
 
+    bool restart = false;
+    std::vector<point> particle_coordinates;
+    std::vector<int> particle_domains;
+    std::vector<int> particle_cell_indices;
+
     void print_options()
     {
 
         std::cout << "[prs] number of agents: " << n_agent << std::endl;
-        std::cout << "[prs] number of iterations: " << n_iteration << std::endl;
         std::cout << "[prs] random seed: " << seed << std::endl;
         std::cout << "[prs] geometric fuzz: " << fuzz << std::endl;
         std::cout << "[prs] subgrid spacing: " << subgrid_dx << std::endl;
@@ -121,6 +123,19 @@ public:
             std::cout << subsurface_data_group[i] << " ";
         std::cout << std::endl;
         std::cout << "[prs] restart file: " << checkfile << std::endl;
+        if (restart) {
+            std::cout << "[prs] start time: " << start_time << " s (";
+            std::cout << start_time / (3600.0 * 24) << " d)" << std::endl;
+            std::cout << "[prs] ";
+            std::cout << particle_coordinates.size() << " particles read in";
+            // for (uint i = 0; i < particle_coordinates.size(); ++ i) {
+            //     std::cout << particle_cell_indices[i] << "(" << particle_domains[i] << ")";
+            //     std::cout << ": {" << particle_coordinates[i][0] << ", ";
+            //     std::cout << particle_coordinates[i][1] << ", " << particle_coordinates[i][2];
+            //     std::cout << "} ";
+            // }
+            std::cout << std::endl;
+        }
         std::cout << std::endl;
 
     }
@@ -135,7 +150,6 @@ public:
         // initialize all read-in values to -999
         fuzz    = -999;
         n_agent = -999;
-        n_iteration = -999;
 
         if (fInStream.is_open()) {
 
@@ -247,10 +261,41 @@ public:
                         pline.value >> subgrid_dx;
                     } else if (!strcmp("seed", pline.key.c_str())) {
                         pline.value >> seed;
-                    } else if (!strcmp("n_iteration", pline.key.c_str())) {
-                        pline.value >> n_iteration;
                     } else if (!strcmp("restartfile", pline.key.c_str())) {
+
                         pline.value >> checkfile;
+
+                        if (strcmp("none", checkfile.c_str())) {
+
+                            restart = true;
+
+                            std::ifstream in(checkfile);
+
+                            if (in.is_open()) {
+
+                                double x, y, z, t;
+                                int cell_index, domain;
+
+                                while (!in.eof()) {
+
+                                    in >> domain >> t >> domain >> cell_index >> x >> y >> z;
+
+                                    if (in.good()) {
+                                        particle_domains.push_back(domain);
+                                        particle_cell_indices.push_back(cell_index);
+                                        particle_coordinates.push_back({ x, y, z });
+                                    }
+
+                                }
+
+                            } else {
+                                std::cerr << "[err] file " << checkfile << " can't be opened" << std::endl;
+                            }
+
+                        }
+
+                    } else if (!strcmp("start_time", pline.key.c_str())) {
+                        pline.value >> start_time;
                     } else {
                         std::cerr << "[err] unrecognised option key" << std::endl;
                     }
diff --git a/particle.hpp b/particle.hpp
index 5f95ff5..8111ebd 100644
--- a/particle.hpp
+++ b/particle.hpp
@@ -73,8 +73,9 @@ public:
         fprintf(fp, "#time\tdomain\tx\ty\tz\n");
 
         for (size_t i = 0; i < time.size(); i ++) {
-            fprintf(fp, "%f\t%d\t%f\t%f\t%f\n", time[i],
+            fprintf(fp, "%f\t%d\t%d\t%f\t%f\t%f\n", time[i],
                     domain[i],
+                    cell_index,
                     coordinates[i][0], coordinates[i][1],
                     coordinates[i][2]);
         }
@@ -116,7 +117,7 @@ public:
         rng = std::default_random_engine(parser.seed);
         dist = std::uniform_real_distribution<double>(0.0, 1.0);
 
-        if (!strcmp("none", parser.checkfile.c_str())) {
+        if (!parser.restart) {
 
             for (int i = 0; i < parser.n_agent; ++ i) {
 
@@ -156,6 +157,27 @@ public:
 
             }
 
+        } else {
+
+            for (uint i = 0; i < parser.particle_coordinates.size(); ++ i) {
+
+                Particle p;
+
+                p.cell_index = parser.particle_cell_indices[i];
+                p.domain.push_back(parser.particle_domains[i]);
+
+                p.x = parser.particle_coordinates[i][0];
+                p.y = parser.particle_coordinates[i][1];
+                p.z = parser.particle_coordinates[i][2];
+                p.t = parser.start_time;
+
+                p.time.push_back(p.t);
+                p.coordinates.push_back(parser.particle_coordinates[i]);
+
+                particles.push_back(p);
+
+            }
+
         }
 
     }
diff --git a/tools/generate_checkpoint.sh b/tools/generate_checkpoint.sh
new file mode 100644
index 0000000..4dcf7b7
--- /dev/null
+++ b/tools/generate_checkpoint.sh
@@ -0,0 +1,16 @@
+#!/usr/bin/env bash
+
+DIRECTORY="../output"
+TARGETDIR="../input"
+COUNTER=0
+
+trash ${TARGETDIR}/checkpoint.input
+touch ${TARGETDIR}/checkpoint.input
+
+for FILE in ${DIRECTORY}/*.csv
+do
+    echo "${COUNTER}	`tail -1 ${FILE}`" | awk '{print $1, $2, $3, $4, $5, $6, $7}' >> ${TARGETDIR}/checkpoint.input
+    let COUNTER=COUNTER+1
+done
+
+echo "[OK] created ${TARGETDIR}/checkpoint.input"
diff --git a/tools/plot-paths-integrated.gp b/tools/plot-paths-integrated.gp
index 723c693..d050d89 100644
--- a/tools/plot-paths-integrated.gp
+++ b/tools/plot-paths-integrated.gp
@@ -21,7 +21,7 @@ set size ratio -1
 
 splot \
 	for [i = 100:999] "../output/p-000".i.".csv" \
-	u ($3/1000):($4/1000):($5/1000)\
+	u ($4/1000):($5/1000):($6/1000)\
 	with lines\
 	notitle
 
diff --git a/tools/plot-paths.gp b/tools/plot-paths.gp
index 6bc2fb1..6c19ced 100644
--- a/tools/plot-paths.gp
+++ b/tools/plot-paths.gp
@@ -1,7 +1,7 @@
 #!/usr/bin/env gnuplot
 
 set term png font "hack"
-set output "output/figure.png"
+set output "../output/figure_surface.png"
 
 set border linecolor rgbcolor "white"
 set key textcolor rgbcolor "white"
@@ -18,8 +18,8 @@ set ytics 1
 set size ratio -1
 
 plot \
-	for [i = 100:999] "output/p-000".i.".csv" \
-	u ($3/1000):($4/1000)\
+	for [i = 100:999] "../output/p-000".i.".csv" \
+	u ($4/1000):($5/1000)\
 	with lines\
 	notitle