diff options
Diffstat (limited to 'src/GRHydro_PPMReconstruct_drv_opt.cc')
-rw-r--r-- | src/GRHydro_PPMReconstruct_drv_opt.cc | 124 |
1 files changed, 95 insertions, 29 deletions
diff --git a/src/GRHydro_PPMReconstruct_drv_opt.cc b/src/GRHydro_PPMReconstruct_drv_opt.cc index 43df1b6..9b70d2c 100644 --- a/src/GRHydro_PPMReconstruct_drv_opt.cc +++ b/src/GRHydro_PPMReconstruct_drv_opt.cc @@ -22,14 +22,14 @@ using namespace std; extern "C" CCTK_INT GRHydro_UseGeneralCoordinates(const cGH * cctkGH); -void oned_slice(int npt,int nstart,int nstride,CCTK_REAL* vec3d,vector<CCTK_REAL> vec1d) +static void oned_slice(int npt,int nstart,int nstride,CCTK_REAL* vec3d,vector<CCTK_REAL> vec1d) { int i; for (i=0; i<npt; i++)vec1d[i] = vec3d[nstart+i*nstride]; return; } -void oned_unslice(int npt,int nstart,int nstride,vector<CCTK_REAL> vec1d,CCTK_REAL* vec3d) +static void oned_unslice(int npt,int nstart,int nstride,vector<CCTK_REAL> vec1d,CCTK_REAL* vec3d) { int i; for (i=0; i<npt; i++)vec3d[nstart+i*nstride]=vec1d[i]; @@ -49,38 +49,18 @@ void oned_unslice(int npt,int nstart,int nstride,vector<CCTK_REAL> vec1d,CCTK_RE * with or without divergence cleaning */ -void GRHydro_PPMReconstruct_drv(CCTK_ARGUMENTS) +extern "C" +void CCTK_FNAME(GRHydro_PPMReconstruct_drv_opt)(cGH const * const restrict & restrict cctkGH) { DECLARE_CCTK_ARGUMENTS; DECLARE_CCTK_PARAMETERS; - CCTK_REAL * restrict g11, * restrict g12, *restrict g13; - CCTK_REAL * restrict g22, * restrict g23, *restrict g33; - CCTK_REAL * restrict beta1, * restrict beta2, * restrict beta3; CCTK_REAL * restrict vup; //Multipatch related pointers if(GRHydro_UseGeneralCoordinates(cctkGH)) { - g11=gaa; - g12=gab; - g13=gac; - g22=gbb; - g23=gbc; - g33=gcc; - beta1=betaa; - beta2=betab; - beta3=betac; vup=lvel; } else { - g11=gxx; - g12=gxy; - g13=gxz; - g22=gyy; - g23=gyz; - g33=gzz; - beta1=betax; - beta2=betay; - beta3=betaz; vup=vel; } @@ -159,7 +139,94 @@ void GRHydro_PPMReconstruct_drv(CCTK_ARGUMENTS) vector<CCTK_REAL> Bvcz1d_minus(nmax_xyz); vector<CCTK_REAL> psidc1d_minus(nmax_xyz); - int npt,nstart,nstride; + int nstart; + + typedef + void (*ppm1d_cxx_ptr_t)(const int nx, \ + const double dx, \ + const double* restrict rho, \ + const double* restrict velx, \ + const double* restrict vely, \ + const double* restrict velz, \ + const double* restrict eps, \ + const double* restrict press, \ + const double* restrict temp, \ + const double* restrict ye, \ + const double* restrict Bvcx, \ + const double* restrict Bvcy, \ + const double* restrict Bvcz, \ + const double* restrict psidc, \ + double* restrict rhominus, \ + double* restrict velxminus, \ + double* restrict velyminus, \ + double* restrict velzminus, \ + double* restrict epsminus, \ + double* restrict tempminus, \ + double* restrict yeminus, \ + double* restrict Bvcxminus, \ + double* restrict Bvcyminus, \ + double* restrict Bvczminus, \ + double* restrict psidcminus, \ + double* restrict rhoplus, \ + double* restrict velxplus, \ + double* restrict velyplus, \ + double* restrict velzplus, \ + double* restrict epsplus, \ + double* restrict tempplus, \ + double* restrict yeplus, \ + double* restrict Bvcxplus, \ + double* restrict Bvcyplus, \ + double* restrict Bvczplus, \ + double* restrict psidcplus); + static ppm1d_cxx_ptr_t ppm1d_cxx_funcs[32]; + ppm1d_cxx_ptr_t ppm1d_cxx_func; + + #define MKINDEX(do_temp,do_ye,do_mhd,dc_flag,do_ppm_detect) ((do_temp<<4)|(do_ye<<3)|(do_mhd<<2)|(dc_flag<<1)|(do_ppm_detect)) + #define SETPOINTER(do_temp,do_ye,do_mhd,dc_flag,do_ppm_detect) \ + ppm1d_cxx_funcs[MKINDEX(do_temp,do_ye,do_mhd,dc_flag,do_ppm_detect)] = \ + GRHydro_ppm1d_cxx<do_temp,do_ye,do_mhd,dc_flag,do_ppm_detect> + + // first stuff without MHD + SETPOINTER(false,false,false,false,false); + SETPOINTER(false,false,false,false,true); + // do_temp can only be true of do_ye is also true + SETPOINTER(true,true,false,false,false); + SETPOINTER(true,true,false,false,true); + // but do_ye can be true and do_temp can be false + SETPOINTER(false,true,false,false,false); + SETPOINTER(false,true,false,false,true); + + // with MHD, but dc_flag false + SETPOINTER(false,false,true,false,false); + SETPOINTER(false,false,true,false,true); + // do_temp can only be true of do_ye is also true + SETPOINTER(true,true,true,false,false); + SETPOINTER(true,true,true,false,true); + // but do_ye can be true and do_temp can be false + SETPOINTER(false,true,true,false,false); + SETPOINTER(false,true,true,false,true); + + // with MHD, dc_flag true + SETPOINTER(false,false,true,true,false); + SETPOINTER(false,false,true,true,true); + // do_temp can only be true of do_ye is also true + SETPOINTER(true,true,true,true,false); + SETPOINTER(true,true,true,true,true); + // but do_ye can be true and do_temp can be false + SETPOINTER(false,true,true,true,false); + SETPOINTER(false,true,true,true,true); + + const bool do_temp = *evolve_temper; + const bool do_ye = *evolve_Y_e; + const bool do_mhd = *evolve_MHD; + const bool dc_flag = clean_divergence; + const bool do_ppm_detect = ppm_detect; + ppm1d_cxx_func = ppm1d_cxx_funcs[MKINDEX(do_temp,do_ye,do_mhd,dc_flag,do_ppm_detect)]; + if(ppm1d_cxx_func == NULL) { + CCTK_VError(__LINE__, __FILE__, CCTK_THORNSTRING, + "Combination do_temp = %d do_ye = %d do_mhd = %d dc_flag = %d do_ppm_detect = %d is not supported for PPM. Please instatate the respective template", + do_temp, do_ye, do_mhd, dc_flag, do_ppm_detect); + } //index = i + j*nx + k*nx*ny; @@ -205,8 +272,7 @@ void GRHydro_PPMReconstruct_drv(CCTK_ARGUMENTS) //WE'll need booleans arguments corresponding to the following //ppm1d_cxx<do_temp,do_ye,do_mhd,clean_divergence,ppm_detect> - ppm1d_cxx<true,true,true,true,true> - (nx,CCTK_DELTA_SPACE(0), + ppm1d_cxx_func(nx,CCTK_DELTA_SPACE(0), &rho1d[0],&velx1d[0],&vely1d[0],&velz1d[0],&eps1d[0],&press1d[0],&temp1d[0], &ye1d[0],&Bvcx1d[0],&Bvcy1d[0],&Bvcz1d[0],&psidc1d[0], &rho1d_minus[0],&velx1d_minus[0],&vely1d_minus[0],&velz1d_minus[0],&eps1d_minus[0],&temp1d_minus[0], @@ -296,7 +362,7 @@ void GRHydro_PPMReconstruct_drv(CCTK_ARGUMENTS) //WE'll need booleans arguments corresponding to the following //ppm1d_cxx<do_temp,do_ye,do_mhd,clean_divergence,ppm_detect> - ppm1d_cxx<true,true,true,true,true> + ppm1d_cxx_func (ny,CCTK_DELTA_SPACE(1), &rho1d[0],&vely1d[0],&velz1d[0],&velx1d[0],&eps1d[0],&press1d[0],&temp1d[0], &ye1d[0],&Bvcy1d[0],&Bvcz1d[0],&Bvcx1d[0],&psidc1d[0], @@ -388,7 +454,7 @@ void GRHydro_PPMReconstruct_drv(CCTK_ARGUMENTS) //WE'll need booleans arguments corresponding to the following //ppm1d_cxx<do_temp,do_ye,do_mhd,clean_divergence,ppm_detect> - ppm1d_cxx<true,true,true,true,true> + ppm1d_cxx_func (ny,CCTK_DELTA_SPACE(1), &rho1d[0],&velz1d[0],&velx1d[0],&vely1d[0],&eps1d[0],&press1d[0],&temp1d[0], &ye1d[0],&Bvcz1d[0],&Bvcx1d[0],&Bvcy1d[0],&psidc1d[0], |