1 #ifndef AMREX_ML_CELL_LINOP_H_
2 #define AMREX_ML_CELL_LINOP_H_
3 #include <AMReX_Config.H>
18 template <
typename MF>
24 using FAB =
typename MF::fab_type;
25 using RT =
typename MF::value_type;
46 void setLevelBC (
int amrlev,
const MF* levelbcdata,
47 const MF* robinbc_a =
nullptr,
48 const MF* robinbc_b =
nullptr,
49 const MF* robinbc_f =
nullptr) final;
63 void updateSolBC (
int amrlev,
const MF& crse_bcdata)
const;
64 void updateCorBC (
int amrlev,
const MF& crse_bcdata)
const;
67 const MLMGBndryT<MF>* bndry=
nullptr,
bool skip_fillboundary=
false)
const;
78 IntVect const& nghost)
const override;
81 const MF& fine_sol,
const MF& fine_rhs)
override;
83 void apply (
int amrlev,
int mglev, MF& out, MF& in,
BCMode bc_mode,
85 void smooth (
int amrlev,
int mglev, MF& sol,
const MF& rhs,
86 bool skip_fillboundary=
false) const final;
89 const MF* crse_bcdata=
nullptr) override;
91 void prepareForFluxes (
int amrlev, const MF* crse_bcdata =
nullptr) override;
94 BCMode bc_mode, const MF* crse_bcdata=
nullptr) final;
97 void reflux (
int crse_amrlev,
98 MF& res, const MF& crse_sol, const MF&,
99 MF&, MF& fine_sol, const MF&) const final;
100 void compFlux (
int amrlev, const
Array<MF*,AMREX_SPACEDIM>& fluxes,
101 MF& sol,
Location loc) const override;
102 void compGrad (
int amrlev, const
Array<MF*,AMREX_SPACEDIM>& grad,
103 MF& sol,
Location loc) const override;
108 MF const& rhs) const override;
114 RT xdoty (
int amrlev,
int mglev, const MF& x, const MF& y,
bool local) const final;
116 virtual
void Fapply (
int amrlev,
int mglev, MF& out, const MF& in) const = 0;
117 virtual
void Fsmooth (
int amrlev,
int mglev, MF& sol, const MF& rhs,
int redblack) const = 0;
119 const
Array<
FAB*,AMREX_SPACEDIM>& flux,
120 const
FAB& sol,
Location loc,
int face_only=0) const = 0;
123 const
Array<MF*,AMREX_SPACEDIM>& ,
127 RT normInf (
int amrlev, MF
const& mf,
bool local)
const override;
131 void avgDownResAmr (
int clev, MF& cres, MF
const& fres)
const override;
140 #ifdef AMREX_USE_HYPRE
170 LinOpBCType crse_fine_bc_type);
179 return bcond[mfi][icomp];
182 return bcloc[mfi][icomp];
219 template <
typename T>
239 template <
typename T>
259 template <
typename T>
264 Array4<int const> mlo;
265 Array4<int const> mhi;
276 Box const& box() const noexcept {
return bx; }
280 template <
typename MF>
287 bctl_dv(bctl.local_size()*ncomp),
292 bcond[mfi].resize(ncomp);
293 bcloc[mfi].resize(ncomp);
299 template <
typename MF>
308 LinOpBCType crse_fine_bc_type)
317 const Box& bx = mfi.validbox();
318 for (
int icomp = 0; icomp < m_ncomp; ++icomp) {
320 BCTuple & bctag = bcond[mfi][icomp];
322 lobc[icomp], hibc[icomp],
323 dx, ratio, interior_bloc,
324 domain_bloc_lo, domain_bloc_hi,
334 for (
int icomp = 0; icomp < m_ncomp; ++icomp) {
336 for (
int m = 0; m < 2*AMREX_SPACEDIM; ++m) {
337 tmp[m].type = bcond[mfi][icomp][m];
338 tmp[m].location = bcloc[mfi][icomp][m];
347 template <
typename MF>
353 template <
typename MF>
366 template <
typename MF>
385 this->
m_dmap[amrlev][mglev],
401 this->
m_dmap[amrlev][mglev],
402 this->
m_geom[amrlev][mglev],
403 face, 0, ngrow, extent, 1,
true);
408 for (
int amrlev = 0; amrlev < this->m_num_amr_levels-1; ++amrlev)
413 this->
m_dmap[amrlev+1][0],
415 this->
m_geom[amrlev+1][0],
417 ratio, amrlev+1, ncomp);
424 #if (AMREX_SPACEDIM != 3)
429 template <
typename MF>
455 const int in_rad = 0;
456 const int out_rad = 1;
457 const int extent_rad = 2;
462 (cba, this->
m_dmap[amrlev][0], in_rad, out_rad, extent_rad, ncomp);
467 const int in_rad = 0;
468 const int out_rad = 1;
469 const int extent_rad = 2;
474 (cba, this->
m_dmap[amrlev][0], in_rad, out_rad, extent_rad, ncomp);
481 m_bndry_cor[amrlev] = std::make_unique<MLMGBndryT<MF>>
483 MF bc_data(this->m_grids[amrlev][0], this->m_dmap[amrlev][0], ncomp, 1);
494 BCType::Dirichlet)}});
504 m_bcondloc[amrlev][mglev] = std::make_unique<BndryCondLoc>(this->
m_grids[amrlev][mglev],
505 this->
m_dmap[amrlev][mglev],
511 template <
typename MF>
514 const MF* robinbc_b,
const MF* robinbc_f)
525 if (a_levelbcdata ==
nullptr) {
531 const MF& bcdata = (a_levelbcdata ==
nullptr) ?
zero : *a_levelbcdata;
543 const int in_rad = 0;
544 const int out_rad = 1;
545 const int extent_rad = 2;
546 const IntVect crse_ratio = br_ref_ratio;
550 (cba, this->
m_dmap[amrlev][0], in_rad, out_rad, extent_rad, ncomp);
561 bcdata, 0, 0, ncomp, br_ref_ratio,
568 m_bndry_sol[amrlev]->setPhysBndryValues(bcdata,0,0,ncomp);
574 m_bndry_sol[amrlev]->setPhysBndryValues(bcdata,0,0,ncomp);
582 const Real* dx = this->
m_geom[amrlev][0].CellSize();
593 AMREX_ASSERT(robinbc_a !=
nullptr && robinbc_b !=
nullptr && robinbc_f !=
nullptr);
597 const Box& domain = this->
m_geom[amrlev][0].Domain();
601 #pragma omp parallel if (Gpu::notInLaunchRegion())
604 Box const& vbx = mfi.validbox();
608 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
611 bool outside_domain_lo = !(domain.
contains(blo));
612 bool outside_domain_hi = !(domain.
contains(bhi));
613 if ((!outside_domain_lo) && (!outside_domain_hi)) {
continue; }
614 for (
int icomp = 0; icomp < ncomp; ++icomp) {
616 if (this->
m_lobc_orig[icomp][idim] == LinOpBCType::Robin && outside_domain_lo)
620 rbc(i,j,k,0) = ra(i,j,k,icomp);
621 rbc(i,j,k,1) = rb(i,j,k,icomp);
622 rbc(i,j,k,2) = rf(i,j,k,icomp);
625 if (this->
m_hibc_orig[icomp][idim] == LinOpBCType::Robin && outside_domain_hi)
629 rbc(i,j,k,0) = ra(i,j,k,icomp);
630 rbc(i,j,k,1) = rb(i,j,k,icomp);
631 rbc(i,j,k,2) = rf(i,j,k,icomp);
640 template <
typename MF>
647 template <
typename MF>
656 this->
m_geom[amrlev-1][0].periodicity());
663 template <
typename MF>
671 this->
m_geom[amrlev-1][0].periodicity());
678 template <
typename MF>
685 BL_ASSERT(mglev == 0 || bc_mode == BCMode::Homogeneous);
686 BL_ASSERT(bndry !=
nullptr || bc_mode == BCMode::Homogeneous);
691 if (!skip_fillboundary) {
692 in.FillBoundary(0, ncomp, this->
m_geom[amrlev][mglev].periodicity(), cross);
695 int flagbc = bc_mode == BCMode::Inhomogeneous;
696 const int imaxorder = this->
maxorder;
698 const Real* dxinv = this->
m_geom[amrlev][mglev].InvCellSize();
699 const RT dxi =
static_cast<RT>(dxinv[0]);
700 const RT dyi = (AMREX_SPACEDIM >= 2) ?
static_cast<RT>(dxinv[1]) :
RT(1.0);
701 const RT dzi = (AMREX_SPACEDIM == 3) ?
static_cast<RT>(dxinv[2]) :
RT(1.0);
703 const auto& maskvals =
m_maskvals[amrlev][mglev];
704 const auto& bcondloc = *
m_bcondloc[amrlev][mglev];
707 const auto& foo = foofab.const_array();
713 "non-cross stencil not support for gpu");
721 tags.reserve(in.local_size()*AMREX_SPACEDIM*ncomp);
724 const Box& vbx = mfi.validbox();
725 const auto& iofab = in.array(mfi);
726 const auto & bdlv = bcondloc.bndryLocs(mfi);
727 const auto & bdcv = bcondloc.bndryConds(mfi);
729 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
730 if (idim != hidden_direction) {
733 const auto& bvlo = (bndry !=
nullptr) ?
735 const auto& bvhi = (bndry !=
nullptr) ?
737 for (
int icomp = 0; icomp < ncomp; ++icomp) {
739 maskvals[olo].const_array(mfi),
740 maskvals[ohi].const_array(mfi),
741 bdlv[icomp][olo], bdlv[icomp][ohi],
743 bdcv[icomp][olo], bdcv[icomp][ohi],
744 vbx.
length(idim), icomp, idim});
755 mllinop_apply_bc_x(0, i, j, k, tag.blen, tag.fab,
756 tag.mask_lo, tag.bctype_lo, tag.bcloc_lo, tag.bcval_lo,
757 imaxorder, dxi, flagbc, tag.comp);
758 mllinop_apply_bc_x(1, i+tag.blen+1, j, k, tag.blen, tag.fab,
759 tag.mask_hi, tag.bctype_hi, tag.bcloc_hi, tag.bcval_hi,
760 imaxorder, dxi, flagbc, tag.comp);
762 #
if (AMREX_SPACEDIM > 1)
764 #
if (AMREX_SPACEDIM > 2)
768 mllinop_apply_bc_y(0, i, j, k, tag.blen, tag.fab,
769 tag.mask_lo, tag.bctype_lo, tag.bcloc_lo, tag.bcval_lo,
770 imaxorder, dyi, flagbc, tag.comp);
771 mllinop_apply_bc_y(1, i, j+tag.blen+1, k, tag.blen, tag.fab,
772 tag.mask_hi, tag.bctype_hi, tag.bcloc_hi, tag.bcval_hi,
773 imaxorder, dyi, flagbc, tag.comp);
775 #
if (AMREX_SPACEDIM > 2)
779 imaxorder, dzi, flagbc, tag.
comp);
782 imaxorder, dzi, flagbc, tag.
comp);
789 if (cross || tensorop)
792 #pragma omp parallel if (Gpu::notInLaunchRegion())
796 const Box& vbx = mfi.validbox();
797 const auto& iofab = in.array(mfi);
799 const auto & bdlv = bcondloc.bndryLocs(mfi);
800 const auto & bdcv = bcondloc.bndryConds(mfi);
802 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim)
804 if (hidden_direction == idim) {
continue; }
809 const int blen = vbx.
length(idim);
810 const auto& mlo = maskvals[olo].array(mfi);
811 const auto& mhi = maskvals[ohi].array(mfi);
812 const auto& bvlo = (bndry !=
nullptr) ? bndry->
bndryValues(olo).const_array(mfi) : foo;
813 const auto& bvhi = (bndry !=
nullptr) ? bndry->
bndryValues(ohi).const_array(mfi) : foo;
814 for (
int icomp = 0; icomp < ncomp; ++icomp) {
815 const BoundCond bctlo = bdcv[icomp][olo];
816 const BoundCond bcthi = bdcv[icomp][ohi];
817 const RT bcllo = bdlv[icomp][olo];
818 const RT bclhi = bdlv[icomp][ohi];
822 imaxorder, dxi, flagbc, icomp);
825 imaxorder, dxi, flagbc, icomp);
826 }
else if (idim == 1) {
829 imaxorder, dyi, flagbc, icomp);
832 imaxorder, dyi, flagbc, icomp);
836 imaxorder, dzi, flagbc, icomp);
839 imaxorder, dzi, flagbc, icomp);
848 amrex::Abort(
"amrex_mllinop_apply_bc not available when BL_NO_FORT=TRUE");
850 if constexpr (std::is_same_v<Real,RT>) {
856 const Box& vbx = mfi.validbox();
858 const auto & bdlv = bcondloc.bndryLocs(mfi);
859 const auto & bdcv = bcondloc.bndryConds(mfi);
872 const auto& fsfab = (bndry !=
nullptr) ? bndry->
bndryValues(ori)[mfi] : foofab;
874 const Mask& m = maskvals[ori][mfi];
877 BL_TO_FORTRAN_ANYD(in[mfi]),
878 BL_TO_FORTRAN_ANYD(m),
880 BL_TO_FORTRAN_ANYD(fsfab),
881 imaxorder, dxinv, flagbc, ncomp, cross);
891 template <
typename MF>
895 const Box& dombx = this->
m_geom[0].back().Domain();
898 const int N = old_ba.
size();
901 for (
int i = 0; i < N; ++i)
904 b.coarsen(grid_size);
909 IntVect big =
b.smallEnd() + grid_size - 1;
912 #if (AMREX_SPACEDIM == 3)
913 for (
int kk = 0; kk < nblks[2]; ++kk) {
915 #if (AMREX_SPACEDIM >= 2)
916 for (
int jj = 0; jj < nblks[1]; ++jj) {
918 for (
int ii = 0; ii < nblks[0]; ++ii)
925 #if (AMREX_SPACEDIM >= 2)
928 #if (AMREX_SPACEDIM == 3)
933 std::sort(bv.begin(), bv.end());
934 bv.erase(std::unique(bv.begin(), bv.end()), bv.end());
941 template <
typename MF>
950 template <
typename MF>
956 Dim3 ratio3 = {2,2,2};
959 ratio3.
y = ratio[1];,
960 ratio3.
z = ratio[2];);
964 auto const& finema =
fine.arrays();
965 auto const& crsema =
crse.const_arrays();
972 finema[box_no](i,j,k,n) += crsema[box_no](ic,jc,kc,n);
979 #pragma omp parallel if (Gpu::notInLaunchRegion())
983 const Box& bx = mfi.tilebox();
991 ffab(i,j,k,n) += cfab(ic,jc,kc,n);
997 template <
typename MF>
1001 const int ncomp = this->
getNComp();
1003 const Geometry& crse_geom = this->
Geom(amrlev,fmglev+1);
1019 cfine.
define(cba,
fine.DistributionMap(), ncomp, ng);
1020 cfine.setVal(
RT(0.0));
1025 bool isEB =
fine.hasEBFabFactory();
1035 #ifdef AMREX_USE_OMP
1036 #pragma omp parallel if (Gpu::notInLaunchRegion())
1040 const Box& bx = mfi.tilebox();
1041 const auto&
ff =
fine.array(mfi);
1042 const auto& cc = cmf->array(mfi);
1047 const auto& flag = (*flags)[mfi];
1054 mlmg_eb_cc_interp_r<2>(tbx,
ff, cc, flg, ncomp);
1065 const bool call_lincc =
true;
1069 #if (AMREX_SPACEDIM == 3)
1090 template <
typename MF>
1095 const int ncomp = this->
getNComp();
1096 const int refratio = this->
AMRRefRatio(famrlev-1);
1105 #ifdef AMREX_USE_OMP
1106 #pragma omp parallel if (Gpu::notInLaunchRegion())
1110 const Box& bx = mfi.tilebox();
1111 auto const&
ff =
fine.array(mfi);
1112 auto const& cc =
crse.const_array(mfi);
1117 const auto& flag = (*flags)[mfi];
1127 mlmg_eb_cc_interp_r<2>(tbx,
ff, cc, flg, ncomp);
1135 mlmg_eb_cc_interp_r<4>(tbx,
ff, cc, flg, ncomp);
1140 amrex::Abort(
"mlmg_eb_cc_interp: only refratio 2 and 4 are supported");
1151 const bool call_lincc =
true;
1173 amrex::Abort(
"mlmg_lin_cc_interp: only refratio 2 and 4 are supported");
1179 template <
typename MF>
1182 const MF& fine_sol,
const MF& fine_rhs)
1185 const int ncomp = this->
getNComp();
1190 template <
typename MF>
1196 applyBC(amrlev, mglev, in, bc_mode, s_mode, bndry);
1197 Fapply(amrlev, mglev, out, in);
1200 template <
typename MF>
1203 bool skip_fillboundary)
const
1206 for (
int redblack = 0; redblack < 2; ++redblack)
1208 applyBC(amrlev, mglev, sol, BCMode::Homogeneous, StateMode::Solution,
1209 nullptr, skip_fillboundary);
1210 Fsmooth(amrlev, mglev, sol, rhs, redblack);
1211 skip_fillboundary =
false;
1215 template <
typename MF>
1218 const MF* crse_bcdata)
1220 BL_PROFILE(
"MLCellLinOp::solutionResidual()");
1221 const int ncomp = this->
getNComp();
1222 if (crse_bcdata !=
nullptr) {
1225 const int mglev = 0;
1226 apply(amrlev, mglev, resid,
x, BCMode::Inhomogeneous, StateMode::Solution,
1233 template <
typename MF>
1237 if (crse_bcdata !=
nullptr) {
1242 template <
typename MF>
1245 BCMode bc_mode,
const MF* crse_bcdata)
1247 BL_PROFILE(
"MLCellLinOp::correctionResidual()");
1248 const int ncomp = this->
getNComp();
1249 if (bc_mode == BCMode::Inhomogeneous)
1256 apply(amrlev, mglev, resid,
x, BCMode::Inhomogeneous, StateMode::Correction,
1262 apply(amrlev, mglev, resid,
x, BCMode::Homogeneous, StateMode::Correction,
nullptr);
1268 template <
typename MF>
1271 MF&, MF& fine_sol,
const MF&)
const
1278 const int ncomp = this->
getNComp();
1280 const int fine_amrlev = crse_amrlev+1;
1282 Real dt = Real(1.0);
1283 const Real* crse_dx = this->
m_geom[crse_amrlev][0].CellSize();
1284 const Real* fine_dx = this->
m_geom[fine_amrlev][0].CellSize();
1286 const int mglev = 0;
1287 applyBC(fine_amrlev, mglev, fine_sol, BCMode::Inhomogeneous, StateMode::Solution,
1293 #ifdef AMREX_USE_OMP
1294 #pragma omp parallel if (Gpu::notInLaunchRegion())
1303 if (fluxreg.CrseHasWork(mfi))
1305 const Box& tbx = mfi.tilebox();
1310 Elixir elify = flux[1].elixir();,
1311 Elixir elifz = flux[2].elixir(););
1312 FFlux(crse_amrlev, mfi, pflux, crse_sol[mfi], Location::FaceCentroid);
1313 fluxreg.CrseAdd(mfi, cpflux, crse_dx, dt,
RunOn::Gpu);
1317 #ifdef AMREX_USE_OMP
1323 if (fluxreg.FineHasWork(mfi))
1325 const Box& tbx = mfi.tilebox();
1326 const int face_only =
true;
1331 Elixir elify = flux[1].elixir();,
1332 Elixir elifz = flux[2].elixir(););
1333 FFlux(fine_amrlev, mfi, pflux, fine_sol[mfi], Location::FaceCentroid, face_only);
1334 fluxreg.FineAdd(mfi, cpflux, fine_dx, dt,
RunOn::Gpu);
1339 fluxreg.Reflux(res);
1342 template <
typename MF>
1349 const int mglev = 0;
1350 const int ncomp = this->
getNComp();
1351 applyBC(amrlev, mglev, sol, BCMode::Inhomogeneous, StateMode::Solution,
1357 #ifdef AMREX_USE_OMP
1358 #pragma omp parallel if (Gpu::notInLaunchRegion())
1365 const Box& tbx = mfi.tilebox();
1370 Elixir elify = flux[1].elixir();,
1371 Elixir elifz = flux[2].elixir(););
1372 FFlux(amrlev, mfi, pflux, sol[mfi], loc);
1373 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
1374 const Box& nbx = mfi.nodaltilebox(idim);
1375 auto const& dst = fluxes[idim]->array(mfi);
1376 auto const& src = pflux[idim]->const_array();
1379 dst(i,j,k,n) = src(i,j,k,n);
1386 template <
typename MF>
1393 if (sol.nComp() > 1) {
1394 amrex::Abort(
"MLCellLinOp::compGrad called, but only works for single-component solves");
1397 const int mglev = 0;
1398 applyBC(amrlev, mglev, sol, BCMode::Inhomogeneous, StateMode::Solution,
1401 const int ncomp = this->
getNComp();
1404 const RT dyi =
static_cast<RT>(this->
m_geom[amrlev][mglev].InvCellSize(1));,
1405 const RT dzi =
static_cast<RT>(this->
m_geom[amrlev][mglev].InvCellSize(2)););
1406 #ifdef AMREX_USE_OMP
1407 #pragma omp parallel if (Gpu::notInLaunchRegion())
1412 const Box& ybx = mfi.nodaltilebox(1);,
1413 const Box& zbx = mfi.nodaltilebox(2););
1414 const auto& s = sol.array(mfi);
1416 const auto& gy = grad[1]->array(mfi);,
1417 const auto& gz = grad[2]->array(mfi););
1421 gx(i,j,k,n) = dxi*(s(i,j,k,n) - s(i-1,j,k,n));
1423 #if (AMREX_SPACEDIM >= 2)
1426 gy(i,j,k,n) = dyi*(s(i,j,k,n) - s(i,j-1,k,n));
1429 #if (AMREX_SPACEDIM == 3)
1432 gz(i,j,k,n) = dzi*(s(i,j,k,n) - s(i,j,k-1,n));
1440 template <
typename MF>
1445 #if (AMREX_SPACEDIM != 3)
1448 const int ncomp = rhs.nComp();
1450 bool cc = rhs.ixType().cellCentered(0);
1454 const RT probxlo =
static_cast<RT>(geom.
ProbLo(0));
1456 #ifdef AMREX_USE_OMP
1457 #pragma omp parallel if (Gpu::notInLaunchRegion())
1461 const Box& tbx = mfi.tilebox();
1462 auto const& rhsarr = rhs.array(mfi);
1463 #if (AMREX_SPACEDIM == 1)
1467 RT rc = probxlo + (i+
RT(0.5))*dx;
1468 rhsarr(i,j,k,n) *= rc*rc;
1473 RT re = probxlo + i*dx;
1474 rhsarr(i,j,k,n) *= re*re;
1477 #elif (AMREX_SPACEDIM == 2)
1481 RT rc = probxlo + (i+
RT(0.5))*dx;
1482 rhsarr(i,j,k,n) *= rc;
1487 RT re = probxlo + i*dx;
1488 rhsarr(i,j,k,n) *= re;
1496 template <
typename MF>
1501 #if (AMREX_SPACEDIM != 3)
1504 const int ncomp = rhs.nComp();
1506 bool cc = rhs.ixType().cellCentered(0);
1510 const RT probxlo =
static_cast<RT>(geom.
ProbLo(0));
1512 #ifdef AMREX_USE_OMP
1513 #pragma omp parallel if (Gpu::notInLaunchRegion())
1517 const Box& tbx = mfi.tilebox();
1518 auto const& rhsarr = rhs.array(mfi);
1519 #if (AMREX_SPACEDIM == 1)
1523 RT rcinv =
RT(1.0)/(probxlo + (i+
RT(0.5))*dx);
1524 rhsarr(i,j,k,n) *= rcinv*rcinv;
1529 RT re = probxlo + i*dx;
1530 RT reinv = (re==
RT(0.0)) ?
RT(0.0) :
RT(1.)/re;
1531 rhsarr(i,j,k,n) *= reinv*reinv;
1534 #elif (AMREX_SPACEDIM == 2)
1538 RT rcinv =
RT(1.0)/(probxlo + (i+
RT(0.5))*dx);
1539 rhsarr(i,j,k,n) *= rcinv;
1544 RT re = probxlo + i*dx;
1545 RT reinv = (re==
RT(0.0)) ?
RT(0.0) :
RT(1.)/re;
1546 rhsarr(i,j,k,n) *= reinv;
1554 template <
typename MF>
1561 const int ncomp = this->
getNComp();
1566 if (factory && !factory->isAllRegular())
1568 if constexpr (std::is_same<MF,MultiFab>()) {
1569 const MultiFab& vfrac = factory->getVolFrac();
1570 for (
int c = 0; c < ncomp; ++c) {
1575 amrex::Abort(
"TODO: MLMG with EB only works with MultiFab");
1581 for (
int c = 0; c < ncomp; ++c) {
1591 template <
typename MF>
1596 const int ncomp = this->
getNComp();
1597 for (
int c = 0; c < ncomp; ++c) {
1598 rhs.plus(-
offset[c], c, 1);
1601 if (!rhs.isAllRegular()) {
1602 if constexpr (std::is_same<MF,MultiFab>()) {
1605 amrex::Abort(
"amrex::EB_set_covered only works with MultiFab");
1611 template <
typename MF>
1615 BL_PROFILE(
"MLCellLinOp::prepareForSolve()");
1617 const int imaxorder = this->
maxorder;
1618 const int ncomp = this->
getNComp();
1624 const auto& bcondloc = *
m_bcondloc[amrlev][mglev];
1625 const auto& maskvals =
m_maskvals[amrlev][mglev];
1627 const RT dxi =
static_cast<RT>(this->
m_geom[amrlev][mglev].InvCellSize(0));
1628 const RT dyi =
static_cast<RT>((AMREX_SPACEDIM >= 2) ? this->
m_geom[amrlev][mglev].InvCellSize(1) : Real(1.0));
1629 const RT dzi =
static_cast<RT>((AMREX_SPACEDIM == 3) ? this->
m_geom[amrlev][mglev].InvCellSize(2) : Real(1.0));
1631 auto& undrrelxr = this->
m_undrrelxr[amrlev][mglev];
1632 MF foo(this->
m_grids[amrlev][mglev], this->
m_dmap[amrlev][mglev], ncomp, 0,
MFInfo().SetAlloc(
false));
1637 (factory) ? &(factory->getMultiEBCellFlagFab()) :
nullptr;
1638 auto area = (factory) ? factory->getAreaFrac()
1643 #ifdef AMREX_USE_GPU
1646 if (factory && !factory->isAllRegular()) {
1647 if constexpr (!std::is_same<MF,MultiFab>()) {
1648 amrex::Abort(
"MLCellLinOp with EB only works with MultiFab");
1651 tags.reserve(foo.local_size()*AMREX_SPACEDIM*ncomp);
1655 const Box& vbx = mfi.validbox();
1657 const auto & bdlv = bcondloc.bndryLocs(mfi);
1658 const auto & bdcv = bcondloc.bndryConds(mfi);
1662 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim)
1669 for (
int icomp = 0; icomp < ncomp; ++icomp) {
1670 tags.emplace_back(MLMGPSEBTag<RT>{undrrelxr[olo].array(mfi),
1671 undrrelxr[ohi].array(mfi),
1673 maskvals[olo].const_array(mfi),
1674 maskvals[ohi].const_array(mfi),
1675 bdlv[icomp][olo], bdlv[icomp][ohi],
1677 bdcv[icomp][olo], bdcv[icomp][ohi],
1678 vbx.
length(idim), icomp, idim});
1685 [=]
AMREX_GPU_DEVICE (
int i,
int j,
int k, MLMGPSEBTag<RT>
const& tag) noexcept
1690 mllinop_comp_interp_coef0_x_eb
1691 (0, i , j, k, tag.blen, tag.flo, tag.mlo, tag.ap,
1692 tag.bctlo, tag.bcllo, imaxorder, dxi, tag.comp);
1693 mllinop_comp_interp_coef0_x_eb
1694 (1, i+tag.blen+1, j, k, tag.blen, tag.fhi, tag.mhi, tag.ap,
1695 tag.bcthi, tag.bclhi, imaxorder, dxi, tag.comp);
1697 #
if (AMREX_SPACEDIM > 1)
1699 #
if (AMREX_SPACEDIM > 2)
1703 mllinop_comp_interp_coef0_y_eb
1704 (0, i, j , k, tag.blen, tag.flo, tag.mlo, tag.ap,
1705 tag.bctlo, tag.bcllo, imaxorder, dyi, tag.comp);
1706 mllinop_comp_interp_coef0_y_eb
1707 (1, i, j+tag.blen+1, k, tag.blen, tag.fhi, tag.mhi, tag.ap,
1708 tag.bcthi, tag.bclhi, imaxorder, dyi, tag.comp);
1710 #
if (AMREX_SPACEDIM > 2)
1712 mllinop_comp_interp_coef0_z_eb
1713 (0, i, j, k , tag.blen, tag.flo, tag.mlo, tag.ap,
1714 tag.bctlo, tag.bcllo, imaxorder, dzi, tag.comp);
1715 mllinop_comp_interp_coef0_z_eb
1716 (1, i, j, k+tag.blen+1, tag.blen, tag.fhi, tag.mhi, tag.ap,
1717 tag.bcthi, tag.bclhi, imaxorder, dzi, tag.comp);
1725 (0, i , j, k, tag.blen, tag.flo, tag.mlo,
1726 tag.bctlo, tag.bcllo, imaxorder, dxi, tag.comp);
1728 (1, i+tag.blen+1, j, k, tag.blen, tag.fhi, tag.mhi,
1729 tag.bcthi, tag.bclhi, imaxorder, dxi, tag.comp);
1731 #if (AMREX_SPACEDIM > 1)
1733 #if (AMREX_SPACEDIM > 2)
1738 (0, i, j , k, tag.blen, tag.flo, tag.mlo,
1739 tag.bctlo, tag.bcllo, imaxorder, dyi, tag.comp);
1741 (1, i, j+tag.blen+1, k, tag.blen, tag.fhi, tag.mhi,
1742 tag.bcthi, tag.bclhi, imaxorder, dyi, tag.comp);
1744 #if (AMREX_SPACEDIM > 2)
1747 (0, i, j, k , tag.blen, tag.flo, tag.mlo,
1748 tag.bctlo, tag.bcllo, imaxorder, dzi, tag.comp);
1750 (1, i, j, k+tag.blen+1, tag.blen, tag.fhi, tag.mhi,
1751 tag.bcthi, tag.bclhi, imaxorder, dzi, tag.comp);
1762 tags.reserve(foo.local_size()*AMREX_SPACEDIM*ncomp);
1766 const Box& vbx = mfi.validbox();
1768 const auto & bdlv = bcondloc.bndryLocs(mfi);
1769 const auto & bdcv = bcondloc.bndryConds(mfi);
1771 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim)
1773 if (idim != hidden_direction) {
1776 for (
int icomp = 0; icomp < ncomp; ++icomp) {
1778 undrrelxr[ohi].array(mfi),
1779 maskvals[olo].const_array(mfi),
1780 maskvals[ohi].const_array(mfi),
1781 bdlv[icomp][olo], bdlv[icomp][ohi],
1783 bdcv[icomp][olo], bdcv[icomp][ohi],
1784 vbx.
length(idim), icomp, idim});
1795 mllinop_comp_interp_coef0_x
1796 (0, i , j, k, tag.blen, tag.flo, tag.mlo,
1797 tag.bctlo, tag.bcllo, imaxorder, dxi, tag.comp);
1798 mllinop_comp_interp_coef0_x
1799 (1, i+tag.blen+1, j, k, tag.blen, tag.fhi, tag.mhi,
1800 tag.bcthi, tag.bclhi, imaxorder, dxi, tag.comp);
1802 #
if (AMREX_SPACEDIM > 1)
1804 #
if (AMREX_SPACEDIM > 2)
1808 mllinop_comp_interp_coef0_y
1809 (0, i, j , k, tag.blen, tag.flo, tag.mlo,
1810 tag.bctlo, tag.bcllo, imaxorder, dyi, tag.comp);
1811 mllinop_comp_interp_coef0_y
1812 (1, i, j+tag.blen+1, k, tag.blen, tag.fhi, tag.mhi,
1813 tag.bcthi, tag.bclhi, imaxorder, dyi, tag.comp);
1815 #
if (AMREX_SPACEDIM > 2)
1831 #ifdef AMREX_USE_OMP
1832 #pragma omp parallel
1836 const Box& vbx = mfi.validbox();
1838 const auto & bdlv = bcondloc.bndryLocs(mfi);
1839 const auto & bdcv = bcondloc.bndryConds(mfi);
1844 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim)
1846 if (idim == hidden_direction) {
continue; }
1851 const int blen = vbx.
length(idim);
1852 const auto& mlo = maskvals[olo].array(mfi);
1853 const auto& mhi = maskvals[ohi].array(mfi);
1854 const auto& flo = undrrelxr[olo].array(mfi);
1855 const auto& fhi = undrrelxr[ohi].array(mfi);
1856 for (
int icomp = 0; icomp < ncomp; ++icomp) {
1857 const BoundCond bctlo = bdcv[icomp][olo];
1858 const BoundCond bcthi = bdcv[icomp][ohi];
1859 const auto bcllo = bdlv[icomp][olo];
1860 const auto bclhi = bdlv[icomp][ohi];
1863 if constexpr (!std::is_same<MF,MultiFab>()) {
1864 amrex::Abort(
"MLCellLinOp with EB only works with MultiFab");
1866 auto const& ap = area[idim]->const_array(mfi);
1868 mllinop_comp_interp_coef0_x_eb
1869 (0, blo, blen, flo, mlo, ap, bctlo, bcllo,
1870 imaxorder, dxi, icomp);
1871 mllinop_comp_interp_coef0_x_eb
1872 (1, bhi, blen, fhi, mhi, ap, bcthi, bclhi,
1873 imaxorder, dxi, icomp);
1874 }
else if (idim == 1) {
1875 mllinop_comp_interp_coef0_y_eb
1876 (0, blo, blen, flo, mlo, ap, bctlo, bcllo,
1877 imaxorder, dyi, icomp);
1878 mllinop_comp_interp_coef0_y_eb
1879 (1, bhi, blen, fhi, mhi, ap, bcthi, bclhi,
1880 imaxorder, dyi, icomp);
1882 mllinop_comp_interp_coef0_z_eb
1883 (0, blo, blen, flo, mlo, ap, bctlo, bcllo,
1884 imaxorder, dzi, icomp);
1885 mllinop_comp_interp_coef0_z_eb
1886 (1, bhi, blen, fhi, mhi, ap, bcthi, bclhi,
1887 imaxorder, dzi, icomp);
1895 (0, blo, blen, flo, mlo, bctlo, bcllo,
1896 imaxorder, dxi, icomp);
1898 (1, bhi, blen, fhi, mhi, bcthi, bclhi,
1899 imaxorder, dxi, icomp);
1900 }
else if (idim == 1) {
1902 (0, blo, blen, flo, mlo, bctlo, bcllo,
1903 imaxorder, dyi, icomp);
1905 (1, bhi, blen, fhi, mhi, bcthi, bclhi,
1906 imaxorder, dyi, icomp);
1909 (0, blo, blen, flo, mlo, bctlo, bcllo,
1910 imaxorder, dzi, icomp);
1912 (1, bhi, blen, fhi, mhi, bcthi, bclhi,
1913 imaxorder, dzi, icomp);
1924 template <
typename MF>
1929 const int ncomp = this->getNComp();
1938 template <
typename MF>
1942 if (!m_volinv.empty()) {
return; }
1944 m_volinv.resize(this->m_num_amr_levels);
1945 for (
int amrlev = 0; amrlev < this->m_num_amr_levels; ++amrlev) {
1946 m_volinv[amrlev].resize(this->NMGLevels(amrlev));
1949 AMREX_ASSERT(this->m_coarse_fine_bc_type == LinOpBCType::Dirichlet || ! this->hasHiddenDimension());
1953 auto f = [&] (
int amrlev,
int mglev) {
1955 const auto *factory =
dynamic_cast<EBFArrayBoxFactory const*
>(this->Factory(amrlev,mglev));
1956 if (factory && !factory->isAllRegular())
1958 if constexpr (std::is_same<MF,MultiFab>()) {
1960 m_volinv[amrlev][mglev] = vfrac.
sum(0,
true);
1962 amrex::Abort(
"MLCellLinOp with EB only works with MultiFab");
1968 if (this->m_coarse_fine_bc_type == LinOpBCType::Dirichlet) {
1969 m_volinv[amrlev][mglev]
1970 =
RT(1.0 / this->compactify(this->Geom(amrlev,mglev).Domain()).d_numPts());
1972 m_volinv[amrlev][mglev]
1973 =
RT(1.0 / this->m_grids[amrlev][mglev].d_numPts());
1981 int mgbottom = this->NMGLevels(0)-1;
1987 if (factory && !factory->isAllRegular())
1989 ParallelAllReduce::Sum<RT>({m_volinv[0][0], m_volinv[0][mgbottom]},
1991 temp1 =
RT(1.0)/m_volinv[0][0];
1992 temp2 =
RT(1.0)/m_volinv[0][mgbottom];
1996 temp1 = m_volinv[0][0];
1997 temp2 = m_volinv[0][mgbottom];
1999 m_volinv[0][0] = temp1;
2000 m_volinv[0][mgbottom] = temp2;
2004 template <
typename MF>
2008 const int ncomp = this->getNComp();
2009 const int finest_level = this->NAMRLevels() - 1;
2012 if (! mf.isAllRegular()) {
2013 if constexpr (!std::is_same<MF,MultiFab>()) {
2014 amrex::Abort(
"MLCellLinOpT with EB only works with MultiFab");
2016 const auto *factory =
dynamic_cast<EBFArrayBoxFactory const*
>(this->Factory(amrlev));
2017 const MultiFab& vfrac = factory->getVolFrac();
2018 if (amrlev == finest_level) {
2019 #ifdef AMREX_USE_GPU
2028 return std::abs(ma[box_no](i,j,k,n)
2029 *vfrac_ma[box_no](i,j,k));
2034 #ifdef AMREX_USE_OMP
2035 #pragma omp parallel reduction(max:norm)
2038 Box const& bx = mfi.tilebox();
2039 auto const& fab = mf.const_array(mfi);
2048 #ifdef AMREX_USE_GPU
2050 auto const& ma = mf.const_arrays();
2051 auto const& mask_ma = m_norm_fine_mask[amrlev]->const_arrays();
2058 if (mask_ma[box_no](i,j,k)) {
2059 return std::abs(ma[box_no](i,j,k,n)
2060 *vfrac_ma[box_no](i,j,k));
2068 #ifdef AMREX_USE_OMP
2069 #pragma omp parallel reduction(max:norm)
2072 Box const& bx = mfi.tilebox();
2073 auto const& fab = mf.const_array(mfi);
2074 auto const&
mask = m_norm_fine_mask[amrlev]->const_array(mfi);
2089 if (amrlev == finest_level) {
2092 norm = mf.norminf(*m_norm_fine_mask[amrlev], 0, ncomp,
IntVect(0),
true);
2100 template <
typename MF>
2104 int ncomp = this->getNComp();
2105 for (
int falev = this->NAMRLevels()-1; falev > 0; --falev)
2108 if (!sol[falev].isAllRegular()) {
2109 if constexpr (std::is_same<MF,MultiFab>()) {
2112 amrex::Abort(
"EB_average_down only works with MultiFab");
2122 template <
typename MF>
2127 if (!fres.isAllRegular()) {
2128 if constexpr (std::is_same<MF,MultiFab>()) {
2130 this->AMRRefRatio(clev));
2132 amrex::Abort(
"EB_average_down only works with MultiFab");
2138 this->AMRRefRatio(clev));
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define AMREX_ALWAYS_ASSERT_WITH_MESSAGE(EX, MSG)
Definition: AMReX_BLassert.H:49
#define BL_ASSERT(EX)
Definition: AMReX_BLassert.H:39
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition: AMReX_BLassert.H:50
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
#define AMREX_LAUNCH_HOST_DEVICE_LAMBDA(...)
Definition: AMReX_GpuLaunch.nolint.H:16
#define AMREX_HOST_DEVICE_PARALLEL_FOR_3D(...)
Definition: AMReX_GpuLaunch.nolint.H:54
#define AMREX_HOST_DEVICE_PARALLEL_FOR_4D(...)
Definition: AMReX_GpuLaunch.nolint.H:55
#define AMREX_GPU_DEVICE
Definition: AMReX_GpuQualifiers.H:18
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
Array4< int const > offset
Definition: AMReX_HypreMLABecLap.cpp:1089
Box cbx
Definition: AMReX_HypreMLABecLap.cpp:1091
Array4< Real > fine
Definition: AMReX_InterpFaceRegister.cpp:90
Array4< int const > mask
Definition: AMReX_InterpFaceRegister.cpp:93
Array4< Real const > crse
Definition: AMReX_InterpFaceRegister.cpp:92
#define AMREX_LOOP_4D(bx, ncomp, i, j, k, n, block)
Definition: AMReX_Loop.nolint.H:16
void amrex_mllinop_apply_bc(const int *lo, const int *hi, amrex_real *phi, const int *philo, const int *phihi, const int *mask, const int *mlo, const int *mhi, int cdir, int bct, amrex_real bcl, const amrex_real *bcval, const int *blo, const int *bhi, int maxorder, const amrex_real *dxinv, int inhomog, int nc, int cross)
#define AMREX_D_TERM(a, b, c)
Definition: AMReX_SPACE.H:129
#define AMREX_D_DECL(a, b, c)
Definition: AMReX_SPACE.H:104
if(!(yy_init))
Definition: amrex_iparser.lex.nolint.H:935
const FabSetT< MF > & bndryValues(Orientation face) const noexcept
Return FabSet on given face.
Definition: AMReX_BndryData.H:77
Maintain an identifier for boundary condition types.
Definition: AMReX_BoundCond.H:20
A collection of Boxes stored in an Array.
Definition: AMReX_BoxArray.H:530
void define(const Box &bx)
Initialize the BoxArray from a single box. It is an error if the BoxArray has already been initialize...
Definition: AMReX_BoxArray.cpp:342
BoxArray & coarsen(int refinement_ratio)
Coarsen each Box in the BoxArray to the specified ratio.
Definition: AMReX_BoxArray.cpp:662
Long size() const noexcept
Return the number of boxes in the BoxArray.
Definition: AMReX_BoxArray.H:577
A class for managing a List of Boxes that share a common IndexType. This class implements operations ...
Definition: AMReX_BoxList.H:52
static AMREX_GPU_HOST_DEVICE BoxND TheUnitBox() noexcept
This static member function returns a constant reference to an object of type BoxND representing the ...
Definition: AMReX_Box.H:737
AMREX_GPU_HOST_DEVICE IntVectND< dim > length() const noexcept
Return the length of the BoxND.
Definition: AMReX_Box.H:146
AMREX_GPU_HOST_DEVICE bool contains(const IntVectND< dim > &p) const noexcept
Returns true if argument is contained within BoxND.
Definition: AMReX_Box.H:204
const Real * CellSize() const noexcept
Returns the cellsize for each coordinate direction.
Definition: AMReX_CoordSys.H:71
Calculates the distribution of FABs to MPI processes.
Definition: AMReX_DistributionMapping.H:41
Definition: AMReX_EBFabFactory.H:22
const MultiFab & getVolFrac() const noexcept
Definition: AMReX_EBFabFactory.H:53
An Array of FortranArrayBox(FAB)-like Objects.
Definition: AMReX_FabArray.H:343
Array4< typename FabArray< FAB >::value_type const > const_array(const MFIter &mfi) const noexcept
Definition: AMReX_FabArray.H:1592
MultiArray4< typename FabArray< FAB >::value_type const > const_arrays() const noexcept
Definition: AMReX_FabArray.H:1674
Definition: AMReX_FabFactory.H:50
Rectangular problem domain geometry.
Definition: AMReX_Geometry.H:73
const Real * ProbLo() const noexcept
Returns the lo end of the problem domain in each dimension.
Definition: AMReX_Geometry.H:178
GpuArray< int, AMREX_SPACEDIM > isPeriodicArray() const noexcept
Definition: AMReX_Geometry.H:347
Periodicity periodicity() const noexcept
Definition: AMReX_Geometry.H:355
const Box & Domain() const noexcept
Returns our rectangular domain.
Definition: AMReX_Geometry.H:210
Definition: AMReX_Tuple.H:93
Definition: AMReX_GpuElixir.H:13
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool allGT(const IntVectND< dim > &rhs) const noexcept
Returns true if this is greater than argument for all components. NOTE: This is NOT a strict weak ord...
Definition: AMReX_IntVect.H:418
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE IntVectND< dim > TheCellVector() noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition: AMReX_IntVect.H:711
An InterpBndryData object adds to a BndryData object the ability to manipulate and set the data store...
Definition: AMReX_InterpBndryData.H:40
a one-thingy-per-box distributed object
Definition: AMReX_LayoutData.H:13
Definition: AMReX_MFIter.H:57
bool isValid() const noexcept
Is the iterator valid i.e. is it associated with a FAB?
Definition: AMReX_MFIter.H:141
Definition: AMReX_MLCellLinOp.H:160
GpuArray< BCTL, 2 *AMREX_SPACEDIM > const * getBCTLPtr(const MFIter &mfi) const noexcept
Definition: AMReX_MLCellLinOp.H:184
const Vector< RealTuple > & bndryLocs(const MFIter &mfi) const noexcept
Definition: AMReX_MLCellLinOp.H:175
Gpu::DeviceVector< GpuArray< BCTL, 2 *AMREX_SPACEDIM > > bctl_dv
Definition: AMReX_MLCellLinOp.H:191
BndryCondLoc(const BoxArray &ba, const DistributionMapping &dm, int ncomp)
Definition: AMReX_MLCellLinOp.H:281
LayoutData< GpuArray< BCTL, 2 *AMREX_SPACEDIM > * > bctl
Definition: AMReX_MLCellLinOp.H:190
const BCTuple & bndryConds(const MFIter &mfi, int icomp) const noexcept
Definition: AMReX_MLCellLinOp.H:178
int m_ncomp
Definition: AMReX_MLCellLinOp.H:192
const RealTuple & bndryLocs(const MFIter &mfi, int icomp) const noexcept
Definition: AMReX_MLCellLinOp.H:181
void setLOBndryConds(const Geometry &geom, const Real *dx, const Vector< Array< BCType, AMREX_SPACEDIM > > &lobc, const Vector< Array< BCType, AMREX_SPACEDIM > > &hibc, IntVect const &ratio, const RealVect &interior_bloc, const Array< Real, AMREX_SPACEDIM > &domain_bloc_lo, const Array< Real, AMREX_SPACEDIM > &domain_bloc_hi, LinOpBCType crse_fine_bc_type)
Definition: AMReX_MLCellLinOp.H:302
LayoutData< Vector< RealTuple > > bcloc
Definition: AMReX_MLCellLinOp.H:189
const Vector< BCTuple > & bndryConds(const MFIter &mfi) const noexcept
Definition: AMReX_MLCellLinOp.H:172
LayoutData< Vector< BCTuple > > bcond
Definition: AMReX_MLCellLinOp.H:188
Definition: AMReX_MLCellLinOp.H:21
void smooth(int amrlev, int mglev, MF &sol, const MF &rhs, bool skip_fillboundary=false) const final
Smooth.
Definition: AMReX_MLCellLinOp.H:1202
virtual void Fsmooth(int amrlev, int mglev, MF &sol, const MF &rhs, int redblack) const =0
Vector< RT > getSolvabilityOffset(int amrlev, int mglev, MF const &rhs) const override
get offset for fixing solvability
Definition: AMReX_MLCellLinOp.H:1556
typename MF::fab_type FAB
Definition: AMReX_MLCellLinOp.H:24
void averageDownSolutionRHS(int camrlev, MF &crse_sol, MF &crse_rhs, const MF &fine_sol, const MF &fine_rhs) override
Average-down data from fine AMR level to coarse AMR level.
Definition: AMReX_MLCellLinOp.H:1181
void averageDownAndSync(Vector< MF > &sol) const override
Definition: AMReX_MLCellLinOp.H:2102
BoxArray makeNGrids(int grid_size) const
Definition: AMReX_MLCellLinOp.H:893
Vector< YAFluxRegisterT< MF > > m_fluxreg
Definition: AMReX_MLCellLinOp.H:204
virtual void applyBC(int amrlev, int mglev, MF &in, BCMode bc_mode, StateMode s_mode, const MLMGBndryT< MF > *bndry=nullptr, bool skip_fillboundary=false) const
Definition: AMReX_MLCellLinOp.H:680
void updateSolBC(int amrlev, const MF &crse_bcdata) const
Definition: AMReX_MLCellLinOp.H:649
Vector< std::unique_ptr< MLMGBndryT< MF > > > m_bndry_sol
Definition: AMReX_MLCellLinOp.H:148
void compGrad(int amrlev, const Array< MF *, AMREX_SPACEDIM > &grad, MF &sol, Location loc) const override
Compute gradients of the solution.
Definition: AMReX_MLCellLinOp.H:1388
void avgDownResAmr(int clev, MF &cres, MF const &fres) const override
Definition: AMReX_MLCellLinOp.H:2124
void reflux(int crse_amrlev, MF &res, const MF &crse_sol, const MF &, MF &, MF &fine_sol, const MF &) const final
Reflux at AMR coarse/fine boundary.
Definition: AMReX_MLCellLinOp.H:1270
virtual void Fapply(int amrlev, int mglev, MF &out, const MF &in) const =0
void defineBC()
Definition: AMReX_MLCellLinOp.H:431
void update() override
Update for reuse.
Definition: AMReX_MLCellLinOp.H:642
typename MF::value_type RT
Definition: AMReX_MLCellLinOp.H:25
LinOpBCType BCType
Definition: AMReX_MLCellLinOp.H:27
void correctionResidual(int amrlev, int mglev, MF &resid, MF &x, const MF &b, BCMode bc_mode, const MF *crse_bcdata=nullptr) final
Compute residual for the residual-correction form, resid = b - L(x)
Definition: AMReX_MLCellLinOp.H:1244
MLCellLinOpT(const MLCellLinOpT< MF > &)=delete
typename MLLinOpT< MF >::BCMode BCMode
Definition: AMReX_MLCellLinOp.H:28
RT normInf(int amrlev, MF const &mf, bool local) const override
Definition: AMReX_MLCellLinOp.H:2006
Array< RT, 2 *BL_SPACEDIM > RealTuple
Definition: AMReX_MLCellLinOp.H:157
virtual bool isCrossStencil() const
Definition: AMReX_MLCellLinOp.H:60
void updateCorBC(int amrlev, const MF &crse_bcdata) const
Definition: AMReX_MLCellLinOp.H:665
Array< BoundCond, 2 *BL_SPACEDIM > BCTuple
Definition: AMReX_MLCellLinOp.H:158
void interpolation(int amrlev, int fmglev, MF &fine, const MF &crse) const override
Add interpolated coarse MG level data to fine MG level data.
Definition: AMReX_MLCellLinOp.H:952
virtual void addInhomogNeumannFlux(int, const Array< MF *, AMREX_SPACEDIM > &, MF const &, bool) const
Definition: AMReX_MLCellLinOp.H:122
RT xdoty(int amrlev, int mglev, const MF &x, const MF &y, bool local) const final
x dot y, used by the bottom solver
Definition: AMReX_MLCellLinOp.H:1926
void prepareForSolve() override
Definition: AMReX_MLCellLinOp.H:1613
void unapplyMetricTerm(int amrlev, int mglev, MF &rhs) const final
unapply metric terms if there are any
Definition: AMReX_MLCellLinOp.H:1498
void apply(int amrlev, int mglev, MF &out, MF &in, BCMode bc_mode, StateMode s_mode, const MLMGBndryT< MF > *bndry=nullptr) const override
Apply the linear operator, out = L(in)
Definition: AMReX_MLCellLinOp.H:1192
void applyMetricTerm(int amrlev, int mglev, MF &rhs) const final
apply metric terms if there are any
Definition: AMReX_MLCellLinOp.H:1442
Vector< std::unique_ptr< BndryRegisterT< MF > > > m_crse_cor_br
Definition: AMReX_MLCellLinOp.H:152
Vector< Vector< std::unique_ptr< BndryCondLoc > > > m_bcondloc
Definition: AMReX_MLCellLinOp.H:194
void setGaussSeidel(bool flag) noexcept
Definition: AMReX_MLCellLinOp.H:58
void defineAuxData()
Definition: AMReX_MLCellLinOp.H:368
Vector< std::unique_ptr< MF > > m_robin_bcval
Definition: AMReX_MLCellLinOp.H:138
Vector< Vector< Array< MultiMask, 2 *AMREX_SPACEDIM > > > m_maskvals
Definition: AMReX_MLCellLinOp.H:200
void define(const Vector< Geometry > &a_geom, const Vector< BoxArray > &a_grids, const Vector< DistributionMapping > &a_dmap, const LPInfo &a_info=LPInfo(), const Vector< FabFactory< FAB > const * > &a_factory={})
Definition: AMReX_MLCellLinOp.H:355
~MLCellLinOpT() override=default
virtual bool isTensorOp() const
Definition: AMReX_MLCellLinOp.H:61
void restriction(int, int, MF &crse, MF &fine) const override
Restriction onto coarse MG level.
Definition: AMReX_MLCellLinOp.H:943
Vector< Vector< BndryRegisterT< MF > > > m_undrrelxr
Definition: AMReX_MLCellLinOp.H:197
int m_interpbndry_halfwidth
Definition: AMReX_MLCellLinOp.H:216
MLCellLinOpT()
Definition: AMReX_MLCellLinOp.H:348
MLCellLinOpT(MLCellLinOpT< MF > &&)=delete
MLCellLinOpT< MF > & operator=(const MLCellLinOpT< MF > &)=delete
void interpAssign(int amrlev, int fmglev, MF &fine, MF &crse) const override
Overwrite fine MG level data with interpolated coarse data.
Definition: AMReX_MLCellLinOp.H:999
void fixSolvabilityByOffset(int amrlev, int mglev, MF &rhs, Vector< RT > const &offset) const override
Definition: AMReX_MLCellLinOp.H:1593
Vector< std::unique_ptr< MLMGBndryT< MF > > > m_bndry_cor
Definition: AMReX_MLCellLinOp.H:151
typename MLLinOpT< MF >::Location Location
Definition: AMReX_MLCellLinOp.H:30
void prepareForFluxes(int amrlev, const MF *crse_bcdata=nullptr) override
Definition: AMReX_MLCellLinOp.H:1235
void solutionResidual(int amrlev, MF &resid, MF &x, const MF &b, const MF *crse_bcdata=nullptr) override
Compute residual for solution.
Definition: AMReX_MLCellLinOp.H:1217
Vector< Vector< RT > > m_volinv
Definition: AMReX_MLCellLinOp.H:214
bool m_has_metric_term
Definition: AMReX_MLCellLinOp.H:146
bool m_use_gauss_seidel
Definition: AMReX_MLCellLinOp.H:206
Vector< std::unique_ptr< BndryRegisterT< MF > > > m_crse_sol_br
Definition: AMReX_MLCellLinOp.H:149
void interpolationAmr(int famrlev, MF &fine, const MF &crse, IntVect const &nghost) const override
Interpolation between AMR levels.
Definition: AMReX_MLCellLinOp.H:1092
void compFlux(int amrlev, const Array< MF *, AMREX_SPACEDIM > &fluxes, MF &sol, Location loc) const override
Compute fluxes.
Definition: AMReX_MLCellLinOp.H:1344
void computeVolInv() const
Definition: AMReX_MLCellLinOp.H:1940
Vector< std::unique_ptr< iMultiFab > > m_norm_fine_mask
Definition: AMReX_MLCellLinOp.H:202
virtual void FFlux(int amrlev, const MFIter &mfi, const Array< FAB *, AMREX_SPACEDIM > &flux, const FAB &sol, Location loc, int face_only=0) const =0
typename MLLinOpT< MF >::StateMode StateMode
Definition: AMReX_MLCellLinOp.H:29
void setLevelBC(int amrlev, const MF *levelbcdata, const MF *robinbc_a=nullptr, const MF *robinbc_b=nullptr, const MF *robinbc_f=nullptr) final
Set boundary conditions for given level. For cell-centered solves only.
Definition: AMReX_MLCellLinOp.H:513
bool needsUpdate() const override
Does it need update if it's reused?
Definition: AMReX_MLCellLinOp.H:53
Definition: AMReX_MLLinOp.H:98
const Geometry & Geom(int amr_lev, int mglev=0) const noexcept
Definition: AMReX_MLLinOp.H:559
const MF * m_coarse_data_for_bc
Definition: AMReX_MLLinOp.H:628
Vector< Vector< std::unique_ptr< FabFactory< FAB > > > > m_factory
Definition: AMReX_MLLinOp.H:601
const Vector< int > & AMRRefRatio() const noexcept
Return AMR refinement ratios.
Definition: AMReX_MLLinOp.H:632
IntVect m_ixtype
Definition: AMReX_MLLinOp.H:589
Array< Real, AMREX_SPACEDIM > m_domain_bloc_hi
Definition: AMReX_MLLinOp.H:622
RealVect m_coarse_bc_loc
Definition: AMReX_MLLinOp.H:627
virtual bool needsUpdate() const
Does it need update if it's reused?
Definition: AMReX_MLLinOp.H:266
Vector< Array< BCType, AMREX_SPACEDIM > > m_hibc_orig
Definition: AMReX_MLLinOp.H:567
FabFactory< FAB > const * Factory(int amr_lev, int mglev=0) const noexcept
Definition: AMReX_MLLinOp.H:637
bool hasRobinBC() const noexcept
Definition: AMReX_MLLinOp.H:1265
Array< Real, AMREX_SPACEDIM > m_domain_bloc_lo
Definition: AMReX_MLLinOp.H:621
Vector< Vector< BoxArray > > m_grids
Definition: AMReX_MLLinOp.H:599
virtual int getNComp() const
Return number of components.
Definition: AMReX_MLLinOp.H:261
Vector< int > m_amr_ref_ratio
Definition: AMReX_MLLinOp.H:584
Vector< int > m_num_mg_levels
Definition: AMReX_MLLinOp.H:586
Vector< Vector< DistributionMapping > > m_dmap
Definition: AMReX_MLLinOp.H:600
IntVect m_coarse_data_crse_ratio
Definition: AMReX_MLLinOp.H:626
Vector< Array< BCType, AMREX_SPACEDIM > > m_lobc_orig
Definition: AMReX_MLLinOp.H:566
bool needsCoarseDataForBC() const noexcept
Needs coarse data for bc?
Definition: AMReX_MLLinOp.H:177
virtual void update()
Update for reuse.
Definition: AMReX_MLLinOp.H:268
bool hasHiddenDimension() const noexcept
Definition: AMReX_MLLinOp.H:695
Vector< Array< BCType, AMREX_SPACEDIM > > m_lobc
Definition: AMReX_MLLinOp.H:562
int hiddenDirection() const noexcept
Definition: AMReX_MLLinOp.H:696
Vector< Vector< Geometry > > m_geom
first Vector is for amr level and second is mg level
Definition: AMReX_MLLinOp.H:598
Box compactify(Box const &b) const noexcept
Definition: AMReX_MLLinOp.H:1272
int maxorder
Definition: AMReX_MLLinOp.H:579
void define(const Vector< Geometry > &a_geom, const Vector< BoxArray > &a_grids, const Vector< DistributionMapping > &a_dmap, const LPInfo &a_info, const Vector< FabFactory< FAB > const * > &a_factory, bool eb_limit_coarsening=true)
Definition: AMReX_MLLinOp.H:757
Vector< IntVect > mg_coarsen_ratio_vec
Definition: AMReX_MLLinOp.H:595
LPInfo info
Definition: AMReX_MLLinOp.H:575
Vector< Array< BCType, AMREX_SPACEDIM > > m_hibc
Definition: AMReX_MLLinOp.H:563
LinOpBCType m_coarse_fine_bc_type
Definition: AMReX_MLLinOp.H:625
int m_num_amr_levels
Definition: AMReX_MLLinOp.H:583
Definition: AMReX_MLMGBndry.H:12
static void setBoxBC(RealTuple &bloc, BCTuple &bctag, const Box &bx, const Box &domain, const Array< LinOpBCType, AMREX_SPACEDIM > &lo, const Array< LinOpBCType, AMREX_SPACEDIM > &hi, const Real *dx, IntVect const &ratio, const RealVect &interior_bloc, const Array< Real, AMREX_SPACEDIM > &domain_bloc_lo, const Array< Real, AMREX_SPACEDIM > &domain_bloc_hi, const GpuArray< int, AMREX_SPACEDIM > &is_periodic, LinOpBCType a_crse_fine_bc_type)
Definition: AMReX_MLMGBndry.H:107
Definition: AMReX_Mask.H:28
A collection (stored as an array) of FArrayBox objects.
Definition: AMReX_MultiFab.H:38
Real sum(int comp=0, bool local=false) const
Returns the sum of component "comp" over the MultiFab – no ghost cells are included.
Definition: AMReX_MultiFab.cpp:1262
An Iterator over the Orientation of Faces of a Box.
Definition: AMReX_Orientation.H:135
Encapsulation of the Orientation of the Faces of a Box.
Definition: AMReX_Orientation.H:29
@ low
Definition: AMReX_Orientation.H:34
@ high
Definition: AMReX_Orientation.H:34
Definition: AMReX_PODVector.H:246
void reserve(size_type a_capacity)
Definition: AMReX_PODVector.H:647
iterator begin() noexcept
Definition: AMReX_PODVector.H:601
iterator end() noexcept
Definition: AMReX_PODVector.H:605
void push_back(const T &a_value)
Definition: AMReX_PODVector.H:556
A Real vector in SpaceDim-dimensional space.
Definition: AMReX_RealVect.H:32
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition: AMReX_Vector.H:27
#define abs(x)
Definition: complex-type.h:85
void copyAsync(HostToDevice, InIter begin, InIter end, OutIter result) noexcept
A host-to-device copy routine. Note this is just a wrapper around memcpy, so it assumes contiguous st...
Definition: AMReX_GpuContainers.H:233
static constexpr HostToDevice hostToDevice
Definition: AMReX_GpuContainers.H:98
void streamSynchronize() noexcept
Definition: AMReX_GpuDevice.H:237
bool inLaunchRegion() noexcept
Definition: AMReX_GpuControl.H:86
bool notInLaunchRegion() noexcept
Definition: AMReX_GpuControl.H:87
void Sum(T &v, MPI_Comm comm)
Definition: AMReX_ParallelReduce.H:204
void Max(KeyValuePair< K, V > &vi, MPI_Comm comm)
Definition: AMReX_ParallelReduce.H:126
MPI_Comm CommunicatorSub() noexcept
sub-communicator for current frame
Definition: AMReX_ParallelContext.H:70
static int f(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:44
static int ff(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:59
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlmg_lin_cc_interp_r2(Box const &bx, Array4< T > const &ff, Array4< T const > const &cc, int nc) noexcept
Definition: AMReX_MLMG_2D_K.H:9
@ max
Definition: AMReX_ParallelReduce.H:17
Definition: AMReX_Amr.cpp:49
std::enable_if_t< std::is_integral_v< T > > ParallelFor(TypeList< CTOs... > ctos, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
Definition: AMReX_CTOParallelForImpl.H:200
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > shift(const BoxND< dim > &b, int dir, int nzones) noexcept
Return a BoxND with indices shifted by nzones in dir direction.
Definition: AMReX_Box.H:1372
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > adjCellLo(const BoxND< dim > &b, int dir, int len=1) noexcept
Returns the cell centered BoxND of length len adjacent to b on the low end along the coordinate direc...
Definition: AMReX_Box.H:1591
void average_down(const MultiFab &S_fine, MultiFab &S_crse, const Geometry &fgeom, const Geometry &cgeom, int scomp, int ncomp, int rr)
Definition: AMReX_MultiFabUtil.cpp:315
BoxND< AMREX_SPACEDIM > Box
Definition: AMReX_BaseFwd.H:27
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlmg_lin_cc_interp_r2(Box const &bx, Array4< T > const &ff, Array4< T const > const &cc, int nc) noexcept
Definition: AMReX_MLMG_1D_K.H:9
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > adjCellHi(const BoxND< dim > &b, int dir, int len=1) noexcept
Similar to adjCellLo but builds an adjacent BoxND on the high end.
Definition: AMReX_Box.H:1612
void EB_set_covered(MultiFab &mf, Real val)
Definition: AMReX_EBMultiFabUtil.cpp:21
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlmg_lin_cc_interp_r4(Box const &bx, Array4< T > const &ff, Array4< T const > const &cc, int nc) noexcept
Definition: AMReX_MLMG_1D_K.H:27
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > grow(const BoxND< dim > &b, int i) noexcept
Grow BoxND in all directions by given amount.
Definition: AMReX_Box.H:1211
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > adjCell(const BoxND< dim > &b, Orientation face, int len=1) noexcept
Similar to adjCellLo and adjCellHi; operates on given face.
Definition: AMReX_Box.H:1634
constexpr AMREX_GPU_HOST_DEVICE GpuTupleElement< I, GpuTuple< Ts... > >::type & get(GpuTuple< Ts... > &tup) noexcept
Definition: AMReX_Tuple.H:179
bool isMFIterSafe(const FabArrayBase &x, const FabArrayBase &y)
Definition: AMReX_MFIter.H:219
ReduceData< Ts... >::Type ParReduce(TypeList< Ops... > operation_list, TypeList< Ts... > type_list, FabArray< FAB > const &fa, IntVect const &nghost, F &&f)
Parallel reduce for MultiFab/FabArray.
Definition: AMReX_ParReduce.H:47
void EB_average_down(const MultiFab &S_fine, MultiFab &S_crse, const MultiFab &vol_fine, const MultiFab &vfrac_fine, int scomp, int ncomp, const IntVect &ratio)
Definition: AMReX_EBMultiFabUtil.cpp:336
IntVectND< AMREX_SPACEDIM > IntVect
Definition: AMReX_BaseFwd.H:30
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition: AMReX.H:111
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE T norm(const GpuComplex< T > &a_z) noexcept
Return the norm (magnitude squared) of a complex number.
Definition: AMReX_GpuComplex.H:344
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > surroundingNodes(const BoxND< dim > &b, int dir) noexcept
Returns a BoxND with NODE based coordinates in direction dir that encloses BoxND b....
Definition: AMReX_Box.H:1399
void Xpay(MF &dst, typename MF::value_type a, MF const &src, int scomp, int dcomp, int ncomp, IntVect const &nghost)
dst = src + a * dst
Definition: AMReX_FabArrayUtility.H:1654
bool TilingIfNotGPU() noexcept
Definition: AMReX_MFIter.H:12
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mllinop_apply_bc_x(int side, Box const &box, int blen, Array4< T > const &phi, Array4< int const > const &mask, BoundCond bct, T bcl, Array4< T const > const &bcval, int maxorder, T dxinv, int inhomog, int icomp) noexcept
Definition: AMReX_MLLinOp_K.H:14
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mllinop_comp_interp_coef0_x(int side, Box const &box, int blen, Array4< T > const &f, Array4< int const > const &mask, BoundCond bct, T bcl, int maxorder, T dxinv, int icomp) noexcept
Definition: AMReX_MLLinOp_K.H:329
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mllinop_apply_bc_z(int side, Box const &box, int blen, Array4< T > const &phi, Array4< int const > const &mask, BoundCond bct, T bcl, Array4< T const > const &bcval, int maxorder, T dzinv, int inhomog, int icomp) noexcept
Definition: AMReX_MLLinOp_K.H:224
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > coarsen(const BoxND< dim > &b, int ref_ratio) noexcept
Coarsen BoxND by given (positive) refinement ratio. NOTE: if type(dir) = CELL centered: lo <- lo/rati...
Definition: AMReX_Box.H:1304
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:221
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mllinop_comp_interp_coef0_y(int side, Box const &box, int blen, Array4< T > const &f, Array4< int const > const &mask, BoundCond bct, T bcl, int maxorder, T dyinv, int icomp) noexcept
Definition: AMReX_MLLinOp_K.H:410
FAB::value_type Dot(FabArray< FAB > const &x, int xcomp, FabArray< FAB > const &y, int ycomp, int ncomp, IntVect const &nghost, bool local=false)
Compute dot products of two FabArrays.
Definition: AMReX_FabArrayUtility.H:1554
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mllinop_comp_interp_coef0_z(int side, Box const &box, int blen, Array4< T > const &f, Array4< int const > const &mask, BoundCond bct, T bcl, int maxorder, T dzinv, int icomp) noexcept
Definition: AMReX_MLLinOp_K.H:491
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mllinop_apply_bc_y(int side, Box const &box, int blen, Array4< T > const &phi, Array4< int const > const &mask, BoundCond bct, T bcl, Array4< T const > const &bcval, int maxorder, T dyinv, int inhomog, int icomp) noexcept
Definition: AMReX_MLLinOp_K.H:119
iMultiFab makeFineMask(const BoxArray &cba, const DistributionMapping &cdm, const BoxArray &fba, const IntVect &ratio, int crse_value, int fine_value)
Definition: AMReX_MultiFabUtil.cpp:602
std::array< T, N > Array
Definition: AMReX_Array.H:23
integer, parameter dp
Definition: AMReX_SDCquadrature.F90:8
Definition: AMReX_Array4.H:61
Definition: AMReX_Dim3.H:12
int x
Definition: AMReX_Dim3.H:12
int z
Definition: AMReX_Dim3.H:12
int y
Definition: AMReX_Dim3.H:12
Definition: AMReX_Array.H:33
Definition: AMReX_MLLinOp.H:35
bool has_metric_term
Definition: AMReX_MLLinOp.H:43
StateMode
Definition: AMReX_MLLinOp.H:86
BCMode
Definition: AMReX_MLLinOp.H:85
Location
Definition: AMReX_MLLinOp.H:87
FabArray memory allocation information.
Definition: AMReX_FabArray.H:65
Definition: AMReX_MFIter.H:20
MFItInfo & EnableTiling(const IntVect &ts=FabArrayBase::mfiter_tile_size) noexcept
Definition: AMReX_MFIter.H:29
MFItInfo & SetDynamic(bool f) noexcept
Definition: AMReX_MFIter.H:34
Definition: AMReX_MLCellLinOp.H:133
BoundCond type
Definition: AMReX_MLCellLinOp.H:134
RT location
Definition: AMReX_MLCellLinOp.H:135
Definition: AMReX_MLCellLinOp.H:220
T bcloc_hi
Definition: AMReX_MLCellLinOp.H:227
Box bx
Definition: AMReX_MLCellLinOp.H:228
Array4< T const > bcval_hi
Definition: AMReX_MLCellLinOp.H:223
Array4< int const > mask_hi
Definition: AMReX_MLCellLinOp.H:225
int comp
Definition: AMReX_MLCellLinOp.H:232
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Box const & box() const noexcept
Definition: AMReX_MLCellLinOp.H:236
int dir
Definition: AMReX_MLCellLinOp.H:233
Array4< T const > bcval_lo
Definition: AMReX_MLCellLinOp.H:222
int blen
Definition: AMReX_MLCellLinOp.H:231
Array4< T > fab
Definition: AMReX_MLCellLinOp.H:221
BoundCond bctype_lo
Definition: AMReX_MLCellLinOp.H:229
T bcloc_lo
Definition: AMReX_MLCellLinOp.H:226
Array4< int const > mask_lo
Definition: AMReX_MLCellLinOp.H:224
BoundCond bctype_hi
Definition: AMReX_MLCellLinOp.H:230
Definition: AMReX_MLCellLinOp.H:240
int blen
Definition: AMReX_MLCellLinOp.H:250
BoundCond bcthi
Definition: AMReX_MLCellLinOp.H:249
int comp
Definition: AMReX_MLCellLinOp.H:251
Array4< T > flo
Definition: AMReX_MLCellLinOp.H:241
Array4< int const > mlo
Definition: AMReX_MLCellLinOp.H:243
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Box const & box() const noexcept
Definition: AMReX_MLCellLinOp.H:255
int dir
Definition: AMReX_MLCellLinOp.H:252
Array4< T > fhi
Definition: AMReX_MLCellLinOp.H:242
T bcllo
Definition: AMReX_MLCellLinOp.H:245
Box bx
Definition: AMReX_MLCellLinOp.H:247
BoundCond bctlo
Definition: AMReX_MLCellLinOp.H:248
Array4< int const > mhi
Definition: AMReX_MLCellLinOp.H:244
T bclhi
Definition: AMReX_MLCellLinOp.H:246
Struct for holding types.
Definition: AMReX_TypeList.H:12