Chombo + EB  3.2
Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
AMRMultiGrid< T > Class Template Reference

#include <AMRMultiGrid.H>

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

Public Member Functions

 AMRMultiGrid ()
 
virtual ~AMRMultiGrid ()
 
void outputAMR (Vector< T *> &a_data, string &a_name, int a_lmax, int a_lbase)
 
virtual void define (const ProblemDomain &a_coarseDomain, AMRLevelOpFactory< T > &a_factory, LinearSolver< T > *a_bottomSolver, int a_numLevels)
 
void addInspector (RefCountedPtr< AMRMultiGridInspector< T > > &a_inspector)
 
virtual void solve (Vector< T *> &a_phi, const Vector< T *> &a_rhs, int l_max, int l_base, bool a_zeroPhi=true, bool forceHomogeneous=false)
 
virtual void solveNoInit (Vector< T *> &a_phi, const Vector< T *> &a_rhs, int l_max, int l_base, bool a_zeroPhi=true, bool forceHomogeneous=false)
 
virtual void solveNoInitResid (Vector< T *> &a_phi, Vector< T *> &a_finalResid, const Vector< T *> &a_rhs, int l_max, int l_base, bool a_zeroPhi=true, bool forceHomogeneous=false)
 use if you want final residual More...
 
void relaxOnlyHomogeneous (Vector< T *> &a_phi, const Vector< T *> &a_rhs, int l_max, int l_base)
 
virtual void AMRVCycle (Vector< T *> &a_correction, Vector< T *> &a_residual, int l, int l_max, int l_base)
 
void setMGCycle (int a_numMG)
 
virtual void init (const Vector< T *> &a_phi, const Vector< T *> &a_rhs, int l_max, int l_base)
 
void revert (const Vector< T *> &a_phi, const Vector< T *> &a_rhs, int l_max, int l_base)
 
AMRLevelOp< T > & levelOp (int level)
 
Real computeAMRResidual (Vector< T *> &a_resid, Vector< T *> &a_phi, const Vector< T *> &a_rhs, int l_max, int l_base, bool a_homogeneousBC=false, bool a_computeNorm=true)
 
Real computeAMRResidual (Vector< T *> &a_phi, const Vector< T *> &a_rhs, int l_max, int l_min)
 
void computeAMROperator (Vector< T *> &a_lph, Vector< T *> &a_phi, int l_max, int l_base, bool a_homogeneousBC=false)
 
Vector< MGLevelOp< T > *> getAllOperators ()
 
Vector< MGLevelOp< T > *> getOperatorsOp ()
 
Vector< Vector< MGLevelOp< T > *> > getOperatorsMG ()
 
Vector< AMRLevelOp< T > *> & getAMROperators ()
 
void setSolverParameters (const int &a_pre, const int &a_post, const int &a_bottom, const int &a_numMG, const int &a_iterMax, const Real &a_eps, const Real &a_hang, const Real &a_normThresh)
 
void setBottomSolver (int l_max, int l_base)
 set up bottom solver for internal MG solver More...
 
void setBottomSolverEpsCushion (Real a_bottomSolverEpsCushion)
 
void getInfo () const
 Dump out the state of the solver. AMR levels, MG levels, bottom solver, box counts, etc. More...
 

Public Attributes

Real m_eps
 
Real m_hang
 
Real m_normThresh
 
bool m_solverParamsSet
 
int m_imin
 
int m_iterMax
 
int m_iterMin
 
int m_verbosity
 
int m_exitStatus
 
int m_pre
 
int m_post
 
int m_bottom
 
int m_numMG
 
int m_maxDepth
 max no. of coarsenings – -1 (default) means coarsen as far as possible More...
 
Real m_convergenceMetric
 
Real m_initial_rnorm
 initial residual from this solve More...
 
Real m_bottomSolverEpsCushion
 

Protected Member Functions

void relax (T &phi, T &R, int depth, int nRelax=2)
 
void computeAMRResidualLevel (Vector< T *> &a_resid, Vector< T *> &a_phi, const Vector< T *> &a_rhs, int l_max, int l_base, int ilev, bool a_homogeneousBC)
 
void clear ()
 

Protected Attributes

Vector< AMRLevelOp< T > * > m_op
 
Vector< MultiGrid< T > * > m_mg
 
Vector< T * > m_correction
 
Vector< T * > m_residual
 
Vector< T * > m_resC
 
Vector< Copierm_resCopier
 
Vector< Copierm_reverseCopier
 
NoOpSolver< T > m_nosolve
 
LinearSolver< T > * m_bottomSolver
 
Vector< char > m_hasInitBeenCalled
 

Private Member Functions

 AMRMultiGrid (const AMRMultiGrid< T > &)
 
AMRMultiGridoperator= (const AMRMultiGrid< T > &)
 

Private Attributes

Vector< RefCountedPtr< AMRMultiGridInspector< T > > > m_inspectors
 

Detailed Description

template<class T>
class AMRMultiGrid< T >

Class to solve elliptic equations using the Martin and Cartwright algorithm.

Constructor & Destructor Documentation

◆ AMRMultiGrid() [1/2]

template<class T >
AMRMultiGrid< T >::AMRMultiGrid ( )

◆ ~AMRMultiGrid()

template<class T >
AMRMultiGrid< T >::~AMRMultiGrid ( )
virtual

◆ AMRMultiGrid() [2/2]

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

Member Function Documentation

◆ outputAMR()

template<class T>
void AMRMultiGrid< T >::outputAMR ( Vector< T *> &  a_data,
string &  a_name,
int  a_lmax,
int  a_lbase 
)

◆ define()

template<class T>
void AMRMultiGrid< T >::define ( const ProblemDomain a_coarseDomain,
AMRLevelOpFactory< T > &  a_factory,
LinearSolver< T > *  a_bottomSolver,
int  a_numLevels 
)
virtual

Define the solver. a_coarseDomain is the index space on the coarsest AMR level. a_factory is the operator factory through which all special information is conveyed. a_bottomSolver is the solver to be used at the termination of multigrid coarsening. It is the client's responsibility to free up the dynamically-allocated memory. a_numLevels is the number of AMR levels.

Reimplemented in AMRFASMultiGrid< T >.

Referenced by AMRFASMultiGrid< T >::define().

◆ addInspector()

template<class T>
void AMRMultiGrid< T >::addInspector ( RefCountedPtr< AMRMultiGridInspector< T > > &  a_inspector)
inline

Add an inspector to the list of inspectors maintained by this AMRMultiGrid instance. It will be given the opportunity to record intermediate data.

◆ solve()

template<class T>
void AMRMultiGrid< T >::solve ( Vector< T *> &  a_phi,
const Vector< T *> &  a_rhs,
int  l_max,
int  l_base,
bool  a_zeroPhi = true,
bool  forceHomogeneous = false 
)
virtual

Solve L(phi) = rho from l_base to l_max. To solve over all levels, l_base = 0 and l_max = max_level = numLevels-1.

Referenced by BaseLevelHeatSolver< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::solveHelm().

◆ solveNoInit()

template<class T>
void AMRMultiGrid< T >::solveNoInit ( Vector< T *> &  a_phi,
const Vector< T *> &  a_rhs,
int  l_max,
int  l_base,
bool  a_zeroPhi = true,
bool  forceHomogeneous = false 
)
virtual

same as "solve" except user has taken the reponsibility of having previously called "init" so solver can allocate temporary holders.

Reimplemented in AMRFASMultiGrid< T >.

Referenced by AMRMultiGrid< LevelData< T > >::solve().

◆ solveNoInitResid()

template<class T>
void AMRMultiGrid< T >::solveNoInitResid ( Vector< T *> &  a_phi,
Vector< T *> &  a_finalResid,
const Vector< T *> &  a_rhs,
int  l_max,
int  l_base,
bool  a_zeroPhi = true,
bool  forceHomogeneous = false 
)
virtual

use if you want final residual

set bottom solver convergence norm and solver tolerance

Referenced by AMRMultiGrid< LevelData< T > >::solveNoInit().

◆ relaxOnlyHomogeneous()

template<class T>
void AMRMultiGrid< T >::relaxOnlyHomogeneous ( Vector< T *> &  a_phi,
const Vector< T *> &  a_rhs,
int  l_max,
int  l_base 
)

◆ AMRVCycle()

template<class T>
void AMRMultiGrid< T >::AMRVCycle ( Vector< T *> &  a_correction,
Vector< T *> &  a_residual,
int  l,
int  l_max,
int  l_base 
)
virtual

◆ setMGCycle()

template<class T >
void AMRMultiGrid< T >::setMGCycle ( int  a_numMG)

◆ init()

template<class T>
void AMRMultiGrid< T >::init ( const Vector< T *> &  a_phi,
const Vector< T *> &  a_rhs,
int  l_max,
int  l_base 
)
virtual

◆ revert()

template<class T>
void AMRMultiGrid< T >::revert ( const Vector< T *> &  a_phi,
const Vector< T *> &  a_rhs,
int  l_max,
int  l_base 
)

◆ levelOp()

template<class T >
AMRLevelOp< T > & AMRMultiGrid< T >::levelOp ( int  level)

◆ computeAMRResidual() [1/2]

template<class T>
Real AMRMultiGrid< T >::computeAMRResidual ( Vector< T *> &  a_resid,
Vector< T *> &  a_phi,
const Vector< T *> &  a_rhs,
int  l_max,
int  l_base,
bool  a_homogeneousBC = false,
bool  a_computeNorm = true 
)

◆ computeAMRResidual() [2/2]

template<class T>
Real AMRMultiGrid< T >::computeAMRResidual ( Vector< T *> &  a_phi,
const Vector< T *> &  a_rhs,
int  l_max,
int  l_min 
)

just return the normed value of computeAMRResidual. used for benchmarking

◆ computeAMROperator()

template<class T>
void AMRMultiGrid< T >::computeAMROperator ( Vector< T *> &  a_lph,
Vector< T *> &  a_phi,
int  l_max,
int  l_base,
bool  a_homogeneousBC = false 
)

lph = L(phi)

◆ getAllOperators()

template<class T >
Vector< MGLevelOp< T > *> AMRMultiGrid< T >::getAllOperators ( )

◆ getOperatorsOp()

template<class T >
Vector< MGLevelOp< T > *> AMRMultiGrid< T >::getOperatorsOp ( )

◆ getOperatorsMG()

template<class T >
Vector< Vector< MGLevelOp< T > *> > AMRMultiGrid< T >::getOperatorsMG ( )

◆ getAMROperators()

template<class T>
Vector< AMRLevelOp<T> * >& AMRMultiGrid< T >::getAMROperators ( )
inline

◆ setSolverParameters()

template<class T >
void AMRMultiGrid< T >::setSolverParameters ( const int &  a_pre,
const int &  a_post,
const int &  a_bottom,
const int &  a_numMG,
const int &  a_iterMax,
const Real a_eps,
const Real a_hang,
const Real a_normThresh 
)

Set parameters of the solve. a_pre is the number of smoothings before averaging. a_post is the number of smoothings after averaging. a_bottom is the number of smoothings at the bottom level. a_numMG = 1 for vcycle, =2 for wcycle (use 1). a_itermax is the max number of v cycles. a_hang is the minimum amount of change per vcycle. a_eps is the solution tolerance. a_normThresh is how close to zero eps*resid is allowed to get.

◆ setBottomSolver()

template<class T >
void AMRMultiGrid< T >::setBottomSolver ( int  l_max,
int  l_base 
)

set up bottom solver for internal MG solver

This function is normally called by the solve(...) function. However, it must be called if solve will not be called (in particular, if only the V-cycle is being used)

Referenced by AMRFASMultiGrid< T >::solveNoInit(), and AMRMultiGrid< LevelData< T > >::solveNoInitResid().

◆ setBottomSolverEpsCushion()

template<class T >
void AMRMultiGrid< T >::setBottomSolverEpsCushion ( Real  a_bottomSolverEpsCushion)

◆ getInfo()

template<class T >
void AMRMultiGrid< T >::getInfo ( ) const

Dump out the state of the solver. AMR levels, MG levels, bottom solver, box counts, etc.

◆ relax()

template<class T>
void AMRMultiGrid< T >::relax ( T &  phi,
T &  R,
int  depth,
int  nRelax = 2 
)
protected

◆ computeAMRResidualLevel()

template<class T>
void AMRMultiGrid< T >::computeAMRResidualLevel ( Vector< T *> &  a_resid,
Vector< T *> &  a_phi,
const Vector< T *> &  a_rhs,
int  l_max,
int  l_base,
int  ilev,
bool  a_homogeneousBC 
)
protected

◆ clear()

template<class T >
void AMRMultiGrid< T >::clear ( )
protected

◆ operator=()

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

Member Data Documentation

◆ m_eps

template<class T>
Real AMRMultiGrid< T >::m_eps

◆ m_hang

template<class T>
Real AMRMultiGrid< T >::m_hang

◆ m_normThresh

template<class T>
Real AMRMultiGrid< T >::m_normThresh

◆ m_solverParamsSet

template<class T>
bool AMRMultiGrid< T >::m_solverParamsSet

◆ m_imin

template<class T>
int AMRMultiGrid< T >::m_imin

◆ m_iterMax

template<class T>
int AMRMultiGrid< T >::m_iterMax

◆ m_iterMin

template<class T>
int AMRMultiGrid< T >::m_iterMin

◆ m_verbosity

template<class T>
int AMRMultiGrid< T >::m_verbosity

◆ m_exitStatus

template<class T>
int AMRMultiGrid< T >::m_exitStatus

◆ m_pre

template<class T>
int AMRMultiGrid< T >::m_pre

◆ m_post

template<class T>
int AMRMultiGrid< T >::m_post

◆ m_bottom

template<class T>
int AMRMultiGrid< T >::m_bottom

◆ m_numMG

template<class T>
int AMRMultiGrid< T >::m_numMG

◆ m_maxDepth

template<class T>
int AMRMultiGrid< T >::m_maxDepth

max no. of coarsenings – -1 (default) means coarsen as far as possible

If using a value besides the default, need to set it before define function is called

Referenced by AMRFASMultiGrid< T >::define(), and AMRMultiGrid< LevelData< T > >::define().

◆ m_convergenceMetric

template<class T>
Real AMRMultiGrid< T >::m_convergenceMetric

◆ m_initial_rnorm

template<class T>
Real AMRMultiGrid< T >::m_initial_rnorm

◆ m_bottomSolverEpsCushion

template<class T>
Real AMRMultiGrid< T >::m_bottomSolverEpsCushion

◆ m_op

template<class T>
Vector<AMRLevelOp<T>*> AMRMultiGrid< T >::m_op
protected

◆ m_mg

template<class T>
Vector<MultiGrid<T> *> AMRMultiGrid< T >::m_mg
protected

◆ m_correction

template<class T>
Vector<T*> AMRMultiGrid< T >::m_correction
protected

◆ m_residual

template<class T>
Vector<T*> AMRMultiGrid< T >::m_residual
protected

◆ m_resC

template<class T>
Vector<T*> AMRMultiGrid< T >::m_resC
protected

◆ m_resCopier

template<class T>
Vector<Copier> AMRMultiGrid< T >::m_resCopier
protected

◆ m_reverseCopier

template<class T>
Vector<Copier> AMRMultiGrid< T >::m_reverseCopier
protected

◆ m_nosolve

template<class T>
NoOpSolver<T> AMRMultiGrid< T >::m_nosolve
protected

◆ m_bottomSolver

template<class T>
LinearSolver<T>* AMRMultiGrid< T >::m_bottomSolver
protected

◆ m_hasInitBeenCalled

template<class T>
Vector<char> AMRMultiGrid< T >::m_hasInitBeenCalled
protected

◆ m_inspectors

template<class T>
Vector<RefCountedPtr<AMRMultiGridInspector<T> > > AMRMultiGrid< T >::m_inspectors
private

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