include/pops/mesh/boundary/physical_bc.hpp Source FileΒΆ

adc_cpp: include/pops/mesh/boundary/physical_bc.hpp Source File
adc_cpp 0.3.0
Model-free C++23 core for coupled hyperbolic-elliptic systems on adaptive (AMR) meshes, with MPI and GPU (Kokkos) backends
physical_bc.hpp
Go to the documentation of this file.
1
12
13#pragma once
14
15#include <pops/core/foundation/types.hpp>
16#include <pops/mesh/index/box2d.hpp>
17#include <pops/mesh/storage/fab2d.hpp>
18#include <pops/mesh/boundary/fill_boundary.hpp>
19#include <pops/mesh/storage/multifab.hpp>
20
21namespace pops {
22
25enum class BCType { Periodic, Foextrap, Dirichlet };
26
34
35namespace detail {
36// NAMED FUNCTORS (not POPS_HD lambdas) for the physical boundary conditions. Same reasons as the
37// rest of the elliptic/mesh path (#93, recipe #64): fill_physical_bc is called from fill_ghosts,
38// itself pulled from the MG V-cycle first-instantiated from an external TU; an extended lambda there
39// trips up device kernel emission under nvcc. Body identical to the old lambdas (Foextrap: copy of
40// the mirror interior cell; Dirichlet: 2 v - reflection) -> bit-identical.
41// x low face: i = lo - k (k = lo - i), Dirichlet mirror at 2 lo - i - 1.
42// Component range [c0, c1): the WHOLE channel for the all-component fill, or a single component for the
43// per-field aux halo override (ADC-369). The loop bounds are the only difference; the all-channel path
44// (c0=0, c1=ncomp) stays bit-identical.
45struct BCFaceXLoKernel {
46 Array4 a;
47 int c0, c1, lo;
48 bool foe;
49 Real val;
50 POPS_HD void operator()(int i, int j) const {
51 for (int c = c0; c < c1; ++c)
52 a(i, j, c) = foe ? a(lo, j, c) : 2 * val - a(2 * lo - i - 1, j, c);
53 }
54};
55struct BCFaceXHiKernel {
56 Array4 a;
57 int c0, c1, hi;
58 bool foe;
59 Real val;
60 POPS_HD void operator()(int i, int j) const {
61 for (int c = c0; c < c1; ++c)
62 a(i, j, c) = foe ? a(hi, j, c) : 2 * val - a(2 * hi - i + 1, j, c);
63 }
64};
65struct BCFaceYLoKernel {
66 Array4 a;
67 int c0, c1, lo;
68 bool foe;
69 Real val;
70 POPS_HD void operator()(int i, int j) const {
71 for (int c = c0; c < c1; ++c)
72 a(i, j, c) = foe ? a(i, lo, c) : 2 * val - a(i, 2 * lo - j - 1, c);
73 }
74};
75struct BCFaceYHiKernel {
76 Array4 a;
77 int c0, c1, hi;
78 bool foe;
79 Real val;
80 POPS_HD void operator()(int i, int j) const {
81 for (int c = c0; c < c1; ++c)
82 a(i, j, c) = foe ? a(i, hi, c) : 2 * val - a(i, 2 * hi - j + 1, c);
83 }
84};
85} // namespace detail
86
96inline void fill_physical_bc_range(MultiFab& mf, const Box2D& domain, const BCRec& bc, int c0,
97 int c1) {
98 const int ng = mf.n_grow();
99 if (ng == 0)
100 return;
101 // All periodic: fill_boundary has already done everything, nothing to read/write here (and we
102 // avoid a useless barrier on the hot path of the periodic multigrid).
103 if (bc.xlo == BCType::Periodic && bc.xhi == BCType::Periodic && bc.ylo == BCType::Periodic &&
104 bc.yhi == BCType::Periodic)
105 return;
106 // Clamp the component range to the channel (defensive; the per-field override passes a single,
107 // facade-validated component). An empty range is a no-op.
108 if (c0 < 0)
109 c0 = 0;
110 if (c1 > mf.ncomp())
111 c1 = mf.ncomp();
112 if (c0 >= c1)
113 return;
114 // Physical edges on DEVICE (for_each_cell -> kernel): ghost = mirror cell (Foextrap: copy of the
115 // 1st interior; Dirichlet: 2 v - reflection). Ghost index <-> layer: for x low, i = lo-k so the
116 // Dirichlet mirror is 2 lo - i - 1 (k = lo - i). No more device_fence nor host access: these
117 // kernels order after copy_shifted (same execution space), and the y-faces (i EXTENDED for the
118 // corners) order after the x-faces on the same stream.
119 for (int li = 0; li < mf.local_size(); ++li) {
120 Fab2D& F = mf.fab(li);
121 const Box2D v = F.box();
122 Array4 a = F.array();
123
124 // --- x-faces, over the EXTENDED j range (j-ghosts included) ---
125 // We extend the j range to the y/theta GHOSTS (j from v.lo[1]-ng to v.hi[1]+ng) instead of the
126 // VALID range alone. Reason (9-point stencil corner, multi-box): a CROSS term (a_rt/a_tr of the
127 // polar operator) reads the DIAGONAL neighbors p(i+-1, j+-1) -> the CORNER ghost (x-physical
128 // CROSSED with y-ghost) must be filled. When y/theta is PERIODIC or borders a NEIGHBOR box,
129 // fill_boundary has already filled the j-ghost row for the INTERIOR x columns; the x-physical
130 // reflection (which reads a(lo, j) / a(2 lo - i - 1, j) at the SAME j) thus correctly extends the
131 // radial edge into the y halo. Without this extension, the corner (x-ghost, y-ghost) stays at 0
132 // and the cross term is WRONG at the box edge (multi-box divergence, cf. test_polar_schur_multibox).
133 // The VALID range alone was enough in 5-point (no diagonal read); that corner was never read.
134 // NOTE: a DOUBLE-physical corner (x AND y non-periodic) is then OVERWRITTEN by the y pass (i
135 // extended, below, which runs AFTER) -> Cartesian behavior unchanged (y wins). In y-physical we
136 // read a(lo, j-ghost) here, possibly not filled, but the result is overwritten: no effect on the
137 // final corner value. Mono-box theta periodic: only the corners (previously at 0, never read in
138 // 5-point) change -> bit-identical for any stencil <= 9-point including the 5-point Cartesian
139 // (the new corner value is only read by a 9-point).
140 const int jglo = v.lo[1] - ng, jghi = v.hi[1] + ng;
141 if (bc.xlo != BCType::Periodic && v.lo[0] == domain.lo[0]) {
142 const int lo = domain.lo[0];
143 const bool foe = bc.xlo == BCType::Foextrap;
144 const Real val = bc.xlo_val;
145 for_each_cell(Box2D{{lo - ng, jglo}, {lo - 1, jghi}},
146 detail::BCFaceXLoKernel{a, c0, c1, lo, foe, val});
147 }
148 if (bc.xhi != BCType::Periodic && v.hi[0] == domain.hi[0]) {
149 const int hi = domain.hi[0];
150 const bool foe = bc.xhi == BCType::Foextrap;
151 const Real val = bc.xhi_val;
152 for_each_cell(Box2D{{hi + 1, jglo}, {hi + ng, jghi}},
153 detail::BCFaceXHiKernel{a, c0, c1, hi, foe, val});
154 }
155
156 // --- y-faces, over the EXTENDED i range (corners via the already-filled x-ghosts) ---
157 const int iglo = v.lo[0] - ng, ighi = v.hi[0] + ng;
158 if (bc.ylo != BCType::Periodic && v.lo[1] == domain.lo[1]) {
159 const int lo = domain.lo[1];
160 const bool foe = bc.ylo == BCType::Foextrap;
161 const Real val = bc.ylo_val;
162 for_each_cell(Box2D{{iglo, lo - ng}, {ighi, lo - 1}},
163 detail::BCFaceYLoKernel{a, c0, c1, lo, foe, val});
164 }
165 if (bc.yhi != BCType::Periodic && v.hi[1] == domain.hi[1]) {
166 const int hi = domain.hi[1];
167 const bool foe = bc.yhi == BCType::Foextrap;
168 const Real val = bc.yhi_val;
169 for_each_cell(Box2D{{iglo, hi + 1}, {ighi, hi + ng}},
170 detail::BCFaceYHiKernel{a, c0, c1, hi, foe, val});
171 }
172 }
173}
174
176inline void fill_physical_bc(MultiFab& mf, const Box2D& domain, const BCRec& bc) {
177 fill_physical_bc_range(mf, domain, bc, 0, mf.ncomp());
178}
179
185inline void fill_physical_bc(MultiFab& mf, const Box2D& domain, const BCRec& bc, int comp) {
186 fill_physical_bc_range(mf, domain, bc, comp, comp + 1);
187}
188
199
204inline BCRec aux_halo_override(const BCRec& shared, const AuxHaloPolicy& p) {
205 BCRec b = shared;
206 if (b.xlo != BCType::Periodic) {
207 b.xlo = p.type;
208 b.xlo_val = p.value;
209 }
210 if (b.xhi != BCType::Periodic) {
211 b.xhi = p.type;
212 b.xhi_val = p.value;
213 }
214 if (b.ylo != BCType::Periodic) {
215 b.ylo = p.type;
216 b.ylo_val = p.value;
217 }
218 if (b.yhi != BCType::Periodic) {
219 b.yhi = p.type;
220 b.yhi_val = p.value;
221 }
222 return b;
223}
224
227inline void fill_ghosts(MultiFab& mf, const Box2D& domain, const BCRec& bc) {
228 Periodicity per{bc.xlo == BCType::Periodic, bc.ylo == BCType::Periodic};
229 fill_boundary(mf, domain, per);
230 fill_physical_bc(mf, domain, bc);
231}
232
233} // namespace pops
box2d.hpp
Box2D: the integer index space of a 2D cell-centered Cartesian grid.
pops::Fab2D
Single-grid data on a Box2D: VALID box + ng ghost layers, ncomp components, component-slow layout.
Definition fab2d.hpp:59
pops::Fab2D::box
const Box2D & box() const
VALID box (without ghosts).
Definition fab2d.hpp:74
pops::Fab2D::array
Array4 array()
WRITE handle (POD device-copyable) over this Fab. Valid as long as the Fab lives.
Definition fab2d.hpp:91
pops::MultiFab
Field distributed over a level: decomposition (BoxArray) + distribution (DistributionMapping) + ncomp...
Definition multifab.hpp:33
pops::MultiFab::fab
Fab2D & fab(int li)
Local fab at index li (0 <= li < local_size()), for writing.
Definition multifab.hpp:67
pops::MultiFab::ncomp
int ncomp() const
Number of components.
Definition multifab.hpp:60
pops::MultiFab::n_grow
int n_grow() const
Number of ghost layers.
Definition multifab.hpp:62
pops::MultiFab::local_size
int local_size() const
Number of fabs OWNED by this rank (bound on local indices).
Definition multifab.hpp:65
fab2d.hpp
Fab2D: single-grid data on a Box2D (in-house equivalent of AMReX's FArrayBox); Array4 / ConstArray4: ...
fill_boundary.hpp
fill_boundary: INTRA-level halo exchange (fills ghosts from neighbors).
multifab.hpp
MultiFab: a field DISTRIBUTED over a level (equivalent of AMReX's MultiFab).
pops
Definition amr_hierarchy.hpp:29
pops::fill_physical_bc_range
void fill_physical_bc_range(MultiFab &mf, const Box2D &domain, const BCRec &bc, int c0, int c1)
Fills the OUT-OF-domain ghosts of the NON-periodic faces of mf according to bc (Foextrap or Dirichlet...
Definition physical_bc.hpp:96
pops::for_each_cell
void for_each_cell(const Box2D &b, F f)
Applies f to EACH cell (i, j) of box b (bounds inclusive), via Kokkos::parallel_for (Serial / OpenMP ...
Definition for_each.hpp:138
pops::Real
double Real
Definition types.hpp:30
pops::aux_halo_override
BCRec aux_halo_override(const BCRec &shared, const AuxHaloPolicy &p)
Builds the effective override BCRec for a per-field aux halo: starts from the SHARED aux BC shared (s...
Definition physical_bc.hpp:204
pops::fill_boundary
void fill_boundary(MultiFab &mf, const Box2D &domain, Periodicity per={})
BLOCKING halo exchange: begin then end immediately (no overlap).
Definition fill_boundary.hpp:333
pops::fill_physical_bc
void fill_physical_bc(MultiFab &mf, const Box2D &domain, const BCRec &bc)
Fills the physical-face ghosts of ALL components per bc (historical entry point, bit-identical).
Definition physical_bc.hpp:176
pops::BCType
BCType
Boundary condition type for a face: Periodic (handled by fill_boundary), Foextrap (zero gradient,...
Definition physical_bc.hpp:25
pops::fill_ghosts
void fill_ghosts(MultiFab &mf, const Box2D &domain, const BCRec &bc)
COMPLETE ghost filling: fill_boundary (interior + periodic, periodicity deduced from bc) THEN fill_ph...
Definition physical_bc.hpp:227
pops::Array4
WRITE POD handle (raw pointer + strides) over a Fab2D buffer, indexed by (i, j, c) IN GLOBAL INDICES ...
Definition fab2d.hpp:29
pops::AuxHaloPolicy
Per-field aux halo policy (ADC-369): a UNIFORM boundary policy for ONE model-named aux component,...
Definition physical_bc.hpp:195
pops::AuxHaloPolicy::value
Real value
Definition physical_bc.hpp:197
pops::AuxHaloPolicy::type
BCType type
Definition physical_bc.hpp:196
pops::BCRec
Boundary conditions for the FOUR faces of the domain (type + associated Dirichlet value).
Definition physical_bc.hpp:29
pops::BCRec::yhi
BCType yhi
Definition physical_bc.hpp:31
pops::BCRec::ylo_val
Real ylo_val
Definition physical_bc.hpp:32
pops::BCRec::xlo
BCType xlo
Definition physical_bc.hpp:30
pops::BCRec::ylo
BCType ylo
Definition physical_bc.hpp:31
pops::BCRec::xlo_val
Real xlo_val
Definition physical_bc.hpp:32
pops::BCRec::yhi_val
Real yhi_val
Definition physical_bc.hpp:32
pops::BCRec::xhi_val
Real xhi_val
Definition physical_bc.hpp:32
pops::BCRec::xhi
BCType xhi
Definition physical_bc.hpp:30
pops::Box2D
2D integer index space, cell-centered.
Definition box2d.hpp:37
pops::Box2D::hi
int hi[2]
Definition box2d.hpp:39
pops::Box2D::lo
int lo[2]
Definition box2d.hpp:38
pops::Periodicity
Per-direction periodicity: halo wrapping in x and/or y during the exchange (false = open edge,...
Definition fill_boundary.hpp:37
pops::detail::BCFaceXHiKernel
Definition physical_bc.hpp:55
pops::detail::BCFaceXHiKernel::foe
bool foe
Definition physical_bc.hpp:58
pops::detail::BCFaceXHiKernel::a
Array4 a
Definition physical_bc.hpp:56
pops::detail::BCFaceXHiKernel::c1
int c1
Definition physical_bc.hpp:57
pops::detail::BCFaceXHiKernel::operator()
POPS_HD void operator()(int i, int j) const
Definition physical_bc.hpp:60
pops::detail::BCFaceXHiKernel::c0
int c0
Definition physical_bc.hpp:57
pops::detail::BCFaceXHiKernel::hi
int hi
Definition physical_bc.hpp:57
pops::detail::BCFaceXHiKernel::val
Real val
Definition physical_bc.hpp:59
pops::detail::BCFaceXLoKernel
Definition physical_bc.hpp:45
pops::detail::BCFaceXLoKernel::operator()
POPS_HD void operator()(int i, int j) const
Definition physical_bc.hpp:50
pops::detail::BCFaceXLoKernel::c1
int c1
Definition physical_bc.hpp:47
pops::detail::BCFaceXLoKernel::foe
bool foe
Definition physical_bc.hpp:48
pops::detail::BCFaceXLoKernel::lo
int lo
Definition physical_bc.hpp:47
pops::detail::BCFaceXLoKernel::val
Real val
Definition physical_bc.hpp:49
pops::detail::BCFaceXLoKernel::a
Array4 a
Definition physical_bc.hpp:46
pops::detail::BCFaceXLoKernel::c0
int c0
Definition physical_bc.hpp:47
pops::detail::BCFaceYHiKernel
Definition physical_bc.hpp:75
pops::detail::BCFaceYHiKernel::hi
int hi
Definition physical_bc.hpp:77
pops::detail::BCFaceYHiKernel::c0
int c0
Definition physical_bc.hpp:77
pops::detail::BCFaceYHiKernel::operator()
POPS_HD void operator()(int i, int j) const
Definition physical_bc.hpp:80
pops::detail::BCFaceYHiKernel::c1
int c1
Definition physical_bc.hpp:77
pops::detail::BCFaceYHiKernel::a
Array4 a
Definition physical_bc.hpp:76
pops::detail::BCFaceYHiKernel::foe
bool foe
Definition physical_bc.hpp:78
pops::detail::BCFaceYHiKernel::val
Real val
Definition physical_bc.hpp:79
pops::detail::BCFaceYLoKernel
Definition physical_bc.hpp:65
pops::detail::BCFaceYLoKernel::c1
int c1
Definition physical_bc.hpp:67
pops::detail::BCFaceYLoKernel::operator()
POPS_HD void operator()(int i, int j) const
Definition physical_bc.hpp:70
pops::detail::BCFaceYLoKernel::c0
int c0
Definition physical_bc.hpp:67
pops::detail::BCFaceYLoKernel::lo
int lo
Definition physical_bc.hpp:67
pops::detail::BCFaceYLoKernel::foe
bool foe
Definition physical_bc.hpp:68
pops::detail::BCFaceYLoKernel::val
Real val
Definition physical_bc.hpp:69
pops::detail::BCFaceYLoKernel::a
Array4 a
Definition physical_bc.hpp:66
types.hpp
Base scalar types and the POPS_HD macro (host+device portability).
POPS_HD
#define POPS_HD
Definition types.hpp:25
Generated by doxygen 1.9.8