From 19c297d49cc8133cdc3e15e80cbb7ac2fb1c34b8 Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Wed, 12 Aug 2020 17:06:34 +0200 Subject: Vectorize 2 --- ML_BSSN/src/ML_BSSN_RHS.cc | 68 +++-- ML_BSSN/src/ML_BSSN_constraints.cc | 104 ++++---- ML_BSSN/src/ML_BSSN_lapse_evol.cc | 41 ++- ML_BSSN/src/eval_timediff.cc | 504 ++++++++++++++++++------------------- cartoon_constr.c | 2 +- 5 files changed, 377 insertions(+), 342 deletions(-) diff --git a/ML_BSSN/src/ML_BSSN_RHS.cc b/ML_BSSN/src/ML_BSSN_RHS.cc index a3bb92e..1152240 100644 --- a/ML_BSSN/src/ML_BSSN_RHS.cc +++ b/ML_BSSN/src/ML_BSSN_RHS.cc @@ -2518,55 +2518,47 @@ static void ML_BSSN_RHS_Body(const cGH* restrict const cctkGH, const int dir, co } #else { - double gu[3][3] = {{ gu11, gu12, gu13 }, - { gu12, gu22, gu23 }, - { gu13, gu23, gu33 }}; - double G[3][3][3] = {{{ G111, G112, G113 }, - { G112, G122, G123 }, - { G113, G123, G133 }}, - {{ G211, G212, G213 }, - { G212, G222, G223 }, - { G213, G223, G233 }}, - {{ G311, G312, G313 }, - { G312, G322, G323 }, - { G313, G323, G333 }}}; + // double gu[3][3] = {{ gu11, gu12, gu13 }, + // { gu12, gu22, gu23 }, + // { gu13, gu23, gu33 }}; + // double G[3][3][3] = {{{ G111, G112, G113 }, + // { G112, G122, G123 }, + // { G113, G123, G133 }}, + // {{ G211, G212, G213 }, + // { G212, G222, G223 }, + // { G213, G223, G233 }}, + // {{ G311, G312, G313 }, + // { G312, G322, G323 }, + // { G313, G323, G333 }}}; - double X1 = 0.0, X3 = 0.0; + // double X1 = 0.0, X3 = 0.0; - double dux_alpha, duz_alpha; + // double dux_alpha, duz_alpha; - for (int fn = 0; fn < 3; fn++) - for (int bz = 0; bz < 3; bz++) { - X1 += gu[fn][bz] * G[0][fn][bz]; - X3 += gu[fn][bz] * G[2][fn][bz]; - } + // for (int fn = 0; fn < 3; fn++) + // for (int bz = 0; bz < 3; bz++) { + // X1 += gu[fn][bz] * G[0][fn][bz]; + // X3 += gu[fn][bz] * G[2][fn][bz]; + // } - dux_alpha = gu[0][0] * PDstandardNth1alpha + gu[0][2] * PDstandardNth3alpha; - duz_alpha = gu[2][0] * PDstandardNth1alpha + gu[2][2] * PDstandardNth3alpha; + // dux_alpha = gu[0][0] * PDstandardNth1alpha + gu[0][2] * PDstandardNth3alpha; + // duz_alpha = gu[2][0] * PDstandardNth1alpha + gu[2][2] * PDstandardNth3alpha; - beta1rhsL = alphaL * alphaL * X1 - alphaL * dux_alpha - alphaL * BetaDriver * beta1L; - beta2rhsL = 0.0; - beta3rhsL = alphaL * alphaL * X3 - alphaL * duz_alpha - alphaL * BetaDriver * beta3L; + // beta1rhsL = alphaL * alphaL * X1 - alphaL * dux_alpha - alphaL * BetaDriver * beta1L; + // beta2rhsL = 0.0; + // beta3rhsL = alphaL * alphaL * X3 - alphaL * duz_alpha - alphaL * BetaDriver * beta3L; } #endif else { - double xx = x[index]; - double zz = z[index]; - double t = cctkGH->cctk_time; - //double rr = sqrt(xx * xx + zz * zz); - //double phi = atan2(zz, xx); + beta1rhsL = theta * ToReal(ShiftGammaCoeff) * + (ToReal(ShiftBCoeff) * B1L + (kone - ToReal(ShiftBCoeff)) * (Xt1L - eta * ToReal(BetaDriver) * beta1L)); - //double gauge_r = -20.0 * (1.0 + 0.0 * cos(2 * phi) * rr * rr) * rr / (1.0 + rr * rr * rr * rr); + beta2rhsL = theta * ToReal(ShiftGammaCoeff) * + (ToReal(ShiftBCoeff) * B2L + (kone - ToReal(ShiftBCoeff)) * (Xt2L - eta * ToReal(BetaDriver) * beta2L)); - beta1rhsL = theta * ShiftGammaCoeff * - (ShiftBCoeff * B1L + (1.0 - ShiftBCoeff) * (Xt1L - eta * BetaDriver * beta1L)); - - beta2rhsL = theta * ShiftGammaCoeff * - (ShiftBCoeff * B2L + (1.0 - ShiftBCoeff) * (Xt2L - eta * BetaDriver * beta2L)); - - beta3rhsL = theta * ShiftGammaCoeff * - (ShiftBCoeff * B3L + (1.0 - ShiftBCoeff) * (Xt3L - eta * BetaDriver * beta3L)); + beta3rhsL = theta * ToReal(ShiftGammaCoeff) * + (ToReal(ShiftBCoeff) * B3L + (kone - ToReal(ShiftBCoeff)) * (Xt3L - eta * ToReal(BetaDriver) * beta3L)); } CCTK_REAL_VEC B1rhsL CCTK_ATTRIBUTE_UNUSED = diff --git a/ML_BSSN/src/ML_BSSN_constraints.cc b/ML_BSSN/src/ML_BSSN_constraints.cc index 62f8953..50c2e48 100644 --- a/ML_BSSN/src/ML_BSSN_constraints.cc +++ b/ML_BSSN/src/ML_BSSN_constraints.cc @@ -240,8 +240,6 @@ static void ML_BSSN_constraints_Body(const cGH* restrict const cctkGH, const int { const ptrdiff_t index CCTK_ATTRIBUTE_UNUSED = di*i + dj*j + dk*k; // vec_iter_counter+=CCTK_REAL_VEC_SIZE; - if (fabs(y[index]) > 1e-8) - continue; /* Assign local copies of grid functions */ @@ -822,51 +820,71 @@ static void ML_BSSN_constraints_Body(const cGH* restrict const cctkGH, const int default: CCTK_BUILTIN_UNREACHABLE(); } - - double xx = x[index]; - int origin = fabs(xx) < 1e-8; + CCTK_REAL_VEC kzero = ToReal(0.0); + CCTK_REAL_VEC kone = ToReal(1.0); + CCTK_REAL_VEC ktwo = ToReal(2.0); + CCTK_REAL_VEC xx = vec_load(x[index]); + CCTK_REAL_VEC x2 = kmul(xx, xx); + + CCTK_REAL_VEC xinv = kdiv(kone, xx); + CCTK_REAL_VEC x2inv = kdiv(kone, x2); - PDstandardNth2gt11 = 0; - PDstandardNth22gt11 = origin ? PDstandardNth11gt22 : PDstandardNth1gt11 / xx - 2 * (gt11[index] - gt22[index]) / (xx * xx); - PDstandardNth12gt11 = 0; - PDstandardNth23gt11 = 0; - PDstandardNth2gt12 = origin ? PDstandardNth1gt11 - PDstandardNth1gt22 : (gt11[index] - gt22[index]) / xx; - PDstandardNth22gt12 = 0; - PDstandardNth12gt12 = origin ? (PDstandardNth11gt11 - PDstandardNth11gt22) / 2 : (PDstandardNth1gt11 - PDstandardNth1gt22) / xx - (gt11[index] - gt22[index]) / (xx * xx); - PDstandardNth23gt12 = origin ? PDstandardNth13gt11 - PDstandardNth13gt22 : (PDstandardNth3gt11 - PDstandardNth3gt22) / xx; - PDstandardNth2gt13 = 0; - PDstandardNth22gt13 = origin ? PDstandardNth11gt13 / 2 : PDstandardNth1gt13 / xx - gt13[index] / (xx * xx); - PDstandardNth12gt13 = 0; - PDstandardNth23gt13 = 0; - PDstandardNth2gt22 = 0; - PDstandardNth22gt22 = origin ? PDstandardNth11gt11 : PDstandardNth1gt22 / xx + 2 * (gt11[index] - gt22[index]) / (xx * xx); - PDstandardNth12gt22 = 0; - PDstandardNth23gt22 = 0; - PDstandardNth2gt23 = origin ? PDstandardNth1gt13 : gt13[index] / xx; - PDstandardNth22gt23 = 0; - PDstandardNth12gt23 = origin ? PDstandardNth11gt13 / 2 : PDstandardNth1gt13 / xx - gt13[index] / (xx * xx); - PDstandardNth23gt23 = origin ? PDstandardNth13gt13 : PDstandardNth3gt13 / xx; - PDstandardNth2gt33 = 0; - PDstandardNth22gt33 = origin ? PDstandardNth11gt33 : PDstandardNth1gt33 / xx; - PDstandardNth12gt33 = 0; - PDstandardNth23gt33 = 0; + CCTK_REAL_VEC absx = kfabs(xx); + CCTK_REAL_VEC eps = ToReal(1e-8); + CCTK_BOOLEAN_VEC origin = kcmplt(absx, eps); - PDstandardNth2phi = 0; - PDstandardNth22phi = origin ? PDstandardNth11phi : PDstandardNth1phi / xx; - PDstandardNth12phi = 0; - PDstandardNth23phi = 0; - PDstandardNth2trK = 0; - PDstandardNth2Xt1 = 0; - PDstandardNth2Xt2 = origin ? PDstandardNth1Xt1 : Xt1[index] / xx; - PDstandardNth2Xt3 = 0; + PDstandardNth2gt11 = kzero; + PDstandardNth22gt11 = kifthen(origin, PDstandardNth11gt22, + ksub(kmul(PDstandardNth1gt11, xinv), kmul(ktwo, kmul(ksub(gt11L, gt22L), x2inv)))); + PDstandardNth12gt11 = kzero; + PDstandardNth23gt11 = kzero; + PDstandardNth2gt12 = kifthen(origin, ksub(PDstandardNth1gt11, PDstandardNth1gt22), + kmul(ksub(gt11L, gt22L), xinv)); + PDstandardNth22gt12 = kzero; + PDstandardNth12gt12 = kifthen(origin, kmul(ksub(PDstandardNth11gt11, PDstandardNth11gt22), khalf), + ksub(kmul(ksub(PDstandardNth1gt11, PDstandardNth1gt22), xinv), kmul(ksub(gt11L, gt22L), x2inv))); + PDstandardNth23gt12 = kifthen(origin, ksub(PDstandardNth13gt11, PDstandardNth13gt22), + kmul(ksub(PDstandardNth3gt11, PDstandardNth3gt22), xinv)); + PDstandardNth2gt13 = kzero; + PDstandardNth22gt13 = kifthen(origin, kmul(PDstandardNth11gt13, khalf), + ksub(kmul(PDstandardNth1gt13, xinv), kmul(gt13L, x2inv))); + PDstandardNth12gt13 = kzero; + PDstandardNth23gt13 = kzero; + PDstandardNth2gt22 = kzero; + PDstandardNth22gt22 = kifthen(origin, PDstandardNth11gt11, + kadd(kmul(PDstandardNth1gt22, xinv), kmul(ktwo, kmul(ksub(gt11L, gt22L), x2inv)))); + PDstandardNth12gt22 = kzero; + PDstandardNth23gt22 = kzero; + PDstandardNth2gt23 = kifthen(origin, PDstandardNth1gt13, kmul(gt13L, xinv)); + PDstandardNth22gt23 = kzero; + PDstandardNth12gt23 = kifthen(origin, kmul(PDstandardNth11gt13, khalf), + ksub(kmul(PDstandardNth1gt13, xinv), kmul(gt13L, x2inv))); + PDstandardNth23gt23 = kifthen(origin, PDstandardNth13gt13, kmul(PDstandardNth3gt13, xinv)); + PDstandardNth2gt33 = kzero; + PDstandardNth22gt33 = kifthen(origin, PDstandardNth11gt33, kmul(PDstandardNth1gt33, xinv)); + PDstandardNth12gt33 = kzero; + PDstandardNth23gt33 = kzero; - PDstandardNth2At11 = 0; //origin ? -2 * PDstandardNth1At12 : -2 * At12[index] / xx; - PDstandardNth2At12 = origin ? PDstandardNth1At11 - PDstandardNth1At22 : (At11[index] - At22[index]) / xx; - PDstandardNth2At13 = 0; //origin ? -PDstandardNth1At23 : - At23[index] / xx; - PDstandardNth2At22 = 0; //origin ? 2 * PDstandardNth1At12 : 2 * At12[index] / xx; - PDstandardNth2At23 = origin ? PDstandardNth1At13 : At13[index] / xx; - PDstandardNth2At33 = 0; + PDstandardNth2phi = kzero; + PDstandardNth22phi = kifthen(origin, PDstandardNth11phi, kmul(PDstandardNth1phi, xinv)); + PDstandardNth12phi = kzero; + PDstandardNth23phi = kzero; + PDstandardNth2Xt1 = kzero; + PDstandardNth2Xt2 = kifthen(origin, PDstandardNth1Xt1, kmul(Xt1L, xinv)); + PDstandardNth2Xt3 = kzero; + PDstandardNth2trK = kzero; + //PDstandardNth2At11 = origin ? -2 * PDstandardNth1At12 : -2 * At12[index] / xx; + //PDstandardNth2At13 = origin ? -PDstandardNth1At23 : - At23[index] / xx; + //PDstandardNth2At22 = origin ? 2 * PDstandardNth1At12 : 2 * At12[index] / xx; + PDstandardNth2At11 = kzero; + PDstandardNth2At12 = kifthen(origin, ksub(PDstandardNth1At11, PDstandardNth1At22), + kmul(ksub(At11L, At22L), xinv)); + PDstandardNth2At13 = kzero; + PDstandardNth2At22 = kzero; + PDstandardNth2At23 = kifthen(origin, PDstandardNth1At13, kmul(At13L, xinv)); + PDstandardNth2At33 = kzero; + /* Calculate temporaries and grid functions */ CCTK_REAL_VEC JacPDstandardNth11gt11 CCTK_ATTRIBUTE_UNUSED; CCTK_REAL_VEC JacPDstandardNth11gt12 CCTK_ATTRIBUTE_UNUSED; diff --git a/ML_BSSN/src/ML_BSSN_lapse_evol.cc b/ML_BSSN/src/ML_BSSN_lapse_evol.cc index e688b5a..0388510 100644 --- a/ML_BSSN/src/ML_BSSN_lapse_evol.cc +++ b/ML_BSSN/src/ML_BSSN_lapse_evol.cc @@ -211,22 +211,47 @@ static void ML_BSSN_lapse_evol_Body(const cGH* restrict const cctkGH, const int break; } - PDupwindNthAnti2alpha = 0.0; - PDupwindNthSymm2alpha = 0.0; + PDupwindNthAnti2alpha = ToReal(0.0); + PDupwindNthSymm2alpha = ToReal(0.0); + + CCTK_REAL_VEC alpharhsL; + +#if 1 if (lapse_gh) { - const double e4phi = conformalMethod == 1 ? 1.0 / (phiL * phiL) : exp(4.0 * phiL); - const double det_gamma = pow(e4phi, 3.0); + CCTK_REAL_VEC e4phi = kifthen(conformalMethod == 1, ToReal(1.0) / (phiL * phiL), kexp(ToReal(4.0) * phiL)); + CCTK_REAL_VEC det_gamma = kpow(e4phi, 3.0); - alpharhs[index] = -(alphaL * alphaL) * trKL + gh_eta_l_bar * pow(alphaL, gh_q) * (alphaL * alphaL) * log(pow(det_gamma, 0.5 * gh_p) / alphaL) - + LapseAdvectionCoeff * (beta1L * PDupwindNthAnti1alpha + beta3L * PDupwindNthAnti3alpha + - fabs(beta1L) * PDupwindNthSymm1alpha + fabs(beta3L) * PDupwindNthSymm3alpha); + alpharhsL = -(alphaL * alphaL) * trKL + ToReal(gh_eta_l_bar) * kpow(alphaL, gh_q) * (alphaL * alphaL) * klog(kpow(det_gamma, 0.5 * gh_p) / alphaL) + + ToReal(LapseAdvectionCoeff) * (beta1L * PDupwindNthAnti1alpha + beta3L * PDupwindNthAnti3alpha + + kfabs(beta1L) * PDupwindNthSymm1alpha + kfabs(beta3L) * PDupwindNthSymm3alpha); } else { - alpharhs[index] = -harmonicF * pow(alphaL, harmonicN) * trKL + WFactor * WL + alpharhsL = -ToReal(harmonicF) * kpow(alphaL, harmonicN) * trKL + ToReal(WFactor) * WL + + ToReal(LapseAdvectionCoeff) * (beta1L * PDupwindNthAnti1alpha + beta3L * PDupwindNthAnti3alpha + + kfabs(beta1L) * PDupwindNthSymm1alpha + kfabs(beta3L) * PDupwindNthSymm3alpha); + } +#else + { + const double t = cctkGH->cctk_time; + const double fact_harmonic = alphaL * alphaL; + const double fact_1plog = 2.0 * alphaL; + double fact; + + if (t < 5.0) + fact = fact_harmonic; + else if (t > 5.5) + fact = fact_1plog; + else + fact = fact_harmonic * exp(-36.0 * pow((t - 5.0) / 0.5, 4)) + fact_1plog * (1.0 - exp(-36.0 * pow((t - 5.0) / 0.5, 4))); + + alpharhs[index] = -fact * trKL + LapseAdvectionCoeff * (beta1L * PDupwindNthAnti1alpha + beta3L * PDupwindNthAnti3alpha + fabs(beta1L) * PDupwindNthSymm1alpha + fabs(beta3L) * PDupwindNthSymm3alpha); } +#endif + vec_store_partial_prepare(i,vecimin,vecimax); + vec_store_nta_partial(alpharhs[index],alpharhsL); } CCTK_ENDLOOP3STR(ML_BSSN_lapse_evol); } diff --git a/ML_BSSN/src/eval_timediff.cc b/ML_BSSN/src/eval_timediff.cc index ec35923..ae4f827 100644 --- a/ML_BSSN/src/eval_timediff.cc +++ b/ML_BSSN/src/eval_timediff.cc @@ -247,7 +247,7 @@ static void eval_timediff_body(const cGH* restrict const cctkGH, const int dir, CCTK_REAL_VEC eTttL, eTtxL, eTtyL, eTtzL, eTxxL, eTxyL, eTxzL, eTyyL, eTyzL, eTzzL CCTK_ATTRIBUTE_UNUSED ; - double beta_term_val; + //double beta_term_val; CCTK_REAL_VEC dJ111L, dJ112L, dJ113L, dJ122L, dJ123L, dJ133L, dJ211L, dJ212L, dJ213L, dJ222L, dJ223L, dJ233L, dJ311L, dJ312L, dJ313L, dJ322L, dJ323L, dJ333L, J11L, J12L, J13L, J21L, J22L, J23L, J31L, J32L, J33L CCTK_ATTRIBUTE_UNUSED ; @@ -2107,65 +2107,65 @@ static void eval_timediff_body(const cGH* restrict const cctkGH, const int dir, CCTK_REAL_VEC G333 CCTK_ATTRIBUTE_UNUSED = kadd(Gt333,kdiv(kmul(ToReal(0.5),kmadd(g33,kmul(kmadd(ddetg1,gu13,kmadd(ddetg2,gu23,kmul(ddetg3,gu33))),ToReal(-0.333333333333333333333333333333)),kmul(ddetg3,ToReal(2)))),detg)); - { - double Gt[3][3][3] = { - { - { Gt111, Gt112, Gt113 }, - { Gt112, Gt122, Gt123 }, - { Gt113, Gt123, Gt133 }, - }, - { - { Gt211, Gt212, Gt213 }, - { Gt212, Gt222, Gt223 }, - { Gt213, Gt223, Gt233 }, - }, - { - { Gt311, Gt312, Gt313 }, - { Gt312, Gt322, Gt323 }, - { Gt313, Gt323, Gt333 }, - } - }; - double gt[3][3] = { - { gt11L, gt12L, gt13L }, - { gt12L, gt22L, gt23L }, - { gt13L, gt23L, gt33L }, - }; - double gtu[3][3] = { - { gtu11, gtu12, gtu13 }, - { gtu12, gtu22, gtu23 }, - { gtu13, gtu23, gtu33 }, - }; + //{ + // double Gt[3][3][3] = { + // { + // { Gt111, Gt112, Gt113 }, + // { Gt112, Gt122, Gt123 }, + // { Gt113, Gt123, Gt133 }, + // }, + // { + // { Gt211, Gt212, Gt213 }, + // { Gt212, Gt222, Gt223 }, + // { Gt213, Gt223, Gt233 }, + // }, + // { + // { Gt311, Gt312, Gt313 }, + // { Gt312, Gt322, Gt323 }, + // { Gt313, Gt323, Gt333 }, + // } + // }; + // double gt[3][3] = { + // { gt11L, gt12L, gt13L }, + // { gt12L, gt22L, gt23L }, + // { gt13L, gt23L, gt33L }, + // }; + // double gtu[3][3] = { + // { gtu11, gtu12, gtu13 }, + // { gtu12, gtu22, gtu23 }, + // { gtu13, gtu23, gtu33 }, + // }; - double dphi[3] = { PDstandardNth1phi, 0.0, PDstandardNth3phi }; + // double dphi[3] = { PDstandardNth1phi, 0.0, PDstandardNth3phi }; - double G[3][3][3]; + // double G[3][3][3]; - for (int idx0 = 0; idx0 < 3; idx0++) - for (int idx1 = 0; idx1 < 3; idx1++) - for (int idx2 = 0; idx2 < 3; idx2++) { - G[idx0][idx1][idx2] = Gt[idx0][idx1][idx2] - (dphi[idx1] * (idx0 == idx2) + dphi[idx2] * (idx0 == idx1) - gt[idx1][idx2] * (gtu[idx0][0] * dphi[0] + gtu[idx0][2] * dphi[2])) / phiL; - } - G111 = G[0][0][0]; - G122 = G[0][1][1]; - G133 = G[0][2][2]; - G112 = G[0][0][1]; - G113 = G[0][0][2]; - G123 = G[0][1][2]; + // for (int idx0 = 0; idx0 < 3; idx0++) + // for (int idx1 = 0; idx1 < 3; idx1++) + // for (int idx2 = 0; idx2 < 3; idx2++) { + // G[idx0][idx1][idx2] = Gt[idx0][idx1][idx2] - (dphi[idx1] * (idx0 == idx2) + dphi[idx2] * (idx0 == idx1) - gt[idx1][idx2] * (gtu[idx0][0] * dphi[0] + gtu[idx0][2] * dphi[2])) / phiL; + // } + // G111 = G[0][0][0]; + // G122 = G[0][1][1]; + // G133 = G[0][2][2]; + // G112 = G[0][0][1]; + // G113 = G[0][0][2]; + // G123 = G[0][1][2]; - G211 = G[1][0][0]; - G222 = G[1][1][1]; - G233 = G[1][2][2]; - G212 = G[1][0][1]; - G213 = G[1][0][2]; - G223 = G[1][1][2]; + // G211 = G[1][0][0]; + // G222 = G[1][1][1]; + // G233 = G[1][2][2]; + // G212 = G[1][0][1]; + // G213 = G[1][0][2]; + // G223 = G[1][1][2]; - G311 = G[2][0][0]; - G322 = G[2][1][1]; - G333 = G[2][2][2]; - G312 = G[2][0][1]; - G313 = G[2][0][2]; - G323 = G[2][1][2]; - } + // G311 = G[2][0][0]; + // G322 = G[2][1][1]; + // G333 = G[2][2][2]; + // G312 = G[2][0][1]; + // G313 = G[2][0][2]; + // G323 = G[2][1][2]; + //} CCTK_REAL_VEC K11 CCTK_ATTRIBUTE_UNUSED = kmadd(At11L,e4phi,kmul(trKL,kmul(g11,ToReal(0.333333333333333333333333333333)))); @@ -2366,237 +2366,237 @@ static void eval_timediff_body(const cGH* restrict const cctkGH, const int dir, CCTK_REAL_VEC Xtdot3L CCTK_ATTRIBUTE_UNUSED = dotXt3; CCTK_REAL_VEC Kdot11L CCTK_ATTRIBUTE_UNUSED = - G111 * JacPDstandardNth1alpha + G211 * JacPDstandardNth2alpha + G311 * JacPDstandardNth3alpha + alphaL * (trKL * K11 + R11 - 2.0 * (K11 * Km11 + K12 * Km12 + K13 * Km13)) - JacPDstandardNth11alpha; + G111 * JacPDstandardNth1alpha + G211 * JacPDstandardNth2alpha + G311 * JacPDstandardNth3alpha + alphaL * (trKL * K11 + R11 - ktwo * (K11 * Km11 + K12 * Km12 + K13 * Km13)) - JacPDstandardNth11alpha; CCTK_REAL_VEC Kdot12L CCTK_ATTRIBUTE_UNUSED = - G112 * JacPDstandardNth1alpha + G212 * JacPDstandardNth2alpha + G312 * JacPDstandardNth3alpha + alphaL * (trKL * K12 + R12 - 2.0 * (K11 * Km21 + K12 * Km22 + K13 * Km23)) - JacPDstandardNth12alpha; + G112 * JacPDstandardNth1alpha + G212 * JacPDstandardNth2alpha + G312 * JacPDstandardNth3alpha + alphaL * (trKL * K12 + R12 - ktwo * (K11 * Km21 + K12 * Km22 + K13 * Km23)) - JacPDstandardNth12alpha; CCTK_REAL_VEC Kdot13L CCTK_ATTRIBUTE_UNUSED = - G113 * JacPDstandardNth1alpha + G213 * JacPDstandardNth2alpha + G313 * JacPDstandardNth3alpha + alphaL * (trKL * K13 + R13 - 2.0 * (K11 * Km31 + K12 * Km32 + K13 * Km33)) - JacPDstandardNth13alpha; + G113 * JacPDstandardNth1alpha + G213 * JacPDstandardNth2alpha + G313 * JacPDstandardNth3alpha + alphaL * (trKL * K13 + R13 - ktwo * (K11 * Km31 + K12 * Km32 + K13 * Km33)) - JacPDstandardNth13alpha; CCTK_REAL_VEC Kdot22L CCTK_ATTRIBUTE_UNUSED = - G122 * JacPDstandardNth1alpha + G222 * JacPDstandardNth2alpha + G322 * JacPDstandardNth3alpha + alphaL * (trKL * K22 + R22 - 2.0 * (K12 * Km21 + K22 * Km22 + K23 * Km23)) - JacPDstandardNth22alpha; + G122 * JacPDstandardNth1alpha + G222 * JacPDstandardNth2alpha + G322 * JacPDstandardNth3alpha + alphaL * (trKL * K22 + R22 - ktwo * (K12 * Km21 + K22 * Km22 + K23 * Km23)) - JacPDstandardNth22alpha; CCTK_REAL_VEC Kdot23L CCTK_ATTRIBUTE_UNUSED = - G123 * JacPDstandardNth1alpha + G223 * JacPDstandardNth2alpha + G323 * JacPDstandardNth3alpha + alphaL * (trKL * K23 + R23 - 2.0 * (K12 * Km31 + K22 * Km32 + K23 * Km33)) - JacPDstandardNth23alpha; + G123 * JacPDstandardNth1alpha + G223 * JacPDstandardNth2alpha + G323 * JacPDstandardNth3alpha + alphaL * (trKL * K23 + R23 - ktwo * (K12 * Km31 + K22 * Km32 + K23 * Km33)) - JacPDstandardNth23alpha; CCTK_REAL_VEC Kdot33L CCTK_ATTRIBUTE_UNUSED = - G133 * JacPDstandardNth1alpha + G233 * JacPDstandardNth2alpha + G333 * JacPDstandardNth3alpha + alphaL * (trKL * K33 + R33 - 2.0 * (K13 * Km31 + K23 * Km32 + K33 * Km33)) - JacPDstandardNth33alpha; + G133 * JacPDstandardNth1alpha + G233 * JacPDstandardNth2alpha + G333 * JacPDstandardNth3alpha + alphaL * (trKL * K33 + R33 - ktwo * (K13 * Km31 + K23 * Km32 + K33 * Km33)) - JacPDstandardNth33alpha; /* compute the beta term */ - { - const double gt[3][3] = { - { gt11L, gt12L, gt13L }, - { gt12L, gt22L, gt23L }, - { gt13L, gt23L, gt33L }, - }; + //{ + // const double gt[3][3] = { + // { gt11L, gt12L, gt13L }, + // { gt12L, gt22L, gt23L }, + // { gt13L, gt23L, gt33L }, + // }; - const double gtu[3][3] = { - { gtu11, gtu12, gtu13 }, - { gtu12, gtu22, gtu23 }, - { gtu13, gtu23, gtu33 }, - }; + // const double gtu[3][3] = { + // { gtu11, gtu12, gtu13 }, + // { gtu12, gtu22, gtu23 }, + // { gtu13, gtu23, gtu33 }, + // }; - const double dgt[3][3][3] = { - { - { PDstandardNth1gt11, PDstandardNth2gt11, PDstandardNth3gt11 }, - { PDstandardNth1gt12, PDstandardNth2gt12, PDstandardNth3gt12 }, - { PDstandardNth1gt13, PDstandardNth2gt13, PDstandardNth3gt13 }, - }, - { - { PDstandardNth1gt12, PDstandardNth2gt12, PDstandardNth3gt12 }, - { PDstandardNth1gt22, PDstandardNth2gt22, PDstandardNth3gt22 }, - { PDstandardNth1gt23, PDstandardNth2gt23, PDstandardNth3gt23 }, - }, - { - { PDstandardNth1gt13, PDstandardNth2gt13, PDstandardNth3gt13 }, - { PDstandardNth1gt23, PDstandardNth2gt23, PDstandardNth3gt23 }, - { PDstandardNth1gt33, PDstandardNth2gt33, PDstandardNth3gt33 }, - }, - }; + // const double dgt[3][3][3] = { + // { + // { PDstandardNth1gt11, PDstandardNth2gt11, PDstandardNth3gt11 }, + // { PDstandardNth1gt12, PDstandardNth2gt12, PDstandardNth3gt12 }, + // { PDstandardNth1gt13, PDstandardNth2gt13, PDstandardNth3gt13 }, + // }, + // { + // { PDstandardNth1gt12, PDstandardNth2gt12, PDstandardNth3gt12 }, + // { PDstandardNth1gt22, PDstandardNth2gt22, PDstandardNth3gt22 }, + // { PDstandardNth1gt23, PDstandardNth2gt23, PDstandardNth3gt23 }, + // }, + // { + // { PDstandardNth1gt13, PDstandardNth2gt13, PDstandardNth3gt13 }, + // { PDstandardNth1gt23, PDstandardNth2gt23, PDstandardNth3gt23 }, + // { PDstandardNth1gt33, PDstandardNth2gt33, PDstandardNth3gt33 }, + // }, + // }; - const double d2gt[3][3][3][3] = { - { - { - { PDstandardNth11gt11, PDstandardNth12gt11, PDstandardNth13gt11 }, - { PDstandardNth12gt11, PDstandardNth22gt11, PDstandardNth23gt11 }, - { PDstandardNth13gt11, PDstandardNth23gt11, PDstandardNth33gt11 }, - }, - { - { PDstandardNth11gt12, PDstandardNth12gt12, PDstandardNth13gt12 }, - { PDstandardNth12gt12, PDstandardNth22gt12, PDstandardNth23gt12 }, - { PDstandardNth13gt12, PDstandardNth23gt12, PDstandardNth33gt12 }, - }, - { - { PDstandardNth11gt13, PDstandardNth12gt13, PDstandardNth13gt13 }, - { PDstandardNth12gt13, PDstandardNth22gt13, PDstandardNth23gt13 }, - { PDstandardNth13gt13, PDstandardNth23gt13, PDstandardNth33gt13 }, - }, - }, - { - { - { PDstandardNth11gt12, PDstandardNth12gt12, PDstandardNth13gt12 }, - { PDstandardNth12gt12, PDstandardNth22gt12, PDstandardNth23gt12 }, - { PDstandardNth13gt12, PDstandardNth23gt12, PDstandardNth33gt12 }, - }, - { - { PDstandardNth11gt22, PDstandardNth12gt22, PDstandardNth13gt22 }, - { PDstandardNth12gt22, PDstandardNth22gt22, PDstandardNth23gt22 }, - { PDstandardNth13gt22, PDstandardNth23gt22, PDstandardNth33gt22 }, - }, - { - { PDstandardNth11gt23, PDstandardNth12gt23, PDstandardNth13gt23 }, - { PDstandardNth12gt23, PDstandardNth22gt23, PDstandardNth23gt23 }, - { PDstandardNth13gt23, PDstandardNth23gt23, PDstandardNth33gt23 }, - }, - }, - { - { - { PDstandardNth11gt13, PDstandardNth12gt13, PDstandardNth13gt13 }, - { PDstandardNth12gt13, PDstandardNth22gt13, PDstandardNth23gt13 }, - { PDstandardNth13gt13, PDstandardNth23gt13, PDstandardNth33gt13 }, - }, - { - { PDstandardNth11gt23, PDstandardNth12gt23, PDstandardNth13gt23 }, - { PDstandardNth12gt23, PDstandardNth22gt23, PDstandardNth23gt23 }, - { PDstandardNth13gt23, PDstandardNth23gt23, PDstandardNth33gt23 }, - }, - { - { PDstandardNth11gt33, PDstandardNth12gt33, PDstandardNth13gt33 }, - { PDstandardNth12gt33, PDstandardNth22gt33, PDstandardNth23gt33 }, - { PDstandardNth13gt33, PDstandardNth23gt33, PDstandardNth33gt33 }, - }, - }, - }; + // const double d2gt[3][3][3][3] = { + // { + // { + // { PDstandardNth11gt11, PDstandardNth12gt11, PDstandardNth13gt11 }, + // { PDstandardNth12gt11, PDstandardNth22gt11, PDstandardNth23gt11 }, + // { PDstandardNth13gt11, PDstandardNth23gt11, PDstandardNth33gt11 }, + // }, + // { + // { PDstandardNth11gt12, PDstandardNth12gt12, PDstandardNth13gt12 }, + // { PDstandardNth12gt12, PDstandardNth22gt12, PDstandardNth23gt12 }, + // { PDstandardNth13gt12, PDstandardNth23gt12, PDstandardNth33gt12 }, + // }, + // { + // { PDstandardNth11gt13, PDstandardNth12gt13, PDstandardNth13gt13 }, + // { PDstandardNth12gt13, PDstandardNth22gt13, PDstandardNth23gt13 }, + // { PDstandardNth13gt13, PDstandardNth23gt13, PDstandardNth33gt13 }, + // }, + // }, + // { + // { + // { PDstandardNth11gt12, PDstandardNth12gt12, PDstandardNth13gt12 }, + // { PDstandardNth12gt12, PDstandardNth22gt12, PDstandardNth23gt12 }, + // { PDstandardNth13gt12, PDstandardNth23gt12, PDstandardNth33gt12 }, + // }, + // { + // { PDstandardNth11gt22, PDstandardNth12gt22, PDstandardNth13gt22 }, + // { PDstandardNth12gt22, PDstandardNth22gt22, PDstandardNth23gt22 }, + // { PDstandardNth13gt22, PDstandardNth23gt22, PDstandardNth33gt22 }, + // }, + // { + // { PDstandardNth11gt23, PDstandardNth12gt23, PDstandardNth13gt23 }, + // { PDstandardNth12gt23, PDstandardNth22gt23, PDstandardNth23gt23 }, + // { PDstandardNth13gt23, PDstandardNth23gt23, PDstandardNth33gt23 }, + // }, + // }, + // { + // { + // { PDstandardNth11gt13, PDstandardNth12gt13, PDstandardNth13gt13 }, + // { PDstandardNth12gt13, PDstandardNth22gt13, PDstandardNth23gt13 }, + // { PDstandardNth13gt13, PDstandardNth23gt13, PDstandardNth33gt13 }, + // }, + // { + // { PDstandardNth11gt23, PDstandardNth12gt23, PDstandardNth13gt23 }, + // { PDstandardNth12gt23, PDstandardNth22gt23, PDstandardNth23gt23 }, + // { PDstandardNth13gt23, PDstandardNth23gt23, PDstandardNth33gt23 }, + // }, + // { + // { PDstandardNth11gt33, PDstandardNth12gt33, PDstandardNth13gt33 }, + // { PDstandardNth12gt33, PDstandardNth22gt33, PDstandardNth23gt33 }, + // { PDstandardNth13gt33, PDstandardNth23gt33, PDstandardNth33gt33 }, + // }, + // }, + // }; - const double Ric[3][3] = { - { R11, R12, R13 }, - { R12, R22, R23 }, - { R13, R23, R33 }, - }; + // const double Ric[3][3] = { + // { R11, R12, R13 }, + // { R12, R22, R23 }, + // { R13, R23, R33 }, + // }; - const double betaval[3] = {beta1L, 0.0, beta3L}; - const double dbeta[3][3] = {{ PDstandardNth1beta1, PDstandardNth2beta1, PDstandardNth3beta1 }, - { PDstandardNth1beta2, PDstandardNth2beta2, PDstandardNth3beta2 }, - { PDstandardNth1beta3, PDstandardNth2beta3, PDstandardNth3beta3 }}; + // const double betaval[3] = {beta1L, 0.0, beta3L}; + // const double dbeta[3][3] = {{ PDstandardNth1beta1, PDstandardNth2beta1, PDstandardNth3beta1 }, + // { PDstandardNth1beta2, PDstandardNth2beta2, PDstandardNth3beta2 }, + // { PDstandardNth1beta3, PDstandardNth2beta3, PDstandardNth3beta3 }}; - const double d2beta[3][3][3] = { - { - { PDstandardNth11beta1, PDstandardNth12beta1, PDstandardNth13beta1 }, - { PDstandardNth12beta1, PDstandardNth22beta1, PDstandardNth23beta1 }, - { PDstandardNth13beta1, PDstandardNth23beta1, PDstandardNth33beta1 }, - }, - { - { PDstandardNth11beta2, PDstandardNth12beta2, PDstandardNth13beta2 }, - { PDstandardNth12beta2, PDstandardNth22beta2, PDstandardNth23beta2 }, - { PDstandardNth13beta2, PDstandardNth23beta2, PDstandardNth33beta2 }, - }, - { - { PDstandardNth11beta3, PDstandardNth12beta3, PDstandardNth13beta3 }, - { PDstandardNth12beta3, PDstandardNth22beta3, PDstandardNth23beta3 }, - { PDstandardNth13beta3, PDstandardNth23beta3, PDstandardNth33beta3 }, - }, - }; + // const double d2beta[3][3][3] = { + // { + // { PDstandardNth11beta1, PDstandardNth12beta1, PDstandardNth13beta1 }, + // { PDstandardNth12beta1, PDstandardNth22beta1, PDstandardNth23beta1 }, + // { PDstandardNth13beta1, PDstandardNth23beta1, PDstandardNth33beta1 }, + // }, + // { + // { PDstandardNth11beta2, PDstandardNth12beta2, PDstandardNth13beta2 }, + // { PDstandardNth12beta2, PDstandardNth22beta2, PDstandardNth23beta2 }, + // { PDstandardNth13beta2, PDstandardNth23beta2, PDstandardNth33beta2 }, + // }, + // { + // { PDstandardNth11beta3, PDstandardNth12beta3, PDstandardNth13beta3 }, + // { PDstandardNth12beta3, PDstandardNth22beta3, PDstandardNth23beta3 }, + // { PDstandardNth13beta3, PDstandardNth23beta3, PDstandardNth33beta3 }, + // }, + // }; - const double G[3][3][3] = { - { - { G111, G112, G113 }, - { G112, G122, G123 }, - { G113, G123, G133 }, - }, - { - { G211, G212, G213 }, - { G212, G222, G223 }, - { G213, G223, G233 }, - }, - { - { G311, G312, G313 }, - { G312, G322, G323 }, - { G313, G323, G333 }, - } - }; + // const double G[3][3][3] = { + // { + // { G111, G112, G113 }, + // { G112, G122, G123 }, + // { G113, G123, G133 }, + // }, + // { + // { G211, G212, G213 }, + // { G212, G222, G223 }, + // { G213, G223, G233 }, + // }, + // { + // { G311, G312, G313 }, + // { G312, G322, G323 }, + // { G313, G323, G333 }, + // } + // }; - const double dalpha[3] = { PDstandardNth1alpha, 0.0, PDstandardNth3alpha }; + // const double dalpha[3] = { PDstandardNth1alpha, 0.0, PDstandardNth3alpha }; - const double dphi[3] = { PDstandardNth1phi, 0.0, PDstandardNth3phi }; - const double d2phi[3][3] = {{ PDstandardNth11phi, PDstandardNth12phi, PDstandardNth13phi }, - { PDstandardNth12phi, PDstandardNth22phi, PDstandardNth23phi }, - { PDstandardNth13phi, PDstandardNth23phi, PDstandardNth33phi }}; + // const double dphi[3] = { PDstandardNth1phi, 0.0, PDstandardNth3phi }; + // const double d2phi[3][3] = {{ PDstandardNth11phi, PDstandardNth12phi, PDstandardNth13phi }, + // { PDstandardNth12phi, PDstandardNth22phi, PDstandardNth23phi }, + // { PDstandardNth13phi, PDstandardNth23phi, PDstandardNth33phi }}; - double gu[3][3], Ricm[3][3], dg[3][3][3], dgu[3][3][3], d2g[3][3][3][3], dG[3][3][3][3]; - double d2beta_term[3] = { 0.0 }; + // double gu[3][3], Ricm[3][3], dg[3][3][3], dgu[3][3][3], d2g[3][3][3][3], dG[3][3][3][3]; + // double d2beta_term[3] = { 0.0 }; - double ric_term; + // double ric_term; - for (int i = 0; i < 3; i++) - for (int j = 0; j < 3; j++) - gu[i][j] = phiL * phiL * gtu[i][j]; + // for (int i = 0; i < 3; i++) + // for (int j = 0; j < 3; j++) + // gu[i][j] = phiL * phiL * gtu[i][j]; - for (int i = 0; i < 3; i++) - for (int j = 0; j < 3; j++) { - double val = 0.0; - for (int k = 0; k < 3; k++) - val += gu[i][k] * Ric[k][j]; - Ricm[i][j] = val; - } + // for (int i = 0; i < 3; i++) + // for (int j = 0; j < 3; j++) { + // double val = 0.0; + // for (int k = 0; k < 3; k++) + // val += gu[i][k] * Ric[k][j]; + // Ricm[i][j] = val; + // } - for (int i = 0; i < 3; i++) - for (int j = 0; j < 3; j++) - for (int k = 0; k < 3; k++) { - dg[j][k][i] = -2.0 * dphi[i] * gt[j][k] / (phiL * SQR(phiL)) + dgt[j][k][i] / SQR(phiL); - } + // for (int i = 0; i < 3; i++) + // for (int j = 0; j < 3; j++) + // for (int k = 0; k < 3; k++) { + // dg[j][k][i] = -2.0 * dphi[i] * gt[j][k] / (phiL * SQR(phiL)) + dgt[j][k][i] / SQR(phiL); + // } - for (int i = 0; i < 3; i++) - for (int k = 0; k < 3; k++) - for (int m = 0; m < 3; m++) { - double val = 0.0; - for (int j = 0; j < 3; j++) - for (int l = 0; l < 3; l++) - val += - gu[k][j] * gu[l][m] * dg[j][l][i]; - dgu[k][m][i] = val; - } + // for (int i = 0; i < 3; i++) + // for (int k = 0; k < 3; k++) + // for (int m = 0; m < 3; m++) { + // double val = 0.0; + // for (int j = 0; j < 3; j++) + // for (int l = 0; l < 3; l++) + // val += - gu[k][j] * gu[l][m] * dg[j][l][i]; + // dgu[k][m][i] = val; + // } - for (int i = 0; i < 3; i++) - for (int j = 0; j < 3; j++) - for (int k = 0; k < 3; k++) - for (int l = 0; l < 3; l++) { - d2g[k][l][i][j] = 6.0 * gt[k][l] * dphi[i] * dphi[j] / SQR(SQR(phiL)) - 2.0 * gt[k][l] * d2phi[i][j] / (phiL * SQR(phiL)) - - 2.0 * dgt[k][l][i] * dphi[j] / (phiL * SQR(phiL)) - 2.0 * dgt[k][l][j] * dphi[i] / (phiL * SQR(phiL)) + d2gt[k][l][i][j] / SQR(phiL); - } + // for (int i = 0; i < 3; i++) + // for (int j = 0; j < 3; j++) + // for (int k = 0; k < 3; k++) + // for (int l = 0; l < 3; l++) { + // d2g[k][l][i][j] = 6.0 * gt[k][l] * dphi[i] * dphi[j] / SQR(SQR(phiL)) - 2.0 * gt[k][l] * d2phi[i][j] / (phiL * SQR(phiL)) + // - 2.0 * dgt[k][l][i] * dphi[j] / (phiL * SQR(phiL)) - 2.0 * dgt[k][l][j] * dphi[i] / (phiL * SQR(phiL)) + d2gt[k][l][i][j] / SQR(phiL); + // } - for (int i = 0; i < 3; i++) - for (int j = 0; j < 3; j++) - for (int k = 0; k < 3; k++) - for (int l = 0; l < 3; l++) { - double val = 0.0; - for (int m = 0; m < 3; m++) { - val += dgu[k][m][i] * (dg[l][m][j] + dg[j][m][l] - dg[j][l][m]) + gu[k][m] * (d2g[l][m][i][j] + d2g[j][m][l][i] - d2g[j][l][m][i]); - } - dG[k][j][l][i] = 0.5 * val; - } + // for (int i = 0; i < 3; i++) + // for (int j = 0; j < 3; j++) + // for (int k = 0; k < 3; k++) + // for (int l = 0; l < 3; l++) { + // double val = 0.0; + // for (int m = 0; m < 3; m++) { + // val += dgu[k][m][i] * (dg[l][m][j] + dg[j][m][l] - dg[j][l][m]) + gu[k][m] * (d2g[l][m][i][j] + d2g[j][m][l][i] - d2g[j][l][m][i]); + // } + // dG[k][j][l][i] = 0.5 * val; + // } - for (int i = 0; i < 3; i++) - for (int j = 0; j < 3; j++) - for (int k = 0; k < 3; k++) { - d2beta_term[k] += gu[i][j] * d2beta[k][i][j]; + // for (int i = 0; i < 3; i++) + // for (int j = 0; j < 3; j++) + // for (int k = 0; k < 3; k++) { + // d2beta_term[k] += gu[i][j] * d2beta[k][i][j]; - for (int l = 0; l < 3; l++) { - d2beta_term[k] += gu[i][j] * (G[k][j][l] * dbeta[l][i] + dG[k][j][l][i] * betaval[l] + G[k][i][l] * dbeta[l][j] - G[l][i][j] * dbeta[k][l]); + // for (int l = 0; l < 3; l++) { + // d2beta_term[k] += gu[i][j] * (G[k][j][l] * dbeta[l][i] + dG[k][j][l][i] * betaval[l] + G[k][i][l] * dbeta[l][j] - G[l][i][j] * dbeta[k][l]); - for (int m = 0; m < 3; m++) { - d2beta_term[k] += gu[i][j] * (G[k][i][l] * G[l][j][m] * betaval[m] - G[l][i][j] * G[k][l][m] * betaval[m]); - } - } - } + // for (int m = 0; m < 3; m++) { + // d2beta_term[k] += gu[i][j] * (G[k][i][l] * G[l][j][m] * betaval[m] - G[l][i][j] * G[k][l][m] * betaval[m]); + // } + // } + // } - ric_term = 0.0; - for (int i = 0; i < 3; i++) - for (int j = 0; j < 3; j++) - ric_term += betaval[j] * Ricm[i][j] * dalpha[i]; + // ric_term = 0.0; + // for (int i = 0; i < 3; i++) + // for (int j = 0; j < 3; j++) + // ric_term += betaval[j] * Ricm[i][j] * dalpha[i]; - beta_term_val = -(d2beta_term[0] * dalpha[0] + d2beta_term[2] * dalpha[2]) - ric_term; - } + // beta_term_val = -(d2beta_term[0] * dalpha[0] + d2beta_term[2] * dalpha[2]) - ric_term; + //} /* Copy local copies back to grid functions */ vec_store_partial_prepare(i,vecimin,vecimax); @@ -2606,7 +2606,7 @@ static void eval_timediff_body(const cGH* restrict const cctkGH, const int dir, vec_store_nta_partial(Kdot22[index],Kdot22L); vec_store_nta_partial(Kdot23[index],Kdot23L); vec_store_nta_partial(Kdot33[index],Kdot33L); - vec_store_nta_partial(beta_term[index], beta_term_val); + //vec_store_nta_partial(beta_term[index], beta_term_val); vec_store_nta_partial(Xtdot1[index],Xtdot1L); vec_store_nta_partial(Xtdot2[index],Xtdot2L); vec_store_nta_partial(Xtdot3[index],Xtdot3L); diff --git a/cartoon_constr.c b/cartoon_constr.c index 91b3445..6f87960 100644 --- a/cartoon_constr.c +++ b/cartoon_constr.c @@ -44,7 +44,7 @@ PDstandardNth23gt33 = kzero; PDstandardNth2phi = kzero; - PDstandardNth22phi = kifthen(origin, PDstandardNth11phi, kmul(PDstandardNth1phi, vinx)); + PDstandardNth22phi = kifthen(origin, PDstandardNth11phi, kmul(PDstandardNth1phi, xinv)); PDstandardNth12phi = kzero; PDstandardNth23phi = kzero; PDstandardNth2Xt1 = kzero; -- cgit v1.2.3