Optimize constraint refresh after regrid

(cherry picked from commit 5c65cea2f0)

AMSS_NCKU_source/bssn_class.C

@@ -372,6 +372,7 @@ bssn_class::bssn_class(double Couranti, double StartTimei, double TotalTimei,
                        : Courant(Couranti), StartTime(StartTimei), TotalTime(TotalTimei), 
                          DumpTime(DumpTimei), d2DumpTime(d2DumpTimei), CheckTime(CheckTimei), AnasTime(AnasTimei),
                          cuda_level0_constraint_cache_valid(false),
+                         ConstraintRefreshLevels(0),
                          Symmetry(Symmetryi), checkrun(checkruni), numepss(numepssi), numepsb(numepsbi), numepsh(numepshi),
 #ifdef With_AHF
                          xc(0), yc(0), zc(0), xr(0), yr(0), zr(0), trigger(0), dTT(0), dumpid(0),
@@ -1014,6 +1015,9 @@ void bssn_class::Initialize()
     }
   }
   GH = new cgh(0, ngfs, Symmetry, pname, checkrun, ErrorMonitor);
+  ConstraintRefreshLevels = new int[GH->levels];
+  for (int il = 0; il < GH->levels; il++)
+    ConstraintRefreshLevels[il] = 0;
   if (checkrun)
     CheckPoint->readcheck_cgh(PhysTime, GH, myrank, nprocs, Symmetry);
   else
@@ -1103,6 +1107,8 @@ bssn_class::~bssn_class()
   DumpList->clearList();
   ConstraintList->clearList();
 
+  delete[] ConstraintRefreshLevels;
+
   delete phio;
   delete trKo;
   delete gxxo;
@@ -2755,7 +2761,11 @@ void bssn_class::RecursiveStep(int lev)
   if (GH->Regrid_Onelevel(lev, Symmetry, BH_num, Porgbr, Porg0,
                       SynchList_cor, OldStateList, StateList, SynchList_pre,
                       fgt(PhysTime - dT_lev, StartTime, dT_lev / 2), ErrorMonitor))
+  {
+  if (ConstraintRefreshLevels)
+    ConstraintRefreshLevels[lev] = 1;
   for (int il = 0; il < GH->levels; il++) { sync_cache_pre[il].invalidate(); sync_cache_cor[il].invalidate(); sync_cache_rp_coarse[il].invalidate(); sync_cache_rp_fine[il].invalidate(); sync_cache_restrict[il].invalidate(); sync_cache_outbd[il].invalidate(); }
+  }
 #endif
 }
 
@@ -7765,7 +7775,8 @@ void bssn_class::Constraint_Out()
     for (int lev = 0; lev < GH->levels; lev++)
     {
       // make sure the data consistent for higher levels
-      if (bssn_constraint_recompute_from_state(lev, cuda_level0_constraint_cache_valid))
+      if ((lev == 0 && bssn_constraint_recompute_from_state(lev, cuda_level0_constraint_cache_valid)) ||
+          (lev > 0 && ConstraintRefreshLevels && ConstraintRefreshLevels[lev]))
       {
         double TRK4 = PhysTime;
         double ndeps = numepsb;
@@ -7978,6 +7989,9 @@ void bssn_class::Constraint_Out()
     Interp_Constraint(false);
 
     LastConsOut = 0;
+    if (ConstraintRefreshLevels)
+      for (int lev = 0; lev < GH->levels; lev++)
+        ConstraintRefreshLevels[lev] = 0;
   }
 }
 

AMSS_NCKU_source/bssn_class.h

@@ -49,6 +49,7 @@ public:
        double AnasTime, DumpTime, d2DumpTime, CheckTime;
        double LastAnas, LastConsOut;
        bool cuda_level0_constraint_cache_valid;
+       int *ConstraintRefreshLevels;
        double Courant;
        double numepss, numepsb, numepsh;
        int Symmetry;