Chombo + EB + MF  3.2
Public Member Functions | Private Member Functions | Private Attributes | List of all members
TGAHelmOp< T > Class Template Referenceabstract

#include <AMRTGA.H>

Inheritance diagram for TGAHelmOp< T >:
Inheritance graph
[legend]

Public Member Functions

 TGAHelmOp ()
 
 TGAHelmOp (bool a_isTimeDependent)
 
virtual ~TGAHelmOp ()
 Destructor. More...
 
virtual void setAlphaAndBeta (const Real &a_alpha, const Real &a_beta)=0
 
virtual void diagonalScale (T &a_rhs, bool a_kappaWeighted)
 
virtual void kappaScale (T &a_rhs)
 for eb only. kappa weight the rhs but do not multiply by the identity coefficient More...
 
virtual void diagonalScale (T &a_rhs)
 
virtual void divideByIdentityCoef (T &a_rhs)=0
 
virtual void applyOpNoBoundary (T &a_ans, const T &a_phi)=0
 
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 T & identityCoef ()
 Allows access to the identity coefficient data for the operator. More...
 
- Public Member Functions inherited from AMRLevelOp< T >
virtual void dumpAMR (Vector< T *> &a_data, string name)
 
virtual void dumpLevel (T &a_data, string name)
 
 AMRLevelOp ()
 Constructor. More...
 
virtual void dumpStuff (Vector< T *> data, string filename)
 
virtual ~AMRLevelOp ()
 Destructor. More...
 
virtual Real AMRNorm (const T &a_coarResid, const T &a_fineResid, const int &a_refRat, const int &a_ord)
 
virtual void outputLevel (T &a_rhs, string &a_name)
 
virtual void outputAMR (Vector< T *> &a_rhs, string &a_name)
 
virtual int refToCoarser ()=0
 
virtual void AMRResidual (T &a_residual, const T &a_phiFine, const T &a_phi, const T &a_phiCoarse, const T &a_rhs, bool a_homogeneousDomBC, AMRLevelOp< T > *a_finerOp)=0
 
virtual void AMRResidualNF (T &a_residual, const T &a_phi, const T &a_phiCoarse, const T &a_rhs, bool a_homogeneousBC)=0
 
virtual void AMRResidualNC (T &a_residual, const T &a_phiFine, const T &a_phi, const T &a_rhs, bool a_homogeneousBC, AMRLevelOp< T > *a_finerOp)=0
 
virtual void AMROperator (T &a_LofPhi, const T &a_phiFine, const T &a_phi, const T &a_phiCoarse, bool a_homogeneousDomBC, AMRLevelOp< T > *a_finerOp)=0
 
virtual void AMROperatorNF (T &a_LofPhi, const T &a_phi, const T &a_phiCoarse, bool a_homogeneousBC)=0
 
virtual void AMROperatorNC (T &a_LofPhi, const T &a_phiFine, const T &a_phi, bool a_homogeneousBC, AMRLevelOp< T > *a_finerOp)=0
 
virtual void AMRRestrict (T &a_resCoarse, const T &a_residual, const T &a_correction, const T &a_coarseCorrection, bool a_skip_res)=0
 
virtual void AMRProlong (T &a_correction, const T &a_coarseCorrection)=0
 
virtual void AMRUpdateResidual (T &a_residual, const T &a_correction, const T &a_coarseCorrection)=0
 
virtual void createCoarsened (T &a_lhs, const T &a_rhs, const int &a_refRat)=0
 
virtual void buildCopier (Copier &a_copier, const T &a_lhs, const T &a_rhs)
 
virtual void assignCopier (T &a_lhs, const T &a_rhs, const Copier &a_copier)
 
virtual void zeroCovered (T &a_lhs, T &a_rhs, const Copier &a_copier)
 
virtual Real localMaxNorm (const T &a_phi)
 
virtual void AMRProlongS (T &a_correction, const T &a_coarseCorrection, T &a_temp, const Copier &a_copier)
 
virtual void AMRProlongS_2 (T &a_correction, const T &a_coarseCorrection, T &a_temp, const Copier &a_copier, const Copier &a_cornerCopier, const AMRLevelOp< LevelData< FArrayBox > > *a_crsOp)
 
virtual void AMRRestrictS (T &a_resCoarse, const T &a_residual, const T &a_correction, const T &a_coarseCorrection, T &scratch, bool a_skip_res=false)
 
virtual unsigned int orderOfAccuracy (void) const
 
virtual void enforceCFConsistency (T &a_coarseCorrection, const T &a_correction)
 This routine is for operators with orderOfAccuracy()>2. More...
 
- Public Member Functions inherited from MGLevelOp< T >
 MGLevelOp ()
 Constructor. More...
 
virtual ~MGLevelOp ()
 Destructor. More...
 
virtual void createCoarser (T &a_coarse, const T &a_fine, bool ghosted)=0
 
virtual void relax (T &a_correction, const T &a_residual, int a_iterations)=0
 
virtual void relaxNF (T &a_phi, const T *a_phiCoarse, const T &a_rhs, int a_iterations)
 specialized no-fine relax function, useful for full-multigrid schemes, defaults to regular relax More...
 
virtual void restrictResidual (T &a_resCoarse, T &a_phiFine, const T &a_rhsFine)=0
 
virtual void restrictResidual (T &a_resCoarse, T &a_phiFine, const T *a_phiCoarse, const T &a_rhsFine, bool homogeneous)
 full-multigrid version of restrictResidual, useful for FAS-type schemes. Defaults to standard restriction More...
 
virtual void restrictR (T &a_phiCoarse, const T &a_phiFine)
 simple restriction function More...
 
virtual void prolongIncrement (T &a_phiThisLevel, const T &a_correctCoarse)=0
 
void addObserver (MGLevelOpObserver< T > *a_observer)
 
virtual void applyOpMg (T &a_lhs, T &a_phi, T *a_phiCoarse, bool a_homogeneous)
 Apply an operator. More...
 
virtual void residualNF (T &a_lhs, T &a_phi, const T *a_phiCoarse, const T &a_rhs, bool a_homogeneous=false)
 
void removeObserver (MGLevelOpObserver< T > *a_observer)
 
void addCoarserObserver (MGLevelOp< T > *a_operator, int a_coarseningFactor)
 
void notifyObserversOfChange ()
 This should be called whenever the operator's data is updated. More...
 
virtual void finerOperatorChanged (const MGLevelOp< T > &a_operator, int a_coarseningFactor)
 
int numObservers () const
 Returns the number of objects observing this operator. More...
 
- Public Member Functions inherited from LinearOp< T >
virtual ~LinearOp ()
 
virtual void residual (T &a_lhs, const T &a_phi, const T &a_rhs, bool a_homogeneous=false)=0
 
virtual void preCond (T &a_cor, const T &a_residual)=0
 
virtual void applyOp (T &a_lhs, const T &a_phi, bool a_homogeneous=false)=0
 
virtual void create (T &a_lhs, const T &a_rhs)=0
 
virtual void clear (T &a_lhs)
 
virtual void assign (T &a_lhs, const T &a_rhs)=0
 
virtual void assignLocal (T &a_lhs, const T &a_rhs)
 
virtual Real dotProduct (const T &a_1, const T &a_2)=0
 
virtual void mDotProduct (const T &a_1, const int a_sz, const T a_2[], Real a_mdots[])
 
virtual void incr (T &a_lhs, const T &a_x, Real a_scale)=0
 
virtual void axby (T &a_lhs, const T &a_x, const T &a_y, Real a_a, Real a_b)=0
 
virtual void scale (T &a_lhs, const Real &a_scale)=0
 
virtual Real norm (const T &a_rhs, int a_ord)=0
 
virtual Real dx () const
 
virtual void setToZero (T &a_lhs)=0
 
virtual void write (const T *a, const char *filename)
 
- Public Member Functions inherited from MGLevelOpObserver< T >
 MGLevelOpObserver ()
 Base level Constructor. Called by all subclasses. More...
 
virtual ~MGLevelOpObserver ()
 Destructor. More...
 
virtual void operatorChanged (const MGLevelOp< T > &a_operator)
 
void setObservee (MGLevelOp< T > *a_observee)
 
void clearObservee ()
 

Private Member Functions

 TGAHelmOp (const TGAHelmOp &)
 
TGAHelmOpoperator= (const TGAHelmOp &)
 

Private Attributes

bool m_isTimeDependent
 
T * m_phonyIdentityCoef
 

Detailed Description

template<class T>
class TGAHelmOp< T >

This operator is meant to represent the general helmholtz operator that AMRTGA will be solving. This operator is of the form alpha A(x) phi + beta div B(x) grad phi = rho. AMRTGA needs to reset the constants alpha and beta often.

Template Parameters
TThe LevelData class that holds solution data for the operator.

Constructor & Destructor Documentation

◆ TGAHelmOp() [1/3]

template<class T>
TGAHelmOp< T >::TGAHelmOp ( )
inline

Base class default constructor. This constructs a time-independent operator.

Referenced by TGAHelmOp< LevelData< MFCellFAB > >::identityCoef().

◆ TGAHelmOp() [2/3]

template<class T>
TGAHelmOp< T >::TGAHelmOp ( bool  a_isTimeDependent)
inlineexplicit

Base class constructor which specifies explicitly whether this operator is time-dependent.

Parameters
a_isTimeDependentIf set to true, this Helmholtz operator will be treated as a time-dependent operator.

◆ ~TGAHelmOp()

template<class T>
virtual TGAHelmOp< T >::~TGAHelmOp ( )
inlinevirtual

Destructor.

◆ TGAHelmOp() [3/3]

template<class T>
TGAHelmOp< T >::TGAHelmOp ( const TGAHelmOp< T > &  )
private

Member Function Documentation

◆ setAlphaAndBeta()

template<class T>
virtual void TGAHelmOp< T >::setAlphaAndBeta ( const Real a_alpha,
const Real a_beta 
)
pure virtual

Sets the scaling constants in the Helmholtz equation.

Parameters
a_alphaThe scaling constant that multiplies the identity term.
a_betaThe scaling constant that multiplies the derivative term.

Implemented in AMRPoissonOp, EBConductivityOp, NWOEBConductivityOp, EBAMRPoissonOp, NWOEBViscousTensorOp, EBViscousTensorOp, VCAMRPoissonOp2, MFPoissonOp, slowEBCO, NWOViscousTensorOp, ViscousTensorOp, and ResistivityOp.

Referenced by AMRTGA< T >::resetAlphaAndBeta(), BaseLevelHeatSolver< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::resetSolverAlphaAndBeta(), and TGAHelmOp< LevelData< MFCellFAB > >::~TGAHelmOp().

◆ diagonalScale() [1/2]

template<class T>
virtual void TGAHelmOp< T >::diagonalScale ( T &  a_rhs,
bool  a_kappaWeighted 
)
inlinevirtual

Scales the right hand side of the Helmholtz equation by the identity term in the operator. If you are solving rho(x) dphi/dt = L(phi), this means multiply by rho (or kappa * rho in embedded boundary calculations.

Parameters
a_rhsThe right hand side of the equation to be scaled.
a_kappaWeightedIf set to true, a_rhs will be scaled by the volume fraction in addition to the identity term.

Reimplemented in AMRPoissonOp, EBConductivityOp, NWOEBViscousTensorOp, NWOEBConductivityOp, EBViscousTensorOp, EBAMRPoissonOp, VCAMRPoissonOp2, slowEBCO, ResistivityOp, NWOViscousTensorOp, and ViscousTensorOp.

Referenced by TGAHelmOp< LevelData< MFCellFAB > >::diagonalScale().

◆ kappaScale()

template<class T>
virtual void TGAHelmOp< T >::kappaScale ( T &  a_rhs)
inlinevirtual

for eb only. kappa weight the rhs but do not multiply by the identity coefficient

Reimplemented in EBAMRPoissonOp, EBConductivityOp, NWOEBViscousTensorOp, NWOEBConductivityOp, and EBViscousTensorOp.

◆ diagonalScale() [2/2]

template<class T>
virtual void TGAHelmOp< T >::diagonalScale ( T &  a_rhs)
inlinevirtual

Scales the right hand side of the Helmholtz equation by the identity term in the operator. This version assumes volume fraction weighting.

Parameters
a_rhsThe right hand side of the equation to be scaled.

Reimplemented in MFPoissonOp, NWOViscousTensorOp, ViscousTensorOp, and ResistivityOp.

◆ divideByIdentityCoef()

template<class T>
virtual void TGAHelmOp< T >::divideByIdentityCoef ( T &  a_rhs)
pure virtual

Divides the right hand side of the Helmholtz equation by the identity coefficient rho(x) in the equation rho(x) dphi/dt = L(phi).

Parameters
a_rhsThe right hand side of the equation to be scaled.

Implemented in AMRPoissonOp, EBConductivityOp, NWOEBViscousTensorOp, EBViscousTensorOp, NWOEBConductivityOp, EBAMRPoissonOp, VCAMRPoissonOp2, MFPoissonOp, NWOViscousTensorOp, ViscousTensorOp, ResistivityOp, and slowEBCO.

Referenced by TGAHelmOp< LevelData< MFCellFAB > >::diagonalScale().

◆ applyOpNoBoundary()

template<class T>
virtual void TGAHelmOp< T >::applyOpNoBoundary ( T &  a_ans,
const T &  a_phi 
)
pure 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.

Implemented in EBAMRPoissonOp, EBConductivityOp, NWOViscousTensorOp, NWOEBViscousTensorOp, ViscousTensorOp, NWOEBConductivityOp, EBViscousTensorOp, slowEBCO, AMRPoissonOp, ResistivityOp, MFPoissonOp, and VCAMRPoissonOp2.

Referenced by TGAHelmOp< LevelData< MFCellFAB > >::diagonalScale().

◆ setTime() [1/2]

template<class T>
virtual void TGAHelmOp< T >::setTime ( Real  a_time)
inlinevirtual

Sets the time-dependent state of the operator. The default implementation does nothing and is appropriate for time-independent operators.

Parameters
a_timeThe time to be used to update the time-dependent operator.

Reimplemented in EBAMRPoissonOp, VCAMRPoissonOp2, and MFPoissonOp.

Referenced by AMRTGA< T >::newOp(), AMRTGA< T >::oneStep(), TGAHelmOp< LevelData< MFCellFAB > >::setTime(), and AMRTGA< T >::setTime().

◆ setTime() [2/2]

template<class T>
virtual void TGAHelmOp< T >::setTime ( Real  a_oldTime,
Real  a_mu,
Real  a_dt 
)
inlinevirtual

Sets the time-dependent state of the operator. This version of setTime allows one to linearly interpolate coefficients across an integration step, since it accepts arguments that define where in the step it is to be updated. The default implementation calls setTime(a_oldTime + a_mu * a_dt).

Parameters
a_oldTimeThe time at the beginning of the current step.
a_muThe fraction of the current step that has elapsed.
a_dtThe size of the current step.

Reimplemented in NWOEBViscousTensorOp, EBConductivityOp, and EBViscousTensorOp.

◆ isTimeDependent()

template<class T>
bool TGAHelmOp< T >::isTimeDependent ( ) const
inline

Returns true if the operator is time-dependent, false otherwise.

◆ identityCoef()

template<class T>
virtual T& TGAHelmOp< T >::identityCoef ( )
inlinevirtual

Allows access to the identity coefficient data for the operator.

Reimplemented in VCAMRPoissonOp2.

◆ operator=()

template<class T>
TGAHelmOp& TGAHelmOp< T >::operator= ( const TGAHelmOp< T > &  )
private

Member Data Documentation

◆ m_isTimeDependent

template<class T>
bool TGAHelmOp< T >::m_isTimeDependent
private

This flag is set if the Helmholtz operator's coefficiens are time-dependent.

Referenced by TGAHelmOp< LevelData< MFCellFAB > >::isTimeDependent().

◆ m_phonyIdentityCoef

template<class T>
T* TGAHelmOp< T >::m_phonyIdentityCoef
private

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