MSBBCXVGD3GRLE5KAI6BKAFRV7SQUWI2SNN43AJAUD3ISRCEXY6QC
PH327KEIT3EWVB6OVULNZDHU4QUMCXBWUDPW6CCBZQJY3JFSPZIAC
722HZ7UFINNE3YKSYKP2NHZ5XEG5QQLQHSKC7PREJZR3EX6RDYUAC
IEFVL4X5BPCLUX5D24ESP7FTFSX2HMCD2GOIFCMSWNUYEHQEJOIAC
KG47IF4CPCUT3BHS34WDRHTH5HYMBTY4OSTB3X7APR2E5ZJ32IYQC
U77PE56ICORZNQW33NXGSEMW7GDHCSSZ4EXB6OHBJSHEG6WHYSSQC
S7LV6VRMBSJPJ765HYMGC5IZ4LCQRJ7VH3D3E36PD6LHVGPQWNYAC
OJZWEAJRLOA5FZAC6FMGBQ6XFQFDNO55IRWLFWDAYWM54MUIU5LAC
QN2UTSQP4IMCMVXZNR24J22N24ASGXF6EWXQ2P3VWVQERUPFAFDQC
2DD222JSYRPHTXKSRXLSOMSCQPZUORNFLLO2P3GMIDELAAMD5MEQC
UUGQGVC4WEKN64WAP7F5QPS2UHGQB5ZLMFRIYNWKMIEBDO3QRX4AC
UZAKARMGORRQG733ZUPJOEGL5FG243I32NCC2SRSFDCZKUQ5A52QC
BVR7DVINVPQG7PA6Z7QYVYNQ43YZL7XCC6AOMSMWMGAAB2Q43STAC
ORAAQXS3UNKXYDJBCTHNEIDJIZDHWMM5EVCZKVMFRKQK2NIQNJGAC
VAF66DTVLDWNG7N2PEQYEH4OH5SPSMFBXKPR2PU67IIM6CVPCJ7AC
PQB3EKQ6MBCXPTW4HB7SGMSTOTYMB3EFZX2573OFCQI6PSOEKCSQC
KQNKYNRSWOY2K7M5PY362RJT4CRUJVECVSCNGKSJJBPYI3NU4GDQC
PUBLIC:
CCTK_REAL regrid_tag TYPE=gf
{
tag
} "Regridding tags"
void CactusAmrMesh::ErrorEst(int lev, TagBoxArray &tags, Real time, int ngrow) {
void CactusAmrMesh::ErrorEst(const int level, TagBoxArray &tags, Real time,
int ngrow) {
// Don't regrid before Cactus is ready to
if (level >= int(ghext->leveldata.size()))
return;
CCTK_VINFO("ErrorEst level %d", level);
return;
const int gi = CCTK_GroupIndex("AMReX::regrid_tag");
assert(gi >= 0);
const int vi = 0;
const int tl = 0;
auto &restrict leveldata = ghext->leveldata.at(level);
auto &restrict groupdata = leveldata.groupdata.at(gi);
MultiFab &mfab = *groupdata.mfab.at(tl);
auto mfitinfo = MFItInfo().SetDynamic(true).EnableTiling({1024000, 16, 32});
#pragma omp parallel
for (MFIter mfi(mfab, mfitinfo); mfi.isValid(); ++mfi) {
const Box &fbx = mfi.fabbox();
const Box &bx = mfi.tilebox();
const Dim3 imin = lbound(bx);
int ash[3], lsh[3];
for (int d = 0; d < dim; ++d)
ash[d] = fbx[orient(d, 1)] - fbx[orient(d, 0)] + 1;
for (int d = 0; d < dim; ++d)
lsh[d] = bx[orient(d, 1)] - bx[orient(d, 0)] + 1;
const Array4<CCTK_REAL> &ctagarr = groupdata.mfab.at(tl)->array(mfi);
const CCTK_REAL *restrict ctags = ctagarr.ptr(imin.x, imin.y, imin.z, vi);
const Array4<char> &atagarr = tags.array(mfi);
char *restrict atags = atagarr.ptr(imin.x, imin.y, imin.z, vi);
constexpr int di = 1;
const int dj = di * ash[0];
const int dk = dj * ash[1];
for (int k = 0; k < lsh[2]; ++k) {
for (int j = 0; j < lsh[1]; ++j) {
#pragma omp simd
for (int i = 0; i < lsh[0]; ++i) {
int idx = di * i + dj * j + dk * k;
atags[idx] = ctags[idx] == 0.0 ? TagBox::CLEAR : TagBox::SET;
}
}
}
}
const int numlevels = maxnumlevels;
ghext->leveldata.resize(numlevels);
for (int level = 0; level < numlevels; ++level) {
GHExt::LevelData &leveldata = ghext->leveldata.at(level);
leveldata.level = level;
// Create coarse grid
const int level = 0;
CCTK_REAL time = 0.0; // dummy time
ghext->amrmesh->MakeNewGrids(time);
SetupLevel(level);
CCTK_REAL time = 0.0; // dummy time
// CCTK_VINFO("BoxArray level %d:", level);
// cout << ghext->amrmesh->boxArray(level) << "\n";
// CCTK_VINFO("DistributionMap level %d:", level);
// cout << ghext->amrmesh->DistributionMap(level) << "\n";
return 0; // unused
} // namespace AMReX
void CreateRefinedGrid(int level) {
CCTK_VINFO("CreateRefinedGrid level %d", level);
// Create grid structure
if (level == 0) {
// Create coarse grid
ghext->amrmesh->MakeNewGrids(time);
} else {
// Create refined grid
int new_finest = -999;
Vector<BoxArray> new_grids;
ghext->amrmesh->MakeNewGrids(0, time, new_finest, new_grids);
// cout << "level=" << level << " new_finest=" << new_finest << "\n";
// assert(new_finest == level);
}
// Create refined grid
CCTK_REAL time = 0.0; // dummy time
int new_finest = -999;
Vector<BoxArray> new_grids;
ghext->amrmesh->MakeNewGrids(0, time, new_finest, new_grids);
cout << "level=" << level << "\n";
cout << "new_finest=" << new_finest << "\n";
assert(new_finest == level - 1 || new_finest == level);
ghext->amrmesh->SetFinestLevel(new_finest);
ghext->amrMesh->SetBoxArray(level, new_grids.at(level));
ghext->amrMesh->SetDistributionMap(...);
const int numgroups = CCTK_NumGroups();
leveldata.groupdata.resize(numgroups);
for (int gi = 0; gi < numgroups; ++gi) {
cGroup group;
int ierr = CCTK_GroupData(gi, &group);
assert(!ierr);
assert(group.grouptype == CCTK_GF);
assert(group.vartype == CCTK_VARIABLE_REAL);
assert(group.disttype == CCTK_DISTRIB_DEFAULT);
assert(group.dim == dim);
void SetupLevel(int level) {
DECLARE_CCTK_PARAMETERS;
CCTK_VINFO("SetupLevel level %d", level);
GHExt::LevelData::GroupData &groupdata = leveldata.groupdata.at(gi);
groupdata.firstvarindex = CCTK_FirstVarIndexI(gi);
groupdata.numvars = group.numvars;
assert(level == int(ghext->leveldata.size()));
ghext->leveldata.resize(level + 1);
GHExt::LevelData &leveldata = ghext->leveldata.at(level);
leveldata.level = level;
// Allocate grid hierarchies
groupdata.mfab.resize(group.numtimelevels);
for (int tl = 0; tl < int(groupdata.mfab.size()); ++tl) {
groupdata.mfab.at(tl) = make_unique<MultiFab>(
ghext->amrmesh->boxArray(leveldata.level),
ghext->amrmesh->DistributionMap(leveldata.level), groupdata.numvars,
ghost_size);
}
}
const int numgroups = CCTK_NumGroups();
leveldata.groupdata.resize(numgroups);
for (int gi = 0; gi < numgroups; ++gi) {
cGroup group;
int ierr = CCTK_GroupData(gi, &group);
assert(!ierr);
assert(group.grouptype == CCTK_GF);
assert(group.vartype == CCTK_VARIABLE_REAL);
assert(group.disttype == CCTK_DISTRIB_DEFAULT);
assert(group.dim == dim);
return 0; // unused
} // namespace AMReX
// Allocate grid hierarchies
groupdata.mfab.resize(group.numtimelevels);
for (int tl = 0; tl < int(groupdata.mfab.size()); ++tl) {
groupdata.mfab.at(tl) = make_unique<MultiFab>(
ghext->amrmesh->boxArray(leveldata.level),
ghext->amrmesh->DistributionMap(leveldata.level), groupdata.numvars,
ghost_size);
}
}
}
const char *recovery_mode = *static_cast<const char *const *>(
CCTK_ParameterGet("recovery_mode", "Cactus", nullptr));
if (!config->recovered || !CCTK_Equals(recovery_mode, "strict")) {
// Set up initial conditions
CCTK_Traverse(cctkGH, "CCTK_INITIAL");
CCTK_Traverse(cctkGH, "CCTK_POSTINITIAL");
CCTK_Traverse(cctkGH, "CCTK_POSTPOSTINITIAL");
CCTK_Traverse(cctkGH, "CCTK_POSTSTEP");
}
const char *recovery_mode = *static_cast<const char *const *>(
CCTK_ParameterGet("recovery_mode", "Cactus", nullptr));
CCTK_Traverse(cctkGH, "CCTK_BASEGRID");
} else {
// Set up initial conditions
for (;;) {
current_level = ghext->amrmesh->finestLevel();
CCTK_Traverse(cctkGH, "CCTK_BASEGRID");
CCTK_Traverse(cctkGH, "CCTK_INITIAL");
CCTK_Traverse(cctkGH, "CCTK_POSTINITIAL");
CCTK_Traverse(cctkGH, "CCTK_POSTPOSTINITIAL");
const int old_numlevels = ghext->amrmesh->finestLevel() + 1;
CreateRefinedGrid(current_level + 1);
const int new_numlevels = ghext->amrmesh->finestLevel() + 1;
assert(new_numlevels == old_numlevels ||
new_numlevels == old_numlevels + 1);
cout << "old_numlevels=" << old_numlevels << "\n";
cout << "new_numlevels=" << new_numlevels << "\n";
cout << "max_numlevels=" << int(ghext->amrmesh->maxLevel() + 1) << "\n";
// Did we create a new level?
if (new_numlevels == old_numlevels)
break;
// Did we create all possible levels?
if (new_numlevels == int(ghext->amrmesh->maxLevel() + 1))
break;
}
current_level = -1;
CCTK_REAL x0[dim], dx[dim];
for (int d = 0; d < dim; ++d) {
x0[d] = CCTK_ORIGIN_SPACE(d);
dx[d] = CCTK_DELTA_SPACE(d);
}
for (int k = 0; k < cctk_lsh[2]; ++k) {
for (int j = 0; j < cctk_lsh[1]; ++j) {
#pragma omp simd
for (int i = 0; i < cctk_lsh[0]; ++i) {
const int idx = CCTK_GFINDEX3D(cctkGH, i, j, k);
CCTK_REAL x = x0[0] + (cctk_lbnd[0] + i) * dx[0];
CCTK_REAL y = x0[1] + (cctk_lbnd[1] + j) * dx[1];
CCTK_REAL z = x0[2] + (cctk_lbnd[2] + k) * dx[2];
CCTK_REAL r = sqrt(sqr(x) + sqr(y) + sqr(z));
tag[idx] = r <= 0.5;
}
}
}
}