Skip to content

Commit

Permalink
CarpetX: correct advancing batches of active levels in Evolve
Browse files Browse the repository at this point in the history
  • Loading branch information
@lwJi
lwJi committed Jul 25, 2024
1 parent 83bd43f commit 579f63d
Showing 1 changed file with 88 additions and 18 deletions.
106 changes: 88 additions & 18 deletions CarpetX/src/schedule.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -84,6 +84,7 @@ optional<active_levels_t> active_levels;
void Reflux(const cGH *cctkGH, int level);
void Restrict(const cGH *cctkGH, int level, const vector<int> &groups);
void Restrict(const cGH *cctkGH, int level);
void SyncAfterRestrict(const cGH *cctkGH, int level);

namespace {
// Convert a (direction, face) pair to an AMReX Orientation
Expand Down Expand Up @@ -493,6 +494,7 @@ void update_cctkGH(cGH *const cctkGH, const cGH *const sourceGH) {
if (cctkGH == sourceGH)
return;
cctkGH->cctk_iteration = sourceGH->cctk_iteration;
cctkGH->cctk_timefac = sourceGH->cctk_timefac;
cctkGH->cctk_time = sourceGH->cctk_time;
cctkGH->cctk_delta_time = sourceGH->cctk_delta_time;
// for (int d = 0; d < dim; ++d)
Expand Down Expand Up @@ -775,7 +777,7 @@ active_levels_t::active_levels_t(const int min_level, const int max_level,
}
active_levels_t::active_levels_t(const int min_level, const int max_level)
: active_levels_t(min_level, max_level, 0, ghext->num_patches()) {}
active_levels_t::active_levels_t() : active_levels_t(0, ghext->num_levels()) {}
active_levels_t::active_levels_t() : active_levels_t(0, 0) {}

void active_levels_t::assert_consistent_iterations() const {
rat64 good_iteration = -1;
Expand Down Expand Up @@ -1095,7 +1097,7 @@ int Initialise(tFleshConfig *config) {
RecoverGridStructure(cctkGH);

assert(!active_levels);
active_levels = make_optional<active_levels_t>();
active_levels = make_optional<active_levels_t>(0, ghext->num_levels());

CCTK_Traverse(cctkGH, "CCTK_BASEGRID");

Expand Down Expand Up @@ -1344,7 +1346,7 @@ int Initialise(tFleshConfig *config) {
CCTK_VINFO("Initialized %d levels", ghext->num_levels());

assert(!active_levels);
active_levels = make_optional<active_levels_t>();
active_levels = make_optional<active_levels_t>(0, ghext->num_levels());

if (!restrict_during_sync) {
// Restrict
Expand All @@ -1353,6 +1355,11 @@ int Initialise(tFleshConfig *config) {
if (leveldata.level != ghext->num_levels() - 1)
Restrict(cctkGH, leveldata.level);
});
// Prolongation
active_levels->loop_coarse_to_fine([&](const auto &leveldata) {
if (leveldata.level != 0)
SyncAfterRestrict(cctkGH, leveldata.level);
});
CCTK_Traverse(cctkGH, "CCTK_POSTRESTRICT");
}

Expand Down Expand Up @@ -1474,8 +1481,6 @@ void CycleTimelevels(cGH *restrict const cctkGH) {
static Timer timer("CycleTimelevels");
Interval interval(timer);

cctkGH->cctk_iteration += 1;
cctkGH->cctk_time += cctkGH->cctk_delta_time;
update_cctkGHs(cctkGH);

// TODO: Parallelize over groups
Expand Down Expand Up @@ -1713,13 +1718,14 @@ int Evolve(tFleshConfig *config) {
for (const auto &leveldata : patchdata.leveldata)
iteration = min(iteration, leveldata.iteration);

cctkGH->cctk_iteration += 1;

// Loop over all levels, in batches that combine levels that don't
// subcycle. The level range is [min_level, max_level).
int min_level = 0;
while (min_level < ghext->num_levels()) {
for (int min_level = 0, max_level = min_level + 1;
min_level < ghext->num_levels();
min_level = max_level, max_level = min_level + 1) {
// Find end of batch
int max_level = min_level + 1;

while (max_level < ghext->num_levels()) {
bool level_is_subcycling_level = false;
for (const auto &patchdata : ghext->patchdata)
Expand All @@ -1734,18 +1740,32 @@ int Evolve(tFleshConfig *config) {
// Skip this batch of levels if it is not active at the current
// iteration
rat64 level_iteration = -1;
for (const auto &patchdata : ghext->patchdata)
if (min_level < int(patchdata.leveldata.size()))
rat64 level_delta_iteration = -1;
for (const auto &patchdata : ghext->patchdata) {
if (min_level < int(patchdata.leveldata.size())) {
level_iteration = patchdata.leveldata.at(min_level).iteration;
level_delta_iteration =
patchdata.leveldata.at(min_level).delta_iteration;
break;
}
}
assert(level_iteration != -1);
assert(level_delta_iteration != -1);
// TODO: if a break is always ok (eg if subcycling factors are strange
// or if a for() loop with a continue is required.
if (level_iteration > iteration)
break;
continue;

// must not terminate loop iteration due to active_levels being reset at
// bootom of loop body
active_levels = make_optional<active_levels_t>(min_level, max_level);

// Advance iteration number on this batch of levels
level_iteration += level_delta_iteration;
active_levels->loop_serially([&](auto &restrict leveldata) {
leveldata.iteration += leveldata.delta_iteration;
assert(level_iteration == leveldata.iteration);
assert(level_delta_iteration == leveldata.delta_iteration);
});

// We cannot invalidate all non-evolved variables. ODESolvers
Expand All @@ -1755,6 +1775,12 @@ int Evolve(tFleshConfig *config) {

CycleTimelevels(cctkGH);

cctkGH->cctk_timefac = (use_subcycling_wip) ? std::pow(2, min_level) : 1;
cctkGH->cctk_time =
(use_subcycling_wip)
? cctkGH->cctk_delta_time * double(level_iteration)
: cctkGH->cctk_time + cctkGH->cctk_delta_time;

CCTK_Traverse(cctkGH, "CCTK_PRESTEP");
CCTK_Traverse(cctkGH, "CCTK_EVOL");

Expand All @@ -1766,12 +1792,35 @@ int Evolve(tFleshConfig *config) {

if (!restrict_during_sync) {
// Restrict
// TODO: These loop bounds are wrong for subcycling
for (int level = ghext->num_levels() - 2; level >= 0; --level)
Restrict(cctkGH, level);
for (int level = min_level - 1; level >= 0; --level) {
rat64 coarse_iteration =
ghext->patchdata.at(0).leveldata.at(level).iteration;
if (coarse_iteration == level_iteration)
Restrict(cctkGH, level);
else
break;
}
// Prolongation
for (int level = 1; level < max_level; ++level) {
rat64 coarse_iteration =
ghext->patchdata.at(0).leveldata.at(level - 1).iteration;
if (coarse_iteration == level_iteration)
SyncAfterRestrict(cctkGH, level);
else
continue;
}
CCTK_Traverse(cctkGH, "CCTK_POSTRESTRICT");
}

for (int level = min_level - 1; level >= 0; --level) {
rat64 coarse_iteration =
ghext->patchdata.at(0).leveldata.at(level).iteration;
if (coarse_iteration == level_iteration)
min_level = level;
else
break;
}
active_levels = make_optional<active_levels_t>(min_level, max_level);
CCTK_Traverse(cctkGH, "CCTK_POSTSTEP");
CCTK_Traverse(cctkGH, "CCTK_CHECKPOINT");
CCTK_Traverse(cctkGH, "CCTK_ANALYSIS");
Expand All @@ -1781,7 +1830,7 @@ int Evolve(tFleshConfig *config) {
total_evolution_output_time += output_finish_time - output_start_time;

active_levels = optional<active_levels_t>();
} // for min_level
} // while min_level

const double finish_time = gettime();
double num_cells = 0;
Expand All @@ -1794,6 +1843,8 @@ int Evolve(tFleshConfig *config) {
total_evolution_time += iteration_time;
const double iterations_per_second = 1 / iteration_time;
const double cell_updates_per_second = num_cells * iterations_per_second;
const double total_evolution_compute_time =
total_evolution_time - total_evolution_output_time;
CCTK_VINFO("Simulation time: %g "
"Iterations per second: %g "
"Simulation time per second: %g",
Expand All @@ -1806,8 +1857,6 @@ int Evolve(tFleshConfig *config) {
"Grid cell updates per second: %g",
num_cells, cell_updates_per_second);

const double total_evolution_compute_time =
total_evolution_time - total_evolution_output_time;
CCTK_VINFO("Performance:");
CCTK_VINFO(" total evolution time: %g sec",
total_evolution_time);
Expand Down Expand Up @@ -2785,6 +2834,27 @@ void Restrict(const cGH *cctkGH, int level) {
Restrict(cctkGH, level, groups);
}

void SyncAfterRestrict(const cGH *cctkGH, int level) {
const int numgroups = CCTK_NumGroups();
vector<int> groups;
groups.reserve(numgroups);
const auto &patchdata0 = ghext->patchdata.at(0);
const auto &leveldata0 = patchdata0.leveldata.at(0);
for (const auto &groupdataptr : leveldata0.groupdata) {
// Restrict only grid functions
if (groupdataptr) {
auto &restrict groupdata = *groupdataptr;
// Restrict only evolved grid functions
if (groupdata.do_checkpoint)
groups.push_back(groupdata.groupindex);
}
}
auto active_levels_bk = active_levels;
active_levels = make_optional<active_levels_t>(level - 1, level + 1);
SyncGroupsByDirI(cctkGH, groups.size(), groups.data(), nullptr);
active_levels = active_levels_bk;
}

// storage handling
namespace {
int GroupStorageCrease(const cGH *cctkGH, int n_groups, const int *groups,
Expand Down

0 comments on commit 579f63d

Please sign in to comment.