Treate Tmunu variables as tensors, instead of accessing them via local helper scalars.

7 files changed, 38 insertions, 208 deletions

diff --git a/ML_ADMConstraints/src/ML_ADMConstraints.c b/ML_ADMConstraints/src/ML_ADMConstraints.c
index 90a8989..30e2457 100644
--- a/ML_ADMConstraints/src/ML_ADMConstraints.c
+++ b/ML_ADMConstraints/src/ML_ADMConstraints.c
@@ -109,8 +109,6 @@ void ML_ADMConstraints_Body(cGH const * restrict const cctkGH, int const dir, in
     // CCTK_REAL R11 = INITVALUE, R12 = INITVALUE, R13 = INITVALUE, R22 = INITVALUE, R23 = INITVALUE, R33 = INITVALUE;
     // CCTK_REAL rho = INITVALUE;
     // CCTK_REAL S1 = INITVALUE, S2 = INITVALUE, S3 = INITVALUE;
-    // CCTK_REAL T00 = INITVALUE, T01 = INITVALUE, T02 = INITVALUE, T03 = INITVALUE, T11 = INITVALUE, T12 = INITVALUE;
-    // CCTK_REAL T13 = INITVALUE, T22 = INITVALUE, T23 = INITVALUE, T33 = INITVALUE;
     // CCTK_REAL trK = INITVALUE;
     // CCTK_REAL trR = INITVALUE;
     
@@ -419,35 +417,15 @@ void ML_ADMConstraints_Body(cGH const * restrict const cctkGH, int const dir, in
     
     CCTK_REAL const trK  =  Km11 + Km22 + Km33;
     
-    CCTK_REAL const T00  =  eTttL;
+    CCTK_REAL const rho  =  pow(alpL,-2)*(eTttL - 2*(betayL*eTtyL + betazL*eTtzL) + 
+          2*(betaxL*(-eTtxL + betayL*eTxyL + betazL*eTxzL) + betayL*betazL*eTyzL) + eTxxL*SQR(betaxL) + eTyyL*SQR(betayL) + 
+          eTzzL*SQR(betazL));
     
-    CCTK_REAL const T01  =  eTtxL;
+    CCTK_REAL const S1  =  (-eTtxL + betaxL*eTxxL + betayL*eTxyL + betazL*eTxzL)*INV(alpL);
     
-    CCTK_REAL const T02  =  eTtyL;
+    CCTK_REAL const S2  =  (-eTtyL + betaxL*eTxyL + betayL*eTyyL + betazL*eTyzL)*INV(alpL);
     
-    CCTK_REAL const T03  =  eTtzL;
-    
-    CCTK_REAL const T11  =  eTxxL;
-    
-    CCTK_REAL const T12  =  eTxyL;
-    
-    CCTK_REAL const T13  =  eTxzL;
-    
-    CCTK_REAL const T22  =  eTyyL;
-    
-    CCTK_REAL const T23  =  eTyzL;
-    
-    CCTK_REAL const T33  =  eTzzL;
-    
-    CCTK_REAL const rho  =  pow(alpL,-2)*(T00 - 2*(betayL*T02 + betazL*T03) + 
-          2*(betaxL*(-T01 + betayL*T12 + betazL*T13) + betayL*betazL*T23) + T11*SQR(betaxL) + T22*SQR(betayL) + 
-          T33*SQR(betazL));
-    
-    CCTK_REAL const S1  =  (-T01 + betaxL*T11 + betayL*T12 + betazL*T13)*INV(alpL);
-    
-    CCTK_REAL const S2  =  (-T02 + betaxL*T12 + betayL*T22 + betazL*T23)*INV(alpL);
-    
-    CCTK_REAL const S3  =  (-T03 + betaxL*T13 + betayL*T23 + betazL*T33)*INV(alpL);
+    CCTK_REAL const S3  =  (-eTtzL + betaxL*eTxzL + betayL*eTyzL + betazL*eTzzL)*INV(alpL);
     
     CCTK_REAL const HL  =  -2*(Km12*Km21 + Km13*Km31 + Km23*Km32) - 50.26548245743669181540229413247204614715*rho + trR - SQR(Km11) - 
         SQR(Km22) - SQR(Km33) + SQR(trK);
diff --git a/ML_ADMConstraints_MP/src/ML_ADMConstraints_MP.c b/ML_ADMConstraints_MP/src/ML_ADMConstraints_MP.c
index 07f750f..4f6597f 100644
--- a/ML_ADMConstraints_MP/src/ML_ADMConstraints_MP.c
+++ b/ML_ADMConstraints_MP/src/ML_ADMConstraints_MP.c
@@ -109,8 +109,6 @@ void ML_ADMConstraints_MP_Body(cGH const * restrict const cctkGH, int const dir,
     // CCTK_REAL R11 = INITVALUE, R12 = INITVALUE, R13 = INITVALUE, R22 = INITVALUE, R23 = INITVALUE, R33 = INITVALUE;
     // CCTK_REAL rho = INITVALUE;
     // CCTK_REAL S1 = INITVALUE, S2 = INITVALUE, S3 = INITVALUE;
-    // CCTK_REAL T00 = INITVALUE, T01 = INITVALUE, T02 = INITVALUE, T03 = INITVALUE, T11 = INITVALUE, T12 = INITVALUE;
-    // CCTK_REAL T13 = INITVALUE, T22 = INITVALUE, T23 = INITVALUE, T33 = INITVALUE;
     // CCTK_REAL trK = INITVALUE;
     // CCTK_REAL trR = INITVALUE;
     
@@ -806,35 +804,15 @@ void ML_ADMConstraints_MP_Body(cGH const * restrict const cctkGH, int const dir,
     
     CCTK_REAL const trK  =  Km11 + Km22 + Km33;
     
-    CCTK_REAL const T00  =  eTttL;
+    CCTK_REAL const rho  =  pow(alpL,-2)*(eTttL - 2*(betayL*eTtyL + betazL*eTtzL) + 
+          2*(betaxL*(-eTtxL + betayL*eTxyL + betazL*eTxzL) + betayL*betazL*eTyzL) + eTxxL*SQR(betaxL) + eTyyL*SQR(betayL) + 
+          eTzzL*SQR(betazL));
     
-    CCTK_REAL const T01  =  eTtxL;
+    CCTK_REAL const S1  =  (-eTtxL + betaxL*eTxxL + betayL*eTxyL + betazL*eTxzL)*INV(alpL);
     
-    CCTK_REAL const T02  =  eTtyL;
+    CCTK_REAL const S2  =  (-eTtyL + betaxL*eTxyL + betayL*eTyyL + betazL*eTyzL)*INV(alpL);
     
-    CCTK_REAL const T03  =  eTtzL;
-    
-    CCTK_REAL const T11  =  eTxxL;
-    
-    CCTK_REAL const T12  =  eTxyL;
-    
-    CCTK_REAL const T13  =  eTxzL;
-    
-    CCTK_REAL const T22  =  eTyyL;
-    
-    CCTK_REAL const T23  =  eTyzL;
-    
-    CCTK_REAL const T33  =  eTzzL;
-    
-    CCTK_REAL const rho  =  pow(alpL,-2)*(T00 - 2*(betayL*T02 + betazL*T03) + 
-          2*(betaxL*(-T01 + betayL*T12 + betazL*T13) + betayL*betazL*T23) + T11*SQR(betaxL) + T22*SQR(betayL) + 
-          T33*SQR(betazL));
-    
-    CCTK_REAL const S1  =  (-T01 + betaxL*T11 + betayL*T12 + betazL*T13)*INV(alpL);
-    
-    CCTK_REAL const S2  =  (-T02 + betaxL*T12 + betayL*T22 + betazL*T23)*INV(alpL);
-    
-    CCTK_REAL const S3  =  (-T03 + betaxL*T13 + betayL*T23 + betazL*T33)*INV(alpL);
+    CCTK_REAL const S3  =  (-eTtzL + betaxL*eTxzL + betayL*eTyzL + betazL*eTzzL)*INV(alpL);
     
     CCTK_REAL const HL  =  -2*(Km12*Km21 + Km13*Km31 + Km23*Km32) - 50.26548245743669181540229413247204614715*rho + trR - SQR(Km11) - 
         SQR(Km22) - SQR(Km33) + SQR(trK);
diff --git a/ML_ADMQuantities/src/ML_ADMQuantities.c b/ML_ADMQuantities/src/ML_ADMQuantities.c
index 51dd888..9253921 100644
--- a/ML_ADMQuantities/src/ML_ADMQuantities.c
+++ b/ML_ADMQuantities/src/ML_ADMQuantities.c
@@ -121,8 +121,6 @@ void ML_ADMQuantities_Body(cGH const * restrict const cctkGH, int const dir, int
     // CCTK_REAL rho = INITVALUE;
     // CCTK_REAL Rt11 = INITVALUE, Rt12 = INITVALUE, Rt13 = INITVALUE, Rt22 = INITVALUE, Rt23 = INITVALUE, Rt33 = INITVALUE;
     // CCTK_REAL S1 = INITVALUE, S2 = INITVALUE, S3 = INITVALUE;
-    // CCTK_REAL T00 = INITVALUE, T01 = INITVALUE, T02 = INITVALUE, T03 = INITVALUE, T11 = INITVALUE, T12 = INITVALUE;
-    // CCTK_REAL T13 = INITVALUE, T22 = INITVALUE, T23 = INITVALUE, T33 = INITVALUE;
     // CCTK_REAL trRt = INITVALUE;
     // CCTK_REAL Xtn1 = INITVALUE, Xtn2 = INITVALUE, Xtn3 = INITVALUE;
     
@@ -611,35 +609,15 @@ void ML_ADMQuantities_Body(cGH const * restrict const cctkGH, int const dir, int
     
     CCTK_REAL const Atm33  =  At13L*gtu31 + At23L*gtu32 + At33L*gtu33;
     
-    CCTK_REAL const T00  =  eTttL;
+    CCTK_REAL const rho  =  pow(alphaL,-2)*(eTttL - 2*(beta2L*eTtyL + beta3L*eTtzL) + 
+          2*(beta1L*(-eTtxL + beta2L*eTxyL + beta3L*eTxzL) + beta2L*beta3L*eTyzL) + eTxxL*SQR(beta1L) + eTyyL*SQR(beta2L) + 
+          eTzzL*SQR(beta3L));
     
-    CCTK_REAL const T01  =  eTtxL;
+    CCTK_REAL const S1  =  (-eTtxL + beta1L*eTxxL + beta2L*eTxyL + beta3L*eTxzL)*INV(alphaL);
     
-    CCTK_REAL const T02  =  eTtyL;
+    CCTK_REAL const S2  =  (-eTtyL + beta1L*eTxyL + beta2L*eTyyL + beta3L*eTyzL)*INV(alphaL);
     
-    CCTK_REAL const T03  =  eTtzL;
-    
-    CCTK_REAL const T11  =  eTxxL;
-    
-    CCTK_REAL const T12  =  eTxyL;
-    
-    CCTK_REAL const T13  =  eTxzL;
-    
-    CCTK_REAL const T22  =  eTyyL;
-    
-    CCTK_REAL const T23  =  eTyzL;
-    
-    CCTK_REAL const T33  =  eTzzL;
-    
-    CCTK_REAL const rho  =  pow(alphaL,-2)*(T00 - 2*(beta2L*T02 + beta3L*T03) + 
-          2*(beta1L*(-T01 + beta2L*T12 + beta3L*T13) + beta2L*beta3L*T23) + T11*SQR(beta1L) + T22*SQR(beta2L) + 
-          T33*SQR(beta3L));
-    
-    CCTK_REAL const S1  =  (-T01 + beta1L*T11 + beta2L*T12 + beta3L*T13)*INV(alphaL);
-    
-    CCTK_REAL const S2  =  (-T02 + beta1L*T12 + beta2L*T22 + beta3L*T23)*INV(alphaL);
-    
-    CCTK_REAL const S3  =  (-T03 + beta1L*T13 + beta2L*T23 + beta3L*T33)*INV(alphaL);
+    CCTK_REAL const S3  =  (-eTtzL + beta1L*eTxzL + beta2L*eTyzL + beta3L*eTzzL)*INV(alphaL);
     
     CCTK_REAL const MadmL  =  0.01989436788648691697111047042156429525431*
         (-((Gt111*Gtlu111 + Gt112*Gtlu112 + Gt113*Gtlu113 + Gt211*Gtlu121 + Gt212*Gtlu122 + Gt213*Gtlu123 + Gt311*Gtlu131 + 
diff --git a/ML_ADMQuantities_MP/src/ML_ADMQuantities_MP.c b/ML_ADMQuantities_MP/src/ML_ADMQuantities_MP.c
index d84a21f..948c148 100644
--- a/ML_ADMQuantities_MP/src/ML_ADMQuantities_MP.c
+++ b/ML_ADMQuantities_MP/src/ML_ADMQuantities_MP.c
@@ -121,8 +121,6 @@ void ML_ADMQuantities_MP_Body(cGH const * restrict const cctkGH, int const dir,
     // CCTK_REAL rho = INITVALUE;
     // CCTK_REAL Rt11 = INITVALUE, Rt12 = INITVALUE, Rt13 = INITVALUE, Rt22 = INITVALUE, Rt23 = INITVALUE, Rt33 = INITVALUE;
     // CCTK_REAL S1 = INITVALUE, S2 = INITVALUE, S3 = INITVALUE;
-    // CCTK_REAL T00 = INITVALUE, T01 = INITVALUE, T02 = INITVALUE, T03 = INITVALUE, T11 = INITVALUE, T12 = INITVALUE;
-    // CCTK_REAL T13 = INITVALUE, T22 = INITVALUE, T23 = INITVALUE, T33 = INITVALUE;
     // CCTK_REAL trRt = INITVALUE;
     // CCTK_REAL Xtn1 = INITVALUE, Xtn2 = INITVALUE, Xtn3 = INITVALUE;
     
@@ -896,35 +894,15 @@ void ML_ADMQuantities_MP_Body(cGH const * restrict const cctkGH, int const dir,
     
     CCTK_REAL const Atm33  =  At13L*gtu31 + At23L*gtu32 + At33L*gtu33;
     
-    CCTK_REAL const T00  =  eTttL;
+    CCTK_REAL const rho  =  pow(alphaL,-2)*(eTttL - 2*(beta2L*eTtyL + beta3L*eTtzL) + 
+          2*(beta1L*(-eTtxL + beta2L*eTxyL + beta3L*eTxzL) + beta2L*beta3L*eTyzL) + eTxxL*SQR(beta1L) + eTyyL*SQR(beta2L) + 
+          eTzzL*SQR(beta3L));
     
-    CCTK_REAL const T01  =  eTtxL;
+    CCTK_REAL const S1  =  (-eTtxL + beta1L*eTxxL + beta2L*eTxyL + beta3L*eTxzL)*INV(alphaL);
     
-    CCTK_REAL const T02  =  eTtyL;
+    CCTK_REAL const S2  =  (-eTtyL + beta1L*eTxyL + beta2L*eTyyL + beta3L*eTyzL)*INV(alphaL);
     
-    CCTK_REAL const T03  =  eTtzL;
-    
-    CCTK_REAL const T11  =  eTxxL;
-    
-    CCTK_REAL const T12  =  eTxyL;
-    
-    CCTK_REAL const T13  =  eTxzL;
-    
-    CCTK_REAL const T22  =  eTyyL;
-    
-    CCTK_REAL const T23  =  eTyzL;
-    
-    CCTK_REAL const T33  =  eTzzL;
-    
-    CCTK_REAL const rho  =  pow(alphaL,-2)*(T00 - 2*(beta2L*T02 + beta3L*T03) + 
-          2*(beta1L*(-T01 + beta2L*T12 + beta3L*T13) + beta2L*beta3L*T23) + T11*SQR(beta1L) + T22*SQR(beta2L) + 
-          T33*SQR(beta3L));
-    
-    CCTK_REAL const S1  =  (-T01 + beta1L*T11 + beta2L*T12 + beta3L*T13)*INV(alphaL);
-    
-    CCTK_REAL const S2  =  (-T02 + beta1L*T12 + beta2L*T22 + beta3L*T23)*INV(alphaL);
-    
-    CCTK_REAL const S3  =  (-T03 + beta1L*T13 + beta2L*T23 + beta3L*T33)*INV(alphaL);
+    CCTK_REAL const S3  =  (-eTtzL + beta1L*eTxzL + beta2L*eTyzL + beta3L*eTzzL)*INV(alphaL);
     
     CCTK_REAL const MadmL  =  0.01989436788648691697111047042156429525431*
         (-((Gt111*Gtlu111 + Gt112*Gtlu112 + Gt113*Gtlu113 + Gt211*Gtlu121 + Gt212*Gtlu122 + Gt213*Gtlu123 + Gt311*Gtlu131 + 
diff --git a/m/McLachlan_ADMConstraints.m b/m/McLachlan_ADMConstraints.m
index a3312a3..2ac7b6c 100644
--- a/m/McLachlan_ADMConstraints.m
+++ b/m/McLachlan_ADMConstraints.m
@@ -155,10 +155,10 @@ K11=kxx; K12=kxy; K22=kyy; K13=kxz; K23=kyz; K33=kzz;
 alpha=alp;
 beta1=betax; beta2=betay; beta3=betaz;
 
-Map [DefineTensor,
-     {eTtt,
-      eTtx, eTty, eTtz,
-      eTxx, eTxy, eTxz, eTyy, eTyz, eTzz}];
+(* Use the TmunuBase variable names *)
+T00=eTtt;
+T01=eTtx; T02=eTty; T03=eTtz;
+T11=eTxx; T12=eTxy; T22=eTyy; T13=eTxz; T23=eTyz; T33=eTzz;
 
 (******************************************************************************)
 (* Expressions *)
@@ -196,13 +196,13 @@ declaredGroupNames = Map [First, declaredGroups];
 
 
 extraGroups =
-  {{"ADMBase::metric",  {gxx, gxy, gxz, gyy, gyz, gzz}},
+  {{"Grid::coordinates", {x, y, z, r}},
+   {"ADMBase::metric",  {gxx, gxy, gxz, gyy, gyz, gzz}},
    {"ADMBase::curv",    {kxx, kxy, kxz, kyy, kyz, kzz}},
    {"ADMBase::lapse",   {alp}},
    {"ADMBase::dtlapse", {dtalp}},
    {"ADMBase::shift",   {betax, betay, betaz}},
    {"ADMBase::dtshift", {dtbetax, dtbetay, dtbetaz}},
-   {"Grid::coordinates", {x, y, z, r}},
    {"TmunuBase::stress_energy_scalar", {eTtt}},
    {"TmunuBase::stress_energy_vector", {eTtx, eTty, eTtz}},
    {"TmunuBase::stress_energy_tensor", {eTxx, eTxy, eTxz, eTyy, eTyz, eTzz}},
@@ -230,7 +230,7 @@ ADMConstraintsCalc =
   Where -> Interior,
   Shorthands -> {detg, gu[ua,ub], G[ua,lb,lc],
                  R[la,lb], trR, Km[la,lb], trK,
-                 T00, T0[la], T[la,lb], rho, S[la]},
+                 rho, S[la]},
   Equations -> 
   {
     detg -> detgExpr,
@@ -247,17 +247,6 @@ ADMConstraintsCalc =
 
     (* Matter terms *)
 
-    T00 -> eTtt,
-    T01 -> eTtx,
-    T02 -> eTty,
-    T03 -> eTtz,
-    T11 -> eTxx,
-    T12 -> eTxy,
-    T13 -> eTxz,
-    T22 -> eTyy,
-    T23 -> eTyz,
-    T33 -> eTzz,
-
     (* rho = n^a n^b T_ab *)
     rho -> 1/alpha^2 (T00 - 2 beta[ui] T0[li] + beta[ui] beta[uj] T[li,lj]),
 
diff --git a/m/McLachlan_ADMQuantities.m b/m/McLachlan_ADMQuantities.m
index fec1d56..932ed0c 100644
--- a/m/McLachlan_ADMQuantities.m
+++ b/m/McLachlan_ADMQuantities.m
@@ -132,15 +132,6 @@ Map [DefineTensor,
      {normal, tangentA, tangentB, dir,
       nn, nu, nlen, nlen2, su, vg,
       J, dJ,
-      (*
-      phi,
-      gt11, gt12, gt13, gt22, gt23, gt33,
-      Xt1, Xt2, Xt3,
-      trK,
-      At11, At12, At13, At22, At23, At33,
-      alpha,
-      beta1, beta2, beta3,
-      *)
       g, K, alpha, beta, detg, gu, G, R, trR, Km, trK,
       phi, gt, At, Xt, Xtn, alpha, A, beta, B, Atm, Atu, trA, Ats, trAts,
       ephi, detgt, gtu, ddetgt, dgtu, ddgtu, Gtl, Gtlu, Gt,
@@ -164,10 +155,10 @@ AssertSymmetricIncreasing [ddgtu[ua,ub,lc,ld], lc, ld];
 (* Use the CartGrid3D variable names *)
 x1=x; x2=y; x3=z;
 
-Map [DefineTensor,
-     {eTtt,
-      eTtx, eTty, eTtz,
-      eTxx, eTxy, eTxz, eTyy, eTyz, eTzz}];
+(* Use the TmunuBase variable names *)
+T00=eTtt;
+T01=eTtx; T02=eTty; T03=eTtz;
+T11=eTxx; T12=eTxy; T22=eTyy; T13=eTxz; T23=eTyz; T33=eTzz;
 
 (******************************************************************************)
 (* Expressions *)
@@ -202,7 +193,8 @@ declaredGroupNames = Map [First, declaredGroups];
 
 
 extraGroups =
-  {{"ADMBase::metric",   {gxx, gxy, gxz, gyy, gyz, gzz}},
+  {{"Grid::coordinates", {x, y, z, r}},
+   {"ADMBase::metric",   {gxx, gxy, gxz, gyy, gyz, gzz}},
    {"ADMBase::curv",     {kxx, kxy, kxz, kyy, kyz, kzz}},
    {"ADMBase::lapse",    {alp}},
    {"ADMBase::dtlapse",  {dtalp}},
@@ -218,7 +210,6 @@ extraGroups =
    {"TmunuBase::stress_energy_scalar", {eTtt}},
    {"TmunuBase::stress_energy_vector", {eTtx, eTty, eTtz}},
    {"TmunuBase::stress_energy_tensor", {eTxx, eTxy, eTxz, eTyy, eTyz, eTzz}},
-   {"Grid::coordinates", {x, y, z, r}},
    {"Coordinates::jacobian",  {J11, J12, J13, J21, J22, J23, J31, J32, J33}},
    {"Coordinates::jacobian2", {dJ111, dJ112, dJ113, dJ122, dJ123, dJ133,
                                dJ211, dJ212, dJ213, dJ222, dJ223, dJ233,
@@ -246,7 +237,7 @@ ADMQuantitiesCalc =
                  Xtn[ua], Rt[la,lb], trRt,
                  Atm[ua,lb],
                  ephi,
-                 T00, T0[la], T[la,lb], rho, S[la]},
+                 rho, S[la]},
   Equations -> 
   {
     detgt          -> 1 (* detgtExpr *),
@@ -279,17 +270,6 @@ ADMQuantitiesCalc =
 
     (* Matter terms *)
 
-    T00 -> eTtt,
-    T01 -> eTtx,
-    T02 -> eTty,
-    T03 -> eTtz,
-    T11 -> eTxx,
-    T12 -> eTxy,
-    T13 -> eTxz,
-    T22 -> eTyy,
-    T23 -> eTyz,
-    T33 -> eTzz,
-
     (* rho = n^a n^b T_ab *)
     rho -> addMatter
            1/alpha^2 (T00 - 2 beta[ui] T0[li] + beta[ui] beta[uj] T[li,lj]),
diff --git a/m/McLachlan_BSSN.m b/m/McLachlan_BSSN.m
index 4c31040..be448fe 100644
--- a/m/McLachlan_BSSN.m
+++ b/m/McLachlan_BSSN.m
@@ -532,19 +532,6 @@ evolCalc =
     
     (* Matter terms *)
     
-    (*
-    T00 -> addMatter eTtt,
-    T01 -> addMatter eTtx,
-    T02 -> addMatter eTty,
-    T03 -> addMatter eTtz,
-    T11 -> addMatter eTxx,
-    T12 -> addMatter eTxy,
-    T13 -> addMatter eTxz,
-    T22 -> addMatter eTyy,
-    T23 -> addMatter eTyz,
-    T33 -> addMatter eTzz,
-    *)
-    
     (* rho = n^a n^b T_ab *)
     rho -> addMatter
            (1/alpha^2 (T00 - 2 beta[ui] T0[li] + beta[ui] beta[uj] T[li,lj])),
@@ -661,18 +648,6 @@ evol1Calc =
     
     (* Matter terms *)
     
-    (*
-    T01 -> addMatter eTtx,
-    T02 -> addMatter eTty,
-    T03 -> addMatter eTtz,
-    T11 -> addMatter eTxx,
-    T12 -> addMatter eTxy,
-    T13 -> addMatter eTxz,
-    T22 -> addMatter eTyy,
-    T23 -> addMatter eTyz,
-    T33 -> addMatter eTzz,
-    *)
-    
     (* S_i = -p^a_i n^b T_ab, where p^a_i = delta^a_i + n^a n_i *)
     S[li] -> addMatter (-1/alpha (T0[li] - beta[uj] T[li,lj])),
     
@@ -786,19 +761,6 @@ evol2Calc =
     
     (* Matter terms *)
     
-    (*
-    T00 -> addMatter eTtt,
-    T01 -> addMatter eTtx,
-    T02 -> addMatter eTty,
-    T03 -> addMatter eTtz,
-    T11 -> addMatter eTxx,
-    T12 -> addMatter eTxy,
-    T13 -> addMatter eTxz,
-    T22 -> addMatter eTyy,
-    T23 -> addMatter eTyz,
-    T33 -> addMatter eTzz,
-    *)
-    
     (* rho = n^a n^b T_ab *)
     rho -> addMatter
            (1/alpha^2 (T00 - 2 beta[ui] T0[li] + beta[ui] beta[uj] T[li,lj])),
@@ -1045,19 +1007,6 @@ constraintsCalc =
     
     (* Matter terms *)
     
-    (*
-    T00 -> eTtt,
-    T01 -> eTtx,
-    T02 -> eTty,
-    T03 -> eTtz,
-    T11 -> eTxx,
-    T12 -> eTxy,
-    T13 -> eTxz,
-    T22 -> eTyy,
-    T23 -> eTyz,
-    T33 -> eTzz,
-    *)
-    
     (* rho = n^a n^b T_ab *)
     rho -> 1/alpha^2 (T00 - 2 beta[ui] T0[li] + beta[ui] beta[uj] T[li,lj]),