diff --git a/Z4co/src/enforce.cxx b/Z4co/src/enforce.cxx
index dbb43f75..7542e015 100644
--- a/Z4co/src/enforce.cxx
+++ b/Z4co/src/enforce.cxx
@@ -44,8 +44,6 @@ extern "C" void Z4co_Enforce(CCTK_ARGUMENTS) {
   typedef simdl<CCTK_REAL> vbool;
   constexpr size_t vsize = tuple_size_v<vreal>;
 
-  const auto delta3 = one<smat<vreal, 3>>()();
-
 #ifdef __CUDACC__
   const nvtxRangeId_t range = nvtxRangeStartA("Z4co_Enforce::enforce");
 #endif
@@ -63,21 +61,20 @@ extern "C" void Z4co_Enforce(CCTK_ARGUMENTS) {
 
         // Enforce floors
 
-        const vreal chi = fmax(vreal(chi_floor - 1), chi_old);
-        const vreal alphaG = fmax(vreal(alphaG_floor - 1), alphaG_old);
+        const vreal chi = fmax(vreal(chi_floor), chi_old);
+        const vreal alphaG = fmax(vreal(alphaG_floor), alphaG_old);
 
         // Enforce algebraic constraints
         // See arXiv:1212.2901 [gr-qc].
 
-        const vreal detgammat_old = calc_det(delta3 + gammat_old);
-        const vreal chi1_old = 1 / cbrt(detgammat_old) - 1;
+        const vreal detgammat_old = calc_det(gammat_old);
+        const vreal chi1_old = 1 / cbrt(detgammat_old);
         const smat<vreal, 3> gammat([&](int a, int b) ARITH_INLINE {
-          return (1 + chi1_old) * (delta3(a, b) + gammat_old(a, b)) -
-                 delta3(a, b);
+          return chi1_old * gammat_old(a, b);
         });
 #ifdef CCTK_DEBUG
-        const vreal detgammat = calc_det(delta3 + gammat);
-        const vreal gammat_norm = maxabs(delta3 + gammat);
+        const vreal detgammat = calc_det(gammat);
+        const vreal gammat_norm = maxabs(gammat);
         const vreal gammat_scale = gammat_norm;
 #if !defined __CUDACC__ && !defined __HIPCC__
         if (!(all(fabs(detgammat - 1) <= 1.0e-12 * gammat_scale))) {
@@ -90,20 +87,19 @@ extern "C" void Z4co_Enforce(CCTK_ARGUMENTS) {
         assert(all(fabs(detgammat - 1) <= 1.0e-12 * gammat_scale));
 #endif
 
-        const smat<vreal, 3> gammatu =
-            calc_inv(delta3 + gammat, vreal(1)) - delta3;
+        const smat<vreal, 3> gammatu = calc_inv(gammat, vreal(1));
 
         const vreal traceAt_old = sum_symm<3>([&](int x, int y) ARITH_INLINE {
-          return (delta3(x, y) + gammatu(x, y)) * At_old(x, y);
+          return gammatu(x, y) * At_old(x, y);
         });
         const smat<vreal, 3> At([&](int a, int b) ARITH_INLINE {
-          return At_old(a, b) - traceAt_old / 3 * (delta3(a, b) + gammat(a, b));
+          return At_old(a, b) - traceAt_old / 3 * gammat(a, b);
         });
 #ifdef CCTK_DEBUG
         const vreal traceAt = sum_symm<3>([&](int x, int y) ARITH_INLINE {
-          return (delta3(x, y) + gammatu(x, y)) * At(x, y);
+          return gammatu(x, y) * At(x, y);
         });
-        const vreal gammatu_norm = maxabs(delta3 + gammatu);
+        const vreal gammatu_norm = maxabs(gammatu);
         const vreal At_norm = maxabs(At);
         const vreal At_scale = fmax(fmax(gammat_norm, gammatu_norm), At_norm);
 #if !defined __CUDACC__ && !defined __HIPCC__
diff --git a/Z4co/src/initial1.cxx b/Z4co/src/initial1.cxx
index 3f6e226a..7fd3d6af 100644
--- a/Z4co/src/initial1.cxx
+++ b/Z4co/src/initial1.cxx
@@ -74,8 +74,6 @@ extern "C" void Z4co_Initial1(CCTK_ARGUMENTS) {
   typedef simdl<CCTK_REAL> vbool;
   constexpr size_t vsize = tuple_size_v<vreal>;
 
-  const auto delta3 = one<smat<vreal, 3>>()();
-
   const Loop::GridDescBaseDevice grid(cctkGH);
 #ifdef __CUDACC__
   const nvtxRangeId_t range = nvtxRangeStartA("Z4co_Initial1::initial1");
@@ -95,11 +93,10 @@ extern "C" void Z4co_Initial1(CCTK_ARGUMENTS) {
         const vreal detg = calc_det(g);
         const smat<vreal, 3> gu = calc_inv(g, detg);
 
-        const vreal chi = 1 / cbrt(detg) - 1;
+        const vreal chi = 1 / cbrt(detg);
 
-        const smat<vreal, 3> gammat([&](int a, int b) ARITH_INLINE {
-          return (1 + chi) * g(a, b) - delta3(a, b);
-        });
+        const smat<vreal, 3> gammat([&](int a, int b)
+                                        ARITH_INLINE { return chi * g(a, b); });
 
         const vreal trK = sum_symm<3>(
             [&](int x, int y) ARITH_INLINE { return gu(x, y) * K(x, y); });
@@ -109,10 +106,10 @@ extern "C" void Z4co_Initial1(CCTK_ARGUMENTS) {
         const vreal Kh = trK - 2 * Theta;
 
         const smat<vreal, 3> At([&](int a, int b) ARITH_INLINE {
-          return (1 + chi) * (K(a, b) - trK / 3 * g(a, b));
+          return chi * (K(a, b) - trK / 3 * g(a, b));
         });
 
-        const vreal alphaG = alp - 1;
+        const vreal alphaG = alp;
 
         const vec<vreal, 3> betaG([&](int a) ARITH_INLINE { return beta(a); });
 
diff --git a/Z4co/src/initial2.cxx b/Z4co/src/initial2.cxx
index 152066c7..fb87e2de 100644
--- a/Z4co/src/initial2.cxx
+++ b/Z4co/src/initial2.cxx
@@ -44,8 +44,6 @@ extern "C" void Z4co_Initial2(CCTK_ARGUMENTS) {
   typedef simdl<CCTK_REAL> vbool;
   constexpr size_t vsize = tuple_size_v<vreal>;
 
-  const auto delta3 = one<smat<vreal, 3>>()();
-
   const Loop::GridDescBaseDevice grid(cctkGH);
 #ifdef __CUDACC__
   const nvtxRangeId_t range = nvtxRangeStartA("Z4co_Initial2::initial2");
@@ -59,19 +57,17 @@ extern "C" void Z4co_Initial2(CCTK_ARGUMENTS) {
         const smat<vreal, 3> gammat = gf_gammat1(mask, index1);
 
         // Calculate Z4c variables (only Gamt)
-        const smat<vreal, 3> gammatu =
-            calc_inv(delta3 + gammat, vreal(1)) - delta3;
+        const smat<vreal, 3> gammatu = calc_inv(gammat, vreal(1));
 
         const smat<vec<vreal, 3>, 3> dgammat([&](int a, int b) {
           return deriv(mask, gf_gammat1(a, b), p.I, dx);
         });
 
         const vec<smat<vreal, 3>, 3> Gammatl = calc_gammal(dgammat);
-        const vec<smat<vreal, 3>, 3> Gammat =
-            calc_gamma(delta3 + gammatu, Gammatl);
+        const vec<smat<vreal, 3>, 3> Gammat = calc_gamma(gammatu, Gammatl);
         const vec<vreal, 3> Gamt([&](int a) ARITH_INLINE {
           return sum_symm<3>([&](int x, int y) ARITH_INLINE {
-            return (delta3(x, y) + gammatu(x, y)) * Gammat(a)(x, y);
+            return gammatu(x, y) * Gammat(a)(x, y);
           });
         });