/* File produced by Kranc */ #define KRANC_C #include #include #include #include #include #include "cctk.h" #include "cctk_Arguments.h" #include "cctk_Parameters.h" #include "GenericFD.h" #include "Differencing.h" #include "cctk_Loop.h" #include "loopcontrol.h" #include "vectors.h" /* Define macros used in calculations */ #define INITVALUE (42) #define ScalarINV(x) ((CCTK_REAL)1.0 / (x)) #define ScalarSQR(x) ((x) * (x)) #define ScalarCUB(x) ((x) * ScalarSQR(x)) #define ScalarQAD(x) (ScalarSQR(ScalarSQR(x))) #define INV(x) (kdiv(ToReal(1.0),x)) #define SQR(x) (kmul(x,x)) #define CUB(x) (kmul(x,SQR(x))) #define QAD(x) (SQR(SQR(x))) extern "C" void psi4_calc_4th_SelectBCs(CCTK_ARGUMENTS) { DECLARE_CCTK_ARGUMENTS; DECLARE_CCTK_PARAMETERS; CCTK_INT ierr = 0; ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, GenericFD_GetBoundaryWidth(cctkGH), -1 /* no table */, "WeylScal4::Psi4i_group","flat"); if (ierr < 0) CCTK_WARN(1, "Failed to register flat BC for WeylScal4::Psi4i_group."); ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, GenericFD_GetBoundaryWidth(cctkGH), -1 /* no table */, "WeylScal4::Psi4r_group","flat"); if (ierr < 0) CCTK_WARN(1, "Failed to register flat BC for WeylScal4::Psi4r_group."); return; } static void psi4_calc_4th_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 imin[3], int const imax[3], int const n_subblock_gfs, CCTK_REAL * restrict const subblock_gfs[]) { DECLARE_CCTK_ARGUMENTS; DECLARE_CCTK_PARAMETERS; /* Include user-supplied include files */ /* Initialise finite differencing variables */ ptrdiff_t const di = 1; ptrdiff_t const dj = CCTK_GFINDEX3D(cctkGH,0,1,0) - CCTK_GFINDEX3D(cctkGH,0,0,0); ptrdiff_t const dk = CCTK_GFINDEX3D(cctkGH,0,0,1) - CCTK_GFINDEX3D(cctkGH,0,0,0); ptrdiff_t const cdi = sizeof(CCTK_REAL) * di; ptrdiff_t const cdj = sizeof(CCTK_REAL) * dj; ptrdiff_t const cdk = sizeof(CCTK_REAL) * dk; CCTK_REAL_VEC const dx = ToReal(CCTK_DELTA_SPACE(0)); CCTK_REAL_VEC const dy = ToReal(CCTK_DELTA_SPACE(1)); CCTK_REAL_VEC const dz = ToReal(CCTK_DELTA_SPACE(2)); CCTK_REAL_VEC const dt = ToReal(CCTK_DELTA_TIME); CCTK_REAL_VEC const t = ToReal(cctk_time); CCTK_REAL_VEC const dxi = INV(dx); CCTK_REAL_VEC const dyi = INV(dy); CCTK_REAL_VEC const dzi = INV(dz); CCTK_REAL_VEC const khalf = ToReal(0.5); CCTK_REAL_VEC const kthird = ToReal(1.0/3.0); CCTK_REAL_VEC const ktwothird = ToReal(2.0/3.0); CCTK_REAL_VEC const kfourthird = ToReal(4.0/3.0); CCTK_REAL_VEC const keightthird = ToReal(8.0/3.0); CCTK_REAL_VEC const hdxi = kmul(ToReal(0.5), dxi); CCTK_REAL_VEC const hdyi = kmul(ToReal(0.5), dyi); CCTK_REAL_VEC const hdzi = kmul(ToReal(0.5), dzi); /* Initialize predefined quantities */ CCTK_REAL_VEC const p1o12dx = kdiv(ToReal(0.0833333333333333333333333333333),dx); CCTK_REAL_VEC const p1o12dy = kdiv(ToReal(0.0833333333333333333333333333333),dy); CCTK_REAL_VEC const p1o12dz = kdiv(ToReal(0.0833333333333333333333333333333),dz); CCTK_REAL_VEC const p1o144dxdy = kdiv(ToReal(0.00694444444444444444444444444444),kmul(dy,dx)); CCTK_REAL_VEC const p1o144dxdz = kdiv(ToReal(0.00694444444444444444444444444444),kmul(dz,dx)); CCTK_REAL_VEC const p1o144dydz = kdiv(ToReal(0.00694444444444444444444444444444),kmul(dz,dy)); CCTK_REAL_VEC const p1o180dx2 = kdiv(ToReal(0.00555555555555555555555555555556),kmul(dx,dx)); CCTK_REAL_VEC const p1o180dy2 = kdiv(ToReal(0.00555555555555555555555555555556),kmul(dy,dy)); CCTK_REAL_VEC const p1o180dz2 = kdiv(ToReal(0.00555555555555555555555555555556),kmul(dz,dz)); CCTK_REAL_VEC const p1o2dx = kdiv(ToReal(0.5),dx); CCTK_REAL_VEC const p1o2dy = kdiv(ToReal(0.5),dy); CCTK_REAL_VEC const p1o2dz = kdiv(ToReal(0.5),dz); CCTK_REAL_VEC const p1o3600dxdy = kdiv(ToReal(0.000277777777777777777777777777778),kmul(dy,dx)); CCTK_REAL_VEC const p1o3600dxdz = kdiv(ToReal(0.000277777777777777777777777777778),kmul(dz,dx)); CCTK_REAL_VEC const p1o3600dydz = kdiv(ToReal(0.000277777777777777777777777777778),kmul(dz,dy)); CCTK_REAL_VEC const p1o4dxdy = kdiv(ToReal(0.25),kmul(dy,dx)); CCTK_REAL_VEC const p1o4dxdz = kdiv(ToReal(0.25),kmul(dz,dx)); CCTK_REAL_VEC const p1o4dydz = kdiv(ToReal(0.25),kmul(dz,dy)); CCTK_REAL_VEC const p1o5040dx2 = kdiv(ToReal(0.000198412698412698412698412698413),kmul(dx,dx)); CCTK_REAL_VEC const p1o5040dy2 = kdiv(ToReal(0.000198412698412698412698412698413),kmul(dy,dy)); CCTK_REAL_VEC const p1o5040dz2 = kdiv(ToReal(0.000198412698412698412698412698413),kmul(dz,dz)); CCTK_REAL_VEC const p1o60dx = kdiv(ToReal(0.0166666666666666666666666666667),dx); CCTK_REAL_VEC const p1o60dy = kdiv(ToReal(0.0166666666666666666666666666667),dy); CCTK_REAL_VEC const p1o60dz = kdiv(ToReal(0.0166666666666666666666666666667),dz); CCTK_REAL_VEC const p1o705600dxdy = kdiv(ToReal(1.41723356009070294784580498866e-6),kmul(dy,dx)); CCTK_REAL_VEC const p1o705600dxdz = kdiv(ToReal(1.41723356009070294784580498866e-6),kmul(dz,dx)); CCTK_REAL_VEC const p1o705600dydz = kdiv(ToReal(1.41723356009070294784580498866e-6),kmul(dz,dy)); CCTK_REAL_VEC const p1o840dx = kdiv(ToReal(0.00119047619047619047619047619048),dx); CCTK_REAL_VEC const p1o840dy = kdiv(ToReal(0.00119047619047619047619047619048),dy); CCTK_REAL_VEC const p1o840dz = kdiv(ToReal(0.00119047619047619047619047619048),dz); CCTK_REAL_VEC const p1odx2 = kdiv(ToReal(1.),kmul(dx,dx)); CCTK_REAL_VEC const p1ody2 = kdiv(ToReal(1.),kmul(dy,dy)); CCTK_REAL_VEC const p1odz2 = kdiv(ToReal(1.),kmul(dz,dz)); CCTK_REAL_VEC const pm1o12dx2 = kdiv(ToReal(-0.0833333333333333333333333333333),kmul(dx,dx)); CCTK_REAL_VEC const pm1o12dy2 = kdiv(ToReal(-0.0833333333333333333333333333333),kmul(dy,dy)); CCTK_REAL_VEC const pm1o12dz2 = kdiv(ToReal(-0.0833333333333333333333333333333),kmul(dz,dz)); /* Jacobian variable pointers */ bool const use_jacobian = (!CCTK_IsFunctionAliased("MultiPatch_GetMap") || MultiPatch_GetMap(cctkGH) != jacobian_identity_map) && strlen(jacobian_group) > 0; if (use_jacobian && strlen(jacobian_derivative_group) == 0) { CCTK_WARN (1, "GenericFD::jacobian_group and GenericFD::jacobian_derivative_group must both be set to valid group names"); } CCTK_REAL const *restrict jacobian_ptrs[9]; if (use_jacobian) GenericFD_GroupDataPointers(cctkGH, jacobian_group, 9, jacobian_ptrs); CCTK_REAL const *restrict const J11 = use_jacobian ? jacobian_ptrs[0] : 0; CCTK_REAL const *restrict const J12 = use_jacobian ? jacobian_ptrs[1] : 0; CCTK_REAL const *restrict const J13 = use_jacobian ? jacobian_ptrs[2] : 0; CCTK_REAL const *restrict const J21 = use_jacobian ? jacobian_ptrs[3] : 0; CCTK_REAL const *restrict const J22 = use_jacobian ? jacobian_ptrs[4] : 0; CCTK_REAL const *restrict const J23 = use_jacobian ? jacobian_ptrs[5] : 0; CCTK_REAL const *restrict const J31 = use_jacobian ? jacobian_ptrs[6] : 0; CCTK_REAL const *restrict const J32 = use_jacobian ? jacobian_ptrs[7] : 0; CCTK_REAL const *restrict const J33 = use_jacobian ? jacobian_ptrs[8] : 0; CCTK_REAL const *restrict jacobian_derivative_ptrs[18]; if (use_jacobian) GenericFD_GroupDataPointers(cctkGH, jacobian_derivative_group, 18, jacobian_derivative_ptrs); CCTK_REAL const *restrict const dJ111 = use_jacobian ? jacobian_derivative_ptrs[0] : 0; CCTK_REAL const *restrict const dJ112 = use_jacobian ? jacobian_derivative_ptrs[1] : 0; CCTK_REAL const *restrict const dJ113 = use_jacobian ? jacobian_derivative_ptrs[2] : 0; CCTK_REAL const *restrict const dJ122 = use_jacobian ? jacobian_derivative_ptrs[3] : 0; CCTK_REAL const *restrict const dJ123 = use_jacobian ? jacobian_derivative_ptrs[4] : 0; CCTK_REAL const *restrict const dJ133 = use_jacobian ? jacobian_derivative_ptrs[5] : 0; CCTK_REAL const *restrict const dJ211 = use_jacobian ? jacobian_derivative_ptrs[6] : 0; CCTK_REAL const *restrict const dJ212 = use_jacobian ? jacobian_derivative_ptrs[7] : 0; CCTK_REAL const *restrict const dJ213 = use_jacobian ? jacobian_derivative_ptrs[8] : 0; CCTK_REAL const *restrict const dJ222 = use_jacobian ? jacobian_derivative_ptrs[9] : 0; CCTK_REAL const *restrict const dJ223 = use_jacobian ? jacobian_derivative_ptrs[10] : 0; CCTK_REAL const *restrict const dJ233 = use_jacobian ? jacobian_derivative_ptrs[11] : 0; CCTK_REAL const *restrict const dJ311 = use_jacobian ? jacobian_derivative_ptrs[12] : 0; CCTK_REAL const *restrict const dJ312 = use_jacobian ? jacobian_derivative_ptrs[13] : 0; CCTK_REAL const *restrict const dJ313 = use_jacobian ? jacobian_derivative_ptrs[14] : 0; CCTK_REAL const *restrict const dJ322 = use_jacobian ? jacobian_derivative_ptrs[15] : 0; CCTK_REAL const *restrict const dJ323 = use_jacobian ? jacobian_derivative_ptrs[16] : 0; CCTK_REAL const *restrict const dJ333 = use_jacobian ? jacobian_derivative_ptrs[17] : 0; /* Assign local copies of arrays functions */ /* Calculate temporaries and arrays functions */ /* Copy local copies back to grid functions */ /* Loop over the grid points */ #pragma omp parallel LC_LOOP3VEC(psi4_calc_4th, i,j,k, imin[0],imin[1],imin[2], imax[0],imax[1],imax[2], cctk_ash[0],cctk_ash[1],cctk_ash[2], CCTK_REAL_VEC_SIZE) { ptrdiff_t const index = di*i + dj*j + dk*k; /* Assign local copies of grid functions */ CCTK_REAL_VEC gxxL = vec_load(gxx[index]); CCTK_REAL_VEC gxyL = vec_load(gxy[index]); CCTK_REAL_VEC gxzL = vec_load(gxz[index]); CCTK_REAL_VEC gyyL = vec_load(gyy[index]); CCTK_REAL_VEC gyzL = vec_load(gyz[index]); CCTK_REAL_VEC gzzL = vec_load(gzz[index]); CCTK_REAL_VEC kxxL = vec_load(kxx[index]); CCTK_REAL_VEC kxyL = vec_load(kxy[index]); CCTK_REAL_VEC kxzL = vec_load(kxz[index]); CCTK_REAL_VEC kyyL = vec_load(kyy[index]); CCTK_REAL_VEC kyzL = vec_load(kyz[index]); CCTK_REAL_VEC kzzL = vec_load(kzz[index]); CCTK_REAL_VEC xL = vec_load(x[index]); CCTK_REAL_VEC yL = vec_load(y[index]); CCTK_REAL_VEC zL = vec_load(z[index]); CCTK_REAL_VEC dJ111L, dJ112L, dJ113L, dJ122L, dJ123L, dJ133L, dJ211L, dJ212L, dJ213L, dJ222L, dJ223L, dJ233L, dJ311L, dJ312L, dJ313L, dJ322L, dJ323L, dJ333L, J11L, J12L, J13L, J21L, J22L, J23L, J31L, J32L, J33L; if (use_jacobian) { dJ111L = vec_load(dJ111[index]); dJ112L = vec_load(dJ112[index]); dJ113L = vec_load(dJ113[index]); dJ122L = vec_load(dJ122[index]); dJ123L = vec_load(dJ123[index]); dJ133L = vec_load(dJ133[index]); dJ211L = vec_load(dJ211[index]); dJ212L = vec_load(dJ212[index]); dJ213L = vec_load(dJ213[index]); dJ222L = vec_load(dJ222[index]); dJ223L = vec_load(dJ223[index]); dJ233L = vec_load(dJ233[index]); dJ311L = vec_load(dJ311[index]); dJ312L = vec_load(dJ312[index]); dJ313L = vec_load(dJ313[index]); dJ322L = vec_load(dJ322[index]); dJ323L = vec_load(dJ323[index]); dJ333L = vec_load(dJ333[index]); J11L = vec_load(J11[index]); J12L = vec_load(J12[index]); J13L = vec_load(J13[index]); J21L = vec_load(J21[index]); J22L = vec_load(J22[index]); J23L = vec_load(J23[index]); J31L = vec_load(J31[index]); J32L = vec_load(J32[index]); J33L = vec_load(J33[index]); } /* Include user supplied include files */ /* Precompute derivatives */ CCTK_REAL_VEC PDstandard4th1gxx; CCTK_REAL_VEC PDstandard4th2gxx; CCTK_REAL_VEC PDstandard4th3gxx; CCTK_REAL_VEC PDstandard4th11gxx; CCTK_REAL_VEC PDstandard4th22gxx; CCTK_REAL_VEC PDstandard4th33gxx; CCTK_REAL_VEC PDstandard4th12gxx; CCTK_REAL_VEC PDstandard4th13gxx; CCTK_REAL_VEC PDstandard4th23gxx; CCTK_REAL_VEC PDstandard4th1gxy; CCTK_REAL_VEC PDstandard4th2gxy; CCTK_REAL_VEC PDstandard4th3gxy; CCTK_REAL_VEC PDstandard4th11gxy; CCTK_REAL_VEC PDstandard4th22gxy; CCTK_REAL_VEC PDstandard4th33gxy; CCTK_REAL_VEC PDstandard4th12gxy; CCTK_REAL_VEC PDstandard4th13gxy; CCTK_REAL_VEC PDstandard4th23gxy; CCTK_REAL_VEC PDstandard4th1gxz; CCTK_REAL_VEC PDstandard4th2gxz; CCTK_REAL_VEC PDstandard4th3gxz; CCTK_REAL_VEC PDstandard4th11gxz; CCTK_REAL_VEC PDstandard4th22gxz; CCTK_REAL_VEC PDstandard4th33gxz; CCTK_REAL_VEC PDstandard4th12gxz; CCTK_REAL_VEC PDstandard4th13gxz; CCTK_REAL_VEC PDstandard4th23gxz; CCTK_REAL_VEC PDstandard4th1gyy; CCTK_REAL_VEC PDstandard4th2gyy; CCTK_REAL_VEC PDstandard4th3gyy; CCTK_REAL_VEC PDstandard4th11gyy; CCTK_REAL_VEC PDstandard4th22gyy; CCTK_REAL_VEC PDstandard4th33gyy; CCTK_REAL_VEC PDstandard4th12gyy; CCTK_REAL_VEC PDstandard4th13gyy; CCTK_REAL_VEC PDstandard4th23gyy; CCTK_REAL_VEC PDstandard4th1gyz; CCTK_REAL_VEC PDstandard4th2gyz; CCTK_REAL_VEC PDstandard4th3gyz; CCTK_REAL_VEC PDstandard4th11gyz; CCTK_REAL_VEC PDstandard4th22gyz; CCTK_REAL_VEC PDstandard4th33gyz; CCTK_REAL_VEC PDstandard4th12gyz; CCTK_REAL_VEC PDstandard4th13gyz; CCTK_REAL_VEC PDstandard4th23gyz; CCTK_REAL_VEC PDstandard4th1gzz; CCTK_REAL_VEC PDstandard4th2gzz; CCTK_REAL_VEC PDstandard4th3gzz; CCTK_REAL_VEC PDstandard4th11gzz; CCTK_REAL_VEC PDstandard4th22gzz; CCTK_REAL_VEC PDstandard4th33gzz; CCTK_REAL_VEC PDstandard4th12gzz; CCTK_REAL_VEC PDstandard4th13gzz; CCTK_REAL_VEC PDstandard4th23gzz; CCTK_REAL_VEC PDstandard4th1kxx; CCTK_REAL_VEC PDstandard4th2kxx; CCTK_REAL_VEC PDstandard4th3kxx; CCTK_REAL_VEC PDstandard4th1kxy; CCTK_REAL_VEC PDstandard4th2kxy; CCTK_REAL_VEC PDstandard4th3kxy; CCTK_REAL_VEC PDstandard4th1kxz; CCTK_REAL_VEC PDstandard4th2kxz; CCTK_REAL_VEC PDstandard4th3kxz; CCTK_REAL_VEC PDstandard4th1kyy; CCTK_REAL_VEC PDstandard4th2kyy; CCTK_REAL_VEC PDstandard4th3kyy; CCTK_REAL_VEC PDstandard4th1kyz; CCTK_REAL_VEC PDstandard4th2kyz; CCTK_REAL_VEC PDstandard4th3kyz; CCTK_REAL_VEC PDstandard4th1kzz; CCTK_REAL_VEC PDstandard4th2kzz; CCTK_REAL_VEC PDstandard4th3kzz; switch(fdOrder) { case 2: PDstandard4th1gxx = PDstandard4th1(&gxx[index]); PDstandard4th2gxx = PDstandard4th2(&gxx[index]); PDstandard4th3gxx = PDstandard4th3(&gxx[index]); PDstandard4th11gxx = PDstandard4th11(&gxx[index]); PDstandard4th22gxx = PDstandard4th22(&gxx[index]); PDstandard4th33gxx = PDstandard4th33(&gxx[index]); PDstandard4th12gxx = PDstandard4th12(&gxx[index]); PDstandard4th13gxx = PDstandard4th13(&gxx[index]); PDstandard4th23gxx = PDstandard4th23(&gxx[index]); PDstandard4th1gxy = PDstandard4th1(&gxy[index]); PDstandard4th2gxy = PDstandard4th2(&gxy[index]); PDstandard4th3gxy = PDstandard4th3(&gxy[index]); PDstandard4th11gxy = PDstandard4th11(&gxy[index]); PDstandard4th22gxy = PDstandard4th22(&gxy[index]); PDstandard4th33gxy = PDstandard4th33(&gxy[index]); PDstandard4th12gxy = PDstandard4th12(&gxy[index]); PDstandard4th13gxy = PDstandard4th13(&gxy[index]); PDstandard4th23gxy = PDstandard4th23(&gxy[index]); PDstandard4th1gxz = PDstandard4th1(&gxz[index]); PDstandard4th2gxz = PDstandard4th2(&gxz[index]); PDstandard4th3gxz = PDstandard4th3(&gxz[index]); PDstandard4th11gxz = PDstandard4th11(&gxz[index]); PDstandard4th22gxz = PDstandard4th22(&gxz[index]); PDstandard4th33gxz = PDstandard4th33(&gxz[index]); PDstandard4th12gxz = PDstandard4th12(&gxz[index]); PDstandard4th13gxz = PDstandard4th13(&gxz[index]); PDstandard4th23gxz = PDstandard4th23(&gxz[index]); PDstandard4th1gyy = PDstandard4th1(&gyy[index]); PDstandard4th2gyy = PDstandard4th2(&gyy[index]); PDstandard4th3gyy = PDstandard4th3(&gyy[index]); PDstandard4th11gyy = PDstandard4th11(&gyy[index]); PDstandard4th22gyy = PDstandard4th22(&gyy[index]); PDstandard4th33gyy = PDstandard4th33(&gyy[index]); PDstandard4th12gyy = PDstandard4th12(&gyy[index]); PDstandard4th13gyy = PDstandard4th13(&gyy[index]); PDstandard4th23gyy = PDstandard4th23(&gyy[index]); PDstandard4th1gyz = PDstandard4th1(&gyz[index]); PDstandard4th2gyz = PDstandard4th2(&gyz[index]); PDstandard4th3gyz = PDstandard4th3(&gyz[index]); PDstandard4th11gyz = PDstandard4th11(&gyz[index]); PDstandard4th22gyz = PDstandard4th22(&gyz[index]); PDstandard4th33gyz = PDstandard4th33(&gyz[index]); PDstandard4th12gyz = PDstandard4th12(&gyz[index]); PDstandard4th13gyz = PDstandard4th13(&gyz[index]); PDstandard4th23gyz = PDstandard4th23(&gyz[index]); PDstandard4th1gzz = PDstandard4th1(&gzz[index]); PDstandard4th2gzz = PDstandard4th2(&gzz[index]); PDstandard4th3gzz = PDstandard4th3(&gzz[index]); PDstandard4th11gzz = PDstandard4th11(&gzz[index]); PDstandard4th22gzz = PDstandard4th22(&gzz[index]); PDstandard4th33gzz = PDstandard4th33(&gzz[index]); PDstandard4th12gzz = PDstandard4th12(&gzz[index]); PDstandard4th13gzz = PDstandard4th13(&gzz[index]); PDstandard4th23gzz = PDstandard4th23(&gzz[index]); PDstandard4th1kxx = PDstandard4th1(&kxx[index]); PDstandard4th2kxx = PDstandard4th2(&kxx[index]); PDstandard4th3kxx = PDstandard4th3(&kxx[index]); PDstandard4th1kxy = PDstandard4th1(&kxy[index]); PDstandard4th2kxy = PDstandard4th2(&kxy[index]); PDstandard4th3kxy = PDstandard4th3(&kxy[index]); PDstandard4th1kxz = PDstandard4th1(&kxz[index]); PDstandard4th2kxz = PDstandard4th2(&kxz[index]); PDstandard4th3kxz = PDstandard4th3(&kxz[index]); PDstandard4th1kyy = PDstandard4th1(&kyy[index]); PDstandard4th2kyy = PDstandard4th2(&kyy[index]); PDstandard4th3kyy = PDstandard4th3(&kyy[index]); PDstandard4th1kyz = PDstandard4th1(&kyz[index]); PDstandard4th2kyz = PDstandard4th2(&kyz[index]); PDstandard4th3kyz = PDstandard4th3(&kyz[index]); PDstandard4th1kzz = PDstandard4th1(&kzz[index]); PDstandard4th2kzz = PDstandard4th2(&kzz[index]); PDstandard4th3kzz = PDstandard4th3(&kzz[index]); break; case 4: PDstandard4th1gxx = PDstandard4th1(&gxx[index]); PDstandard4th2gxx = PDstandard4th2(&gxx[index]); PDstandard4th3gxx = PDstandard4th3(&gxx[index]); PDstandard4th11gxx = PDstandard4th11(&gxx[index]); PDstandard4th22gxx = PDstandard4th22(&gxx[index]); PDstandard4th33gxx = PDstandard4th33(&gxx[index]); PDstandard4th12gxx = PDstandard4th12(&gxx[index]); PDstandard4th13gxx = PDstandard4th13(&gxx[index]); PDstandard4th23gxx = PDstandard4th23(&gxx[index]); PDstandard4th1gxy = PDstandard4th1(&gxy[index]); PDstandard4th2gxy = PDstandard4th2(&gxy[index]); PDstandard4th3gxy = PDstandard4th3(&gxy[index]); PDstandard4th11gxy = PDstandard4th11(&gxy[index]); PDstandard4th22gxy = PDstandard4th22(&gxy[index]); PDstandard4th33gxy = PDstandard4th33(&gxy[index]); PDstandard4th12gxy = PDstandard4th12(&gxy[index]); PDstandard4th13gxy = PDstandard4th13(&gxy[index]); PDstandard4th23gxy = PDstandard4th23(&gxy[index]); PDstandard4th1gxz = PDstandard4th1(&gxz[index]); PDstandard4th2gxz = PDstandard4th2(&gxz[index]); PDstandard4th3gxz = PDstandard4th3(&gxz[index]); PDstandard4th11gxz = PDstandard4th11(&gxz[index]); PDstandard4th22gxz = PDstandard4th22(&gxz[index]); PDstandard4th33gxz = PDstandard4th33(&gxz[index]); PDstandard4th12gxz = PDstandard4th12(&gxz[index]); PDstandard4th13gxz = PDstandard4th13(&gxz[index]); PDstandard4th23gxz = PDstandard4th23(&gxz[index]); PDstandard4th1gyy = PDstandard4th1(&gyy[index]); PDstandard4th2gyy = PDstandard4th2(&gyy[index]); PDstandard4th3gyy = PDstandard4th3(&gyy[index]); PDstandard4th11gyy = PDstandard4th11(&gyy[index]); PDstandard4th22gyy = PDstandard4th22(&gyy[index]); PDstandard4th33gyy = PDstandard4th33(&gyy[index]); PDstandard4th12gyy = PDstandard4th12(&gyy[index]); PDstandard4th13gyy = PDstandard4th13(&gyy[index]); PDstandard4th23gyy = PDstandard4th23(&gyy[index]); PDstandard4th1gyz = PDstandard4th1(&gyz[index]); PDstandard4th2gyz = PDstandard4th2(&gyz[index]); PDstandard4th3gyz = PDstandard4th3(&gyz[index]); PDstandard4th11gyz = PDstandard4th11(&gyz[index]); PDstandard4th22gyz = PDstandard4th22(&gyz[index]); PDstandard4th33gyz = PDstandard4th33(&gyz[index]); PDstandard4th12gyz = PDstandard4th12(&gyz[index]); PDstandard4th13gyz = PDstandard4th13(&gyz[index]); PDstandard4th23gyz = PDstandard4th23(&gyz[index]); PDstandard4th1gzz = PDstandard4th1(&gzz[index]); PDstandard4th2gzz = PDstandard4th2(&gzz[index]); PDstandard4th3gzz = PDstandard4th3(&gzz[index]); PDstandard4th11gzz = PDstandard4th11(&gzz[index]); PDstandard4th22gzz = PDstandard4th22(&gzz[index]); PDstandard4th33gzz = PDstandard4th33(&gzz[index]); PDstandard4th12gzz = PDstandard4th12(&gzz[index]); PDstandard4th13gzz = PDstandard4th13(&gzz[index]); PDstandard4th23gzz = PDstandard4th23(&gzz[index]); PDstandard4th1kxx = PDstandard4th1(&kxx[index]); PDstandard4th2kxx = PDstandard4th2(&kxx[index]); PDstandard4th3kxx = PDstandard4th3(&kxx[index]); PDstandard4th1kxy = PDstandard4th1(&kxy[index]); PDstandard4th2kxy = PDstandard4th2(&kxy[index]); PDstandard4th3kxy = PDstandard4th3(&kxy[index]); PDstandard4th1kxz = PDstandard4th1(&kxz[index]); PDstandard4th2kxz = PDstandard4th2(&kxz[index]); PDstandard4th3kxz = PDstandard4th3(&kxz[index]); PDstandard4th1kyy = PDstandard4th1(&kyy[index]); PDstandard4th2kyy = PDstandard4th2(&kyy[index]); PDstandard4th3kyy = PDstandard4th3(&kyy[index]); PDstandard4th1kyz = PDstandard4th1(&kyz[index]); PDstandard4th2kyz = PDstandard4th2(&kyz[index]); PDstandard4th3kyz = PDstandard4th3(&kyz[index]); PDstandard4th1kzz = PDstandard4th1(&kzz[index]); PDstandard4th2kzz = PDstandard4th2(&kzz[index]); PDstandard4th3kzz = PDstandard4th3(&kzz[index]); break; case 6: PDstandard4th1gxx = PDstandard4th1(&gxx[index]); PDstandard4th2gxx = PDstandard4th2(&gxx[index]); PDstandard4th3gxx = PDstandard4th3(&gxx[index]); PDstandard4th11gxx = PDstandard4th11(&gxx[index]); PDstandard4th22gxx = PDstandard4th22(&gxx[index]); PDstandard4th33gxx = PDstandard4th33(&gxx[index]); PDstandard4th12gxx = PDstandard4th12(&gxx[index]); PDstandard4th13gxx = PDstandard4th13(&gxx[index]); PDstandard4th23gxx = PDstandard4th23(&gxx[index]); PDstandard4th1gxy = PDstandard4th1(&gxy[index]); PDstandard4th2gxy = PDstandard4th2(&gxy[index]); PDstandard4th3gxy = PDstandard4th3(&gxy[index]); PDstandard4th11gxy = PDstandard4th11(&gxy[index]); PDstandard4th22gxy = PDstandard4th22(&gxy[index]); PDstandard4th33gxy = PDstandard4th33(&gxy[index]); PDstandard4th12gxy = PDstandard4th12(&gxy[index]); PDstandard4th13gxy = PDstandard4th13(&gxy[index]); PDstandard4th23gxy = PDstandard4th23(&gxy[index]); PDstandard4th1gxz = PDstandard4th1(&gxz[index]); PDstandard4th2gxz = PDstandard4th2(&gxz[index]); PDstandard4th3gxz = PDstandard4th3(&gxz[index]); PDstandard4th11gxz = PDstandard4th11(&gxz[index]); PDstandard4th22gxz = PDstandard4th22(&gxz[index]); PDstandard4th33gxz = PDstandard4th33(&gxz[index]); PDstandard4th12gxz = PDstandard4th12(&gxz[index]); PDstandard4th13gxz = PDstandard4th13(&gxz[index]); PDstandard4th23gxz = PDstandard4th23(&gxz[index]); PDstandard4th1gyy = PDstandard4th1(&gyy[index]); PDstandard4th2gyy = PDstandard4th2(&gyy[index]); PDstandard4th3gyy = PDstandard4th3(&gyy[index]); PDstandard4th11gyy = PDstandard4th11(&gyy[index]); PDstandard4th22gyy = PDstandard4th22(&gyy[index]); PDstandard4th33gyy = PDstandard4th33(&gyy[index]); PDstandard4th12gyy = PDstandard4th12(&gyy[index]); PDstandard4th13gyy = PDstandard4th13(&gyy[index]); PDstandard4th23gyy = PDstandard4th23(&gyy[index]); PDstandard4th1gyz = PDstandard4th1(&gyz[index]); PDstandard4th2gyz = PDstandard4th2(&gyz[index]); PDstandard4th3gyz = PDstandard4th3(&gyz[index]); PDstandard4th11gyz = PDstandard4th11(&gyz[index]); PDstandard4th22gyz = PDstandard4th22(&gyz[index]); PDstandard4th33gyz = PDstandard4th33(&gyz[index]); PDstandard4th12gyz = PDstandard4th12(&gyz[index]); PDstandard4th13gyz = PDstandard4th13(&gyz[index]); PDstandard4th23gyz = PDstandard4th23(&gyz[index]); PDstandard4th1gzz = PDstandard4th1(&gzz[index]); PDstandard4th2gzz = PDstandard4th2(&gzz[index]); PDstandard4th3gzz = PDstandard4th3(&gzz[index]); PDstandard4th11gzz = PDstandard4th11(&gzz[index]); PDstandard4th22gzz = PDstandard4th22(&gzz[index]); PDstandard4th33gzz = PDstandard4th33(&gzz[index]); PDstandard4th12gzz = PDstandard4th12(&gzz[index]); PDstandard4th13gzz = PDstandard4th13(&gzz[index]); PDstandard4th23gzz = PDstandard4th23(&gzz[index]); PDstandard4th1kxx = PDstandard4th1(&kxx[index]); PDstandard4th2kxx = PDstandard4th2(&kxx[index]); PDstandard4th3kxx = PDstandard4th3(&kxx[index]); PDstandard4th1kxy = PDstandard4th1(&kxy[index]); PDstandard4th2kxy = PDstandard4th2(&kxy[index]); PDstandard4th3kxy = PDstandard4th3(&kxy[index]); PDstandard4th1kxz = PDstandard4th1(&kxz[index]); PDstandard4th2kxz = PDstandard4th2(&kxz[index]); PDstandard4th3kxz = PDstandard4th3(&kxz[index]); PDstandard4th1kyy = PDstandard4th1(&kyy[index]); PDstandard4th2kyy = PDstandard4th2(&kyy[index]); PDstandard4th3kyy = PDstandard4th3(&kyy[index]); PDstandard4th1kyz = PDstandard4th1(&kyz[index]); PDstandard4th2kyz = PDstandard4th2(&kyz[index]); PDstandard4th3kyz = PDstandard4th3(&kyz[index]); PDstandard4th1kzz = PDstandard4th1(&kzz[index]); PDstandard4th2kzz = PDstandard4th2(&kzz[index]); PDstandard4th3kzz = PDstandard4th3(&kzz[index]); break; case 8: PDstandard4th1gxx = PDstandard4th1(&gxx[index]); PDstandard4th2gxx = PDstandard4th2(&gxx[index]); PDstandard4th3gxx = PDstandard4th3(&gxx[index]); PDstandard4th11gxx = PDstandard4th11(&gxx[index]); PDstandard4th22gxx = PDstandard4th22(&gxx[index]); PDstandard4th33gxx = PDstandard4th33(&gxx[index]); PDstandard4th12gxx = PDstandard4th12(&gxx[index]); PDstandard4th13gxx = PDstandard4th13(&gxx[index]); PDstandard4th23gxx = PDstandard4th23(&gxx[index]); PDstandard4th1gxy = PDstandard4th1(&gxy[index]); PDstandard4th2gxy = PDstandard4th2(&gxy[index]); PDstandard4th3gxy = PDstandard4th3(&gxy[index]); PDstandard4th11gxy = PDstandard4th11(&gxy[index]); PDstandard4th22gxy = PDstandard4th22(&gxy[index]); PDstandard4th33gxy = PDstandard4th33(&gxy[index]); PDstandard4th12gxy = PDstandard4th12(&gxy[index]); PDstandard4th13gxy = PDstandard4th13(&gxy[index]); PDstandard4th23gxy = PDstandard4th23(&gxy[index]); PDstandard4th1gxz = PDstandard4th1(&gxz[index]); PDstandard4th2gxz = PDstandard4th2(&gxz[index]); PDstandard4th3gxz = PDstandard4th3(&gxz[index]); PDstandard4th11gxz = PDstandard4th11(&gxz[index]); PDstandard4th22gxz = PDstandard4th22(&gxz[index]); PDstandard4th33gxz = PDstandard4th33(&gxz[index]); PDstandard4th12gxz = PDstandard4th12(&gxz[index]); PDstandard4th13gxz = PDstandard4th13(&gxz[index]); PDstandard4th23gxz = PDstandard4th23(&gxz[index]); PDstandard4th1gyy = PDstandard4th1(&gyy[index]); PDstandard4th2gyy = PDstandard4th2(&gyy[index]); PDstandard4th3gyy = PDstandard4th3(&gyy[index]); PDstandard4th11gyy = PDstandard4th11(&gyy[index]); PDstandard4th22gyy = PDstandard4th22(&gyy[index]); PDstandard4th33gyy = PDstandard4th33(&gyy[index]); PDstandard4th12gyy = PDstandard4th12(&gyy[index]); PDstandard4th13gyy = PDstandard4th13(&gyy[index]); PDstandard4th23gyy = PDstandard4th23(&gyy[index]); PDstandard4th1gyz = PDstandard4th1(&gyz[index]); PDstandard4th2gyz = PDstandard4th2(&gyz[index]); PDstandard4th3gyz = PDstandard4th3(&gyz[index]); PDstandard4th11gyz = PDstandard4th11(&gyz[index]); PDstandard4th22gyz = PDstandard4th22(&gyz[index]); PDstandard4th33gyz = PDstandard4th33(&gyz[index]); PDstandard4th12gyz = PDstandard4th12(&gyz[index]); PDstandard4th13gyz = PDstandard4th13(&gyz[index]); PDstandard4th23gyz = PDstandard4th23(&gyz[index]); PDstandard4th1gzz = PDstandard4th1(&gzz[index]); PDstandard4th2gzz = PDstandard4th2(&gzz[index]); PDstandard4th3gzz = PDstandard4th3(&gzz[index]); PDstandard4th11gzz = PDstandard4th11(&gzz[index]); PDstandard4th22gzz = PDstandard4th22(&gzz[index]); PDstandard4th33gzz = PDstandard4th33(&gzz[index]); PDstandard4th12gzz = PDstandard4th12(&gzz[index]); PDstandard4th13gzz = PDstandard4th13(&gzz[index]); PDstandard4th23gzz = PDstandard4th23(&gzz[index]); PDstandard4th1kxx = PDstandard4th1(&kxx[index]); PDstandard4th2kxx = PDstandard4th2(&kxx[index]); PDstandard4th3kxx = PDstandard4th3(&kxx[index]); PDstandard4th1kxy = PDstandard4th1(&kxy[index]); PDstandard4th2kxy = PDstandard4th2(&kxy[index]); PDstandard4th3kxy = PDstandard4th3(&kxy[index]); PDstandard4th1kxz = PDstandard4th1(&kxz[index]); PDstandard4th2kxz = PDstandard4th2(&kxz[index]); PDstandard4th3kxz = PDstandard4th3(&kxz[index]); PDstandard4th1kyy = PDstandard4th1(&kyy[index]); PDstandard4th2kyy = PDstandard4th2(&kyy[index]); PDstandard4th3kyy = PDstandard4th3(&kyy[index]); PDstandard4th1kyz = PDstandard4th1(&kyz[index]); PDstandard4th2kyz = PDstandard4th2(&kyz[index]); PDstandard4th3kyz = PDstandard4th3(&kyz[index]); PDstandard4th1kzz = PDstandard4th1(&kzz[index]); PDstandard4th2kzz = PDstandard4th2(&kzz[index]); PDstandard4th3kzz = PDstandard4th3(&kzz[index]); break; } /* Calculate temporaries and grid functions */ CCTK_REAL_VEC JacPDstandard4th11gyy; CCTK_REAL_VEC JacPDstandard4th11gyz; CCTK_REAL_VEC JacPDstandard4th11gzz; CCTK_REAL_VEC JacPDstandard4th12gxy; CCTK_REAL_VEC JacPDstandard4th12gxz; CCTK_REAL_VEC JacPDstandard4th12gyz; CCTK_REAL_VEC JacPDstandard4th12gzz; CCTK_REAL_VEC JacPDstandard4th13gxz; CCTK_REAL_VEC JacPDstandard4th1gxx; CCTK_REAL_VEC JacPDstandard4th1gxy; CCTK_REAL_VEC JacPDstandard4th1gxz; CCTK_REAL_VEC JacPDstandard4th1gyy; CCTK_REAL_VEC JacPDstandard4th1gyz; CCTK_REAL_VEC JacPDstandard4th1gzz; CCTK_REAL_VEC JacPDstandard4th1kxy; CCTK_REAL_VEC JacPDstandard4th1kxz; CCTK_REAL_VEC JacPDstandard4th1kyy; CCTK_REAL_VEC JacPDstandard4th1kyz; CCTK_REAL_VEC JacPDstandard4th1kzz; CCTK_REAL_VEC JacPDstandard4th21gxy; CCTK_REAL_VEC JacPDstandard4th22gxx; CCTK_REAL_VEC JacPDstandard4th22gxz; CCTK_REAL_VEC JacPDstandard4th22gzz; CCTK_REAL_VEC JacPDstandard4th23gxx; CCTK_REAL_VEC JacPDstandard4th23gxy; CCTK_REAL_VEC JacPDstandard4th23gxz; CCTK_REAL_VEC JacPDstandard4th23gyz; CCTK_REAL_VEC JacPDstandard4th2gxx; CCTK_REAL_VEC JacPDstandard4th2gxy; CCTK_REAL_VEC JacPDstandard4th2gxz; CCTK_REAL_VEC JacPDstandard4th2gyy; CCTK_REAL_VEC JacPDstandard4th2gyz; CCTK_REAL_VEC JacPDstandard4th2gzz; CCTK_REAL_VEC JacPDstandard4th2kxx; CCTK_REAL_VEC JacPDstandard4th2kxy; CCTK_REAL_VEC JacPDstandard4th2kxz; CCTK_REAL_VEC JacPDstandard4th2kyz; CCTK_REAL_VEC JacPDstandard4th2kzz; CCTK_REAL_VEC JacPDstandard4th31gxy; CCTK_REAL_VEC JacPDstandard4th31gxz; CCTK_REAL_VEC JacPDstandard4th31gyy; CCTK_REAL_VEC JacPDstandard4th31gyz; CCTK_REAL_VEC JacPDstandard4th32gyz; CCTK_REAL_VEC JacPDstandard4th33gxx; CCTK_REAL_VEC JacPDstandard4th33gxy; CCTK_REAL_VEC JacPDstandard4th33gyy; CCTK_REAL_VEC JacPDstandard4th3gxx; CCTK_REAL_VEC JacPDstandard4th3gxy; CCTK_REAL_VEC JacPDstandard4th3gxz; CCTK_REAL_VEC JacPDstandard4th3gyy; CCTK_REAL_VEC JacPDstandard4th3gyz; CCTK_REAL_VEC JacPDstandard4th3gzz; CCTK_REAL_VEC JacPDstandard4th3kxx; CCTK_REAL_VEC JacPDstandard4th3kxy; CCTK_REAL_VEC JacPDstandard4th3kxz; CCTK_REAL_VEC JacPDstandard4th3kyy; CCTK_REAL_VEC JacPDstandard4th3kyz; if (use_jacobian) { JacPDstandard4th1gxx = kmadd(J11L,PDstandard4th1gxx,kmadd(J21L,PDstandard4th2gxx,kmul(J31L,PDstandard4th3gxx))); JacPDstandard4th1gxy = kmadd(J11L,PDstandard4th1gxy,kmadd(J21L,PDstandard4th2gxy,kmul(J31L,PDstandard4th3gxy))); JacPDstandard4th1gxz = kmadd(J11L,PDstandard4th1gxz,kmadd(J21L,PDstandard4th2gxz,kmul(J31L,PDstandard4th3gxz))); JacPDstandard4th1gyy = kmadd(J11L,PDstandard4th1gyy,kmadd(J21L,PDstandard4th2gyy,kmul(J31L,PDstandard4th3gyy))); JacPDstandard4th1gyz = kmadd(J11L,PDstandard4th1gyz,kmadd(J21L,PDstandard4th2gyz,kmul(J31L,PDstandard4th3gyz))); JacPDstandard4th1gzz = kmadd(J11L,PDstandard4th1gzz,kmadd(J21L,PDstandard4th2gzz,kmul(J31L,PDstandard4th3gzz))); JacPDstandard4th1kxy = kmadd(J11L,PDstandard4th1kxy,kmadd(J21L,PDstandard4th2kxy,kmul(J31L,PDstandard4th3kxy))); JacPDstandard4th1kxz = kmadd(J11L,PDstandard4th1kxz,kmadd(J21L,PDstandard4th2kxz,kmul(J31L,PDstandard4th3kxz))); JacPDstandard4th1kyy = kmadd(J11L,PDstandard4th1kyy,kmadd(J21L,PDstandard4th2kyy,kmul(J31L,PDstandard4th3kyy))); JacPDstandard4th1kyz = kmadd(J11L,PDstandard4th1kyz,kmadd(J21L,PDstandard4th2kyz,kmul(J31L,PDstandard4th3kyz))); JacPDstandard4th1kzz = kmadd(J11L,PDstandard4th1kzz,kmadd(J21L,PDstandard4th2kzz,kmul(J31L,PDstandard4th3kzz))); JacPDstandard4th2gxx = kmadd(J12L,PDstandard4th1gxx,kmadd(J22L,PDstandard4th2gxx,kmul(J32L,PDstandard4th3gxx))); JacPDstandard4th2gxy = kmadd(J12L,PDstandard4th1gxy,kmadd(J22L,PDstandard4th2gxy,kmul(J32L,PDstandard4th3gxy))); JacPDstandard4th2gxz = kmadd(J12L,PDstandard4th1gxz,kmadd(J22L,PDstandard4th2gxz,kmul(J32L,PDstandard4th3gxz))); JacPDstandard4th2gyy = kmadd(J12L,PDstandard4th1gyy,kmadd(J22L,PDstandard4th2gyy,kmul(J32L,PDstandard4th3gyy))); JacPDstandard4th2gyz = kmadd(J12L,PDstandard4th1gyz,kmadd(J22L,PDstandard4th2gyz,kmul(J32L,PDstandard4th3gyz))); JacPDstandard4th2gzz = kmadd(J12L,PDstandard4th1gzz,kmadd(J22L,PDstandard4th2gzz,kmul(J32L,PDstandard4th3gzz))); JacPDstandard4th2kxx = kmadd(J12L,PDstandard4th1kxx,kmadd(J22L,PDstandard4th2kxx,kmul(J32L,PDstandard4th3kxx))); JacPDstandard4th2kxy = kmadd(J12L,PDstandard4th1kxy,kmadd(J22L,PDstandard4th2kxy,kmul(J32L,PDstandard4th3kxy))); JacPDstandard4th2kxz = kmadd(J12L,PDstandard4th1kxz,kmadd(J22L,PDstandard4th2kxz,kmul(J32L,PDstandard4th3kxz))); JacPDstandard4th2kyz = kmadd(J12L,PDstandard4th1kyz,kmadd(J22L,PDstandard4th2kyz,kmul(J32L,PDstandard4th3kyz))); JacPDstandard4th2kzz = kmadd(J12L,PDstandard4th1kzz,kmadd(J22L,PDstandard4th2kzz,kmul(J32L,PDstandard4th3kzz))); JacPDstandard4th3gxx = kmadd(J13L,PDstandard4th1gxx,kmadd(J23L,PDstandard4th2gxx,kmul(J33L,PDstandard4th3gxx))); JacPDstandard4th3gxy = kmadd(J13L,PDstandard4th1gxy,kmadd(J23L,PDstandard4th2gxy,kmul(J33L,PDstandard4th3gxy))); JacPDstandard4th3gxz = kmadd(J13L,PDstandard4th1gxz,kmadd(J23L,PDstandard4th2gxz,kmul(J33L,PDstandard4th3gxz))); JacPDstandard4th3gyy = kmadd(J13L,PDstandard4th1gyy,kmadd(J23L,PDstandard4th2gyy,kmul(J33L,PDstandard4th3gyy))); JacPDstandard4th3gyz = kmadd(J13L,PDstandard4th1gyz,kmadd(J23L,PDstandard4th2gyz,kmul(J33L,PDstandard4th3gyz))); JacPDstandard4th3gzz = kmadd(J13L,PDstandard4th1gzz,kmadd(J23L,PDstandard4th2gzz,kmul(J33L,PDstandard4th3gzz))); JacPDstandard4th3kxx = kmadd(J13L,PDstandard4th1kxx,kmadd(J23L,PDstandard4th2kxx,kmul(J33L,PDstandard4th3kxx))); JacPDstandard4th3kxy = kmadd(J13L,PDstandard4th1kxy,kmadd(J23L,PDstandard4th2kxy,kmul(J33L,PDstandard4th3kxy))); JacPDstandard4th3kxz = kmadd(J13L,PDstandard4th1kxz,kmadd(J23L,PDstandard4th2kxz,kmul(J33L,PDstandard4th3kxz))); JacPDstandard4th3kyy = kmadd(J13L,PDstandard4th1kyy,kmadd(J23L,PDstandard4th2kyy,kmul(J33L,PDstandard4th3kyy))); JacPDstandard4th3kyz = kmadd(J13L,PDstandard4th1kyz,kmadd(J23L,PDstandard4th2kyz,kmul(J33L,PDstandard4th3kyz))); JacPDstandard4th11gyy = kmadd(dJ111L,PDstandard4th1gyy,kmadd(dJ211L,PDstandard4th2gyy,kmadd(dJ311L,PDstandard4th3gyy,kmadd(PDstandard4th11gyy,kmul(J11L,J11L),kmadd(PDstandard4th22gyy,kmul(J21L,J21L),kmadd(PDstandard4th33gyy,kmul(J31L,J31L),kmul(kmadd(J11L,kmadd(J21L,PDstandard4th12gyy,kmul(J31L,PDstandard4th13gyy)),kmul(J21L,kmul(J31L,PDstandard4th23gyy))),ToReal(2.)))))))); JacPDstandard4th11gyz = kmadd(dJ111L,PDstandard4th1gyz,kmadd(dJ211L,PDstandard4th2gyz,kmadd(dJ311L,PDstandard4th3gyz,kmadd(PDstandard4th11gyz,kmul(J11L,J11L),kmadd(PDstandard4th22gyz,kmul(J21L,J21L),kmadd(PDstandard4th33gyz,kmul(J31L,J31L),kmul(kmadd(J11L,kmadd(J21L,PDstandard4th12gyz,kmul(J31L,PDstandard4th13gyz)),kmul(J21L,kmul(J31L,PDstandard4th23gyz))),ToReal(2.)))))))); JacPDstandard4th11gzz = kmadd(dJ111L,PDstandard4th1gzz,kmadd(dJ211L,PDstandard4th2gzz,kmadd(dJ311L,PDstandard4th3gzz,kmadd(PDstandard4th11gzz,kmul(J11L,J11L),kmadd(PDstandard4th22gzz,kmul(J21L,J21L),kmadd(PDstandard4th33gzz,kmul(J31L,J31L),kmul(kmadd(J11L,kmadd(J21L,PDstandard4th12gzz,kmul(J31L,PDstandard4th13gzz)),kmul(J21L,kmul(J31L,PDstandard4th23gzz))),ToReal(2.)))))))); JacPDstandard4th22gxx = kmadd(dJ122L,PDstandard4th1gxx,kmadd(dJ222L,PDstandard4th2gxx,kmadd(dJ322L,PDstandard4th3gxx,kmadd(PDstandard4th11gxx,kmul(J12L,J12L),kmadd(PDstandard4th22gxx,kmul(J22L,J22L),kmadd(PDstandard4th33gxx,kmul(J32L,J32L),kmul(kmadd(J12L,kmadd(J22L,PDstandard4th12gxx,kmul(J32L,PDstandard4th13gxx)),kmul(J22L,kmul(J32L,PDstandard4th23gxx))),ToReal(2.)))))))); JacPDstandard4th22gxz = kmadd(dJ122L,PDstandard4th1gxz,kmadd(dJ222L,PDstandard4th2gxz,kmadd(dJ322L,PDstandard4th3gxz,kmadd(PDstandard4th11gxz,kmul(J12L,J12L),kmadd(PDstandard4th22gxz,kmul(J22L,J22L),kmadd(PDstandard4th33gxz,kmul(J32L,J32L),kmul(kmadd(J12L,kmadd(J22L,PDstandard4th12gxz,kmul(J32L,PDstandard4th13gxz)),kmul(J22L,kmul(J32L,PDstandard4th23gxz))),ToReal(2.)))))))); JacPDstandard4th22gzz = kmadd(dJ122L,PDstandard4th1gzz,kmadd(dJ222L,PDstandard4th2gzz,kmadd(dJ322L,PDstandard4th3gzz,kmadd(PDstandard4th11gzz,kmul(J12L,J12L),kmadd(PDstandard4th22gzz,kmul(J22L,J22L),kmadd(PDstandard4th33gzz,kmul(J32L,J32L),kmul(kmadd(J12L,kmadd(J22L,PDstandard4th12gzz,kmul(J32L,PDstandard4th13gzz)),kmul(J22L,kmul(J32L,PDstandard4th23gzz))),ToReal(2.)))))))); JacPDstandard4th33gxx = kmadd(dJ133L,PDstandard4th1gxx,kmadd(dJ233L,PDstandard4th2gxx,kmadd(dJ333L,PDstandard4th3gxx,kmadd(PDstandard4th11gxx,kmul(J13L,J13L),kmadd(PDstandard4th22gxx,kmul(J23L,J23L),kmadd(PDstandard4th33gxx,kmul(J33L,J33L),kmul(kmadd(J13L,kmadd(J23L,PDstandard4th12gxx,kmul(J33L,PDstandard4th13gxx)),kmul(J23L,kmul(J33L,PDstandard4th23gxx))),ToReal(2.)))))))); JacPDstandard4th33gxy = kmadd(dJ133L,PDstandard4th1gxy,kmadd(dJ233L,PDstandard4th2gxy,kmadd(dJ333L,PDstandard4th3gxy,kmadd(PDstandard4th11gxy,kmul(J13L,J13L),kmadd(PDstandard4th22gxy,kmul(J23L,J23L),kmadd(PDstandard4th33gxy,kmul(J33L,J33L),kmul(kmadd(J13L,kmadd(J23L,PDstandard4th12gxy,kmul(J33L,PDstandard4th13gxy)),kmul(J23L,kmul(J33L,PDstandard4th23gxy))),ToReal(2.)))))))); JacPDstandard4th33gyy = kmadd(dJ133L,PDstandard4th1gyy,kmadd(dJ233L,PDstandard4th2gyy,kmadd(dJ333L,PDstandard4th3gyy,kmadd(PDstandard4th11gyy,kmul(J13L,J13L),kmadd(PDstandard4th22gyy,kmul(J23L,J23L),kmadd(PDstandard4th33gyy,kmul(J33L,J33L),kmul(kmadd(J13L,kmadd(J23L,PDstandard4th12gyy,kmul(J33L,PDstandard4th13gyy)),kmul(J23L,kmul(J33L,PDstandard4th23gyy))),ToReal(2.)))))))); JacPDstandard4th12gxy = kmadd(J12L,kmadd(J11L,PDstandard4th11gxy,kmadd(J21L,PDstandard4th12gxy,kmul(J31L,PDstandard4th13gxy))),kmadd(J11L,kmadd(J22L,PDstandard4th12gxy,kmul(J32L,PDstandard4th13gxy)),kmadd(dJ112L,PDstandard4th1gxy,kmadd(J22L,kmadd(J21L,PDstandard4th22gxy,kmul(J31L,PDstandard4th23gxy)),kmadd(dJ212L,PDstandard4th2gxy,kmadd(J32L,kmadd(J21L,PDstandard4th23gxy,kmul(J31L,PDstandard4th33gxy)),kmul(dJ312L,PDstandard4th3gxy))))))); JacPDstandard4th12gxz = kmadd(J12L,kmadd(J11L,PDstandard4th11gxz,kmadd(J21L,PDstandard4th12gxz,kmul(J31L,PDstandard4th13gxz))),kmadd(J11L,kmadd(J22L,PDstandard4th12gxz,kmul(J32L,PDstandard4th13gxz)),kmadd(dJ112L,PDstandard4th1gxz,kmadd(J22L,kmadd(J21L,PDstandard4th22gxz,kmul(J31L,PDstandard4th23gxz)),kmadd(dJ212L,PDstandard4th2gxz,kmadd(J32L,kmadd(J21L,PDstandard4th23gxz,kmul(J31L,PDstandard4th33gxz)),kmul(dJ312L,PDstandard4th3gxz))))))); JacPDstandard4th12gyz = kmadd(J12L,kmadd(J11L,PDstandard4th11gyz,kmadd(J21L,PDstandard4th12gyz,kmul(J31L,PDstandard4th13gyz))),kmadd(J11L,kmadd(J22L,PDstandard4th12gyz,kmul(J32L,PDstandard4th13gyz)),kmadd(dJ112L,PDstandard4th1gyz,kmadd(J22L,kmadd(J21L,PDstandard4th22gyz,kmul(J31L,PDstandard4th23gyz)),kmadd(dJ212L,PDstandard4th2gyz,kmadd(J32L,kmadd(J21L,PDstandard4th23gyz,kmul(J31L,PDstandard4th33gyz)),kmul(dJ312L,PDstandard4th3gyz))))))); JacPDstandard4th12gzz = kmadd(J12L,kmadd(J11L,PDstandard4th11gzz,kmadd(J21L,PDstandard4th12gzz,kmul(J31L,PDstandard4th13gzz))),kmadd(J11L,kmadd(J22L,PDstandard4th12gzz,kmul(J32L,PDstandard4th13gzz)),kmadd(dJ112L,PDstandard4th1gzz,kmadd(J22L,kmadd(J21L,PDstandard4th22gzz,kmul(J31L,PDstandard4th23gzz)),kmadd(dJ212L,PDstandard4th2gzz,kmadd(J32L,kmadd(J21L,PDstandard4th23gzz,kmul(J31L,PDstandard4th33gzz)),kmul(dJ312L,PDstandard4th3gzz))))))); JacPDstandard4th13gxz = kmadd(J13L,kmadd(J11L,PDstandard4th11gxz,kmadd(J21L,PDstandard4th12gxz,kmul(J31L,PDstandard4th13gxz))),kmadd(J11L,kmadd(J23L,PDstandard4th12gxz,kmul(J33L,PDstandard4th13gxz)),kmadd(dJ113L,PDstandard4th1gxz,kmadd(J23L,kmadd(J21L,PDstandard4th22gxz,kmul(J31L,PDstandard4th23gxz)),kmadd(dJ213L,PDstandard4th2gxz,kmadd(J33L,kmadd(J21L,PDstandard4th23gxz,kmul(J31L,PDstandard4th33gxz)),kmul(dJ313L,PDstandard4th3gxz))))))); JacPDstandard4th21gxy = kmadd(J12L,kmadd(J11L,PDstandard4th11gxy,kmadd(J21L,PDstandard4th12gxy,kmul(J31L,PDstandard4th13gxy))),kmadd(J11L,kmadd(J22L,PDstandard4th12gxy,kmul(J32L,PDstandard4th13gxy)),kmadd(dJ112L,PDstandard4th1gxy,kmadd(J22L,kmadd(J21L,PDstandard4th22gxy,kmul(J31L,PDstandard4th23gxy)),kmadd(dJ212L,PDstandard4th2gxy,kmadd(J32L,kmadd(J21L,PDstandard4th23gxy,kmul(J31L,PDstandard4th33gxy)),kmul(dJ312L,PDstandard4th3gxy))))))); JacPDstandard4th23gxx = kmadd(J13L,kmadd(J12L,PDstandard4th11gxx,kmadd(J22L,PDstandard4th12gxx,kmul(J32L,PDstandard4th13gxx))),kmadd(J12L,kmadd(J23L,PDstandard4th12gxx,kmul(J33L,PDstandard4th13gxx)),kmadd(dJ123L,PDstandard4th1gxx,kmadd(J23L,kmadd(J22L,PDstandard4th22gxx,kmul(J32L,PDstandard4th23gxx)),kmadd(dJ223L,PDstandard4th2gxx,kmadd(J33L,kmadd(J22L,PDstandard4th23gxx,kmul(J32L,PDstandard4th33gxx)),kmul(dJ323L,PDstandard4th3gxx))))))); JacPDstandard4th23gxy = kmadd(J13L,kmadd(J12L,PDstandard4th11gxy,kmadd(J22L,PDstandard4th12gxy,kmul(J32L,PDstandard4th13gxy))),kmadd(J12L,kmadd(J23L,PDstandard4th12gxy,kmul(J33L,PDstandard4th13gxy)),kmadd(dJ123L,PDstandard4th1gxy,kmadd(J23L,kmadd(J22L,PDstandard4th22gxy,kmul(J32L,PDstandard4th23gxy)),kmadd(dJ223L,PDstandard4th2gxy,kmadd(J33L,kmadd(J22L,PDstandard4th23gxy,kmul(J32L,PDstandard4th33gxy)),kmul(dJ323L,PDstandard4th3gxy))))))); JacPDstandard4th23gxz = kmadd(J13L,kmadd(J12L,PDstandard4th11gxz,kmadd(J22L,PDstandard4th12gxz,kmul(J32L,PDstandard4th13gxz))),kmadd(J12L,kmadd(J23L,PDstandard4th12gxz,kmul(J33L,PDstandard4th13gxz)),kmadd(dJ123L,PDstandard4th1gxz,kmadd(J23L,kmadd(J22L,PDstandard4th22gxz,kmul(J32L,PDstandard4th23gxz)),kmadd(dJ223L,PDstandard4th2gxz,kmadd(J33L,kmadd(J22L,PDstandard4th23gxz,kmul(J32L,PDstandard4th33gxz)),kmul(dJ323L,PDstandard4th3gxz))))))); JacPDstandard4th23gyz = kmadd(J13L,kmadd(J12L,PDstandard4th11gyz,kmadd(J22L,PDstandard4th12gyz,kmul(J32L,PDstandard4th13gyz))),kmadd(J12L,kmadd(J23L,PDstandard4th12gyz,kmul(J33L,PDstandard4th13gyz)),kmadd(dJ123L,PDstandard4th1gyz,kmadd(J23L,kmadd(J22L,PDstandard4th22gyz,kmul(J32L,PDstandard4th23gyz)),kmadd(dJ223L,PDstandard4th2gyz,kmadd(J33L,kmadd(J22L,PDstandard4th23gyz,kmul(J32L,PDstandard4th33gyz)),kmul(dJ323L,PDstandard4th3gyz))))))); JacPDstandard4th31gxy = kmadd(J13L,kmadd(J11L,PDstandard4th11gxy,kmadd(J21L,PDstandard4th12gxy,kmul(J31L,PDstandard4th13gxy))),kmadd(J11L,kmadd(J23L,PDstandard4th12gxy,kmul(J33L,PDstandard4th13gxy)),kmadd(dJ113L,PDstandard4th1gxy,kmadd(J23L,kmadd(J21L,PDstandard4th22gxy,kmul(J31L,PDstandard4th23gxy)),kmadd(dJ213L,PDstandard4th2gxy,kmadd(J33L,kmadd(J21L,PDstandard4th23gxy,kmul(J31L,PDstandard4th33gxy)),kmul(dJ313L,PDstandard4th3gxy))))))); JacPDstandard4th31gxz = kmadd(J13L,kmadd(J11L,PDstandard4th11gxz,kmadd(J21L,PDstandard4th12gxz,kmul(J31L,PDstandard4th13gxz))),kmadd(J11L,kmadd(J23L,PDstandard4th12gxz,kmul(J33L,PDstandard4th13gxz)),kmadd(dJ113L,PDstandard4th1gxz,kmadd(J23L,kmadd(J21L,PDstandard4th22gxz,kmul(J31L,PDstandard4th23gxz)),kmadd(dJ213L,PDstandard4th2gxz,kmadd(J33L,kmadd(J21L,PDstandard4th23gxz,kmul(J31L,PDstandard4th33gxz)),kmul(dJ313L,PDstandard4th3gxz))))))); JacPDstandard4th31gyy = kmadd(J13L,kmadd(J11L,PDstandard4th11gyy,kmadd(J21L,PDstandard4th12gyy,kmul(J31L,PDstandard4th13gyy))),kmadd(J11L,kmadd(J23L,PDstandard4th12gyy,kmul(J33L,PDstandard4th13gyy)),kmadd(dJ113L,PDstandard4th1gyy,kmadd(J23L,kmadd(J21L,PDstandard4th22gyy,kmul(J31L,PDstandard4th23gyy)),kmadd(dJ213L,PDstandard4th2gyy,kmadd(J33L,kmadd(J21L,PDstandard4th23gyy,kmul(J31L,PDstandard4th33gyy)),kmul(dJ313L,PDstandard4th3gyy))))))); JacPDstandard4th31gyz = kmadd(J13L,kmadd(J11L,PDstandard4th11gyz,kmadd(J21L,PDstandard4th12gyz,kmul(J31L,PDstandard4th13gyz))),kmadd(J11L,kmadd(J23L,PDstandard4th12gyz,kmul(J33L,PDstandard4th13gyz)),kmadd(dJ113L,PDstandard4th1gyz,kmadd(J23L,kmadd(J21L,PDstandard4th22gyz,kmul(J31L,PDstandard4th23gyz)),kmadd(dJ213L,PDstandard4th2gyz,kmadd(J33L,kmadd(J21L,PDstandard4th23gyz,kmul(J31L,PDstandard4th33gyz)),kmul(dJ313L,PDstandard4th3gyz))))))); JacPDstandard4th32gyz = kmadd(J13L,kmadd(J12L,PDstandard4th11gyz,kmadd(J22L,PDstandard4th12gyz,kmul(J32L,PDstandard4th13gyz))),kmadd(J12L,kmadd(J23L,PDstandard4th12gyz,kmul(J33L,PDstandard4th13gyz)),kmadd(dJ123L,PDstandard4th1gyz,kmadd(J23L,kmadd(J22L,PDstandard4th22gyz,kmul(J32L,PDstandard4th23gyz)),kmadd(dJ223L,PDstandard4th2gyz,kmadd(J33L,kmadd(J22L,PDstandard4th23gyz,kmul(J32L,PDstandard4th33gyz)),kmul(dJ323L,PDstandard4th3gyz))))))); } else { JacPDstandard4th1gxx = PDstandard4th1gxx; JacPDstandard4th1gxy = PDstandard4th1gxy; JacPDstandard4th1gxz = PDstandard4th1gxz; JacPDstandard4th1gyy = PDstandard4th1gyy; JacPDstandard4th1gyz = PDstandard4th1gyz; JacPDstandard4th1gzz = PDstandard4th1gzz; JacPDstandard4th1kxy = PDstandard4th1kxy; JacPDstandard4th1kxz = PDstandard4th1kxz; JacPDstandard4th1kyy = PDstandard4th1kyy; JacPDstandard4th1kyz = PDstandard4th1kyz; JacPDstandard4th1kzz = PDstandard4th1kzz; JacPDstandard4th2gxx = PDstandard4th2gxx; JacPDstandard4th2gxy = PDstandard4th2gxy; JacPDstandard4th2gxz = PDstandard4th2gxz; JacPDstandard4th2gyy = PDstandard4th2gyy; JacPDstandard4th2gyz = PDstandard4th2gyz; JacPDstandard4th2gzz = PDstandard4th2gzz; JacPDstandard4th2kxx = PDstandard4th2kxx; JacPDstandard4th2kxy = PDstandard4th2kxy; JacPDstandard4th2kxz = PDstandard4th2kxz; JacPDstandard4th2kyz = PDstandard4th2kyz; JacPDstandard4th2kzz = PDstandard4th2kzz; JacPDstandard4th3gxx = PDstandard4th3gxx; JacPDstandard4th3gxy = PDstandard4th3gxy; JacPDstandard4th3gxz = PDstandard4th3gxz; JacPDstandard4th3gyy = PDstandard4th3gyy; JacPDstandard4th3gyz = PDstandard4th3gyz; JacPDstandard4th3gzz = PDstandard4th3gzz; JacPDstandard4th3kxx = PDstandard4th3kxx; JacPDstandard4th3kxy = PDstandard4th3kxy; JacPDstandard4th3kxz = PDstandard4th3kxz; JacPDstandard4th3kyy = PDstandard4th3kyy; JacPDstandard4th3kyz = PDstandard4th3kyz; JacPDstandard4th11gyy = PDstandard4th11gyy; JacPDstandard4th11gyz = PDstandard4th11gyz; JacPDstandard4th11gzz = PDstandard4th11gzz; JacPDstandard4th22gxx = PDstandard4th22gxx; JacPDstandard4th22gxz = PDstandard4th22gxz; JacPDstandard4th22gzz = PDstandard4th22gzz; JacPDstandard4th33gxx = PDstandard4th33gxx; JacPDstandard4th33gxy = PDstandard4th33gxy; JacPDstandard4th33gyy = PDstandard4th33gyy; JacPDstandard4th12gxy = PDstandard4th12gxy; JacPDstandard4th12gxz = PDstandard4th12gxz; JacPDstandard4th12gyz = PDstandard4th12gyz; JacPDstandard4th12gzz = PDstandard4th12gzz; JacPDstandard4th13gxz = PDstandard4th13gxz; JacPDstandard4th21gxy = PDstandard4th12gxy; JacPDstandard4th23gxx = PDstandard4th23gxx; JacPDstandard4th23gxy = PDstandard4th23gxy; JacPDstandard4th23gxz = PDstandard4th23gxz; JacPDstandard4th23gyz = PDstandard4th23gyz; JacPDstandard4th31gxy = PDstandard4th13gxy; JacPDstandard4th31gxz = PDstandard4th13gxz; JacPDstandard4th31gyy = PDstandard4th13gyy; JacPDstandard4th31gyz = PDstandard4th13gyz; JacPDstandard4th32gyz = PDstandard4th23gyz; } CCTK_REAL_VEC detg = knmsub(gyyL,kmul(gxzL,gxzL),knmsub(gxxL,kmul(gyzL,gyzL),kmadd(gzzL,kmsub(gxxL,gyyL,kmul(gxyL,gxyL)),kmul(gxyL,kmul(gxzL,kmul(gyzL,ToReal(2.))))))); CCTK_REAL_VEC invdetg = kdiv(ToReal(1.),detg); CCTK_REAL_VEC gInv11 = kmul(invdetg,kmsub(gyyL,gzzL,kmul(gyzL,gyzL))); CCTK_REAL_VEC gInv12 = kmul(invdetg,kmsub(gxzL,gyzL,kmul(gxyL,gzzL))); CCTK_REAL_VEC gInv13 = kmul(invdetg,kmsub(gxyL,gyzL,kmul(gxzL,gyyL))); CCTK_REAL_VEC gInv21 = kmul(invdetg,kmsub(gxzL,gyzL,kmul(gxyL,gzzL))); CCTK_REAL_VEC gInv22 = kmul(invdetg,kmsub(gxxL,gzzL,kmul(gxzL,gxzL))); CCTK_REAL_VEC gInv23 = kmul(invdetg,kmsub(gxyL,gxzL,kmul(gxxL,gyzL))); CCTK_REAL_VEC gInv31 = kmul(invdetg,kmsub(gxyL,gyzL,kmul(gxzL,gyyL))); CCTK_REAL_VEC gInv32 = kmul(invdetg,kmsub(gxyL,gxzL,kmul(gxxL,gyzL))); CCTK_REAL_VEC gInv33 = kmul(invdetg,kmsub(gxxL,gyyL,kmul(gxyL,gxyL))); CCTK_REAL_VEC gamma111 = kmul(ToReal(0.5),kmadd(gInv11,JacPDstandard4th1gxx,knmsub(gInv12,JacPDstandard4th2gxx,kmsub(kmadd(gInv12,JacPDstandard4th1gxy,kmul(gInv13,JacPDstandard4th1gxz)),ToReal(2.),kmul(gInv13,JacPDstandard4th3gxx))))); CCTK_REAL_VEC gamma211 = kmul(ToReal(0.5),kmadd(gInv21,JacPDstandard4th1gxx,knmsub(gInv22,JacPDstandard4th2gxx,kmsub(kmadd(gInv22,JacPDstandard4th1gxy,kmul(gInv23,JacPDstandard4th1gxz)),ToReal(2.),kmul(gInv23,JacPDstandard4th3gxx))))); CCTK_REAL_VEC gamma311 = kmul(ToReal(0.5),kmadd(gInv31,JacPDstandard4th1gxx,knmsub(gInv32,JacPDstandard4th2gxx,kmsub(kmadd(gInv32,JacPDstandard4th1gxy,kmul(gInv33,JacPDstandard4th1gxz)),ToReal(2.),kmul(gInv33,JacPDstandard4th3gxx))))); CCTK_REAL_VEC gamma121 = kmul(kmadd(gInv12,JacPDstandard4th1gyy,kmadd(gInv11,JacPDstandard4th2gxx,kmul(gInv13,kadd(JacPDstandard4th1gyz,ksub(JacPDstandard4th2gxz,JacPDstandard4th3gxy))))),ToReal(0.5)); CCTK_REAL_VEC gamma221 = kmul(kmadd(gInv22,JacPDstandard4th1gyy,kmadd(gInv21,JacPDstandard4th2gxx,kmul(gInv23,kadd(JacPDstandard4th1gyz,ksub(JacPDstandard4th2gxz,JacPDstandard4th3gxy))))),ToReal(0.5)); CCTK_REAL_VEC gamma321 = kmul(kmadd(gInv32,JacPDstandard4th1gyy,kmadd(gInv31,JacPDstandard4th2gxx,kmul(gInv33,kadd(JacPDstandard4th1gyz,ksub(JacPDstandard4th2gxz,JacPDstandard4th3gxy))))),ToReal(0.5)); CCTK_REAL_VEC gamma131 = kmul(kmadd(gInv13,JacPDstandard4th1gzz,kmadd(gInv11,JacPDstandard4th3gxx,kmul(gInv12,kadd(JacPDstandard4th1gyz,ksub(JacPDstandard4th3gxy,JacPDstandard4th2gxz))))),ToReal(0.5)); CCTK_REAL_VEC gamma231 = kmul(kmadd(gInv23,JacPDstandard4th1gzz,kmadd(gInv21,JacPDstandard4th3gxx,kmul(gInv22,kadd(JacPDstandard4th1gyz,ksub(JacPDstandard4th3gxy,JacPDstandard4th2gxz))))),ToReal(0.5)); CCTK_REAL_VEC gamma331 = kmul(kmadd(gInv33,JacPDstandard4th1gzz,kmadd(gInv31,JacPDstandard4th3gxx,kmul(gInv32,kadd(JacPDstandard4th1gyz,ksub(JacPDstandard4th3gxy,JacPDstandard4th2gxz))))),ToReal(0.5)); CCTK_REAL_VEC gamma122 = kmul(ToReal(0.5),kmadd(gInv12,JacPDstandard4th2gyy,kmadd(gInv11,kmsub(JacPDstandard4th2gxy,ToReal(2.),JacPDstandard4th1gyy),kmul(gInv13,kmsub(JacPDstandard4th2gyz,ToReal(2.),JacPDstandard4th3gyy))))); CCTK_REAL_VEC gamma222 = kmul(ToReal(0.5),kmadd(gInv22,JacPDstandard4th2gyy,kmadd(gInv21,kmsub(JacPDstandard4th2gxy,ToReal(2.),JacPDstandard4th1gyy),kmul(gInv23,kmsub(JacPDstandard4th2gyz,ToReal(2.),JacPDstandard4th3gyy))))); CCTK_REAL_VEC gamma322 = kmul(ToReal(0.5),kmadd(gInv32,JacPDstandard4th2gyy,kmadd(gInv31,kmsub(JacPDstandard4th2gxy,ToReal(2.),JacPDstandard4th1gyy),kmul(gInv33,kmsub(JacPDstandard4th2gyz,ToReal(2.),JacPDstandard4th3gyy))))); CCTK_REAL_VEC gamma132 = kmul(kmadd(gInv13,JacPDstandard4th2gzz,kmadd(gInv12,JacPDstandard4th3gyy,kmul(gInv11,kadd(JacPDstandard4th2gxz,ksub(JacPDstandard4th3gxy,JacPDstandard4th1gyz))))),ToReal(0.5)); CCTK_REAL_VEC gamma232 = kmul(kmadd(gInv23,JacPDstandard4th2gzz,kmadd(gInv22,JacPDstandard4th3gyy,kmul(gInv21,kadd(JacPDstandard4th2gxz,ksub(JacPDstandard4th3gxy,JacPDstandard4th1gyz))))),ToReal(0.5)); CCTK_REAL_VEC gamma332 = kmul(kmadd(gInv33,JacPDstandard4th2gzz,kmadd(gInv32,JacPDstandard4th3gyy,kmul(gInv31,kadd(JacPDstandard4th2gxz,ksub(JacPDstandard4th3gxy,JacPDstandard4th1gyz))))),ToReal(0.5)); CCTK_REAL_VEC gamma133 = kmul(ToReal(0.5),kmadd(gInv13,JacPDstandard4th3gzz,kmadd(gInv11,kmsub(JacPDstandard4th3gxz,ToReal(2.),JacPDstandard4th1gzz),kmul(gInv12,kmsub(JacPDstandard4th3gyz,ToReal(2.),JacPDstandard4th2gzz))))); CCTK_REAL_VEC gamma233 = kmul(ToReal(0.5),kmadd(gInv23,JacPDstandard4th3gzz,kmadd(gInv21,kmsub(JacPDstandard4th3gxz,ToReal(2.),JacPDstandard4th1gzz),kmul(gInv22,kmsub(JacPDstandard4th3gyz,ToReal(2.),JacPDstandard4th2gzz))))); CCTK_REAL_VEC gamma333 = kmul(ToReal(0.5),kmadd(gInv33,JacPDstandard4th3gzz,kmadd(gInv31,kmsub(JacPDstandard4th3gxz,ToReal(2.),JacPDstandard4th1gzz),kmul(gInv32,kmsub(JacPDstandard4th3gyz,ToReal(2.),JacPDstandard4th2gzz))))); CCTK_REAL_VEC xmoved = ksub(xL,ToReal(xorig)); CCTK_REAL_VEC ymoved = ksub(yL,ToReal(yorig)); CCTK_REAL_VEC zmoved = ksub(zL,ToReal(zorig)); CCTK_REAL_VEC va1 = kneg(ymoved); CCTK_REAL_VEC va2 = kadd(xmoved,ToReal(offset)); CCTK_REAL_VEC va3 = ToReal(0.); CCTK_REAL_VEC vb1 = kadd(xmoved,ToReal(offset)); CCTK_REAL_VEC vb2 = ymoved; CCTK_REAL_VEC vb3 = zmoved; CCTK_REAL_VEC vc1 = kmul(ksqrt(detg),kmadd(vb3,kmsub(gInv11,va2,kmul(gInv12,va1)),kmadd(vb1,kmsub(gInv12,va3,kmul(gInv13,va2)),kmul(vb2,kmsub(gInv13,va1,kmul(gInv11,va3)))))); CCTK_REAL_VEC vc2 = kmul(ksqrt(detg),kmadd(vb3,kmsub(gInv21,va2,kmul(gInv22,va1)),kmadd(vb1,kmsub(gInv22,va3,kmul(gInv23,va2)),kmul(vb2,kmsub(gInv23,va1,kmul(gInv21,va3)))))); CCTK_REAL_VEC vc3 = kmul(ksqrt(detg),kmadd(vb3,kmsub(gInv31,va2,kmul(gInv32,va1)),kmadd(vb1,kmsub(gInv32,va3,kmul(gInv33,va2)),kmul(vb2,kmsub(gInv33,va1,kmul(gInv31,va3)))))); CCTK_REAL_VEC wa1 = va1; CCTK_REAL_VEC wa2 = va2; CCTK_REAL_VEC wa3 = va3; CCTK_REAL_VEC omega11 = kmadd(gxxL,kmul(wa1,wa1),kmadd(gyyL,kmul(wa2,wa2),kmadd(gzzL,kmul(wa3,wa3),kmul(kmadd(gyzL,kmul(wa2,wa3),kmul(wa1,kmadd(gxyL,wa2,kmul(gxzL,wa3)))),ToReal(2.))))); CCTK_REAL_VEC ea1 = kdiv(wa1,ksqrt(omega11)); CCTK_REAL_VEC ea2 = kdiv(wa2,ksqrt(omega11)); CCTK_REAL_VEC ea3 = kdiv(wa3,ksqrt(omega11)); CCTK_REAL_VEC omega12 = kmadd(ea1,kmadd(gxxL,vb1,kmadd(gxyL,vb2,kmul(gxzL,vb3))),kmadd(ea2,kmadd(gxyL,vb1,kmadd(gyyL,vb2,kmul(gyzL,vb3))),kmul(ea3,kmadd(gxzL,vb1,kmadd(gyzL,vb2,kmul(gzzL,vb3)))))); CCTK_REAL_VEC wb1 = knmsub(ea1,omega12,vb1); CCTK_REAL_VEC wb2 = knmsub(ea2,omega12,vb2); CCTK_REAL_VEC wb3 = knmsub(ea3,omega12,vb3); CCTK_REAL_VEC omega22 = kmadd(gxxL,kmul(wb1,wb1),kmadd(gyyL,kmul(wb2,wb2),kmadd(gzzL,kmul(wb3,wb3),kmul(kmadd(gyzL,kmul(wb2,wb3),kmul(wb1,kmadd(gxyL,wb2,kmul(gxzL,wb3)))),ToReal(2.))))); CCTK_REAL_VEC eb1 = kdiv(wb1,ksqrt(omega22)); CCTK_REAL_VEC eb2 = kdiv(wb2,ksqrt(omega22)); CCTK_REAL_VEC eb3 = kdiv(wb3,ksqrt(omega22)); CCTK_REAL_VEC omega13 = kmadd(ea1,kmadd(gxxL,vc1,kmadd(gxyL,vc2,kmul(gxzL,vc3))),kmadd(ea2,kmadd(gxyL,vc1,kmadd(gyyL,vc2,kmul(gyzL,vc3))),kmul(ea3,kmadd(gxzL,vc1,kmadd(gyzL,vc2,kmul(gzzL,vc3)))))); CCTK_REAL_VEC omega23 = kmadd(eb1,kmadd(gxxL,vc1,kmadd(gxyL,vc2,kmul(gxzL,vc3))),kmadd(eb2,kmadd(gxyL,vc1,kmadd(gyyL,vc2,kmul(gyzL,vc3))),kmul(eb3,kmadd(gxzL,vc1,kmadd(gyzL,vc2,kmul(gzzL,vc3)))))); CCTK_REAL_VEC wc1 = ksub(vc1,kmadd(eb1,omega23,kmul(ea1,omega13))); CCTK_REAL_VEC wc2 = ksub(vc2,kmadd(eb2,omega23,kmul(ea2,omega13))); CCTK_REAL_VEC wc3 = ksub(vc3,kmadd(eb3,omega23,kmul(ea3,omega13))); CCTK_REAL_VEC omega33 = kmadd(gxxL,kmul(wc1,wc1),kmadd(gyyL,kmul(wc2,wc2),kmadd(gzzL,kmul(wc3,wc3),kmul(kmadd(gyzL,kmul(wc2,wc3),kmul(wc1,kmadd(gxyL,wc2,kmul(gxzL,wc3)))),ToReal(2.))))); CCTK_REAL_VEC ec1 = kdiv(wc1,ksqrt(omega33)); CCTK_REAL_VEC ec2 = kdiv(wc2,ksqrt(omega33)); CCTK_REAL_VEC ec3 = kdiv(wc3,ksqrt(omega33)); CCTK_REAL_VEC isqrt2 = ToReal(0.707106781186547524); CCTK_REAL_VEC n1 = kneg(kmul(eb1,isqrt2)); CCTK_REAL_VEC n2 = kneg(kmul(eb2,isqrt2)); CCTK_REAL_VEC n3 = kneg(kmul(eb3,isqrt2)); CCTK_REAL_VEC rm1 = kmul(ec1,isqrt2); CCTK_REAL_VEC rm2 = kmul(ec2,isqrt2); CCTK_REAL_VEC rm3 = kmul(ec3,isqrt2); CCTK_REAL_VEC im1 = kmul(ea1,isqrt2); CCTK_REAL_VEC im2 = kmul(ea2,isqrt2); CCTK_REAL_VEC im3 = kmul(ea3,isqrt2); CCTK_REAL_VEC rmbar1 = kmul(ec1,isqrt2); CCTK_REAL_VEC rmbar2 = kmul(ec2,isqrt2); CCTK_REAL_VEC rmbar3 = kmul(ec3,isqrt2); CCTK_REAL_VEC imbar1 = kneg(kmul(ea1,isqrt2)); CCTK_REAL_VEC imbar2 = kneg(kmul(ea2,isqrt2)); CCTK_REAL_VEC imbar3 = kneg(kmul(ea3,isqrt2)); CCTK_REAL_VEC nn = isqrt2; CCTK_REAL_VEC R1212 = kmul(ToReal(0.5),kadd(JacPDstandard4th12gxy,kadd(JacPDstandard4th21gxy,kmadd(kmadd(gxyL,kmul(gamma122,gamma211),kmadd(gyyL,kmul(gamma211,gamma222),kmadd(gxzL,kmul(gamma122,gamma311),kmadd(gyzL,kmul(gamma222,gamma311),kmadd(gyzL,kmul(gamma211,gamma322),kmadd(gzzL,kmul(gamma311,gamma322),kmul(gamma111,kmadd(gxxL,gamma122,kmadd(gxyL,gamma222,kmul(gxzL,gamma322)))))))))),ToReal(-2.),ksub(ksub(kmadd(gxxL,kmul(kmul(gamma121,gamma121),ToReal(2.)),kmadd(gyyL,kmul(kmul(gamma221,gamma221),ToReal(2.)),kmadd(gzzL,kmul(kmul(gamma321,gamma321),ToReal(2.)),kmul(kmadd(gyzL,kmul(gamma221,gamma321),kmul(gamma121,kmadd(gxyL,gamma221,kmul(gxzL,gamma321)))),ToReal(4.))))),JacPDstandard4th22gxx),JacPDstandard4th11gyy))))); CCTK_REAL_VEC R1213 = kmadd(gamma121,kmadd(gxxL,gamma131,kmadd(gxyL,gamma231,kmul(gxzL,gamma331))),kmul(ToReal(0.5),kadd(JacPDstandard4th12gxz,kadd(JacPDstandard4th31gxy,kmadd(kmadd(gyyL,kmul(gamma211,gamma232),kmadd(gyzL,kmul(gamma232,gamma311),kmadd(gamma132,kmadd(gxyL,gamma211,kmul(gxzL,gamma311)),kmadd(gyzL,kmul(gamma211,gamma332),kmadd(gzzL,kmul(gamma311,gamma332),kmul(gamma111,kmadd(gxxL,gamma132,kmadd(gxyL,gamma232,kmul(gxzL,gamma332))))))))),ToReal(-2.),ksub(kmsub(kmadd(gxyL,kmul(gamma131,gamma221),kmadd(gyyL,kmul(gamma221,gamma231),kmadd(gxzL,kmul(gamma131,gamma321),kmadd(gyzL,kmul(gamma231,gamma321),kmadd(gyzL,kmul(gamma221,gamma331),kmul(gzzL,kmul(gamma321,gamma331))))))),ToReal(2.),JacPDstandard4th23gxx),JacPDstandard4th11gyz)))))); CCTK_REAL_VEC R1223 = kmadd(gamma122,kmadd(gxxL,gamma131,kmadd(gxyL,gamma231,kmul(gxzL,gamma331))),kmul(ToReal(0.5),kadd(JacPDstandard4th22gxz,kadd(JacPDstandard4th31gyy,kmadd(kmadd(gyyL,kmul(gamma221,gamma232),kmadd(gyzL,kmul(gamma232,gamma321),kmadd(gamma132,kmadd(gxyL,gamma221,kmul(gxzL,gamma321)),kmadd(gyzL,kmul(gamma221,gamma332),kmadd(gzzL,kmul(gamma321,gamma332),kmul(gamma121,kmadd(gxxL,gamma132,kmadd(gxyL,gamma232,kmul(gxzL,gamma332))))))))),ToReal(-2.),ksub(kmsub(kmadd(gxyL,kmul(gamma131,gamma222),kmadd(gyyL,kmul(gamma222,gamma231),kmadd(gxzL,kmul(gamma131,gamma322),kmadd(gyzL,kmul(gamma231,gamma322),kmadd(gyzL,kmul(gamma222,gamma331),kmul(gzzL,kmul(gamma322,gamma331))))))),ToReal(2.),JacPDstandard4th23gxy),JacPDstandard4th12gyz)))))); CCTK_REAL_VEC R1313 = kmul(ToReal(0.5),kadd(JacPDstandard4th13gxz,kadd(JacPDstandard4th31gxz,kmadd(kmadd(gxyL,kmul(gamma133,gamma211),kmadd(gyyL,kmul(gamma211,gamma233),kmadd(gxzL,kmul(gamma133,gamma311),kmadd(gyzL,kmul(gamma233,gamma311),kmadd(gyzL,kmul(gamma211,gamma333),kmadd(gzzL,kmul(gamma311,gamma333),kmul(gamma111,kmadd(gxxL,gamma133,kmadd(gxyL,gamma233,kmul(gxzL,gamma333)))))))))),ToReal(-2.),ksub(ksub(kmadd(gxxL,kmul(kmul(gamma131,gamma131),ToReal(2.)),kmadd(gyyL,kmul(kmul(gamma231,gamma231),ToReal(2.)),kmadd(gzzL,kmul(kmul(gamma331,gamma331),ToReal(2.)),kmul(kmadd(gyzL,kmul(gamma231,gamma331),kmul(gamma131,kmadd(gxyL,gamma231,kmul(gxzL,gamma331)))),ToReal(4.))))),JacPDstandard4th33gxx),JacPDstandard4th11gzz))))); CCTK_REAL_VEC R1323 = kmadd(gamma131,kmadd(gxxL,gamma132,kmadd(gxyL,gamma232,kmul(gxzL,gamma332))),kmul(ToReal(0.5),kadd(JacPDstandard4th23gxz,kadd(JacPDstandard4th31gyz,kmadd(kmadd(gyyL,kmul(gamma221,gamma233),kmadd(gyzL,kmul(gamma233,gamma321),kmadd(gamma133,kmadd(gxyL,gamma221,kmul(gxzL,gamma321)),kmadd(gyzL,kmul(gamma221,gamma333),kmadd(gzzL,kmul(gamma321,gamma333),kmul(gamma121,kmadd(gxxL,gamma133,kmadd(gxyL,gamma233,kmul(gxzL,gamma333))))))))),ToReal(-2.),ksub(kmsub(kmadd(gxyL,kmul(gamma132,gamma231),kmadd(gyyL,kmul(gamma231,gamma232),kmadd(gxzL,kmul(gamma132,gamma331),kmadd(gyzL,kmul(gamma232,gamma331),kmadd(gyzL,kmul(gamma231,gamma332),kmul(gzzL,kmul(gamma331,gamma332))))))),ToReal(2.),JacPDstandard4th33gxy),JacPDstandard4th12gzz)))))); CCTK_REAL_VEC R2323 = kmul(ToReal(0.5),kadd(JacPDstandard4th23gyz,kadd(JacPDstandard4th32gyz,kmadd(kmadd(gxyL,kmul(gamma133,gamma222),kmadd(gyyL,kmul(gamma222,gamma233),kmadd(gxzL,kmul(gamma133,gamma322),kmadd(gyzL,kmul(gamma233,gamma322),kmadd(gyzL,kmul(gamma222,gamma333),kmadd(gzzL,kmul(gamma322,gamma333),kmul(gamma122,kmadd(gxxL,gamma133,kmadd(gxyL,gamma233,kmul(gxzL,gamma333)))))))))),ToReal(-2.),ksub(ksub(kmadd(gxxL,kmul(kmul(gamma132,gamma132),ToReal(2.)),kmadd(gyyL,kmul(kmul(gamma232,gamma232),ToReal(2.)),kmadd(gzzL,kmul(kmul(gamma332,gamma332),ToReal(2.)),kmul(kmadd(gyzL,kmul(gamma232,gamma332),kmul(gamma132,kmadd(gxyL,gamma232,kmul(gxzL,gamma332)))),ToReal(4.))))),JacPDstandard4th33gyy),JacPDstandard4th22gzz))))); CCTK_REAL_VEC R4p1212 = kmadd(kxxL,kyyL,knmsub(kxyL,kxyL,R1212)); CCTK_REAL_VEC R4p1213 = kmadd(kxxL,kyzL,knmsub(kxyL,kxzL,R1213)); CCTK_REAL_VEC R4p1223 = kmadd(kxyL,kyzL,knmsub(kxzL,kyyL,R1223)); CCTK_REAL_VEC R4p1313 = kmadd(kxxL,kzzL,knmsub(kxzL,kxzL,R1313)); CCTK_REAL_VEC R4p1323 = kmadd(kxyL,kzzL,knmsub(kxzL,kyzL,R1323)); CCTK_REAL_VEC R4p2323 = kmadd(kyyL,kzzL,knmsub(kyzL,kyzL,R2323)); CCTK_REAL_VEC Ro111 = ToReal(0.); CCTK_REAL_VEC Ro112 = kmadd(kxxL,gamma121,kmadd(kxzL,gamma321,kadd(JacPDstandard4th1kxy,knmsub(kyyL,gamma211,knmsub(kyzL,gamma311,kmsub(kxyL,ksub(gamma221,gamma111),JacPDstandard4th2kxx)))))); CCTK_REAL_VEC Ro113 = kmadd(kxxL,gamma131,kmadd(kxyL,gamma231,kadd(JacPDstandard4th1kxz,knmsub(kyzL,gamma211,knmsub(kzzL,gamma311,kmsub(kxzL,ksub(gamma331,gamma111),JacPDstandard4th3kxx)))))); CCTK_REAL_VEC Ro121 = kmadd(kyyL,gamma211,kmadd(kyzL,gamma311,kadd(JacPDstandard4th2kxx,knmsub(kxxL,gamma121,knmsub(kxzL,gamma321,kmsub(kxyL,ksub(gamma111,gamma221),JacPDstandard4th1kxy)))))); CCTK_REAL_VEC Ro122 = ToReal(0.); CCTK_REAL_VEC Ro123 = kmadd(kxyL,gamma131,kmadd(kyyL,gamma231,kadd(JacPDstandard4th2kxz,knmsub(kxzL,gamma121,knmsub(kzzL,gamma321,kmsub(kyzL,ksub(gamma331,gamma221),JacPDstandard4th3kxy)))))); CCTK_REAL_VEC Ro131 = kmadd(kyzL,gamma211,kmadd(kzzL,gamma311,kadd(JacPDstandard4th3kxx,knmsub(kxxL,gamma131,knmsub(kxyL,gamma231,kmsub(kxzL,ksub(gamma111,gamma331),JacPDstandard4th1kxz)))))); CCTK_REAL_VEC Ro132 = kmadd(kxzL,gamma121,kmadd(kzzL,gamma321,kadd(JacPDstandard4th3kxy,knmsub(kxyL,gamma131,knmsub(kyyL,gamma231,kmsub(kyzL,ksub(gamma221,gamma331),JacPDstandard4th2kxz)))))); CCTK_REAL_VEC Ro133 = ToReal(0.); CCTK_REAL_VEC Ro211 = ToReal(0.); CCTK_REAL_VEC Ro212 = kmadd(kxxL,gamma122,kmadd(kxzL,gamma322,kadd(JacPDstandard4th1kyy,knmsub(kyyL,gamma221,knmsub(kyzL,gamma321,kmsub(kxyL,ksub(gamma222,gamma121),JacPDstandard4th2kxy)))))); CCTK_REAL_VEC Ro213 = kmadd(kxxL,gamma132,kmadd(kxyL,gamma232,kadd(JacPDstandard4th1kyz,knmsub(kyzL,gamma221,knmsub(kzzL,gamma321,kmsub(kxzL,ksub(gamma332,gamma121),JacPDstandard4th3kxy)))))); CCTK_REAL_VEC Ro221 = kmadd(kyyL,gamma221,kmadd(kyzL,gamma321,kadd(JacPDstandard4th2kxy,knmsub(kxxL,gamma122,knmsub(kxzL,gamma322,kmsub(kxyL,ksub(gamma121,gamma222),JacPDstandard4th1kyy)))))); CCTK_REAL_VEC Ro222 = ToReal(0.); CCTK_REAL_VEC Ro223 = kmadd(kxyL,gamma132,kmadd(kyyL,gamma232,kadd(JacPDstandard4th2kyz,knmsub(kxzL,gamma122,knmsub(kzzL,gamma322,kmsub(kyzL,ksub(gamma332,gamma222),JacPDstandard4th3kyy)))))); CCTK_REAL_VEC Ro231 = kmadd(kyzL,gamma221,kmadd(kzzL,gamma321,kadd(JacPDstandard4th3kxy,knmsub(kxxL,gamma132,knmsub(kxyL,gamma232,kmsub(kxzL,ksub(gamma121,gamma332),JacPDstandard4th1kyz)))))); CCTK_REAL_VEC Ro232 = kmadd(kxzL,gamma122,kmadd(kzzL,gamma322,kadd(JacPDstandard4th3kyy,knmsub(kxyL,gamma132,knmsub(kyyL,gamma232,kmsub(kyzL,ksub(gamma222,gamma332),JacPDstandard4th2kyz)))))); CCTK_REAL_VEC Ro233 = ToReal(0.); CCTK_REAL_VEC Ro311 = ToReal(0.); CCTK_REAL_VEC Ro312 = kmadd(kxxL,gamma132,kmadd(kxzL,gamma332,kadd(JacPDstandard4th1kyz,knmsub(kyyL,gamma231,knmsub(kyzL,gamma331,kmsub(kxyL,ksub(gamma232,gamma131),JacPDstandard4th2kxz)))))); CCTK_REAL_VEC Ro313 = kmadd(kxxL,gamma133,kmadd(kxyL,gamma233,kadd(JacPDstandard4th1kzz,knmsub(kyzL,gamma231,knmsub(kzzL,gamma331,kmsub(kxzL,ksub(gamma333,gamma131),JacPDstandard4th3kxz)))))); CCTK_REAL_VEC Ro321 = kmadd(kyyL,gamma231,kmadd(kyzL,gamma331,kadd(JacPDstandard4th2kxz,knmsub(kxxL,gamma132,knmsub(kxzL,gamma332,kmsub(kxyL,ksub(gamma131,gamma232),JacPDstandard4th1kyz)))))); CCTK_REAL_VEC Ro322 = ToReal(0.); CCTK_REAL_VEC Ro323 = kmadd(kxyL,gamma133,kmadd(kyyL,gamma233,kadd(JacPDstandard4th2kzz,knmsub(kxzL,gamma132,knmsub(kzzL,gamma332,kmsub(kyzL,ksub(gamma333,gamma232),JacPDstandard4th3kyz)))))); CCTK_REAL_VEC Ro331 = kmadd(kyzL,gamma231,kmadd(kzzL,gamma331,kadd(JacPDstandard4th3kxz,knmsub(kxxL,gamma133,knmsub(kxyL,gamma233,kmsub(kxzL,ksub(gamma131,gamma333),JacPDstandard4th1kzz)))))); CCTK_REAL_VEC Ro332 = kmadd(kxzL,gamma132,kmadd(kzzL,gamma332,kadd(JacPDstandard4th3kyz,knmsub(kxyL,gamma133,knmsub(kyyL,gamma233,kmsub(kyzL,ksub(gamma232,gamma333),JacPDstandard4th2kzz)))))); CCTK_REAL_VEC Ro333 = ToReal(0.); CCTK_REAL_VEC Rojo11 = kmadd(kadd(gInv23,gInv32),kmadd(kxxL,kyzL,knmsub(kxyL,kxzL,R1213)),kmadd(gInv22,kmadd(kxxL,kyyL,knmsub(kxyL,kxyL,R1212)),kmul(gInv33,kmadd(kxxL,kzzL,knmsub(kxzL,kxzL,R1313))))); CCTK_REAL_VEC Rojo12 = kmadd(gInv23,R1223,knmsub(gInv21,R1212,knmsub(gInv31,R1213,kmadd(gInv12,kmsub(kxyL,kxyL,kmul(kxxL,kyyL)),kmadd(gInv32,kmsub(kxyL,kyzL,kmul(kxzL,kyyL)),kmadd(gInv13,kmsub(kxyL,kxzL,kmul(kxxL,kyzL)),kmul(gInv33,kmadd(kxyL,kzzL,knmsub(kxzL,kyzL,R1323))))))))); CCTK_REAL_VEC Rojo13 = kmadd(kxzL,kmul(kyzL,gInv23),kmadd(gInv13,kmul(kxzL,kxzL),knmsub(gInv21,R1213,knmsub(gInv31,R1313,knmsub(gInv32,R1323,kmadd(gInv12,kmsub(kxyL,kxzL,kmul(kxxL,kyzL)),kmsub(gInv22,kmsub(kxzL,kyyL,kmadd(kxyL,kyzL,R1223)),kmul(kzzL,kmadd(kxyL,gInv23,kmul(kxxL,gInv13)))))))))); CCTK_REAL_VEC Rojo21 = kmadd(gInv32,R1223,knmsub(gInv12,R1212,knmsub(gInv13,R1213,kmadd(gInv21,kmsub(kxyL,kxyL,kmul(kxxL,kyyL)),kmadd(gInv23,kmsub(kxyL,kyzL,kmul(kxzL,kyyL)),kmadd(gInv31,kmsub(kxyL,kxzL,kmul(kxxL,kyzL)),kmul(gInv33,kmadd(kxyL,kzzL,knmsub(kxzL,kyzL,R1323))))))))); CCTK_REAL_VEC Rojo22 = kmadd(kadd(gInv13,gInv31),kmsub(kxzL,kyyL,kmadd(kxyL,kyzL,R1223)),kmadd(gInv11,kmadd(kxxL,kyyL,knmsub(kxyL,kxyL,R1212)),kmul(gInv33,kmadd(kyyL,kzzL,knmsub(kyzL,kyzL,R2323))))); CCTK_REAL_VEC Rojo23 = kmadd(gInv12,R1223,knmsub(gInv31,R1323,knmsub(gInv32,R2323,kmadd(gInv11,kmadd(kxxL,kyzL,knmsub(kxyL,kxzL,R1213)),kmadd(gInv21,kmsub(kxyL,kyzL,kmul(kxzL,kyyL)),kmadd(gInv13,kmsub(kxzL,kyzL,kmul(kxyL,kzzL)),kmul(gInv23,kmsub(kyzL,kyzL,kmul(kyyL,kzzL))))))))); CCTK_REAL_VEC Rojo31 = kmadd(kxzL,kmul(kyzL,gInv32),kmadd(gInv31,kmul(kxzL,kxzL),knmsub(gInv12,R1213,knmsub(gInv13,R1313,knmsub(gInv23,R1323,kmadd(gInv21,kmsub(kxyL,kxzL,kmul(kxxL,kyzL)),kmsub(gInv22,kmsub(kxzL,kyyL,kmadd(kxyL,kyzL,R1223)),kmul(kzzL,kmadd(kxyL,gInv32,kmul(kxxL,gInv31)))))))))); CCTK_REAL_VEC Rojo32 = kmadd(gInv21,R1223,knmsub(gInv13,R1323,knmsub(gInv23,R2323,kmadd(gInv11,kmadd(kxxL,kyzL,knmsub(kxyL,kxzL,R1213)),kmadd(gInv12,kmsub(kxyL,kyzL,kmul(kxzL,kyyL)),kmadd(gInv31,kmsub(kxzL,kyzL,kmul(kxyL,kzzL)),kmul(gInv32,kmsub(kyzL,kyzL,kmul(kyyL,kzzL))))))))); CCTK_REAL_VEC Rojo33 = kmadd(kadd(gInv12,gInv21),kmadd(kxyL,kzzL,knmsub(kxzL,kyzL,R1323)),kmadd(gInv11,kmadd(kxxL,kzzL,knmsub(kxzL,kxzL,R1313)),kmul(gInv22,kmadd(kyyL,kzzL,knmsub(kyzL,kyzL,R2323))))); CCTK_REAL_VEC Psi4rL = kmadd(rmbar1,kmul(rmbar2,kmul(Rojo12,kmul(nn,nn))),kmadd(rmbar1,kmul(rmbar3,kmul(Rojo13,kmul(nn,nn))),kmadd(rmbar1,kmul(rmbar2,kmul(Rojo21,kmul(nn,nn))),kmadd(rmbar2,kmul(rmbar3,kmul(Rojo23,kmul(nn,nn))),kmadd(rmbar1,kmul(rmbar3,kmul(Rojo31,kmul(nn,nn))),kmadd(rmbar2,kmul(rmbar3,kmul(Rojo32,kmul(nn,nn))),kmadd(R4p1212,kmul(kmul(n2,n2),kmul(rmbar1,rmbar1)),kmadd(Rojo11,kmul(kmul(nn,nn),kmul(rmbar1,rmbar1)),kmadd(R4p1212,kmul(kmul(n1,n1),kmul(rmbar2,rmbar2)),kmadd(Rojo22,kmul(kmul(nn,nn),kmul(rmbar2,rmbar2)),kmadd(R4p1313,kmul(kmul(n1,n1),kmul(rmbar3,rmbar3)),kmadd(R4p2323,kmul(kmul(n2,n2),kmul(rmbar3,rmbar3)),kmadd(kmadd(n1,kmul(n2,kmul(R4p1212,kmul(rmbar1,rmbar2))),kmadd(n1,kmul(n3,kmul(R4p1213,kmul(rmbar1,rmbar2))),kmadd(n1,kmul(n3,kmul(R4p1313,kmul(rmbar1,rmbar3))),kmadd(n1,kmul(n3,kmul(R4p1323,kmul(rmbar2,rmbar3))),kmadd(n2,kmul(n3,kmul(R4p2323,kmul(rmbar2,rmbar3))),kmadd(n1,kmul(n2,kmul(R4p1323,kmul(imbar3,imbar3))),kmadd(n1,kmul(nn,kmul(Ro313,kmul(imbar3,imbar3))),kmadd(n3,kmul(nn,kmul(Ro333,kmul(imbar3,imbar3))),kmul(R4p1223,kmul(rmbar1,kmul(rmbar3,kmul(n2,n2)))))))))))),ToReal(-2.),kmadd(n1,kmul(n3,kmul(R4p1223,kmul(kmul(rmbar2,rmbar2),ToReal(-2.)))),kmadd(R4p1313,kmsub(kmul(n3,n3),kmul(rmbar1,rmbar1),kmul(kmul(imbar3,imbar3),kmul(n1,n1))),kmadd(R4p2323,kmsub(kmul(n3,n3),kmul(rmbar2,rmbar2),kmul(kmul(imbar3,imbar3),kmul(n2,n2))),kmadd(Rojo33,kmsub(kmul(nn,nn),kmul(rmbar3,rmbar3),kmul(kmul(imbar3,imbar3),kmul(nn,nn))),kmadd(kmadd(n1,kmul(nn,kmul(rmbar1,kmul(rmbar2,Ro112))),kmadd(n1,kmul(nn,kmul(rmbar1,kmul(rmbar3,Ro113))),kmadd(n2,kmul(nn,kmul(rmbar1,kmul(rmbar3,Ro123))),kmadd(n3,kmul(nn,kmul(rmbar1,kmul(rmbar2,Ro132))),kmadd(n3,kmul(nn,kmul(rmbar1,kmul(rmbar3,Ro133))),kmadd(n1,kmul(nn,kmul(rmbar1,kmul(rmbar2,Ro211))),kmadd(n1,kmul(nn,kmul(rmbar2,kmul(rmbar3,Ro213))),kmadd(n2,kmul(nn,kmul(rmbar2,kmul(rmbar3,Ro223))),kmadd(n3,kmul(nn,kmul(rmbar1,kmul(rmbar2,Ro231))),kmadd(n3,kmul(nn,kmul(rmbar2,kmul(rmbar3,Ro233))),kmadd(n1,kmul(nn,kmul(rmbar1,kmul(rmbar3,Ro311))),kmadd(n1,kmul(nn,kmul(rmbar2,kmul(rmbar3,Ro312))),kmadd(n2,kmul(nn,kmul(rmbar2,kmul(rmbar3,Ro322))),kmadd(n3,kmul(nn,kmul(rmbar1,kmul(rmbar3,Ro331))),kmadd(n3,kmul(nn,kmul(rmbar2,kmul(rmbar3,Ro332))),kmadd(R4p1213,kmul(rmbar2,kmul(rmbar3,kmul(n1,n1))),kmul(R4p1323,kmul(rmbar1,kmul(rmbar2,kmul(n3,n3)))))))))))))))))))),ToReal(2.),kmadd(n1,kmul(nn,kmul(Ro111,kmul(kmul(rmbar1,rmbar1),ToReal(2.)))),kmadd(n2,kmul(nn,kmul(Ro121,kmul(kmul(rmbar1,rmbar1),ToReal(2.)))),kmadd(n3,kmul(nn,kmul(Ro131,kmul(kmul(rmbar1,rmbar1),ToReal(2.)))),kmadd(n1,kmul(nn,kmul(Ro212,kmul(kmul(rmbar2,rmbar2),ToReal(2.)))),kmadd(n2,kmul(nn,kmul(Ro222,kmul(kmul(rmbar2,rmbar2),ToReal(2.)))),kmadd(n3,kmul(nn,kmul(Ro232,kmul(kmul(rmbar2,rmbar2),ToReal(2.)))),kmadd(n1,kmul(n2,kmul(R4p1323,kmul(kmul(rmbar3,rmbar3),ToReal(2.)))),kmadd(n1,kmul(nn,kmul(Ro313,kmul(kmul(rmbar3,rmbar3),ToReal(2.)))),kmadd(n2,kmul(nn,kmul(Ro323,kmul(kmul(rmbar3,rmbar3),ToReal(2.)))),kmadd(n3,kmul(nn,kmul(Ro333,kmul(kmul(rmbar3,rmbar3),ToReal(2.)))),knmsub(kmul(imbar2,imbar2),kmadd(R4p1212,kmul(n1,n1),kmadd(R4p2323,kmul(n3,n3),kmadd(Rojo22,kmul(nn,nn),kmadd(n1,kmul(n3,kmul(R4p1223,ToReal(-2.))),kmadd(n1,kmul(nn,kmul(Ro212,ToReal(2.))),kmadd(n2,kmul(nn,kmul(Ro222,ToReal(2.))),kmul(n3,kmul(nn,kmul(Ro232,ToReal(2.)))))))))),kmadd(n2,kmadd(n1,kmul(R4p1213,kmul(rmbar1,kmul(rmbar3,ToReal(-2.)))),kmadd(n3,kmul(R4p1323,kmul(rmbar1,kmul(rmbar3,ToReal(-2.)))),kmadd(nn,kmul(Ro323,kmul(kmul(imbar3,imbar3),ToReal(-2.))),kmadd(n3,kmul(R4p1223,kmul(rmbar1,kmul(rmbar2,ToReal(2.)))),kmadd(n1,kmul(R4p1223,kmul(rmbar2,kmul(rmbar3,ToReal(2.)))),kmadd(nn,kmul(rmbar1,kmul(rmbar2,kmul(Ro122,ToReal(2.)))),kmadd(nn,kmul(rmbar1,kmul(rmbar2,kmul(Ro221,ToReal(2.)))),kmadd(nn,kmul(rmbar1,kmul(rmbar3,kmul(Ro321,ToReal(2.)))),kmul(n3,kmul(R4p1213,kmul(kmul(rmbar1,rmbar1),ToReal(2.)))))))))))),knmsub(kmul(imbar1,imbar1),kmadd(R4p1212,kmul(n2,n2),kmadd(R4p1313,kmul(n3,n3),kmadd(n2,kmul(kmadd(n3,R4p1213,kmul(nn,Ro121)),ToReal(2.)),kmadd(n3,kmul(nn,kmul(Ro131,ToReal(2.))),kmul(nn,kmadd(nn,Rojo11,kmul(n1,kmul(Ro111,ToReal(2.))))))))),kmsub(imbar1,kmul(imbar3,kmadd(R4p1223,kmul(kmul(n2,n2),ToReal(2.)),kmadd(n1,kmul(kmadd(n2,R4p1213,kmsub(n3,R4p1313,kmul(nn,kadd(Ro113,Ro311)))),ToReal(2.)),kmsub(n2,kmul(kmsub(n3,R4p1323,kmul(nn,kadd(Ro123,Ro321))),ToReal(2.)),kmul(nn,kmadd(nn,kadd(Rojo13,Rojo31),kmul(n3,kmul(kadd(Ro133,Ro331),ToReal(2.))))))))),kmul(imbar2,kmadd(imbar3,kmadd(R4p1213,kmul(kmul(n1,n1),ToReal(2.)),kmadd(n1,kmul(kmadd(n2,R4p1223,kmsub(nn,kadd(Ro213,Ro312),kmul(n3,R4p1323))),ToReal(2.)),kmadd(n2,kmadd(n3,kmul(R4p2323,ToReal(-2.)),kmul(nn,kmul(kadd(Ro223,Ro322),ToReal(2.)))),kmul(nn,kmadd(nn,kadd(Rojo23,Rojo32),kmul(n3,kmul(kadd(Ro233,Ro332),ToReal(2.)))))))),kmul(imbar1,kmadd(Rojo12,kmul(nn,nn),kmadd(Rojo21,kmul(nn,nn),kmadd(n1,kmul(ToReal(-2.),kmadd(n2,R4p1212,kmsub(n3,R4p1213,kmul(nn,kadd(Ro112,Ro211))))),kmadd(n3,kmul(nn,kmul(Ro132,ToReal(2.))),kmadd(n2,kmul(kmadd(n3,R4p1223,kmul(nn,kadd(Ro122,Ro221))),ToReal(2.)),kmadd(n3,kmul(nn,kmul(Ro231,ToReal(2.))),kmul(R4p1323,kmul(kmul(n3,n3),ToReal(2.)))))))))))))))))))))))))))))))))))))))))))); CCTK_REAL_VEC Psi4iL = knmadd(im1,kmadd(kmadd(R4p1313,rm1,kmul(R4p1323,rm2)),kmul(kmul(n3,n3),ToReal(2.)),kmadd(kmul(n2,n2),kmul(kmsub(R4p1212,rm1,kmul(R4p1223,rm3)),ToReal(2.)),kmadd(n2,kmul(ToReal(2.),knmsub(n1,kmadd(R4p1212,rm2,kmul(R4p1213,rm3)),kmadd(n3,kmadd(R4p1223,rm2,kmsub(R4p1213,kmul(rm1,ToReal(2.)),kmul(R4p1323,rm3))),kmul(nn,kmadd(rm2,kadd(Ro122,Ro221),kmadd(rm3,kadd(Ro123,Ro321),kmul(rm1,kmul(Ro121,ToReal(2.))))))))),kmadd(n3,kmul(ToReal(2.),kmsub(nn,kmadd(rm2,kadd(Ro132,Ro231),kmadd(rm3,kadd(Ro133,Ro331),kmul(rm1,kmul(Ro131,ToReal(2.))))),kmul(n1,kmadd(R4p1213,rm2,kmul(R4p1313,rm3))))),kmul(nn,kmadd(n1,kmul(ToReal(2.),kmadd(rm2,kadd(Ro112,Ro211),kmadd(rm3,kadd(Ro113,Ro311),kmul(rm1,kmul(Ro111,ToReal(2.)))))),kmul(nn,kmadd(rm2,kadd(Rojo12,Rojo21),kmadd(rm3,kadd(Rojo13,Rojo31),kmul(rm1,kmul(Rojo11,ToReal(2.)))))))))))),kmadd(im2,kmadd(rm1,kmul(Rojo12,kmul(nn,nn)),kmadd(rm1,kmul(Rojo21,kmul(nn,nn)),kmadd(rm3,kmul(Rojo23,kmul(nn,nn)),kmadd(rm3,kmul(Rojo32,kmul(nn,nn)),kmadd(n3,kmul(nn,kmul(rm1,kmul(Ro132,ToReal(2.)))),kmadd(n3,kmul(nn,kmul(rm1,kmul(Ro231,ToReal(2.)))),kmadd(n3,kmul(nn,kmul(rm3,kmul(Ro233,ToReal(2.)))),kmadd(n3,kmul(nn,kmul(rm3,kmul(Ro332,ToReal(2.)))),kmadd(kmadd(R4p1212,rm2,kmul(R4p1213,rm3)),kmul(kmul(n1,n1),ToReal(2.)),kmadd(R4p1323,kmul(rm1,kmul(kmul(n3,n3),ToReal(2.))),kmadd(R4p2323,kmul(rm2,kmul(kmul(n3,n3),ToReal(2.))),kmadd(rm2,kmul(Rojo22,kmul(kmul(nn,nn),ToReal(2.))),kmadd(n1,kmul(ToReal(2.),kmadd(n2,kmsub(R4p1223,rm3,kmul(R4p1212,rm1)),kmsub(nn,kmadd(rm1,kadd(Ro112,Ro211),kmadd(rm3,kadd(Ro213,Ro312),kmul(rm2,kmul(Ro212,ToReal(2.))))),kmul(n3,kmadd(R4p1213,rm1,kmadd(R4p1323,rm3,kmul(R4p1223,kmul(rm2,ToReal(2.))))))))),kmadd(n2,kmul(ToReal(2.),kmadd(n3,kmsub(R4p1223,rm1,kmul(R4p2323,rm3)),kmul(nn,kmadd(rm1,kadd(Ro122,Ro221),kmadd(rm3,kadd(Ro223,Ro322),kmul(rm2,kmul(Ro222,ToReal(2.)))))))),kmul(n3,kmul(nn,kmul(rm2,kmul(Ro232,ToReal(4.))))))))))))))))))),kmul(im3,kmadd(kmadd(R4p1213,rm2,kmul(R4p1313,rm3)),kmul(kmul(n1,n1),ToReal(2.)),kmadd(kmul(n2,n2),kmadd(R4p1223,kmul(rm1,ToReal(-2.)),kmul(R4p2323,kmul(rm3,ToReal(2.)))),kmadd(n1,kmul(ToReal(2.),knmsub(n3,kmadd(R4p1313,rm1,kmul(R4p1323,rm2)),kmadd(n2,kmadd(R4p1223,rm2,kmsub(R4p1323,kmul(rm3,ToReal(2.)),kmul(R4p1213,rm1))),kmul(nn,kmadd(rm1,kadd(Ro113,Ro311),kmadd(rm2,kadd(Ro213,Ro312),kmul(rm3,kmul(Ro313,ToReal(2.))))))))),kmadd(n2,kmul(ToReal(2.),kmsub(nn,kmadd(rm1,kadd(Ro123,Ro321),kmadd(rm2,kadd(Ro223,Ro322),kmul(rm3,kmul(Ro323,ToReal(2.))))),kmul(n3,kmadd(R4p1323,rm1,kmul(R4p2323,rm2))))),kmul(nn,kmadd(n3,kmul(ToReal(2.),kmadd(rm1,kadd(Ro133,Ro331),kmadd(rm2,kadd(Ro233,Ro332),kmul(rm3,kmul(Ro333,ToReal(2.)))))),kmul(nn,kmadd(rm1,kadd(Rojo13,Rojo31),kmadd(rm2,kadd(Rojo23,Rojo32),kmul(rm3,kmul(Rojo33,ToReal(2.))))))))))))))); /* Copy local copies back to grid functions */ vec_store_partial_prepare(i,lc_imin,lc_imax); vec_store_nta_partial(Psi4i[index],Psi4iL); vec_store_nta_partial(Psi4r[index],Psi4rL); } LC_ENDLOOP3VEC(psi4_calc_4th); } extern "C" void psi4_calc_4th(CCTK_ARGUMENTS) { DECLARE_CCTK_ARGUMENTS; DECLARE_CCTK_PARAMETERS; if (verbose > 1) { CCTK_VInfo(CCTK_THORNSTRING,"Entering psi4_calc_4th_Body"); } if (cctk_iteration % psi4_calc_4th_calc_every != psi4_calc_4th_calc_offset) { return; } const char *const groups[] = { "admbase::curv", "admbase::metric", "grid::coordinates", "WeylScal4::Psi4i_group", "WeylScal4::Psi4r_group"}; GenericFD_AssertGroupStorage(cctkGH, "psi4_calc_4th", 5, groups); switch(fdOrder) { case 2: GenericFD_EnsureStencilFits(cctkGH, "psi4_calc_4th", 2, 2, 2); break; case 4: GenericFD_EnsureStencilFits(cctkGH, "psi4_calc_4th", 2, 2, 2); break; case 6: GenericFD_EnsureStencilFits(cctkGH, "psi4_calc_4th", 2, 2, 2); break; case 8: GenericFD_EnsureStencilFits(cctkGH, "psi4_calc_4th", 2, 2, 2); break; } GenericFD_LoopOverInterior(cctkGH, psi4_calc_4th_Body); if (verbose > 1) { CCTK_VInfo(CCTK_THORNSTRING,"Leaving psi4_calc_4th_Body"); } }