#include <AMReX_MLLinOp.H>
Classes | |
struct | CommContainer |
Public Types | |
using | MFType = MF |
using | FAB = typename FabDataType< MF >::fab_type |
using | RT = typename FabDataType< MF >::value_type |
using | BCType = LinOpBCType |
using | BCMode = LinOpEnumType::BCMode |
using | StateMode = LinOpEnumType::StateMode |
using | Location = LinOpEnumType::Location |
Public Member Functions | |
MLLinOpT ()=default | |
virtual | ~MLLinOpT ()=default |
MLLinOpT (const MLLinOpT< MF > &)=delete | |
MLLinOpT (MLLinOpT< MF > &&)=delete | |
MLLinOpT< MF > & | operator= (const MLLinOpT< MF > &)=delete |
MLLinOpT< MF > & | operator= (MLLinOpT< MF > &&)=delete |
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) |
virtual std::string | name () const |
void | setDomainBC (const Array< BCType, AMREX_SPACEDIM > &lobc, const Array< BCType, AMREX_SPACEDIM > &hibc) noexcept |
Boundary of the whole domain. More... | |
void | setDomainBC (const Vector< Array< BCType, AMREX_SPACEDIM > > &lobc, const Vector< Array< BCType, AMREX_SPACEDIM > > &hibc) noexcept |
Boundary of the whole domain. More... | |
void | setDomainBCLoc (const Array< Real, AMREX_SPACEDIM > &lo_bcloc, const Array< Real, AMREX_SPACEDIM > &hi_bcloc) noexcept |
Set location of domain boundaries. More... | |
bool | needsCoarseDataForBC () const noexcept |
Needs coarse data for bc? More... | |
void | setCoarseFineBC (const MF *crse, int crse_ratio, LinOpBCType bc_type=LinOpBCType::Dirichlet) noexcept |
Set coarse/fine boundary conditions. For cell-centered solves only. More... | |
void | setCoarseFineBC (const MF *crse, IntVect const &crse_ratio, LinOpBCType bc_type=LinOpBCType::Dirichlet) noexcept |
template<typename AMF , std::enable_if_t<!std::is_same_v< MF, AMF >, int > = 0> | |
void | setCoarseFineBC (const AMF *crse, int crse_ratio, LinOpBCType bc_type=LinOpBCType::Dirichlet) noexcept |
template<typename AMF , std::enable_if_t<!std::is_same_v< MF, AMF >, int > = 0> | |
void | setCoarseFineBC (const AMF *crse, IntVect const &crse_ratio, LinOpBCType bc_type=LinOpBCType::Dirichlet) noexcept |
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. More... | |
template<typename AMF , std::enable_if_t<!std::is_same_v< MF, AMF >, int > = 0> | |
void | setLevelBC (int amrlev, const AMF *levelbcdata, const AMF *robinbc_a=nullptr, const AMF *robinbc_b=nullptr, const AMF *robinbc_f=nullptr) |
void | setVerbose (int v) noexcept |
Set verbosity. More... | |
void | setMaxOrder (int o) noexcept |
Set order of interpolation at coarse/fine boundary. More... | |
int | getMaxOrder () const noexcept |
Get order of interpolation at coarse/fine boundary. More... | |
void | setEnforceSingularSolvable (bool o) noexcept |
bool | getEnforceSingularSolvable () const noexcept |
virtual BottomSolver | getDefaultBottomSolver () const |
virtual int | getNComp () const |
Return number of components. More... | |
virtual int | getNGrow (int=0, int=0) const |
virtual bool | needsUpdate () const |
Does it need update if it's reused? More... | |
virtual void | update () |
Update for reuse. More... | |
virtual void | restriction (int amrlev, int cmglev, MF &crse, MF &fine) const =0 |
Restriction onto coarse MG level. More... | |
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. More... | |
virtual void | interpAssign (int amrlev, int fmglev, MF &fine, MF &crse) const |
Overwrite fine MG level data with interpolated coarse data. More... | |
virtual void | interpolationAmr (int famrlev, MF &fine, const MF &crse, IntVect const &nghost) const |
Interpolation between AMR levels. More... | |
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. More... | |
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) More... | |
virtual void | smooth (int amrlev, int mglev, MF &sol, const MF &rhs, bool skip_fillboundary=false) const =0 |
Smooth. More... | |
virtual void | normalize (int amrlev, int mglev, MF &mf) const |
Divide mf by the diagonal component of the operator. Used by bicgstab. More... | |
virtual void | solutionResidual (int amrlev, MF &resid, MF &x, const MF &b, const MF *crse_bcdata=nullptr)=0 |
Compute residual for solution. More... | |
virtual void | prepareForFluxes (int, const MF *=nullptr) |
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) More... | |
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. More... | |
virtual void | compFlux (int, const Array< MF *, AMREX_SPACEDIM > &, MF &, Location) const |
Compute fluxes. More... | |
virtual void | compGrad (int, const Array< MF *, AMREX_SPACEDIM > &, MF &, Location) const |
Compute gradients of the solution. More... | |
virtual void | applyMetricTerm (int, int, MF &) const |
apply metric terms if there are any More... | |
virtual void | unapplyMetricTerm (int, int, MF &) const |
unapply metric terms if there are any More... | |
virtual void | unimposeNeumannBC (int, MF &) const |
This is needed for our nodal projection solver. More... | |
virtual void | applyInhomogNeumannTerm (int, MF &) const |
Extra terms introduced when we treat inhomogeneous Nuemann BC as homogeneous. More... | |
virtual void | applyOverset (int, MF &) const |
for overset solver only More... | |
virtual void | scaleRHS (int, MF &) const |
scale RHS to fix solvability More... | |
virtual Vector< RT > | getSolvabilityOffset (int, int, MF const &) const |
get offset for fixing solvability More... | |
virtual void | fixSolvabilityByOffset (int, int, MF &, Vector< RT > const &) const |
fix solvability by subtracting offset from RHS More... | |
virtual void | prepareForSolve ()=0 |
virtual void | prepareForGMRES () |
virtual void | setDirichletNodesToZero (int, int, MF &) const |
virtual bool | isSingular (int amrlev) const =0 |
Is it singular on given AMR level? More... | |
virtual bool | isBottomSingular () const =0 |
Is the bottom of MG singular? More... | |
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 More... | |
virtual std::unique_ptr< MLLinOpT< MF > > | makeNLinOp (int) const |
virtual void | getFluxes (const Vector< Array< MF *, AMREX_SPACEDIM > > &, const Vector< MF * > &, Location) const |
virtual void | getFluxes (const Vector< MF * > &, const Vector< MF * > &) const |
virtual bool | supportNSolve () const |
virtual void | copyNSolveSolution (MF &, MF const &) const |
virtual void | postSolve (Vector< MF > &) const |
virtual RT | normInf (int amrlev, MF const &mf, bool local) const =0 |
virtual void | averageDownAndSync (Vector< MF > &sol) const =0 |
virtual void | avgDownResAmr (int clev, MF &cres, MF const &fres) const |
virtual void | avgDownResMG (int clev, MF &cres, MF const &fres) const |
bool | isMFIterSafe (int amrlev, int mglev1, int mglev2) const |
int | NAMRLevels () const noexcept |
Return the number of AMR levels. More... | |
int | NMGLevels (int amrlev) const noexcept |
Return the number of MG levels at given AMR level. More... | |
const Geometry & | Geom (int amr_lev, int mglev=0) const noexcept |
Public Attributes | |
Vector< Array< BCType, AMREX_SPACEDIM > > | m_lobc |
Vector< Array< BCType, AMREX_SPACEDIM > > | m_hibc |
Vector< Array< BCType, AMREX_SPACEDIM > > | m_lobc_orig |
Vector< Array< BCType, AMREX_SPACEDIM > > | m_hibc_orig |
Protected Member Functions | |
const Vector< int > & | AMRRefRatio () const noexcept |
Return AMR refinement ratios. More... | |
int | AMRRefRatio (int amr_lev) const noexcept |
Return AMR refinement ratio at given AMR level. More... | |
FabFactory< FAB > const * | Factory (int amr_lev, int mglev=0) const noexcept |
GpuArray< BCType, AMREX_SPACEDIM > | LoBC (int icomp=0) const noexcept |
GpuArray< BCType, AMREX_SPACEDIM > | HiBC (int icomp=0) const noexcept |
bool | hasBC (BCType bct) const noexcept |
bool | hasInhomogNeumannBC () const noexcept |
bool | hasRobinBC () const noexcept |
virtual bool | supportRobinBC () const noexcept |
virtual bool | supportInhomogNeumannBC () const noexcept |
bool | isBottomActive () const noexcept |
MPI_Comm | BottomCommunicator () const noexcept |
MPI_Comm | Communicator () const noexcept |
void | setCoarseFineBCLocation (const RealVect &cloc) noexcept |
bool | doAgglomeration () const noexcept |
bool | doConsolidation () const noexcept |
bool | doSemicoarsening () const noexcept |
bool | isCellCentered () const noexcept |
virtual IntVect | getNGrowVectRestriction () const |
virtual void | make (Vector< Vector< MF > > &mf, IntVect const &ng) const |
virtual MF | make (int amrlev, int mglev, IntVect const &ng) const |
virtual MF | makeAlias (MF const &mf) const |
virtual MF | makeCoarseMG (int amrlev, int mglev, IntVect const &ng) const |
virtual MF | makeCoarseAmr (int famrlev, IntVect const &ng) const |
virtual std::unique_ptr< FabFactory< FAB > > | makeFactory (int, int) const |
virtual void | resizeMultiGrid (int new_size) |
bool | hasHiddenDimension () const noexcept |
int | hiddenDirection () const noexcept |
Box | compactify (Box const &b) const noexcept |
template<typename T > | |
Array4< T > | compactify (Array4< T > const &a) const noexcept |
template<typename T > | |
T | get_d0 (T const &dx, T const &dy, T const &) const noexcept |
template<typename T > | |
T | get_d1 (T const &, T const &dy, T const &dz) const noexcept |
Static Protected Attributes | |
static constexpr int | mg_coarsen_ratio = 2 |
static constexpr int | mg_box_min_width = 2 |
Private Member Functions | |
void | defineGrids (const Vector< Geometry > &a_geom, const Vector< BoxArray > &a_grids, const Vector< DistributionMapping > &a_dmap, const Vector< FabFactory< FAB > const * > &a_factory) |
void | defineBC () |
MPI_Comm | makeSubCommunicator (const DistributionMapping &dm) |
virtual void | checkPoint (std::string const &) const |
Static Private Member Functions | |
static void | makeAgglomeratedDMap (const Vector< BoxArray > &ba, Vector< DistributionMapping > &dm) |
static void | makeConsolidatedDMap (const Vector< BoxArray > &ba, Vector< DistributionMapping > &dm, int ratio, int strategy) |
Private Attributes | |
Vector< std::unique_ptr< MF > > | levelbc_raii |
Vector< std::unique_ptr< MF > > | robin_a_raii |
Vector< std::unique_ptr< MF > > | robin_b_raii |
Vector< std::unique_ptr< MF > > | robin_f_raii |
Friends | |
template<typename T > | |
class | MLMGT |
template<typename T > | |
class | MLCGSolverT |
template<typename T > | |
class | MLPoissonT |
template<typename T > | |
class | MLABecLaplacianT |
template<typename T > | |
class | GMRESMLMGT |
using amrex::MLLinOpT< MF >::BCMode = LinOpEnumType::BCMode |
using amrex::MLLinOpT< MF >::BCType = LinOpBCType |
using amrex::MLLinOpT< MF >::FAB = typename FabDataType<MF>::fab_type |
using amrex::MLLinOpT< MF >::Location = LinOpEnumType::Location |
using amrex::MLLinOpT< MF >::MFType = MF |
using amrex::MLLinOpT< MF >::RT = typename FabDataType<MF>::value_type |
using amrex::MLLinOpT< MF >::StateMode = LinOpEnumType::StateMode |
|
default |
|
virtualdefault |
|
delete |
|
delete |
|
inlineprotectednoexcept |
Return AMR refinement ratios.
|
inlineprotectednoexcept |
Return AMR refinement ratio at given AMR level.
|
pure virtual |
Apply the linear operator, out = L(in)
amrlev | AMR level |
mglev | MG level |
out | output |
in | input |
bc_mode | Is the BC homogeneous or inhomogeneous? |
s_mode | Are data data solution or correction? |
bndry | object for handling coarse/fine and physical boundaries |
Implemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
Extra terms introduced when we treat inhomogeneous Nuemann BC as homogeneous.
Reimplemented in amrex::MLCellABecLapT< MF >.
|
inlinevirtual |
apply metric terms if there are any
Reimplemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
for overset solver only
Reimplemented in amrex::MLCellABecLapT< MF >.
|
pure virtual |
Implemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
Average-down data from fine AMR level to coarse AMR level.
camrlev | coarse AMR level |
crse_sol | solutoin on coarse AMR level |
crse_rhs | RHS on coarse AMR level |
fine_sol | solution on fine AMR level |
fine_rhs | RHS on fine AMR level |
Reimplemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
Reimplemented in amrex::MLCellLinOpT< MF >.
|
virtual |
|
inlineprotectednoexcept |
|
inlineprivatevirtual |
Reimplemented in amrex::MLNodeLaplacian.
|
inlineprotectednoexcept |
|
inlineprotectednoexcept |
|
protectednoexcept |
|
inlinevirtual |
Compute fluxes.
amrlev | AMR level |
fluxes | fluxes |
sol | solution |
loc | location of the fluxes |
Reimplemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
Compute gradients of the solution.
amrlev | AMR level |
grad | grad(sol) |
sol | solution |
loc | location of the gradients |
Reimplemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
Reimplemented in amrex::MLPoissonT< MF >, and amrex::MLABecLaplacianT< MF >.
|
pure virtual |
Compute residual for the residual-correction form, resid = b - L(x)
amrlev | AMR level |
mglev | MG level |
resid | residual |
x | unknown in the residual-correction form |
b | RHS in the residual-correction form |
bc_mode | Is the BC homogeneous or inhomogeneous? |
crse_bc_data | optional argument providing BC at coarse/fine boundary. |
Implemented in amrex::MLCellLinOpT< MF >.
void amrex::MLLinOpT< MF >::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 |
||
) |
|
private |
|
private |
|
inlineprotectednoexcept |
|
inlineprotectednoexcept |
|
inlineprotectednoexcept |
|
inlineprotectednoexcept |
|
inlinevirtual |
fix solvability by subtracting offset from RHS
|
inlinenoexcept |
|
inlineprotectednoexcept |
|
inlineprotectednoexcept |
|
inlinevirtual |
Reimplemented in amrex::MLNodeLaplacian.
|
inlinenoexcept |
Get the flag for whether the solver should try to make singular problem solvable.
|
inlinevirtual |
Reimplemented in amrex::MLCellABecLapT< MF >.
|
inlinevirtual |
Reimplemented in amrex::MLCellABecLapT< MF >.
|
inlinenoexcept |
Get order of interpolation at coarse/fine boundary.
|
inlinevirtual |
Return number of components.
Reimplemented in amrex::MLEBABecLap, amrex::MLALaplacianT< MF >, amrex::MLABecLaplacianT< MF >, amrex::MLTensorOp, and amrex::MLEBTensorOp.
|
inlinevirtual |
|
inlineprotectedvirtual |
Reimplemented in amrex::MLCurlCurl.
|
inlinevirtual |
get offset for fixing solvability
Reimplemented in amrex::MLCellLinOpT< MF >.
|
protectednoexcept |
|
inlineprotectednoexcept |
|
protectednoexcept |
|
protectednoexcept |
|
inlineprotectednoexcept |
|
inlineprotectednoexcept |
|
inlinevirtual |
Overwrite fine MG level data with interpolated coarse data.
amrlev | AMR level |
fmglev | fine MG level |
fine | fine MG level data |
crse | coarse MG level data |
Reimplemented in amrex::MLCellLinOpT< MF >.
|
pure virtual |
Add interpolated coarse MG level data to fine MG level data.
amrlev | AMR level |
fmglev | fine MG level |
crse | coarse data. |
fine | fine data. |
Implemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
Interpolation between AMR levels.
famrlev | fine AMR level |
fine | fine level data |
crse | coarse level data |
nghost | number of ghost cells |
Reimplemented in amrex::MLCellLinOpT< MF >.
|
inlineprotectednoexcept |
|
pure virtual |
Is the bottom of MG singular?
Implemented in amrex::MLNodeLinOp, amrex::MLEBABecLap, amrex::MLCurlCurl, amrex::MLABecLaplacianT< MF >, amrex::MLTensorOp, amrex::MLPoissonT< MF >, amrex::MLNodeABecLaplacian, amrex::MLEBTensorOp, amrex::MLEBNodeFDLaplacian, and amrex::MLALaplacianT< MF >.
|
inlineprotectednoexcept |
bool amrex::MLLinOpT< MF >::isMFIterSafe | ( | int | amrlev, |
int | mglev1, | ||
int | mglev2 | ||
) | const |
|
pure virtual |
Is it singular on given AMR level?
Implemented in amrex::MLCurlCurl, amrex::MLTensorOp, amrex::MLNodeABecLaplacian, amrex::MLEBTensorOp, amrex::MLEBNodeFDLaplacian, amrex::MLNodeLinOp, amrex::MLEBABecLap, amrex::MLABecLaplacianT< MF >, amrex::MLPoissonT< MF >, and amrex::MLALaplacianT< MF >.
|
inlineprotectednoexcept |
|
protectedvirtual |
Reimplemented in amrex::MLCurlCurl.
|
protectedvirtual |
|
staticprivate |
|
protectedvirtual |
|
protectedvirtual |
Reimplemented in amrex::MLCurlCurl.
|
protectedvirtual |
Reimplemented in amrex::MLCurlCurl.
|
staticprivate |
|
inlineprotectedvirtual |
Reimplemented in amrex::MLEBABecLap.
|
inlinevirtual |
Reimplemented in amrex::MLEBABecLap, amrex::MLALaplacianT< MF >, amrex::MLABecLaplacianT< MF >, and amrex::MLPoissonT< MF >.
|
private |
|
inlinevirtual |
Reimplemented in amrex::MLNodeTensorLaplacian, amrex::MLNodeLaplacian, amrex::MLNodeABecLaplacian, amrex::MLEBNodeFDLaplacian, and amrex::MLCurlCurl.
|
inlinenoexcept |
Return the number of AMR levels.
|
inlinenoexcept |
Needs coarse data for bc?
If the lowest level grids does not cover the entire domain, coarse level data are needed for supplying Dirichlet bc at coarse/fine boundary, even when the domain bc is not Dirichlet.
|
inlinevirtual |
Does it need update if it's reused?
Reimplemented in amrex::MLEBABecLap, amrex::MLCellLinOpT< MF >, amrex::MLCellABecLapT< MF >, amrex::MLALaplacianT< MF >, amrex::MLABecLaplacianT< MF >, amrex::MLTensorOp, amrex::MLNodeABecLaplacian, and amrex::MLEBTensorOp.
|
inlinenoexcept |
Return the number of MG levels at given AMR level.
|
inlinevirtual |
Divide mf by the diagonal component of the operator. Used by bicgstab.
Reimplemented in amrex::MLPoissonT< MF >, amrex::MLALaplacianT< MF >, and amrex::MLABecLaplacianT< MF >.
|
pure virtual |
Implemented in amrex::MLCellLinOpT< MF >.
|
delete |
|
delete |
|
inlinevirtual |
|
inlinevirtual |
Reimplemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
Reimplemented in amrex::MLNodeLinOp.
|
pure virtual |
Implemented in amrex::MLNodeLinOp, amrex::MLEBABecLap, amrex::MLCurlCurl, amrex::MLCellLinOpT< MF >, amrex::MLCellABecLapT< MF >, amrex::MLABecLaplacianT< MF >, amrex::MLTensorOp, amrex::MLPoissonT< MF >, amrex::MLNodeTensorLaplacian, amrex::MLNodeLaplacian, amrex::MLNodeABecLaplacian, amrex::MLEBTensorOp, amrex::MLEBNodeFDLaplacian, and amrex::MLALaplacianT< MF >.
|
inlinevirtual |
Reflux at AMR coarse/fine boundary.
crse_amrlev | coarse AMR level |
res | coarse level residual |
crse_sol | coarse level solution |
crse_rhs | coarse level RHS |
fine_res | fine level residual |
fine_sol | fine level solution |
fine_rhs | fine level RHS |
Reimplemented in amrex::MLCellLinOpT< MF >.
|
protectedvirtual |
Reimplemented in amrex::MLNodeLinOp, and amrex::MLNodeLaplacian.
|
pure virtual |
Restriction onto coarse MG level.
amrlev | AMR level |
cmglev | coarse MG level |
crse | coarse data. This is the output. |
fine | fine data. This is the input. Some operators might need to fill ghost cells. This is why it's not a const reference. |
Implemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
scale RHS to fix solvability
|
noexcept |
|
noexcept |
|
noexcept |
Set coarse/fine boundary conditions. For cell-centered solves only.
If we want to do a linear solve where the boundary conditions on the coarsest AMR level of the solve come from a coarser level (e.g. the base AMR level of the solve is > 0 and does not cover the entire domain), we must explicitly provide the coarser data. Boundary conditions from a coarser level are Dirichlet by default. The MultiFab crse does not need to have ghost cells and is at a coarser resolution than the coarsest AMR level of the solve; it is used to supply (interpolated) boundary conditions for the solve. NOTE: If this is called, it must be called before setLevelBC
. If crse is nullptr, then the bc values are assumed to be zero. The coarse/fine BC type can be changed to homogeneous Neumann by the bc_type argument. In that case, use nullptr for the crse argument.
crse | the coarse AMR level data |
crse_ratio | the coarsening ratio between fine and coarse AMR levels. |
bc_type | optional. It's Dirichlet by default, and can be Neumann. |
|
noexcept |
|
inlineprotectednoexcept |
|
inlinevirtual |
For GMRES to work, this might need to be implemented to mask out Dirichlet nodes or cells (e.g., EB covered cells or overset cells)
Reimplemented in amrex::MLCellABecLapT< MF >.
|
noexcept |
Boundary of the whole domain.
This functions must be called, and must be called before other bc functions. This version is for single-component solve or when all the components have the same BC types.
lobc | lower boundaries |
hibc | upper boundaries |
|
noexcept |
Boundary of the whole domain.
This functions must be called, and must be called before other bc functions. This version is for multi-component solve.
lobc | lower boundaries |
hibc | upper boundaries |
|
noexcept |
Set location of domain boundaries.
By default, domain BC is on the domain face. If that's the case, this function doesn't need to be called. However, one could use this function to set non-zero domain BC locations. Note all values should be >= 0. If this function is called, it MUST be called before setLevelBC.
|
inlinenoexcept |
Set the flag for whether the solver should try to make singular problem solvable, which is on by default.
void amrex::MLLinOpT< MF >::setLevelBC | ( | int | amrlev, |
const AMF * | levelbcdata, | ||
const AMF * | robinbc_a = nullptr , |
||
const AMF * | robinbc_b = nullptr , |
||
const AMF * | robinbc_f = nullptr |
||
) |
|
pure virtual |
Set boundary conditions for given level. For cell-centered solves only.
This must be called for each level. Argument levelbcdata
is used to supply Dirichlet or Neumann bc at the physical domain; if those data are homogeneous we can pass nullptr instead of levelbcdata. Regardless, this function must be called. If used, the MultiFab levelbcdata must have one ghost cell. Only the data outside the physical domain will be used. It is assumed that the data in those ghost cells outside the domain live exactly on the face of the physical domain. Argument amrlev
is relative level such that the lowest to the solver is always 0. The optional arguments robinbc_[a|b|f] provide Robin boundary condition a*phi + b*dphi/dn = f
. Note that d./dn
is d./dx
at the upper boundary and -d./dx
at the lower boundary, for Robin BC. However, for inhomogeneous Neumann BC, the value in leveldata is assumed to be d./dx
.
Implemented in amrex::MLCellLinOpT< MF >.
|
inlinenoexcept |
Set order of interpolation at coarse/fine boundary.
|
inlinenoexcept |
Set verbosity.
|
pure virtual |
Smooth.
amrlev | AMR level |
mglev | MG level |
sol | unknowns |
rhs | RHS |
skip_fillboundary | flag controlling whether ghost cell filling can be skipped. |
Implemented in amrex::MLCellLinOpT< MF >.
|
pure virtual |
Compute residual for solution.
amrlev | AMR level |
resid | residual |
x | solution |
b | RHS |
crse_bc_data | optional argument providing BC at coarse/fine boundary. |
Implemented in amrex::MLCellLinOpT< MF >.
|
inlineprotectedvirtualnoexcept |
Reimplemented in amrex::MLCellABecLapT< MF >.
|
inlinevirtual |
Reimplemented in amrex::MLPoissonT< MF >, and amrex::MLABecLaplacianT< MF >.
|
inlineprotectedvirtualnoexcept |
Reimplemented in amrex::MLEBABecLap, and amrex::MLABecLaplacianT< MF >.
|
inlinevirtual |
unapply metric terms if there are any
Reimplemented in amrex::MLCellLinOpT< MF >.
|
inlinevirtual |
This is needed for our nodal projection solver.
|
inlinevirtual |
Update for reuse.
Reimplemented in amrex::MLEBABecLap, amrex::MLCellLinOpT< MF >, amrex::MLCellABecLapT< MF >, amrex::MLALaplacianT< MF >, amrex::MLABecLaplacianT< MF >, amrex::MLTensorOp, amrex::MLNodeABecLaplacian, and amrex::MLEBTensorOp.
|
pure virtual |
x dot y, used by the bottom solver
Implemented in amrex::MLCellLinOpT< MF >.
|
protected |
|
protected |
|
private |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
first Vector is for amr level and second is mg level
|
protected |
Vector<Array<BCType, AMREX_SPACEDIM> > amrex::MLLinOpT< MF >::m_hibc |
Vector<Array<BCType, AMREX_SPACEDIM> > amrex::MLLinOpT< MF >::m_hibc_orig |
|
protected |
Vector<Array<BCType, AMREX_SPACEDIM> > amrex::MLLinOpT< MF >::m_lobc |
Vector<Array<BCType, AMREX_SPACEDIM> > amrex::MLLinOpT< MF >::m_lobc_orig |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
staticconstexprprotected |
|
staticconstexprprotected |
|
protected |
|
protected |
|
private |
|
private |
|
private |
|
protected |