From 5c65cea2f03a90058f297a1d661b8ab5152bcc2f Mon Sep 17 00:00:00 2001
From: CGH0S7 <776459475@qq.com>
Date: Mon, 13 Apr 2026 11:39:50 +0800
Subject: [PATCH] Optimize constraint refresh after regrid

---
 AMSS_NCKU_source/bssn_class.C | 57 +++++++++++++++++++++--------------
 AMSS_NCKU_source/bssn_class.h |  9 +++---
 2 files changed, 40 insertions(+), 26 deletions(-)

diff --git a/AMSS_NCKU_source/bssn_class.C b/AMSS_NCKU_source/bssn_class.C
index bcd9117..6a12b14 100644
--- a/AMSS_NCKU_source/bssn_class.C
+++ b/AMSS_NCKU_source/bssn_class.C
@@ -64,8 +64,9 @@ bssn_class::bssn_class(double Couranti, double StartTimei, double TotalTimei,
 #ifdef With_AHF
                          xc(0), yc(0), zc(0), xr(0), yr(0), zr(0), trigger(0), dTT(0), dumpid(0),
 #endif
-                         a_lev(a_levi), maxl(maxli), decn(decni), maxrex(maxrexi), drex(drexi),
-                         CheckPoint(0)
+                         a_lev(a_levi), maxl(maxli), decn(decni), maxrex(maxrexi), drex(drexi),
+                         ConstraintRefreshLevels(0),
+                         CheckPoint(0)
                          // CheckPoint(0)
 {
   MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
@@ -701,9 +702,12 @@ void bssn_class::Initialize()
       exit(0);
     }
   }
-  GH = new cgh(0, ngfs, Symmetry, pname, checkrun, ErrorMonitor);
-  if (checkrun)
-    CheckPoint->readcheck_cgh(PhysTime, GH, myrank, nprocs, Symmetry);
+  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
     GH->compose_cgh(nprocs);
 #ifdef WithShell
@@ -789,9 +793,11 @@ bssn_class::~bssn_class()
   SynchList_pre->clearList();
   SynchList_cor->clearList();
   DumpList->clearList();
-  ConstraintList->clearList();
-
-  delete phio;
+  ConstraintList->clearList();
+
+  delete[] ConstraintRefreshLevels;
+
+  delete phio;
   delete trKo;
   delete gxxo;
   delete gxyo;
@@ -2437,13 +2443,17 @@ void bssn_class::RecursiveStep(int lev)
     if(lev>0) Parallel::Restrict_after(GH->PatL[lev-1],GH->PatL[lev],StateList,StateList,Symmetry);
 #endif
 
-#if (REGLEV == 0)
-  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))
-  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
-}
+#if (REGLEV == 0)
+  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
+}
 
 //================================================================================================
 
@@ -7261,8 +7271,8 @@ void bssn_class::Constraint_Out()
     for (int lev = 0; lev < GH->levels; lev++)
     {
       // make sure the data consistent for higher levels
-      if (lev > 0) // if the constrait quantities can be reused from the step rhs calculation
-      {
+      if (lev > 0 && ConstraintRefreshLevels && ConstraintRefreshLevels[lev]) // only refresh levels whose grid layout changed after evolution
+      {
         double TRK4 = PhysTime;
         double ndeps = numepsb;
         int pre = 0;
@@ -7471,11 +7481,14 @@ void bssn_class::Constraint_Out()
       */
     }
 
-    Interp_Constraint(false);
-
-    LastConsOut = 0;
-  }
-}
+    Interp_Constraint(false);
+
+    LastConsOut = 0;
+    if (ConstraintRefreshLevels)
+      for (int lev = 0; lev < GH->levels; lev++)
+        ConstraintRefreshLevels[lev] = 0;
+  }
+}
 
 //================================================================================================
 
diff --git a/AMSS_NCKU_source/bssn_class.h b/AMSS_NCKU_source/bssn_class.h
index 5a8eb2a..c067902 100644
--- a/AMSS_NCKU_source/bssn_class.h
+++ b/AMSS_NCKU_source/bssn_class.h
@@ -45,10 +45,11 @@ public:
        int checkrun;
        char checkfilename[50];
        int Steps;
-       double StartTime, TotalTime;
-       double AnasTime, DumpTime, d2DumpTime, CheckTime;
-       double LastAnas, LastConsOut;
-       double Courant;
+       double StartTime, TotalTime;
+       double AnasTime, DumpTime, d2DumpTime, CheckTime;
+       double LastAnas, LastConsOut;
+       int *ConstraintRefreshLevels;
+       double Courant;
        double numepss, numepsb, numepsh;
        int Symmetry;
        int maxl, decn;