diff options
33 files changed, 3324 insertions, 301 deletions
diff --git a/ML_BSSN/param.ccl b/ML_BSSN/param.ccl index 2324c01..d105c58 100644 --- a/ML_BSSN/param.ccl +++ b/ML_BSSN/param.ccl @@ -213,7 +213,13 @@ CCTK_INT ML_BSSN_convertFromADMBaseGamma_calc_every "ML_BSSN_convertFromADMBaseG } 1 restricted: -CCTK_INT ML_BSSN_setBetaDriver_calc_every "ML_BSSN_setBetaDriver_calc_every" +CCTK_INT ML_BSSN_setBetaDriverConstant_calc_every "ML_BSSN_setBetaDriverConstant_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_setBetaDriverSpatial_calc_every "ML_BSSN_setBetaDriverSpatial_calc_every" { *:* :: "" } 1 @@ -321,7 +327,13 @@ CCTK_INT ML_BSSN_convertFromADMBaseGamma_calc_offset "ML_BSSN_convertFromADMBase } 0 restricted: -CCTK_INT ML_BSSN_setBetaDriver_calc_offset "ML_BSSN_setBetaDriver_calc_offset" +CCTK_INT ML_BSSN_setBetaDriverConstant_calc_offset "ML_BSSN_setBetaDriverConstant_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_setBetaDriverSpatial_calc_offset "ML_BSSN_setBetaDriverSpatial_calc_offset" { *:* :: "" } 0 diff --git a/ML_BSSN/schedule.ccl b/ML_BSSN/schedule.ccl index d8422cd..bd0f54d 100644 --- a/ML_BSSN/schedule.ccl +++ b/ML_BSSN/schedule.ccl @@ -368,12 +368,21 @@ if (CCTK_EQUALS(my_initial_data, "ADMBase")) } +if (CCTK_EQUALS(UseSpatialBetaDriver, "no")) +{ + schedule ML_BSSN_setBetaDriverConstant IN ML_BSSN_InitEta + { + LANG: C + } "ML_BSSN_setBetaDriverConstant" +} + + if (CCTK_EQUALS(UseSpatialBetaDriver, "yes")) { - schedule ML_BSSN_setBetaDriver AT initial AFTER ADMBase_PostInitial AFTER ML_BSSN_convertFromADMBase + schedule ML_BSSN_setBetaDriverSpatial IN ML_BSSN_InitEta { LANG: C - } "ML_BSSN_setBetaDriver" + } "ML_BSSN_setBetaDriverSpatial" } schedule ML_BSSN_RHS IN NoSuchGroup diff --git a/ML_BSSN/src/ML_BSSN_Minkowski.c b/ML_BSSN/src/ML_BSSN_Minkowski.c index c467955..29622aa 100644 --- a/ML_BSSN/src/ML_BSSN_Minkowski.c +++ b/ML_BSSN/src/ML_BSSN_Minkowski.c @@ -106,7 +106,6 @@ void ML_BSSN_Minkowski_Body(cGH const * restrict const cctkGH, int const dir, in // CCTK_REAL At11L = INITVALUE, At12L = INITVALUE, At13L = INITVALUE, At22L = INITVALUE, At23L = INITVALUE, At33L = INITVALUE; // CCTK_REAL B1L = INITVALUE, B2L = INITVALUE, B3L = INITVALUE; // CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE; - // CCTK_REAL etaL = INITVALUE; // CCTK_REAL gt11L = INITVALUE, gt12L = INITVALUE, gt13L = INITVALUE, gt22L = INITVALUE, gt23L = INITVALUE, gt33L = INITVALUE; // CCTK_REAL phiL = INITVALUE; // CCTK_REAL trKL = INITVALUE; @@ -176,8 +175,6 @@ void ML_BSSN_Minkowski_Body(cGH const * restrict const cctkGH, int const dir, in CCTK_REAL const B3L = 0; - CCTK_REAL const etaL = BetaDriver; - /* Copy local copies back to grid functions */ A[index] = AL; @@ -194,7 +191,6 @@ void ML_BSSN_Minkowski_Body(cGH const * restrict const cctkGH, int const dir, in beta1[index] = beta1L; beta2[index] = beta2L; beta3[index] = beta3L; - eta[index] = etaL; gt11[index] = gt11L; gt12[index] = gt12L; gt13[index] = gt13L; diff --git a/ML_BSSN/src/ML_BSSN_RHS.c b/ML_BSSN/src/ML_BSSN_RHS.c index f207bcc..ad6e249 100644 --- a/ML_BSSN/src/ML_BSSN_RHS.c +++ b/ML_BSSN/src/ML_BSSN_RHS.c @@ -132,8 +132,6 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons // CCTK_REAL Rphi11 = INITVALUE, Rphi12 = INITVALUE, Rphi13 = INITVALUE, Rphi22 = INITVALUE, Rphi23 = INITVALUE, Rphi33 = 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 trAts = INITVALUE; // CCTK_REAL trS = INITVALUE; // CCTK_REAL Xtn1 = INITVALUE, Xtn2 = INITVALUE, Xtn3 = INITVALUE; @@ -922,37 +920,17 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons CCTK_REAL const R33 = Rphi33 + Rt33; - 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 S3 = (-eTtzL + beta1L*eTxzL + beta2L*eTyzL + beta3L*eTzzL)*INV(alphaL); - 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 trS = gu11*T11 + gu22*T22 + 2*(gu21*T12 + gu31*T13 + gu32*T23) + gu33*T33; + CCTK_REAL const trS = eTxxL*gu11 + eTyyL*gu22 + 2*(eTxyL*gu21 + eTxzL*gu31 + eTyzL*gu32) + eTzzL*gu33; CCTK_REAL const phirhsL = PDupwindNth1(phi, i, j, k)*beta1L + PDupwindNth2(phi, i, j, k)*beta2L + PDupwindNth3(phi, i, j, k)*beta3L + (PDstandardNth1beta1 + PDstandardNth2beta2 + PDstandardNth3beta3)* @@ -1062,7 +1040,7 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons At11L*(1.333333333333333333333333333333333333333*PDstandardNth1beta1 - 0.6666666666666666666666666666666666666667*(PDstandardNth2beta2 + PDstandardNth3beta3) + alphaL*trKL) + em4phi*(Ats11 - 0.3333333333333333333333333333333333333333*g11*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T11 + 8.377580409572781969233715688745341024526*g11*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxxL + 8.377580409572781969233715688745341024526*g11*trS)); CCTK_REAL const At12rhsL = -2.*alphaL*(At11L*Atm12 + At12L*Atm22 + At13L*Atm32) + PDupwindNth1(At12, i, j, k)*beta1L + PDupwindNth2(At12, i, j, k)*beta2L + PDupwindNth3(At12, i, j, k)*beta3L + At22L*PDstandardNth1beta2 + @@ -1070,7 +1048,7 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons At12L*(0.3333333333333333333333333333333333333333*(PDstandardNth1beta1 + PDstandardNth2beta2) - 0.6666666666666666666666666666666666666667*PDstandardNth3beta3 + alphaL*trKL) + em4phi*(Ats12 - 0.3333333333333333333333333333333333333333*g12*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T12 + 8.377580409572781969233715688745341024526*g12*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxyL + 8.377580409572781969233715688745341024526*g12*trS)); CCTK_REAL const At13rhsL = -2.*alphaL*(At11L*Atm13 + At12L*Atm23 + At13L*Atm33) + PDupwindNth1(At13, i, j, k)*beta1L + PDupwindNth2(At13, i, j, k)*beta2L + PDupwindNth3(At13, i, j, k)*beta3L + At23L*PDstandardNth1beta2 + @@ -1078,7 +1056,7 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons At13L*(-0.6666666666666666666666666666666666666667*PDstandardNth2beta2 + 0.3333333333333333333333333333333333333333*(PDstandardNth1beta1 + PDstandardNth3beta3) + alphaL*trKL) + em4phi*(Ats13 - 0.3333333333333333333333333333333333333333*g13*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T13 + 8.377580409572781969233715688745341024526*g13*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxzL + 8.377580409572781969233715688745341024526*g13*trS)); CCTK_REAL const At22rhsL = -2.*alphaL*(At12L*Atm12 + At22L*Atm22 + At23L*Atm32) + PDupwindNth1(At22, i, j, k)*beta1L + PDupwindNth2(At22, i, j, k)*beta2L + PDupwindNth3(At22, i, j, k)*beta3L + @@ -1086,7 +1064,7 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons At22L*(1.333333333333333333333333333333333333333*PDstandardNth2beta2 - 0.6666666666666666666666666666666666666667*(PDstandardNth1beta1 + PDstandardNth3beta3) + alphaL*trKL) + em4phi*(Ats22 - 0.3333333333333333333333333333333333333333*g22*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T22 + 8.377580409572781969233715688745341024526*g22*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTyyL + 8.377580409572781969233715688745341024526*g22*trS)); CCTK_REAL const At23rhsL = -2.*alphaL*(At12L*Atm13 + At22L*Atm23 + At23L*Atm33) + PDupwindNth1(At23, i, j, k)*beta1L + PDupwindNth2(At23, i, j, k)*beta2L + PDupwindNth3(At23, i, j, k)*beta3L + At13L*PDstandardNth2beta1 + @@ -1094,7 +1072,7 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons At23L*(-0.6666666666666666666666666666666666666667*PDstandardNth1beta1 + 0.3333333333333333333333333333333333333333*(PDstandardNth2beta2 + PDstandardNth3beta3) + alphaL*trKL) + em4phi*(Ats23 - 0.3333333333333333333333333333333333333333*g23*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T23 + 8.377580409572781969233715688745341024526*g23*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTyzL + 8.377580409572781969233715688745341024526*g23*trS)); CCTK_REAL const At33rhsL = -2.*alphaL*(At13L*Atm13 + At23L*Atm23 + At33L*Atm33) + PDupwindNth1(At33, i, j, k)*beta1L + PDupwindNth2(At33, i, j, k)*beta2L + PDupwindNth3(At33, i, j, k)*beta3L + @@ -1102,7 +1080,7 @@ void ML_BSSN_RHS_Body(cGH const * restrict const cctkGH, int const dir, int cons At33L*(-0.6666666666666666666666666666666666666667*(PDstandardNth1beta1 + PDstandardNth2beta2) + 1.333333333333333333333333333333333333333*PDstandardNth3beta3 + alphaL*trKL) + em4phi*(Ats33 - 0.3333333333333333333333333333333333333333*g33*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T33 + 8.377580409572781969233715688745341024526*g33*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTzzL + 8.377580409572781969233715688745341024526*g33*trS)); CCTK_REAL const alpharhsL = (PDupwindNth1(alpha, i, j, k)*beta1L + PDupwindNth2(alpha, i, j, k)*beta2L + PDupwindNth3(alpha, i, j, k)*beta3L)*LapseAdvectionCoeff + diff --git a/ML_BSSN/src/ML_BSSN_RHS1.c b/ML_BSSN/src/ML_BSSN_RHS1.c index 71e93dc..63d25bb 100644 --- a/ML_BSSN/src/ML_BSSN_RHS1.c +++ b/ML_BSSN/src/ML_BSSN_RHS1.c @@ -111,8 +111,6 @@ void ML_BSSN_RHS1_Body(cGH const * restrict const cctkGH, int const dir, int con // CCTK_REAL Gt311 = INITVALUE, Gt312 = INITVALUE, Gt313 = INITVALUE, Gt322 = INITVALUE, Gt323 = INITVALUE, Gt333 = INITVALUE; // CCTK_REAL gtu11 = INITVALUE, gtu21 = INITVALUE, gtu22 = INITVALUE, gtu31 = INITVALUE, gtu32 = INITVALUE, gtu33 = INITVALUE; // CCTK_REAL S1 = INITVALUE, S2 = INITVALUE, S3 = INITVALUE; - // CCTK_REAL T01 = INITVALUE, T02 = INITVALUE, T03 = INITVALUE, T11 = INITVALUE, T12 = INITVALUE, T13 = INITVALUE; - // CCTK_REAL T22 = INITVALUE, T23 = INITVALUE, T33 = INITVALUE; // CCTK_REAL Xtn1 = INITVALUE, Xtn2 = INITVALUE, Xtn3 = INITVALUE; /* Declare local copies of grid functions */ @@ -410,29 +408,11 @@ void ML_BSSN_RHS1_Body(cGH const * restrict const cctkGH, int const dir, int con CCTK_REAL const Atu33 = Atm31*gtu31 + Atm32*gtu32 + Atm33*gtu33; - 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 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 phirhsL = PDupwindNth1(phi, i, j, k)*beta1L + PDupwindNth2(phi, i, j, k)*beta2L + PDupwindNth3(phi, i, j, k)*beta3L + (PDstandardNth1beta1 + PDstandardNth2beta2 + PDstandardNth3beta3)* diff --git a/ML_BSSN/src/ML_BSSN_RHS2.c b/ML_BSSN/src/ML_BSSN_RHS2.c index a8ff2dc..3c6b439 100644 --- a/ML_BSSN/src/ML_BSSN_RHS2.c +++ b/ML_BSSN/src/ML_BSSN_RHS2.c @@ -138,8 +138,6 @@ void ML_BSSN_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int con // CCTK_REAL rho = INITVALUE; // CCTK_REAL Rphi11 = INITVALUE, Rphi12 = INITVALUE, Rphi13 = INITVALUE, Rphi22 = INITVALUE, Rphi23 = INITVALUE, Rphi33 = INITVALUE; // CCTK_REAL Rt11 = INITVALUE, Rt12 = INITVALUE, Rt13 = INITVALUE, Rt22 = INITVALUE, Rt23 = INITVALUE, Rt33 = 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 trAts = INITVALUE; // CCTK_REAL trS = INITVALUE; // CCTK_REAL Xtn1 = INITVALUE, Xtn2 = INITVALUE, Xtn3 = INITVALUE; @@ -742,31 +740,11 @@ void ML_BSSN_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int con CCTK_REAL const R33 = Rphi33 + Rt33; - 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 T02 = eTtyL; - - 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 trS = gu11*T11 + gu22*T22 + 2*(gu21*T12 + gu31*T13 + gu32*T23) + gu33*T33; + CCTK_REAL const trS = eTxxL*gu11 + eTyyL*gu22 + 2*(eTxyL*gu21 + eTxzL*gu31 + eTyzL*gu32) + eTzzL*gu33; CCTK_REAL const trKrhsL = PDupwindNth1(trK, i, j, k)*beta1L + PDupwindNth2(trK, i, j, k)*beta2L + PDupwindNth3(trK, i, j, k)*beta3L + (G111*gu11 + G122*gu22 + 2.*(G112*gu21 + G113*gu31 + G123*gu32) + G133*gu33)* @@ -804,7 +782,7 @@ void ML_BSSN_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int con At11L*(1.333333333333333333333333333333333333333*PDstandardNth1beta1 - 0.6666666666666666666666666666666666666667*(PDstandardNth2beta2 + PDstandardNth3beta3) + alphaL*trKL) + em4phi*(Ats11 - 0.3333333333333333333333333333333333333333*g11*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T11 + 8.377580409572781969233715688745341024526*g11*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxxL + 8.377580409572781969233715688745341024526*g11*trS)); CCTK_REAL const At12rhsL = -2.*alphaL*(At11L*Atm12 + At12L*Atm22 + At13L*Atm32) + PDupwindNth1(At12, i, j, k)*beta1L + PDupwindNth2(At12, i, j, k)*beta2L + PDupwindNth3(At12, i, j, k)*beta3L + At22L*PDstandardNth1beta2 + @@ -812,7 +790,7 @@ void ML_BSSN_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int con At12L*(0.3333333333333333333333333333333333333333*(PDstandardNth1beta1 + PDstandardNth2beta2) - 0.6666666666666666666666666666666666666667*PDstandardNth3beta3 + alphaL*trKL) + em4phi*(Ats12 - 0.3333333333333333333333333333333333333333*g12*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T12 + 8.377580409572781969233715688745341024526*g12*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxyL + 8.377580409572781969233715688745341024526*g12*trS)); CCTK_REAL const At13rhsL = -2.*alphaL*(At11L*Atm13 + At12L*Atm23 + At13L*Atm33) + PDupwindNth1(At13, i, j, k)*beta1L + PDupwindNth2(At13, i, j, k)*beta2L + PDupwindNth3(At13, i, j, k)*beta3L + At23L*PDstandardNth1beta2 + @@ -820,7 +798,7 @@ void ML_BSSN_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int con At13L*(-0.6666666666666666666666666666666666666667*PDstandardNth2beta2 + 0.3333333333333333333333333333333333333333*(PDstandardNth1beta1 + PDstandardNth3beta3) + alphaL*trKL) + em4phi*(Ats13 - 0.3333333333333333333333333333333333333333*g13*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T13 + 8.377580409572781969233715688745341024526*g13*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxzL + 8.377580409572781969233715688745341024526*g13*trS)); CCTK_REAL const At22rhsL = -2.*alphaL*(At12L*Atm12 + At22L*Atm22 + At23L*Atm32) + PDupwindNth1(At22, i, j, k)*beta1L + PDupwindNth2(At22, i, j, k)*beta2L + PDupwindNth3(At22, i, j, k)*beta3L + @@ -828,7 +806,7 @@ void ML_BSSN_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int con At22L*(1.333333333333333333333333333333333333333*PDstandardNth2beta2 - 0.6666666666666666666666666666666666666667*(PDstandardNth1beta1 + PDstandardNth3beta3) + alphaL*trKL) + em4phi*(Ats22 - 0.3333333333333333333333333333333333333333*g22*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T22 + 8.377580409572781969233715688745341024526*g22*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTyyL + 8.377580409572781969233715688745341024526*g22*trS)); CCTK_REAL const At23rhsL = -2.*alphaL*(At12L*Atm13 + At22L*Atm23 + At23L*Atm33) + PDupwindNth1(At23, i, j, k)*beta1L + PDupwindNth2(At23, i, j, k)*beta2L + PDupwindNth3(At23, i, j, k)*beta3L + At13L*PDstandardNth2beta1 + @@ -836,7 +814,7 @@ void ML_BSSN_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int con At23L*(-0.6666666666666666666666666666666666666667*PDstandardNth1beta1 + 0.3333333333333333333333333333333333333333*(PDstandardNth2beta2 + PDstandardNth3beta3) + alphaL*trKL) + em4phi*(Ats23 - 0.3333333333333333333333333333333333333333*g23*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T23 + 8.377580409572781969233715688745341024526*g23*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTyzL + 8.377580409572781969233715688745341024526*g23*trS)); CCTK_REAL const At33rhsL = -2.*alphaL*(At13L*Atm13 + At23L*Atm23 + At33L*Atm33) + PDupwindNth1(At33, i, j, k)*beta1L + PDupwindNth2(At33, i, j, k)*beta2L + PDupwindNth3(At33, i, j, k)*beta3L + @@ -844,7 +822,7 @@ void ML_BSSN_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int con At33L*(-0.6666666666666666666666666666666666666667*(PDstandardNth1beta1 + PDstandardNth2beta2) + 1.333333333333333333333333333333333333333*PDstandardNth3beta3 + alphaL*trKL) + em4phi*(Ats33 - 0.3333333333333333333333333333333333333333*g33*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T33 + 8.377580409572781969233715688745341024526*g33*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTzzL + 8.377580409572781969233715688745341024526*g33*trS)); CCTK_REAL const alpharhsL = (PDupwindNth1(alpha, i, j, k)*beta1L + PDupwindNth2(alpha, i, j, k)*beta2L + PDupwindNth3(alpha, i, j, k)*beta3L)*LapseAdvectionCoeff + diff --git a/ML_BSSN/src/ML_BSSN_constraints.c b/ML_BSSN/src/ML_BSSN_constraints.c index f63c3fd..48db1b3 100644 --- a/ML_BSSN/src/ML_BSSN_constraints.c +++ b/ML_BSSN/src/ML_BSSN_constraints.c @@ -117,8 +117,6 @@ void ML_BSSN_constraints_Body(cGH const * restrict const cctkGH, int const dir, // CCTK_REAL Rphi11 = INITVALUE, Rphi12 = INITVALUE, Rphi13 = INITVALUE, Rphi22 = INITVALUE, Rphi23 = INITVALUE, Rphi33 = 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 trR = INITVALUE; /* Declare local copies of grid functions */ @@ -793,35 +791,15 @@ void ML_BSSN_constraints_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 HL = -2.*(Atm12*Atm21 + Atm13*Atm31 + Atm23*Atm32) - 50.26548245743669181540229413247204614715*rho + trR - 1.*(SQR(Atm11) + SQR(Atm22) + SQR(Atm33)) + 0.6666666666666666666666666666666666666667*SQR(trKL); diff --git a/ML_BSSN/src/ML_BSSN_convertFromADMBase.c b/ML_BSSN/src/ML_BSSN_convertFromADMBase.c index 95c7d4c..1b93cb2 100644 --- a/ML_BSSN/src/ML_BSSN_convertFromADMBase.c +++ b/ML_BSSN/src/ML_BSSN_convertFromADMBase.c @@ -112,7 +112,6 @@ void ML_BSSN_convertFromADMBase_Body(cGH const * restrict const cctkGH, int cons // CCTK_REAL betaxL = INITVALUE; // CCTK_REAL betayL = INITVALUE; // CCTK_REAL betazL = INITVALUE; - // CCTK_REAL etaL = INITVALUE; // CCTK_REAL gt11L = INITVALUE, gt12L = INITVALUE, gt13L = INITVALUE, gt22L = INITVALUE, gt23L = INITVALUE, gt33L = INITVALUE; // CCTK_REAL gxxL = INITVALUE; // CCTK_REAL gxyL = INITVALUE; @@ -223,8 +222,6 @@ void ML_BSSN_convertFromADMBase_Body(cGH const * restrict const cctkGH, int cons CCTK_REAL const beta3L = betazL; - CCTK_REAL const etaL = BetaDriver; - /* Copy local copies back to grid functions */ alpha[index] = alphaL; @@ -237,7 +234,6 @@ void ML_BSSN_convertFromADMBase_Body(cGH const * restrict const cctkGH, int cons beta1[index] = beta1L; beta2[index] = beta2L; beta3[index] = beta3L; - eta[index] = etaL; gt11[index] = gt11L; gt12[index] = gt12L; gt13[index] = gt13L; diff --git a/ML_BSSN/src/ML_BSSN_setBetaDriver.c b/ML_BSSN/src/ML_BSSN_setBetaDriverSpatial.c index 2f6a35f..1e726e2 100644 --- a/ML_BSSN/src/ML_BSSN_setBetaDriver.c +++ b/ML_BSSN/src/ML_BSSN_setBetaDriverSpatial.c @@ -20,7 +20,7 @@ #define CUB(x) ((x) * (x) * (x)) #define QAD(x) ((x) * (x) * (x) * (x)) -void ML_BSSN_setBetaDriver_Body(cGH const * restrict const cctkGH, int const dir, int const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], int const min[3], int const max[3], int const n_subblock_gfs, CCTK_REAL * restrict const subblock_gfs[]) +void ML_BSSN_setBetaDriverSpatial_Body(cGH const * restrict const cctkGH, int const dir, int const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], int const min[3], int const max[3], int const n_subblock_gfs, CCTK_REAL * restrict const subblock_gfs[]) { DECLARE_CCTK_ARGUMENTS; DECLARE_CCTK_PARAMETERS; @@ -44,10 +44,10 @@ void ML_BSSN_setBetaDriver_Body(cGH const * restrict const cctkGH, int const dir if (verbose > 1) { - CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_BSSN_setBetaDriver_Body"); + CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_BSSN_setBetaDriverSpatial_Body"); } - if (cctk_iteration % ML_BSSN_setBetaDriver_calc_every != ML_BSSN_setBetaDriver_calc_offset) + if (cctk_iteration % ML_BSSN_setBetaDriverSpatial_calc_every != ML_BSSN_setBetaDriverSpatial_calc_offset) { return; } @@ -89,7 +89,7 @@ void ML_BSSN_setBetaDriver_Body(cGH const * restrict const cctkGH, int const dir /* Loop over the grid points */ #pragma omp parallel - LC_LOOP3 (ML_BSSN_setBetaDriver, + LC_LOOP3 (ML_BSSN_setBetaDriverSpatial, i,j,k, min[0],min[1],min[2], max[0],max[1],max[2], cctk_lsh[0],cctk_lsh[1],cctk_lsh[2]) { @@ -127,13 +127,13 @@ void ML_BSSN_setBetaDriver_Body(cGH const * restrict const cctkGH, int const dir /* Copy local copies back to subblock grid functions */ } - LC_ENDLOOP3 (ML_BSSN_setBetaDriver); + LC_ENDLOOP3 (ML_BSSN_setBetaDriverSpatial); } -void ML_BSSN_setBetaDriver(CCTK_ARGUMENTS) +void ML_BSSN_setBetaDriverSpatial(CCTK_ARGUMENTS) { DECLARE_CCTK_ARGUMENTS; DECLARE_CCTK_PARAMETERS; - GenericFD_LoopOverEverything(cctkGH, &ML_BSSN_setBetaDriver_Body); + GenericFD_LoopOverEverything(cctkGH, &ML_BSSN_setBetaDriverSpatial_Body); } diff --git a/ML_BSSN/src/make.code.defn b/ML_BSSN/src/make.code.defn index dcab216..eed6fb1 100644 --- a/ML_BSSN/src/make.code.defn +++ b/ML_BSSN/src/make.code.defn @@ -1,3 +1,3 @@ # File produced by Kranc -SRCS = Startup.c RegisterMoL.c RegisterSymmetries.c ML_BSSN_Minkowski.c ML_BSSN_convertFromADMBase.c ML_BSSN_convertFromADMBaseGamma.c ML_BSSN_setBetaDriver.c ML_BSSN_RHS.c ML_BSSN_RHS1.c ML_BSSN_RHS2.c ML_BSSN_RHSStaticBoundary.c ML_BSSN_RHSRadiativeBoundary.c ML_BSSN_enforce.c ML_BSSN_enforce2.c ML_BSSN_boundary.c ML_BSSN_convertToADMBase.c ML_BSSN_convertToADMBaseDtLapseShift.c ML_BSSN_convertToADMBaseDtLapseShiftBoundary.c ML_BSSN_convertToADMBaseFakeDtLapseShift.c ML_BSSN_constraints.c ML_BSSN_constraints_boundary.c Boundaries.c +SRCS = Startup.c RegisterMoL.c RegisterSymmetries.c ML_BSSN_Minkowski.c ML_BSSN_convertFromADMBase.c ML_BSSN_convertFromADMBaseGamma.c ML_BSSN_setBetaDriverConstant.c ML_BSSN_setBetaDriverSpatial.c ML_BSSN_RHS.c ML_BSSN_RHS1.c ML_BSSN_RHS2.c ML_BSSN_RHSStaticBoundary.c ML_BSSN_RHSRadiativeBoundary.c ML_BSSN_enforce.c ML_BSSN_enforce2.c ML_BSSN_boundary.c ML_BSSN_convertToADMBase.c ML_BSSN_convertToADMBaseDtLapseShift.c ML_BSSN_convertToADMBaseDtLapseShiftBoundary.c ML_BSSN_convertToADMBaseFakeDtLapseShift.c ML_BSSN_constraints.c ML_BSSN_constraints_boundary.c Boundaries.c diff --git a/ML_BSSN_Helper/schedule.ccl b/ML_BSSN_Helper/schedule.ccl index 48ebedc..6b6c4e8 100644 --- a/ML_BSSN_Helper/schedule.ccl +++ b/ML_BSSN_Helper/schedule.ccl @@ -56,6 +56,20 @@ if (CCTK_EQUALS (evolution_method, "ML_BSSN")) { + SCHEDULE GROUP ML_BSSN_InitEta AT initial + { + } "Initialize BetaDriver" + + SCHEDULE GROUP ML_BSSN_InitEta AT postregrid + { + } "Initialize BetaDriver" + + SCHEDULE GROUP ML_BSSN_InitEta AT post_recover_variables + { + } "Initialize BetaDriver" + + + #SCHEDULE GROUP ML_BSSN_evolCalcGroup AT postinitial AFTER MoL_PostStep #{ #} "Calculate BSSN RHS" @@ -144,7 +158,11 @@ if (CCTK_EQUALS (evolution_method, "ML_BSSN")) { TRIGGERS: ML_BSSN::ML_mom } "Calculate ADM variables" } - + + SCHEDULE GROUP ML_BSSN_convertToADMBaseGroupWrapper AT CCTK_POST_RECOVER_VARIABLES + { + } "Calculate ADM variables" + SCHEDULE ML_BSSN_SelectBCsADMBase IN ML_BSSN_convertToADMBaseGroupWrapper AFTER ML_BSSN_convertToADMBaseGroup { LANG: C @@ -162,5 +180,4 @@ if (CCTK_EQUALS (evolution_method, "ML_BSSN")) { TRIGGERS: ML_BSSN::ML_Ham TRIGGERS: ML_BSSN::ML_mom } "Calculate BSSN constraints" - } diff --git a/ML_BSSN_M/param.ccl b/ML_BSSN_M/param.ccl new file mode 100644 index 0000000..588df7b --- /dev/null +++ b/ML_BSSN_M/param.ccl @@ -0,0 +1,1408 @@ +# File produced by Kranc + + +shares: ADMBase + + +EXTENDS CCTK_KEYWORD evolution_method "evolution_method" +{ + ML_BSSN_M :: "" +} + + +EXTENDS CCTK_KEYWORD lapse_evolution_method "lapse_evolution_method" +{ + ML_BSSN_M :: "" +} + + +EXTENDS CCTK_KEYWORD shift_evolution_method "shift_evolution_method" +{ + ML_BSSN_M :: "" +} + + +EXTENDS CCTK_KEYWORD dtlapse_evolution_method "dtlapse_evolution_method" +{ + ML_BSSN_M :: "" +} + + +EXTENDS CCTK_KEYWORD dtshift_evolution_method "dtshift_evolution_method" +{ + ML_BSSN_M :: "" +} + + + +shares: GenericFD + +USES CCTK_INT stencil_width +USES CCTK_INT stencil_width_x +USES CCTK_INT stencil_width_y +USES CCTK_INT stencil_width_z +USES CCTK_INT boundary_width + + +shares: MethodOfLines + +USES CCTK_INT MoL_Num_Evolved_Vars + +restricted: +CCTK_INT verbose "verbose" +{ + *:* :: "" +} 0 + +restricted: +CCTK_REAL harmonicF "d/dt alpha = - f alpha^n K (harmonic=1, 1+log=2)" +{ + "*:*" :: "" +} 1 + +restricted: +CCTK_REAL AlphaDriver "AlphaDriver" +{ + "*:*" :: "" +} 0 + +restricted: +CCTK_REAL ShiftGammaCoeff "ShiftGammaCoeff" +{ + "*:*" :: "" +} 0 + +restricted: +CCTK_REAL BetaDriver "BetaDriver" +{ + "*:*" :: "" +} 0 + +restricted: +CCTK_REAL LapseAdvectionCoeff "Factor in front of the shift advection terms in 1+log" +{ + "*:*" :: "" +} 1 + +restricted: +CCTK_REAL ShiftAdvectionCoeff "Factor in front of the shift advection terms in gamma driver" +{ + "*:*" :: "" +} 1 + +restricted: +CCTK_REAL MinimumLapse "Minimum value of the lapse function" +{ + "*:*" :: "" +} -1 + +restricted: +CCTK_REAL SpatialBetaDriverRadius "Radius at which the BetaDriver starts to be reduced" +{ + "*:*" :: "" +} 1000000000000 + +restricted: +CCTK_INT harmonicN "d/dt alpha = - f alpha^n K (harmonic=2, 1+log=1)" +{ + *:* :: "" +} 2 + +restricted: +CCTK_INT ShiftAlphaPower "ShiftAlphaPower" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT conformalMethod "Treatment of conformal factor" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT useMatter "Add matter terms" +{ + *:* :: "" +} 1 + +private: +KEYWORD my_initial_data "my_initial_data" +{ + "ADMBase" :: "ADMBase" + "Minkowski" :: "Minkowski" +} "ADMBase" + +restricted: +KEYWORD my_initial_boundary_condition "my_initial_boundary_condition" +{ + "none" :: "none" +} "none" + +restricted: +KEYWORD my_rhs_boundary_condition "my_rhs_boundary_condition" +{ + "none" :: "none" + "static" :: "static" + "radiative" :: "radiative" +} "none" + +private: +KEYWORD my_boundary_condition "my_boundary_condition" +{ + "none" :: "none" + "Minkowski" :: "Minkowski" +} "none" + +restricted: +KEYWORD calculate_ADMBase_variables_at "calculate_ADMBase_variables_at" +{ + "MoL_PostStep" :: "MoL_PostStep" + "CCTK_EVOL" :: "CCTK_EVOL" + "CCTK_ANALYSIS" :: "CCTK_ANALYSIS" +} "MoL_PostStep" + +restricted: +KEYWORD UseSpatialBetaDriver "UseSpatialBetaDriver" +{ + "no" :: "no" + "yes" :: "yes" +} "no" + +private: +KEYWORD dt_lapse_shift_method "Treatment of ADMBase dtlapse and dtshift" +{ + "correct" :: "correct" + "noLapseShiftAdvection" :: "noLapseShiftAdvection" +} "correct" + +restricted: +CCTK_INT ML_BSSN_M_MaxNumEvolvedVars "Number of evolved variables used by this thorn" ACCUMULATOR-BASE=MethodofLines::MoL_Num_Evolved_Vars +{ + 25:25 :: "Number of evolved variables used by this thorn" +} 25 + +restricted: +CCTK_INT timelevels "Number of active timelevels" +{ + 0:4 :: "" +} 4 + +restricted: +CCTK_INT rhs_timelevels "Number of active RHS timelevels" +{ + 0:4 :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_Minkowski_calc_every "ML_BSSN_M_Minkowski_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_convertFromADMBase_calc_every "ML_BSSN_M_convertFromADMBase_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_convertFromADMBaseGamma_calc_every "ML_BSSN_M_convertFromADMBaseGamma_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_setBetaDriverSpatial_calc_every "ML_BSSN_M_setBetaDriverSpatial_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_setBetaDriverConstant_calc_every "ML_BSSN_M_setBetaDriverConstant_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_RHS_calc_every "ML_BSSN_M_RHS_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_RHSStaticBoundary_calc_every "ML_BSSN_M_RHSStaticBoundary_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_RHSRadiativeBoundary_calc_every "ML_BSSN_M_RHSRadiativeBoundary_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_enforce_calc_every "ML_BSSN_M_enforce_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_boundary_calc_every "ML_BSSN_M_boundary_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_convertToADMBase_calc_every "ML_BSSN_M_convertToADMBase_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_convertToADMBaseDtLapseShift_calc_every "ML_BSSN_M_convertToADMBaseDtLapseShift_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_convertToADMBaseDtLapseShiftBoundary_calc_every "ML_BSSN_M_convertToADMBaseDtLapseShiftBoundary_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_convertToADMBaseFakeDtLapseShift_calc_every "ML_BSSN_M_convertToADMBaseFakeDtLapseShift_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_constraints_calc_every "ML_BSSN_M_constraints_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_constraints_boundary_calc_every "ML_BSSN_M_constraints_boundary_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_M_Minkowski_calc_offset "ML_BSSN_M_Minkowski_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_convertFromADMBase_calc_offset "ML_BSSN_M_convertFromADMBase_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_convertFromADMBaseGamma_calc_offset "ML_BSSN_M_convertFromADMBaseGamma_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_setBetaDriverSpatial_calc_offset "ML_BSSN_M_setBetaDriverSpatial_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_setBetaDriverConstant_calc_offset "ML_BSSN_M_setBetaDriverConstant_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_RHS_calc_offset "ML_BSSN_M_RHS_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_RHSStaticBoundary_calc_offset "ML_BSSN_M_RHSStaticBoundary_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_RHSRadiativeBoundary_calc_offset "ML_BSSN_M_RHSRadiativeBoundary_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_enforce_calc_offset "ML_BSSN_M_enforce_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_boundary_calc_offset "ML_BSSN_M_boundary_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_convertToADMBase_calc_offset "ML_BSSN_M_convertToADMBase_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_convertToADMBaseDtLapseShift_calc_offset "ML_BSSN_M_convertToADMBaseDtLapseShift_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_convertToADMBaseDtLapseShiftBoundary_calc_offset "ML_BSSN_M_convertToADMBaseDtLapseShiftBoundary_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_convertToADMBaseFakeDtLapseShift_calc_offset "ML_BSSN_M_convertToADMBaseFakeDtLapseShift_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_constraints_calc_offset "ML_BSSN_M_constraints_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_M_constraints_boundary_calc_offset "ML_BSSN_M_constraints_boundary_calc_offset" +{ + *:* :: "" +} 0 + +private: +KEYWORD At11_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD At12_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD At13_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD At22_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD At23_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD At33_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD A_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD B1_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD B2_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD B3_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD Xt1_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD Xt2_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD Xt3_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD alpha_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD phi_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD gt11_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD gt12_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD gt13_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD gt22_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD gt23_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD gt33_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD beta1_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD beta2_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD beta3_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD trK_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD ML_curv_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD ML_dtlapse_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD ML_dtshift_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD ML_Gamma_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD ML_lapse_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD ML_log_confac_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD ML_metric_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD ML_shift_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +KEYWORD ML_trace_curv_bound "Boundary condition to implement" +{ + "flat" :: "Flat boundary condition" + "none" :: "No boundary condition" + "static" :: "Boundaries held fixed" + "radiative" :: "Radiation boundary condition" + "scalar" :: "Dirichlet boundary condition" + "newrad" :: "Improved radiative boundary condition" + "skip" :: "skip boundary condition code" +} "skip" + +private: +CCTK_REAL At11_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL At12_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL At13_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL At22_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL At23_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL At33_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL A_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL B1_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL B2_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL B3_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL Xt1_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL Xt2_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL Xt3_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL alpha_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL phi_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL gt11_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL gt12_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL gt13_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL gt22_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL gt23_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL gt33_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL beta1_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL beta2_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL beta3_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL trK_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL ML_curv_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL ML_dtlapse_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL ML_dtshift_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL ML_Gamma_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL ML_lapse_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL ML_log_confac_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL ML_metric_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL ML_shift_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL ML_trace_curv_bound_speed "characteristic speed at boundary" +{ + "0:*" :: "outgoing characteristic speed > 0" +} 1. + +private: +CCTK_REAL At11_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL At12_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL At13_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL At22_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL At23_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL At33_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL A_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL B1_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL B2_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL B3_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL Xt1_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL Xt2_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL Xt3_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL alpha_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL phi_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL gt11_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL gt12_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL gt13_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL gt22_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL gt23_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL gt33_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL beta1_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL beta2_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL beta3_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL trK_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL ML_curv_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL ML_dtlapse_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL ML_dtshift_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL ML_Gamma_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL ML_lapse_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL ML_log_confac_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL ML_metric_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL ML_shift_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL ML_trace_curv_bound_limit "limit value for r -> infinity" +{ + "*:*" :: "value of limit value is unrestricted" +} 0. + +private: +CCTK_REAL At11_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL At12_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL At13_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL At22_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL At23_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL At33_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL A_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL B1_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL B2_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL B3_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL Xt1_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL Xt2_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL Xt3_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL alpha_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL phi_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL gt11_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL gt12_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL gt13_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL gt22_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL gt23_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL gt33_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL beta1_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL beta2_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL beta3_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL trK_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL ML_curv_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL ML_dtlapse_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL ML_dtshift_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL ML_Gamma_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL ML_lapse_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL ML_log_confac_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL ML_metric_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL ML_shift_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + +private: +CCTK_REAL ML_trace_curv_bound_scalar "Dirichlet boundary value" +{ + "*:*" :: "unrestricted" +} 0. + diff --git a/ML_BSSN_M/schedule.ccl b/ML_BSSN_M/schedule.ccl new file mode 100644 index 0000000..41cd90a --- /dev/null +++ b/ML_BSSN_M/schedule.ccl @@ -0,0 +1,498 @@ +# File produced by Kranc + + +STORAGE: ML_BetaDriver[1] + +STORAGE: ML_cons_detg[1] + +STORAGE: ML_cons_Gamma[1] + +STORAGE: ML_cons_traceA[1] + +STORAGE: ML_Ham[1] + +STORAGE: ML_mom[1] + +if (timelevels == 1) +{ + STORAGE: ML_curv[1] +} +if (timelevels == 2) +{ + STORAGE: ML_curv[2] +} +if (timelevels == 3) +{ + STORAGE: ML_curv[3] +} +if (timelevels == 4) +{ + STORAGE: ML_curv[4] +} + +if (timelevels == 1) +{ + STORAGE: ML_dtlapse[1] +} +if (timelevels == 2) +{ + STORAGE: ML_dtlapse[2] +} +if (timelevels == 3) +{ + STORAGE: ML_dtlapse[3] +} +if (timelevels == 4) +{ + STORAGE: ML_dtlapse[4] +} + +if (timelevels == 1) +{ + STORAGE: ML_dtshift[1] +} +if (timelevels == 2) +{ + STORAGE: ML_dtshift[2] +} +if (timelevels == 3) +{ + STORAGE: ML_dtshift[3] +} +if (timelevels == 4) +{ + STORAGE: ML_dtshift[4] +} + +if (timelevels == 1) +{ + STORAGE: ML_Gamma[1] +} +if (timelevels == 2) +{ + STORAGE: ML_Gamma[2] +} +if (timelevels == 3) +{ + STORAGE: ML_Gamma[3] +} +if (timelevels == 4) +{ + STORAGE: ML_Gamma[4] +} + +if (timelevels == 1) +{ + STORAGE: ML_lapse[1] +} +if (timelevels == 2) +{ + STORAGE: ML_lapse[2] +} +if (timelevels == 3) +{ + STORAGE: ML_lapse[3] +} +if (timelevels == 4) +{ + STORAGE: ML_lapse[4] +} + +if (timelevels == 1) +{ + STORAGE: ML_log_confac[1] +} +if (timelevels == 2) +{ + STORAGE: ML_log_confac[2] +} +if (timelevels == 3) +{ + STORAGE: ML_log_confac[3] +} +if (timelevels == 4) +{ + STORAGE: ML_log_confac[4] +} + +if (timelevels == 1) +{ + STORAGE: ML_metric[1] +} +if (timelevels == 2) +{ + STORAGE: ML_metric[2] +} +if (timelevels == 3) +{ + STORAGE: ML_metric[3] +} +if (timelevels == 4) +{ + STORAGE: ML_metric[4] +} + +if (timelevels == 1) +{ + STORAGE: ML_shift[1] +} +if (timelevels == 2) +{ + STORAGE: ML_shift[2] +} +if (timelevels == 3) +{ + STORAGE: ML_shift[3] +} +if (timelevels == 4) +{ + STORAGE: ML_shift[4] +} + +if (timelevels == 1) +{ + STORAGE: ML_trace_curv[1] +} +if (timelevels == 2) +{ + STORAGE: ML_trace_curv[2] +} +if (timelevels == 3) +{ + STORAGE: ML_trace_curv[3] +} +if (timelevels == 4) +{ + STORAGE: ML_trace_curv[4] +} + +if (rhs_timelevels == 1) +{ + STORAGE: ML_curvrhs[1] +} +if (rhs_timelevels == 2) +{ + STORAGE: ML_curvrhs[2] +} +if (rhs_timelevels == 3) +{ + STORAGE: ML_curvrhs[3] +} +if (rhs_timelevels == 4) +{ + STORAGE: ML_curvrhs[4] +} + +if (rhs_timelevels == 1) +{ + STORAGE: ML_dtlapserhs[1] +} +if (rhs_timelevels == 2) +{ + STORAGE: ML_dtlapserhs[2] +} +if (rhs_timelevels == 3) +{ + STORAGE: ML_dtlapserhs[3] +} +if (rhs_timelevels == 4) +{ + STORAGE: ML_dtlapserhs[4] +} + +if (rhs_timelevels == 1) +{ + STORAGE: ML_dtshiftrhs[1] +} +if (rhs_timelevels == 2) +{ + STORAGE: ML_dtshiftrhs[2] +} +if (rhs_timelevels == 3) +{ + STORAGE: ML_dtshiftrhs[3] +} +if (rhs_timelevels == 4) +{ + STORAGE: ML_dtshiftrhs[4] +} + +if (rhs_timelevels == 1) +{ + STORAGE: ML_Gammarhs[1] +} +if (rhs_timelevels == 2) +{ + STORAGE: ML_Gammarhs[2] +} +if (rhs_timelevels == 3) +{ + STORAGE: ML_Gammarhs[3] +} +if (rhs_timelevels == 4) +{ + STORAGE: ML_Gammarhs[4] +} + +if (rhs_timelevels == 1) +{ + STORAGE: ML_lapserhs[1] +} +if (rhs_timelevels == 2) +{ + STORAGE: ML_lapserhs[2] +} +if (rhs_timelevels == 3) +{ + STORAGE: ML_lapserhs[3] +} +if (rhs_timelevels == 4) +{ + STORAGE: ML_lapserhs[4] +} + +if (rhs_timelevels == 1) +{ + STORAGE: ML_log_confacrhs[1] +} +if (rhs_timelevels == 2) +{ + STORAGE: ML_log_confacrhs[2] +} +if (rhs_timelevels == 3) +{ + STORAGE: ML_log_confacrhs[3] +} +if (rhs_timelevels == 4) +{ + STORAGE: ML_log_confacrhs[4] +} + +if (rhs_timelevels == 1) +{ + STORAGE: ML_metricrhs[1] +} +if (rhs_timelevels == 2) +{ + STORAGE: ML_metricrhs[2] +} +if (rhs_timelevels == 3) +{ + STORAGE: ML_metricrhs[3] +} +if (rhs_timelevels == 4) +{ + STORAGE: ML_metricrhs[4] +} + +if (rhs_timelevels == 1) +{ + STORAGE: ML_shiftrhs[1] +} +if (rhs_timelevels == 2) +{ + STORAGE: ML_shiftrhs[2] +} +if (rhs_timelevels == 3) +{ + STORAGE: ML_shiftrhs[3] +} +if (rhs_timelevels == 4) +{ + STORAGE: ML_shiftrhs[4] +} + +if (rhs_timelevels == 1) +{ + STORAGE: ML_trace_curvrhs[1] +} +if (rhs_timelevels == 2) +{ + STORAGE: ML_trace_curvrhs[2] +} +if (rhs_timelevels == 3) +{ + STORAGE: ML_trace_curvrhs[3] +} +if (rhs_timelevels == 4) +{ + STORAGE: ML_trace_curvrhs[4] +} + +schedule ML_BSSN_M_Startup at STARTUP +{ + LANG: C + OPTIONS: meta +} "create banner" + +schedule ML_BSSN_M_RegisterVars in MoL_Register +{ + LANG: C + OPTIONS: meta +} "Register Variables for MoL" + +schedule ML_BSSN_M_RegisterSymmetries in SymmetryRegister +{ + LANG: C + OPTIONS: meta +} "register symmetries" + + +if (CCTK_EQUALS(my_initial_data, "Minkowski")) +{ + schedule ML_BSSN_M_Minkowski IN ADMBase_InitialData + { + LANG: C + } "ML_BSSN_M_Minkowski" +} + + +if (CCTK_EQUALS(my_initial_data, "ADMBase")) +{ + schedule ML_BSSN_M_convertFromADMBase AT initial AFTER ADMBase_PostInitial + { + LANG: C + } "ML_BSSN_M_convertFromADMBase" +} + + +if (CCTK_EQUALS(my_initial_data, "ADMBase")) +{ + schedule ML_BSSN_M_convertFromADMBaseGamma AT initial AFTER ML_BSSN_M_convertFromADMBase + { + LANG: C + SYNC: ML_dtlapse + SYNC: ML_dtshift + SYNC: ML_Gamma + } "ML_BSSN_M_convertFromADMBaseGamma" +} + + +if (CCTK_EQUALS(UseSpatialBetaDriver, "yes")) +{ + schedule ML_BSSN_M_setBetaDriverSpatial IN ML_BSSN_M_init_eta + { + LANG: C + } "ML_BSSN_M_setBetaDriverSpatial" +} + + +if (CCTK_EQUALS(UseSpatialBetaDriver, "no")) +{ + schedule ML_BSSN_M_setBetaDriverConstant IN ML_BSSN_M_init_eta + { + LANG: C + } "ML_BSSN_M_setBetaDriverConstant" +} + +schedule ML_BSSN_M_RHS IN ML_BSSN_M_evolCalcGroup +{ + LANG: C +} "ML_BSSN_M_RHS" + + +if (CCTK_EQUALS(my_rhs_boundary_condition, "static")) +{ + schedule ML_BSSN_M_RHSStaticBoundary IN MoL_CalcRHS + { + LANG: C + } "ML_BSSN_M_RHSStaticBoundary" +} + + +if (CCTK_EQUALS(my_rhs_boundary_condition, "radiative")) +{ + schedule ML_BSSN_M_RHSRadiativeBoundary IN MoL_CalcRHS + { + LANG: C + } "ML_BSSN_M_RHSRadiativeBoundary" +} + +schedule ML_BSSN_M_enforce IN MoL_PostStep BEFORE ML_BSSN_M_BoundConds +{ + LANG: C +} "ML_BSSN_M_enforce" + + +if (CCTK_EQUALS(my_boundary_condition, "Minkowski")) +{ + schedule ML_BSSN_M_boundary IN MoL_PostStep + { + LANG: C + } "ML_BSSN_M_boundary" +} + +schedule ML_BSSN_M_convertToADMBase IN ML_BSSN_M_convertToADMBaseGroup +{ + LANG: C +} "ML_BSSN_M_convertToADMBase" + + +if (CCTK_EQUALS(dt_lapse_shift_method, "correct")) +{ + schedule ML_BSSN_M_convertToADMBaseDtLapseShift IN ML_BSSN_M_convertToADMBaseGroup + { + LANG: C + SYNC: ADMBase::dtlapse + SYNC: ADMBase::dtshift + } "ML_BSSN_M_convertToADMBaseDtLapseShift" +} + + +if (CCTK_EQUALS(dt_lapse_shift_method, "correct")) +{ + schedule ML_BSSN_M_convertToADMBaseDtLapseShiftBoundary IN ML_BSSN_M_convertToADMBaseGroup + { + LANG: C + } "ML_BSSN_M_convertToADMBaseDtLapseShiftBoundary" +} + + +if (CCTK_EQUALS(dt_lapse_shift_method, "noLapseShiftAdvection")) +{ + schedule ML_BSSN_M_convertToADMBaseFakeDtLapseShift IN ML_BSSN_M_convertToADMBaseGroup + { + LANG: C + } "ML_BSSN_M_convertToADMBaseFakeDtLapseShift" +} + +schedule ML_BSSN_M_constraints IN ML_BSSN_M_constraintsCalcGroup +{ + LANG: C + SYNC: ML_cons_detg + SYNC: ML_cons_Gamma + SYNC: ML_cons_traceA + SYNC: ML_Ham + SYNC: ML_mom +} "ML_BSSN_M_constraints" + +schedule ML_BSSN_M_constraints_boundary IN ML_BSSN_M_constraintsCalcGroup AFTER ML_BSSN_M_constraints +{ + LANG: C +} "ML_BSSN_M_constraints_boundary" + +schedule ML_BSSN_M_SelectBoundConds in MoL_PostStep +{ + LANG: C + OPTIONS: level + SYNC: ML_curv + SYNC: ML_dtlapse + SYNC: ML_dtshift + SYNC: ML_Gamma + SYNC: ML_lapse + SYNC: ML_log_confac + SYNC: ML_metric + SYNC: ML_shift + SYNC: ML_trace_curv +} "select boundary conditions" + +schedule ML_BSSN_M_CheckBoundaries at BASEGRID +{ + LANG: C + OPTIONS: meta +} "check boundaries treatment" + +schedule group ApplyBCs as ML_BSSN_M_ApplyBCs in MoL_PostStep after ML_BSSN_M_SelectBoundConds +{ + # no language specified +} "Apply boundary conditions controlled by thorn Boundary" diff --git a/ML_BSSN_M/src/ML_BSSN_M_Minkowski.c b/ML_BSSN_M/src/ML_BSSN_M_Minkowski.c new file mode 100644 index 0000000..a4163e4 --- /dev/null +++ b/ML_BSSN_M/src/ML_BSSN_M_Minkowski.c @@ -0,0 +1,219 @@ +/* File produced by Kranc */ + +#define KRANC_C + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include "cctk.h" +#include "cctk_Arguments.h" +#include "cctk_Parameters.h" +#include "GenericFD.h" +#include "Differencing.h" +#include "loopcontrol.h" + +/* Define macros used in calculations */ +#define INITVALUE (42) +#define INV(x) ((1.0) / (x)) +#define SQR(x) ((x) * (x)) +#define CUB(x) ((x) * (x) * (x)) +#define QAD(x) ((x) * (x) * (x) * (x)) + +void ML_BSSN_M_Minkowski_Body(cGH const * const cctkGH, CCTK_INT const dir, CCTK_INT const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], CCTK_INT const min[3], CCTK_INT const max[3], CCTK_INT const n_subblock_gfs, CCTK_REAL * const subblock_gfs[]) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + + /* Declare finite differencing variables */ + CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE; + CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE; + CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE; + CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE; + + + /* Declare predefined quantities */ + CCTK_REAL p1o12dx = INITVALUE; + CCTK_REAL p1o12dy = INITVALUE; + CCTK_REAL p1o12dz = INITVALUE; + CCTK_REAL p1o144dxdy = INITVALUE; + CCTK_REAL p1o144dxdz = INITVALUE; + CCTK_REAL p1o144dydz = INITVALUE; + CCTK_REAL p1odx = INITVALUE; + CCTK_REAL p1ody = INITVALUE; + CCTK_REAL p1odz = INITVALUE; + CCTK_REAL pm1o12dx2 = INITVALUE; + CCTK_REAL pm1o12dy2 = INITVALUE; + CCTK_REAL pm1o12dz2 = INITVALUE; + + if (verbose > 1) + { + CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_BSSN_M_Minkowski_Body"); + } + + if (cctk_iteration % ML_BSSN_M_Minkowski_calc_every != ML_BSSN_M_Minkowski_calc_offset) + { + return; + } + + /* Include user-supplied include files */ + + /* Initialise finite differencing variables */ + dx = CCTK_DELTA_SPACE(0); + dy = CCTK_DELTA_SPACE(1); + dz = CCTK_DELTA_SPACE(2); + dxi = 1.0 / dx; + dyi = 1.0 / dy; + dzi = 1.0 / dz; + khalf = 0.5; + kthird = 1/3.0; + ktwothird = 2.0/3.0; + kfourthird = 4.0/3.0; + keightthird = 8.0/3.0; + hdxi = 0.5 * dxi; + hdyi = 0.5 * dyi; + hdzi = 0.5 * dzi; + + /* Initialize predefined quantities */ + p1o12dx = INV(dx)/12.; + p1o12dy = INV(dy)/12.; + p1o12dz = INV(dz)/12.; + p1o144dxdy = (INV(dx)*INV(dy))/144.; + p1o144dxdz = (INV(dx)*INV(dz))/144.; + p1o144dydz = (INV(dy)*INV(dz))/144.; + p1odx = INV(dx); + p1ody = INV(dy); + p1odz = INV(dz); + pm1o12dx2 = -pow(dx,-2)/12.; + pm1o12dy2 = -pow(dy,-2)/12.; + pm1o12dz2 = -pow(dz,-2)/12.; + + /* Loop over the grid points */ + #pragma omp parallel + LC_LOOP3 (ML_BSSN_M_Minkowski, + i,j,k, min[0],min[1],min[2], max[0],max[1],max[2], + cctk_lsh[0],cctk_lsh[1],cctk_lsh[2]) + { + int index = INITVALUE; + int subblock_index = INITVALUE; + index = CCTK_GFINDEX3D(cctkGH,i,j,k); + subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2])); + + /* Declare shorthands */ + + /* Declare local copies of grid functions */ + CCTK_REAL AL = INITVALUE; + CCTK_REAL alphaL = INITVALUE; + CCTK_REAL At11L = INITVALUE, At12L = INITVALUE, At13L = INITVALUE, At22L = INITVALUE, At23L = INITVALUE, At33L = INITVALUE; + CCTK_REAL B1L = INITVALUE, B2L = INITVALUE, B3L = INITVALUE; + CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE; + CCTK_REAL gt11L = INITVALUE, gt12L = INITVALUE, gt13L = INITVALUE, gt22L = INITVALUE, gt23L = INITVALUE, gt33L = INITVALUE; + CCTK_REAL phiL = INITVALUE; + CCTK_REAL trKL = INITVALUE; + CCTK_REAL Xt1L = INITVALUE, Xt2L = INITVALUE, Xt3L = INITVALUE; + /* Declare precomputed derivatives*/ + + /* Declare derivatives */ + + /* Assign local copies of grid functions */ + + /* Assign local copies of subblock grid functions */ + + /* Include user supplied include files */ + + /* Precompute derivatives (new style) */ + + /* Precompute derivatives (old style) */ + + /* Calculate temporaries and grid functions */ + phiL = IfThen(conformalMethod,1,0); + + gt11L = 1; + + gt12L = 0; + + gt13L = 0; + + gt22L = 1; + + gt23L = 0; + + gt33L = 1; + + trKL = 0; + + At11L = 0; + + At12L = 0; + + At13L = 0; + + At22L = 0; + + At23L = 0; + + At33L = 0; + + Xt1L = 0; + + Xt2L = 0; + + Xt3L = 0; + + alphaL = 1; + + AL = 0; + + beta1L = 0; + + beta2L = 0; + + beta3L = 0; + + B1L = 0; + + B2L = 0; + + B3L = 0; + + + /* Copy local copies back to grid functions */ + A[index] = AL; + alpha[index] = alphaL; + At11[index] = At11L; + At12[index] = At12L; + At13[index] = At13L; + At22[index] = At22L; + At23[index] = At23L; + At33[index] = At33L; + B1[index] = B1L; + B2[index] = B2L; + B3[index] = B3L; + beta1[index] = beta1L; + beta2[index] = beta2L; + beta3[index] = beta3L; + gt11[index] = gt11L; + gt12[index] = gt12L; + gt13[index] = gt13L; + gt22[index] = gt22L; + gt23[index] = gt23L; + gt33[index] = gt33L; + phi[index] = phiL; + trK[index] = trKL; + Xt1[index] = Xt1L; + Xt2[index] = Xt2L; + Xt3[index] = Xt3L; + + /* Copy local copies back to subblock grid functions */ + } + LC_ENDLOOP3 (ML_BSSN_M_Minkowski); +} + +void ML_BSSN_M_Minkowski(CCTK_ARGUMENTS) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + GenericFD_LoopOverEverything(cctkGH, &ML_BSSN_M_Minkowski_Body); +} diff --git a/ML_BSSN_M/src/ML_BSSN_M_convertFromADMBase.c b/ML_BSSN_M/src/ML_BSSN_M_convertFromADMBase.c new file mode 100644 index 0000000..5033a9c --- /dev/null +++ b/ML_BSSN_M/src/ML_BSSN_M_convertFromADMBase.c @@ -0,0 +1,274 @@ +/* File produced by Kranc */ + +#define KRANC_C + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include "cctk.h" +#include "cctk_Arguments.h" +#include "cctk_Parameters.h" +#include "GenericFD.h" +#include "Differencing.h" +#include "loopcontrol.h" + +/* Define macros used in calculations */ +#define INITVALUE (42) +#define INV(x) ((1.0) / (x)) +#define SQR(x) ((x) * (x)) +#define CUB(x) ((x) * (x) * (x)) +#define QAD(x) ((x) * (x) * (x) * (x)) + +void ML_BSSN_M_convertFromADMBase_Body(cGH const * const cctkGH, CCTK_INT const dir, CCTK_INT const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], CCTK_INT const min[3], CCTK_INT const max[3], CCTK_INT const n_subblock_gfs, CCTK_REAL * const subblock_gfs[]) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + + /* Declare finite differencing variables */ + CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE; + CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE; + CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE; + CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE; + + + /* Declare predefined quantities */ + CCTK_REAL p1o12dx = INITVALUE; + CCTK_REAL p1o12dy = INITVALUE; + CCTK_REAL p1o12dz = INITVALUE; + CCTK_REAL p1o144dxdy = INITVALUE; + CCTK_REAL p1o144dxdz = INITVALUE; + CCTK_REAL p1o144dydz = INITVALUE; + CCTK_REAL p1odx = INITVALUE; + CCTK_REAL p1ody = INITVALUE; + CCTK_REAL p1odz = INITVALUE; + CCTK_REAL pm1o12dx2 = INITVALUE; + CCTK_REAL pm1o12dy2 = INITVALUE; + CCTK_REAL pm1o12dz2 = INITVALUE; + + if (verbose > 1) + { + CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_BSSN_M_convertFromADMBase_Body"); + } + + if (cctk_iteration % ML_BSSN_M_convertFromADMBase_calc_every != ML_BSSN_M_convertFromADMBase_calc_offset) + { + return; + } + + /* Include user-supplied include files */ + + /* Initialise finite differencing variables */ + dx = CCTK_DELTA_SPACE(0); + dy = CCTK_DELTA_SPACE(1); + dz = CCTK_DELTA_SPACE(2); + dxi = 1.0 / dx; + dyi = 1.0 / dy; + dzi = 1.0 / dz; + khalf = 0.5; + kthird = 1/3.0; + ktwothird = 2.0/3.0; + kfourthird = 4.0/3.0; + keightthird = 8.0/3.0; + hdxi = 0.5 * dxi; + hdyi = 0.5 * dyi; + hdzi = 0.5 * dzi; + + /* Initialize predefined quantities */ + p1o12dx = INV(dx)/12.; + p1o12dy = INV(dy)/12.; + p1o12dz = INV(dz)/12.; + p1o144dxdy = (INV(dx)*INV(dy))/144.; + p1o144dxdz = (INV(dx)*INV(dz))/144.; + p1o144dydz = (INV(dy)*INV(dz))/144.; + p1odx = INV(dx); + p1ody = INV(dy); + p1odz = INV(dz); + pm1o12dx2 = -pow(dx,-2)/12.; + pm1o12dy2 = -pow(dy,-2)/12.; + pm1o12dz2 = -pow(dz,-2)/12.; + + /* Loop over the grid points */ + #pragma omp parallel + LC_LOOP3 (ML_BSSN_M_convertFromADMBase, + i,j,k, min[0],min[1],min[2], max[0],max[1],max[2], + cctk_lsh[0],cctk_lsh[1],cctk_lsh[2]) + { + int index = INITVALUE; + int subblock_index = INITVALUE; + index = CCTK_GFINDEX3D(cctkGH,i,j,k); + subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2])); + + /* Declare shorthands */ + CCTK_REAL detg = INITVALUE; + CCTK_REAL em4phi = INITVALUE; + CCTK_REAL g11 = INITVALUE, g12 = INITVALUE, g13 = INITVALUE, g22 = INITVALUE, g23 = INITVALUE, g33 = INITVALUE; + CCTK_REAL gu11 = INITVALUE, gu21 = INITVALUE, gu22 = INITVALUE, gu31 = INITVALUE, gu32 = INITVALUE, gu33 = INITVALUE; + CCTK_REAL K11 = INITVALUE, K12 = INITVALUE, K13 = INITVALUE, K22 = INITVALUE, K23 = INITVALUE, K33 = INITVALUE; + + /* Declare local copies of grid functions */ + CCTK_REAL alpL = INITVALUE; + CCTK_REAL alphaL = INITVALUE; + CCTK_REAL At11L = INITVALUE, At12L = INITVALUE, At13L = INITVALUE, At22L = INITVALUE, At23L = INITVALUE, At33L = INITVALUE; + CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE; + CCTK_REAL betaxL = INITVALUE; + CCTK_REAL betayL = INITVALUE; + CCTK_REAL betazL = INITVALUE; + CCTK_REAL gt11L = INITVALUE, gt12L = INITVALUE, gt13L = INITVALUE, gt22L = INITVALUE, gt23L = INITVALUE, gt33L = INITVALUE; + CCTK_REAL gxxL = INITVALUE; + CCTK_REAL gxyL = INITVALUE; + CCTK_REAL gxzL = INITVALUE; + CCTK_REAL gyyL = INITVALUE; + CCTK_REAL gyzL = INITVALUE; + CCTK_REAL gzzL = INITVALUE; + CCTK_REAL kxxL = INITVALUE; + CCTK_REAL kxyL = INITVALUE; + CCTK_REAL kxzL = INITVALUE; + CCTK_REAL kyyL = INITVALUE; + CCTK_REAL kyzL = INITVALUE; + CCTK_REAL kzzL = INITVALUE; + CCTK_REAL phiL = INITVALUE; + CCTK_REAL trKL = INITVALUE; + /* Declare precomputed derivatives*/ + + /* Declare derivatives */ + + /* Assign local copies of grid functions */ + alpL = alp[index]; + betaxL = betax[index]; + betayL = betay[index]; + betazL = betaz[index]; + gxxL = gxx[index]; + gxyL = gxy[index]; + gxzL = gxz[index]; + gyyL = gyy[index]; + gyzL = gyz[index]; + gzzL = gzz[index]; + kxxL = kxx[index]; + kxyL = kxy[index]; + kxzL = kxz[index]; + kyyL = kyy[index]; + kyzL = kyz[index]; + kzzL = kzz[index]; + phiL = phi[index]; + trKL = trK[index]; + + /* Assign local copies of subblock grid functions */ + + /* Include user supplied include files */ + + /* Precompute derivatives (new style) */ + + /* Precompute derivatives (old style) */ + + /* Calculate temporaries and grid functions */ + g11 = gxxL; + + g12 = gxyL; + + g13 = gxzL; + + g22 = gyyL; + + g23 = gyzL; + + g33 = gzzL; + + detg = 2*g12*g13*g23 + g33*(g11*g22 - SQR(g12)) - g22*SQR(g13) - g11*SQR(g23); + + gu11 = INV(detg)*(g22*g33 - SQR(g23)); + + gu21 = (g13*g23 - g12*g33)*INV(detg); + + gu31 = (-(g13*g22) + g12*g23)*INV(detg); + + gu22 = INV(detg)*(g11*g33 - SQR(g13)); + + gu32 = (g12*g13 - g11*g23)*INV(detg); + + gu33 = INV(detg)*(g11*g22 - SQR(g12)); + + phiL = IfThen(conformalMethod,pow(detg,-0.16666666666666666),Log(detg)/12.); + + em4phi = IfThen(conformalMethod,SQR(phiL),exp(-4*phiL)); + + gt11L = em4phi*g11; + + gt12L = em4phi*g12; + + gt13L = em4phi*g13; + + gt22L = em4phi*g22; + + gt23L = em4phi*g23; + + gt33L = em4phi*g33; + + K11 = kxxL; + + K12 = kxyL; + + K13 = kxzL; + + K22 = kyyL; + + K23 = kyzL; + + K33 = kzzL; + + trKL = gu11*K11 + gu22*K22 + 2*(gu21*K12 + gu31*K13 + gu32*K23) + gu33*K33; + + At11L = em4phi*(K11 - g11*kthird*trKL); + + At12L = em4phi*(K12 - g12*kthird*trKL); + + At13L = em4phi*(K13 - g13*kthird*trKL); + + At22L = em4phi*(K22 - g22*kthird*trKL); + + At23L = em4phi*(K23 - g23*kthird*trKL); + + At33L = em4phi*(K33 - g33*kthird*trKL); + + alphaL = alpL; + + beta1L = betaxL; + + beta2L = betayL; + + beta3L = betazL; + + + /* Copy local copies back to grid functions */ + alpha[index] = alphaL; + At11[index] = At11L; + At12[index] = At12L; + At13[index] = At13L; + At22[index] = At22L; + At23[index] = At23L; + At33[index] = At33L; + beta1[index] = beta1L; + beta2[index] = beta2L; + beta3[index] = beta3L; + gt11[index] = gt11L; + gt12[index] = gt12L; + gt13[index] = gt13L; + gt22[index] = gt22L; + gt23[index] = gt23L; + gt33[index] = gt33L; + phi[index] = phiL; + trK[index] = trKL; + + /* Copy local copies back to subblock grid functions */ + } + LC_ENDLOOP3 (ML_BSSN_M_convertFromADMBase); +} + +void ML_BSSN_M_convertFromADMBase(CCTK_ARGUMENTS) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + GenericFD_LoopOverEverything(cctkGH, &ML_BSSN_M_convertFromADMBase_Body); +} diff --git a/ML_BSSN_M/src/ML_BSSN_M_convertFromADMBaseGamma.c b/ML_BSSN_M/src/ML_BSSN_M_convertFromADMBaseGamma.c new file mode 100644 index 0000000..1b93415 --- /dev/null +++ b/ML_BSSN_M/src/ML_BSSN_M_convertFromADMBaseGamma.c @@ -0,0 +1,308 @@ +/* File produced by Kranc */ + +#define KRANC_C + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include "cctk.h" +#include "cctk_Arguments.h" +#include "cctk_Parameters.h" +#include "GenericFD.h" +#include "Differencing.h" +#include "loopcontrol.h" + +/* Define macros used in calculations */ +#define INITVALUE (42) +#define INV(x) ((1.0) / (x)) +#define SQR(x) ((x) * (x)) +#define CUB(x) ((x) * (x) * (x)) +#define QAD(x) ((x) * (x) * (x) * (x)) + +void ML_BSSN_M_convertFromADMBaseGamma_Body(cGH const * const cctkGH, CCTK_INT const dir, CCTK_INT const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], CCTK_INT const min[3], CCTK_INT const max[3], CCTK_INT const n_subblock_gfs, CCTK_REAL * const subblock_gfs[]) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + + /* Declare finite differencing variables */ + CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE; + CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE; + CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE; + CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE; + + + /* Declare predefined quantities */ + CCTK_REAL p1o12dx = INITVALUE; + CCTK_REAL p1o12dy = INITVALUE; + CCTK_REAL p1o12dz = INITVALUE; + CCTK_REAL p1o144dxdy = INITVALUE; + CCTK_REAL p1o144dxdz = INITVALUE; + CCTK_REAL p1o144dydz = INITVALUE; + CCTK_REAL p1odx = INITVALUE; + CCTK_REAL p1ody = INITVALUE; + CCTK_REAL p1odz = INITVALUE; + CCTK_REAL pm1o12dx2 = INITVALUE; + CCTK_REAL pm1o12dy2 = INITVALUE; + CCTK_REAL pm1o12dz2 = INITVALUE; + + if (verbose > 1) + { + CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_BSSN_M_convertFromADMBaseGamma_Body"); + } + + if (cctk_iteration % ML_BSSN_M_convertFromADMBaseGamma_calc_every != ML_BSSN_M_convertFromADMBaseGamma_calc_offset) + { + return; + } + + /* Include user-supplied include files */ + + /* Initialise finite differencing variables */ + dx = CCTK_DELTA_SPACE(0); + dy = CCTK_DELTA_SPACE(1); + dz = CCTK_DELTA_SPACE(2); + dxi = 1.0 / dx; + dyi = 1.0 / dy; + dzi = 1.0 / dz; + khalf = 0.5; + kthird = 1/3.0; + ktwothird = 2.0/3.0; + kfourthird = 4.0/3.0; + keightthird = 8.0/3.0; + hdxi = 0.5 * dxi; + hdyi = 0.5 * dyi; + hdzi = 0.5 * dzi; + + /* Initialize predefined quantities */ + p1o12dx = INV(dx)/12.; + p1o12dy = INV(dy)/12.; + p1o12dz = INV(dz)/12.; + p1o144dxdy = (INV(dx)*INV(dy))/144.; + p1o144dxdz = (INV(dx)*INV(dz))/144.; + p1o144dydz = (INV(dy)*INV(dz))/144.; + p1odx = INV(dx); + p1ody = INV(dy); + p1odz = INV(dz); + pm1o12dx2 = -pow(dx,-2)/12.; + pm1o12dy2 = -pow(dy,-2)/12.; + pm1o12dz2 = -pow(dz,-2)/12.; + + /* Loop over the grid points */ + #pragma omp parallel + LC_LOOP3 (ML_BSSN_M_convertFromADMBaseGamma, + i,j,k, min[0],min[1],min[2], max[0],max[1],max[2], + cctk_lsh[0],cctk_lsh[1],cctk_lsh[2]) + { + int index = INITVALUE; + int subblock_index = INITVALUE; + index = CCTK_GFINDEX3D(cctkGH,i,j,k); + subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2])); + + /* Declare shorthands */ + CCTK_REAL detgt = INITVALUE; + CCTK_REAL dir1 = INITVALUE, dir2 = INITVALUE, dir3 = INITVALUE; + CCTK_REAL Gt111 = INITVALUE, Gt112 = INITVALUE, Gt113 = INITVALUE, Gt122 = INITVALUE, Gt123 = INITVALUE, Gt133 = INITVALUE; + CCTK_REAL Gt211 = INITVALUE, Gt212 = INITVALUE, Gt213 = INITVALUE, Gt222 = INITVALUE, Gt223 = INITVALUE, Gt233 = INITVALUE; + CCTK_REAL Gt311 = INITVALUE, Gt312 = INITVALUE, Gt313 = INITVALUE, Gt322 = INITVALUE, Gt323 = INITVALUE, Gt333 = INITVALUE; + CCTK_REAL gtu11 = INITVALUE, gtu21 = INITVALUE, gtu22 = INITVALUE, gtu31 = INITVALUE, gtu32 = INITVALUE, gtu33 = INITVALUE; + + /* Declare local copies of grid functions */ + CCTK_REAL AL = INITVALUE; + CCTK_REAL alphaL = INITVALUE; + CCTK_REAL B1L = INITVALUE, B2L = INITVALUE, B3L = INITVALUE; + CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE; + CCTK_REAL dtalpL = INITVALUE; + CCTK_REAL dtbetaxL = INITVALUE; + CCTK_REAL dtbetayL = INITVALUE; + CCTK_REAL dtbetazL = INITVALUE; + CCTK_REAL gt11L = INITVALUE, gt12L = INITVALUE, gt13L = INITVALUE, gt22L = INITVALUE, gt23L = INITVALUE, gt33L = INITVALUE; + CCTK_REAL Xt1L = INITVALUE, Xt2L = INITVALUE, Xt3L = INITVALUE; + /* Declare precomputed derivatives*/ + + /* Declare derivatives */ + CCTK_REAL PDstandardNth1gt11 = INITVALUE; + CCTK_REAL PDstandardNth2gt11 = INITVALUE; + CCTK_REAL PDstandardNth3gt11 = INITVALUE; + CCTK_REAL PDstandardNth1gt12 = INITVALUE; + CCTK_REAL PDstandardNth2gt12 = INITVALUE; + CCTK_REAL PDstandardNth3gt12 = INITVALUE; + CCTK_REAL PDstandardNth1gt13 = INITVALUE; + CCTK_REAL PDstandardNth2gt13 = INITVALUE; + CCTK_REAL PDstandardNth3gt13 = INITVALUE; + CCTK_REAL PDstandardNth1gt22 = INITVALUE; + CCTK_REAL PDstandardNth2gt22 = INITVALUE; + CCTK_REAL PDstandardNth3gt22 = INITVALUE; + CCTK_REAL PDstandardNth1gt23 = INITVALUE; + CCTK_REAL PDstandardNth2gt23 = INITVALUE; + CCTK_REAL PDstandardNth3gt23 = INITVALUE; + CCTK_REAL PDstandardNth1gt33 = INITVALUE; + CCTK_REAL PDstandardNth2gt33 = INITVALUE; + CCTK_REAL PDstandardNth3gt33 = INITVALUE; + + /* Assign local copies of grid functions */ + alphaL = alpha[index]; + beta1L = beta1[index]; + beta2L = beta2[index]; + beta3L = beta3[index]; + dtalpL = dtalp[index]; + dtbetaxL = dtbetax[index]; + dtbetayL = dtbetay[index]; + dtbetazL = dtbetaz[index]; + gt11L = gt11[index]; + gt12L = gt12[index]; + gt13L = gt13[index]; + gt22L = gt22[index]; + gt23L = gt23[index]; + gt33L = gt33[index]; + + /* Assign local copies of subblock grid functions */ + + /* Include user supplied include files */ + + /* Precompute derivatives (new style) */ + PDstandardNth1gt11 = PDstandardNth1(gt11, i, j, k); + PDstandardNth2gt11 = PDstandardNth2(gt11, i, j, k); + PDstandardNth3gt11 = PDstandardNth3(gt11, i, j, k); + PDstandardNth1gt12 = PDstandardNth1(gt12, i, j, k); + PDstandardNth2gt12 = PDstandardNth2(gt12, i, j, k); + PDstandardNth3gt12 = PDstandardNth3(gt12, i, j, k); + PDstandardNth1gt13 = PDstandardNth1(gt13, i, j, k); + PDstandardNth2gt13 = PDstandardNth2(gt13, i, j, k); + PDstandardNth3gt13 = PDstandardNth3(gt13, i, j, k); + PDstandardNth1gt22 = PDstandardNth1(gt22, i, j, k); + PDstandardNth2gt22 = PDstandardNth2(gt22, i, j, k); + PDstandardNth3gt22 = PDstandardNth3(gt22, i, j, k); + PDstandardNth1gt23 = PDstandardNth1(gt23, i, j, k); + PDstandardNth2gt23 = PDstandardNth2(gt23, i, j, k); + PDstandardNth3gt23 = PDstandardNth3(gt23, i, j, k); + PDstandardNth1gt33 = PDstandardNth1(gt33, i, j, k); + PDstandardNth2gt33 = PDstandardNth2(gt33, i, j, k); + PDstandardNth3gt33 = PDstandardNth3(gt33, i, j, k); + + /* Precompute derivatives (old style) */ + + /* Calculate temporaries and grid functions */ + dir1 = Sign(beta1L); + + dir2 = Sign(beta2L); + + dir3 = Sign(beta3L); + + detgt = 1; + + gtu11 = INV(detgt)*(gt22L*gt33L - SQR(gt23L)); + + gtu21 = (gt13L*gt23L - gt12L*gt33L)*INV(detgt); + + gtu31 = (-(gt13L*gt22L) + gt12L*gt23L)*INV(detgt); + + gtu22 = INV(detgt)*(gt11L*gt33L - SQR(gt13L)); + + gtu32 = (gt12L*gt13L - gt11L*gt23L)*INV(detgt); + + gtu33 = INV(detgt)*(gt11L*gt22L - SQR(gt12L)); + + Gt111 = khalf*(gtu11*PDstandardNth1gt11 + 2*(gtu21*PDstandardNth1gt12 + gtu31*PDstandardNth1gt13) - + gtu21*PDstandardNth2gt11 - gtu31*PDstandardNth3gt11); + + Gt211 = khalf*(gtu21*PDstandardNth1gt11 + 2*(gtu22*PDstandardNth1gt12 + gtu32*PDstandardNth1gt13) - + gtu22*PDstandardNth2gt11 - gtu32*PDstandardNth3gt11); + + Gt311 = khalf*(gtu31*PDstandardNth1gt11 + 2*(gtu32*PDstandardNth1gt12 + gtu33*PDstandardNth1gt13) - + gtu32*PDstandardNth2gt11 - gtu33*PDstandardNth3gt11); + + Gt112 = khalf*(gtu21*PDstandardNth1gt22 + gtu11*PDstandardNth2gt11 + + gtu31*(PDstandardNth1gt23 + PDstandardNth2gt13 - PDstandardNth3gt12)); + + Gt212 = khalf*(gtu22*PDstandardNth1gt22 + gtu21*PDstandardNth2gt11 + + gtu32*(PDstandardNth1gt23 + PDstandardNth2gt13 - PDstandardNth3gt12)); + + Gt312 = khalf*(gtu32*PDstandardNth1gt22 + gtu31*PDstandardNth2gt11 + + gtu33*(PDstandardNth1gt23 + PDstandardNth2gt13 - PDstandardNth3gt12)); + + Gt113 = khalf*(gtu31*PDstandardNth1gt33 + gtu11*PDstandardNth3gt11 + + gtu21*(PDstandardNth1gt23 - PDstandardNth2gt13 + PDstandardNth3gt12)); + + Gt213 = khalf*(gtu32*PDstandardNth1gt33 + gtu21*PDstandardNth3gt11 + + gtu22*(PDstandardNth1gt23 - PDstandardNth2gt13 + PDstandardNth3gt12)); + + Gt313 = khalf*(gtu33*PDstandardNth1gt33 + gtu31*PDstandardNth3gt11 + + gtu32*(PDstandardNth1gt23 - PDstandardNth2gt13 + PDstandardNth3gt12)); + + Gt122 = khalf*(gtu11*(-PDstandardNth1gt22 + 2*PDstandardNth2gt12) + gtu21*PDstandardNth2gt22 + + gtu31*(2*PDstandardNth2gt23 - PDstandardNth3gt22)); + + Gt222 = khalf*(gtu21*(-PDstandardNth1gt22 + 2*PDstandardNth2gt12) + gtu22*PDstandardNth2gt22 + + gtu32*(2*PDstandardNth2gt23 - PDstandardNth3gt22)); + + Gt322 = khalf*(gtu31*(-PDstandardNth1gt22 + 2*PDstandardNth2gt12) + gtu32*PDstandardNth2gt22 + + gtu33*(2*PDstandardNth2gt23 - PDstandardNth3gt22)); + + Gt123 = khalf*(gtu31*PDstandardNth2gt33 + gtu11*(-PDstandardNth1gt23 + PDstandardNth2gt13 + PDstandardNth3gt12) + + gtu21*PDstandardNth3gt22); + + Gt223 = khalf*(gtu32*PDstandardNth2gt33 + gtu21*(-PDstandardNth1gt23 + PDstandardNth2gt13 + PDstandardNth3gt12) + + gtu22*PDstandardNth3gt22); + + Gt323 = khalf*(gtu33*PDstandardNth2gt33 + gtu31*(-PDstandardNth1gt23 + PDstandardNth2gt13 + PDstandardNth3gt12) + + gtu32*PDstandardNth3gt22); + + Gt133 = khalf*(-(gtu11*PDstandardNth1gt33) - gtu21*PDstandardNth2gt33 + 2*gtu11*PDstandardNth3gt13 + + 2*gtu21*PDstandardNth3gt23 + gtu31*PDstandardNth3gt33); + + Gt233 = khalf*(-(gtu21*PDstandardNth1gt33) - gtu22*PDstandardNth2gt33 + 2*gtu21*PDstandardNth3gt13 + + 2*gtu22*PDstandardNth3gt23 + gtu32*PDstandardNth3gt33); + + Gt333 = khalf*(-(gtu31*PDstandardNth1gt33) - gtu32*PDstandardNth2gt33 + 2*gtu31*PDstandardNth3gt13 + + 2*gtu32*PDstandardNth3gt23 + gtu33*PDstandardNth3gt33); + + Xt1L = Gt111*gtu11 + Gt122*gtu22 + 2*(Gt112*gtu21 + Gt113*gtu31 + Gt123*gtu32) + Gt133*gtu33; + + Xt2L = Gt211*gtu11 + Gt222*gtu22 + 2*(Gt212*gtu21 + Gt213*gtu31 + Gt223*gtu32) + Gt233*gtu33; + + Xt3L = Gt311*gtu11 + Gt322*gtu22 + 2*(Gt312*gtu21 + Gt313*gtu31 + Gt323*gtu32) + Gt333*gtu33; + + AL = -(dtalpL*(-1 + LapseAdvectionCoeff)*INV(harmonicF)*pow(alphaL,-harmonicN)); + + B1L = 6*IfThen(ShiftGammaCoeff,dtbetaxL*INV(ShiftGammaCoeff),0) + + IfThen(ShiftGammaCoeff,(dtbetaxL - PDupwindNth1(beta1, i, j, k)*beta1L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff), + 0) + IfThen(ShiftGammaCoeff,(dtbetaxL - PDupwindNth2(beta1, i, j, k)*beta2L*ShiftAdvectionCoeff)* + INV(ShiftGammaCoeff),0) + IfThen(ShiftGammaCoeff, + (dtbetaxL - PDupwindNth3(beta1, i, j, k)*beta3L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0); + + B2L = 6*IfThen(ShiftGammaCoeff,dtbetayL*INV(ShiftGammaCoeff),0) + + IfThen(ShiftGammaCoeff,(dtbetayL - PDupwindNth1(beta2, i, j, k)*beta1L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff), + 0) + IfThen(ShiftGammaCoeff,(dtbetayL - PDupwindNth2(beta2, i, j, k)*beta2L*ShiftAdvectionCoeff)* + INV(ShiftGammaCoeff),0) + IfThen(ShiftGammaCoeff, + (dtbetayL - PDupwindNth3(beta2, i, j, k)*beta3L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0); + + B3L = 6*IfThen(ShiftGammaCoeff,dtbetazL*INV(ShiftGammaCoeff),0) + + IfThen(ShiftGammaCoeff,(dtbetazL - PDupwindNth1(beta3, i, j, k)*beta1L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff), + 0) + IfThen(ShiftGammaCoeff,(dtbetazL - PDupwindNth2(beta3, i, j, k)*beta2L*ShiftAdvectionCoeff)* + INV(ShiftGammaCoeff),0) + IfThen(ShiftGammaCoeff, + (dtbetazL - PDupwindNth3(beta3, i, j, k)*beta3L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0); + + + /* Copy local copies back to grid functions */ + A[index] = AL; + B1[index] = B1L; + B2[index] = B2L; + B3[index] = B3L; + Xt1[index] = Xt1L; + Xt2[index] = Xt2L; + Xt3[index] = Xt3L; + + /* Copy local copies back to subblock grid functions */ + } + LC_ENDLOOP3 (ML_BSSN_M_convertFromADMBaseGamma); +} + +void ML_BSSN_M_convertFromADMBaseGamma(CCTK_ARGUMENTS) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + GenericFD_LoopOverInterior(cctkGH, &ML_BSSN_M_convertFromADMBaseGamma_Body); +} diff --git a/ML_BSSN_M/src/ML_BSSN_M_setBetaDriverSpatial.c b/ML_BSSN_M/src/ML_BSSN_M_setBetaDriverSpatial.c new file mode 100644 index 0000000..afc6635 --- /dev/null +++ b/ML_BSSN_M/src/ML_BSSN_M_setBetaDriverSpatial.c @@ -0,0 +1,141 @@ +/* File produced by Kranc */ + +#define KRANC_C + +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include "cctk.h" +#include "cctk_Arguments.h" +#include "cctk_Parameters.h" +#include "GenericFD.h" +#include "Differencing.h" +#include "loopcontrol.h" + +/* Define macros used in calculations */ +#define INITVALUE (42) +#define INV(x) ((1.0) / (x)) +#define SQR(x) ((x) * (x)) +#define CUB(x) ((x) * (x) * (x)) +#define QAD(x) ((x) * (x) * (x) * (x)) + +void ML_BSSN_M_setBetaDriverSpatial_Body(cGH const * const cctkGH, CCTK_INT const dir, CCTK_INT const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], CCTK_INT const min[3], CCTK_INT const max[3], CCTK_INT const n_subblock_gfs, CCTK_REAL * const subblock_gfs[]) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + + /* Declare finite differencing variables */ + CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE; + CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE; + CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE; + CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE; + + + /* Declare predefined quantities */ + CCTK_REAL p1o12dx = INITVALUE; + CCTK_REAL p1o12dy = INITVALUE; + CCTK_REAL p1o12dz = INITVALUE; + CCTK_REAL p1o144dxdy = INITVALUE; + CCTK_REAL p1o144dxdz = INITVALUE; + CCTK_REAL p1o144dydz = INITVALUE; + CCTK_REAL p1odx = INITVALUE; + CCTK_REAL p1ody = INITVALUE; + CCTK_REAL p1odz = INITVALUE; + CCTK_REAL pm1o12dx2 = INITVALUE; + CCTK_REAL pm1o12dy2 = INITVALUE; + CCTK_REAL pm1o12dz2 = INITVALUE; + + if (verbose > 1) + { + CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_BSSN_M_setBetaDriverSpatial_Body"); + } + + if (cctk_iteration % ML_BSSN_M_setBetaDriverSpatial_calc_every != ML_BSSN_M_setBetaDriverSpatial_calc_offset) + { + return; + } + + /* Include user-supplied include files */ + + /* Initialise finite differencing variables */ + dx = CCTK_DELTA_SPACE(0); + dy = CCTK_DELTA_SPACE(1); + dz = CCTK_DELTA_SPACE(2); + dxi = 1.0 / dx; + dyi = 1.0 / dy; + dzi = 1.0 / dz; + khalf = 0.5; + kthird = 1/3.0; + ktwothird = 2.0/3.0; + kfourthird = 4.0/3.0; + keightthird = 8.0/3.0; + hdxi = 0.5 * dxi; + hdyi = 0.5 * dyi; + hdzi = 0.5 * dzi; + + /* Initialize predefined quantities */ + p1o12dx = INV(dx)/12.; + p1o12dy = INV(dy)/12.; + p1o12dz = INV(dz)/12.; + p1o144dxdy = (INV(dx)*INV(dy))/144.; + p1o144dxdz = (INV(dx)*INV(dz))/144.; + p1o144dydz = (INV(dy)*INV(dz))/144.; + p1odx = INV(dx); + p1ody = INV(dy); + p1odz = INV(dz); + pm1o12dx2 = -pow(dx,-2)/12.; + pm1o12dy2 = -pow(dy,-2)/12.; + pm1o12dz2 = -pow(dz,-2)/12.; + + /* Loop over the grid points */ + #pragma omp parallel + LC_LOOP3 (ML_BSSN_M_setBetaDriverSpatial, + i,j,k, min[0],min[1],min[2], max[0],max[1],max[2], + cctk_lsh[0],cctk_lsh[1],cctk_lsh[2]) + { + int index = INITVALUE; + int subblock_index = INITVALUE; + index = CCTK_GFINDEX3D(cctkGH,i,j,k); + subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2])); + + /* Declare shorthands */ + + /* Declare local copies of grid functions */ + CCTK_REAL etaL = INITVALUE; + CCTK_REAL rL = INITVALUE; + /* Declare precomputed derivatives*/ + + /* Declare derivatives */ + + /* Assign local copies of grid functions */ + rL = r[index]; + + /* Assign local copies of subblock grid functions */ + + /* Include user supplied include files */ + + /* Precompute derivatives (new style) */ + + /* Precompute derivatives (old style) */ + + /* Calculate temporaries and grid functions */ + etaL = BetaDriver*IfThen(rL > SpatialBetaDriverRadius,SpatialBetaDriverRadius*INV(rL),1); + + + /* Copy local copies back to grid functions */ + eta[index] = etaL; + + /* Copy local copies back to subblock grid functions */ + } + LC_ENDLOOP3 (ML_BSSN_M_setBetaDriverSpatial); +} + +void ML_BSSN_M_setBetaDriverSpatial(CCTK_ARGUMENTS) +{ + DECLARE_CCTK_ARGUMENTS; + DECLARE_CCTK_PARAMETERS; + + GenericFD_LoopOverEverything(cctkGH, &ML_BSSN_M_setBetaDriverSpatial_Body); +} diff --git a/ML_BSSN_M/src/make.code.defn b/ML_BSSN_M/src/make.code.defn new file mode 100644 index 0000000..2c89d2a --- /dev/null +++ b/ML_BSSN_M/src/make.code.defn @@ -0,0 +1,3 @@ +# File produced by Kranc + +SRCS = Startup.c RegisterMoL.c RegisterSymmetries.c ML_BSSN_M_Minkowski.c ML_BSSN_M_convertFromADMBase.c ML_BSSN_M_convertFromADMBaseGamma.c ML_BSSN_M_setBetaDriverSpatial.c ML_BSSN_M_setBetaDriverConstant.c ML_BSSN_M_RHS.c ML_BSSN_M_RHSStaticBoundary.c ML_BSSN_M_RHSRadiativeBoundary.c ML_BSSN_M_enforce.c ML_BSSN_M_boundary.c ML_BSSN_M_convertToADMBase.c ML_BSSN_M_convertToADMBaseDtLapseShift.c ML_BSSN_M_convertToADMBaseDtLapseShiftBoundary.c ML_BSSN_M_convertToADMBaseFakeDtLapseShift.c ML_BSSN_M_constraints.c ML_BSSN_M_constraints_boundary.c Boundaries.c diff --git a/ML_BSSN_MP/param.ccl b/ML_BSSN_MP/param.ccl index f4beb67..f991fd1 100644 --- a/ML_BSSN_MP/param.ccl +++ b/ML_BSSN_MP/param.ccl @@ -213,7 +213,13 @@ CCTK_INT ML_BSSN_MP_convertFromADMBaseGamma_calc_every "ML_BSSN_MP_convertFromAD } 1 restricted: -CCTK_INT ML_BSSN_MP_setBetaDriver_calc_every "ML_BSSN_MP_setBetaDriver_calc_every" +CCTK_INT ML_BSSN_MP_setBetaDriverConstant_calc_every "ML_BSSN_MP_setBetaDriverConstant_calc_every" +{ + *:* :: "" +} 1 + +restricted: +CCTK_INT ML_BSSN_MP_setBetaDriverSpatial_calc_every "ML_BSSN_MP_setBetaDriverSpatial_calc_every" { *:* :: "" } 1 @@ -321,7 +327,13 @@ CCTK_INT ML_BSSN_MP_convertFromADMBaseGamma_calc_offset "ML_BSSN_MP_convertFromA } 0 restricted: -CCTK_INT ML_BSSN_MP_setBetaDriver_calc_offset "ML_BSSN_MP_setBetaDriver_calc_offset" +CCTK_INT ML_BSSN_MP_setBetaDriverConstant_calc_offset "ML_BSSN_MP_setBetaDriverConstant_calc_offset" +{ + *:* :: "" +} 0 + +restricted: +CCTK_INT ML_BSSN_MP_setBetaDriverSpatial_calc_offset "ML_BSSN_MP_setBetaDriverSpatial_calc_offset" { *:* :: "" } 0 diff --git a/ML_BSSN_MP/schedule.ccl b/ML_BSSN_MP/schedule.ccl index 3b38d75..e5cf2ae 100644 --- a/ML_BSSN_MP/schedule.ccl +++ b/ML_BSSN_MP/schedule.ccl @@ -368,12 +368,21 @@ if (CCTK_EQUALS(my_initial_data, "ADMBase")) } +if (CCTK_EQUALS(UseSpatialBetaDriver, "no")) +{ + schedule ML_BSSN_MP_setBetaDriverConstant IN ML_BSSN_MP_InitEta + { + LANG: C + } "ML_BSSN_MP_setBetaDriverConstant" +} + + if (CCTK_EQUALS(UseSpatialBetaDriver, "yes")) { - schedule ML_BSSN_MP_setBetaDriver AT initial AFTER ADMBase_PostInitial AFTER ML_BSSN_MP_convertFromADMBase + schedule ML_BSSN_MP_setBetaDriverSpatial IN ML_BSSN_MP_InitEta { LANG: C - } "ML_BSSN_MP_setBetaDriver" + } "ML_BSSN_MP_setBetaDriverSpatial" } schedule ML_BSSN_MP_RHS IN NoSuchGroup diff --git a/ML_BSSN_MP/src/ML_BSSN_MP_Minkowski.c b/ML_BSSN_MP/src/ML_BSSN_MP_Minkowski.c index 74efd8b..12429c4 100644 --- a/ML_BSSN_MP/src/ML_BSSN_MP_Minkowski.c +++ b/ML_BSSN_MP/src/ML_BSSN_MP_Minkowski.c @@ -106,7 +106,6 @@ void ML_BSSN_MP_Minkowski_Body(cGH const * restrict const cctkGH, int const dir, // CCTK_REAL At11L = INITVALUE, At12L = INITVALUE, At13L = INITVALUE, At22L = INITVALUE, At23L = INITVALUE, At33L = INITVALUE; // CCTK_REAL B1L = INITVALUE, B2L = INITVALUE, B3L = INITVALUE; // CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE; - // CCTK_REAL etaL = INITVALUE; // CCTK_REAL gt11L = INITVALUE, gt12L = INITVALUE, gt13L = INITVALUE, gt22L = INITVALUE, gt23L = INITVALUE, gt33L = INITVALUE; // CCTK_REAL phiL = INITVALUE; // CCTK_REAL trKL = INITVALUE; @@ -176,8 +175,6 @@ void ML_BSSN_MP_Minkowski_Body(cGH const * restrict const cctkGH, int const dir, CCTK_REAL const B3L = 0; - CCTK_REAL const etaL = BetaDriver; - /* Copy local copies back to grid functions */ A[index] = AL; @@ -194,7 +191,6 @@ void ML_BSSN_MP_Minkowski_Body(cGH const * restrict const cctkGH, int const dir, beta1[index] = beta1L; beta2[index] = beta2L; beta3[index] = beta3L; - eta[index] = etaL; gt11[index] = gt11L; gt12[index] = gt12L; gt13[index] = gt13L; diff --git a/ML_BSSN_MP/src/ML_BSSN_MP_RHS.c b/ML_BSSN_MP/src/ML_BSSN_MP_RHS.c index 83e567c..a8d8cee 100644 --- a/ML_BSSN_MP/src/ML_BSSN_MP_RHS.c +++ b/ML_BSSN_MP/src/ML_BSSN_MP_RHS.c @@ -132,8 +132,6 @@ void ML_BSSN_MP_RHS_Body(cGH const * restrict const cctkGH, int const dir, int c // CCTK_REAL Rphi11 = INITVALUE, Rphi12 = INITVALUE, Rphi13 = INITVALUE, Rphi22 = INITVALUE, Rphi23 = INITVALUE, Rphi33 = 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 trAts = INITVALUE; // CCTK_REAL trS = INITVALUE; // CCTK_REAL Xtn1 = INITVALUE, Xtn2 = INITVALUE, Xtn3 = INITVALUE; @@ -1164,37 +1162,17 @@ void ML_BSSN_MP_RHS_Body(cGH const * restrict const cctkGH, int const dir, int c CCTK_REAL const R33 = Rphi33 + Rt33; - 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 S3 = (-eTtzL + beta1L*eTxzL + beta2L*eTyzL + beta3L*eTzzL)*INV(alphaL); - 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 trS = gu11*T11 + gu22*T22 + 2*(gu21*T12 + gu31*T13 + gu32*T23) + gu33*T33; + CCTK_REAL const trS = eTxxL*gu11 + eTyyL*gu22 + 2*(eTxyL*gu21 + eTxzL*gu31 + eTyzL*gu32) + eTzzL*gu33; CCTK_REAL const phirhsL = PDupwindNth1(phi, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(phi, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1581,7 +1559,7 @@ void ML_BSSN_MP_RHS_Body(cGH const * restrict const cctkGH, int const dir, int c J21L*(At11L*PDstandardNth2beta1 + At12L*PDstandardNth2beta2 + At13L*PDstandardNth2beta3) + J31L*(At11L*PDstandardNth3beta1 + At12L*PDstandardNth3beta2 + At13L*PDstandardNth3beta3)) + alphaL*(-2*(At11L*Atm11 + At12L*Atm21 + At13L*Atm31) + At11L*trKL) + - em4phi*(Ats11 - g11*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(T11 - g11*kthird*trS)); + em4phi*(Ats11 - g11*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(eTxxL - g11*kthird*trS)); CCTK_REAL const At12rhsL = PDupwindNth1(At12, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At12, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1597,7 +1575,7 @@ void ML_BSSN_MP_RHS_Body(cGH const * restrict const cctkGH, int const dir, int c (At23L*J31L + At13L*J32L - 0.6666666666666666666666666666666666666667*At12L*J33L)*PDstandardNth3beta3 + alphaL*(-2.*(At11L*Atm12 + At12L*Atm22 + At13L*Atm32) + At12L*trKL) + em4phi*(Ats12 - 0.3333333333333333333333333333333333333333*g12*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T12 + 8.377580409572781969233715688745341024526*g12*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxyL + 8.377580409572781969233715688745341024526*g12*trS)); CCTK_REAL const At13rhsL = PDupwindNth1(At13, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At13, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1613,7 +1591,7 @@ void ML_BSSN_MP_RHS_Body(cGH const * restrict const cctkGH, int const dir, int c (At33L*J31L + 0.3333333333333333333333333333333333333333*At13L*J33L)*PDstandardNth3beta3 + alphaL*(-2.*(At11L*Atm13 + At12L*Atm23 + At13L*Atm33) + At13L*trKL) + em4phi*(Ats13 - 0.3333333333333333333333333333333333333333*g13*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T13 + 8.377580409572781969233715688745341024526*g13*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxzL + 8.377580409572781969233715688745341024526*g13*trS)); CCTK_REAL const At22rhsL = PDupwindNth1(At22, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At22, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1625,7 +1603,7 @@ void ML_BSSN_MP_RHS_Body(cGH const * restrict const cctkGH, int const dir, int c J22L*(At12L*PDstandardNth2beta1 + At22L*PDstandardNth2beta2 + At23L*PDstandardNth2beta3) + J32L*(At12L*PDstandardNth3beta1 + At22L*PDstandardNth3beta2 + At23L*PDstandardNth3beta3)) + alphaL*(-2*(At12L*Atm12 + At22L*Atm22 + At23L*Atm32) + At22L*trKL) + - em4phi*(Ats22 - g22*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(T22 - g22*kthird*trS)); + em4phi*(Ats22 - g22*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(eTyyL - g22*kthird*trS)); CCTK_REAL const At23rhsL = PDupwindNth1(At23, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At23, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1641,7 +1619,7 @@ void ML_BSSN_MP_RHS_Body(cGH const * restrict const cctkGH, int const dir, int c (At33L*J32L + 0.3333333333333333333333333333333333333333*At23L*J33L)*PDstandardNth3beta3 + alphaL*(-2.*(At12L*Atm13 + At22L*Atm23 + At23L*Atm33) + At23L*trKL) + em4phi*(Ats23 - 0.3333333333333333333333333333333333333333*g23*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T23 + 8.377580409572781969233715688745341024526*g23*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTyzL + 8.377580409572781969233715688745341024526*g23*trS)); CCTK_REAL const At33rhsL = PDupwindNth1(At33, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At33, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1653,7 +1631,7 @@ void ML_BSSN_MP_RHS_Body(cGH const * restrict const cctkGH, int const dir, int c J23L*(At13L*PDstandardNth2beta1 + At23L*PDstandardNth2beta2 + At33L*PDstandardNth2beta3) + J33L*(At13L*PDstandardNth3beta1 + At23L*PDstandardNth3beta2 + At33L*PDstandardNth3beta3)) + alphaL*(-2*(At13L*Atm13 + At23L*Atm23 + At33L*Atm33) + At33L*trKL) + - em4phi*(Ats33 - g33*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(T33 - g33*kthird*trS)); + em4phi*(Ats33 - g33*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(eTzzL - g33*kthird*trS)); CCTK_REAL const alpharhsL = (PDupwindNth1(alpha, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(alpha, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + diff --git a/ML_BSSN_MP/src/ML_BSSN_MP_RHS1.c b/ML_BSSN_MP/src/ML_BSSN_MP_RHS1.c index 6895e5c..b5d6f9e 100644 --- a/ML_BSSN_MP/src/ML_BSSN_MP_RHS1.c +++ b/ML_BSSN_MP/src/ML_BSSN_MP_RHS1.c @@ -111,8 +111,6 @@ void ML_BSSN_MP_RHS1_Body(cGH const * restrict const cctkGH, int const dir, int // CCTK_REAL Gt311 = INITVALUE, Gt312 = INITVALUE, Gt313 = INITVALUE, Gt322 = INITVALUE, Gt323 = INITVALUE, Gt333 = INITVALUE; // CCTK_REAL gtu11 = INITVALUE, gtu21 = INITVALUE, gtu22 = INITVALUE, gtu31 = INITVALUE, gtu32 = INITVALUE, gtu33 = INITVALUE; // CCTK_REAL S1 = INITVALUE, S2 = INITVALUE, S3 = INITVALUE; - // CCTK_REAL T01 = INITVALUE, T02 = INITVALUE, T03 = INITVALUE, T11 = INITVALUE, T12 = INITVALUE, T13 = INITVALUE; - // CCTK_REAL T22 = INITVALUE, T23 = INITVALUE, T33 = INITVALUE; // CCTK_REAL Xtn1 = INITVALUE, Xtn2 = INITVALUE, Xtn3 = INITVALUE; /* Declare local copies of grid functions */ @@ -493,29 +491,11 @@ void ML_BSSN_MP_RHS1_Body(cGH const * restrict const cctkGH, int const dir, int CCTK_REAL const Atu33 = Atm31*gtu31 + Atm32*gtu32 + Atm33*gtu33; - 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 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 phirhsL = PDupwindNth1(phi, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(phi, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + diff --git a/ML_BSSN_MP/src/ML_BSSN_MP_RHS2.c b/ML_BSSN_MP/src/ML_BSSN_MP_RHS2.c index d54bf6a..3745c32 100644 --- a/ML_BSSN_MP/src/ML_BSSN_MP_RHS2.c +++ b/ML_BSSN_MP/src/ML_BSSN_MP_RHS2.c @@ -138,8 +138,6 @@ void ML_BSSN_MP_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int // CCTK_REAL rho = INITVALUE; // CCTK_REAL Rphi11 = INITVALUE, Rphi12 = INITVALUE, Rphi13 = INITVALUE, Rphi22 = INITVALUE, Rphi23 = INITVALUE, Rphi33 = INITVALUE; // CCTK_REAL Rt11 = INITVALUE, Rt12 = INITVALUE, Rt13 = INITVALUE, Rt22 = INITVALUE, Rt23 = INITVALUE, Rt33 = 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 trAts = INITVALUE; // CCTK_REAL trS = INITVALUE; // CCTK_REAL Xtn1 = INITVALUE, Xtn2 = INITVALUE, Xtn3 = INITVALUE; @@ -952,31 +950,11 @@ void ML_BSSN_MP_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int CCTK_REAL const R33 = Rphi33 + Rt33; - 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 T02 = eTtyL; - - 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 trS = gu11*T11 + gu22*T22 + 2*(gu21*T12 + gu31*T13 + gu32*T23) + gu33*T33; + CCTK_REAL const trS = eTxxL*gu11 + eTyyL*gu22 + 2*(eTxyL*gu21 + eTxzL*gu31 + eTyzL*gu32) + eTzzL*gu33; CCTK_REAL const trKrhsL = PDupwindNth1(trK, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(trK, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1075,7 +1053,7 @@ void ML_BSSN_MP_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int J21L*(At11L*PDstandardNth2beta1 + At12L*PDstandardNth2beta2 + At13L*PDstandardNth2beta3) + J31L*(At11L*PDstandardNth3beta1 + At12L*PDstandardNth3beta2 + At13L*PDstandardNth3beta3)) + alphaL*(-2*(At11L*Atm11 + At12L*Atm21 + At13L*Atm31) + At11L*trKL) + - em4phi*(Ats11 - g11*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(T11 - g11*kthird*trS)); + em4phi*(Ats11 - g11*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(eTxxL - g11*kthird*trS)); CCTK_REAL const At12rhsL = PDupwindNth1(At12, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At12, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1091,7 +1069,7 @@ void ML_BSSN_MP_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int (At23L*J31L + At13L*J32L - 0.6666666666666666666666666666666666666667*At12L*J33L)*PDstandardNth3beta3 + alphaL*(-2.*(At11L*Atm12 + At12L*Atm22 + At13L*Atm32) + At12L*trKL) + em4phi*(Ats12 - 0.3333333333333333333333333333333333333333*g12*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T12 + 8.377580409572781969233715688745341024526*g12*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxyL + 8.377580409572781969233715688745341024526*g12*trS)); CCTK_REAL const At13rhsL = PDupwindNth1(At13, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At13, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1107,7 +1085,7 @@ void ML_BSSN_MP_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int (At33L*J31L + 0.3333333333333333333333333333333333333333*At13L*J33L)*PDstandardNth3beta3 + alphaL*(-2.*(At11L*Atm13 + At12L*Atm23 + At13L*Atm33) + At13L*trKL) + em4phi*(Ats13 - 0.3333333333333333333333333333333333333333*g13*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T13 + 8.377580409572781969233715688745341024526*g13*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTxzL + 8.377580409572781969233715688745341024526*g13*trS)); CCTK_REAL const At22rhsL = PDupwindNth1(At22, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At22, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1119,7 +1097,7 @@ void ML_BSSN_MP_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int J22L*(At12L*PDstandardNth2beta1 + At22L*PDstandardNth2beta2 + At23L*PDstandardNth2beta3) + J32L*(At12L*PDstandardNth3beta1 + At22L*PDstandardNth3beta2 + At23L*PDstandardNth3beta3)) + alphaL*(-2*(At12L*Atm12 + At22L*Atm22 + At23L*Atm32) + At22L*trKL) + - em4phi*(Ats22 - g22*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(T22 - g22*kthird*trS)); + em4phi*(Ats22 - g22*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(eTyyL - g22*kthird*trS)); CCTK_REAL const At23rhsL = PDupwindNth1(At23, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At23, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1135,7 +1113,7 @@ void ML_BSSN_MP_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int (At33L*J32L + 0.3333333333333333333333333333333333333333*At23L*J33L)*PDstandardNth3beta3 + alphaL*(-2.*(At12L*Atm13 + At22L*Atm23 + At23L*Atm33) + At23L*trKL) + em4phi*(Ats23 - 0.3333333333333333333333333333333333333333*g23*trAts + - alphaL*(-25.13274122871834590770114706623602307358*T23 + 8.377580409572781969233715688745341024526*g23*trS)); + alphaL*(-25.13274122871834590770114706623602307358*eTyzL + 8.377580409572781969233715688745341024526*g23*trS)); CCTK_REAL const At33rhsL = PDupwindNth1(At33, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(At33, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + @@ -1147,7 +1125,7 @@ void ML_BSSN_MP_RHS2_Body(cGH const * restrict const cctkGH, int const dir, int J23L*(At13L*PDstandardNth2beta1 + At23L*PDstandardNth2beta2 + At33L*PDstandardNth2beta3) + J33L*(At13L*PDstandardNth3beta1 + At23L*PDstandardNth3beta2 + At33L*PDstandardNth3beta3)) + alphaL*(-2*(At13L*Atm13 + At23L*Atm23 + At33L*Atm33) + At33L*trKL) + - em4phi*(Ats33 - g33*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(T33 - g33*kthird*trS)); + em4phi*(Ats33 - g33*kthird*trAts - 25.13274122871834590770114706623602307358*alphaL*(eTzzL - g33*kthird*trS)); CCTK_REAL const alpharhsL = (PDupwindNth1(alpha, i, j, k)*(beta1L*J11L + beta2L*J12L + beta3L*J13L) + PDupwindNth2(alpha, i, j, k)*(beta1L*J21L + beta2L*J22L + beta3L*J23L) + diff --git a/ML_BSSN_MP/src/ML_BSSN_MP_constraints.c b/ML_BSSN_MP/src/ML_BSSN_MP_constraints.c index 1c9b39b..3974724 100644 --- a/ML_BSSN_MP/src/ML_BSSN_MP_constraints.c +++ b/ML_BSSN_MP/src/ML_BSSN_MP_constraints.c @@ -117,8 +117,6 @@ void ML_BSSN_MP_constraints_Body(cGH const * restrict const cctkGH, int const di // CCTK_REAL Rphi11 = INITVALUE, Rphi12 = INITVALUE, Rphi13 = INITVALUE, Rphi22 = INITVALUE, Rphi23 = INITVALUE, Rphi33 = 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 trR = INITVALUE; /* Declare local copies of grid functions */ @@ -1035,35 +1033,15 @@ void ML_BSSN_MP_constraints_Body(cGH const * restrict const cctkGH, int const di 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 HL = -2.*(Atm12*Atm21 + Atm13*Atm31 + Atm23*Atm32) - 50.26548245743669181540229413247204614715*rho + trR - 1.*(SQR(Atm11) + SQR(Atm22) + SQR(Atm33)) + 0.6666666666666666666666666666666666666667*SQR(trKL); diff --git a/ML_BSSN_MP/src/ML_BSSN_MP_convertFromADMBase.c b/ML_BSSN_MP/src/ML_BSSN_MP_convertFromADMBase.c index 2601733..764bc4a 100644 --- a/ML_BSSN_MP/src/ML_BSSN_MP_convertFromADMBase.c +++ b/ML_BSSN_MP/src/ML_BSSN_MP_convertFromADMBase.c @@ -112,7 +112,6 @@ void ML_BSSN_MP_convertFromADMBase_Body(cGH const * restrict const cctkGH, int c // CCTK_REAL betaxL = INITVALUE; // CCTK_REAL betayL = INITVALUE; // CCTK_REAL betazL = INITVALUE; - // CCTK_REAL etaL = INITVALUE; // CCTK_REAL gt11L = INITVALUE, gt12L = INITVALUE, gt13L = INITVALUE, gt22L = INITVALUE, gt23L = INITVALUE, gt33L = INITVALUE; // CCTK_REAL gxxL = INITVALUE; // CCTK_REAL gxyL = INITVALUE; @@ -223,8 +222,6 @@ void ML_BSSN_MP_convertFromADMBase_Body(cGH const * restrict const cctkGH, int c CCTK_REAL const beta3L = betazL; - CCTK_REAL const etaL = BetaDriver; - /* Copy local copies back to grid functions */ alpha[index] = alphaL; @@ -237,7 +234,6 @@ void ML_BSSN_MP_convertFromADMBase_Body(cGH const * restrict const cctkGH, int c beta1[index] = beta1L; beta2[index] = beta2L; beta3[index] = beta3L; - eta[index] = etaL; gt11[index] = gt11L; gt12[index] = gt12L; gt13[index] = gt13L; diff --git a/ML_BSSN_MP/src/ML_BSSN_MP_setBetaDriver.c b/ML_BSSN_MP/src/ML_BSSN_MP_setBetaDriverSpatial.c index ec8e6c4..d517184 100644 --- a/ML_BSSN_MP/src/ML_BSSN_MP_setBetaDriver.c +++ b/ML_BSSN_MP/src/ML_BSSN_MP_setBetaDriverSpatial.c @@ -20,7 +20,7 @@ #define CUB(x) ((x) * (x) * (x)) #define QAD(x) ((x) * (x) * (x) * (x)) -void ML_BSSN_MP_setBetaDriver_Body(cGH const * restrict const cctkGH, int const dir, int const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], int const min[3], int const max[3], int const n_subblock_gfs, CCTK_REAL * restrict const subblock_gfs[]) +void ML_BSSN_MP_setBetaDriverSpatial_Body(cGH const * restrict const cctkGH, int const dir, int const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], int const min[3], int const max[3], int const n_subblock_gfs, CCTK_REAL * restrict const subblock_gfs[]) { DECLARE_CCTK_ARGUMENTS; DECLARE_CCTK_PARAMETERS; @@ -44,10 +44,10 @@ void ML_BSSN_MP_setBetaDriver_Body(cGH const * restrict const cctkGH, int const if (verbose > 1) { - CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_BSSN_MP_setBetaDriver_Body"); + CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_BSSN_MP_setBetaDriverSpatial_Body"); } - if (cctk_iteration % ML_BSSN_MP_setBetaDriver_calc_every != ML_BSSN_MP_setBetaDriver_calc_offset) + if (cctk_iteration % ML_BSSN_MP_setBetaDriverSpatial_calc_every != ML_BSSN_MP_setBetaDriverSpatial_calc_offset) { return; } @@ -89,7 +89,7 @@ void ML_BSSN_MP_setBetaDriver_Body(cGH const * restrict const cctkGH, int const /* Loop over the grid points */ #pragma omp parallel - LC_LOOP3 (ML_BSSN_MP_setBetaDriver, + LC_LOOP3 (ML_BSSN_MP_setBetaDriverSpatial, i,j,k, min[0],min[1],min[2], max[0],max[1],max[2], cctk_lsh[0],cctk_lsh[1],cctk_lsh[2]) { @@ -127,13 +127,13 @@ void ML_BSSN_MP_setBetaDriver_Body(cGH const * restrict const cctkGH, int const /* Copy local copies back to subblock grid functions */ } - LC_ENDLOOP3 (ML_BSSN_MP_setBetaDriver); + LC_ENDLOOP3 (ML_BSSN_MP_setBetaDriverSpatial); } -void ML_BSSN_MP_setBetaDriver(CCTK_ARGUMENTS) +void ML_BSSN_MP_setBetaDriverSpatial(CCTK_ARGUMENTS) { DECLARE_CCTK_ARGUMENTS; DECLARE_CCTK_PARAMETERS; - GenericFD_LoopOverEverything(cctkGH, &ML_BSSN_MP_setBetaDriver_Body); + GenericFD_LoopOverEverything(cctkGH, &ML_BSSN_MP_setBetaDriverSpatial_Body); } diff --git a/ML_BSSN_MP/src/make.code.defn b/ML_BSSN_MP/src/make.code.defn index d89c879..792fefe 100644 --- a/ML_BSSN_MP/src/make.code.defn +++ b/ML_BSSN_MP/src/make.code.defn @@ -1,3 +1,3 @@ # File produced by Kranc -SRCS = Startup.c RegisterMoL.c RegisterSymmetries.c ML_BSSN_MP_Minkowski.c ML_BSSN_MP_convertFromADMBase.c ML_BSSN_MP_convertFromADMBaseGamma.c ML_BSSN_MP_setBetaDriver.c ML_BSSN_MP_RHS.c ML_BSSN_MP_RHS1.c ML_BSSN_MP_RHS2.c ML_BSSN_MP_RHSStaticBoundary.c ML_BSSN_MP_RHSRadiativeBoundary.c ML_BSSN_MP_enforce.c ML_BSSN_MP_enforce2.c ML_BSSN_MP_boundary.c ML_BSSN_MP_convertToADMBase.c ML_BSSN_MP_convertToADMBaseDtLapseShift.c ML_BSSN_MP_convertToADMBaseDtLapseShiftBoundary.c ML_BSSN_MP_convertToADMBaseFakeDtLapseShift.c ML_BSSN_MP_constraints.c ML_BSSN_MP_constraints_boundary.c Boundaries.c +SRCS = Startup.c RegisterMoL.c RegisterSymmetries.c ML_BSSN_MP_Minkowski.c ML_BSSN_MP_convertFromADMBase.c ML_BSSN_MP_convertFromADMBaseGamma.c ML_BSSN_MP_setBetaDriverConstant.c ML_BSSN_MP_setBetaDriverSpatial.c ML_BSSN_MP_RHS.c ML_BSSN_MP_RHS1.c ML_BSSN_MP_RHS2.c ML_BSSN_MP_RHSStaticBoundary.c ML_BSSN_MP_RHSRadiativeBoundary.c ML_BSSN_MP_enforce.c ML_BSSN_MP_enforce2.c ML_BSSN_MP_boundary.c ML_BSSN_MP_convertToADMBase.c ML_BSSN_MP_convertToADMBaseDtLapseShift.c ML_BSSN_MP_convertToADMBaseDtLapseShiftBoundary.c ML_BSSN_MP_convertToADMBaseFakeDtLapseShift.c ML_BSSN_MP_constraints.c ML_BSSN_MP_constraints_boundary.c Boundaries.c diff --git a/ML_BSSN_MP_Helper/schedule.ccl b/ML_BSSN_MP_Helper/schedule.ccl index cfec9cd..ee198a6 100644 --- a/ML_BSSN_MP_Helper/schedule.ccl +++ b/ML_BSSN_MP_Helper/schedule.ccl @@ -56,6 +56,20 @@ if (CCTK_EQUALS (evolution_method, "ML_BSSN_MP")) { + SCHEDULE GROUP ML_BSSN_MP_InitEta AT initial + { + } "Initialize BetaDriver" + + SCHEDULE GROUP ML_BSSN_MP_InitEta AT postregrid + { + } "Initialize BetaDriver" + + SCHEDULE GROUP ML_BSSN_MP_InitEta AT post_recover_variables + { + } "Initialize BetaDriver" + + + #SCHEDULE GROUP ML_BSSN_MP_evolCalcGroup AT postinitial AFTER MoL_PostStep #{ #} "Calculate BSSN RHS" @@ -144,7 +158,11 @@ if (CCTK_EQUALS (evolution_method, "ML_BSSN_MP")) { TRIGGERS: ML_BSSN_MP::ML_mom } "Calculate ADM variables" } - + + SCHEDULE GROUP ML_BSSN_MP_convertToADMBaseGroupWrapper AT CCTK_POST_RECOVER_VARIABLES + { + } "Calculate ADM variables" + SCHEDULE ML_BSSN_MP_SelectBCsADMBase IN ML_BSSN_MP_convertToADMBaseGroupWrapper AFTER ML_BSSN_MP_convertToADMBaseGroup { LANG: C @@ -162,5 +180,4 @@ if (CCTK_EQUALS (evolution_method, "ML_BSSN_MP")) { TRIGGERS: ML_BSSN_MP::ML_Ham TRIGGERS: ML_BSSN_MP::ML_mom } "Calculate BSSN constraints" - } diff --git a/ML_BSSN_M_Helper/schedule.ccl b/ML_BSSN_M_Helper/schedule.ccl new file mode 100644 index 0000000..6600156 --- /dev/null +++ b/ML_BSSN_M_Helper/schedule.ccl @@ -0,0 +1,180 @@ +if (CCTK_EQUALS (evolution_method, "ML_BSSN_M")) { + + if (timelevels == 1) { + STORAGE: ADMBase::metric[1] + STORAGE: ADMBase::curv[1] + STORAGE: ADMBase::lapse[1] + STORAGE: ADMBase::shift[1] + STORAGE: ADMBase::dtlapse[1] + STORAGE: ADMBase::dtshift[1] + } else if (timelevels == 2) { + STORAGE: ADMBase::metric[2] + STORAGE: ADMBase::curv[2] + STORAGE: ADMBase::lapse[2] + STORAGE: ADMBase::shift[2] + STORAGE: ADMBase::dtlapse[2] + STORAGE: ADMBase::dtshift[2] + } else if (timelevels == 3) { + STORAGE: ADMBase::metric[3] + STORAGE: ADMBase::curv[3] + STORAGE: ADMBase::lapse[3] + STORAGE: ADMBase::shift[3] + STORAGE: ADMBase::dtlapse[3] + STORAGE: ADMBase::dtshift[3] + } else if (timelevels == 4) { + #STORAGE: ADMBase::metric[4] + #STORAGE: ADMBase::curv[4] + #STORAGE: ADMBase::lapse[4] + #STORAGE: ADMBase::shift[4] + #STORAGE: ADMBase::dtlapse[4] + #STORAGE: ADMBase::dtshift[4] + STORAGE: ADMBase::metric[3] + STORAGE: ADMBase::curv[3] + STORAGE: ADMBase::lapse[3] + STORAGE: ADMBase::shift[3] + STORAGE: ADMBase::dtlapse[3] + STORAGE: ADMBase::dtshift[3] + } + + SCHEDULE ML_BSSN_M_RegisterSlicing AT startup + { + LANG: C + OPTIONS: meta + } "Register slicing" + + SCHEDULE ML_BSSN_M_SetGroupTags AT startup BEFORE Driver_Startup + { + LANG: C + OPTIONS: meta + } "Set checkpointing and prolongation group tags" + + SCHEDULE ML_BSSN_M_RegisterConstrained IN MoL_Register + { + LANG: C + OPTIONS: meta + } "Register ADMBase variables as constrained" + + + + #SCHEDULE GROUP ML_BSSN_M_evolCalcGroup AT postinitial AFTER MoL_PostStep + #{ + #} "Calculate BSSN RHS" + SCHEDULE GROUP MoL_CalcRHS AT postinitial AFTER MoL_PostStep + { + } "Evaluate RHS" + + SCHEDULE GROUP ML_BSSN_M_evolCalcGroup IN MoL_CalcRHS + { + } "Calculate BSSN RHS" + + SCHEDULE GROUP ML_BSSN_M_evolCalcGroup AT analysis + { + TRIGGERS: ML_BSSN_M::ML_log_confacrhs + TRIGGERS: ML_BSSN_M::ML_metricrhs + TRIGGERS: ML_BSSN_M::ML_Gammarhs + TRIGGERS: ML_BSSN_M::ML_trace_curvrhs + TRIGGERS: ML_BSSN_M::ML_curvrhs + TRIGGERS: ML_BSSN_M::ML_lapserhs + TRIGGERS: ML_BSSN_M::ML_dtlapserhs + TRIGGERS: ML_BSSN_M::ML_shiftrhs + TRIGGERS: ML_BSSN_M::ML_dtshiftrhs + } "Calculate BSSN RHS" + + + + if (CCTK_EQUALS (my_initial_boundary_condition, "extrapolate-gammas")) + { + SCHEDULE ML_BSSN_M_ExtrapolateGammas AT initial AFTER ML_BSSN_M_convertFromADMBaseGamma + { + LANG: C + SYNC: ML_Gamma + SYNC: ML_dtlapse + SYNC: ML_dtshift + } "Extrapolate Gammas and time derivatives of lapse and shift" + } + + if (CCTK_EQUALS (my_rhs_boundary_condition, "NewRad")) + { + SCHEDULE ML_BSSN_M_NewRad IN ML_BSSN_M_evolCalcGroup AFTER ML_BSSN_M_RHS + { + LANG: C + #SYNC: ML_curvrhs + #SYNC: ML_dtlapserhs + #SYNC: ML_dtshiftrhs + #SYNC: ML_Gammarhs + #SYNC: ML_lapserhs + #SYNC: ML_log_confacrhs + #SYNC: ML_metricrhs + #SYNC: ML_shiftrhs + #SYNC: ML_trace_curvrhs + } "Apply NewRad boundary conditions to RHS" + } + + + + SCHEDULE GROUP ML_BSSN_M_convertToADMBaseGroup IN ML_BSSN_M_convertToADMBaseGroupWrapper + { + } "Calculate ADM variables" + + if (CCTK_EQUALS (calculate_ADMBase_variables_at, "MoL_PostStep")) + { + # if (timelevels > 1) + # { + # SCHEDULE ML_BSSN_M_CopyADMBase AT prestep + # { + # LANG: C + # } "Copy ADMBase variables to current time level" + # } + SCHEDULE GROUP ML_BSSN_M_convertToADMBaseGroupWrapper IN MoL_PostStep AFTER (ML_BSSN_M_ApplyBCs ML_BSSN_M_enforce) BEFORE (ADMBase_SetADMVars Whisky_PostStep) + { + } "Calculate ADM variables" + } + else if (CCTK_EQUALS (calculate_ADMBase_variables_at, "CCTK_EVOL")) + { + SCHEDULE GROUP ML_BSSN_M_convertToADMBaseGroupWrapper AT evol AFTER MoL_Evolution BEFORE (ADMBase_SetADMVars Whisky_PostStep) + { + } "Calculate ADM variables" + } + else if (CCTK_EQUALS (calculate_ADMBase_variables_at, "CCTK_ANALYSIS")) + { + SCHEDULE GROUP ML_BSSN_M_convertToADMBaseGroupWrapper AT analysis BEFORE (ADMBase_SetADMVars Whisky_PostStep) + { + TRIGGERS: ML_BSSN_M::ML_Ham + TRIGGERS: ML_BSSN_M::ML_mom + } "Calculate ADM variables" + } + + SCHEDULE GROUP ML_BSSN_M_convertToADMBaseGroupWrapper AT CCTK_POST_RECOVER_VARIABLES + { + } "Calculate ADM variables" + + SCHEDULE ML_BSSN_M_SelectBCsADMBase IN ML_BSSN_M_convertToADMBaseGroupWrapper AFTER ML_BSSN_M_convertToADMBaseGroup + { + LANG: C + OPTIONS: level + } "Apply boundary conditions to ADMBase variables" + + SCHEDULE GROUP ApplyBCs AS ML_BSSN_M_ApplyBCsADMBase IN ML_BSSN_M_convertToADMBaseGroupWrapper AFTER ML_BSSN_M_SelectBCsADMBase + { + } "Apply boundary conditions to ADMBase variables" + + + + SCHEDULE GROUP ML_BSSN_M_constraintsCalcGroup AT analysis + { + TRIGGERS: ML_BSSN_M::ML_Ham + TRIGGERS: ML_BSSN_M::ML_mom + } "Calculate BSSN constraints" + + SCHEDULE GROUP ML_BSSN_M_init_eta AT initial + { + } "Initialize BetaDriver" + + SCHEDULE GROUP ML_BSSN_M_init_eta AT postregrid + { + } "Initialize BetaDriver" + + SCHEDULE GROUP ML_BSSN_M_init_eta AT post_recover_variables + { + } "Initialize BetaDriver" +} diff --git a/ML_BSSN_M_Helper/src/SetGroupTags.c b/ML_BSSN_M_Helper/src/SetGroupTags.c new file mode 100644 index 0000000..79e02fb --- /dev/null +++ b/ML_BSSN_M_Helper/src/SetGroupTags.c @@ -0,0 +1,70 @@ +#include <cctk.h> +#include <cctk_Parameters.h> +#include <util_Table.h> + +#include <assert.h> + +static void +set_group_tags (int const checkpoint, int const prolongate, + char const * restrict const gn); + +int +ML_BSSN_M_SetGroupTags (void) +{ + DECLARE_CCTK_PARAMETERS; + + set_group_tags (0, 1, "ADMBase::metric"); + set_group_tags (0, 1, "ADMBase::curv"); + set_group_tags (0, 1, "ADMBase::lapse"); + set_group_tags (0, 1, "ADMBase::shift"); + set_group_tags (0, 1, "ADMBase::dtlapse"); + set_group_tags (0, 1, "ADMBase::dtshift"); + + set_group_tags (0, 0, "ML_BSSN_M::ML_cons_detg"); + set_group_tags (0, 0, "ML_BSSN_M::ML_cons_Gamma"); + set_group_tags (0, 0, "ML_BSSN_M::ML_cons_traceA"); + set_group_tags (0, 0, "ML_BSSN_M::ML_Ham"); + set_group_tags (0, 0, "ML_BSSN_M::ML_mom"); + set_group_tags (0, 0, "ML_BSSN_M::ML_curvrhs"); + set_group_tags (0, 0, "ML_BSSN_M::ML_BetaDriver"); + + int const checkpoint = rhs_timelevels > 1; + set_group_tags (checkpoint, 0, "ML_BSSN_M::ML_dtlapserhs"); + set_group_tags (checkpoint, 0, "ML_BSSN_M::ML_dtshiftrhs"); + set_group_tags (checkpoint, 0, "ML_BSSN_M::ML_Gammarhs"); + set_group_tags (checkpoint, 0, "ML_BSSN_M::ML_lapserhs"); + set_group_tags (checkpoint, 0, "ML_BSSN_M::ML_log_confacrhs"); + set_group_tags (checkpoint, 0, "ML_BSSN_M::ML_metricrhs"); + set_group_tags (checkpoint, 0, "ML_BSSN_M::ML_shiftrhs"); + set_group_tags (checkpoint, 0, "ML_BSSN_M::ML_trace_curvrhs"); + + return 0; +} + +static void +set_group_tags (int const checkpoint, int const prolongate, + char const * restrict const gn) +{ + assert (gn); + + int const gi = CCTK_GroupIndex (gn); + assert (gi >= 0); + + int const table = CCTK_GroupTagsTableI (gi); + assert (table >= 0); + + if (! checkpoint) { + int ierr; + ierr = Util_TableSetString (table, "no", "Checkpoint"); + assert (! ierr); + + ierr = Util_TableSetString (table, "no", "Persistent"); + assert (! ierr); + } + + if (! prolongate) { + int ierr; + ierr = Util_TableSetString (table, "none", "Prolongation"); + assert (! ierr); + } +} diff --git a/m/McLachlan_BSSN.m b/m/McLachlan_BSSN.m index 1b4d19f..4c31040 100644 --- a/m/McLachlan_BSSN.m +++ b/m/McLachlan_BSSN.m @@ -184,10 +184,10 @@ admdtalpha=dtalp; admbeta1=betax; admbeta2=betay; admbeta3=betaz; admdtbeta1=dtbetax; admdtbeta2=dtbetay; admdtbeta3=dtbetaz; -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 *) @@ -235,7 +235,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}}, @@ -244,7 +245,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, @@ -274,8 +274,7 @@ initialCalc = alpha -> 1, A -> 0, beta[ua] -> 0, - B[ua] -> 0, - eta -> BetaDriver + B[ua] -> 0 } }; @@ -318,9 +317,7 @@ convertFromADMBaseCalc = alpha -> admalpha, - beta[ua] -> admbeta[ua], - - eta -> BetaDriver + beta[ua] -> admbeta[ua] } }; @@ -364,10 +361,21 @@ convertFromADMBaseGammaCalc = } }; -setBetaDriverCalc = +setBetaDriverConstantCalc = +{ + Name -> BSSN <> "_setBetaDriverConstant", + Schedule -> {"IN "<> BSSN <> "_InitEta"}, + ConditionalOnKeyword -> {"UseSpatialBetaDriver", "no"}, + Equations -> + { + eta -> BetaDriver + } +}; + +setBetaDriverSpatialCalc = { - Name -> BSSN <> "_setBetaDriver", - Schedule -> {"AT initial AFTER ADMBase_PostInitial AFTER " <> BSSN <> "_convertFromADMBase"}, + Name -> BSSN <> "_setBetaDriverSpatial", + Schedule -> {"IN "<> BSSN <> "_InitEta"}, ConditionalOnKeyword -> {"UseSpatialBetaDriver", "yes"}, Equations -> { @@ -470,7 +478,7 @@ evolCalc = Atm[ua,lb], Atu[ua,ub], e4phi, em4phi, cdphi[la], cdphi2[la,lb], g[la,lb], detg, gu[ua,ub], G[ua,lb,lc], Ats[la,lb], trAts, - T00, T0[la], T[la,lb], rho, S[la], trS, fac1, fac2}, + rho, S[la], trS, fac1, fac2}, Equations -> { dir[ua] -> Sign[beta[ua]], @@ -524,6 +532,7 @@ evolCalc = (* Matter terms *) + (* T00 -> addMatter eTtt, T01 -> addMatter eTtx, T02 -> addMatter eTty, @@ -534,6 +543,7 @@ evolCalc = T22 -> addMatter eTyy, T23 -> addMatter eTyz, T33 -> addMatter eTzz, + *) (* rho = n^a n^b T_ab *) rho -> addMatter @@ -627,7 +637,7 @@ evol1Calc = Gt[ua,lb,lc], Xtn[ua], Rt[la,lb], Rphi[la,lb], R[la,lb], Atm[ua,lb], Atu[ua,ub], e4phi, em4phi, cdphi[la], cdphi2[la,lb], g[la,lb], detg, - T00, T0[la], T[la,lb], rho, S[la], trS, fac1, fac2}, + rho, S[la], trS, fac1, fac2}, Equations -> { dir[ua] -> Sign[beta[ua]], @@ -651,6 +661,7 @@ evol1Calc = (* Matter terms *) + (* T01 -> addMatter eTtx, T02 -> addMatter eTty, T03 -> addMatter eTtz, @@ -660,6 +671,7 @@ evol1Calc = 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])), @@ -720,7 +732,7 @@ evol2Calc = Atm[ua,lb], e4phi, em4phi, cdphi[la], cdphi2[la,lb], g[la,lb], detg, gu[ua,ub], G[ua,lb,lc], Ats[la,lb], trAts, - T00, T0[la], T[la,lb], rho, S[la], trS, fac1, fac2}, + rho, S[la], trS, fac1, fac2}, Equations -> { dir[ua] -> Sign[beta[ua]], @@ -774,6 +786,7 @@ evol2Calc = (* Matter terms *) + (* T00 -> addMatter eTtt, T01 -> addMatter eTtx, T02 -> addMatter eTty, @@ -784,6 +797,7 @@ evol2Calc = T22 -> addMatter eTyy, T23 -> addMatter eTyz, T33 -> addMatter eTzz, + *) (* rho = n^a n^b T_ab *) rho -> addMatter @@ -952,7 +966,7 @@ constraintsCalc = g[la,lb], detg, gu[ua,ub], ddetg[la], G[ua,lb,lc], Rt[la,lb], Rphi[la,lb], R[la,lb], trR, Atm[la,lb], gK[la,lb,lc], cdphi[la], cdphi2[la,lb], - T00, T0[la], T[la,lb], rho, S[la], fac1, fac2}, + rho, S[la], fac1, fac2}, Equations -> { detgt -> 1 (* detgtExpr *), @@ -1031,6 +1045,7 @@ constraintsCalc = (* Matter terms *) + (* T00 -> eTtt, T01 -> eTtx, T02 -> eTty, @@ -1041,6 +1056,7 @@ constraintsCalc = 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]), @@ -1259,7 +1275,8 @@ calculations = initialCalc, convertFromADMBaseCalc, convertFromADMBaseGammaCalc, - setBetaDriverCalc, + setBetaDriverConstantCalc, + setBetaDriverSpatialCalc, evolCalc, evol1Calc, evol2Calc, RHSStaticBoundaryCalc, diff --git a/m/prototype/ML_BSSN_Helper/schedule.ccl b/m/prototype/ML_BSSN_Helper/schedule.ccl index 48ebedc..6b6c4e8 100644 --- a/m/prototype/ML_BSSN_Helper/schedule.ccl +++ b/m/prototype/ML_BSSN_Helper/schedule.ccl @@ -56,6 +56,20 @@ if (CCTK_EQUALS (evolution_method, "ML_BSSN")) { + SCHEDULE GROUP ML_BSSN_InitEta AT initial + { + } "Initialize BetaDriver" + + SCHEDULE GROUP ML_BSSN_InitEta AT postregrid + { + } "Initialize BetaDriver" + + SCHEDULE GROUP ML_BSSN_InitEta AT post_recover_variables + { + } "Initialize BetaDriver" + + + #SCHEDULE GROUP ML_BSSN_evolCalcGroup AT postinitial AFTER MoL_PostStep #{ #} "Calculate BSSN RHS" @@ -144,7 +158,11 @@ if (CCTK_EQUALS (evolution_method, "ML_BSSN")) { TRIGGERS: ML_BSSN::ML_mom } "Calculate ADM variables" } - + + SCHEDULE GROUP ML_BSSN_convertToADMBaseGroupWrapper AT CCTK_POST_RECOVER_VARIABLES + { + } "Calculate ADM variables" + SCHEDULE ML_BSSN_SelectBCsADMBase IN ML_BSSN_convertToADMBaseGroupWrapper AFTER ML_BSSN_convertToADMBaseGroup { LANG: C @@ -162,5 +180,4 @@ if (CCTK_EQUALS (evolution_method, "ML_BSSN")) { TRIGGERS: ML_BSSN::ML_Ham TRIGGERS: ML_BSSN::ML_mom } "Calculate BSSN constraints" - } |