Correct schedule, in particular for checkpoint/recovery

Do not mark ADMBase variables for non-checkpointing if they have
multiple timelevels. (Variables with multiple timelevels must always
be checkpointed, because the past timelevels cannot be regenerated
after recovery.)

Finally remove all perl post-processing of the auto-generated code;
instead, use proper Kranc mechanisms.

Schedule the ADM constraints and ADM quantities after MoL_PostStep,
since this is where the ADMBase variables are set.

Schedule enforcing the BSSN constraints in the new schedule group
MoL_PostStepModify, since they should not be enforced after recovery.
(This would lead to inconsistencies at floating-point round-off
level.)

Regenerate all thorns.

5 files changed, 67 insertions, 42 deletions

diff --git a/m/Makefile b/m/Makefile
index 48c7197..6c6bde5 100644
--- a/m/Makefile
+++ b/m/Makefile
@@ -8,54 +8,52 @@ all: McLachlan_ADM.out McLachlan_BSSN.out McLachlan_ADMConstraints.out McLachlan
 McLachlan_ADM.out: McLachlan_ADM.m
 	rm -rf ML_ADM*
 	./runmath.sh $^
-	for thorn in ML_ADM*; do \
-		./copy-if-changed.sh $$thorn ../$$thorn; \
+	for thorn in ML_ADM*; do				\
+		./copy-if-changed.sh $$thorn ../$$thorn;	\
 	done
 
 McLachlan_BSSN.out: McLachlan_BSSN.m
 	rm -rf ML_BSSN*
 	./runmath.sh $^
-	for thorn in ML_BSSN*; do \
-		./copy-if-changed.sh $$thorn ../$$thorn; \
-		./create-helper-thorn.sh $$thorn; \
+	for thorn in ML_BSSN*; do					   \
+		./copy-if-changed.sh $$thorn ../$$thorn &&		   \
+		./create-helper-thorn.sh $$thorn &&			   \
 		./copy-if-changed.sh $${thorn}_Helper ../$${thorn}_Helper; \
 	done
 
 McLachlan_ADMConstraints.out: McLachlan_ADMConstraints.m
 	rm -rf ML_ADMConstraints*
 	./runmath.sh $^
-	perl -pi -e 's/MoL_PostStep/MoL_PseudoEvolution/g' ML_ADMConstraints*/schedule.ccl
-	for thorn in ML_ADMConstraints*; do \
-		./copy-if-changed.sh $$thorn ../$$thorn; \
+	for thorn in ML_ADMConstraints*; do			\
+		./copy-if-changed.sh $$thorn ../$$thorn;	\
 	done
 
 McLachlan_ADMQuantities.out: McLachlan_ADMQuantities.m
 	rm -rf ML_ADMQuantities*
 	./runmath.sh $^
-	perl -pi -e 's/MoL_PostStep/MoL_PseudoEvolution/g' ML_ADMQuantities*/schedule.ccl
-	for thorn in ML_ADMQuantities*; do \
-		./copy-if-changed.sh $$thorn ../$$thorn; \
+	for thorn in ML_ADMQuantities*; do			\
+		./copy-if-changed.sh $$thorn ../$$thorn;	\
 	done
 
 WaveToy.out: WaveToy.m
 	rm -rf ML_WaveToy
 	./runmath.sh $^
-	for thorn in ML_WaveToy; do \
-		./copy-if-changed.sh $$thorn ../$$thorn; \
+	for thorn in ML_WaveToy; do				\
+		./copy-if-changed.sh $$thorn ../$$thorn;	\
 	done
 
 WaveToyFO.out: WaveToyFO.m
 	rm -rf ML_WaveToyFO
 	./runmath.sh $^
-	for thorn in ML_WaveToyFO; do \
-		./copy-if-changed.sh $$thorn ../$$thorn; \
+	for thorn in ML_WaveToyFO; do				\
+		./copy-if-changed.sh $$thorn ../$$thorn;	\
 	done
 
 hydro.out: hydro.m
 	rm -rf ML_hydro
 	./runmath.sh $^
-	for thorn in ML_hydro; do \
-		./copy-if-changed.sh $$thorn ../$$thorn; \
+	for thorn in ML_hydro; do				\
+		./copy-if-changed.sh $$thorn ../$$thorn;	\
 	done
 
 clean:
diff --git a/m/McLachlan_ADMConstraints.m b/m/McLachlan_ADMConstraints.m
index 1abb5ae..bb9447b 100644
--- a/m/McLachlan_ADMConstraints.m
+++ b/m/McLachlan_ADMConstraints.m
@@ -138,6 +138,7 @@ ADMConstraintsCalc =
   Name -> ADMConstraints,
   Schedule -> Automatic,
   Where -> Interior,
+  After -> "MoL_PostStep",
   Shorthands -> {detg, gu[ua,ub], G[ua,lb,lc],
                  R[la,lb], trR, Km[la,lb], trK,
                  rho, S[la]},
diff --git a/m/McLachlan_ADMQuantities.m b/m/McLachlan_ADMQuantities.m
index dfa5d22..14131f8 100644
--- a/m/McLachlan_ADMQuantities.m
+++ b/m/McLachlan_ADMQuantities.m
@@ -212,6 +212,7 @@ ADMQuantitiesCalc =
   Name -> ADMQuantities,
   Schedule -> Automatic,
   Where -> Interior,
+  After -> "MoL_PostStep",
   Shorthands -> {detgt, gtu[ua,ub], dgtu[ua,ub,lc],
                  Gtl[la,lb,lc], Gtlu[la,lb,uc], Gt[ua,lb,lc],
                  Xtn[ua], Rt[la,lb], trRt,
diff --git a/m/McLachlan_BSSN.m b/m/McLachlan_BSSN.m
index b05c8aa..927b894 100644
--- a/m/McLachlan_BSSN.m
+++ b/m/McLachlan_BSSN.m
@@ -808,7 +808,6 @@ evolCalc2 = PartialCalculation[evolCalc, "2",
     dot[At[la,lb]]
   }];
 
-
 dissCalc =
 {
   Name -> BSSN <> "_Dissipation",
@@ -826,6 +825,24 @@ dissCalc =
   }
 };
 
+dissCalcs =
+Table[
+{
+  Name -> BSSN <> "_Dissipation_" <> ToString[var /. {Tensor[n_,__] -> n}],
+  Schedule -> {"IN " <> BSSN <> "_evolCalcGroup " <>
+               "AFTER (" <> BSSN <> "_RHS1 " <> BSSN <> "_RHS2)"},
+  ConditionalOnKeyword -> {"apply_dissipation", "always"},
+  Where -> InteriorNoSync,
+  Shorthands -> {epsdiss[ua]},
+  Equations ->
+  {
+    epsdiss[ua] -> EpsDiss,
+    dot[var]    -> dot[var] + epsdiss[ux] PDdiss[var,lx]
+  }
+},
+  {var, {phi, gt[la,lb], Xt[ui], trK, At[la,lb], alpha, A, beta[ua], B[ua]}}
+];
+
 RHSStaticBoundaryCalc =
 {
   Name -> BSSN <> "_RHSStaticBoundary",
@@ -910,11 +927,14 @@ RHSRadiativeBoundaryCalc =
 enforceCalc =
 {
   Name -> BSSN <> "_enforce",
-  Schedule -> {"IN MoL_PostStep BEFORE " <> BSSN <> "_SelectBoundConds"},
+  Schedule -> {"IN MoL_PostStepModify"},
   Shorthands -> {detgt, gtu[ua,ub], trAt},
   Equations -> 
   {
-    (* Enforcing the constraints needs to be a projection, because it
+    (* The following comment is still interesting, but is not correct
+       any more since it is now scheduled in MoL_PostStepModify instead:
+
+       Enforcing the constraints needs to be a projection, because it
        is applied in MoL_PostStep and may thus be applied multiple
        times, not only during time evolution. Therefore detgt has to
        be calculated correctly, without assuming that det gt_ij = 1,
@@ -963,6 +983,7 @@ constraintsCalc =
 {
   Name -> BSSN <> "_constraints",
   Schedule -> Automatic,
+  After -> "MoL_PostStep",
   Where -> Interior,
   Shorthands -> {detgt, ddetgt[la], gtu[ua,ub],
                  Gt[ua,lb,lc], Gtl[la,lb,lc], Gtlu[la,lb,uc], Xtn[ua],
@@ -1293,6 +1314,7 @@ realParameters =
 (******************************************************************************)
 
 calculations =
+Join[
 {
   initialCalc,
   convertFromADMBaseCalc,
@@ -1314,7 +1336,9 @@ calculations =
   convertToADMBaseFakeDtLapseShiftCalc,
   (* constraintsCalc, *)
   constraintsCalc1, constraintsCalc2
-};
+},
+  {} (*dissCalcs*)
+];
 
 CreateKrancThornTT [groups, ".", BSSN,
   Calculations -> calculations,
diff --git a/m/prototype/ML_BSSN_Helper/src/SetGroupTags.c b/m/prototype/ML_BSSN_Helper/src/SetGroupTags.c
index 4b24e56..7ae24bb 100644
--- a/m/prototype/ML_BSSN_Helper/src/SetGroupTags.c
+++ b/m/prototype/ML_BSSN_Helper/src/SetGroupTags.c
@@ -15,29 +15,30 @@ ML_BSSN_SetGroupTags (void)
 {
   DECLARE_CCTK_PARAMETERS;
   
-  set_group_tags (0, 0, 1, "ADMBase::metric");
-  set_group_tags (0, 0, 1, "ADMBase::curv");
-  set_group_tags (0, 0, 1, "ADMBase::lapse");
-  set_group_tags (0, 0, 1, "ADMBase::shift");
-  set_group_tags (0, 0, 1, "ADMBase::dtlapse");
-  set_group_tags (0, 0, 1, "ADMBase::dtshift");
+  int const checkpoint = timelevels > 1;
+  set_group_tags (checkpoint, checkpoint, 1, "ADMBase::metric");
+  set_group_tags (checkpoint, checkpoint, 1, "ADMBase::curv");
+  set_group_tags (checkpoint, checkpoint, 1, "ADMBase::lapse");
+  set_group_tags (checkpoint, checkpoint, 1, "ADMBase::shift");
+  set_group_tags (checkpoint, checkpoint, 1, "ADMBase::dtlapse");
+  set_group_tags (checkpoint, checkpoint, 1, "ADMBase::dtshift");
   
-  set_group_tags (0, 0, 0, "ML_BSSN::ML_cons_detg");
-  set_group_tags (0, 0, 0, "ML_BSSN::ML_cons_Gamma");
-  set_group_tags (0, 0, 0, "ML_BSSN::ML_cons_traceA");
-  set_group_tags (0, 0, 0, "ML_BSSN::ML_Ham");
-  set_group_tags (0, 0, 0, "ML_BSSN::ML_mom");
+  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_cons_detg");
+  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_cons_Gamma");
+  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_cons_traceA");
+  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_Ham");
+  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_mom");
   
-  int const checkpoint = rhs_timelevels > 1;
-  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_log_confacrhs");
-  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_metricrhs");
-  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_Gammarhs");
-  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_trace_curvrhs");
-  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_curvrhs");
-  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_lapserhs");
-  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_dtlapserhs");
-  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_shiftrhs");
-  set_group_tags (checkpoint, checkpoint, 0, "ML_BSSN::ML_dtshiftrhs");
+  int const rhs_checkpoint = rhs_timelevels > 1;
+  set_group_tags (rhs_checkpoint, rhs_checkpoint, 0, "ML_BSSN::ML_log_confacrhs");
+  set_group_tags (rhs_checkpoint, rhs_checkpoint, 0, "ML_BSSN::ML_metricrhs");
+  set_group_tags (rhs_checkpoint, rhs_checkpoint, 0, "ML_BSSN::ML_Gammarhs");
+  set_group_tags (rhs_checkpoint, rhs_checkpoint, 0, "ML_BSSN::ML_trace_curvrhs");
+  set_group_tags (rhs_checkpoint, rhs_checkpoint, 0, "ML_BSSN::ML_curvrhs");
+  set_group_tags (rhs_checkpoint, rhs_checkpoint, 0, "ML_BSSN::ML_lapserhs");
+  set_group_tags (rhs_checkpoint, rhs_checkpoint, 0, "ML_BSSN::ML_dtlapserhs");
+  set_group_tags (rhs_checkpoint, rhs_checkpoint, 0, "ML_BSSN::ML_shiftrhs");
+  set_group_tags (rhs_checkpoint, rhs_checkpoint, 0, "ML_BSSN::ML_dtshiftrhs");
   
   return 0;
 }