diff options
Diffstat (limited to 'ML_BSSN/src/ML_BSSN_RHS.c')
| -rw-r--r-- | ML_BSSN/src/ML_BSSN_RHS.c | 13 |
1 files changed, 8 insertions, 5 deletions
diff --git a/ML_BSSN/src/ML_BSSN_RHS.c b/ML_BSSN/src/ML_BSSN_RHS.c index 241c214..70d492b 100644 --- a/ML_BSSN/src/ML_BSSN_RHS.c +++ b/ML_BSSN/src/ML_BSSN_RHS.c @@ -109,6 +109,7 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons // CCTK_REAL dir1 = INITVALUE, dir2 = INITVALUE, dir3 = INITVALUE; // CCTK_REAL e4phi = INITVALUE; // CCTK_REAL em4phi = INITVALUE; + // CCTK_REAL eta = INITVALUE; // CCTK_REAL fac1 = INITVALUE, fac2 = INITVALUE; // CCTK_REAL g11 = INITVALUE; // CCTK_REAL G111 = INITVALUE, G112 = INITVALUE, G113 = INITVALUE; @@ -144,7 +145,6 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons // CCTK_REAL At22L = INITVALUE, At22rhsL = INITVALUE, At23L = INITVALUE, At23rhsL = INITVALUE, At33L = INITVALUE, At33rhsL = INITVALUE; // CCTK_REAL B1L = INITVALUE, B1rhsL = INITVALUE, B2L = INITVALUE, B2rhsL = INITVALUE, B3L = INITVALUE, B3rhsL = INITVALUE; // CCTK_REAL beta1L = INITVALUE, beta1rhsL = INITVALUE, beta2L = INITVALUE, beta2rhsL = INITVALUE, beta3L = INITVALUE, beta3rhsL = INITVALUE; - // CCTK_REAL etaL = INITVALUE; // CCTK_REAL eTttL = INITVALUE; // CCTK_REAL eTtxL = INITVALUE; // CCTK_REAL eTtyL = INITVALUE; @@ -158,6 +158,7 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons // CCTK_REAL gt11L = INITVALUE, gt11rhsL = INITVALUE, gt12L = INITVALUE, gt12rhsL = INITVALUE, gt13L = INITVALUE, gt13rhsL = INITVALUE; // CCTK_REAL gt22L = INITVALUE, gt22rhsL = INITVALUE, gt23L = INITVALUE, gt23rhsL = INITVALUE, gt33L = INITVALUE, gt33rhsL = INITVALUE; // CCTK_REAL phiL = INITVALUE, phirhsL = INITVALUE; + // CCTK_REAL rL = INITVALUE; // CCTK_REAL trKL = INITVALUE, trKrhsL = INITVALUE; // CCTK_REAL Xt1L = INITVALUE, Xt1rhsL = INITVALUE, Xt2L = INITVALUE, Xt2rhsL = INITVALUE, Xt3L = INITVALUE, Xt3rhsL = INITVALUE; /* Declare precomputed derivatives*/ @@ -290,7 +291,6 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons CCTK_REAL const beta1L = beta1[index]; CCTK_REAL const beta2L = beta2[index]; CCTK_REAL const beta3L = beta3[index]; - CCTK_REAL const etaL = eta[index]; CCTK_REAL const eTttL = (*stress_energy_state) ? (eTtt[index]) : 0.0; CCTK_REAL const eTtxL = (*stress_energy_state) ? (eTtx[index]) : 0.0; CCTK_REAL const eTtyL = (*stress_energy_state) ? (eTty[index]) : 0.0; @@ -308,6 +308,7 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons CCTK_REAL const gt23L = gt23[index]; CCTK_REAL const gt33L = gt33[index]; CCTK_REAL const phiL = phi[index]; + CCTK_REAL const rL = r[index]; CCTK_REAL const trKL = trK[index]; CCTK_REAL const Xt1L = Xt1[index]; CCTK_REAL const Xt2L = Xt2[index]; @@ -1088,6 +1089,8 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons CCTK_REAL const ArhsL = (-1 + LapseAdvectionCoeff)*(AL*AlphaDriver - trKrhsL); + CCTK_REAL const eta = BetaDriver*IfThen(rL > SpatialBetaDriverRadius,SpatialBetaDriverRadius*INV(rL),1); + CCTK_REAL const beta1rhsL = (PDupwindNth1(beta1, i, j, k)*beta1L + PDupwindNth2(beta1, i, j, k)*beta2L + PDupwindNth3(beta1, i, j, k)*beta3L)*ShiftAdvectionCoeff + B1L*ShiftGammaCoeff; @@ -1097,15 +1100,15 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons CCTK_REAL const beta3rhsL = (PDupwindNth1(beta3, i, j, k)*beta1L + PDupwindNth2(beta3, i, j, k)*beta2L + PDupwindNth3(beta3, i, j, k)*beta3L)*ShiftAdvectionCoeff + B3L*ShiftGammaCoeff; - CCTK_REAL const B1rhsL = -(B1L*etaL) + ((PDupwindNth1(B1, i, j, k) - PDupwindNth1(Xt1, i, j, k))*beta1L + + CCTK_REAL const B1rhsL = -(B1L*eta) + ((PDupwindNth1(B1, i, j, k) - PDupwindNth1(Xt1, i, j, k))*beta1L + (PDupwindNth2(B1, i, j, k) - PDupwindNth2(Xt1, i, j, k))*beta2L + (PDupwindNth3(B1, i, j, k) - PDupwindNth3(Xt1, i, j, k))*beta3L)*ShiftAdvectionCoeff + Xt1rhsL; - CCTK_REAL const B2rhsL = -(B2L*etaL) + ((PDupwindNth1(B2, i, j, k) - PDupwindNth1(Xt2, i, j, k))*beta1L + + CCTK_REAL const B2rhsL = -(B2L*eta) + ((PDupwindNth1(B2, i, j, k) - PDupwindNth1(Xt2, i, j, k))*beta1L + (PDupwindNth2(B2, i, j, k) - PDupwindNth2(Xt2, i, j, k))*beta2L + (PDupwindNth3(B2, i, j, k) - PDupwindNth3(Xt2, i, j, k))*beta3L)*ShiftAdvectionCoeff + Xt2rhsL; - CCTK_REAL const B3rhsL = -(B3L*etaL) + ((PDupwindNth1(B3, i, j, k) - PDupwindNth1(Xt3, i, j, k))*beta1L + + CCTK_REAL const B3rhsL = -(B3L*eta) + ((PDupwindNth1(B3, i, j, k) - PDupwindNth1(Xt3, i, j, k))*beta1L + (PDupwindNth2(B3, i, j, k) - PDupwindNth2(Xt3, i, j, k))*beta2L + (PDupwindNth3(B3, i, j, k) - PDupwindNth3(Xt3, i, j, k))*beta3L)*ShiftAdvectionCoeff + Xt3rhsL; |