E3MBKFT4GEFDAGZQQW4OROY5F6FWC46G6MRH54GDYTGO7O5YSRIAC VIK5E6DBUCP5HGDVHDP6SWTHH7ODEUIGKOGSXR7VLVK46LDD4W6QC BSMJ4V7GV3EOGY4KCSTOJQUOFE2OOCIKQETE4WC2WRNLWBQIBW3QC GSGI6HIWST43XFEORCMKH6VPEGXHEZG5EK4JUUNWP3XFYUKGNA3AC KCIWCVZOHG44WBOLKI2XK33WPHPRI5FWCETF4AOGTPZISKCW3CLQC 722HZ7UFINNE3YKSYKP2NHZ5XEG5QQLQHSKC7PREJZR3EX6RDYUAC BPRNUTY7MHK7LK4EY5MY5OFFG3ABOL7LWXD574L35M74YSQPULFAC BVR7DVINVPQG7PA6Z7QYVYNQ43YZL7XCC6AOMSMWMGAAB2Q43STAC TOBGHRPKEPSXDN56WGSZNWOMCBVJ4KUSLWYWI56MC2RR3MM3KLZAC PQB3EKQ6MBCXPTW4HB7SGMSTOTYMB3EFZX2573OFCQI6PSOEKCSQC 2DD222JSYRPHTXKSRXLSOMSCQPZUORNFLLO2P3GMIDELAAMD5MEQC JD6PQOJ6YYNQYEEWEXO2NM7NVYNBUI6V7ZU6Q3FNHGAT2VYOF5WAC UUGQGVC4WEKN64WAP7F5QPS2UHGQB5ZLMFRIYNWKMIEBDO3QRX4AC WE6MDRN5SPK3THM4COLQFE3IUWBCQ5ZYUIAUCBJAZVEMMOVTNBOAC WASO7G5FJXRXWNH2U2FLUNEKU6VE63OI3HUYP64BVD4LMD6KE7OQC UTHNLH3J3VG7BBJPOKGE6DY7Z2QHDLY72TR2VMEDQYP73SIWZGKAC ULVYXPG2IOJ724MTFVKJLQMPMVL5VAJUOFH7BIXPXVJ2SUGOFRFAC NGD6QWRDHWBMQXL2YXZVF3BFQALSII3LH6BMIOFL7VFMGPLHCZEQC PG2P3JQPLWZRLCAJ7JY2B7G7AM5DHOE2CJUKIPWREI3NAUKZF3TAC YIQN7NJTGEVKW7JZHL6CTH6EPCIXCNBYNURIGXPYZAOUX3VAJQMAC for (auto &restrict leveldata : ghext->leveldata) {
const int min_level = current_level == -1 ? 0 : current_level;const int max_level =current_level == -1 ? ghext->leveldata.size() : current_level + 1;for (int level = min_level; level < max_level; ++level) {auto &restrict leveldata = ghext->leveldata.at(level);
# SCHEDULE HydroToyAMReX_EstimateError AT postinitial# {# LANG: C# WRITES: AMReX::regrid_error(interior)# } "Estimate local error for regridding initial conditions"
SCHEDULE HydroToyAMReX_Pressure AT initial AFTER HydroToyAMReX_Boundaries{LANG: CREADS: conserved(everywhere)WRITES: primitive(everywhere)} "Calculate pressure"SCHEDULE HydroToyAMReX_Fluxes AT initial AFTER HydroToyAMReX_Pressure{LANG: CREADS: conserved(everywhere)READS: primitive(everywhere)WRITES: flux_x(interior) flux_y(interior) flux_z(interior)SYNC: flux_x flux_y flux_z} "Calculate the hydro fluxes"SCHEDULE HydroToyAMReX_EstimateError AT postinitial AFTER HydroToyAMReX_Boundaries{LANG: CREADS: conserved(everywhere)WRITES: AMReX::regrid_error(interior)} "Estimate local error for regridding initial conditions"
# SCHEDULE HydroToyAMReX_EstimateError AT poststep
SCHEDULE HydroToyAMReX_Fluxes AT poststep AFTER HydroToyAMReX_Pressure{LANG: CREADS: conserved(everywhere)READS: primitive(everywhere)WRITES: flux_x(interior) flux_y(interior) flux_z(interior)SYNC: flux_x flux_y flux_z} "Calculate the hydro fluxes"# SCHEDULE HydroToyAMReX_Output AT poststep
# WRITES: AMReX::regrid_error(interior)# } "Estimate local error for regridding during evolution"
# READS: conserved(everywhere)# } "Output grid data"SCHEDULE HydroToyAMReX_EstimateError AT poststep{LANG: CREADS: conserved(everywhere)WRITES: AMReX::regrid_error(interior)} "Estimate local error for regridding during evolution"
# SCHEDULE HydroToyAMReX_Evolve1 AT evol# {# LANG: C# READS: conserved_p(everywhere)# READS: flux_x(everywhere) flux_y(everywhere) flux_z(everywhere)# WRITES: conserved(everywhere)# } "Evolve the hydro equations, step 1"## SCHEDULE HydroToyAMReX_Pressure AT evol AFTER HydroToyAMReX_Evolve1# {# LANG: C# READS: conserved(everywhere)# WRITES: primitive(everywhere)# } "Calculate pressure"## SCHEDULE HydroToyAMReX_Fluxes AT evol AFTER HydroToyAMReX_Pressure# {# LANG: C# READS: conserved(everywhere)# READS: primitive(everywhere)# WRITES: flux_x(interior) flux_y(interior) flux_z(interior)# SYNC: flux_x flux_y flux_z# } "Calculate the hydro fluxes"## SCHEDULE HydroToyAMReX_Evolve AT evol AFTER HydroToyAMReX_Fluxes# {# LANG: C# READS: conserved_p(everywhere)# READS: flux_x(everywhere) flux_y(everywhere) flux_z(everywhere)# WRITES: conserved(everywhere)# } "Evolve the hydro equations, step 2"
# RestrictionSCHEDULE HydroToyAMReX_Boundaries AT postrestrict{LANG: CSYNC: conserved} "Apply boundary conditions after restricting"
fxmomx_(p.I) =mkflux(momx_, [&](auto I) { return momx_(I) * velx_(I) + press_(I); });fxmomy_(p.I) = mkflux(momy_, [&](auto I) { return momy_(I) * velx_(I); });fxmomz_(p.I) = mkflux(momz_, [&](auto I) { return momz_(I) * velx_(I); });
fxmomx_(p.I) = calcflux(momx_, [&](auto I) { return momx_(I) * velx_(I) + press_(I); });fxmomy_(p.I) = calcflux(momy_, [&](auto I) { return momy_(I) * velx_(I); });fxmomz_(p.I) = calcflux(momz_, [&](auto I) { return momz_(I) * velx_(I); });
fymomx_(p.I) = mkflux(momx_, [&](auto I) { return momx_(I) * vely_(I); });fymomy_(p.I) =mkflux(momy_, [&](auto I) { return momy_(I) * vely_(I) + press_(I); });fymomz_(p.I) = mkflux(momz_, [&](auto I) { return momz_(I) * vely_(I); });
fymomx_(p.I) = calcflux(momx_, [&](auto I) { return momx_(I) * vely_(I); });fymomy_(p.I) = calcflux(momy_, [&](auto I) { return momy_(I) * vely_(I) + press_(I); });fymomz_(p.I) = calcflux(momz_, [&](auto I) { return momz_(I) * vely_(I); });
fzmomx_(p.I) = mkflux(momx_, [&](auto I) { return momx_(I) * velz_(I); });fzmomy_(p.I) = mkflux(momy_, [&](auto I) { return momy_(I) * velz_(I); });fzmomz_(p.I) =mkflux(momz_, [&](auto I) { return momz_(I) * velz_(I) + press_(I); });
fzmomx_(p.I) = calcflux(momx_, [&](auto I) { return momx_(I) * velz_(I); });fzmomy_(p.I) = calcflux(momy_, [&](auto I) { return momy_(I) * velz_(I); });fzmomz_(p.I) = calcflux(momz_, [&](auto I) { return momz_(I) * velz_(I) + press_(I); });
auto mkupdate{[&](auto &fx, auto &fy, auto &fz) {return 0.5 * dt *((fx(p.I + p.DI(0)) - fx(p.I)) / dx +(fy(p.I + p.DI(1)) - fy(p.I)) / dy +(fz(p.I + p.DI(2)) - fz(p.I)) / dz);
auto calcupdate{[&](auto &fx, auto &fy, auto &fz) {return dt * ((fx(p.I + p.DI(0)) - fx(p.I)) / dx +(fy(p.I + p.DI(1)) - fy(p.I)) / dy +(fz(p.I + p.DI(2)) - fz(p.I)) / dz);
momx_(p.I) = momx_p_(p.I) - mkupdate(fxmomx_, fymomx_, fzmomx_);momy_(p.I) = momy_p_(p.I) - mkupdate(fxmomy_, fymomy_, fzmomy_);momz_(p.I) = momz_p_(p.I) - mkupdate(fxmomz_, fymomz_, fzmomz_);
momx_(p.I) = momx_p_(p.I) - calcupdate(fxmomx_, fymomx_, fzmomx_);momy_(p.I) = momy_p_(p.I) - calcupdate(fxmomy_, fymomy_, fzmomy_);momz_(p.I) = momz_p_(p.I) - calcupdate(fxmomz_, fymomz_, fzmomz_);
const CCTK_REAL dt = CCTK_DELTA_TIME;const CCTK_REAL dx = CCTK_DELTA_SPACE(0);const CCTK_REAL dy = CCTK_DELTA_SPACE(1);const CCTK_REAL dz = CCTK_DELTA_SPACE(2);
const Loop::GF3D<const CCTK_REAL, 1, 1, 1> rho_(cctkGH, rho);const Loop::GF3D<const CCTK_REAL, 1, 1, 1> momx_(cctkGH, momx);const Loop::GF3D<const CCTK_REAL, 1, 1, 1> momy_(cctkGH, momy);const Loop::GF3D<const CCTK_REAL, 1, 1, 1> momz_(cctkGH, momz);const Loop::GF3D<const CCTK_REAL, 1, 1, 1> etot_(cctkGH, etot);
const Loop::GF3D<const CCTK_REAL, 1, 1, 1> rho_p_(cctkGH, rho_p);const Loop::GF3D<const CCTK_REAL, 1, 1, 1> momx_p_(cctkGH, momx_p);const Loop::GF3D<const CCTK_REAL, 1, 1, 1> momy_p_(cctkGH, momy_p);const Loop::GF3D<const CCTK_REAL, 1, 1, 1> momz_p_(cctkGH, momz_p);const Loop::GF3D<const CCTK_REAL, 1, 1, 1> etot_p_(cctkGH, etot_p);const Loop::GF3D<const CCTK_REAL, 0, 1, 1> fxrho_(cctkGH, fxrho);const Loop::GF3D<const CCTK_REAL, 0, 1, 1> fxmomx_(cctkGH, fxmomx);const Loop::GF3D<const CCTK_REAL, 0, 1, 1> fxmomy_(cctkGH, fxmomy);const Loop::GF3D<const CCTK_REAL, 0, 1, 1> fxmomz_(cctkGH, fxmomz);const Loop::GF3D<const CCTK_REAL, 0, 1, 1> fxetot_(cctkGH, fxetot);const Loop::GF3D<const CCTK_REAL, 1, 0, 1> fyrho_(cctkGH, fyrho);const Loop::GF3D<const CCTK_REAL, 1, 0, 1> fymomx_(cctkGH, fymomx);const Loop::GF3D<const CCTK_REAL, 1, 0, 1> fymomy_(cctkGH, fymomy);const Loop::GF3D<const CCTK_REAL, 1, 0, 1> fymomz_(cctkGH, fymomz);const Loop::GF3D<const CCTK_REAL, 1, 0, 1> fyetot_(cctkGH, fyetot);
const Loop::GF3D<CCTK_REAL, 1, 1, 1> regrid_error_(cctkGH, regrid_error);
const Loop::GF3D<const CCTK_REAL, 1, 1, 0> fzrho_(cctkGH, fzrho);const Loop::GF3D<const CCTK_REAL, 1, 1, 0> fzmomx_(cctkGH, fzmomx);const Loop::GF3D<const CCTK_REAL, 1, 1, 0> fzmomy_(cctkGH, fzmomy);const Loop::GF3D<const CCTK_REAL, 1, 1, 0> fzmomz_(cctkGH, fzmomz);const Loop::GF3D<const CCTK_REAL, 1, 1, 0> fzetot_(cctkGH, fzetot);const Loop::GF3D<CCTK_REAL, 1, 1, 1> rho_(cctkGH, rho);const Loop::GF3D<CCTK_REAL, 1, 1, 1> momx_(cctkGH, momx);const Loop::GF3D<CCTK_REAL, 1, 1, 1> momy_(cctkGH, momy);const Loop::GF3D<CCTK_REAL, 1, 1, 1> momz_(cctkGH, momz);const Loop::GF3D<CCTK_REAL, 1, 1, 1> etot_(cctkGH, etot);// Transport// dt rho + d_i (rho vel^i) = 0// dt mom_j + d_i (mom_j vel^i) = 0// dt etot + d_i (etot vel^i) = 0Loop::loop_all<1, 1, 1>(cctkGH, [&](const Loop::PointDesc &p) {auto mkupdate{[&](auto &fx, auto &fy, auto &fz) {return dt * ((fx(p.I + p.DI(0)) - fx(p.I)) / dx +(fy(p.I + p.DI(1)) - fy(p.I)) / dy +(fz(p.I + p.DI(2)) - fz(p.I)) / dz);
Loop::loop_int<1, 1, 1>(cctkGH, [&](const Loop::PointDesc &p) {auto calcerr{[&](auto &var_) {CCTK_REAL err{0};for (int d = 0; d < dim; ++d) {auto varm = var_(p.I - p.DI(d));auto var0 = var_(p.I);auto varp = var_(p.I + p.DI(d));err = fmax(err, fabs(varm - 2 * var0 + varp));}return err;
rho_(p.I) = rho_p_(p.I) - mkupdate(fxrho_, fyrho_, fzrho_);momx_(p.I) = momx_p_(p.I) - mkupdate(fxmomx_, fymomx_, fzmomx_);momy_(p.I) = momy_p_(p.I) - mkupdate(fxmomy_, fymomy_, fzmomy_);momz_(p.I) = momz_p_(p.I) - mkupdate(fxmomz_, fymomz_, fzmomz_);etot_(p.I) = etot_p_(p.I) - mkupdate(fxetot_, fyetot_, fzetot_);
regrid_error_(p.I) =fmax(fmax(fmax(fmax(calcerr(rho_), calcerr(momx_)), calcerr(momy_)),calcerr(momz_)),calcerr(etot_));