#include <AMRNodeSolver.H>
Collaboration diagram for AMRNodeSolver:
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< Real > | m_dxLevel |
Vector< DisjointBoxLayout > | m_gridsLevel |
Vector< ProblemDomain > | m_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 &) | |
AMRNodeSolver & | operator= (const AMRNodeSolver &) |
Friends | |
class | AMRNodeLevelMG |
|
Creates an AMRNodeSolver whose internal state is undefined. Need to call define() function to use any of the functionality of the class. |
|
Creates a fully-defined AMRNodeSolver. Calls define() function with identical arguments. |
|
Creates a fully-defined AMRNodeSolver. Calls define() function with identical arguments. |
|
Destructor. |
|
|
|
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.
|
|
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.
|
|
Returns |
|
Set number of multigrid smoothings on way up V-cycle; Default is 4. |
|
Set number of multigrid smoothings on way down V-cycle; Default is 4. |
|
Set tolerance of iterative solution. Default is 1.0e-10. |
|
Set whether bottom solver does smoothing at bottom level. Default is set in LevelNodeSolver. |
|
Set tolerance of bottom solver. Default is set in LevelNodeSolver. |
|
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. |
|
Set max number of iterations. Default is 42. |
|
Set min number of iterations. Only relevant when residual is not decreasing fast enough to satisfy the "operator tolerance". Default is 4. |
|
Set max number of iterations for bottom solver. Default is set in LevelNodeSolver. |
|
Set whether the solver gives output while solving. Default is |
|
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.
|
|
Does one relaxation V-cycle using a multigrid solver.
|
|
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.
|
|
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:
This function also stores the residual at level a_ilev in m_amrmgLevel[a_ilev]->m_resid. |
|
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:
|
|
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:
|
|
Set data to default values. Not for external use. |
|
Delete m_amrmgLevel Vector. |
|
|
|
|
|
solver tolerance |
|
operator tolerance: iteration will stop if (new_residual/old_residual) > 1-m_operatorTolerance |
|
index of coarsest level on which solution is to be computed |
|
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 |
|
maximum number of solver iterations |
|
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 |
|
whether the full define() function has been called |
|
multigrid solver for coarsest level |
|
whether to write verbose output; default is |