diff options
Diffstat (limited to 'ML_BSSN/src/ML_BSSN_convertToADMBase.cc')
| -rw-r--r-- | ML_BSSN/src/ML_BSSN_convertToADMBase.cc | 150 |
1 files changed, 44 insertions, 106 deletions
diff --git a/ML_BSSN/src/ML_BSSN_convertToADMBase.cc b/ML_BSSN/src/ML_BSSN_convertToADMBase.cc index 6afda7b..de1d9df 100644 --- a/ML_BSSN/src/ML_BSSN_convertToADMBase.cc +++ b/ML_BSSN/src/ML_BSSN_convertToADMBase.cc @@ -29,8 +29,6 @@ static void ML_BSSN_convertToADMBase_Body(cGH const * restrict const cctkGH, int DECLARE_CCTK_PARAMETERS; - /* Declare finite differencing variables */ - /* Include user-supplied include files */ /* Initialise finite differencing variables */ @@ -67,9 +65,9 @@ static void ML_BSSN_convertToADMBase_Body(cGH const * restrict const cctkGH, int CCTK_REAL_VEC const p1o12dx = kmul(INV(dx),ToReal(0.0833333333333333333333333333333)); CCTK_REAL_VEC const p1o12dy = kmul(INV(dy),ToReal(0.0833333333333333333333333333333)); CCTK_REAL_VEC const p1o12dz = kmul(INV(dz),ToReal(0.0833333333333333333333333333333)); - CCTK_REAL_VEC const p1o144dxdy = kmul(INV(dx),kmul(INV(dy),ToReal(0.00694444444444444444444444444444))); - CCTK_REAL_VEC const p1o144dxdz = kmul(INV(dx),kmul(INV(dz),ToReal(0.00694444444444444444444444444444))); - CCTK_REAL_VEC const p1o144dydz = kmul(INV(dy),kmul(INV(dz),ToReal(0.00694444444444444444444444444444))); + CCTK_REAL_VEC const p1o144dxdy = kmul(INV(kmul(dx,dy)),ToReal(0.00694444444444444444444444444444)); + CCTK_REAL_VEC const p1o144dxdz = kmul(INV(kmul(dx,dz)),ToReal(0.00694444444444444444444444444444)); + CCTK_REAL_VEC const p1o144dydz = kmul(INV(kmul(dy,dz)),ToReal(0.00694444444444444444444444444444)); CCTK_REAL_VEC const p1o1680dx = kmul(INV(dx),ToReal(0.000595238095238095238095238095238)); CCTK_REAL_VEC const p1o1680dy = kmul(INV(dy),ToReal(0.000595238095238095238095238095238)); CCTK_REAL_VEC const p1o1680dz = kmul(INV(dz),ToReal(0.000595238095238095238095238095238)); @@ -88,14 +86,14 @@ static void ML_BSSN_convertToADMBase_Body(cGH const * restrict const cctkGH, int CCTK_REAL_VEC const p1o2dx = kmul(INV(dx),ToReal(0.5)); CCTK_REAL_VEC const p1o2dy = kmul(INV(dy),ToReal(0.5)); CCTK_REAL_VEC const p1o2dz = kmul(INV(dz),ToReal(0.5)); - CCTK_REAL_VEC const p1o3600dxdy = kmul(INV(dx),kmul(INV(dy),ToReal(0.000277777777777777777777777777778))); - CCTK_REAL_VEC const p1o3600dxdz = kmul(INV(dx),kmul(INV(dz),ToReal(0.000277777777777777777777777777778))); - CCTK_REAL_VEC const p1o3600dydz = kmul(INV(dy),kmul(INV(dz),ToReal(0.000277777777777777777777777777778))); + CCTK_REAL_VEC const p1o3600dxdy = kmul(INV(kmul(dx,dy)),ToReal(0.000277777777777777777777777777778)); + CCTK_REAL_VEC const p1o3600dxdz = kmul(INV(kmul(dx,dz)),ToReal(0.000277777777777777777777777777778)); + CCTK_REAL_VEC const p1o3600dydz = kmul(INV(kmul(dy,dz)),ToReal(0.000277777777777777777777777777778)); CCTK_REAL_VEC const p1o4dx = kmul(INV(dx),ToReal(0.25)); - CCTK_REAL_VEC const p1o4dxdy = kmul(INV(dx),kmul(INV(dy),ToReal(0.25))); - CCTK_REAL_VEC const p1o4dxdz = kmul(INV(dx),kmul(INV(dz),ToReal(0.25))); + CCTK_REAL_VEC const p1o4dxdy = kmul(INV(kmul(dx,dy)),ToReal(0.25)); + CCTK_REAL_VEC const p1o4dxdz = kmul(INV(kmul(dx,dz)),ToReal(0.25)); CCTK_REAL_VEC const p1o4dy = kmul(INV(dy),ToReal(0.25)); - CCTK_REAL_VEC const p1o4dydz = kmul(INV(dy),kmul(INV(dz),ToReal(0.25))); + CCTK_REAL_VEC const p1o4dydz = kmul(INV(kmul(dy,dz)),ToReal(0.25)); CCTK_REAL_VEC const p1o4dz = kmul(INV(dz),ToReal(0.25)); CCTK_REAL_VEC const p1o5040dx2 = kmul(INV(SQR(dx)),ToReal(0.000198412698412698412698412698413)); CCTK_REAL_VEC const p1o5040dy2 = kmul(INV(SQR(dy)),ToReal(0.000198412698412698412698412698413)); @@ -109,9 +107,9 @@ static void ML_BSSN_convertToADMBase_Body(cGH const * restrict const cctkGH, int CCTK_REAL_VEC const p1o64dx = kmul(INV(dx),ToReal(0.015625)); CCTK_REAL_VEC const p1o64dy = kmul(INV(dy),ToReal(0.015625)); CCTK_REAL_VEC const p1o64dz = kmul(INV(dz),ToReal(0.015625)); - CCTK_REAL_VEC const p1o705600dxdy = kmul(INV(dx),kmul(INV(dy),ToReal(1.41723356009070294784580498866e-6))); - CCTK_REAL_VEC const p1o705600dxdz = kmul(INV(dx),kmul(INV(dz),ToReal(1.41723356009070294784580498866e-6))); - CCTK_REAL_VEC const p1o705600dydz = kmul(INV(dy),kmul(INV(dz),ToReal(1.41723356009070294784580498866e-6))); + CCTK_REAL_VEC const p1o705600dxdy = kmul(INV(kmul(dx,dy)),ToReal(1.41723356009070294784580498866e-6)); + CCTK_REAL_VEC const p1o705600dxdz = kmul(INV(kmul(dx,dz)),ToReal(1.41723356009070294784580498866e-6)); + CCTK_REAL_VEC const p1o705600dydz = kmul(INV(kmul(dy,dz)),ToReal(1.41723356009070294784580498866e-6)); CCTK_REAL_VEC const p1o840dx = kmul(INV(dx),ToReal(0.00119047619047619047619047619048)); CCTK_REAL_VEC const p1o840dy = kmul(INV(dy),ToReal(0.00119047619047619047619047619048)); CCTK_REAL_VEC const p1o840dz = kmul(INV(dz),ToReal(0.00119047619047619047619047619048)); @@ -192,7 +190,7 @@ static void ML_BSSN_convertToADMBase_Body(cGH const * restrict const cctkGH, int /* Loop over the grid points */ #pragma omp parallel - LC_LOOP3VEC (ML_BSSN_convertToADMBase, + LC_LOOP3VEC(ML_BSSN_convertToADMBase, i,j,k, imin[0],imin[1],imin[2], imax[0],imax[1],imax[2], cctk_lsh[0],cctk_lsh[1],cctk_lsh[2], CCTK_REAL_VEC_SIZE) @@ -289,96 +287,26 @@ static void ML_BSSN_convertToADMBase_Body(cGH const * restrict const cctkGH, int CCTK_REAL_VEC betazL = beta3L; - /* If necessary, store only partial vectors after the first iteration */ - - if (CCTK_REAL_VEC_SIZE > 2 && CCTK_BUILTIN_EXPECT(i < lc_imin && i+CCTK_REAL_VEC_SIZE > lc_imax, 0)) - { - ptrdiff_t const elt_count_lo = lc_imin-i; - ptrdiff_t const elt_count_hi = lc_imax-i; - vec_store_nta_partial_mid(alp[index],alpL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(betax[index],betaxL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(betay[index],betayL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(betaz[index],betazL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(gxx[index],gxxL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(gxy[index],gxyL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(gxz[index],gxzL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(gyy[index],gyyL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(gyz[index],gyzL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(gzz[index],gzzL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(kxx[index],kxxL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(kxy[index],kxyL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(kxz[index],kxzL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(kyy[index],kyyL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(kyz[index],kyzL,elt_count_lo,elt_count_hi); - vec_store_nta_partial_mid(kzz[index],kzzL,elt_count_lo,elt_count_hi); - break; - } - - /* If necessary, store only partial vectors after the first iteration */ - - if (CCTK_REAL_VEC_SIZE > 1 && CCTK_BUILTIN_EXPECT(i < lc_imin, 0)) - { - ptrdiff_t const elt_count = lc_imin-i; - vec_store_nta_partial_hi(alp[index],alpL,elt_count); - vec_store_nta_partial_hi(betax[index],betaxL,elt_count); - vec_store_nta_partial_hi(betay[index],betayL,elt_count); - vec_store_nta_partial_hi(betaz[index],betazL,elt_count); - vec_store_nta_partial_hi(gxx[index],gxxL,elt_count); - vec_store_nta_partial_hi(gxy[index],gxyL,elt_count); - vec_store_nta_partial_hi(gxz[index],gxzL,elt_count); - vec_store_nta_partial_hi(gyy[index],gyyL,elt_count); - vec_store_nta_partial_hi(gyz[index],gyzL,elt_count); - vec_store_nta_partial_hi(gzz[index],gzzL,elt_count); - vec_store_nta_partial_hi(kxx[index],kxxL,elt_count); - vec_store_nta_partial_hi(kxy[index],kxyL,elt_count); - vec_store_nta_partial_hi(kxz[index],kxzL,elt_count); - vec_store_nta_partial_hi(kyy[index],kyyL,elt_count); - vec_store_nta_partial_hi(kyz[index],kyzL,elt_count); - vec_store_nta_partial_hi(kzz[index],kzzL,elt_count); - continue; - } - - /* If necessary, store only partial vectors after the last iteration */ - - if (CCTK_REAL_VEC_SIZE > 1 && CCTK_BUILTIN_EXPECT(i+CCTK_REAL_VEC_SIZE > lc_imax, 0)) - { - ptrdiff_t const elt_count = lc_imax-i; - vec_store_nta_partial_lo(alp[index],alpL,elt_count); - vec_store_nta_partial_lo(betax[index],betaxL,elt_count); - vec_store_nta_partial_lo(betay[index],betayL,elt_count); - vec_store_nta_partial_lo(betaz[index],betazL,elt_count); - vec_store_nta_partial_lo(gxx[index],gxxL,elt_count); - vec_store_nta_partial_lo(gxy[index],gxyL,elt_count); - vec_store_nta_partial_lo(gxz[index],gxzL,elt_count); - vec_store_nta_partial_lo(gyy[index],gyyL,elt_count); - vec_store_nta_partial_lo(gyz[index],gyzL,elt_count); - vec_store_nta_partial_lo(gzz[index],gzzL,elt_count); - vec_store_nta_partial_lo(kxx[index],kxxL,elt_count); - vec_store_nta_partial_lo(kxy[index],kxyL,elt_count); - vec_store_nta_partial_lo(kxz[index],kxzL,elt_count); - vec_store_nta_partial_lo(kyy[index],kyyL,elt_count); - vec_store_nta_partial_lo(kyz[index],kyzL,elt_count); - vec_store_nta_partial_lo(kzz[index],kzzL,elt_count); - break; - } - vec_store_nta(alp[index],alpL); - vec_store_nta(betax[index],betaxL); - vec_store_nta(betay[index],betayL); - vec_store_nta(betaz[index],betazL); - vec_store_nta(gxx[index],gxxL); - vec_store_nta(gxy[index],gxyL); - vec_store_nta(gxz[index],gxzL); - vec_store_nta(gyy[index],gyyL); - vec_store_nta(gyz[index],gyzL); - vec_store_nta(gzz[index],gzzL); - vec_store_nta(kxx[index],kxxL); - vec_store_nta(kxy[index],kxyL); - vec_store_nta(kxz[index],kxzL); - vec_store_nta(kyy[index],kyyL); - vec_store_nta(kyz[index],kyzL); - vec_store_nta(kzz[index],kzzL); + /* Copy local copies back to grid functions */ + vec_store_partial_prepare(i,lc_imin,lc_imax); + vec_store_nta_partial(alp[index],alpL); + vec_store_nta_partial(betax[index],betaxL); + vec_store_nta_partial(betay[index],betayL); + vec_store_nta_partial(betaz[index],betazL); + vec_store_nta_partial(gxx[index],gxxL); + vec_store_nta_partial(gxy[index],gxyL); + vec_store_nta_partial(gxz[index],gxzL); + vec_store_nta_partial(gyy[index],gyyL); + vec_store_nta_partial(gyz[index],gyzL); + vec_store_nta_partial(gzz[index],gzzL); + vec_store_nta_partial(kxx[index],kxxL); + vec_store_nta_partial(kxy[index],kxyL); + vec_store_nta_partial(kxz[index],kxzL); + vec_store_nta_partial(kyy[index],kyyL); + vec_store_nta_partial(kyz[index],kyzL); + vec_store_nta_partial(kzz[index],kzzL); } - LC_ENDLOOP3VEC (ML_BSSN_convertToADMBase); + LC_ENDLOOP3VEC(ML_BSSN_convertToADMBase); } extern "C" void ML_BSSN_convertToADMBase(CCTK_ARGUMENTS) @@ -397,7 +325,17 @@ extern "C" void ML_BSSN_convertToADMBase(CCTK_ARGUMENTS) return; } - const char *groups[] = {"ADMBase::curv","ADMBase::lapse","ADMBase::metric","ADMBase::shift","ML_BSSN::ML_curv","ML_BSSN::ML_lapse","ML_BSSN::ML_log_confac","ML_BSSN::ML_metric","ML_BSSN::ML_shift","ML_BSSN::ML_trace_curv"}; + const char *const groups[] = { + "ADMBase::curv", + "ADMBase::lapse", + "ADMBase::metric", + "ADMBase::shift", + "ML_BSSN::ML_curv", + "ML_BSSN::ML_lapse", + "ML_BSSN::ML_log_confac", + "ML_BSSN::ML_metric", + "ML_BSSN::ML_shift", + "ML_BSSN::ML_trace_curv"}; GenericFD_AssertGroupStorage(cctkGH, "ML_BSSN_convertToADMBase", 10, groups); switch(fdOrder) @@ -415,7 +353,7 @@ extern "C" void ML_BSSN_convertToADMBase(CCTK_ARGUMENTS) break; } - GenericFD_LoopOverEverything(cctkGH, &ML_BSSN_convertToADMBase_Body); + GenericFD_LoopOverEverything(cctkGH, ML_BSSN_convertToADMBase_Body); if (verbose > 1) { |