From 2c6053350188da6f8b7e63675323157a25a980be Mon Sep 17 00:00:00 2001
From: CGH0S7 <776459475@qq.com>
Date: Sat, 28 Feb 2026 10:57:40 +0800
Subject: [PATCH] =?UTF-8?q?bssn=5Frhs=5Fc:=20=E8=9E=8D=E5=90=88=E9=80=86?=
 =?UTF-8?q?=E5=BA=A6=E8=A7=84+Gamma=E7=BA=A6=E6=9D=9F+Christoffel=E4=B8=89?=
 =?UTF-8?q?=E5=BE=AA=E7=8E=AF=E4=B8=BA=E4=B8=80=E5=BE=AA=E7=8E=AF=20(Modif?=
 =?UTF-8?q?y=202)?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

逆度规计算结果保留在局部标量中直接复用，减少对gupxx..gupzz数组的重复读取，每步加速0.01秒

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
---
 AMSS_NCKU_source/bssn_rhs_c.C | 217 +++++++++++++++-------------------
 1 file changed, 94 insertions(+), 123 deletions(-)

diff --git a/AMSS_NCKU_source/bssn_rhs_c.C b/AMSS_NCKU_source/bssn_rhs_c.C
index 1086026..26f8486 100644
--- a/AMSS_NCKU_source/bssn_rhs_c.C
+++ b/AMSS_NCKU_source/bssn_rhs_c.C
@@ -145,138 +145,109 @@ int f_compute_rhs_bssn(int *ex, double &T,
                         gxx[i] * betaxz[i] + gxy[i] * betayz[i] + gyz[i] * betayx[i] 
                         + gzz[i] * betazx[i] - gxz[i] * betayy[i];                  
         }
-        // 1ms //
+        // Fused: inverse metric + Gamma constraint + Christoffel (3 loops -> 1)
         for(int i=0;i<all;i+=1){
-            double det = gxx[i] * gyy[i] * gzz[i] + gxy[i] * gyz[i] * gxz[i] + gxz[i] * gxy[i] * gyz[i] - 
-                        gxz[i] * gyy[i] * gxz[i] - gxy[i] * gxy[i] * gzz[i] - gxx[i] * gyz[i] * gyz[i];  
-            gupxx[i] = (gyy[i] * gzz[i] - gyz[i] * gyz[i]) / det;
-            gupxy[i] = -(gxy[i] * gzz[i] - gyz[i] * gxz[i]) / det;
-            gupxz[i] = (gxy[i] * gyz[i] - gyy[i] * gxz[i]) / det;
-            gupyy[i] = (gxx[i] * gzz[i] - gxz[i] * gxz[i]) / det;
-            gupyz[i] = -(gxx[i] * gyz[i] - gxy[i] * gxz[i]) / det;
-            gupzz[i] = (gxx[i] * gyy[i] - gxy[i] * gxy[i]) / det;   
-        }
-        // 2.2ms //
-        if(co==0){
-            for (int i=0;i<all;i+=1) {
+            double det = gxx[i] * gyy[i] * gzz[i] + gxy[i] * gyz[i] * gxz[i] + gxz[i] * gxy[i] * gyz[i] -
+                        gxz[i] * gyy[i] * gxz[i] - gxy[i] * gxy[i] * gzz[i] - gxx[i] * gyz[i] * gyz[i];
+            double lg_xx = (gyy[i] * gzz[i] - gyz[i] * gyz[i]) / det;
+            double lg_xy = -(gxy[i] * gzz[i] - gyz[i] * gxz[i]) / det;
+            double lg_xz = (gxy[i] * gyz[i] - gyy[i] * gxz[i]) / det;
+            double lg_yy = (gxx[i] * gzz[i] - gxz[i] * gxz[i]) / det;
+            double lg_yz = -(gxx[i] * gyz[i] - gxy[i] * gxz[i]) / det;
+            double lg_zz = (gxx[i] * gyy[i] - gxy[i] * gxy[i]) / det;
+            gupxx[i] = lg_xx; gupxy[i] = lg_xy; gupxz[i] = lg_xz;
+            gupyy[i] = lg_yy; gupyz[i] = lg_yz; gupzz[i] = lg_zz;
+
+            if(co==0){
                 Gmx_Res[i] = Gamx[i] - (
-                    gupxx[i] * (gupxx[i]*gxxx[i] + gupxy[i]*gxyx[i] + gupxz[i]*gxzx[i]) +
-                    gupxy[i] * (gupxx[i]*gxyx[i] + gupxy[i]*gyyx[i] + gupxz[i]*gyzx[i]) +
-                    gupxz[i] * (gupxx[i]*gxzx[i] + gupxy[i]*gyzx[i] + gupxz[i]*gzzx[i]) +
-
-                    gupxx[i] * (gupxy[i]*gxxy[i] + gupyy[i]*gxyy[i] + gupyz[i]*gxzy[i]) +
-                    gupxy[i] * (gupxy[i]*gxyy[i] + gupyy[i]*gyyy[i] + gupyz[i]*gyzy[i]) +
-                    gupxz[i] * (gupxy[i]*gxzy[i] + gupyy[i]*gyzy[i] + gupyz[i]*gzzy[i]) +
-
-                    gupxx[i] * (gupxz[i]*gxxz[i] + gupyz[i]*gxyz[i] + gupzz[i]*gxzz[i]) +
-                    gupxy[i] * (gupxz[i]*gxyz[i] + gupyz[i]*gyyz[i] + gupzz[i]*gyzz[i]) +
-                    gupxz[i] * (gupxz[i]*gxzz[i] + gupyz[i]*gyzz[i] + gupzz[i]*gzzz[i])
+                    lg_xx * (lg_xx*gxxx[i] + lg_xy*gxyx[i] + lg_xz*gxzx[i]) +
+                    lg_xy * (lg_xx*gxyx[i] + lg_xy*gyyx[i] + lg_xz*gyzx[i]) +
+                    lg_xz * (lg_xx*gxzx[i] + lg_xy*gyzx[i] + lg_xz*gzzx[i]) +
+                    lg_xx * (lg_xy*gxxy[i] + lg_yy*gxyy[i] + lg_yz*gxzy[i]) +
+                    lg_xy * (lg_xy*gxyy[i] + lg_yy*gyyy[i] + lg_yz*gyzy[i]) +
+                    lg_xz * (lg_xy*gxzy[i] + lg_yy*gyzy[i] + lg_yz*gzzy[i]) +
+                    lg_xx * (lg_xz*gxxz[i] + lg_yz*gxyz[i] + lg_zz*gxzz[i]) +
+                    lg_xy * (lg_xz*gxyz[i] + lg_yz*gyyz[i] + lg_zz*gyzz[i]) +
+                    lg_xz * (lg_xz*gxzz[i] + lg_yz*gyzz[i] + lg_zz*gzzz[i])
                 );
-
                 Gmy_Res[i] = Gamy[i] - (
-                    gupxx[i] * (gupxy[i]*gxxx[i] + gupyy[i]*gxyx[i] + gupyz[i]*gxzx[i]) +
-                    gupxy[i] * (gupxy[i]*gxyx[i] + gupyy[i]*gyyx[i] + gupyz[i]*gyzx[i]) +
-                    gupxz[i] * (gupxy[i]*gxzx[i] + gupyy[i]*gyzx[i] + gupyz[i]*gzzx[i]) +
-
-                    gupxy[i] * (gupxy[i]*gxxy[i] + gupyy[i]*gxyy[i] + gupyz[i]*gxzy[i]) +
-                    gupyy[i] * (gupxy[i]*gxyy[i] + gupyy[i]*gyyy[i] + gupyz[i]*gyzy[i]) +
-                    gupyz[i] * (gupxy[i]*gxzy[i] + gupyy[i]*gyzy[i] + gupyz[i]*gzzy[i]) +
-
-                    gupxy[i] * (gupxz[i]*gxxz[i] + gupyz[i]*gxyz[i] + gupzz[i]*gxzz[i]) +
-                    gupyy[i] * (gupxz[i]*gxyz[i] + gupyz[i]*gyyz[i] + gupzz[i]*gyzz[i]) +
-                    gupyz[i] * (gupxz[i]*gxzz[i] + gupyz[i]*gyzz[i] + gupzz[i]*gzzz[i])
+                    lg_xx * (lg_xy*gxxx[i] + lg_yy*gxyx[i] + lg_yz*gxzx[i]) +
+                    lg_xy * (lg_xy*gxyx[i] + lg_yy*gyyx[i] + lg_yz*gyzx[i]) +
+                    lg_xz * (lg_xy*gxzx[i] + lg_yy*gyzx[i] + lg_yz*gzzx[i]) +
+                    lg_xy * (lg_xy*gxxy[i] + lg_yy*gxyy[i] + lg_yz*gxzy[i]) +
+                    lg_yy * (lg_xy*gxyy[i] + lg_yy*gyyy[i] + lg_yz*gyzy[i]) +
+                    lg_yz * (lg_xy*gxzy[i] + lg_yy*gyzy[i] + lg_yz*gzzy[i]) +
+                    lg_xy * (lg_xz*gxxz[i] + lg_yz*gxyz[i] + lg_zz*gxzz[i]) +
+                    lg_yy * (lg_xz*gxyz[i] + lg_yz*gyyz[i] + lg_zz*gyzz[i]) +
+                    lg_yz * (lg_xz*gxzz[i] + lg_yz*gyzz[i] + lg_zz*gzzz[i])
                 );
-
                 Gmz_Res[i] = Gamz[i] - (
-                    gupxx[i] * (gupxz[i]*gxxx[i] + gupyz[i]*gxyx[i] + gupzz[i]*gxzx[i]) +
-                    gupxy[i] * (gupxz[i]*gxyx[i] + gupyz[i]*gyyx[i] + gupzz[i]*gyzx[i]) +
-                    gupxz[i] * (gupxz[i]*gxzx[i] + gupyz[i]*gyzx[i] + gupzz[i]*gzzx[i]) +
-
-                    gupxy[i] * (gupxz[i]*gxxy[i] + gupyz[i]*gxyy[i] + gupzz[i]*gxzy[i]) +
-                    gupyy[i] * (gupxz[i]*gxyy[i] + gupyz[i]*gyyy[i] + gupzz[i]*gyzy[i]) +
-                    gupyz[i] * (gupxz[i]*gxzy[i] + gupyz[i]*gyzy[i] + gupzz[i]*gzzy[i]) +
-
-                    gupxz[i] * (gupxz[i]*gxxz[i] + gupyz[i]*gxyz[i] + gupzz[i]*gxzz[i]) +
-                    gupyz[i] * (gupxz[i]*gxyz[i] + gupyz[i]*gyyz[i] + gupzz[i]*gyzz[i]) +
-                    gupzz[i] * (gupxz[i]*gxzz[i] + gupyz[i]*gyzz[i] + gupzz[i]*gzzz[i])
+                    lg_xx * (lg_xz*gxxx[i] + lg_yz*gxyx[i] + lg_zz*gxzx[i]) +
+                    lg_xy * (lg_xz*gxyx[i] + lg_yz*gyyx[i] + lg_zz*gyzx[i]) +
+                    lg_xz * (lg_xz*gxzx[i] + lg_yz*gyzx[i] + lg_zz*gzzx[i]) +
+                    lg_xy * (lg_xz*gxxy[i] + lg_yz*gxyy[i] + lg_zz*gxzy[i]) +
+                    lg_yy * (lg_xz*gxyy[i] + lg_yz*gyyy[i] + lg_zz*gyzy[i]) +
+                    lg_yz * (lg_xz*gxzy[i] + lg_yz*gyzy[i] + lg_zz*gzzy[i]) +
+                    lg_xz * (lg_xz*gxxz[i] + lg_yz*gxyz[i] + lg_zz*gxzz[i]) +
+                    lg_yz * (lg_xz*gxyz[i] + lg_yz*gyyz[i] + lg_zz*gyzz[i]) +
+                    lg_zz * (lg_xz*gxzz[i] + lg_yz*gyzz[i] + lg_zz*gzzz[i])
                 );
             }
-        }
-        // 5ms //
-        for (int i=0;i<all;i+=1) {
-
-            Gamxxx[i] = HALF * ( gupxx[i]*gxxx[i]
-                            + gupxy[i]*(TWO*gxyx[i] - gxxy[i])
-                            + gupxz[i]*(TWO*gxzx[i] - gxxz[i]) );
-
-            Gamyxx[i] = HALF * ( gupxy[i]*gxxx[i]
-                            + gupyy[i]*(TWO*gxyx[i] - gxxy[i])
-                            + gupyz[i]*(TWO*gxzx[i] - gxxz[i]) );
-
-            Gamzxx[i] = HALF * ( gupxz[i]*gxxx[i]
-                            + gupyz[i]*(TWO*gxyx[i] - gxxy[i])
-                            + gupzz[i]*(TWO*gxzx[i] - gxxz[i]) );
-
-            Gamxyy[i] = HALF * ( gupxx[i]*(TWO*gxyy[i] - gyyx[i])
-                            + gupxy[i]*gyyy[i]
-                            + gupxz[i]*(TWO*gyzy[i] - gyyz[i]) );
-
-            Gamyyy[i] = HALF * ( gupxy[i]*(TWO*gxyy[i] - gyyx[i])
-                            + gupyy[i]*gyyy[i]
-                            + gupyz[i]*(TWO*gyzy[i] - gyyz[i]) );
-
-            Gamzyy[i] = HALF * ( gupxz[i]*(TWO*gxyy[i] - gyyx[i])
-                            + gupyz[i]*gyyy[i]
-                            + gupzz[i]*(TWO*gyzy[i] - gyyz[i]) );
-
-            Gamxzz[i] = HALF * ( gupxx[i]*(TWO*gxzz[i] - gzzx[i])
-                            + gupxy[i]*(TWO*gyzz[i] - gzzy[i])
-                            + gupxz[i]*gzzz[i] );
-
-            Gamyzz[i] = HALF * ( gupxy[i]*(TWO*gxzz[i] - gzzx[i])
-                            + gupyy[i]*(TWO*gyzz[i] - gzzy[i])
-                            + gupyz[i]*gzzz[i] );
-
-            Gamzzz[i] = HALF * ( gupxz[i]*(TWO*gxzz[i] - gzzx[i])
-                            + gupyz[i]*(TWO*gyzz[i] - gzzy[i])
-                            + gupzz[i]*gzzz[i] );
-
-            Gamxxy[i] = HALF * ( gupxx[i]*gxxy[i]
-                            + gupxy[i]*gyyx[i]
-                            + gupxz[i]*(gxzy[i] + gyzx[i] - gxyz[i]) );
-
-            Gamyxy[i] = HALF * ( gupxy[i]*gxxy[i]
-                            + gupyy[i]*gyyx[i]
-                            + gupyz[i]*(gxzy[i] + gyzx[i] - gxyz[i]) );
-
-            Gamzxy[i] = HALF * ( gupxz[i]*gxxy[i]
-                            + gupyz[i]*gyyx[i]
-                            + gupzz[i]*(gxzy[i] + gyzx[i] - gxyz[i]) );
-
-            Gamxxz[i] = HALF * ( gupxx[i]*gxxz[i]
-                            + gupxy[i]*(gxyz[i] + gyzx[i] - gxzy[i])
-                            + gupxz[i]*gzzx[i] );
-
-            Gamyxz[i] = HALF * ( gupxy[i]*gxxz[i]
-                            + gupyy[i]*(gxyz[i] + gyzx[i] - gxzy[i])
-                            + gupyz[i]*gzzx[i] );
-
-            Gamzxz[i] = HALF * ( gupxz[i]*gxxz[i]
-                            + gupyz[i]*(gxyz[i] + gyzx[i] - gxzy[i])
-                            + gupzz[i]*gzzx[i] );
-
-            Gamxyz[i] = HALF * ( gupxx[i]*(gxyz[i] + gxzy[i] - gyzx[i])
-                            + gupxy[i]*gyyz[i]
-                            + gupxz[i]*gzzy[i] );
-
-            Gamyyz[i] = HALF * ( gupxy[i]*(gxyz[i] + gxzy[i] - gyzx[i])
-                            + gupyy[i]*gyyz[i]
-                            + gupyz[i]*gzzy[i] );
-
-            Gamzyz[i] = HALF * ( gupxz[i]*(gxyz[i] + gxzy[i] - gyzx[i])
-                            + gupyz[i]*gyyz[i]
-                            + gupzz[i]*gzzy[i] );
 
+            Gamxxx[i] = HALF * ( lg_xx*gxxx[i]
+                            + lg_xy*(TWO*gxyx[i] - gxxy[i])
+                            + lg_xz*(TWO*gxzx[i] - gxxz[i]) );
+            Gamyxx[i] = HALF * ( lg_xy*gxxx[i]
+                            + lg_yy*(TWO*gxyx[i] - gxxy[i])
+                            + lg_yz*(TWO*gxzx[i] - gxxz[i]) );
+            Gamzxx[i] = HALF * ( lg_xz*gxxx[i]
+                            + lg_yz*(TWO*gxyx[i] - gxxy[i])
+                            + lg_zz*(TWO*gxzx[i] - gxxz[i]) );
+            Gamxyy[i] = HALF * ( lg_xx*(TWO*gxyy[i] - gyyx[i])
+                            + lg_xy*gyyy[i]
+                            + lg_xz*(TWO*gyzy[i] - gyyz[i]) );
+            Gamyyy[i] = HALF * ( lg_xy*(TWO*gxyy[i] - gyyx[i])
+                            + lg_yy*gyyy[i]
+                            + lg_yz*(TWO*gyzy[i] - gyyz[i]) );
+            Gamzyy[i] = HALF * ( lg_xz*(TWO*gxyy[i] - gyyx[i])
+                            + lg_yz*gyyy[i]
+                            + lg_zz*(TWO*gyzy[i] - gyyz[i]) );
+            Gamxzz[i] = HALF * ( lg_xx*(TWO*gxzz[i] - gzzx[i])
+                            + lg_xy*(TWO*gyzz[i] - gzzy[i])
+                            + lg_xz*gzzz[i] );
+            Gamyzz[i] = HALF * ( lg_xy*(TWO*gxzz[i] - gzzx[i])
+                            + lg_yy*(TWO*gyzz[i] - gzzy[i])
+                            + lg_yz*gzzz[i] );
+            Gamzzz[i] = HALF * ( lg_xz*(TWO*gxzz[i] - gzzx[i])
+                            + lg_yz*(TWO*gyzz[i] - gzzy[i])
+                            + lg_zz*gzzz[i] );
+            Gamxxy[i] = HALF * ( lg_xx*gxxy[i]
+                            + lg_xy*gyyx[i]
+                            + lg_xz*(gxzy[i] + gyzx[i] - gxyz[i]) );
+            Gamyxy[i] = HALF * ( lg_xy*gxxy[i]
+                            + lg_yy*gyyx[i]
+                            + lg_yz*(gxzy[i] + gyzx[i] - gxyz[i]) );
+            Gamzxy[i] = HALF * ( lg_xz*gxxy[i]
+                            + lg_yz*gyyx[i]
+                            + lg_zz*(gxzy[i] + gyzx[i] - gxyz[i]) );
+            Gamxxz[i] = HALF * ( lg_xx*gxxz[i]
+                            + lg_xy*(gxyz[i] + gyzx[i] - gxzy[i])
+                            + lg_xz*gzzx[i] );
+            Gamyxz[i] = HALF * ( lg_xy*gxxz[i]
+                            + lg_yy*(gxyz[i] + gyzx[i] - gxzy[i])
+                            + lg_yz*gzzx[i] );
+            Gamzxz[i] = HALF * ( lg_xz*gxxz[i]
+                            + lg_yz*(gxyz[i] + gyzx[i] - gxzy[i])
+                            + lg_zz*gzzx[i] );
+            Gamxyz[i] = HALF * ( lg_xx*(gxyz[i] + gxzy[i] - gyzx[i])
+                            + lg_xy*gyyz[i]
+                            + lg_xz*gzzy[i] );
+            Gamyyz[i] = HALF * ( lg_xy*(gxyz[i] + gxzy[i] - gyzx[i])
+                            + lg_yy*gyyz[i]
+                            + lg_yz*gzzy[i] );
+            Gamzyz[i] = HALF * ( lg_xz*(gxyz[i] + gxzy[i] - gyzx[i])
+                            + lg_yz*gyyz[i]
+                            + lg_zz*gzzy[i] );
         }
         // 1.8ms //
         for (int i=0;i<all;i+=1) {