diff options
Diffstat (limited to 'ML_BSSN/src/ML_BSSN_Advect.cc')
| -rw-r--r-- | ML_BSSN/src/ML_BSSN_Advect.cc | 79 |
1 files changed, 78 insertions, 1 deletions
diff --git a/ML_BSSN/src/ML_BSSN_Advect.cc b/ML_BSSN/src/ML_BSSN_Advect.cc index 2543f21..15053d6 100644 --- a/ML_BSSN/src/ML_BSSN_Advect.cc +++ b/ML_BSSN/src/ML_BSSN_Advect.cc @@ -53,6 +53,9 @@ extern "C" void ML_BSSN_Advect_SelectBCs(CCTK_ARGUMENTS) ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, GenericFD_GetBoundaryWidth(cctkGH), -1 /* no table */, "ML_BSSN::ML_shiftrhs","flat"); if (ierr < 0) CCTK_WARN(1, "Failed to register flat BC for ML_BSSN::ML_shiftrhs."); + ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, GenericFD_GetBoundaryWidth(cctkGH), -1 /* no table */, "ML_BSSN::ML_Thetarhs","flat"); + if (ierr < 0) + CCTK_WARN(1, "Failed to register flat BC for ML_BSSN::ML_Thetarhs."); ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, GenericFD_GetBoundaryWidth(cctkGH), -1 /* no table */, "ML_BSSN::ML_trace_curvrhs","flat"); if (ierr < 0) CCTK_WARN(1, "Failed to register flat BC for ML_BSSN::ML_trace_curvrhs."); @@ -277,6 +280,8 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir CCTK_REAL_VEC gt33rhsL = vec_load(gt33rhs[index]); CCTK_REAL_VEC phiL = vec_load(phi[index]); CCTK_REAL_VEC phirhsL = vec_load(phirhs[index]); + CCTK_REAL_VEC ThetaL = vec_load(Theta[index]); + CCTK_REAL_VEC ThetarhsL = vec_load(Thetarhs[index]); CCTK_REAL_VEC trKL = vec_load(trK[index]); CCTK_REAL_VEC trKrhsL = vec_load(trKrhs[index]); CCTK_REAL_VEC Xt1L = vec_load(Xt1[index]); @@ -431,6 +436,12 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir CCTK_REAL_VEC PDupwindNthSymm2phi; CCTK_REAL_VEC PDupwindNthAnti3phi; CCTK_REAL_VEC PDupwindNthSymm3phi; + CCTK_REAL_VEC PDupwindNthAnti1Theta; + CCTK_REAL_VEC PDupwindNthSymm1Theta; + CCTK_REAL_VEC PDupwindNthAnti2Theta; + CCTK_REAL_VEC PDupwindNthSymm2Theta; + CCTK_REAL_VEC PDupwindNthAnti3Theta; + CCTK_REAL_VEC PDupwindNthSymm3Theta; CCTK_REAL_VEC PDupwindNthAnti1trK; CCTK_REAL_VEC PDupwindNthSymm1trK; CCTK_REAL_VEC PDupwindNthAnti2trK; @@ -585,6 +596,12 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir PDupwindNthSymm2phi = PDupwindNthSymmfdOrder22(&phi[index]); PDupwindNthAnti3phi = PDupwindNthAntifdOrder23(&phi[index]); PDupwindNthSymm3phi = PDupwindNthSymmfdOrder23(&phi[index]); + PDupwindNthAnti1Theta = PDupwindNthAntifdOrder21(&Theta[index]); + PDupwindNthSymm1Theta = PDupwindNthSymmfdOrder21(&Theta[index]); + PDupwindNthAnti2Theta = PDupwindNthAntifdOrder22(&Theta[index]); + PDupwindNthSymm2Theta = PDupwindNthSymmfdOrder22(&Theta[index]); + PDupwindNthAnti3Theta = PDupwindNthAntifdOrder23(&Theta[index]); + PDupwindNthSymm3Theta = PDupwindNthSymmfdOrder23(&Theta[index]); PDupwindNthAnti1trK = PDupwindNthAntifdOrder21(&trK[index]); PDupwindNthSymm1trK = PDupwindNthSymmfdOrder21(&trK[index]); PDupwindNthAnti2trK = PDupwindNthAntifdOrder22(&trK[index]); @@ -738,6 +755,12 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir PDupwindNthSymm2phi = PDupwindNthSymmfdOrder42(&phi[index]); PDupwindNthAnti3phi = PDupwindNthAntifdOrder43(&phi[index]); PDupwindNthSymm3phi = PDupwindNthSymmfdOrder43(&phi[index]); + PDupwindNthAnti1Theta = PDupwindNthAntifdOrder41(&Theta[index]); + PDupwindNthSymm1Theta = PDupwindNthSymmfdOrder41(&Theta[index]); + PDupwindNthAnti2Theta = PDupwindNthAntifdOrder42(&Theta[index]); + PDupwindNthSymm2Theta = PDupwindNthSymmfdOrder42(&Theta[index]); + PDupwindNthAnti3Theta = PDupwindNthAntifdOrder43(&Theta[index]); + PDupwindNthSymm3Theta = PDupwindNthSymmfdOrder43(&Theta[index]); PDupwindNthAnti1trK = PDupwindNthAntifdOrder41(&trK[index]); PDupwindNthSymm1trK = PDupwindNthSymmfdOrder41(&trK[index]); PDupwindNthAnti2trK = PDupwindNthAntifdOrder42(&trK[index]); @@ -891,6 +914,12 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir PDupwindNthSymm2phi = PDupwindNthSymmfdOrder62(&phi[index]); PDupwindNthAnti3phi = PDupwindNthAntifdOrder63(&phi[index]); PDupwindNthSymm3phi = PDupwindNthSymmfdOrder63(&phi[index]); + PDupwindNthAnti1Theta = PDupwindNthAntifdOrder61(&Theta[index]); + PDupwindNthSymm1Theta = PDupwindNthSymmfdOrder61(&Theta[index]); + PDupwindNthAnti2Theta = PDupwindNthAntifdOrder62(&Theta[index]); + PDupwindNthSymm2Theta = PDupwindNthSymmfdOrder62(&Theta[index]); + PDupwindNthAnti3Theta = PDupwindNthAntifdOrder63(&Theta[index]); + PDupwindNthSymm3Theta = PDupwindNthSymmfdOrder63(&Theta[index]); PDupwindNthAnti1trK = PDupwindNthAntifdOrder61(&trK[index]); PDupwindNthSymm1trK = PDupwindNthSymmfdOrder61(&trK[index]); PDupwindNthAnti2trK = PDupwindNthAntifdOrder62(&trK[index]); @@ -1044,6 +1073,12 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir PDupwindNthSymm2phi = PDupwindNthSymmfdOrder82(&phi[index]); PDupwindNthAnti3phi = PDupwindNthAntifdOrder83(&phi[index]); PDupwindNthSymm3phi = PDupwindNthSymmfdOrder83(&phi[index]); + PDupwindNthAnti1Theta = PDupwindNthAntifdOrder81(&Theta[index]); + PDupwindNthSymm1Theta = PDupwindNthSymmfdOrder81(&Theta[index]); + PDupwindNthAnti2Theta = PDupwindNthAntifdOrder82(&Theta[index]); + PDupwindNthSymm2Theta = PDupwindNthSymmfdOrder82(&Theta[index]); + PDupwindNthAnti3Theta = PDupwindNthAntifdOrder83(&Theta[index]); + PDupwindNthSymm3Theta = PDupwindNthSymmfdOrder83(&Theta[index]); PDupwindNthAnti1trK = PDupwindNthAntifdOrder81(&trK[index]); PDupwindNthSymm1trK = PDupwindNthSymmfdOrder81(&trK[index]); PDupwindNthAnti2trK = PDupwindNthAntifdOrder82(&trK[index]); @@ -1099,6 +1134,7 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir CCTK_REAL_VEC JacPDupwindNthAnti1gt23; CCTK_REAL_VEC JacPDupwindNthAnti1gt33; CCTK_REAL_VEC JacPDupwindNthAnti1phi; + CCTK_REAL_VEC JacPDupwindNthAnti1Theta; CCTK_REAL_VEC JacPDupwindNthAnti1trK; CCTK_REAL_VEC JacPDupwindNthAnti1Xt1; CCTK_REAL_VEC JacPDupwindNthAnti1Xt2; @@ -1124,6 +1160,7 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir CCTK_REAL_VEC JacPDupwindNthAnti2gt23; CCTK_REAL_VEC JacPDupwindNthAnti2gt33; CCTK_REAL_VEC JacPDupwindNthAnti2phi; + CCTK_REAL_VEC JacPDupwindNthAnti2Theta; CCTK_REAL_VEC JacPDupwindNthAnti2trK; CCTK_REAL_VEC JacPDupwindNthAnti2Xt1; CCTK_REAL_VEC JacPDupwindNthAnti2Xt2; @@ -1149,6 +1186,7 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir CCTK_REAL_VEC JacPDupwindNthAnti3gt23; CCTK_REAL_VEC JacPDupwindNthAnti3gt33; CCTK_REAL_VEC JacPDupwindNthAnti3phi; + CCTK_REAL_VEC JacPDupwindNthAnti3Theta; CCTK_REAL_VEC JacPDupwindNthAnti3trK; CCTK_REAL_VEC JacPDupwindNthAnti3Xt1; CCTK_REAL_VEC JacPDupwindNthAnti3Xt2; @@ -1174,6 +1212,7 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir CCTK_REAL_VEC JacPDupwindNthSymm1gt23; CCTK_REAL_VEC JacPDupwindNthSymm1gt33; CCTK_REAL_VEC JacPDupwindNthSymm1phi; + CCTK_REAL_VEC JacPDupwindNthSymm1Theta; CCTK_REAL_VEC JacPDupwindNthSymm1trK; CCTK_REAL_VEC JacPDupwindNthSymm1Xt1; CCTK_REAL_VEC JacPDupwindNthSymm1Xt2; @@ -1199,6 +1238,7 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir CCTK_REAL_VEC JacPDupwindNthSymm2gt23; CCTK_REAL_VEC JacPDupwindNthSymm2gt33; CCTK_REAL_VEC JacPDupwindNthSymm2phi; + CCTK_REAL_VEC JacPDupwindNthSymm2Theta; CCTK_REAL_VEC JacPDupwindNthSymm2trK; CCTK_REAL_VEC JacPDupwindNthSymm2Xt1; CCTK_REAL_VEC JacPDupwindNthSymm2Xt2; @@ -1224,6 +1264,7 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir CCTK_REAL_VEC JacPDupwindNthSymm3gt23; CCTK_REAL_VEC JacPDupwindNthSymm3gt33; CCTK_REAL_VEC JacPDupwindNthSymm3phi; + CCTK_REAL_VEC JacPDupwindNthSymm3Theta; CCTK_REAL_VEC JacPDupwindNthSymm3trK; CCTK_REAL_VEC JacPDupwindNthSymm3Xt1; CCTK_REAL_VEC JacPDupwindNthSymm3Xt2; @@ -1294,6 +1335,9 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthAnti1phi = kmadd(J11L,PDupwindNthAnti1phi,kmadd(J21L,PDupwindNthAnti2phi,kmul(J31L,PDupwindNthAnti3phi))); + JacPDupwindNthAnti1Theta = + kmadd(J11L,PDupwindNthAnti1Theta,kmadd(J21L,PDupwindNthAnti2Theta,kmul(J31L,PDupwindNthAnti3Theta))); + JacPDupwindNthAnti1trK = kmadd(J11L,PDupwindNthAnti1trK,kmadd(J21L,PDupwindNthAnti2trK,kmul(J31L,PDupwindNthAnti3trK))); @@ -1369,6 +1413,9 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthSymm1phi = kmadd(J11L,PDupwindNthSymm1phi,kmadd(J21L,PDupwindNthSymm2phi,kmul(J31L,PDupwindNthSymm3phi))); + JacPDupwindNthSymm1Theta = + kmadd(J11L,PDupwindNthSymm1Theta,kmadd(J21L,PDupwindNthSymm2Theta,kmul(J31L,PDupwindNthSymm3Theta))); + JacPDupwindNthSymm1trK = kmadd(J11L,PDupwindNthSymm1trK,kmadd(J21L,PDupwindNthSymm2trK,kmul(J31L,PDupwindNthSymm3trK))); @@ -1444,6 +1491,9 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthAnti2phi = kmadd(J12L,PDupwindNthAnti1phi,kmadd(J22L,PDupwindNthAnti2phi,kmul(J32L,PDupwindNthAnti3phi))); + JacPDupwindNthAnti2Theta = + kmadd(J12L,PDupwindNthAnti1Theta,kmadd(J22L,PDupwindNthAnti2Theta,kmul(J32L,PDupwindNthAnti3Theta))); + JacPDupwindNthAnti2trK = kmadd(J12L,PDupwindNthAnti1trK,kmadd(J22L,PDupwindNthAnti2trK,kmul(J32L,PDupwindNthAnti3trK))); @@ -1519,6 +1569,9 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthSymm2phi = kmadd(J12L,PDupwindNthSymm1phi,kmadd(J22L,PDupwindNthSymm2phi,kmul(J32L,PDupwindNthSymm3phi))); + JacPDupwindNthSymm2Theta = + kmadd(J12L,PDupwindNthSymm1Theta,kmadd(J22L,PDupwindNthSymm2Theta,kmul(J32L,PDupwindNthSymm3Theta))); + JacPDupwindNthSymm2trK = kmadd(J12L,PDupwindNthSymm1trK,kmadd(J22L,PDupwindNthSymm2trK,kmul(J32L,PDupwindNthSymm3trK))); @@ -1594,6 +1647,9 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthAnti3phi = kmadd(J13L,PDupwindNthAnti1phi,kmadd(J23L,PDupwindNthAnti2phi,kmul(J33L,PDupwindNthAnti3phi))); + JacPDupwindNthAnti3Theta = + kmadd(J13L,PDupwindNthAnti1Theta,kmadd(J23L,PDupwindNthAnti2Theta,kmul(J33L,PDupwindNthAnti3Theta))); + JacPDupwindNthAnti3trK = kmadd(J13L,PDupwindNthAnti1trK,kmadd(J23L,PDupwindNthAnti2trK,kmul(J33L,PDupwindNthAnti3trK))); @@ -1669,6 +1725,9 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthSymm3phi = kmadd(J13L,PDupwindNthSymm1phi,kmadd(J23L,PDupwindNthSymm2phi,kmul(J33L,PDupwindNthSymm3phi))); + JacPDupwindNthSymm3Theta = + kmadd(J13L,PDupwindNthSymm1Theta,kmadd(J23L,PDupwindNthSymm2Theta,kmul(J33L,PDupwindNthSymm3Theta))); + JacPDupwindNthSymm3trK = kmadd(J13L,PDupwindNthSymm1trK,kmadd(J23L,PDupwindNthSymm2trK,kmul(J33L,PDupwindNthSymm3trK))); @@ -1725,6 +1784,8 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthAnti1phi = PDupwindNthAnti1phi; + JacPDupwindNthAnti1Theta = PDupwindNthAnti1Theta; + JacPDupwindNthAnti1trK = PDupwindNthAnti1trK; JacPDupwindNthAnti1Xt1 = PDupwindNthAnti1Xt1; @@ -1775,6 +1836,8 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthSymm1phi = PDupwindNthSymm1phi; + JacPDupwindNthSymm1Theta = PDupwindNthSymm1Theta; + JacPDupwindNthSymm1trK = PDupwindNthSymm1trK; JacPDupwindNthSymm1Xt1 = PDupwindNthSymm1Xt1; @@ -1825,6 +1888,8 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthAnti2phi = PDupwindNthAnti2phi; + JacPDupwindNthAnti2Theta = PDupwindNthAnti2Theta; + JacPDupwindNthAnti2trK = PDupwindNthAnti2trK; JacPDupwindNthAnti2Xt1 = PDupwindNthAnti2Xt1; @@ -1875,6 +1940,8 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthSymm2phi = PDupwindNthSymm2phi; + JacPDupwindNthSymm2Theta = PDupwindNthSymm2Theta; + JacPDupwindNthSymm2trK = PDupwindNthSymm2trK; JacPDupwindNthSymm2Xt1 = PDupwindNthSymm2Xt1; @@ -1925,6 +1992,8 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthAnti3phi = PDupwindNthAnti3phi; + JacPDupwindNthAnti3Theta = PDupwindNthAnti3Theta; + JacPDupwindNthAnti3trK = PDupwindNthAnti3trK; JacPDupwindNthAnti3Xt1 = PDupwindNthAnti3Xt1; @@ -1975,6 +2044,8 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir JacPDupwindNthSymm3phi = PDupwindNthSymm3phi; + JacPDupwindNthSymm3Theta = PDupwindNthSymm3Theta; + JacPDupwindNthSymm3trK = PDupwindNthSymm3trK; JacPDupwindNthSymm3Xt1 = PDupwindNthSymm3Xt1; @@ -2014,6 +2085,9 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir Xt3rhsL = kadd(Xt3rhsL,kmadd(beta1L,JacPDupwindNthAnti1Xt3,kmadd(beta2L,JacPDupwindNthAnti2Xt3,kmadd(beta3L,JacPDupwindNthAnti3Xt3,kmadd(JacPDupwindNthSymm1Xt3,kfabs(beta1L),kmadd(JacPDupwindNthSymm2Xt3,kfabs(beta2L),kmul(JacPDupwindNthSymm3Xt3,kfabs(beta3L)))))))); + ThetarhsL = kadd(ThetarhsL,IfThen(conformalMethod == + 2,kmadd(beta1L,JacPDupwindNthAnti1Theta,kmadd(beta2L,JacPDupwindNthAnti2Theta,kmadd(beta3L,JacPDupwindNthAnti3Theta,kmadd(JacPDupwindNthSymm1Theta,kfabs(beta1L),kmadd(JacPDupwindNthSymm2Theta,kfabs(beta2L),kmul(JacPDupwindNthSymm3Theta,kfabs(beta3L))))))),ToReal(0))); + trKrhsL = kadd(trKrhsL,kmadd(beta1L,JacPDupwindNthAnti1trK,kmadd(beta2L,JacPDupwindNthAnti2trK,kmadd(beta3L,JacPDupwindNthAnti3trK,kmadd(JacPDupwindNthSymm1trK,kfabs(beta1L),kmadd(JacPDupwindNthSymm2trK,kfabs(beta2L),kmul(JacPDupwindNthSymm3trK,kfabs(beta3L)))))))); @@ -2082,6 +2156,7 @@ static void ML_BSSN_Advect_Body(cGH const * restrict const cctkGH, int const dir vec_store_nta_partial(gt23rhs[index],gt23rhsL); vec_store_nta_partial(gt33rhs[index],gt33rhsL); vec_store_nta_partial(phirhs[index],phirhsL); + vec_store_nta_partial(Thetarhs[index],ThetarhsL); vec_store_nta_partial(trKrhs[index],trKrhsL); vec_store_nta_partial(Xt1rhs[index],Xt1rhsL); vec_store_nta_partial(Xt2rhs[index],Xt2rhsL); @@ -2123,9 +2198,11 @@ extern "C" void ML_BSSN_Advect(CCTK_ARGUMENTS) "ML_BSSN::ML_metricrhs", "ML_BSSN::ML_shift", "ML_BSSN::ML_shiftrhs", + "ML_BSSN::ML_Theta", + "ML_BSSN::ML_Thetarhs", "ML_BSSN::ML_trace_curv", "ML_BSSN::ML_trace_curvrhs"}; - GenericFD_AssertGroupStorage(cctkGH, "ML_BSSN_Advect", 18, groups); + GenericFD_AssertGroupStorage(cctkGH, "ML_BSSN_Advect", 20, groups); switch(fdOrder) { |