Chombo + EB  3.2
Public Member Functions | Public Attributes | Static Public Attributes | Protected Member Functions | Protected Attributes | List of all members
AMRPoissonOp Class Reference

#include <AMRPoissonOp.H>

Inheritance diagram for AMRPoissonOp:
Inheritance graph
[legend]

Public Member Functions

virtual Real dx () const
 
AMRPoissonOp functions
 AMRPoissonOp ()
 
virtual ~AMRPoissonOp ()
 
void define (const DisjointBoxLayout &a_grids, const DisjointBoxLayout &a_gridsFiner, const DisjointBoxLayout &a_gridsCoarser, const Real &a_dxLevel, int a_refRatio, int a_refRatioFiner, const ProblemDomain &a_domain, BCHolder a_bc, const Copier &a_exchange, const CFRegion &a_cfregion, const int a_nComp=1)
 
void define (const DisjointBoxLayout &a_grids, const DisjointBoxLayout &a_gridsFiner, const Real &a_dxLevel, int a_refRatio, int a_refRatioFiner, const ProblemDomain &a_domain, BCHolder a_bc, const Copier &a_exchange, const CFRegion &a_cfregion, const int a_nComp=1)
 
void define (const DisjointBoxLayout &a_grids, const DisjointBoxLayout &a_baseBAPtr, const Real &a_dxLevel, int a_refRatio, const ProblemDomain &a_domain, BCHolder a_bc, const Copier &a_exchange, const CFRegion &a_cfregion, int a_numComp=1)
 
void define (const DisjointBoxLayout &a_grids, const Real &a_dx, const ProblemDomain &a_domain, BCHolder a_bc, const Copier &a_exchange, const CFRegion &a_cfregion)
 
void define (const DisjointBoxLayout &a_grids, const Real &a_dx, const ProblemDomain &a_domain, BCHolder a_bc)
 
void define (const DisjointBoxLayout &a_grids, const DisjointBoxLayout *a_baseBAPtr, Real a_dxLevel, int a_refRatio, const ProblemDomain &a_domain, BCHolder a_bc)
 Full define function that mimics the old PoissonOp. More...
 
virtual void residual (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs, bool a_homogeneous=false)
 
virtual void residualNF (LevelData< FArrayBox > &a_lhs, LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > *a_phiCoarse, const LevelData< FArrayBox > &a_rhs, bool a_homogeneous=false)
 Residual which uses coarse-level boundary conditions, but ignores finer levels. More...
 
virtual void residualI (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs, bool a_homogeneous=false)
 despite what you might think, the "I" here means "Ignore the coarse-fine boundary" More...
 
virtual void preCond (LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_residual)
 
virtual void applyOp (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, bool a_homogeneous=false)
 
virtual void applyOpI (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, bool a_homogeneous=false)
 despite what you might think, the "I" here means "Ignore the coarse-fine boundary" More...
 
virtual void applyOpNoBoundary (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
 
virtual void create (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
 
virtual void createCoarsened (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs, const int &a_refRat)
 
virtual void assign (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
 
virtual void assignLocal (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
 
virtual void buildCopier (Copier &a_copier, const LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
 
virtual void assignCopier (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs, const Copier &a_copier)
 
virtual void zeroCovered (LevelData< FArrayBox > &a_lhs, LevelData< FArrayBox > &a_rhs, const Copier &a_copier)
 
virtual Real dotProduct (const LevelData< FArrayBox > &a_1, const LevelData< FArrayBox > &a_2)
 
virtual void mDotProduct (const LevelData< FArrayBox > &a_1, const int a_sz, const LevelData< FArrayBox > a_2[], Real a_mdots[])
 
virtual void incr (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_x, Real a_scale)
 
virtual void axby (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_x, const LevelData< FArrayBox > &a_y, Real a_a, Real a_b)
 
virtual void scale (LevelData< FArrayBox > &a_lhs, const Real &a_scale)
 
virtual Real norm (const LevelData< FArrayBox > &a_x, int a_ord)
 
virtual Real localMaxNorm (const LevelData< FArrayBox > &a_x)
 
virtual void setToZero (LevelData< FArrayBox > &a_x)
 
MGLevelOp functions
virtual void relax (LevelData< FArrayBox > &a_e, const LevelData< FArrayBox > &a_residual, int a_iterations)
 
virtual void relaxNF (LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > *a_phiCoarse, const LevelData< FArrayBox > &a_rhs, int a_iterations)
 Do CF interpolation then relax as normal. More...
 
virtual void createCoarser (LevelData< FArrayBox > &a_coarse, const LevelData< FArrayBox > &a_fine, bool a_ghosted)
 
virtual void restrictResidual (LevelData< FArrayBox > &a_resCoarse, LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_rhsFine)
 
virtual void prolongIncrement (LevelData< FArrayBox > &a_phiThisLevel, const LevelData< FArrayBox > &a_correctCoarse)
 
AMRLevelOp functions
virtual int refToCoarser ()
 
virtual void AMRResidual (LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, const LevelData< FArrayBox > &a_rhs, bool a_homogeneousPhysBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
virtual void AMRResidualNC (LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs, bool a_homogeneousPhysBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
virtual void AMRResidualNF (LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, const LevelData< FArrayBox > &a_rhs, bool a_homogeneousPhysBC)
 
virtual void AMROperator (LevelData< FArrayBox > &a_LofPhi, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, bool a_homogeneousDomBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
virtual void AMROperatorNC (LevelData< FArrayBox > &a_LofPhi, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, bool a_homogeneousBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
virtual void AMROperatorNF (LevelData< FArrayBox > &a_LofPhi, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, bool a_homogeneousBC)
 
virtual void AMRRestrict (LevelData< FArrayBox > &a_resCoarse, const LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection, bool a_skip_res=false)
 
virtual void AMRRestrictS (LevelData< FArrayBox > &a_resCoarse, const LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection, LevelData< FArrayBox > &a_scratch, bool a_skip_res=false)
 
virtual void AMRProlong (LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection)
 
virtual void AMRProlongS (LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection, LevelData< FArrayBox > &a_temp, const Copier &a_copier)
 
virtual void AMRProlongS_2 (LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection, LevelData< FArrayBox > &a_temp, const Copier &a_copier, const Copier &a_cornerCopier, const AMRLevelOp< LevelData< FArrayBox > > *a_crsOp)
 
virtual void AMRUpdateResidual (LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection)
 
virtual Real AMRNorm (const LevelData< FArrayBox > &a_coarseResid, const LevelData< FArrayBox > &a_fineResid, const int &a_refRat, const int &a_ord)
 
virtual void setAlphaAndBeta (const Real &a_alpha, const Real &a_beta)
 For tga to reset stuff. More...
 
virtual void setBC (const BCHolder &a_bc)
 Change boundary conditions. More...
 
virtual void diagonalScale (LevelData< FArrayBox > &a_rhs, bool a_kappaWeighted)
 For tga stuff—in this case a noop. More...
 
virtual void divideByIdentityCoef (LevelData< FArrayBox > &a_rhs)
 For tga stuff—in this case a noop. More...
 
virtual void fillGrad (const LevelData< FArrayBox > &a_phi)
 These functions are part of the LevelTGA interface...... More...
 
virtual void reflux (const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, LevelData< FArrayBox > &a_residual, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
virtual void getFlux (FluxBox &a_flux, const LevelData< FArrayBox > &a_data, const Box &a_grid, const DataIndex &a_dit, Real a_scale)
 
virtual void write (const LevelData< FArrayBox > *a_data, const char *a_filename)
 
- Public Member Functions inherited from LevelTGAHelmOp< LevelData< FArrayBox >, FluxBox >
 LevelTGAHelmOp ()
 
 LevelTGAHelmOp (bool a_isTimeIndependent)
 
virtual ~LevelTGAHelmOp ()
 Destructor. More...
 
- Public Member Functions inherited from TGAHelmOp< LevelData< FArrayBox > >
 TGAHelmOp ()
 
 TGAHelmOp (bool a_isTimeDependent)
 
virtual ~TGAHelmOp ()
 Destructor. More...
 
virtual void diagonalScale (LevelData< FArrayBox > &a_rhs)
 
virtual void kappaScale (LevelData< FArrayBox > &a_rhs)
 for eb only. kappa weight the rhs but do not multiply by the identity coefficient More...
 
virtual void setTime (Real a_time)
 
virtual void setTime (Real a_oldTime, Real a_mu, Real a_dt)
 
bool isTimeDependent () const
 Returns true if the operator is time-dependent, false otherwise. More...
 
virtual LevelData< FArrayBox > & identityCoef ()
 Allows access to the identity coefficient data for the operator. More...
 
- Public Member Functions inherited from AMRLevelOp< LevelData< FArrayBox > >
virtual void dumpAMR (Vector< LevelData< FArrayBox > * > &a_data, string name)
 
virtual void dumpLevel (LevelData< FArrayBox > &a_data, string name)
 
 AMRLevelOp ()
 Constructor. More...
 
virtual void dumpStuff (Vector< LevelData< FArrayBox > * > data, string filename)
 
virtual ~AMRLevelOp ()
 Destructor. More...
 
virtual void outputLevel (LevelData< FArrayBox > &a_rhs, string &a_name)
 
virtual void outputAMR (Vector< LevelData< FArrayBox > * > &a_rhs, string &a_name)
 
virtual unsigned int orderOfAccuracy (void) const
 
virtual void enforceCFConsistency (LevelData< FArrayBox > &a_coarseCorrection, const LevelData< FArrayBox > &a_correction)
 This routine is for operators with orderOfAccuracy()>2. More...
 
- Public Member Functions inherited from MGLevelOp< LevelData< FArrayBox > >
 MGLevelOp ()
 Constructor. More...
 
virtual ~MGLevelOp ()
 Destructor. More...
 
virtual void restrictResidual (LevelData< FArrayBox > &a_resCoarse, LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > *a_phiCoarse, const LevelData< FArrayBox > &a_rhsFine, bool homogeneous)
 full-multigrid version of restrictResidual, useful for FAS-type schemes. Defaults to standard restriction More...
 
virtual void restrictR (LevelData< FArrayBox > &a_phiCoarse, const LevelData< FArrayBox > &a_phiFine)
 simple restriction function More...
 
void addObserver (MGLevelOpObserver< LevelData< FArrayBox > > *a_observer)
 
virtual void applyOpMg (LevelData< FArrayBox > &a_lhs, LevelData< FArrayBox > &a_phi, LevelData< FArrayBox > *a_phiCoarse, bool a_homogeneous)
 Apply an operator. More...
 
void removeObserver (MGLevelOpObserver< LevelData< FArrayBox > > *a_observer)
 
void addCoarserObserver (MGLevelOp< LevelData< FArrayBox > > *a_operator, int a_coarseningFactor)
 
void notifyObserversOfChange ()
 This should be called whenever the operator's data is updated. More...
 
virtual void finerOperatorChanged (const MGLevelOp< LevelData< FArrayBox > > &a_operator, int a_coarseningFactor)
 
int numObservers () const
 Returns the number of objects observing this operator. More...
 
- Public Member Functions inherited from LinearOp< LevelData< FArrayBox > >
virtual ~LinearOp ()
 
virtual void clear (LevelData< FArrayBox > &a_lhs)
 
- Public Member Functions inherited from MGLevelOpObserver< LevelData< FArrayBox > >
 MGLevelOpObserver ()
 Base level Constructor. Called by all subclasses. More...
 
virtual ~MGLevelOpObserver ()
 Destructor. More...
 
virtual void operatorChanged (const MGLevelOp< LevelData< FArrayBox > > &a_operator)
 
void setObservee (MGLevelOp< LevelData< FArrayBox > > *a_observee)
 
void clearObservee ()
 

Public Attributes

Real m_alpha
 public constants More...
 
Real m_beta
 
Real m_aCoef
 
Real m_bCoef
 
Real m_dxCrse
 
Vector< IntVectm_colors
 

Static Public Attributes

static int s_exchangeMode
 
static int s_relaxMode
 
static int s_maxCoarse
 
static int s_prolongType
 

Protected Member Functions

virtual void levelGSRB (LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs)
 
virtual void levelMultiColor (LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs)
 
virtual void looseGSRB (LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs)
 
virtual void overlapGSRB (LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs)
 
virtual void levelGSRBLazy (LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs)
 
virtual void levelJacobi (LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs)
 
virtual void homogeneousCFInterp (LevelData< FArrayBox > &a_phif)
 
virtual void homogeneousCFInterp (LevelData< FArrayBox > &a_phif, const DataIndex &a_datInd, int a_idir, Side::LoHiSide a_hiorlo)
 
virtual void singleBoxCFInterp (FArrayBox &a_phi)
 
virtual void interpOnIVSHomo (LevelData< FArrayBox > &a_phif, const DataIndex &a_datInd, const int a_idir, const Side::LoHiSide a_hiorlo, const IntVectSet &a_interpIVS)
 
virtual void getFlux (FArrayBox &a_flux, const FArrayBox &a_data, const Box &a_edgebox, int a_dir, int a_ref=1) const
 
virtual void getFlux (FArrayBox &a_flux, const FArrayBox &a_data, int a_dir, int a_ref=1) const
 

Protected Attributes

Real m_dx
 
ProblemDomain m_domain
 
LevelDataOps< FArrayBoxm_levelOps
 
BCHolder m_bc
 
CFRegion m_cfregion
 
Copier m_exchangeCopier
 
QuadCFInterp m_interpWithCoarser
 
LevelFluxRegister m_levfluxreg
 
DisjointBoxLayout m_coarsenedMGrids
 
int m_refToCoarser
 
int m_refToFiner
 

Detailed Description

Operator for solving (alpha I + beta*Laplacian)(phi) = rho over an AMR hierarchy.

Constructor & Destructor Documentation

◆ AMRPoissonOp()

AMRPoissonOp::AMRPoissonOp ( )
inline

◆ ~AMRPoissonOp()

virtual AMRPoissonOp::~AMRPoissonOp ( )
inlinevirtual

Member Function Documentation

◆ define() [1/6]

void AMRPoissonOp::define ( const DisjointBoxLayout a_grids,
const DisjointBoxLayout a_gridsFiner,
const DisjointBoxLayout a_gridsCoarser,
const Real a_dxLevel,
int  a_refRatio,
int  a_refRatioFiner,
const ProblemDomain a_domain,
BCHolder  a_bc,
const Copier a_exchange,
const CFRegion a_cfregion,
const int  a_nComp = 1 
)

full define function for AMRLevelOp with both coarser and finer levels

Referenced by ~AMRPoissonOp(), AMRPoissonOpFactory::~AMRPoissonOpFactory(), and VCAMRPoissonOp2Factory::~VCAMRPoissonOp2Factory().

◆ define() [2/6]

void AMRPoissonOp::define ( const DisjointBoxLayout a_grids,
const DisjointBoxLayout a_gridsFiner,
const Real a_dxLevel,
int  a_refRatio,
int  a_refRatioFiner,
const ProblemDomain a_domain,
BCHolder  a_bc,
const Copier a_exchange,
const CFRegion a_cfregion,
const int  a_nComp = 1 
)

full define function for AMRLevelOp with finer levels, but no coarser

◆ define() [3/6]

void AMRPoissonOp::define ( const DisjointBoxLayout a_grids,
const DisjointBoxLayout a_baseBAPtr,
const Real a_dxLevel,
int  a_refRatio,
const ProblemDomain a_domain,
BCHolder  a_bc,
const Copier a_exchange,
const CFRegion a_cfregion,
int  a_numComp = 1 
)

define function for AMRLevelOp which has no finer AMR level Jamie Parkinson added a_numComp parameter so that we define CF interp objects correctly

◆ define() [4/6]

void AMRPoissonOp::define ( const DisjointBoxLayout a_grids,
const Real a_dx,
const ProblemDomain a_domain,
BCHolder  a_bc,
const Copier a_exchange,
const CFRegion a_cfregion 
)

define function for AMRLevelOp which has no finer or coarser AMR level

◆ define() [5/6]

void AMRPoissonOp::define ( const DisjointBoxLayout a_grids,
const Real a_dx,
const ProblemDomain a_domain,
BCHolder  a_bc 
)

define function for AMRLevelOp which has no finer or coarser AMR level

◆ define() [6/6]

void AMRPoissonOp::define ( const DisjointBoxLayout a_grids,
const DisjointBoxLayout a_baseBAPtr,
Real  a_dxLevel,
int  a_refRatio,
const ProblemDomain a_domain,
BCHolder  a_bc 
)

Full define function that mimics the old PoissonOp.

Makes all coarse-fine information and sets internal variables

◆ residual()

virtual void AMRPoissonOp::residual ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneous = false 
)
virtual

Say you are solving L(phi) = rhs. Make a_lhs = L(a_phi) - a_rhs. If a_homogeneous is true, evaluate the operator using homogeneous boundary conditions.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ residualNF()

virtual void AMRPoissonOp::residualNF ( LevelData< FArrayBox > &  a_lhs,
LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > *  a_phiCoarse,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneous = false 
)
virtual

Residual which uses coarse-level boundary conditions, but ignores finer levels.

useful for full-multigrid types of schemes

Reimplemented from MGLevelOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ residualI()

virtual void AMRPoissonOp::residualI ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneous = false 
)
virtual

despite what you might think, the "I" here means "Ignore the coarse-fine boundary"

Reimplemented in VCAMRPoissonOp2.

Referenced by ~AMRPoissonOp().

◆ preCond()

virtual void AMRPoissonOp::preCond ( LevelData< FArrayBox > &  a_cor,
const LevelData< FArrayBox > &  a_residual 
)
virtual

Given the current state of the residual the correction, apply your preconditioner to a_cor.

Implements LinearOp< LevelData< FArrayBox > >.

Reimplemented in VCAMRPoissonOp2.

Referenced by ~AMRPoissonOp().

◆ applyOp()

virtual void AMRPoissonOp::applyOp ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_phi,
bool  a_homogeneous = false 
)
virtual

In the context of solving L(phi) = rhs, set a_lhs = L(a_phi). If a_homogeneous is true, evaluate the operator using homogeneous boundary conditions.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ applyOpI()

virtual void AMRPoissonOp::applyOpI ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_phi,
bool  a_homogeneous = false 
)
virtual

despite what you might think, the "I" here means "Ignore the coarse-fine boundary"

Reimplemented in VCAMRPoissonOp2.

Referenced by ~AMRPoissonOp().

◆ applyOpNoBoundary()

virtual void AMRPoissonOp::applyOpNoBoundary ( LevelData< FArrayBox > &  a_ans,
const LevelData< FArrayBox > &  a_phi 
)
virtual

Apply the differential operator without any boundary or coarse-fine boundary conditions and no finer level

Parameters
a_ansThe result of the application of the operator to a_phi.
a_phiThe data to which the operator is applied.

Implements TGAHelmOp< LevelData< FArrayBox > >.

Reimplemented in VCAMRPoissonOp2.

Referenced by ~AMRPoissonOp().

◆ create()

virtual void AMRPoissonOp::create ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs 
)
virtual

Creat data holder a_lhs that mirrors a_rhs. You do not need to copy the data of a_rhs, just make a holder the same size.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ createCoarsened()

virtual void AMRPoissonOp::createCoarsened ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs,
const int &  a_refRat 
)
virtual

◆ assign()

virtual void AMRPoissonOp::assign ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs 
)
virtual

Set a_lhs equal to a_rhs.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ assignLocal()

virtual void AMRPoissonOp::assignLocal ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs 
)
virtual

Reimplemented from LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ buildCopier()

virtual void AMRPoissonOp::buildCopier ( Copier a_copier,
const LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs 
)
virtual

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ assignCopier()

virtual void AMRPoissonOp::assignCopier ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs,
const Copier a_copier 
)
virtual

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ zeroCovered()

virtual void AMRPoissonOp::zeroCovered ( LevelData< FArrayBox > &  a_lhs,
LevelData< FArrayBox > &  a_rhs,
const Copier a_copier 
)
virtual

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ dotProduct()

virtual Real AMRPoissonOp::dotProduct ( const LevelData< FArrayBox > &  a_1,
const LevelData< FArrayBox > &  a_2 
)
virtual

Compute and return the dot product of a_1 and a_2. In most contexts, this means return the sum over all data points of a_1*a_2.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ mDotProduct()

virtual void AMRPoissonOp::mDotProduct ( const LevelData< FArrayBox > &  a_1,
const int  a_sz,
const LevelData< FArrayBox a_2[],
Real  a_mdots[] 
)
virtual

Reimplemented from LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ incr()

virtual void AMRPoissonOp::incr ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_x,
Real  a_scale 
)
virtual

Increment by scaled amount (a_lhs += a_scale*a_x).

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ axby()

virtual void AMRPoissonOp::axby ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_x,
const LevelData< FArrayBox > &  a_y,
Real  a_a,
Real  a_b 
)
virtual

Set input to a scaled sum (a_lhs = a_a*a_x + a_b*a_y).

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ scale()

virtual void AMRPoissonOp::scale ( LevelData< FArrayBox > &  a_lhs,
const Real a_scale 
)
virtual

Multiply the input by a given scale (a_lhs *= a_scale).

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ norm()

virtual Real AMRPoissonOp::norm ( const LevelData< FArrayBox > &  a_rhs,
int  a_ord 
)
virtual

Return the norm of a_rhs. a_ord == 0 max norm, a_ord == 1 sum(abs(a_rhs)), else, L(a_ord) norm.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ localMaxNorm()

virtual Real AMRPoissonOp::localMaxNorm ( const LevelData< FArrayBox > &  a_x)
virtual

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ setToZero()

virtual void AMRPoissonOp::setToZero ( LevelData< FArrayBox > &  a_lhs)
virtual

Set a_lhs to zero.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ relax()

virtual void AMRPoissonOp::relax ( LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_residual,
int  a_iterations 
)
virtual

Use your relaxation operator to remove the high frequency wave numbers from the correction so that it may be averaged to a coarser refinement. A point relaxtion scheme, for example takes the form a_correction -= lambda*(L(a_correction) - a_residual).

Implements MGLevelOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ relaxNF()

virtual void AMRPoissonOp::relaxNF ( LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > *  a_phiCoarse,
const LevelData< FArrayBox > &  a_rhs,
int  a_iterations 
)
virtual

Do CF interpolation then relax as normal.

useful for full-multigrid/FAS schemes

Reimplemented from MGLevelOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ createCoarser()

virtual void AMRPoissonOp::createCoarser ( LevelData< FArrayBox > &  a_coarse,
const LevelData< FArrayBox > &  a_fine,
bool  ghosted 
)
virtual

Create a coarsened (by two) version of the input data. This does not include averaging the data. So if a_fine is over a Box of (0, 0, 0) (63, 63, 63), a_fine should be over a Box (0, 0, 0) (31, 31, 31).

Implements MGLevelOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ restrictResidual()

virtual void AMRPoissonOp::restrictResidual ( LevelData< FArrayBox > &  a_resCoarse,
LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_rhsFine 
)
virtual

calculate restricted residual a_resCoarse[2h] = I[h->2h] (rhsFine[h] - L[h](phiFine[h])

Implements MGLevelOp< LevelData< FArrayBox > >.

Reimplemented in VCAMRPoissonOp2.

Referenced by ~AMRPoissonOp().

◆ prolongIncrement()

virtual void AMRPoissonOp::prolongIncrement ( LevelData< FArrayBox > &  a_phiThisLevel,
const LevelData< FArrayBox > &  a_correctCoarse 
)
virtual

correct the fine solution based on coarse correction a_phiThisLevel += I[2h->h](a_correctCoarse)

Implements MGLevelOp< LevelData< FArrayBox > >.

Referenced by ~AMRPoissonOp().

◆ refToCoarser()

virtual int AMRPoissonOp::refToCoarser ( )
inlinevirtual

◆ AMRResidual()

virtual void AMRPoissonOp::AMRResidual ( LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_phiCoarse,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneousPhysBC,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)
virtual

a_residual = a_rhs - L(a_phi, a_phiFine, a_phiCoarse)

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRResidualNC()

virtual void AMRPoissonOp::AMRResidualNC ( LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneousPhysBC,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)
virtual

residual assuming no more coarser AMR levels

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRResidualNF()

virtual void AMRPoissonOp::AMRResidualNF ( LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_phiCoarse,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneousPhysBC 
)
virtual

a_residual = a_rhs - L(a_phi, a_phiCoarse)

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMROperator()

virtual void AMRPoissonOp::AMROperator ( LevelData< FArrayBox > &  a_LofPhi,
const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_phiCoarse,
bool  a_homogeneousDomBC,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)
virtual

Apply the AMR operator, including coarse-fine matching

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMROperatorNC()

virtual void AMRPoissonOp::AMROperatorNC ( LevelData< FArrayBox > &  a_LofPhi,
const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
bool  a_homogeneousBC,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)
virtual

Apply the AMR operator, including coarse-fine matching assume no coarser AMR level

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMROperatorNF()

virtual void AMRPoissonOp::AMROperatorNF ( LevelData< FArrayBox > &  a_LofPhi,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_phiCoarse,
bool  a_homogeneousBC 
)
virtual

Apply the AMR operator, including coarse-fine matching. assume no finer AMR level

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRRestrict()

virtual void AMRPoissonOp::AMRRestrict ( LevelData< FArrayBox > &  a_resCoarse,
const LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection,
bool  a_skip_res = false 
)
virtual

a_resCoarse = I[h-2h]( a_residual - L(a_correction, a_coarseCorrection)) it is assumed that a_resCoarse has already been filled in with the coarse version of AMRResidualNF and that this operation is free to overwrite in the overlap regions.

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRRestrictS()

virtual void AMRPoissonOp::AMRRestrictS ( LevelData< FArrayBox > &  a_resCoarse,
const LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection,
LevelData< FArrayBox > &  a_scratch,
bool  a_skip_res = false 
)
virtual

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRProlong()

virtual void AMRPoissonOp::AMRProlong ( LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection 
)
virtual

a_correction += I[h->h](a_coarseCorrection)

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRProlongS()

virtual void AMRPoissonOp::AMRProlongS ( LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection,
LevelData< FArrayBox > &  a_temp,
const Copier a_copier 
)
virtual

optimization of AMRProlong that sends in the existing temporary and copier

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRProlongS_2()

virtual void AMRPoissonOp::AMRProlongS_2 ( LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection,
LevelData< FArrayBox > &  a_temp,
const Copier a_copier,
const Copier a_cornerCopier,
const AMRLevelOp< LevelData< FArrayBox > > *  a_crsOp 
)
virtual

optimization of AMRProlong that sends in the existing temporary and copier – higher order

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRUpdateResidual()

virtual void AMRPoissonOp::AMRUpdateResidual ( LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection 
)
virtual

a_residual = a_residual - L(a_correction, a_coarseCorrection)

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRNorm()

virtual Real AMRPoissonOp::AMRNorm ( const LevelData< FArrayBox > &  a_coarseResid,
const LevelData< FArrayBox > &  a_fineResid,
const int &  a_refRat,
const int &  a_ord 
)
virtual

compute norm over all cells on coarse not covered by finer

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ setAlphaAndBeta()

virtual void AMRPoissonOp::setAlphaAndBeta ( const Real a_alpha,
const Real a_beta 
)
virtual

For tga to reset stuff.

Implements TGAHelmOp< LevelData< FArrayBox > >.

Reimplemented in VCAMRPoissonOp2.

Referenced by refToCoarser().

◆ setBC()

virtual void AMRPoissonOp::setBC ( const BCHolder a_bc)
virtual

Change boundary conditions.

Referenced by refToCoarser().

◆ diagonalScale()

virtual void AMRPoissonOp::diagonalScale ( LevelData< FArrayBox > &  a_rhs,
bool  a_kappaWeighted 
)
inlinevirtual

For tga stuff—in this case a noop.

Reimplemented from TGAHelmOp< LevelData< FArrayBox > >.

Reimplemented in VCAMRPoissonOp2.

◆ divideByIdentityCoef()

virtual void AMRPoissonOp::divideByIdentityCoef ( LevelData< FArrayBox > &  a_rhs)
inlinevirtual

For tga stuff—in this case a noop.

Implements TGAHelmOp< LevelData< FArrayBox > >.

Reimplemented in VCAMRPoissonOp2.

◆ fillGrad()

virtual void AMRPoissonOp::fillGrad ( const LevelData< FArrayBox > &  a_phi)
inlinevirtual

These functions are part of the LevelTGA interface......

Implements LevelTGAHelmOp< LevelData< FArrayBox >, FluxBox >.

References reflux().

◆ reflux()

virtual void AMRPoissonOp::reflux ( const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
LevelData< FArrayBox > &  a_residual,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)
virtual

Reimplemented in VCAMRPoissonOp2.

Referenced by fillGrad().

◆ getFlux() [1/3]

virtual void AMRPoissonOp::getFlux ( FluxBox a_flux,
const LevelData< FArrayBox > &  a_data,
const Box a_grid,
const DataIndex a_dit,
Real  a_scale 
)
inlinevirtual

◆ write()

virtual void AMRPoissonOp::write ( const LevelData< FArrayBox > *  a,
const char *  filename 
)
virtual

Debugging aid for solvers. Print out a "T" to a file named "filename" default implementation is to print out a message saying "LinearOp::write not implemented"

Reimplemented from LinearOp< LevelData< FArrayBox > >.

Referenced by getFlux().

◆ dx()

virtual Real AMRPoissonOp::dx ( ) const
inlinevirtual

Return dx at this level of refinement

Reimplemented from LinearOp< LevelData< FArrayBox > >.

References m_dx.

◆ levelGSRB()

virtual void AMRPoissonOp::levelGSRB ( LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs 
)
protectedvirtual

Reimplemented in VCAMRPoissonOp2.

◆ levelMultiColor()

virtual void AMRPoissonOp::levelMultiColor ( LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs 
)
protectedvirtual

Reimplemented in VCAMRPoissonOp2.

◆ looseGSRB()

virtual void AMRPoissonOp::looseGSRB ( LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs 
)
protectedvirtual

Reimplemented in VCAMRPoissonOp2.

◆ overlapGSRB()

virtual void AMRPoissonOp::overlapGSRB ( LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs 
)
protectedvirtual

Reimplemented in VCAMRPoissonOp2.

◆ levelGSRBLazy()

virtual void AMRPoissonOp::levelGSRBLazy ( LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs 
)
protectedvirtual

Reimplemented in VCAMRPoissonOp2.

◆ levelJacobi()

virtual void AMRPoissonOp::levelJacobi ( LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs 
)
protectedvirtual

Reimplemented in VCAMRPoissonOp2.

◆ homogeneousCFInterp() [1/2]

virtual void AMRPoissonOp::homogeneousCFInterp ( LevelData< FArrayBox > &  a_phif)
protectedvirtual

◆ homogeneousCFInterp() [2/2]

virtual void AMRPoissonOp::homogeneousCFInterp ( LevelData< FArrayBox > &  a_phif,
const DataIndex a_datInd,
int  a_idir,
Side::LoHiSide  a_hiorlo 
)
protectedvirtual

◆ singleBoxCFInterp()

virtual void AMRPoissonOp::singleBoxCFInterp ( FArrayBox a_phi)
protectedvirtual

◆ interpOnIVSHomo()

virtual void AMRPoissonOp::interpOnIVSHomo ( LevelData< FArrayBox > &  a_phif,
const DataIndex a_datInd,
const int  a_idir,
const Side::LoHiSide  a_hiorlo,
const IntVectSet a_interpIVS 
)
protectedvirtual

◆ getFlux() [2/3]

virtual void AMRPoissonOp::getFlux ( FArrayBox a_flux,
const FArrayBox a_data,
const Box a_edgebox,
int  a_dir,
int  a_ref = 1 
) const
protectedvirtual

◆ getFlux() [3/3]

virtual void AMRPoissonOp::getFlux ( FArrayBox a_flux,
const FArrayBox a_data,
int  a_dir,
int  a_ref = 1 
) const
protectedvirtual

Member Data Documentation

◆ m_alpha

Real AMRPoissonOp::m_alpha

public constants

◆ m_beta

Real AMRPoissonOp::m_beta

◆ m_aCoef

Real AMRPoissonOp::m_aCoef

◆ m_bCoef

Real AMRPoissonOp::m_bCoef

◆ m_dxCrse

Real AMRPoissonOp::m_dxCrse

◆ m_colors

Vector<IntVect> AMRPoissonOp::m_colors

◆ s_exchangeMode

int AMRPoissonOp::s_exchangeMode
static

◆ s_relaxMode

int AMRPoissonOp::s_relaxMode
static

◆ s_maxCoarse

int AMRPoissonOp::s_maxCoarse
static

◆ s_prolongType

int AMRPoissonOp::s_prolongType
static

◆ m_dx

Real AMRPoissonOp::m_dx
protected

Referenced by dx().

◆ m_domain

ProblemDomain AMRPoissonOp::m_domain
protected

◆ m_levelOps

LevelDataOps<FArrayBox> AMRPoissonOp::m_levelOps
protected

◆ m_bc

BCHolder AMRPoissonOp::m_bc
protected

◆ m_cfregion

CFRegion AMRPoissonOp::m_cfregion
protected

◆ m_exchangeCopier

Copier AMRPoissonOp::m_exchangeCopier
protected

◆ m_interpWithCoarser

QuadCFInterp AMRPoissonOp::m_interpWithCoarser
protected

◆ m_levfluxreg

LevelFluxRegister AMRPoissonOp::m_levfluxreg
protected

◆ m_coarsenedMGrids

DisjointBoxLayout AMRPoissonOp::m_coarsenedMGrids
protected

◆ m_refToCoarser

int AMRPoissonOp::m_refToCoarser
protected

Referenced by refToCoarser().

◆ m_refToFiner

int AMRPoissonOp::m_refToFiner
protected

The documentation for this class was generated from the following file: