1#ifndef AMREX_ML_LINOP_H_
2#define AMREX_ML_LINOP_H_
3#include <AMReX_Config.H>
5#if defined(AMREX_USE_HYPRE) && (AMREX_SPACEDIM > 1)
10#if defined(AMREX_USE_PETSC) && (AMREX_SPACEDIM > 1)
90template <
typename T>
class MLMGT;
91template <
typename T>
class MLCGSolverT;
92template <
typename T>
class MLPoissonT;
93template <
typename T>
class MLABecLaplacianT;
94template <
typename T>
class GMRESMLMGT;
101 template <
typename T>
friend class MLMGT;
129 bool eb_limit_coarsening =
true);
131 [[nodiscard]]
virtual std::string
name ()
const {
return std::string(
"Unspecified"); }
201 LinOpBCType bc_type = LinOpBCType::Dirichlet)
noexcept;
204 LinOpBCType bc_type = LinOpBCType::Dirichlet)
noexcept;
206 template <
typename AMF, std::enable_if_t<!std::is_same_v<MF,AMF>,
int> = 0>
208 LinOpBCType bc_type = LinOpBCType::Dirichlet)
noexcept;
210 template <
typename AMF, std::enable_if_t<!std::is_same_v<MF,AMF>,
int> = 0>
212 LinOpBCType bc_type = LinOpBCType::Dirichlet)
noexcept;
235 const MF* =
nullptr) = 0;
237 template <
typename AMF, std::enable_if_t<!std::is_same_v<MF,AMF>,
int> = 0>
239 const AMF* robinbc_a =
nullptr,
240 const AMF* robinbc_b =
nullptr,
241 const AMF* robinbc_f =
nullptr);
261 [[nodiscard]]
virtual int getNComp ()
const {
return 1; }
263 [[nodiscard]]
virtual int getNGrow (
int = 0,
int = 0)
const {
return 0; }
302 amrex::Abort(
"MLLinOpT::interpAssign: Must be implemented for FMG cycle");
317 amrex::Abort(
"MLLinOpT::interpolationAmr: Must be implemented for composite solves across multiple AMR levels");
330 const MF& fine_sol,
const MF& fine_rhs)
333 amrex::Abort(
"MLLinOpT::averageDownSolutionRHS: Must be implemented for composite solves across multiple AMR levels");
347 virtual void apply (
int amrlev,
int mglev, MF& out, MF& in,
BCMode bc_mode,
359 virtual void smooth (
int amrlev,
int mglev, MF& sol,
const MF& rhs,
360 bool skip_fillboundary=
false)
const = 0;
363 virtual void normalize (
int amrlev,
int mglev, MF& mf)
const {
377 const MF* crse_bcdata=
nullptr) = 0;
393 BCMode bc_mode,
const MF* crse_bcdata=
nullptr) = 0;
407 MF& res,
const MF& crse_sol,
const MF& crse_rhs,
408 MF& fine_res, MF& fine_sol,
const MF& fine_rhs)
const
412 amrex::Abort(
"MLLinOpT::reflux: Must be implemented for composite solves across multiple AMR levels");
426 amrex::Abort(
"AMReX_MLLinOp::compFlux::How did we get here?");
440 amrex::Abort(
"AMReX_MLLinOp::compGrad::How did we get here?");
458 [[nodiscard]]
virtual bool scaleRHS (
int , MF* )
const {
464 MF
const& )
const {
return {}; }
479 amrex::Warning(
"This function might need to be implemented for GMRES to work with this LinOp.");
483 [[nodiscard]]
virtual bool isSingular (
int amrlev)
const = 0;
488 virtual RT xdoty (
int amrlev,
int mglev,
const MF&
x,
const MF& y,
bool local)
const = 0;
494 virtual std::unique_ptr<MLLinOpT<MF>>
makeNLinOp (
int )
const
496 amrex::Abort(
"MLLinOp::makeNLinOp: N-Solve not supported");
503 amrex::Abort(
"MLLinOp::getFluxes: How did we get here?");
507 amrex::Abort(
"MLLinOp::getFluxes: How did we get here?");
513 amrex::Abort(
"MLLinOp::getEBFluxes: How did we get here?");
517#if defined(AMREX_USE_HYPRE) && (AMREX_SPACEDIM > 1)
519 amrex::Abort(
"MLLinOp::makeHypre: How did we get here?");
522 [[nodiscard]]
virtual std::unique_ptr<HypreNodeLap> makeHypreNodeLap(
524 const std::string& )
const
526 amrex::Abort(
"MLLinOp::makeHypreNodeLap: How did we get here?");
531#if defined(AMREX_USE_PETSC) && (AMREX_SPACEDIM > 1)
532 [[nodiscard]]
virtual std::unique_ptr<PETScABecLap> makePETSc ()
const {
533 amrex::Abort(
"MLLinOp::makePETSc: How did we get here?");
544 [[nodiscard]]
virtual RT normInf (
int amrlev, MF
const& mf,
bool local)
const = 0;
551 amrex::Abort(
"MLLinOpT::avgDownResAmr: Must be implemented for composite solves across multiple AMR levels");
560 [[nodiscard]]
bool isMFIterSafe (
int amrlev,
int mglev1,
int mglev2)
const;
568 [[nodiscard]]
const Geometry&
Geom (
int amr_lev,
int mglev=0) const noexcept {
return m_geom[amr_lev][mglev]; }
692 [[nodiscard]]
virtual MF
make (
int amrlev,
int mglev,
IntVect const& ng)
const;
694 [[nodiscard]]
virtual MF
makeAlias (MF
const& mf)
const;
700 [[nodiscard]]
virtual std::unique_ptr<FabFactory<FAB> >
makeFactory (
int ,
int )
const {
701 return std::make_unique<DefaultFabFactory<FAB>>();
710 template <
typename T>
714 return Array4<T>(a.dataPtr(), {a.begin.y,a.begin.z,0}, {a.end.y,a.end.z,1}, a.nComp());
716 return Array4<T>(a.dataPtr(), {a.begin.x,a.begin.z,0}, {a.end.x,a.end.z,1}, a.nComp());
718 return Array4<T>(a.dataPtr(), {a.begin.x,a.begin.y,0}, {a.end.x,a.end.y,1}, a.nComp());
724 template <
typename T>
725 [[nodiscard]] T
get_d0 (T
const& dx, T
const& dy, T
const&)
const noexcept
734 template <
typename T>
735 [[nodiscard]] T
get_d1 (T
const&, T
const& dy, T
const& dz)
const noexcept
753 int ratio,
int strategy);
766template <
typename MF>
773 [[maybe_unused]]
bool eb_limit_coarsening)
782 if (info.con_grid_size <= 0) { info.con_grid_size =
AMREX_D_PICK(32, 16, 8); }
792 if (!a_factory.empty() && eb_limit_coarsening) {
795 info.max_coarsening_level = std::min(info.max_coarsening_level,
796 f->maxCoarseningLevel());
800 defineGrids(a_geom, a_grids, a_dmap, a_factory);
804template <
typename MF>
814 if ( ! a_factory.empty() ) {
816 if (ebf && !(ebf->isAllRegular())) {
817 mg_domain_min_width = 4;
822 m_num_amr_levels = 0;
823 for (
int amrlev = 0; amrlev < a_geom.
size(); amrlev++) {
824 if (!a_grids[amrlev].empty()) {
829 m_amr_ref_ratio.resize(m_num_amr_levels);
830 m_num_mg_levels.resize(m_num_amr_levels);
832 m_geom.resize(m_num_amr_levels);
833 m_grids.resize(m_num_amr_levels);
834 m_dmap.resize(m_num_amr_levels);
835 m_factory.resize(m_num_amr_levels);
839 const RealBox& rb = a_geom[0].ProbDomain();
840 const int coord = a_geom[0].Coord();
843 IntVect mg_coarsen_ratio_v(mg_coarsen_ratio);
844 IntVect mg_box_min_width_v(mg_box_min_width);
845 IntVect mg_domain_min_width_v(mg_domain_min_width);
846 if (hasHiddenDimension()) {
848 "Hidden direction only supported for 3d level solve");
849 mg_coarsen_ratio_v[info.hidden_direction] = 1;
850 mg_box_min_width_v[info.hidden_direction] = 0;
851 mg_domain_min_width_v[info.hidden_direction] = 0;
855 for (
int amrlev = m_num_amr_levels-1; amrlev > 0; --amrlev)
857 m_num_mg_levels[amrlev] = 1;
858 m_geom[amrlev].push_back(a_geom[amrlev]);
859 m_grids[amrlev].push_back(a_grids[amrlev]);
860 m_dmap[amrlev].push_back(a_dmap[amrlev]);
861 if (amrlev < a_factory.size()) {
862 m_factory[amrlev].emplace_back(a_factory[amrlev]->clone());
867 IntVect rr = mg_coarsen_ratio_v;
868 const Box& dom = a_geom[amrlev].Domain();
869 for (
int i = 0; i < 2; ++i)
874 if (cdom == a_geom[amrlev-1].Domain()) {
break; }
876 ++(m_num_mg_levels[amrlev]);
878 m_geom[amrlev].emplace_back(cdom, rb, coord, is_per);
880 m_grids[amrlev].push_back(a_grids[amrlev]);
882 m_grids[amrlev].back().coarsen(rr);
884 m_dmap[amrlev].push_back(a_dmap[amrlev]);
886 rr *= mg_coarsen_ratio_v;
889#if (AMREX_SPACEDIM > 1)
890 if (hasHiddenDimension()) {
891 m_amr_ref_ratio[amrlev-1] = rr[AMREX_SPACEDIM-info.hidden_direction];
895 m_amr_ref_ratio[amrlev-1] = rr[0];
900 m_num_mg_levels[0] = 1;
901 m_geom[0].push_back(a_geom[0]);
902 m_grids[0].push_back(a_grids[0]);
903 m_dmap[0].push_back(a_dmap[0]);
904 if (!a_factory.empty()) {
905 m_factory[0].emplace_back(a_factory[0]->clone());
910 m_domain_covered.
resize(m_num_amr_levels,
false);
911 auto npts0 = m_grids[0][0].numPts();
912 m_domain_covered[0] = (npts0 == compactify(m_geom[0][0].Domain()).numPts());
913 for (
int amrlev = 1; amrlev < m_num_amr_levels; ++amrlev)
915 if (!m_domain_covered[amrlev-1]) {
break; }
916 m_domain_covered[amrlev] = (m_grids[amrlev][0].numPts() ==
917 compactify(m_geom[amrlev][0].Domain()).numPts());
921 bool aggable =
false;
923 if (m_grids[0][0].size() > 1 && info.do_agglomeration)
925 if (m_domain_covered[0])
927 aggbox = m_geom[0][0].Domain();
928 if (hasHiddenDimension()) {
929 aggbox.
makeSlab(hiddenDirection(), m_grids[0][0][0].smallEnd(hiddenDirection()));
935 aggbox = m_grids[0][0].minimalBox();
936 aggable = (aggbox.
numPts() == npts0);
942 int agg_lev = 0, con_lev = 0;
945 && info.semicoarsening_direction >= -1
946 && info.semicoarsening_direction < AMREX_SPACEDIM );
948 if (info.do_agglomeration && aggable)
950 Box dbx = m_geom[0][0].Domain();
952 Real
const nbxs =
static_cast<Real
>(m_grids[0][0].size());
953 Real
const threshold_npts =
static_cast<Real
>(
AMREX_D_TERM(info.agg_grid_size,
955 *info.agg_grid_size));
962 for (
int lev = 0; lev < info.max_coarsening_level; ++lev)
964 IntVect rr_level = mg_coarsen_ratio_v;
965 bool const do_semicoarsening_level = info.do_semicoarsening
966 && numsclevs < info.max_semicoarsening_level;
967 if (do_semicoarsening_level
968 && info.semicoarsening_direction != -1)
970 rr_level[info.semicoarsening_direction] = 1;
973 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
975 rr_dir[idim] = rr_level[idim];
976 is_coarsenable[idim] = dbx.
coarsenable(rr_dir, mg_domain_min_width_v)
978 if (!is_coarsenable[idim] && do_semicoarsening_level
979 && info.semicoarsening_direction == -1)
981 is_coarsenable[idim] =
true;
988 if (do_semicoarsening_level && info.semicoarsening_direction == -1) {
990 int n_ones =
AMREX_D_TERM(
static_cast<int>(rr_level[0] == 1),
991 +
static_cast<int>(rr_level[1] == 1),
992 +
static_cast<int>(rr_level[2] == 1));
993 if (n_ones > 1) {
break; }
995 if (rr_level != mg_coarsen_ratio_v) {
999 accum_coarsen_ratio.push_back(accum_coarsen_ratio.back()*rr_level);
1000 domainboxes.push_back(dbx.
coarsen(rr_level));
1001 boundboxes.push_back(bbx.
coarsen(rr_level));
1002 bool to_agg = (bbx.
d_numPts() / nbxs) < 0.999*threshold_npts;
1003 agg_flag.push_back(to_agg);
1006 for (
int lev = 1, nlevs =
static_cast<int>(domainboxes.size()); lev < nlevs; ++lev) {
1007 if (!agged && !agg_flag[lev] &&
1008 a_grids[0].coarsenable(accum_coarsen_ratio[lev], mg_box_min_width_v))
1010 m_grids[0].push_back(
amrex::coarsen(a_grids[0], accum_coarsen_ratio[lev]));
1011 m_dmap[0].push_back(a_dmap[0]);
1013 IntVect cr = domainboxes[lev-1].length() / domainboxes[lev].length();
1014 if (!m_grids[0].back().coarsenable(cr)) {
1017 m_grids[0].emplace_back(boundboxes[lev]);
1018 IntVect max_grid_size(info.agg_grid_size);
1019 if (info.do_semicoarsening && info.max_semicoarsening_level >= lev
1020 && info.semicoarsening_direction != -1)
1023 AMREX_D_TERM(
int mgs_0 = (max_grid_size[0]+blen[0]-1) / blen[0];,
1024 int mgs_1 = (max_grid_size[1]+blen[1]-1) / blen[1];,
1025 int mgs_2 = (max_grid_size[2]+blen[2]-1) / blen[2]);
1026 max_grid_size[info.semicoarsening_direction]
1029 m_grids[0].back().maxSize(max_grid_size);
1036 m_geom[0].emplace_back(domainboxes[lev],rb,coord,is_per);
1041 Long consolidation_threshold = 0;
1042 Real avg_npts = 0.0;
1043 if (info.do_consolidation) {
1045 consolidation_threshold =
AMREX_D_TERM(Long(info.con_grid_size),
1046 *info.con_grid_size,
1047 *info.con_grid_size);
1050 Box const& dom0 = a_geom[0].Domain();
1053 for (
int lev = 0; lev < info.max_coarsening_level; ++lev)
1055 IntVect rr_level = mg_coarsen_ratio_v;
1056 bool do_semicoarsening_level = info.do_semicoarsening
1057 && numsclevs < info.max_semicoarsening_level;
1058 if (do_semicoarsening_level
1059 && info.semicoarsening_direction != -1)
1061 rr_level[info.semicoarsening_direction] = 1;
1064 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
1066 rr_dir[idim] = rr_vec[idim] * rr_level[idim];
1067 is_coarsenable[idim] = dom0.
coarsenable(rr_dir, mg_domain_min_width_v)
1068 && a_grids[0].coarsenable(rr_dir, mg_box_min_width_v);
1069 if (!is_coarsenable[idim] && do_semicoarsening_level
1070 && info.semicoarsening_direction == -1)
1072 is_coarsenable[idim] =
true;
1079 if (do_semicoarsening_level && info.semicoarsening_direction == -1) {
1081 int n_ones =
AMREX_D_TERM(
static_cast<int>(rr_level[0] == 1),
1082 +
static_cast<int>(rr_level[1] == 1),
1083 +
static_cast<int>(rr_level[2] == 1));
1084 if (n_ones > 1) {
break; }
1086 if (rr_level != mg_coarsen_ratio_v) {
1091 m_geom[0].emplace_back(
amrex::coarsen(dom0, rr_vec), rb, coord, is_per);
1094 if (info.do_consolidation)
1096 if (avg_npts/
static_cast<Real
>(
AMREX_D_TERM(rr_vec[0], *rr_vec[1], *rr_vec[2]))
1097 < Real(0.999)*
static_cast<Real
>(consolidation_threshold))
1100 con_lev = m_dmap[0].
size();
1105 m_dmap[0].push_back(m_dmap[0].back());
1110 m_dmap[0].push_back(a_dmap[0]);
1115 m_num_mg_levels[0] = m_grids[0].size();
1117 for (
int mglev = 0; mglev < m_num_mg_levels[0] - 1; mglev++){
1118 const Box& fine_domain = m_geom[0][mglev].Domain();
1119 const Box& crse_domain = m_geom[0][mglev+1].Domain();
1120 mg_coarsen_ratio_vec.push_back(fine_domain.
length()/crse_domain.
length());
1123 for (
int amrlev = 0; amrlev < m_num_amr_levels; ++amrlev) {
1124 if (AMRRefRatio(amrlev) == 4 && mg_coarsen_ratio_vec.empty()) {
1125 mg_coarsen_ratio_vec.push_back(
IntVect(2));
1131 makeAgglomeratedDMap(m_grids[0], m_dmap[0]);
1135 makeConsolidatedDMap(m_grids[0], m_dmap[0], info.con_ratio, info.con_strategy);
1140 m_bottom_comm = makeSubCommunicator(m_dmap[0].back());
1144 m_bottom_comm = m_default_comm;
1147 m_do_agglomeration = agged;
1148 m_do_consolidation = coned;
1152 Print() <<
"MLLinOp::defineGrids(): agglomerated AMR level 0 starting at MG level "
1153 << agg_lev <<
" of " << m_num_mg_levels[0] <<
"\n";
1155 Print() <<
"MLLinOp::defineGrids(): consolidated AMR level 0 starting at MG level "
1156 << con_lev <<
" of " << m_num_mg_levels[0]
1157 <<
" (ratio = " << info.con_ratio <<
")" <<
"\n";
1159 Print() <<
"MLLinOp::defineGrids(): no agglomeration or consolidation of AMR level 0\n";
1163 for (
int amrlev = 0; amrlev < m_num_amr_levels; ++amrlev)
1165 for (
int mglev = 1; mglev < m_num_mg_levels[amrlev]; ++mglev)
1167 m_factory[amrlev].emplace_back(makeFactory(amrlev,mglev));
1171 for (
int amrlev = 1; amrlev < m_num_amr_levels; ++amrlev)
1174 "MLLinOp: grids not coarsenable between AMR levels");
1178template <
typename MF>
1182 m_needs_coarse_data_for_bc = !m_domain_covered[0];
1184 levelbc_raii.resize(m_num_amr_levels);
1185 robin_a_raii.resize(m_num_amr_levels);
1186 robin_b_raii.resize(m_num_amr_levels);
1187 robin_f_raii.resize(m_num_amr_levels);
1190template <
typename MF>
1195 const int ncomp = getNComp();
1200template <
typename MF>
1205 const int ncomp = getNComp();
1207 "MLLinOp::setDomainBC: wrong size");
1210 m_lobc_orig = m_lobc;
1211 m_hibc_orig = m_hibc;
1212 for (
int icomp = 0; icomp < ncomp; ++icomp) {
1213 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
1214 if (m_geom[0][0].isPeriodic(idim)) {
1216 m_hibc[icomp][idim] == BCType::Periodic);
1219 m_hibc[icomp][idim] != BCType::Periodic);
1222 if (m_lobc[icomp][idim] == LinOpBCType::inhomogNeumann ||
1223 m_lobc[icomp][idim] == LinOpBCType::Robin)
1225 m_lobc[icomp][idim] = LinOpBCType::Neumann;
1228 if (m_hibc[icomp][idim] == LinOpBCType::inhomogNeumann ||
1229 m_hibc[icomp][idim] == LinOpBCType::Robin)
1231 m_hibc[icomp][idim] = LinOpBCType::Neumann;
1236 if (hasHiddenDimension()) {
1237 const int hd = hiddenDirection();
1238 for (
int n = 0; n < ncomp; ++n) {
1239 m_lobc[n][hd] = LinOpBCType::Neumann;
1240 m_hibc[n][hd] = LinOpBCType::Neumann;
1244 if (hasInhomogNeumannBC() && !supportInhomogNeumannBC()) {
1245 amrex::Abort(
"Inhomogeneous Neumann BC not supported");
1247 if (hasRobinBC() && !supportRobinBC()) {
1252template <
typename MF>
1256 int ncomp = m_lobc_orig.size();
1257 for (
int n = 0; n < ncomp; ++n) {
1258 for (
int idim = 0; idim <AMREX_SPACEDIM; ++idim) {
1259 if (m_lobc_orig[n][idim] == bct || m_hibc_orig[n][idim] == bct) {
1267template <
typename MF>
1271 return hasBC(BCType::inhomogNeumann);
1274template <
typename MF>
1278 return hasBC(BCType::Robin);
1281template <
typename MF>
1285#if (AMREX_SPACEDIM == 3)
1286 if (info.hasHiddenDimension()) {
1288 const auto& hi =
b.bigEnd();
1289 if (info.hidden_direction == 0) {
1291 }
else if (info.hidden_direction == 1) {
1303template <
typename MF>
1311 for (
int i = 1, N=
static_cast<int>(ba.
size()); i < N; ++i)
1321 for (
int iproc = 0; iproc < nprocs; ++iproc) {
1323 for (
int ibox : sfc[iproc]) {
1327 dm[i].define(std::move(pmap));
1332template <
typename MF>
1336 int ratio,
int strategy)
1338 BL_PROFILE(
"MLLinOp::makeConsolidatedDMap()");
1342 for (
int i = 1, N=
static_cast<int>(ba.
size()); i < N; ++i)
1349 const auto& pmap_fine = dm[i-1].ProcessorMap();
1352 if (strategy == 1) {
1353 for (
auto&
x: pmap) {
1356 }
else if (strategy == 2) {
1357 int nprocs_con =
static_cast<int>(std::ceil(
static_cast<Real
>(nprocs)
1358 /
static_cast<Real
>(factor)));
1359 for (
auto&
x: pmap) {
1360 auto d = std::div(
x,nprocs_con);
1363 }
else if (strategy == 3) {
1364 if (factor == ratio) {
1366 for (
int iproc = 0; iproc < nprocs; ++iproc) {
1367 for (
int ibox : sfc[iproc]) {
1372 for (
auto&
x: pmap) {
1378 dm[i].define(std::move(pmap));
1382 dm[i].define(std::move(pmap_g));
1388template <
typename MF>
1392 BL_PROFILE(
"MLLinOp::makeSubCommunicator()");
1397 std::sort(newgrp_ranks.begin(), newgrp_ranks.end());
1398 auto last = std::unique(newgrp_ranks.begin(), newgrp_ranks.end());
1399 newgrp_ranks.erase(last, newgrp_ranks.end());
1403 MPI_Comm_group(m_default_comm, &defgrp);
1405 MPI_Group_incl(defgrp,
static_cast<int>(newgrp_ranks.
size()), newgrp_ranks.data(), &newgrp);
1409 newgrp_ranks.data(),
static_cast<int>(newgrp_ranks.
size()));
1410 MPI_Group_incl(defgrp,
static_cast<int>(local_newgrp_ranks.
size()), local_newgrp_ranks.data(), &newgrp);
1413 MPI_Comm_create(m_default_comm, newgrp, &newcomm);
1415 m_raii_comm = std::make_unique<CommContainer>(newcomm);
1417 MPI_Group_free(&defgrp);
1418 MPI_Group_free(&newgrp);
1423 return m_default_comm;
1427template <
typename MF>
1432 m_domain_bloc_lo = lo_bcloc;
1433 m_domain_bloc_hi = hi_bcloc;
1436template <
typename MF>
1439 LinOpBCType bc_type)
noexcept
1441 setCoarseFineBC(
crse,
IntVect(crse_ratio), bc_type);
1444template <
typename MF>
1447 LinOpBCType bc_type)
noexcept
1449 m_coarse_data_for_bc =
crse;
1450 m_coarse_data_crse_ratio = crse_ratio;
1451 m_coarse_fine_bc_type = bc_type;
1454template <
typename MF>
1455template <
typename AMF, std::enable_if_t<!std::is_same_v<MF,AMF>,
int>>
1458 LinOpBCType bc_type)
noexcept
1460 setCoarseFineBC(
crse,
IntVect(crse_ratio), bc_type);
1463template <
typename MF>
1464template <
typename AMF, std::enable_if_t<!std::is_same_v<MF,AMF>,
int>>
1467 LinOpBCType bc_type)
noexcept
1469 m_coarse_data_for_bc_raii = MF(
crse->boxArray(),
crse->DistributionMap(),
1471 m_coarse_data_for_bc_raii.LocalCopy(*
crse, 0, 0,
crse->nComp(),
1473 m_coarse_data_for_bc = &m_coarse_data_for_bc_raii;
1474 m_coarse_data_crse_ratio = crse_ratio;
1475 m_coarse_fine_bc_type = bc_type;
1478template <
typename MF>
1483 mf.resize(m_num_amr_levels);
1484 for (
int alev = 0; alev < m_num_amr_levels; ++alev) {
1485 mf[alev].resize(m_num_mg_levels[alev]);
1486 for (
int mlev = 0; mlev < m_num_mg_levels[alev]; ++mlev) {
1487 mf[alev][mlev] = make(alev, mlev, ng);
1492template <
typename MF>
1496 if constexpr (IsMultiFabLike_v<MF>) {
1498 m_dmap[amrlev][mglev], getNComp(), ng,
MFInfo(),
1499 *m_factory[amrlev][mglev]);
1507template <
typename MF>
1511 if constexpr (IsMultiFabLike_v<MF>) {
1515 amrex::Abort(
"MLLinOpT::makeAlias: how did we get here?");
1520template <
typename MF>
1524 if constexpr (IsMultiFabLike_v<MF>) {
1525 BoxArray cba = m_grids[amrlev][mglev];
1526 IntVect ratio = (amrlev > 0) ?
IntVect(2) : mg_coarsen_ratio_vec[mglev];
1529 return MF(cba, m_dmap[amrlev][mglev], getNComp(), ng);
1532 amrex::Abort(
"MLLinOpT::makeCoarseMG: how did we get here?");
1537template <
typename MF>
1541 if constexpr (IsMultiFabLike_v<MF>) {
1542 BoxArray cba = m_grids[famrlev][0];
1543 IntVect ratio(AMRRefRatio(famrlev-1));
1546 return MF(cba, m_dmap[famrlev][0], getNComp(), ng);
1549 amrex::Abort(
"MLLinOpT::makeCoarseAmr: how did we get here?");
1554template <
typename MF>
1558 if (new_size <= 0 || new_size >= m_num_mg_levels[0]) {
return; }
1560 m_num_mg_levels[0] = new_size;
1562 m_geom[0].resize(new_size);
1563 m_grids[0].resize(new_size);
1564 m_dmap[0].resize(new_size);
1565 m_factory[0].resize(new_size);
1567 if (m_bottom_comm != m_default_comm) {
1568 m_bottom_comm = makeSubCommunicator(m_dmap[0].back());
1572template <
typename MF>
1578 const int ncomp = this->getNComp();
1580 if (!fres.isAllRegular()) {
1581 if constexpr (std::is_same<MF,MultiFab>()) {
1583 mg_coarsen_ratio_vec[clev-1]);
1585 amrex::Abort(
"EB_average_down only works with MultiFab");
1593 amrex::Abort(
"For non-FabArray, MLLinOpT<MF>::avgDownResMG should be overridden.");
1597template <
typename MF>
1601 return m_dmap[amrlev][mglev1] == m_dmap[amrlev][mglev2]
1605template <
typename MF>
1606template <
typename AMF, std::enable_if_t<!std::is_same_v<MF,AMF>,
int>>
1609 const AMF* robinbc_a,
const AMF* robinbc_b,
1610 const AMF* robinbc_f)
1612 const int ncomp = this->getNComp();
1614 levelbc_raii[amrlev] = std::make_unique<MF>(levelbcdata->boxArray(),
1615 levelbcdata->DistributionMap(),
1616 ncomp, levelbcdata->nGrowVect());
1617 levelbc_raii[amrlev]->LocalCopy(*levelbcdata, 0, 0, ncomp,
1618 levelbcdata->nGrowVect());
1620 levelbc_raii[amrlev].reset();
1624 robin_a_raii[amrlev] = std::make_unique<MF>(robinbc_a->boxArray(),
1625 robinbc_a->DistributionMap(),
1626 ncomp, robinbc_a->nGrowVect());
1627 robin_a_raii[amrlev]->LocalCopy(*robinbc_a, 0, 0, ncomp,
1628 robinbc_a->nGrowVect());
1630 robin_a_raii[amrlev].reset();
1634 robin_b_raii[amrlev] = std::make_unique<MF>(robinbc_b->boxArray(),
1635 robinbc_b->DistributionMap(),
1636 ncomp, robinbc_b->nGrowVect());
1637 robin_b_raii[amrlev]->LocalCopy(*robinbc_b, 0, 0, ncomp,
1638 robinbc_b->nGrowVect());
1640 robin_b_raii[amrlev].reset();
1644 robin_f_raii[amrlev] = std::make_unique<MF>(robinbc_f->boxArray(),
1645 robinbc_f->DistributionMap(),
1646 ncomp, robinbc_f->nGrowVect());
1647 robin_f_raii[amrlev]->LocalCopy(*robinbc_f, 0, 0, ncomp,
1648 robinbc_f->nGrowVect());
1650 robin_f_raii[amrlev].reset();
1653 this->setLevelBC(amrlev, levelbc_raii[amrlev].
get(), robin_a_raii[amrlev].
get(),
1654 robin_b_raii[amrlev].
get(), robin_f_raii[amrlev].
get());
1657template <
typename MF>
1662 return xdoty(0,0,*
x[0],*y[0],
false);
1665template <
typename MF>
1670 auto r = xdoty(0,0,*
x[0],*
x[0],
false);
1671 return std::sqrt(
r);
#define BL_PROFILE(a)
Definition AMReX_BLProfiler.H:551
#define BL_ASSERT(EX)
Definition AMReX_BLassert.H:39
#define AMREX_ASSERT_WITH_MESSAGE(EX, MSG)
Definition AMReX_BLassert.H:37
#define AMREX_ASSERT(EX)
Definition AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition AMReX_BLassert.H:50
Array4< Real > fine
Definition AMReX_InterpFaceRegister.cpp:90
Array4< Real const > crse
Definition AMReX_InterpFaceRegister.cpp:92
#define AMREX_D_TERM(a, b, c)
Definition AMReX_SPACE.H:129
#define AMREX_D_PICK(a, b, c)
Definition AMReX_SPACE.H:151
int MPI_Comm
Definition AMReX_ccse-mpi.H:47
int MPI_Group
Definition AMReX_ccse-mpi.H:48
static constexpr int MPI_COMM_NULL
Definition AMReX_ccse-mpi.H:55
A collection of Boxes stored in an Array.
Definition AMReX_BoxArray.H:550
static bool SameRefs(const BoxArray &lhs, const BoxArray &rhs)
whether two BoxArrays share the same data
Definition AMReX_BoxArray.H:820
BoxArray & coarsen(int refinement_ratio)
Coarsen each Box in the BoxArray to the specified ratio.
BoxArray & convert(IndexType typ)
Apply Box::convert(IndexType) to each Box in the BoxArray.
AMREX_GPU_HOST_DEVICE BoxND & coarsen(int ref_ratio) noexcept
Coarsen BoxND by given (positive) refinement ratio. NOTE: if type(dir) = CELL centered: lo <- lo/rati...
Definition AMReX_Box.H:708
AMREX_GPU_HOST_DEVICE BoxND & makeSlab(int direction, int slab_index) noexcept
Definition AMReX_Box.H:779
AMREX_GPU_HOST_DEVICE double d_numPts() const noexcept
Returns the number of points contained in the BoxND. This is intended for use only in diagnostic mess...
Definition AMReX_Box.H:366
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 coarsenable(const IntVectND< dim > &refrat, const IntVectND< dim > &min_width) const noexcept
Definition AMReX_Box.H:745
AMREX_GPU_HOST_DEVICE const IntVectND< dim > & smallEnd() const &noexcept
Get the smallend of the BoxND.
Definition AMReX_Box.H:105
AMREX_GPU_HOST_DEVICE Long numPts() const noexcept
Returns the number of points contained in the BoxND.
Definition AMReX_Box.H:346
Definition AMReX_FabFactory.H:76
Calculates the distribution of FABs to MPI processes.
Definition AMReX_DistributionMapping.H:41
static DistributionMapping makeSFC(const MultiFab &weight, bool sort=true)
Definition AMReX_DistributionMapping.cpp:1764
const Vector< int > & ProcessorMap() const noexcept
Returns a constant reference to the mapping of boxes in the underlying BoxArray to the CPU that holds...
Definition AMReX_DistributionMapping.cpp:47
Definition AMReX_EBFabFactory.H:24
Definition AMReX_FabFactory.H:50
Solve using GMRES with multigrid as preconditioner.
Definition AMReX_GMRES_MLMG.H:21
Rectangular problem domain geometry.
Definition AMReX_Geometry.H:73
Interface
Definition AMReX_Hypre.H:21
AMREX_GPU_HOST_DEVICE static AMREX_FORCE_INLINE constexpr std::size_t size() noexcept
Definition AMReX_IntVect.H:725
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE IntVectND< new_dim > resize(int fill_extra=0) const noexcept
Returns a new IntVectND of size new_dim by either shrinking or expanding this IntVectND.
Definition AMReX_IntVect.H:768
Definition AMReX_MLABecLaplacian.H:15
Definition AMReX_MLCGSolver.H:12
Definition AMReX_MLLinOp.H:98
const MF * m_coarse_data_for_bc
Definition AMReX_MLLinOp.H:637
virtual void smooth(int amrlev, int mglev, MF &sol, const MF &rhs, bool skip_fillboundary=false) const =0
Smooth.
void setDomainBC(const Array< BCType, AMREX_SPACEDIM > &lobc, const Array< BCType, AMREX_SPACEDIM > &hibc) noexcept
Boundary of the whole domain.
Definition AMReX_MLLinOp.H:1192
Vector< Vector< std::unique_ptr< FabFactory< FAB > > > > m_factory
Definition AMReX_MLLinOp.H:610
virtual void avgDownResMG(int clev, MF &cres, MF const &fres) const
Definition AMReX_MLLinOp.H:1574
int NAMRLevels() const noexcept
Return the number of AMR levels.
Definition AMReX_MLLinOp.H:563
bool m_do_consolidation
Definition AMReX_MLLinOp.H:601
virtual void compGrad(int, const Array< MF *, AMREX_SPACEDIM > &, MF &, Location) const
Compute gradients of the solution.
Definition AMReX_MLLinOp.H:437
bool isCellCentered() const noexcept
Definition AMReX_MLLinOp.H:684
IntVect m_ixtype
Definition AMReX_MLLinOp.H:598
void setVerbose(int v) noexcept
Set verbosity.
Definition AMReX_MLLinOp.H:244
Array< Real, AMREX_SPACEDIM > m_domain_bloc_hi
Definition AMReX_MLLinOp.H:631
bool isMFIterSafe(int amrlev, int mglev1, int mglev2) const
Definition AMReX_MLLinOp.H:1599
RealVect m_coarse_bc_loc
Definition AMReX_MLLinOp.H:636
virtual bool needsUpdate() const
Does it need update if it's reused?
Definition AMReX_MLLinOp.H:266
virtual void interpolation(int amrlev, int fmglev, MF &fine, const MF &crse) const =0
Add interpolated coarse MG level data to fine MG level data.
virtual void setLevelBC(int, const MF *, const MF *=nullptr, const MF *=nullptr, const MF *=nullptr)=0
Set boundary conditions for given level. For cell-centered solves only.
virtual MF make(int amrlev, int mglev, IntVect const &ng) const
Definition AMReX_MLLinOp.H:1494
FabFactory< FAB > const * Factory(int amr_lev, int mglev=0) const noexcept
Definition AMReX_MLLinOp.H:648
T get_d0(T const &dx, T const &dy, T const &) const noexcept
Definition AMReX_MLLinOp.H:725
MPI_Comm BottomCommunicator() const noexcept
Definition AMReX_MLLinOp.H:675
void setEnforceSingularSolvable(bool o) noexcept
Definition AMReX_MLLinOp.H:253
MPI_Comm Communicator() const noexcept
Definition AMReX_MLLinOp.H:676
virtual void getFluxes(const Vector< Array< MF *, AMREX_SPACEDIM > > &, const Vector< MF * > &, Location) const
Definition AMReX_MLLinOp.H:500
int mg_domain_min_width
Definition AMReX_MLLinOp.H:582
void setMaxOrder(int o) noexcept
Set order of interpolation at coarse/fine boundary.
Definition AMReX_MLLinOp.H:247
Vector< Array< BCType, AMREX_SPACEDIM > > m_hibc_orig
Definition AMReX_MLLinOp.H:576
virtual void interpAssign(int amrlev, int fmglev, MF &fine, MF &crse) const
Overwrite fine MG level data with interpolated coarse data.
Definition AMReX_MLLinOp.H:299
virtual std::string name() const
Definition AMReX_MLLinOp.H:131
void setCoarseFineBC(const AMF *crse, IntVect const &crse_ratio, LinOpBCType bc_type=LinOpBCType::Dirichlet) noexcept
Definition AMReX_MLLinOp.H:1466
bool doAgglomeration() const noexcept
Definition AMReX_MLLinOp.H:680
Vector< std::unique_ptr< MF > > levelbc_raii
Definition AMReX_MLLinOp.H:760
MF m_coarse_data_for_bc_raii
Definition AMReX_MLLinOp.H:638
virtual void setDirichletNodesToZero(int, int, MF &) const
Definition AMReX_MLLinOp.H:476
MLLinOpT< MF > & operator=(const MLLinOpT< MF > &)=delete
std::unique_ptr< CommContainer > m_raii_comm
Definition AMReX_MLLinOp.H:628
bool m_do_semicoarsening
Definition AMReX_MLLinOp.H:603
bool hasRobinBC() const noexcept
Definition AMReX_MLLinOp.H:1276
Array< Real, AMREX_SPACEDIM > m_domain_bloc_lo
Definition AMReX_MLLinOp.H:630
virtual void resizeMultiGrid(int new_size)
Definition AMReX_MLLinOp.H:1556
Vector< Vector< BoxArray > > m_grids
Definition AMReX_MLLinOp.H:608
virtual MF makeCoarseAmr(int famrlev, IntVect const &ng) const
Definition AMReX_MLLinOp.H:1539
Vector< std::unique_ptr< MF > > robin_f_raii
Definition AMReX_MLLinOp.H:763
bool m_do_agglomeration
Definition AMReX_MLLinOp.H:600
void setCoarseFineBC(const AMF *crse, int crse_ratio, LinOpBCType bc_type=LinOpBCType::Dirichlet) noexcept
Definition AMReX_MLLinOp.H:1457
virtual MF makeAlias(MF const &mf) const
Definition AMReX_MLLinOp.H:1509
void setDomainBC(const Vector< Array< BCType, AMREX_SPACEDIM > > &lobc, const Vector< Array< BCType, AMREX_SPACEDIM > > &hibc) noexcept
Boundary of the whole domain.
Definition AMReX_MLLinOp.H:1202
Array4< T > compactify(Array4< T > const &a) const noexcept
Definition AMReX_MLLinOp.H:711
static constexpr int mg_coarsen_ratio
Definition AMReX_MLLinOp.H:580
virtual void copyNSolveSolution(MF &, MF const &) const
Definition AMReX_MLLinOp.H:540
virtual void solutionResidual(int amrlev, MF &resid, MF &x, const MF &b, const MF *crse_bcdata=nullptr)=0
Compute residual for solution.
int getMaxOrder() const noexcept
Get order of interpolation at coarse/fine boundary.
Definition AMReX_MLLinOp.H:249
virtual int getNComp() const
Return number of components.
Definition AMReX_MLLinOp.H:261
virtual void getFluxes(const Vector< MF * > &, const Vector< MF * > &) const
Definition AMReX_MLLinOp.H:505
void defineBC()
Definition AMReX_MLLinOp.H:1180
void setCoarseFineBCLocation(const RealVect &cloc) noexcept
Definition AMReX_MLLinOp.H:678
Vector< int > m_amr_ref_ratio
Definition AMReX_MLLinOp.H:593
MPI_Comm m_default_comm
Definition AMReX_MLLinOp.H:613
static void makeAgglomeratedDMap(const Vector< BoxArray > &ba, Vector< DistributionMapping > &dm)
Definition AMReX_MLLinOp.H:1305
bool isBottomActive() const noexcept
Definition AMReX_MLLinOp.H:673
virtual Vector< RT > getSolvabilityOffset(int, int, MF const &) const
get offset for fixing solvability
Definition AMReX_MLLinOp.H:463
void defineGrids(const Vector< Geometry > &a_geom, const Vector< BoxArray > &a_grids, const Vector< DistributionMapping > &a_dmap, const Vector< FabFactory< FAB > const * > &a_factory)
Definition AMReX_MLLinOp.H:806
virtual BottomSolver getDefaultBottomSolver() const
Definition AMReX_MLLinOp.H:258
typename FabDataType< MF >::fab_type FAB
Definition AMReX_MLLinOp.H:108
virtual RT normInf(int amrlev, MF const &mf, bool local) const =0
Vector< int > m_num_mg_levels
Definition AMReX_MLLinOp.H:595
bool hasBC(BCType bct) const noexcept
Definition AMReX_MLLinOp.H:1254
static void makeConsolidatedDMap(const Vector< BoxArray > &ba, Vector< DistributionMapping > &dm, int ratio, int strategy)
Definition AMReX_MLLinOp.H:1334
Vector< Vector< DistributionMapping > > m_dmap
Definition AMReX_MLLinOp.H:609
int verbose
Definition AMReX_MLLinOp.H:586
virtual void checkPoint(std::string const &) const
Definition AMReX_MLLinOp.H:756
IntVect m_coarse_data_crse_ratio
Definition AMReX_MLLinOp.H:635
LinOpBCType BCType
Definition AMReX_MLLinOp.H:111
virtual void correctionResidual(int amrlev, int mglev, MF &resid, MF &x, const MF &b, BCMode bc_mode, const MF *crse_bcdata=nullptr)=0
Compute residual for the residual-correction form, resid = b - L(x)
virtual void compFlux(int, const Array< MF *, AMREX_SPACEDIM > &, MF &, Location) const
Compute fluxes.
Definition AMReX_MLLinOp.H:423
const Vector< int > & AMRRefRatio() const noexcept
Return AMR refinement ratios.
Definition AMReX_MLLinOp.H:643
MLLinOpT(MLLinOpT< MF > &&)=delete
Vector< Array< BCType, AMREX_SPACEDIM > > m_lobc_orig
Definition AMReX_MLLinOp.H:575
void setCoarseFineBC(const MF *crse, int crse_ratio, LinOpBCType bc_type=LinOpBCType::Dirichlet) noexcept
Set coarse/fine boundary conditions. For cell-centered solves only.
Definition AMReX_MLLinOp.H:1438
virtual void apply(int amrlev, int mglev, MF &out, MF &in, BCMode bc_mode, StateMode s_mode, const MLMGBndryT< MF > *bndry=nullptr) const =0
Apply the linear operator, out = L(in)
GpuArray< BCType, AMREX_SPACEDIM > HiBC(int icomp=0) const noexcept
Definition AMReX_MLLinOp.H:657
void setDomainBCLoc(const Array< Real, AMREX_SPACEDIM > &lo_bcloc, const Array< Real, AMREX_SPACEDIM > &hi_bcloc) noexcept
Set location of domain boundaries.
Definition AMReX_MLLinOp.H:1429
bool needsCoarseDataForBC() const noexcept
Needs coarse data for bc?
Definition AMReX_MLLinOp.H:177
typename FabDataType< MF >::value_type RT
Definition AMReX_MLLinOp.H:109
virtual void update()
Update for reuse.
Definition AMReX_MLLinOp.H:268
bool m_precond_mode
Definition AMReX_MLLinOp.H:640
virtual std::unique_ptr< FabFactory< FAB > > makeFactory(int, int) const
Definition AMReX_MLLinOp.H:700
virtual void applyMetricTerm(int, int, MF &) const
apply metric terms if there are any
Definition AMReX_MLLinOp.H:444
virtual void unapplyMetricTerm(int, int, MF &) const
unapply metric terms if there are any
Definition AMReX_MLLinOp.H:446
bool hasHiddenDimension() const noexcept
Definition AMReX_MLLinOp.H:706
virtual bool isBottomSingular() const =0
Is the bottom of MG singular?
virtual void reflux(int crse_amrlev, MF &res, const MF &crse_sol, const MF &crse_rhs, MF &fine_res, MF &fine_sol, const MF &fine_rhs) const
Reflux at AMR coarse/fine boundary.
Definition AMReX_MLLinOp.H:406
virtual IntVect getNGrowVectRestriction() const
Definition AMReX_MLLinOp.H:686
Vector< std::unique_ptr< MF > > robin_b_raii
Definition AMReX_MLLinOp.H:762
static constexpr int mg_box_min_width
Definition AMReX_MLLinOp.H:581
virtual RT dotProductPrecond(Vector< MF const * > const &x, Vector< MF const * > const &y) const
Definition AMReX_MLLinOp.H:1659
virtual MF makeCoarseMG(int amrlev, int mglev, IntVect const &ng) const
Definition AMReX_MLLinOp.H:1522
virtual void fixSolvabilityByOffset(int, int, MF &, Vector< RT > const &) const
fix solvability by subtracting offset from RHS
Definition AMReX_MLLinOp.H:467
const Geometry & Geom(int amr_lev, int mglev=0) const noexcept
Definition AMReX_MLLinOp.H:568
virtual void endPrecondBC()
Definition AMReX_MLLinOp.H:558
Vector< Array< BCType, AMREX_SPACEDIM > > m_lobc
Definition AMReX_MLLinOp.H:571
int hiddenDirection() const noexcept
Definition AMReX_MLLinOp.H:707
Vector< int > m_domain_covered
Definition AMReX_MLLinOp.H:611
const MLLinOpT< MF > * m_parent
Definition AMReX_MLLinOp.H:596
bool doSemicoarsening() const noexcept
Definition AMReX_MLLinOp.H:682
virtual bool supportNSolve() const
Definition AMReX_MLLinOp.H:538
virtual bool supportRobinBC() const noexcept
Definition AMReX_MLLinOp.H:667
virtual void normalize(int amrlev, int mglev, MF &mf) const
Divide mf by the diagonal component of the operator. Used by bicgstab.
Definition AMReX_MLLinOp.H:363
virtual void applyInhomogNeumannTerm(int, MF &) const
Extra terms introduced when we treat inhomogeneous Nuemann BC as homogeneous.
Definition AMReX_MLLinOp.H:452
virtual void avgDownResAmr(int clev, MF &cres, MF const &fres) const
Definition AMReX_MLLinOp.H:548
Vector< Vector< Geometry > > m_geom
first Vector is for amr level and second is mg level
Definition AMReX_MLLinOp.H:607
virtual RT norm2Precond(Vector< MF const * > const &x) const
Definition AMReX_MLLinOp.H:1667
MLLinOpT(const MLLinOpT< MF > &)=delete
virtual bool scaleRHS(int, MF *) const
scale RHS to fix solvability
Definition AMReX_MLLinOp.H:458
virtual void averageDownAndSync(Vector< MF > &sol) const =0
void setCoarseFineBC(const MF *crse, IntVect const &crse_ratio, LinOpBCType bc_type=LinOpBCType::Dirichlet) noexcept
Definition AMReX_MLLinOp.H:1446
GpuArray< BCType, AMREX_SPACEDIM > LoBC(int icomp=0) const noexcept
Definition AMReX_MLLinOp.H:652
Box compactify(Box const &b) const noexcept
Definition AMReX_MLLinOp.H:1283
bool m_needs_coarse_data_for_bc
Definition AMReX_MLLinOp.H:633
int maxorder
Definition AMReX_MLLinOp.H:588
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:768
Vector< IntVect > mg_coarsen_ratio_vec
Definition AMReX_MLLinOp.H:604
MF MFType
Definition AMReX_MLLinOp.H:107
virtual void preparePrecond()
Definition AMReX_MLLinOp.H:472
virtual ~MLLinOpT()=default
virtual void averageDownSolutionRHS(int camrlev, MF &crse_sol, MF &crse_rhs, const MF &fine_sol, const MF &fine_rhs)
Average-down data from fine AMR level to coarse AMR level.
Definition AMReX_MLLinOp.H:329
LPInfo info
Definition AMReX_MLLinOp.H:584
MPI_Comm makeSubCommunicator(const DistributionMapping &dm)
Definition AMReX_MLLinOp.H:1390
Vector< Array< BCType, AMREX_SPACEDIM > > m_hibc
Definition AMReX_MLLinOp.H:572
int NMGLevels(int amrlev) const noexcept
Return the number of MG levels at given AMR level.
Definition AMReX_MLLinOp.H:566
T get_d1(T const &, T const &dy, T const &dz) const noexcept
Definition AMReX_MLLinOp.H:735
bool enforceSingularSolvable
Definition AMReX_MLLinOp.H:590
void setLevelBC(int amrlev, const AMF *levelbcdata, const AMF *robinbc_a=nullptr, const AMF *robinbc_b=nullptr, const AMF *robinbc_f=nullptr)
Definition AMReX_MLLinOp.H:1608
virtual RT xdoty(int amrlev, int mglev, const MF &x, const MF &y, bool local) const =0
x dot y, used by the bottom solver
virtual void interpolationAmr(int famrlev, MF &fine, const MF &crse, IntVect const &nghost) const
Interpolation between AMR levels.
Definition AMReX_MLLinOp.H:313
bool doConsolidation() const noexcept
Definition AMReX_MLLinOp.H:681
virtual std::unique_ptr< MLLinOpT< MF > > makeNLinOp(int) const
Definition AMReX_MLLinOp.H:494
virtual bool supportInhomogNeumannBC() const noexcept
Definition AMReX_MLLinOp.H:668
virtual void prepareForFluxes(int, const MF *=nullptr)
Definition AMReX_MLLinOp.H:379
LinOpBCType m_coarse_fine_bc_type
Definition AMReX_MLLinOp.H:634
virtual bool isSingular(int amrlev) const =0
Is it singular on given AMR level?
virtual void postSolve(Vector< MF * > const &) const
Definition AMReX_MLLinOp.H:542
int AMRRefRatio(int amr_lev) const noexcept
Return AMR refinement ratio at given AMR level.
Definition AMReX_MLLinOp.H:646
MPI_Comm m_bottom_comm
Definition AMReX_MLLinOp.H:614
virtual void restriction(int amrlev, int cmglev, MF &crse, MF &fine) const =0
Restriction onto coarse MG level.
virtual void beginPrecondBC()
Definition AMReX_MLLinOp.H:557
virtual void applyOverset(int, MF &) const
for overset solver only
Definition AMReX_MLLinOp.H:455
virtual void unimposeNeumannBC(int, MF &) const
This is needed for our nodal projection solver.
Definition AMReX_MLLinOp.H:449
bool getEnforceSingularSolvable() const noexcept
Definition AMReX_MLLinOp.H:256
virtual void prepareForSolve()=0
virtual void make(Vector< Vector< MF > > &mf, IntVect const &ng) const
Definition AMReX_MLLinOp.H:1480
bool hasInhomogNeumannBC() const noexcept
Definition AMReX_MLLinOp.H:1269
Vector< std::unique_ptr< MF > > robin_a_raii
Definition AMReX_MLLinOp.H:761
virtual int getNGrow(int=0, int=0) const
Definition AMReX_MLLinOp.H:263
int m_num_amr_levels
Definition AMReX_MLLinOp.H:592
Definition AMReX_MLMGBndry.H:12
Definition AMReX_MLMG.H:12
Definition AMReX_MLPoisson.H:16
This class provides the user with a few print options.
Definition AMReX_Print.H:35
A Box with real dimensions. A RealBox is OK iff volume >= 0.
Definition AMReX_RealBox.H:21
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
Long size() const noexcept
Definition AMReX_Vector.H:50
bool notInLaunchRegion() noexcept
Definition AMReX_GpuControl.H:87
MPI_Comm CommunicatorSub() noexcept
sub-communicator for current frame
Definition AMReX_ParallelContext.H:70
int local_to_global_rank(int rank) noexcept
translate between local rank and global rank
Definition AMReX_ParallelContext.H:95
int global_to_local_rank(int rank) noexcept
Definition AMReX_ParallelContext.H:98
int NProcsSub() noexcept
number of ranks in current frame
Definition AMReX_ParallelContext.H:74
MPI_Comm Communicator() noexcept
Definition AMReX_ParallelDescriptor.H:210
Definition AMReX_Amr.cpp:49
@ make_alias
Definition AMReX_MakeType.H:7
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:314
BoxND< AMREX_SPACEDIM > Box
Definition AMReX_BaseFwd.H:27
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > enclosedCells(const BoxND< dim > &b, int dir) noexcept
Returns a BoxND with CELL based coordinates in direction dir that is enclosed by b....
Definition AMReX_Box.H:1463
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > convert(const BoxND< dim > &b, const IntVectND< dim > &typ) noexcept
Returns a BoxND with different type.
Definition AMReX_Box.H:1435
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 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 constexpr GpuTupleElement< I, GpuTuple< Ts... > >::type & get(GpuTuple< Ts... > &tup) noexcept
Definition AMReX_Tuple.H:179
BottomSolver
Definition AMReX_MLLinOp.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
void Warning(const std::string &msg)
Print out warning message to cerr.
Definition AMReX.cpp:231
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition AMReX.cpp:225
int verbose
Definition AMReX_DistributionMapping.cpp:36
std::unique_ptr< Hypre > makeHypre(const BoxArray &grids, const DistributionMapping &dmap, const Geometry &geom, MPI_Comm comm_, Hypre::Interface interface, const iMultiFab *overset_mask)
Definition AMReX_Hypre.cpp:12
std::array< T, N > Array
Definition AMReX_Array.H:24
Definition AMReX_Array4.H:61
Definition AMReX_FabDataType.H:9
Definition AMReX_Array.H:34
Definition AMReX_TypeTraits.H:29
Definition AMReX_MLLinOp.H:35
LPInfo & setConsolidationRatio(int x) noexcept
Definition AMReX_MLLinOp.H:54
int con_strategy
Definition AMReX_MLLinOp.H:42
bool do_semicoarsening
Definition AMReX_MLLinOp.H:38
bool has_metric_term
Definition AMReX_MLLinOp.H:43
LPInfo & setConsolidationGridSize(int x) noexcept
Definition AMReX_MLLinOp.H:53
int max_semicoarsening_level
Definition AMReX_MLLinOp.H:45
bool hasHiddenDimension() const noexcept
Definition AMReX_MLLinOp.H:62
int con_ratio
Definition AMReX_MLLinOp.H:41
bool do_consolidation
Definition AMReX_MLLinOp.H:37
int con_grid_size
Definition AMReX_MLLinOp.H:40
LPInfo & setSemicoarsening(bool x) noexcept
Definition AMReX_MLLinOp.H:51
LPInfo & setHiddenDirection(int n) noexcept
Definition AMReX_MLLinOp.H:60
LPInfo & setConsolidation(bool x) noexcept
Definition AMReX_MLLinOp.H:50
LPInfo & setSemicoarseningDirection(int n) noexcept
Definition AMReX_MLLinOp.H:59
LPInfo & setMaxSemicoarseningLevel(int n) noexcept
Definition AMReX_MLLinOp.H:58
bool do_agglomeration
Definition AMReX_MLLinOp.H:36
LPInfo & setMetricTerm(bool x) noexcept
Definition AMReX_MLLinOp.H:56
static constexpr int getDefaultConsolidationGridSize()
Definition AMReX_MLLinOp.H:74
LPInfo & setConsolidationStrategy(int x) noexcept
Definition AMReX_MLLinOp.H:55
int max_coarsening_level
Definition AMReX_MLLinOp.H:44
int agg_grid_size
Definition AMReX_MLLinOp.H:39
static constexpr int getDefaultAgglomerationGridSize()
Definition AMReX_MLLinOp.H:66
int hidden_direction
Definition AMReX_MLLinOp.H:47
LPInfo & setAgglomerationGridSize(int x) noexcept
Definition AMReX_MLLinOp.H:52
LPInfo & setMaxCoarseningLevel(int n) noexcept
Definition AMReX_MLLinOp.H:57
LPInfo & setAgglomeration(bool x) noexcept
Definition AMReX_MLLinOp.H:49
int semicoarsening_direction
Definition AMReX_MLLinOp.H:46
Definition AMReX_MLLinOp.H:84
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:66
Definition AMReX_MLLinOp.H:615
void operator=(const CommContainer &)=delete
~CommContainer()
Definition AMReX_MLLinOp.H:622
MPI_Comm comm
Definition AMReX_MLLinOp.H:616
CommContainer(const CommContainer &)=delete
CommContainer(CommContainer &&)=delete
CommContainer(MPI_Comm m) noexcept
Definition AMReX_MLLinOp.H:617