Main Page   Namespace List   Class Hierarchy   Alphabetical List   Compound List   File List   Compound Members   File Members  

AMRSolver Class Reference

Class which manages grid hierarchy and composite elliptic solution. More...

#include <AMRSolver.H>

Collaboration diagram for AMRSolver:

Collaboration graph
[legend]
List of all members.

Public Methods

 AMRSolver ()
 Creates a AMRSolver whose internal state is undefined.

 AMRSolver (const Vector< DisjointBoxLayout > &a_gridsLevel, const Vector< Box > &a_domainLevel, const Vector< Real > &a_dxLevel, const Vector< int > &a_refRatio, int a_numLevels, int a_lBase, const LevelOp *const a_opin, int a_ncomp=1, bool a_limitCoarsening=false)
 Creates a fully-defined AMRSolver.

 AMRSolver (const Vector< DisjointBoxLayout > &a_gridsLevel, const Vector< ProblemDomain > &a_domainLevel, const Vector< Real > &a_dxLevel, const Vector< int > &a_refRatio, int a_numLevels, int a_lBase, const LevelOp *const a_opin, int a_ncomp=1, bool a_limitCoarsening=false)
 Creates a fully-defined AMRSolver.

virtual void define (const Vector< DisjointBoxLayout > &a_gridsLevel, const Vector< Box > &a_domainLevel, const Vector< Real > &a_dxLevel, const Vector< int > &a_refRatio, int a_numLevels, int a_lBase, const LevelOp *const a_opin, int a_ncomp=1, bool a_limitCoarsening=false)
 Defines AMRSolver object's internal state.

virtual void define (const Vector< DisjointBoxLayout > &a_gridsLevel, const Vector< ProblemDomain > &a_domainLevel, const Vector< Real > &a_dxLevel, const Vector< int > &a_refRatio, int a_numLevels, int a_lBase, const LevelOp *const a_opin, int a_ncomp=1, bool a_limitCoarsening=false)
 Defines AMRSolver object's internal state.

virtual ~AMRSolver ()
 Destructor.

bool isDefined () const
 Returns true if full define function has been called.

void setNumSmoothUp (int a_numSmoothUp)
 Set number of multigrid smoothings on way up v-cycle. Default is 4.

void setNumSmoothDown (int a_numSmoothDown)
 Set number of multigrid smoothings on way down v-cycle; Default is 4.

void setTolerance (Real a_tolerance)
 Set tolerance of iterative solution. Default is 1.0e-10.

void setOperatorTolerance (Real a_operatorTolerance)
 Set "operator tolerance" of iterative solution.

void setMaxIter (int a_maxIter)
 Set max number of iterations. Default is 42.

void setMinIter (int a_minIter)
 Set min number of iterations.

void setNumVCyclesBottom (int a_numVCycleBottom)
 Set the number of V-cycles performed at the bottom of the AMR v-cycle.

void setConvergenceMetric (Real a_metric, int a_comp)
 set the convergence metrics (if appropriate) for this problem

void setNormType (int a_normType)
 set norm type to use for evaluating convergence

void solveAMR (Vector< LevelData< FArrayBox > * > &a_phiLevel, const Vector< LevelData< FArrayBox > * > &a_rhsLevel, bool a_initializePhiToZero=true)
 Solves on hierarchy to tolerance m_tolerance.

virtual void AMRVCycleMG (Vector< LevelData< FArrayBox > * > &a_corrLevel, const Vector< LevelData< FArrayBox > * > &a_residLevel)
 Does one relaxation V-cycle using a MG solver.

Vector< RealcomputeResidualNorm (Vector< LevelData< FArrayBox > * > &a_phiLevel, const Vector< LevelData< FArrayBox > * > &a_rhsLevel, int a_normType)
 Calculate norm of multilevel residual on levels lBase to lmax.

Vector< RealcomputeResidualNorm (Vector< LevelData< FArrayBox > * > &a_residLevel, Vector< LevelData< FArrayBox > * > &a_phiLevel, const Vector< LevelData< FArrayBox > * > &a_rhsLevel, int a_normType)
 Calculate norm of multilevel residual on levels lBase to lmax.

void computeAMRResidual (Vector< LevelData< FArrayBox > * > &a_phiLevel, const Vector< LevelData< FArrayBox > * > &a_rhsLevel, LevelData< FArrayBox > &a_res, int a_ilev)
 Calculate multilevel residual on level ilev.

void applyAMROperator (Vector< LevelData< FArrayBox > * > &a_phiLevel, LevelData< FArrayBox > &a_LofPhi, int a_ilev)
 Calculate multilevel L(phi). includes refluxing and all that.

void applyAMROperatorHphys (Vector< LevelData< FArrayBox > * > &a_phiLevel, LevelData< FArrayBox > &a_LofPhi, int a_ilev)
 Calculate multilevel L(phi) with homogeneous physical BCs.

void setVerbose (bool a_verbose)
 set whether the solver does i/o while solving. default is true


Protected Methods

virtual void setDefaultValues ()
virtual void clear ()
Vector< RealcomputeNorm (const Vector< LevelData< FArrayBox > * > &a_phiVect, int a_normType)
LevelData< FArrayBox > & getResid (int a_level)
LevelData< FArrayBox > & getCorr (int a_level)

Protected Attributes

Real m_errorTolerance
Real m_tolerance
Vector< Realm_convergenceMetrics
Real m_operatorTolerance
int m_numLevels
int m_finestLevel
Vector< int > m_refRatio
Vector< Realm_dxLevel
Vector< DisjointBoxLayoutm_gridsLevel
Vector< ProblemDomainm_domainLevel
int m_maxIter
int m_minIter
int m_numSmoothUp
int m_numVCyclesBottom
int m_numSmoothDown
int m_normType
Vector< AMRLevelMG * > m_amrmgLevel
bool m_isDefined
LevelSolver m_levelSolver
int m_lBase
int m_ncomp
bool m_verbose

Friends

class AMRLevelMG

Detailed Description

Class which manages grid hierarchy and composite elliptic solution.

AMRSolver manages the AMR/multigrid solution to the elliptic equation on a multiple level grid that satisfies certain proper nesting conditions. It can be used either as a solver, or to apply the AMR / multigrid V-cycle as a preconditioner for some other iterative method, such as a Krylov method.


Constructor & Destructor Documentation

AMRSolver::AMRSolver  
 

Creates a AMRSolver whose internal state is undefined.

Need to call define(...) function to use any of the functionality of the class

AMRSolver::AMRSolver const Vector< DisjointBoxLayout > &    a_gridsLevel,
const Vector< Box > &    a_domainLevel,
const Vector< Real > &    a_dxLevel,
const Vector< int > &    a_refRatio,
int    a_numLevels,
int    a_lBase,
const LevelOp *const    a_opin,
int    a_ncomp = 1,
bool    a_limitCoarsening = false
 

Creates a fully-defined AMRSolver.

Calls define function with identical arguments

AMRSolver::AMRSolver const Vector< DisjointBoxLayout > &    a_gridsLevel,
const Vector< ProblemDomain > &    a_domainLevel,
const Vector< Real > &    a_dxLevel,
const Vector< int > &    a_refRatio,
int    a_numLevels,
int    a_lBase,
const LevelOp *const    a_opin,
int    a_ncomp = 1,
bool    a_limitCoarsening = false
 

Creates a fully-defined AMRSolver.

Calls define function with identical arguments

virtual AMRSolver::~AMRSolver   [virtual]
 

Destructor.


Member Function Documentation

virtual void AMRSolver::AMRVCycleMG Vector< LevelData< FArrayBox > * > &    a_corrLevel,
const Vector< LevelData< FArrayBox > * > &    a_residLevel
[virtual]
 

Does one relaxation V-cycle using a MG solver.

Problem is assumed to already be in residual-correction form.

Inputs:

  • corrLevel - pointers to current guess at the correction values for levels lMin = max(lBase-1,0) ... lMax. Vector index corresponds to level number.
  • residLevel - pointers to AMR residual for levels lMin ... lMax. Vector index corresponds to level number.
Outputs:
  • corrLevel - LevelData<FArrayBox>s pointed to for levels lMin,...,lMax are updated in place.

void AMRSolver::applyAMROperator Vector< LevelData< FArrayBox > * > &    a_phiLevel,
LevelData< FArrayBox > &    a_LofPhi,
int    a_ilev
 

Calculate multilevel L(phi). includes refluxing and all that.

Includes refluxing and all that. This is the three-level operator.

void AMRSolver::applyAMROperatorHphys Vector< LevelData< FArrayBox > * > &    a_phiLevel,
LevelData< FArrayBox > &    a_LofPhi,
int    a_ilev
 

Calculate multilevel L(phi) with homogeneous physical BCs.

Calculate multilevel L(phi) with homogeneous physical BCs (but with inhomogeneous C/F BCs). includes refluxing and all that This is the three-level operator.

virtual void AMRSolver::clear   [protected, virtual]
 

Clear internal storage

void AMRSolver::computeAMRResidual Vector< LevelData< FArrayBox > * > &    a_phiLevel,
const Vector< LevelData< FArrayBox > * > &    a_rhsLevel,
LevelData< FArrayBox > &    a_res,
int    a_ilev
 

Calculate multilevel residual on level ilev.

Vector<Real> AMRSolver::computeNorm const Vector< LevelData< FArrayBox > * > &    a_phiVect,
int    a_normType
[protected]
 

compute multilevel norm on levels lBase to lMax for all components in a_phiVect. Multilevel norm does not include covered regions.

Vector<Real> AMRSolver::computeResidualNorm Vector< LevelData< FArrayBox > * > &    a_residLevel,
Vector< LevelData< FArrayBox > * > &    a_phiLevel,
const Vector< LevelData< FArrayBox > * > &    a_rhsLevel,
int    a_normType
 

Calculate norm of multilevel residual on levels lBase to lmax.

Does not include data covered by finer levels. Also returns residual in argument a_residLevel.

Vector<Real> AMRSolver::computeResidualNorm Vector< LevelData< FArrayBox > * > &    a_phiLevel,
const Vector< LevelData< FArrayBox > * > &    a_rhsLevel,
int    a_normType
 

Calculate norm of multilevel residual on levels lBase to lmax.

Does not include data covered by finer levels.

virtual void AMRSolver::define const Vector< DisjointBoxLayout > &    a_gridsLevel,
const Vector< ProblemDomain > &    a_domainLevel,
const Vector< Real > &    a_dxLevel,
const Vector< int > &    a_refRatio,
int    a_numLevels,
int    a_lBase,
const LevelOp *const    a_opin,
int    a_ncomp = 1,
bool    a_limitCoarsening = false
[virtual]
 

Defines AMRSolver object's internal state.

Inputs:

  • gridsLevel -- The grids at all levels. Each element in the vector is a level's worth of grids. gridsLevel[0] are grids for level 0 and so forth. Vector index corresponds to level number.
  • domainLevel -- The domains at all levels. Each element in the vector is a level's domain. domainLevel[0] is the domain for level 0 and so forth. Vector index corresponds to level number.
  • dxLevel -- The grid spacing at all levels. Each element in the vector is a level's $\dx$. dxLevel[0] -- is $\dx$ for level 0 and so forth. Vector index corresponds to level number.
  • refRatio -- The refinement ratio between all levels. refRatio[0] is the refinement ratio between level 0 and level 1; Vector index corresponds to level number.
  • numlevels -- The number of AMR levels in the calculation. The length of the Vector arguments has to be at least numlevels.
  • lBase -- coarsest level on which solution is to be computed. This needs to be set at the time of definition, in order to build the bottom LevelSolver.
  • a_opin -- The levelop to use in the solution.
  • a_limitCoarsening: if true, only do multigrid coarsening down to next coarser AMR level (only coarsen by a_nRefCrse). If false, coarsen as far as possible. Only relevant when lBase > 0.

virtual void AMRSolver::define const Vector< DisjointBoxLayout > &    a_gridsLevel,
const Vector< Box > &    a_domainLevel,
const Vector< Real > &    a_dxLevel,
const Vector< int > &    a_refRatio,
int    a_numLevels,
int    a_lBase,
const LevelOp *const    a_opin,
int    a_ncomp = 1,
bool    a_limitCoarsening = false
[virtual]
 

Defines AMRSolver object's internal state.

Inputs:

  • gridsLevel -- The grids at all levels. Each element in the vector is a level's worth of grids. gridsLevel[0] are grids for level 0 and so forth. Vector index corresponds to level number.
  • domainLevel -- The domains at all levels. Each element in the vector is a level's domain. domainLevel[0] is the domain for level 0 and so forth. Vector index corresponds to level number.
  • dxLevel -- The grid spacing at all levels. Each element in the vector is a level's $\dx$. dxLevel[0] -- is $\dx$ for level 0 and so forth. Vector index corresponds to level number.
  • refRatio -- The refinement ratio between all levels. refRatio[0] is the refinement ratio between level 0 and level 1; Vector index corresponds to level number.
  • numlevels -- The number of AMR levels in the calculation. The length of the Vector arguments has to be at least numlevels.
  • lBase -- coarsest level on which solution is to be computed. This needs to be set at the time of definition, in order to build the bottom LevelSolver.
  • a_opin -- The levelop to use in the solution.
  • a_limitCoarsening: if true, only do multigrid coarsening down to next coarser AMR level (only coarsen by a_nRefCrse). If false, coarsen as far as possible. Only relevant when lBase > 0.

LevelData<FArrayBox>& AMRSolver::getCorr int    a_level [protected]
 

LevelData<FArrayBox>& AMRSolver::getResid int    a_level [protected]
 

bool AMRSolver::isDefined   const
 

Returns true if full define function has been called.

void AMRSolver::setConvergenceMetric Real    a_metric,
int    a_comp
 

set the convergence metrics (if appropriate) for this problem

sets metric against which we define convergence (for example, the norm of the RHS). Also passes a_metric through to the LevelMG (and eventually to a bottom solver.)

virtual void AMRSolver::setDefaultValues   [protected, virtual]
 

Set data to default values. not for external use

void AMRSolver::setMaxIter int    a_maxIter
 

Set max number of iterations. Default is 42.

void AMRSolver::setMinIter int    a_minIter
 

Set min number of iterations.

Only relevant when residual is not decreasing fast enough to satisfy the "operator tolerance". Default is 4.

void AMRSolver::setNormType int    a_normType
 

set norm type to use for evaluating convergence

void AMRSolver::setNumSmoothDown int    a_numSmoothDown
 

Set number of multigrid smoothings on way down v-cycle; Default is 4.

void AMRSolver::setNumSmoothUp int    a_numSmoothUp
 

Set number of multigrid smoothings on way up v-cycle. Default is 4.

void AMRSolver::setNumVCyclesBottom int    a_numVCycleBottom
 

Set the number of V-cycles performed at the bottom of the AMR v-cycle.

Set the number of v-cycles performed from the base level to the maximum coarsening of the base level during one overall v-cycle on the AMR hierarchy. Default is 1.

void AMRSolver::setOperatorTolerance Real    a_operatorTolerance
 

Set "operator tolerance" of iterative solution.

Iteration will stop if (new_residual/old_residual) > 1-operatorTolerance (and at least minIter iterations have been performed). Default is 1.0e-5.

void AMRSolver::setTolerance Real    a_tolerance
 

Set tolerance of iterative solution. Default is 1.0e-10.

void AMRSolver::setVerbose bool    a_verbose
 

set whether the solver does i/o while solving. default is true

void AMRSolver::solveAMR Vector< LevelData< FArrayBox > * > &    a_phiLevel,
const Vector< LevelData< FArrayBox > * > &    a_rhsLevel,
bool    a_initializePhiToZero = true
 

Solves on hierarchy to tolerance m_tolerance.

Solves the elliptic equation over the hierarchy of levels lBase ... lMax where $lMax = numlevels-1$. If lBase > 0, then the data at level lBase - 1 is used to interpolate boundary conditions at boundary cells that are not adjacent to the domain boundary.

Inputs:

  • phiLevel - pointers to current guess at the solution values for levels (lMin = max(lBase-1,0)) ... lMax. Vector index corresponds to level number.
  • rhsLevel - pointers to right-hand side for levels lmin ... lMax. Vector index corresponds to level number.
Outputs:
  • phiLevel - LevelData<FArrayBox>s pointed to for levels lMin,...,lMax are updated in place.


Friends And Related Function Documentation

friend class AMRLevelMG [friend]
 


Member Data Documentation

Vector<AMRLevelMG*> AMRSolver::m_amrmgLevel [protected]
 

Array of AMRLevelMGs containing data at each refinement level

Vector<Real> AMRSolver::m_convergenceMetrics [protected]
 

metrics used to judge convergence

Vector<ProblemDomain> AMRSolver::m_domainLevel [protected]
 

Domains.

Vector<Real> AMRSolver::m_dxLevel [protected]
 

Grid spacing.

Real AMRSolver::m_errorTolerance [protected]
 

Error tolerance for error estimation

int AMRSolver::m_finestLevel [protected]
 

Maximum number of working level is m_numLevels - 1

Vector<DisjointBoxLayout> AMRSolver::m_gridsLevel [protected]
 

Grids.

bool AMRSolver::m_isDefined [protected]
 

Has this solver been defined?

int AMRSolver::m_lBase [protected]
 

LevelSolver AMRSolver::m_levelSolver [protected]
 

int AMRSolver::m_maxIter [protected]
 

Maximum number of solver iterations

int AMRSolver::m_minIter [protected]
 

Minimum number of solver iterations before "hung converence" criteria is valid

int AMRSolver::m_ncomp [protected]
 

int AMRSolver::m_normType [protected]
 

Norm type to use for evaluatiing convergence

int AMRSolver::m_numLevels [protected]
 

Number of levels allowed

int AMRSolver::m_numSmoothDown [protected]
 

Number of smoother iterations on down cycle

int AMRSolver::m_numSmoothUp [protected]
 

Number of smoother iterations on up cycle

int AMRSolver::m_numVCyclesBottom [protected]
 

Number of VCycles for max iter of levelsolver

Real AMRSolver::m_operatorTolerance [protected]
 

Tolerance on maxOperator

Vector<int> AMRSolver::m_refRatio [protected]
 

Refinement ratio.

Real AMRSolver::m_tolerance [protected]
 

Solver tolerance

bool AMRSolver::m_verbose [protected]
 

default is true


The documentation for this class was generated from the following file:
Generated on Wed Jan 19 17:55:27 2005 for Chombo&INSwithParticles by doxygen1.2.16