Carpet: Set times of past timelevels correctly

Redesign the way in which the times of the past timelevels are set. If cctk_delta_time changes during initialisation, re-set the times of the past timelevels accordingly.

5 files changed, 74 insertions, 12 deletions

diff --git a/Carpet/Carpet/src/CallFunction.cc b/Carpet/Carpet/src/CallFunction.cc
index ffd695c45..73c34d73c 100644
--- a/Carpet/Carpet/src/CallFunction.cc
+++ b/Carpet/Carpet/src/CallFunction.cc
@@ -1,6 +1,7 @@
 #include <algorithm>
 #include <cassert>
 #include <cstdlib>
+#include <cstring>
 #include <map>
 #include <string>
 
@@ -342,6 +343,9 @@ namespace Carpet {
       }
       Timer & timer = * ti->second;
       
+      // Save the time step size
+      CCTK_REAL const saved_cctk_delta_time = cctkGH->cctk_delta_time;
+      
       user_timer.start();
       timer.start();
       int const res = CCTK_CallFunction (function, attribute, data);
@@ -349,6 +353,32 @@ namespace Carpet {
       timer.stop();
       user_timer.stop();
       
+      // Manage the time step size. If the time step size changes
+      // during initialisation, assume it is thorn Time, and update
+      // the time hierarchy. If it changes during evolution, assume it
+      // is MoL, and do nothing.
+      if (cctkGH->cctk_iteration == 0 and
+          cctkGH->cctk_delta_time != saved_cctk_delta_time)
+      {
+        // The user changed cctk_delta_time during initialisation --
+        // update our internals and the time hierarchy
+        delta_time = cctkGH->cctk_delta_time * timereflevelfact / mglevelfact;
+        for (int ml=0; ml<mglevels; ++ml) {
+          for (int rl=0; rl<reflevels; ++rl) {
+            // Update the time delta
+            CCTK_REAL const dt =
+              delta_time / timereffacts.AT(rl) * ipow(mgfact, ml);
+            tt->set_delta(ml,rl,dt);
+            CCTK_REAL const t0 = tt->get_time(ml,rl,0);
+            // Update the times of the past timelevels
+            for (int tl=1; tl<timelevels; ++tl) {
+              CCTK_REAL const t = t0 - tl * dt;
+              tt->set_time(ml,rl,tl,t);
+            }
+          }
+        }
+      }
+      
     }
     CallAfterRoutines (cctkGH, function, attribute, data);
   }
diff --git a/Carpet/Carpet/src/Evolve.cc b/Carpet/Carpet/src/Evolve.cc
index 1747df301..8337e98aa 100644
--- a/Carpet/Carpet/src/Evolve.cc
+++ b/Carpet/Carpet/src/Evolve.cc
@@ -298,7 +298,7 @@ namespace Carpet {
                 do_global_mode = do_late_global_mode;
                 do_meta_mode = do_late_meta_mode;
                 
-                BEGIN_TIMELEVEL_LOOP (cctkGH) {
+                BEGIN_TIMELEVEL_LOOP(cctkGH) {
                   
                   Waypoint ("Postregrid at iteration %d time %g timelevel %d%s%s",
                             cctkGH->cctk_iteration,
diff --git a/Carpet/Carpet/src/Initialise.cc b/Carpet/Carpet/src/Initialise.cc
index 557c81230..9a82d2c23 100644
--- a/Carpet/Carpet/src/Initialise.cc
+++ b/Carpet/Carpet/src/Initialise.cc
@@ -273,7 +273,7 @@ namespace Carpet {
           do_global_mode = do_early_global_mode; // on first iteration, coarsest grid
           do_meta_mode = do_early_meta_mode; // on first iteration, coarsest grid
           
-          BEGIN_TIMELEVEL_LOOP (cctkGH) {
+          BEGIN_TIMELEVEL_LOOP(cctkGH) {
             
             Waypoint ("Recovering II at iteration %d time %g timelevel %d%s%s",
                       cctkGH->cctk_iteration,
@@ -285,11 +285,11 @@ namespace Carpet {
             // Post recover variables
             ScheduleTraverse (where, "CCTK_POST_RECOVER_VARIABLES", cctkGH);
             
-            // Checking
-            PoisonCheck (cctkGH, currenttime);
-            
           } END_TIMELEVEL_LOOP;
           
+          // Checking
+          PoisonCheck (cctkGH, currenttime);
+          
           CheckChecksums (cctkGH, allbutcurrenttime);
           
           EndTimingLevel (cctkGH);
@@ -360,7 +360,7 @@ namespace Carpet {
           
           if (init_each_timelevel) {
             
-            BEGIN_TIMELEVEL_LOOP (cctkGH) {
+            BEGIN_TIMELEVEL_LOOP(cctkGH) {
               
               // Set up the initial data
               ScheduleTraverse (where, "CCTK_INITIAL", cctkGH);
@@ -841,8 +841,8 @@ namespace Carpet {
                 do_global_mode = do_late_global_mode;
                 do_meta_mode = do_late_meta_mode;
                 
-                BEGIN_TIMELEVEL_LOOP (cctkGH) {
-                
+                BEGIN_TIMELEVEL_LOOP(cctkGH) {
+                  
                   Waypoint ("Postregrid at iteration %d time %g timelevel %d%s%s",
                             cctkGH->cctk_iteration,
                             (double)cctkGH->cctk_time,
@@ -852,7 +852,7 @@ namespace Carpet {
                   
                   // Postregrid
                   ScheduleTraverse (where, "CCTK_POSTREGRID", cctkGH);
-                  
+                
                 } END_TIMELEVEL_LOOP;
                 
                 EndTimingLevel (cctkGH);
@@ -1009,7 +1009,7 @@ namespace Carpet {
                 
                 if (init_each_timelevel) {
                   
-                  BEGIN_TIMELEVEL_LOOP (cctkGH) {
+                  BEGIN_TIMELEVEL_LOOP(cctkGH) {
                     
                     // Postregrid
                     ScheduleTraverse (where, "CCTK_POSTREGRIDINITIAL", cctkGH);
diff --git a/Carpet/Carpet/src/modes.cc b/Carpet/Carpet/src/modes.cc
index d101a509b..ec26ada6a 100644
--- a/Carpet/Carpet/src/modes.cc
+++ b/Carpet/Carpet/src/modes.cc
@@ -223,8 +223,7 @@ namespace Carpet {
     
     Checkpoint ("Leaving global mode");
 
-    // Save and unset time delta
-    delta_time = cctkGH->cctk_delta_time / mglevelfact;
+    // Unset time delta
     cctkGH->cctk_delta_time = 0.0;
     if (maps == 1) {
       // Save and unset space delta
diff --git a/Carpet/Carpet/src/modes.hh b/Carpet/Carpet/src/modes.hh
index 37d3ec646..fe322c7be 100644
--- a/Carpet/Carpet/src/modes.hh
+++ b/Carpet/Carpet/src/modes.hh
@@ -269,6 +269,39 @@ namespace Carpet {
   
   
   
+#if 0
+// TODO: Introduce such maybe-loops for the other loops as well
+
+#define BEGIN_MAYBE_TIMELEVEL_LOOP(cctkGH)                              \
+  do {                                                                  \
+    bool timelevel_maybe_loop_ = true;                                  \
+    assert (do_allow_past_timelevels);                                  \
+    int const min_tl = timelevel == -1 ? 0 : timelevel;                 \
+    int const max_tl = timelevel == -1 ? timelevels : timelevel+1;      \
+    int const old_timelevel = timelevel;                                \
+    bool const old_do_allow_past_timelevels = do_allow_past_timelevels; \
+    do_allow_past_timelevels = timelevel != -1;                         \
+    for (timelevel = min_tl; timelevel < max_tl; ++ timelevel) {        \
+      if (mglevel != -1 and reflevel != -1) {                           \
+        cctkGH->cctk_time = tt->get_time (mglevel, reflevel, timelevel); \
+      }                                                                 \
+      {
+#define END_MAYBE_TIMELEVEL_LOOP                                \
+      }                                                         \
+    }                                                           \
+    assert (timelevel_maybe_loop_);                             \
+    timelevel_maybe_loop_ = false;                              \
+    timelevel = old_timelevel;                                  \
+    if (mglevel != -1 and reflevel != -1) {                     \
+      int const tl = timelevel == -1 ? 0 : timelevel;           \
+      cctkGH->cctk_time = tt->get_time (mglevel, reflevel, tl); \
+    }                                                           \
+    do_allow_past_timelevels = old_do_allow_past_timelevels;    \
+  } while (false)
+#endif
+
+
+
   // Mode escapes
   
   class singlemap_escape {