AMRNodeSolver Class Reference

#include <AMRNodeSolver.H>

Collaboration diagram for AMRNodeSolver:

Collaboration graph
[legend]

List of all members.


Detailed Description

Class which manages grid hierarchy and composite elliptic solution.

AMRNodeSolver manages the AMR / multigrid solution to the elliptic equation on a set of grids on multiple levels satisfying 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.

Public Member Functions

Constructors, destructor and defines
 AMRNodeSolver ()
 AMRNodeSolver (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 NodeLevelOp *const a_opin, int a_minLength=1)
 AMRNodeSolver (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 NodeLevelOp *const a_opin, int a_minLength=1)
 ~AMRNodeSolver ()
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 NodeLevelOp *const a_opin, int a_minLength=1)
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 NodeLevelOp *const a_opin, int a_minLength=1)
Access functions
bool isDefined () const
Parameter-setting functions
void setNumSmoothUp (int a_numSmoothUp)
void setNumSmoothDown (int a_numSmoothDown)
void setTolerance (Real a_tolerance)
void setBottomSmoothing (bool a_doBottomSmooth)
void setBottomTolerance (Real a_tolerance)
void setOperatorTolerance (Real a_operatorTolerance)
void setMaxIter (int a_maxIter)
void setMinIter (int a_minIter)
void setBottomMaxIter (int a_maxIter)
void setVerbose (bool a_verbose)
Data modification functions
void solveAMR (Vector< LevelData< NodeFArrayBox > * > &a_phiLevel, const Vector< LevelData< NodeFArrayBox > * > &a_rhsLevel)
void AMRVCycleMG (Vector< LevelData< NodeFArrayBox > * > &a_phiLevel, const Vector< LevelData< NodeFArrayBox > * > &a_rhsLevel)
Real computeResidualNorm (int a_normType)
void computeAMRResidual (LevelData< NodeFArrayBox > &a_res, Vector< LevelData< NodeFArrayBox > * > &a_phiLevel, const Vector< LevelData< NodeFArrayBox > * > &a_rhsLevel, int a_ilev)
void applyAMROperator (LevelData< NodeFArrayBox > &a_lofPhi, Vector< LevelData< NodeFArrayBox > * > &a_phiLevel, int a_ilev)
void applyAMRGradient (LevelData< NodeFArrayBox > &a_gradPhi, Vector< LevelData< NodeFArrayBox > * > &a_phiLevel, int a_ilev)

Protected Member Functions

void setDefaultValues ()
void clear ()

Protected Attributes

Real m_tolerance
Real m_operatorTolerance
int m_lBase
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_numSmoothDown
Vector< AMRNodeLevelMG * > m_amrmgLevel
bool m_isDefined
LevelNodeSolver m_levelSolver
bool m_verbose

Private Member Functions

 AMRNodeSolver (const AMRNodeSolver &)
AMRNodeSolveroperator= (const AMRNodeSolver &)

Friends

class AMRNodeLevelMG


Constructor & Destructor Documentation

AMRNodeSolver::AMRNodeSolver (  ) 

Creates an AMRNodeSolver whose internal state is undefined. Need to call define() function to use any of the functionality of the class.

AMRNodeSolver::AMRNodeSolver ( 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 NodeLevelOp *const   a_opin,
int  a_minLength = 1 
)

Creates a fully-defined AMRNodeSolver. Calls define() function with identical arguments.

AMRNodeSolver::AMRNodeSolver ( 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 NodeLevelOp *const   a_opin,
int  a_minLength = 1 
)

Creates a fully-defined AMRNodeSolver. Calls define() function with identical arguments.

AMRNodeSolver::~AMRNodeSolver (  ) 

Destructor.

AMRNodeSolver::AMRNodeSolver ( const AMRNodeSolver  )  [private]


Member Function Documentation

void AMRNodeSolver::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 NodeLevelOp *const   a_opin,
int  a_minLength = 1 
)

Defines AMRNodeSolver's internal state.

Except for a_refRatio, each Vector argument contains a component for each level, and Vector index corresponds to level number.

Parameters:
a_gridsLevel  the grids (CELL-centered) at all levels; Vector length at least a_numLevels
a_domainLevel  the CELL-centered physical domains at all levels; Vector length at least a_numLevels
a_dxLevel  the grid spacing at all levels; Vector length at least a_numLevels
a_refRatio  refinement ratios between adjacent levels, starting with a_refRatio[0], the refinement ratio between levels 0 and 1; Vector length at least a_numLevels-1
a_numLevels  the number of AMR levels in the calculation
a_lBase  index of coarsest level on which solution is to be computed; must be set at time of definition in order to build the bottom LevelNodeSolver
a_opin  pointer to the NodeLevelOp to use in the solution
a_minLength  minimum length of maximally coarsened box in LevelNodeSolver, or 0 for no coarsening

void AMRNodeSolver::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 NodeLevelOp *const   a_opin,
int  a_minLength = 1 
)

Defines AMRNodeSolver's internal state.

Except for a_refRatio, each Vector argument contains a component for each level, and Vector index corresponds to level number.

Parameters:
a_gridsLevel  the grids (CELL-centered) at all levels; Vector length at least a_numLevels
a_domainLevel  the CELL-centered physical domains at all levels; Vector length at least a_numLevels
a_dxLevel  the grid spacing at all levels; Vector length at least a_numLevels
a_refRatio  refinement ratios between adjacent levels, starting with a_refRatio[0], the refinement ratio between levels 0 and 1; Vector length at least a_numLevels-1
a_numLevels  the number of AMR levels in the calculation
a_lBase  index of coarsest level on which solution is to be computed; must be set at time of definition in order to build the bottom LevelNodeSolver
a_opin  pointer to the NodeLevelOp to use in the solution
a_minLength  minimum length of maximally coarsened box in LevelNodeSolver, or 0 for no coarsening

bool AMRNodeSolver::isDefined (  )  const

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

void AMRNodeSolver::setNumSmoothUp ( int  a_numSmoothUp  ) 

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

void AMRNodeSolver::setNumSmoothDown ( int  a_numSmoothDown  ) 

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

void AMRNodeSolver::setTolerance ( Real  a_tolerance  ) 

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

void AMRNodeSolver::setBottomSmoothing ( bool  a_doBottomSmooth  ) 

Set whether bottom solver does smoothing at bottom level. Default is set in LevelNodeSolver.

void AMRNodeSolver::setBottomTolerance ( Real  a_tolerance  ) 

Set tolerance of bottom solver. Default is set in LevelNodeSolver.

void AMRNodeSolver::setOperatorTolerance ( Real  a_operatorTolerance  ) 

Set "operator tolerance" of iterative solution. Iteration will stop if (new_residual/old_residual) > 1-a_operatorTolerance (and at least m_minIter iterations have been performed). Default is 1.0e-5.

void AMRNodeSolver::setMaxIter ( int  a_maxIter  ) 

Set max number of iterations. Default is 42.

void AMRNodeSolver::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 AMRNodeSolver::setBottomMaxIter ( int  a_maxIter  ) 

Set max number of iterations for bottom solver. Default is set in LevelNodeSolver.

void AMRNodeSolver::setVerbose ( bool  a_verbose  ) 

Set whether the solver gives output while solving. Default is true.

void AMRNodeSolver::solveAMR ( Vector< LevelData< NodeFArrayBox > * > &  a_phiLevel,
const Vector< LevelData< NodeFArrayBox > * > &  a_rhsLevel 
)

Solves the elliptic equation over the hierarchy of levels m_lBase ... m_finestLevel where m_finestLevel = m_numLevels-1.

If m_lBase > 0, then the data at level m_lBase - 1 is used to interpolate boundary conditions at boundary nodes that are not adjacent to the domain boundary.

Solves to tolerance m_tolerance.

Parameters:
a_phiLevel  pointers to current guess at the solution values for levels (lMin = max(m_lBase-1,0)) ... m_finestLevel. Vector index corresponds to level number. These values are updated in place.
a_rhsLevel  pointers to right-hand side for levels lMin ... m_finestLevel. Not modified.

void AMRNodeSolver::AMRVCycleMG ( Vector< LevelData< NodeFArrayBox > * > &  a_phiLevel,
const Vector< LevelData< NodeFArrayBox > * > &  a_rhsLevel 
)

Does one relaxation V-cycle using a multigrid solver.

Parameters:
a_phiLevel  pointers to current guess at the solution values for levels (lMin = max(m_lBase-1,0)) ... m_finestLevel. Vector index corresponds to level number. These values are updated in place.
a_rhsLevel  pointers to right-hand side for levels lMin ... m_finestLevel. Not modified.

Real AMRNodeSolver::computeResidualNorm ( int  a_normType  ) 

Calculate norm of multilevel residual on levels m_lBase to m_finestLevel. Does not include data covered by finer levels (not valid data). Assumes that the residual has already been computed in m_amrmgLevel[ilev]->m_resid for all levels ilev from m_lBase to m_finestLevel.

Parameters:
a_normType  which norm to take, or 0 for max norm

void AMRNodeSolver::computeAMRResidual ( LevelData< NodeFArrayBox > &  a_res,
Vector< LevelData< NodeFArrayBox > * > &  a_phiLevel,
const Vector< LevelData< NodeFArrayBox > * > &  a_rhsLevel,
int  a_ilev 
)

Calculate residual on level a_ilev using the multilevel operator. The residual is rhs - L(phi).

As a side effect, this function modifies a_phiLevel[a_ilev] by:

  • projecting down from interior coarse nodes of a_phiLevel[a_ilev+1], if this finer level exists;
  • interpolating from a_phiLevel[a_ilev-1] at the coarse/fine boundary, if this coarser level exists;
  • applying the physical boundary conditions.

This function also stores the residual at level a_ilev in m_amrmgLevel[a_ilev]->m_resid.

Parameters:
a_res  residual, rhs - L(phi), at level a_ilev
a_phiLevel  pointers to solution, indexed by level; max Vector index must be at least a_ilev
a_rhsLevel  pointers to right-hand side, indexed by level; max Vector index must be at least a_ilev
a_ilev  level at which the residual is evaluated

void AMRNodeSolver::applyAMROperator ( LevelData< NodeFArrayBox > &  a_lofPhi,
Vector< LevelData< NodeFArrayBox > * > &  a_phiLevel,
int  a_ilev 
)

Calculate multilevel L(phi) at level a_ilev. This also depends on data from the next coarser level and the next finer level.

As a side effect, this function modifies a_phiLevel[a_ilev] by:

  • projecting down from interior coarse nodes of a_phiLevel[a_ilev+1], if this finer level exists;
  • interpolating from a_phiLevel[a_ilev-1] at the coarse/fine boundary, if this coarser level exists;
  • applying the physical boundary conditions.
Parameters:
a_lofPhi  L(phi) result of operator evaluation at level a_ilev
a_phiLevel  pointers to solution, indexed by level; max Vector index must be at least a_ilev
a_ilev  level at which operator is evaluated

void AMRNodeSolver::applyAMRGradient ( LevelData< NodeFArrayBox > &  a_gradPhi,
Vector< LevelData< NodeFArrayBox > * > &  a_phiLevel,
int  a_ilev 
)

Calculate multilevel gradient at level a_ilev. This also depends on data from the next coarser level and the next finer level.

As a side effect, this function modifies a_phiLevel[a_ilev] by:

  • projecting down from interior coarse nodes of a_phiLevel[a_ilev+1], if this finer level exists;
  • interpolating from a_phiLevel[a_ilev-1] at the coarse/fine boundary, if this coarser level exists;
  • applying the physical boundary conditions.
Parameters:
a_gradPhi  gradient of phi at level a_ilev
a_phiLevel  pointers to solution, indexed by level; max Vector index must be at least a_ilev
a_ilev  level at which gradient is evaluated

void AMRNodeSolver::setDefaultValues (  )  [protected]

Set data to default values. Not for external use.

void AMRNodeSolver::clear (  )  [protected]

Delete m_amrmgLevel Vector.

AMRNodeSolver& AMRNodeSolver::operator= ( const AMRNodeSolver  )  [private]


Friends And Related Function Documentation

friend class AMRNodeLevelMG [friend]


Member Data Documentation

solver tolerance

operator tolerance: iteration will stop if (new_residual/old_residual) > 1-m_operatorTolerance

int AMRNodeSolver::m_lBase [protected]

index of coarsest level on which solution is to be computed

int AMRNodeSolver::m_numLevels [protected]

number of levels allowed

index of finest level used, equal to m_numLevels-1

refinement ratio between levels; m_refRatio[ilev] is refinement ratio between levels ilev and ilev+1

grid spacing at each level

grids, CELL-centered (disjoint) at each level

CELL-centered physical domain at each level

int AMRNodeSolver::m_maxIter [protected]

maximum number of solver iterations

int AMRNodeSolver::m_minIter [protected]

minimum number of solver iterations before "hung convergence" criterion applies

number of smoother iterations on sweep up

number of smoother iterations on sweep down

data at each refinement level

bool AMRNodeSolver::m_isDefined [protected]

whether the full define() function has been called

multigrid solver for coarsest level

bool AMRNodeSolver::m_verbose [protected]

whether to write verbose output; default is true


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

Generated on Tue Apr 14 14:22:50 2009 for Chombo + EB by  doxygen 1.5.5