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

#include <EBPoissonOp.H>

Inheritance diagram for EBPoissonOp:
Inheritance graph
[legend]

Public Member Functions

 ~EBPoissonOp ()
 
 EBPoissonOp ()
 
 EBPoissonOp (const EBLevelGrid &a_eblg, const EBLevelGrid &a_eblgCoarMG, const RefCountedPtr< BaseDomainBC > &a_domainBC, const RefCountedPtr< BaseEBBC > &a_ebBC, const RealVect &a_dx, const RealVect &a_origin, const bool &a_hasMGObjects, const int &a_numPreCondIters, const int &a_relaxType, const int &a_orderEB, const Real &a_alpha, const Real &a_beta, const IntVect &a_ghostCellsPhi, const IntVect &a_ghostCellsRHS)
 
virtual void residual (LevelData< EBCellFAB > &a_residual, const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs, bool a_homogeneousPhysBC=false)
 
virtual void preCond (LevelData< EBCellFAB > &a_opPhi, const LevelData< EBCellFAB > &a_phi)
 
virtual void GSColorAllRegular (LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs, const IntVect &a_color, const Real &a_weight, const bool &a_homogeneousPhysBC)
 
virtual void GSColorAllIrregular (LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs, const IntVect &a_color, const bool &a_homogeneousPhysBC, int icolor)
 
virtual void applyOp (LevelData< EBCellFAB > &a_opPhi, const LevelData< EBCellFAB > &a_phi, bool a_homogeneousPhysBC, DataIterator &a_dit, bool do_exchange=true)
 
virtual void applyOp (LevelData< EBCellFAB > &a_opPhi, const LevelData< EBCellFAB > &a_phi, bool a_homogeneousPhysBC)
 
virtual void create (LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs)
 
virtual void assign (LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs)
 
virtual Real dotProduct (const LevelData< EBCellFAB > &a_1, const LevelData< EBCellFAB > &a_2)
 
virtual void incr (LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_x, Real a_scale)
 
virtual void axby (LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_x, const LevelData< EBCellFAB > &a_y, Real a_a, Real a_b)
 
virtual void scale (LevelData< EBCellFAB > &a_lhs, const Real &a_scale)
 
virtual Real norm (const LevelData< EBCellFAB > &a_rhs, int a_ord)
 
virtual void setToZero (LevelData< EBCellFAB > &a_lhs)
 
virtual void setVal (LevelData< EBCellFAB > &a_lhs, const Real &a_value)
 
virtual void createCoarser (LevelData< EBCellFAB > &a_coarse, const LevelData< EBCellFAB > &a_fine, bool a_ghosted)
 
virtual void relax (LevelData< EBCellFAB > &a_e, const LevelData< EBCellFAB > &a_residual, int a_iterations)
 
virtual void restrictResidual (LevelData< EBCellFAB > &a_resCoarse, LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_rhsFine)
 
virtual void prolongIncrement (LevelData< EBCellFAB > &a_phiThisLevel, const LevelData< EBCellFAB > &a_correctCoarse)
 
void levelMulticolorGS (LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs)
 
void colorGS (LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs, const IntVect &color, int icolor)
 
- Public Member Functions inherited from MGLevelOp< LevelData< EBCellFAB > >
 MGLevelOp ()
 Constructor. More...
 
virtual ~MGLevelOp ()
 Destructor. More...
 
virtual void relaxNF (LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > *a_phiCoarse, const LevelData< EBCellFAB > &a_rhs, int a_iterations)
 specialized no-fine relax function, useful for full-multigrid schemes, defaults to regular relax More...
 
virtual void restrictResidual (LevelData< EBCellFAB > &a_resCoarse, LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > *a_phiCoarse, const LevelData< EBCellFAB > &a_rhsFine, bool homogeneous)
 full-multigrid version of restrictResidual, useful for FAS-type schemes. Defaults to standard restriction More...
 
virtual void restrictR (LevelData< EBCellFAB > &a_phiCoarse, const LevelData< EBCellFAB > &a_phiFine)
 simple restriction function More...
 
void addObserver (MGLevelOpObserver< LevelData< EBCellFAB > > *a_observer)
 
virtual void applyOpMg (LevelData< EBCellFAB > &a_lhs, LevelData< EBCellFAB > &a_phi, LevelData< EBCellFAB > *a_phiCoarse, bool a_homogeneous)
 Apply an operator. More...
 
virtual void residualNF (LevelData< EBCellFAB > &a_lhs, LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > *a_phiCoarse, const LevelData< EBCellFAB > &a_rhs, bool a_homogeneous=false)
 
void removeObserver (MGLevelOpObserver< LevelData< EBCellFAB > > *a_observer)
 
void addCoarserObserver (MGLevelOp< LevelData< EBCellFAB > > *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< EBCellFAB > > &a_operator, int a_coarseningFactor)
 
int numObservers () const
 Returns the number of objects observing this operator. More...
 
- Public Member Functions inherited from LinearOp< LevelData< EBCellFAB > >
virtual ~LinearOp ()
 
virtual void clear (LevelData< EBCellFAB > &a_lhs)
 
virtual void assignLocal (LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs)
 
virtual void mDotProduct (const LevelData< EBCellFAB > &a_1, const int a_sz, const LevelData< EBCellFAB > a_2[], Real a_mdots[])
 
virtual Real dx () const
 
virtual void write (const LevelData< EBCellFAB > *a, const char *filename)
 
- Public Member Functions inherited from MGLevelOpObserver< LevelData< EBCellFAB > >
 MGLevelOpObserver ()
 Base level Constructor. Called by all subclasses. More...
 
virtual ~MGLevelOpObserver ()
 Destructor. More...
 
virtual void operatorChanged (const MGLevelOp< LevelData< EBCellFAB > > &a_operator)
 
void setObservee (MGLevelOp< LevelData< EBCellFAB > > *a_observee)
 
void clearObservee ()
 

Static Public Member Functions

static bool nextColor (IntVect &color, const IntVect &limit)
 
static int getIndex (const IntVect &a_color)
 
static IntVect getColor (const int &a_icolor)
 

Protected Member Functions

void defineStencils ()
 
void getJacobiRelaxCoeff (LevelData< EBCellFAB > &a_relaxCoeff)
 

Protected Attributes

const IntVect m_ghostCellsPhi
 
const IntVect m_ghostCellsRHS
 
Vector< IntVectm_colors
 
EBLevelGrid m_eblg
 
EBLevelGrid m_eblgCoarMG
 
RefCountedPtr< BaseDomainBCm_domainBC
 
RefCountedPtr< BaseEBBCm_ebBC
 
RealVect m_dx
 
RealVect m_invDx
 
RealVect m_invDx2
 
Real m_dxScale
 
Real m_alpha
 
Real m_beta
 
Real m_time
 
RealVect m_origin
 
int m_orderEB
 
int m_numPreCondIters
 
int m_relaxType
 
LayoutData< RefCountedPtr< EBStencil > > m_opEBStencil
 
LayoutData< RefCountedPtr< EBStencil > > m_colorEBStencil [EBPO_NUMSTEN]
 
LayoutData< RefCountedPtr< EBStencil > > m_rhsColorEBStencil [EBPO_NUMSTEN]
 
LayoutData< BaseIVFAB< Real > > m_alphaDiagWeight
 
LayoutData< BaseIVFAB< Real > > m_betaDiagWeight
 
LayoutData< VoFIteratorm_vofItIrreg
 
LayoutData< VoFIteratorm_vofItIrregColor [EBPO_NUMSTEN]
 
LayoutData< VoFIteratorm_vofItIrregDomLo [SpaceDim]
 
LayoutData< VoFIteratorm_vofItIrregDomHi [SpaceDim]
 
LayoutData< VoFIteratorm_vofItIrregColorDomLo [EBPO_NUMSTEN][SpaceDim]
 
LayoutData< VoFIteratorm_vofItIrregColorDomHi [EBPO_NUMSTEN][SpaceDim]
 
bool m_hasMGObjects
 
EBMGAverage m_ebAverageMG
 
EBMGInterp m_ebInterpMG
 
DisjointBoxLayout m_dblCoarMG
 
EBISLayout m_ebislCoarMG
 
ProblemDomain m_domainCoarMG
 

Private Member Functions

void getOpVoFStencil (VoFStencil &a_stencil, const EBISBox &a_ebisbox, const VolIndex &a_vof)
 
void getOpVoFStencil (VoFStencil &a_stencil, const int &a_dir, const Vector< VolIndex > &a_allMonotoneVoFs, const EBISBox &a_ebisbox, const VolIndex &a_vof, const bool &a_lowOrder)
 
void getOpFaceStencil (VoFStencil &a_stencil, const Vector< VolIndex > &a_allMonotoneVofs, const EBISBox &a_ebisbox, const VolIndex &a_vof, int a_dir, const Side::LoHiSide &a_side, const FaceIndex &a_face, const bool &a_lowOrder)
 
void levelJacobi (LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs)
 
void applyOpRegular (int idir, Box *a_loBox, Box *a_hiBox, int *a_hasLo, int *a_hasHi, Box &a_curOpPhiBox, Box &a_curPhiBox, int a_nComps, BaseFab< Real > &a_curOpPhiFAB, const BaseFab< Real > &a_curPhiFAB, bool a_homogeneousPhysBC, const DataIndex &a_dit, const Real &a_beta)
 
void applyOpRegularAllDirs (Box *a_loBox, Box *a_hiBox, int *a_hasLo, int *a_hasHi, Box &a_curOpPhiBox, Box &a_curPhiBox, int a_nComps, BaseFab< Real > &a_curOpPhiFAB, const BaseFab< Real > &a_curPhiFAB, bool a_homogeneousPhysBC, const DataIndex &a_dit, const Real &a_beta)
 
void applyDomainFlux (Box *a_loBox, Box *a_hiBox, int *a_hasLo, int *a_hasHi, Box &a_curPhiBox, int a_nComps, BaseFab< Real > &a_phiFAB, bool a_homogeneousPhysBC, const DataIndex &a_dit, const Real &a_beta)
 
void getInvDiagRHS (LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs)
 
 EBPoissonOp (const EBPoissonOp &a_opin)
 
void operator= (const EBPoissonOp &a_opin)
 

Detailed Description

Operator to solve (alpha + beta lapl)phi = rhs. This follows the AMRLevelOp interface.

Constructor & Destructor Documentation

◆ ~EBPoissonOp()

EBPoissonOp::~EBPoissonOp ( )

◆ EBPoissonOp() [1/3]

EBPoissonOp::EBPoissonOp ( )

◆ EBPoissonOp() [2/3]

EBPoissonOp::EBPoissonOp ( const EBLevelGrid a_eblg,
const EBLevelGrid a_eblgCoarMG,
const RefCountedPtr< BaseDomainBC > &  a_domainBC,
const RefCountedPtr< BaseEBBC > &  a_ebBC,
const RealVect a_dx,
const RealVect a_origin,
const bool &  a_hasMGObjects,
const int &  a_numPreCondIters,
const int &  a_relaxType,
const int &  a_orderEB,
const Real a_alpha,
const Real a_beta,
const IntVect a_ghostCellsPhi,
const IntVect a_ghostCellsRHS 
)

If you are approaching this operator from this interface, consider backing away and using EBPoissonOpFactory to generate these objects. Really. a_eblg, : grid at this level \ a_eblgCoarMG, : grid at intermediate multigrid level \ a_domainBC, : domain boundary conditions at this level \ a_ebBC: eb boundary conditions at this level \ a_dx: grid spacing at this level \ a_origin: offset to lowest corner of the domain \ a_hasCoarserMG: true if there is a coarser MultiGrid level. \ a_preCondIters: number of iterations to do for pre-conditioning \ a_relaxType: 0 means point Jacobi, 1 is Gauss-Seidel. \ a_orderEB: 0 to not do flux interpolation at cut faces. \ a_alpha: coefficent of identity \ a_beta: coefficient of laplacian.\ a_ghostCellsPhi: Number of ghost cells in phi, correction (typically one)\ a_ghostCellsRhs: Number of ghost cells in RHS, residual, lphi (typically zero)\ Ghost cell arguments are there for caching reasons. Once you set them, an error is thrown if you send in data that does not match.

◆ EBPoissonOp() [3/3]

EBPoissonOp::EBPoissonOp ( const EBPoissonOp a_opin)
inlineprivate

References MayDay::Error().

Member Function Documentation

◆ residual()

virtual void EBPoissonOp::residual ( LevelData< EBCellFAB > &  a_lhs,
const LevelData< EBCellFAB > &  a_phi,
const LevelData< EBCellFAB > &  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< EBCellFAB > >.

◆ preCond()

virtual void EBPoissonOp::preCond ( LevelData< EBCellFAB > &  a_cor,
const LevelData< EBCellFAB > &  a_residual 
)
virtual

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

Implements LinearOp< LevelData< EBCellFAB > >.

◆ GSColorAllRegular()

virtual void EBPoissonOp::GSColorAllRegular ( LevelData< EBCellFAB > &  a_phi,
const LevelData< EBCellFAB > &  a_rhs,
const IntVect a_color,
const Real a_weight,
const bool &  a_homogeneousPhysBC 
)
virtual

◆ GSColorAllIrregular()

virtual void EBPoissonOp::GSColorAllIrregular ( LevelData< EBCellFAB > &  a_phi,
const LevelData< EBCellFAB > &  a_rhs,
const IntVect a_color,
const bool &  a_homogeneousPhysBC,
int  icolor 
)
virtual

◆ applyOp() [1/2]

virtual void EBPoissonOp::applyOp ( LevelData< EBCellFAB > &  a_opPhi,
const LevelData< EBCellFAB > &  a_phi,
bool  a_homogeneousPhysBC,
DataIterator a_dit,
bool  do_exchange = true 
)
virtual

this is the linearop function.

Referenced by applyOp().

◆ applyOp() [2/2]

virtual void EBPoissonOp::applyOp ( LevelData< EBCellFAB > &  a_lhs,
const LevelData< EBCellFAB > &  a_phi,
bool  a_homogeneous 
)
inlinevirtual

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< EBCellFAB > >.

References applyOp(), assign(), axby(), colorGS(), create(), createCoarser(), LayoutData< T >::dataIterator(), defineStencils(), dotProduct(), getColor(), getIndex(), incr(), levelMulticolorGS(), nextColor(), norm(), prolongIncrement(), relax(), restrictResidual(), scale(), setToZero(), and setVal().

◆ create()

virtual void EBPoissonOp::create ( LevelData< EBCellFAB > &  a_lhs,
const LevelData< EBCellFAB > &  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< EBCellFAB > >.

Referenced by applyOp().

◆ assign()

virtual void EBPoissonOp::assign ( LevelData< EBCellFAB > &  a_lhs,
const LevelData< EBCellFAB > &  a_rhs 
)
virtual

Set a_lhs equal to a_rhs.

Implements LinearOp< LevelData< EBCellFAB > >.

Referenced by applyOp().

◆ dotProduct()

virtual Real EBPoissonOp::dotProduct ( const LevelData< EBCellFAB > &  a_1,
const LevelData< EBCellFAB > &  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< EBCellFAB > >.

Referenced by applyOp().

◆ incr()

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

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

Implements LinearOp< LevelData< EBCellFAB > >.

Referenced by applyOp().

◆ axby()

virtual void EBPoissonOp::axby ( LevelData< EBCellFAB > &  a_lhs,
const LevelData< EBCellFAB > &  a_x,
const LevelData< EBCellFAB > &  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< EBCellFAB > >.

Referenced by applyOp().

◆ scale()

virtual void EBPoissonOp::scale ( LevelData< EBCellFAB > &  a_lhs,
const Real a_scale 
)
virtual

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

Implements LinearOp< LevelData< EBCellFAB > >.

Referenced by applyOp().

◆ norm()

virtual Real EBPoissonOp::norm ( const LevelData< EBCellFAB > &  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< EBCellFAB > >.

Referenced by applyOp().

◆ setToZero()

virtual void EBPoissonOp::setToZero ( LevelData< EBCellFAB > &  a_lhs)
virtual

Set a_lhs to zero.

Implements LinearOp< LevelData< EBCellFAB > >.

Referenced by applyOp().

◆ setVal()

virtual void EBPoissonOp::setVal ( LevelData< EBCellFAB > &  a_lhs,
const Real a_value 
)
virtual

Referenced by applyOp().

◆ createCoarser()

virtual void EBPoissonOp::createCoarser ( LevelData< EBCellFAB > &  a_coarse,
const LevelData< EBCellFAB > &  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< EBCellFAB > >.

Referenced by applyOp().

◆ relax()

virtual void EBPoissonOp::relax ( LevelData< EBCellFAB > &  a_correction,
const LevelData< EBCellFAB > &  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< EBCellFAB > >.

Referenced by applyOp().

◆ restrictResidual()

virtual void EBPoissonOp::restrictResidual ( LevelData< EBCellFAB > &  a_resCoarse,
LevelData< EBCellFAB > &  a_phiFine,
const LevelData< EBCellFAB > &  a_rhsFine 
)
virtual

Calculate restricted residual: a_resCoarse[2h] = I[h->2h] (a_rhsFine[h] - L[h](a_phiFine[h]))

Implements MGLevelOp< LevelData< EBCellFAB > >.

Referenced by applyOp().

◆ prolongIncrement()

virtual void EBPoissonOp::prolongIncrement ( LevelData< EBCellFAB > &  a_phiThisLevel,
const LevelData< EBCellFAB > &  a_correctCoarse 
)
virtual

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

Implements MGLevelOp< LevelData< EBCellFAB > >.

Referenced by applyOp().

◆ nextColor()

static bool EBPoissonOp::nextColor ( IntVect color,
const IntVect limit 
)
static

Referenced by applyOp().

◆ getIndex()

static int EBPoissonOp::getIndex ( const IntVect a_color)
static

Referenced by applyOp().

◆ getColor()

static IntVect EBPoissonOp::getColor ( const int &  a_icolor)
static

Referenced by applyOp().

◆ levelMulticolorGS()

void EBPoissonOp::levelMulticolorGS ( LevelData< EBCellFAB > &  a_phi,
const LevelData< EBCellFAB > &  a_rhs 
)

Referenced by applyOp().

◆ colorGS()

void EBPoissonOp::colorGS ( LevelData< EBCellFAB > &  a_phi,
const LevelData< EBCellFAB > &  a_rhs,
const IntVect color,
int  icolor 
)

Referenced by applyOp().

◆ defineStencils()

void EBPoissonOp::defineStencils ( )
protected

Referenced by applyOp().

◆ getJacobiRelaxCoeff()

void EBPoissonOp::getJacobiRelaxCoeff ( LevelData< EBCellFAB > &  a_relaxCoeff)
protected

Get the coefficients (for all cells) to multiply the residual by before adding it to the previous iterate for Jacobi iteration.

◆ getOpVoFStencil() [1/2]

void EBPoissonOp::getOpVoFStencil ( VoFStencil a_stencil,
const EBISBox a_ebisbox,
const VolIndex a_vof 
)
private

◆ getOpVoFStencil() [2/2]

void EBPoissonOp::getOpVoFStencil ( VoFStencil a_stencil,
const int &  a_dir,
const Vector< VolIndex > &  a_allMonotoneVoFs,
const EBISBox a_ebisbox,
const VolIndex a_vof,
const bool &  a_lowOrder 
)
private

◆ getOpFaceStencil()

void EBPoissonOp::getOpFaceStencil ( VoFStencil a_stencil,
const Vector< VolIndex > &  a_allMonotoneVofs,
const EBISBox a_ebisbox,
const VolIndex a_vof,
int  a_dir,
const Side::LoHiSide a_side,
const FaceIndex a_face,
const bool &  a_lowOrder 
)
private

◆ levelJacobi()

void EBPoissonOp::levelJacobi ( LevelData< EBCellFAB > &  a_phi,
const LevelData< EBCellFAB > &  a_rhs 
)
private

◆ applyOpRegular()

void EBPoissonOp::applyOpRegular ( int  idir,
Box a_loBox,
Box a_hiBox,
int *  a_hasLo,
int *  a_hasHi,
Box a_curOpPhiBox,
Box a_curPhiBox,
int  a_nComps,
BaseFab< Real > &  a_curOpPhiFAB,
const BaseFab< Real > &  a_curPhiFAB,
bool  a_homogeneousPhysBC,
const DataIndex a_dit,
const Real a_beta 
)
private

applyOp() on the regular cells. Only used in line relaxation stuff opPhi comes in holding alpha*phi. this adds on beta*lapl phi

◆ applyOpRegularAllDirs()

void EBPoissonOp::applyOpRegularAllDirs ( Box a_loBox,
Box a_hiBox,
int *  a_hasLo,
int *  a_hasHi,
Box a_curOpPhiBox,
Box a_curPhiBox,
int  a_nComps,
BaseFab< Real > &  a_curOpPhiFAB,
const BaseFab< Real > &  a_curPhiFAB,
bool  a_homogeneousPhysBC,
const DataIndex a_dit,
const Real a_beta 
)
private

applyOp() on the regular cells for all directions opPhi comes in holding alpha*phi. this adds on beta*lapl phi

◆ applyDomainFlux()

void EBPoissonOp::applyDomainFlux ( Box a_loBox,
Box a_hiBox,
int *  a_hasLo,
int *  a_hasHi,
Box a_curPhiBox,
int  a_nComps,
BaseFab< Real > &  a_phiFAB,
bool  a_homogeneousPhysBC,
const DataIndex a_dit,
const Real a_beta 
)
private

◆ getInvDiagRHS()

void EBPoissonOp::getInvDiagRHS ( LevelData< EBCellFAB > &  a_lhs,
const LevelData< EBCellFAB > &  a_rhs 
)
private

This function computes: a_lhs = (1/diagonal)*a_rhs Use this function to initialize the preconditioner

◆ operator=()

void EBPoissonOp::operator= ( const EBPoissonOp a_opin)
inlineprivate

References MayDay::Error().

Member Data Documentation

◆ m_ghostCellsPhi

const IntVect EBPoissonOp::m_ghostCellsPhi
protected

◆ m_ghostCellsRHS

const IntVect EBPoissonOp::m_ghostCellsRHS
protected

◆ m_colors

Vector<IntVect> EBPoissonOp::m_colors
protected

◆ m_eblg

EBLevelGrid EBPoissonOp::m_eblg
protected

◆ m_eblgCoarMG

EBLevelGrid EBPoissonOp::m_eblgCoarMG
protected

◆ m_domainBC

RefCountedPtr<BaseDomainBC> EBPoissonOp::m_domainBC
protected

◆ m_ebBC

RefCountedPtr<BaseEBBC> EBPoissonOp::m_ebBC
protected

◆ m_dx

RealVect EBPoissonOp::m_dx
protected

◆ m_invDx

RealVect EBPoissonOp::m_invDx
protected

◆ m_invDx2

RealVect EBPoissonOp::m_invDx2
protected

◆ m_dxScale

Real EBPoissonOp::m_dxScale
protected

◆ m_alpha

Real EBPoissonOp::m_alpha
protected

◆ m_beta

Real EBPoissonOp::m_beta
protected

◆ m_time

Real EBPoissonOp::m_time
protected

◆ m_origin

RealVect EBPoissonOp::m_origin
protected

◆ m_orderEB

int EBPoissonOp::m_orderEB
protected

◆ m_numPreCondIters

int EBPoissonOp::m_numPreCondIters
protected

◆ m_relaxType

int EBPoissonOp::m_relaxType
protected

◆ m_opEBStencil

LayoutData<RefCountedPtr<EBStencil> > EBPoissonOp::m_opEBStencil
protected

◆ m_colorEBStencil

LayoutData<RefCountedPtr<EBStencil> > EBPoissonOp::m_colorEBStencil[EBPO_NUMSTEN]
protected

◆ m_rhsColorEBStencil

LayoutData<RefCountedPtr<EBStencil> > EBPoissonOp::m_rhsColorEBStencil[EBPO_NUMSTEN]
protected

◆ m_alphaDiagWeight

LayoutData<BaseIVFAB<Real> > EBPoissonOp::m_alphaDiagWeight
protected

◆ m_betaDiagWeight

LayoutData<BaseIVFAB<Real> > EBPoissonOp::m_betaDiagWeight
protected

◆ m_vofItIrreg

LayoutData<VoFIterator > EBPoissonOp::m_vofItIrreg
protected

◆ m_vofItIrregColor

LayoutData<VoFIterator > EBPoissonOp::m_vofItIrregColor[EBPO_NUMSTEN]
protected

◆ m_vofItIrregDomLo

LayoutData<VoFIterator > EBPoissonOp::m_vofItIrregDomLo[SpaceDim]
protected

◆ m_vofItIrregDomHi

LayoutData<VoFIterator > EBPoissonOp::m_vofItIrregDomHi[SpaceDim]
protected

◆ m_vofItIrregColorDomLo

LayoutData<VoFIterator > EBPoissonOp::m_vofItIrregColorDomLo[EBPO_NUMSTEN][SpaceDim]
protected

◆ m_vofItIrregColorDomHi

LayoutData<VoFIterator > EBPoissonOp::m_vofItIrregColorDomHi[EBPO_NUMSTEN][SpaceDim]
protected

◆ m_hasMGObjects

bool EBPoissonOp::m_hasMGObjects
protected

◆ m_ebAverageMG

EBMGAverage EBPoissonOp::m_ebAverageMG
protected

◆ m_ebInterpMG

EBMGInterp EBPoissonOp::m_ebInterpMG
protected

◆ m_dblCoarMG

DisjointBoxLayout EBPoissonOp::m_dblCoarMG
protected

◆ m_ebislCoarMG

EBISLayout EBPoissonOp::m_ebislCoarMG
protected

◆ m_domainCoarMG

ProblemDomain EBPoissonOp::m_domainCoarMG
protected

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