Skip to content

Commit

Permalink
remove MG prolongation operations since they are now user defined
Browse files Browse the repository at this point in the history
  • Loading branch information
@lroberts36
lroberts36 committed Sep 25, 2024
1 parent 4fd9883 commit bca3a1b
Show file tree
Hide file tree
Showing 2 changed files with 3 additions and 120 deletions.
98 changes: 0 additions & 98 deletions src/prolong_restrict/pr_ops.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -288,104 +288,6 @@ using ProlongateSharedMinMod = ProlongateSharedGeneral<true, false>;
using ProlongateSharedLinear = ProlongateSharedGeneral<false, false>;
using ProlongatePiecewiseConstant = ProlongateSharedGeneral<false, true>;

enum class MGProlongationType { Constant, Linear, Quadratic, Kwak };

template <MGProlongationType type>
struct ProlongateSharedMG {
static constexpr bool OperationRequired(TopologicalElement fel,
TopologicalElement cel) {
if (fel != TopologicalElement::CC) return false;
return fel == cel;
}

KOKKOS_FORCEINLINE_FUNCTION
static Real QuadraticFactor(int d) {
if (d == 0) return 1.0; // Indicates this dimension is not included
if (d == 1 || d == -1) return 30.0 / 32.0;
if (d == 3 || d == -3) return 5.0 / 32.0;
if (d == 5 || d == -5) return -3.0 / 32.0;
return 0.0;
}

KOKKOS_FORCEINLINE_FUNCTION
static Real LinearFactor(int d, bool up_bound, bool lo_bound) {
if (d == 0) return 1.0; // Indicates this dimension is not included
if (d == 1) return (2.0 + !up_bound) / 4.0;
if (d == -1) return (2.0 + !lo_bound) / 4.0;
if (d == 3) return !up_bound / 4.0;
if (d == -3) return !lo_bound / 4.0;
return 0.0;
}

KOKKOS_FORCEINLINE_FUNCTION
static Real ConstantFactor(int d) {
if (d == 0) return 1.0; // Indicates this dimension is not included
if (d == 1 || d == -1) return 1.0;
return 0.0;
}

template <int DIM, TopologicalElement el = TopologicalElement::CC,
TopologicalElement /*cel*/ = TopologicalElement::CC>
KOKKOS_FORCEINLINE_FUNCTION static void
Do(const int l, const int m, const int n, const int k, const int j, const int i,
const IndexRange &ckb, const IndexRange &cjb, const IndexRange &cib,
const IndexRange &kb, const IndexRange &jb, const IndexRange &ib,
const Coordinates_t &coords, const Coordinates_t &coarse_coords,
const ParArrayND<Real, VariableState> *pcoarse,
const ParArrayND<Real, VariableState> *pfine) {
using namespace util;
auto &coarse = *pcoarse;
auto &fine = *pfine;

constexpr int element_idx = static_cast<int>(el) % 3;

const int fi = (DIM > 0) ? (i - cib.s) * 2 + ib.s : ib.s;
const int fj = (DIM > 1) ? (j - cjb.s) * 2 + jb.s : jb.s;
const int fk = (DIM > 2) ? (k - ckb.s) * 2 + kb.s : kb.s;

for (int fok = 0; fok < 1 + (DIM > 2); ++fok) {
for (int foj = 0; foj < 1 + (DIM > 1); ++foj) {
for (int foi = 0; foi < 1 + (DIM > 0); ++foi) {
auto &f = fine(element_idx, l, m, n, fk + fok, fj + foj, fi + foi);
f = 0.0;
const bool lo_bound_x = ((fi + foi) == ib.s);
const bool up_bound_x = ((fi + foi) == ib.e);
const bool lo_bound_y = ((fj + foj) == jb.s);
const bool up_bound_y = ((fj + foj) == jb.e);
const bool lo_bound_z = ((fk + fok) == kb.s);
const bool up_bound_z = ((fk + fok) == kb.e);
for (int ok = -(DIM > 2); ok < 1 + (DIM > 2); ++ok) {
for (int oj = -(DIM > 1); oj < 1 + (DIM > 1); ++oj) {
for (int oi = -(DIM > 0); oi < 1 + (DIM > 0); ++oi) {
const int dx = 4 * oi - foi + 1;
const int dy = (DIM > 1) ? 4 * oj - foj + 1 : 0;
const int dz = (DIM > 2) ? 4 * ok - fok + 1 : 0;
if constexpr (MGProlongationType::Linear == type) {
f += LinearFactor(dx, lo_bound_x, up_bound_x) *
LinearFactor(dy, lo_bound_y, up_bound_y) *
LinearFactor(dz, lo_bound_z, up_bound_z) *
coarse(element_idx, l, m, n, k + ok, j + oj, i + oi);
} else if constexpr (MGProlongationType::Kwak == type) {
const Real fac =
((dx <= 1) + (dy <= 1 && DIM > 1) + (dz <= 1 && DIM > 2)) /
(2.0 * DIM);
f += fac * coarse(element_idx, l, m, n, k + ok, j + oj, i + oi);
} else if constexpr (MGProlongationType::Quadratic == type) {
f += QuadraticFactor(dx) * QuadraticFactor(dy) * QuadraticFactor(dz) *
coarse(element_idx, l, m, n, k + ok, j + oj, i + oi);
} else {
f += ConstantFactor(dx) * ConstantFactor(dy) * ConstantFactor(dz) *
coarse(element_idx, l, m, n, k + ok, j + oj, i + oi);
}
}
}
}
}
}
}
}
};

struct ProlongateInternalAverage {
static constexpr bool OperationRequired(TopologicalElement fel,
TopologicalElement cel) {
Expand Down
25 changes: 3 additions & 22 deletions src/solvers/mg_solver.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -101,27 +101,7 @@ class MGSolver {
Metadata::GMGProlongate, Metadata::OneCopy},
shape);

if (params_.prolongation == "Linear") {
mres_err.RegisterRefinementOps<ProlongateSharedMG<MGProlongationType::Linear>,
RestrictAverage>();
} else if (params_.prolongation == "Kwak") {
mres_err.RegisterRefinementOps<ProlongateSharedMG<MGProlongationType::Kwak>,
RestrictAverage>();
} else if (params_.prolongation == "Quadratic") {
mres_err.RegisterRefinementOps<ProlongateSharedMG<MGProlongationType::Quadratic>,
RestrictAverage>();
} else if (params_.prolongation == "Constant") {
mres_err.RegisterRefinementOps<ProlongateSharedMG<MGProlongationType::Constant>,
RestrictAverage>();
} else if (params_.prolongation == "OldLinear") {
mres_err.RegisterRefinementOps<ProlongateSharedLinear, RestrictAverage>();
} else if (params_.prolongation == "User") {
mres_err.RegisterRefinementOps<ProlongateSharedMG<MGProlongationType::Constant>,
RestrictAverage>();
} else {
printf("Requested prolongation type: %s\n", params_.prolongation.c_str());
PARTHENON_FAIL("Unknown multi-grid prolongation type.");
}
mres_err.RegisterRefinementOps<ProlongateSharedLinear, RestrictAverage>();
pkg->AddField(res_err::name(), mres_err);

auto mtemp =
Expand Down Expand Up @@ -508,7 +488,8 @@ class MGSolver {
if (params_.prolongation == "User") {
prolongate = eqs_.template Prolongate<res_err>(tl, set_from_coarser, md_comm);
} else {
prolongate = tl.AddTask(set_from_coarser,
prolongate =
tl.AddTask(set_from_coarser,
BTF(ProlongateBounds<BoundaryType::gmg_prolongate_recv>), md_comm);
}

Expand Down

0 comments on commit bca3a1b

Please sign in to comment.