Public Member Functions | Private Member Functions | Private Attributes | Friends
GAMGSolver Class Reference

Geometric agglomerated algebraic multigrid solver. More...

Inheritance diagram for GAMGSolver:
Inheritance graph
[legend]
Collaboration diagram for GAMGSolver:
Collaboration graph
[legend]

Public Member Functions

 TypeName ("GAMG")
 Runtime type information. More...
 
 GAMGSolver (const word &fieldName, const lduMatrix &matrix, const FieldField< Field, scalar > &interfaceBouCoeffs, const FieldField< Field, scalar > &interfaceIntCoeffs, const lduInterfaceFieldPtrsList &interfaces, const dictionary &solverControls)
 Construct from lduMatrix and solver controls. More...
 
virtual ~GAMGSolver ()
 Destructor. More...
 
virtual solverPerformance solve (scalarField &psi, const scalarField &source, const direction cmpt=0) const
 Solve. More...
 
- Public Member Functions inherited from lduMatrix::solver
virtual const wordtype () const =0
 Runtime type information. More...
 
 declareRunTimeSelectionTable (autoPtr, solver, symMatrix,(const word &fieldName, const lduMatrix &matrix, const FieldField< Field, scalar > &interfaceBouCoeffs, const FieldField< Field, scalar > &interfaceIntCoeffs, const lduInterfaceFieldPtrsList &interfaces, const dictionary &solverControls),(fieldName, matrix, interfaceBouCoeffs, interfaceIntCoeffs, interfaces, solverControls))
 
 declareRunTimeSelectionTable (autoPtr, solver, asymMatrix,(const word &fieldName, const lduMatrix &matrix, const FieldField< Field, scalar > &interfaceBouCoeffs, const FieldField< Field, scalar > &interfaceIntCoeffs, const lduInterfaceFieldPtrsList &interfaces, const dictionary &solverControls),(fieldName, matrix, interfaceBouCoeffs, interfaceIntCoeffs, interfaces, solverControls))
 
 solver (const word &fieldName, const lduMatrix &matrix, const FieldField< Field, scalar > &interfaceBouCoeffs, const FieldField< Field, scalar > &interfaceIntCoeffs, const lduInterfaceFieldPtrsList &interfaces, const dictionary &solverControls)
 
virtual ~solver ()
 Destructor. More...
 
const wordfieldName () const
 
const lduMatrixmatrix () const
 
const FieldField< Field, scalar > & interfaceBouCoeffs () const
 
const FieldField< Field, scalar > & interfaceIntCoeffs () const
 
const lduInterfaceFieldPtrsListinterfaces () const
 
virtual void read (const dictionary &)
 Read and reset the solver parameters from the given stream. More...
 
scalar normFactor (const scalarField &psi, const scalarField &source, const scalarField &Apsi, scalarField &tmpField) const
 Return the matrix norm used to normalise the residual for the. More...
 

Private Member Functions

virtual void readControls ()
 Read control parameters from the control dictionary. More...
 
const lduInterfaceFieldPtrsListinterfaceLevel (const label i) const
 Simplified access to interface level. More...
 
const lduMatrixmatrixLevel (const label i) const
 Simplified access to matrix level. More...
 
const FieldField< Field, scalar > & interfaceBouCoeffsLevel (const label i) const
 Simplified access to interface boundary coeffs level. More...
 
const FieldField< Field, scalar > & interfaceIntCoeffsLevel (const label i) const
 Simplified access to interface internal coeffs level. More...
 
void agglomerateMatrix (const label fineLevelIndex, const lduMesh &coarseMesh, const lduInterfacePtrsList &coarseMeshInterfaces)
 Agglomerate coarse matrix. Supply mesh to use - so we can. More...
 
void agglomerateInterfaceCoefficients (const label fineLevelIndex, const lduInterfacePtrsList &coarseMeshInterfaces, PtrList< lduInterfaceField > &coarsePrimInterfaces, lduInterfaceFieldPtrsList &coarseInterfaces, FieldField< Field, scalar > &coarseInterfaceBouCoeffs, FieldField< Field, scalar > &coarseInterfaceIntCoeffs) const
 Agglomerate coarse interface coefficients. More...
 
void gatherMatrices (const labelList &procIDs, const lduMesh &dummyMesh, const label meshComm, const lduMatrix &mat, const FieldField< Field, scalar > &interfaceBouCoeffs, const FieldField< Field, scalar > &interfaceIntCoeffs, const lduInterfaceFieldPtrsList &interfaces, PtrList< lduMatrix > &otherMats, PtrList< FieldField< Field, scalar > > &otherBouCoeffs, PtrList< FieldField< Field, scalar > > &otherIntCoeffs, List< boolList > &otherTransforms, List< List< label > > &otherRanks) const
 Collect matrices from other processors. More...
 
void procAgglomerateMatrix (const labelList &procAgglomMap, const List< label > &agglomProcIDs, const label levelI, autoPtr< lduMatrix > &allMatrixPtr, FieldField< Field, scalar > &allInterfaceBouCoeffs, FieldField< Field, scalar > &allInterfaceIntCoeffs, PtrList< lduInterfaceField > &allPrimitiveInterfaces, lduInterfaceFieldPtrsList &allInterfaces) const
 Agglomerate processor matrices. More...
 
void procAgglomerateMatrix (const labelList &procAgglomMap, const List< label > &agglomProcIDs, const label levelI)
 Agglomerate processor matrices. More...
 
void interpolate (scalarField &psi, scalarField &Apsi, const lduMatrix &m, const FieldField< Field, scalar > &interfaceBouCoeffs, const lduInterfaceFieldPtrsList &interfaces, const direction cmpt) const
 Interpolate the correction after injected prolongation. More...
 
void interpolate (scalarField &psi, scalarField &Apsi, const lduMatrix &m, const FieldField< Field, scalar > &interfaceBouCoeffs, const lduInterfaceFieldPtrsList &interfaces, const labelList &restrictAddressing, const scalarField &psiC, const direction cmpt) const
 Interpolate the correction after injected prolongation and. More...
 
void scale (scalarField &field, scalarField &Acf, const lduMatrix &A, const FieldField< Field, scalar > &interfaceLevelBouCoeffs, const lduInterfaceFieldPtrsList &interfaceLevel, const scalarField &source, const direction cmpt) const
 Calculate and apply the scaling factor from Acf, coarseSource. More...
 
void initVcycle (PtrList< scalarField > &coarseCorrFields, PtrList< scalarField > &coarseSources, PtrList< lduMatrix::smoother > &smoothers, scalarField &scratch1, scalarField &scratch2) const
 Initialise the data structures for the V-cycle. More...
 
void Vcycle (const PtrList< lduMatrix::smoother > &smoothers, scalarField &psi, const scalarField &source, scalarField &Apsi, scalarField &finestCorrection, scalarField &finestResidual, scalarField &scratch1, scalarField &scratch2, PtrList< scalarField > &coarseCorrFields, PtrList< scalarField > &coarseSources, const direction cmpt=0) const
 Perform a single GAMG V-cycle with pre, post and finest smoothing. More...
 
void solveCoarsestLevel (scalarField &coarsestCorrField, const scalarField &coarsestSource) const
 Solve the coarsest level with either an iterative or direct solver. More...
 

Private Attributes

bool cacheAgglomeration_
 
label nPreSweeps_
 Number of pre-smoothing sweeps. More...
 
label preSweepsLevelMultiplier_
 Lever multiplier for the number of pre-smoothing sweeps. More...
 
label maxPreSweeps_
 Maximum number of pre-smoothing sweeps. More...
 
label nPostSweeps_
 Number of post-smoothing sweeps. More...
 
label postSweepsLevelMultiplier_
 Lever multiplier for the number of post-smoothing sweeps. More...
 
label maxPostSweeps_
 Maximum number of post-smoothing sweeps. More...
 
label nFinestSweeps_
 Number of smoothing sweeps on finest mesh. More...
 
bool interpolateCorrection_
 Choose if the corrections should be interpolated after injection. More...
 
bool scaleCorrection_
 Choose if the corrections should be scaled. More...
 
bool directSolveCoarsest_
 Direct or iteratively solve the coarsest level. More...
 
const GAMGAgglomerationagglomeration_
 The agglomeration. More...
 
PtrList< lduMatrixmatrixLevels_
 Hierarchy of matrix levels. More...
 
PtrList< PtrList< lduInterfaceField > > primitiveInterfaceLevels_
 Hierarchy of interfaces. More...
 
PtrList< lduInterfaceFieldPtrsListinterfaceLevels_
 Hierarchy of interfaces in lduInterfaceFieldPtrs form. More...
 
PtrList< FieldField< Field, scalar > > interfaceLevelsBouCoeffs_
 Hierarchy of interface boundary coefficients. More...
 
PtrList< FieldField< Field, scalar > > interfaceLevelsIntCoeffs_
 Hierarchy of interface internal coefficients. More...
 
autoPtr< LUscalarMatrixcoarsestLUMatrixPtr_
 LU decompsed coarsest matrix. More...
 

Friends

class GAMGPreconditioner
 

Additional Inherited Members

- Static Public Member Functions inherited from lduMatrix::solver
static autoPtr< solverNew (const word &fieldName, const lduMatrix &matrix, const FieldField< Field, scalar > &interfaceBouCoeffs, const FieldField< Field, scalar > &interfaceIntCoeffs, const lduInterfaceFieldPtrsList &interfaces, const dictionary &solverControls)
 Return a new solver. More...
 
- Protected Attributes inherited from lduMatrix::solver
word fieldName_
 
const lduMatrixmatrix_
 
const FieldField< Field, scalar > & interfaceBouCoeffs_
 
const FieldField< Field, scalar > & interfaceIntCoeffs_
 
lduInterfaceFieldPtrsList interfaces_
 
dictionary controlDict_
 Dictionary of controls. More...
 
label maxIter_
 Maximum number of iterations in the solver. More...
 
label minIter_
 Minimum number of iterations in the solver. More...
 
scalar tolerance_
 Final convergence tolerance. More...
 
scalar relTol_
 Convergence tolerance relative to the initial. More...
 

Detailed Description

Geometric agglomerated algebraic multigrid solver.

Characteristics:

Source files

Definition at line 70 of file GAMGSolver.H.

Constructor & Destructor Documentation

◆ GAMGSolver()

GAMGSolver ( const word fieldName,
const lduMatrix matrix,
const FieldField< Field, scalar > &  interfaceBouCoeffs,
const FieldField< Field, scalar > &  interfaceIntCoeffs,
const lduInterfaceFieldPtrsList interfaces,
const dictionary solverControls 
)

◆ ~GAMGSolver()

~GAMGSolver ( )
virtual

Destructor.

Definition at line 291 of file GAMGSolver.C.

References GAMGSolver::agglomeration_, and GAMGSolver::cacheAgglomeration_.

Member Function Documentation

◆ readControls()

void readControls ( )
privatevirtual

Read control parameters from the control dictionary.

Reimplemented from lduMatrix::solver.

Reimplemented in GAMGPreconditioner.

Definition at line 302 of file GAMGSolver.C.

References Foam::endl(), Foam::Pout, and lduMatrix::solver::readControls().

Referenced by GAMGPreconditioner::readControls().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ interfaceLevel()

const Foam::lduInterfaceFieldPtrsList & interfaceLevel ( const label  i) const
private

Simplified access to interface level.

Definition at line 359 of file GAMGSolver.C.

◆ matrixLevel()

const Foam::lduMatrix & matrixLevel ( const label  i) const
private

Simplified access to matrix level.

Definition at line 345 of file GAMGSolver.C.

◆ interfaceBouCoeffsLevel()

const Foam::FieldField< Foam::Field, Foam::scalar > & interfaceBouCoeffsLevel ( const label  i) const
private

Simplified access to interface boundary coeffs level.

Definition at line 376 of file GAMGSolver.C.

◆ interfaceIntCoeffsLevel()

const Foam::FieldField< Foam::Field, Foam::scalar > & interfaceIntCoeffsLevel ( const label  i) const
private

Simplified access to interface internal coeffs level.

Definition at line 393 of file GAMGSolver.C.

◆ agglomerateMatrix()

void agglomerateMatrix ( const label  fineLevelIndex,
const lduMesh coarseMesh,
const lduInterfacePtrsList coarseMeshInterfaces 
)
private

Agglomerate coarse matrix. Supply mesh to use - so we can.

construct temporary matrix on the fine mesh (instead of the coarse mesh)

Definition at line 34 of file GAMGSolverAgglomerateMatrix.C.

References lduMesh::comm(), lduMatrix::diag(), forAll, lduMatrix::hasLower(), lduMatrix::lower(), lduMatrix::mesh(), UPtrList::set(), PtrList::set(), UPtrList::size(), and lduMatrix::upper().

Here is the call graph for this function:

◆ agglomerateInterfaceCoefficients()

void agglomerateInterfaceCoefficients ( const label  fineLevelIndex,
const lduInterfacePtrsList coarseMeshInterfaces,
PtrList< lduInterfaceField > &  coarsePrimInterfaces,
lduInterfaceFieldPtrsList coarseInterfaces,
FieldField< Field, scalar > &  coarseInterfaceBouCoeffs,
FieldField< Field, scalar > &  coarseInterfaceIntCoeffs 
) const
private

Agglomerate coarse interface coefficients.

Definition at line 198 of file GAMGSolverAgglomerateMatrix.C.

References forAll, Foam::compressible::New(), scalarField(), UPtrList::set(), and PtrList::set().

Here is the call graph for this function:

◆ gatherMatrices()

void gatherMatrices ( const labelList procIDs,
const lduMesh dummyMesh,
const label  meshComm,
const lduMatrix mat,
const FieldField< Field, scalar > &  interfaceBouCoeffs,
const FieldField< Field, scalar > &  interfaceIntCoeffs,
const lduInterfaceFieldPtrsList interfaces,
PtrList< lduMatrix > &  otherMats,
PtrList< FieldField< Field, scalar > > &  otherBouCoeffs,
PtrList< FieldField< Field, scalar > > &  otherIntCoeffs,
List< boolList > &  otherTransforms,
List< List< label > > &  otherRanks 
) const
private

Collect matrices from other processors.

Definition at line 285 of file GAMGSolverAgglomerateMatrix.C.

References Foam::endl(), forAll, interface(), Foam::Pout, Foam::refCast(), scalarField(), UPtrList::set(), PtrList::set(), List::setSize(), PtrList::setSize(), and List::size().

Here is the call graph for this function:

◆ procAgglomerateMatrix() [1/2]

void procAgglomerateMatrix ( const labelList procAgglomMap,
const List< label > &  agglomProcIDs,
const label  levelI,
autoPtr< lduMatrix > &  allMatrixPtr,
FieldField< Field, scalar > &  allInterfaceBouCoeffs,
FieldField< Field, scalar > &  allInterfaceIntCoeffs,
PtrList< lduInterfaceField > &  allPrimitiveInterfaces,
lduInterfaceFieldPtrsList allInterfaces 
) const
private

◆ procAgglomerateMatrix() [2/2]

void procAgglomerateMatrix ( const labelList procAgglomMap,
const List< label > &  agglomProcIDs,
const label  levelI 
)
private

Agglomerate processor matrices.

Definition at line 746 of file GAMGSolverAgglomerateMatrix.C.

◆ interpolate() [1/2]

void interpolate ( scalarField psi,
scalarField Apsi,
const lduMatrix m,
const FieldField< Field, scalar > &  interfaceBouCoeffs,
const lduInterfaceFieldPtrsList interfaces,
const direction  cmpt 
) const
private

Interpolate the correction after injected prolongation.

Definition at line 31 of file GAMGSolverInterpolate.C.

References UList::begin(), lduMatrix::diag(), lduMatrix::initMatrixInterfaces(), lduMatrix::lduAddr(), lduMatrix::lower(), lduAddressing::lowerAddr(), psi, lduMatrix::updateMatrixInterfaces(), lduMatrix::upper(), and lduAddressing::upperAddr().

Here is the call graph for this function:

◆ interpolate() [2/2]

void interpolate ( scalarField psi,
scalarField Apsi,
const lduMatrix m,
const FieldField< Field, scalar > &  interfaceBouCoeffs,
const lduInterfaceFieldPtrsList interfaces,
const labelList restrictAddressing,
const scalarField psiC,
const direction  cmpt 
) const
private

Interpolate the correction after injected prolongation and.

/ re-normalise

Definition at line 86 of file GAMGSolverInterpolate.C.

References lduMatrix::diag(), Foam::MULES::interpolate(), and psi.

Here is the call graph for this function:

◆ scale()

void scale ( scalarField field,
scalarField Acf,
const lduMatrix A,
const FieldField< Field, scalar > &  interfaceLevelBouCoeffs,
const lduInterfaceFieldPtrsList interfaceLevel,
const scalarField source,
const direction  cmpt 
) const
private

Calculate and apply the scaling factor from Acf, coarseSource.

and coarseField. At the same time do a Jacobi iteration on the coarseField using the Acf provided after the coarseField values are used for the scaling factor.

Definition at line 32 of file GAMGSolverScale.C.

References A(), forAll, Foam::Pout, sf(), Foam::stabilise(), Vector2D< Cmpt >::x(), and Vector2D< Cmpt >::y().

Here is the call graph for this function:

◆ initVcycle()

void initVcycle ( PtrList< scalarField > &  coarseCorrFields,
PtrList< scalarField > &  coarseSources,
PtrList< lduMatrix::smoother > &  smoothers,
scalarField scratch1,
scalarField scratch2 
) const
private

Initialise the data structures for the V-cycle.

Definition at line 449 of file GAMGSolverSolve.C.

References lduMatrix::diag(), forAll, Foam::max(), Foam::compressible::New(), scalarField(), PtrList::set(), and PtrList::setSize().

Here is the call graph for this function:

◆ Vcycle()

void Vcycle ( const PtrList< lduMatrix::smoother > &  smoothers,
scalarField psi,
const scalarField source,
scalarField Apsi,
scalarField finestCorrection,
scalarField finestResidual,
scalarField scratch1,
scalarField scratch2,
PtrList< scalarField > &  coarseCorrFields,
PtrList< scalarField > &  coarseSources,
const direction  cmpt = 0 
) const
private

Perform a single GAMG V-cycle with pre, post and finest smoothing.

Definition at line 151 of file GAMGSolverSolve.C.

References Foam::endl(), forAll, Foam::MULES::interpolate(), Foam::min(), Foam::Pout, psi, and PtrList::set().

Here is the call graph for this function:

◆ solveCoarsestLevel()

void solveCoarsestLevel ( scalarField coarsestCorrField,
const scalarField coarsestSource 
) const
private

Solve the coarsest level with either an iterative or direct solver.

Definition at line 523 of file GAMGSolverSolve.C.

References Foam::Info, messageStream::masterStream(), SolverPerformance::print(), PBiCG::solve(), and PCG::solve().

Here is the call graph for this function:

◆ TypeName()

TypeName ( "GAMG"  )

Runtime type information.

◆ solve()

Foam::solverPerformance solve ( scalarField psi,
const scalarField source,
const direction  cmpt = 0 
) const
virtual

Friends And Related Function Documentation

◆ GAMGPreconditioner

friend class GAMGPreconditioner
friend

Definition at line 305 of file GAMGSolver.H.

Field Documentation

◆ cacheAgglomeration_

bool cacheAgglomeration_
private

Definition at line 76 of file GAMGSolver.H.

Referenced by GAMGSolver::~GAMGSolver().

◆ nPreSweeps_

label nPreSweeps_
private

Number of pre-smoothing sweeps.

Definition at line 79 of file GAMGSolver.H.

◆ preSweepsLevelMultiplier_

label preSweepsLevelMultiplier_
private

Lever multiplier for the number of pre-smoothing sweeps.

Definition at line 82 of file GAMGSolver.H.

◆ maxPreSweeps_

label maxPreSweeps_
private

Maximum number of pre-smoothing sweeps.

Definition at line 85 of file GAMGSolver.H.

◆ nPostSweeps_

label nPostSweeps_
private

Number of post-smoothing sweeps.

Definition at line 88 of file GAMGSolver.H.

◆ postSweepsLevelMultiplier_

label postSweepsLevelMultiplier_
private

Lever multiplier for the number of post-smoothing sweeps.

Definition at line 91 of file GAMGSolver.H.

◆ maxPostSweeps_

label maxPostSweeps_
private

Maximum number of post-smoothing sweeps.

Definition at line 94 of file GAMGSolver.H.

◆ nFinestSweeps_

label nFinestSweeps_
private

Number of smoothing sweeps on finest mesh.

Definition at line 97 of file GAMGSolver.H.

◆ interpolateCorrection_

bool interpolateCorrection_
private

Choose if the corrections should be interpolated after injection.

By default corrections are not interpolated.

Definition at line 101 of file GAMGSolver.H.

◆ scaleCorrection_

bool scaleCorrection_
private

Choose if the corrections should be scaled.

By default corrections for symmetric matrices are scaled but not for asymmetric matrices.

Definition at line 106 of file GAMGSolver.H.

◆ directSolveCoarsest_

bool directSolveCoarsest_
private

Direct or iteratively solve the coarsest level.

Definition at line 109 of file GAMGSolver.H.

◆ agglomeration_

const GAMGAgglomeration& agglomeration_
private

The agglomeration.

Definition at line 112 of file GAMGSolver.H.

Referenced by GAMGSolver::~GAMGSolver().

◆ matrixLevels_

PtrList<lduMatrix> matrixLevels_
private

Hierarchy of matrix levels.

Definition at line 115 of file GAMGSolver.H.

◆ primitiveInterfaceLevels_

PtrList<PtrList<lduInterfaceField> > primitiveInterfaceLevels_
private

Hierarchy of interfaces.

Definition at line 118 of file GAMGSolver.H.

◆ interfaceLevels_

PtrList<lduInterfaceFieldPtrsList> interfaceLevels_
private

Hierarchy of interfaces in lduInterfaceFieldPtrs form.

Definition at line 121 of file GAMGSolver.H.

◆ interfaceLevelsBouCoeffs_

PtrList<FieldField<Field, scalar> > interfaceLevelsBouCoeffs_
private

Hierarchy of interface boundary coefficients.

Definition at line 124 of file GAMGSolver.H.

◆ interfaceLevelsIntCoeffs_

PtrList<FieldField<Field, scalar> > interfaceLevelsIntCoeffs_
private

Hierarchy of interface internal coefficients.

Definition at line 127 of file GAMGSolver.H.

◆ coarsestLUMatrixPtr_

autoPtr<LUscalarMatrix> coarsestLUMatrixPtr_
private

LU decompsed coarsest matrix.

Definition at line 130 of file GAMGSolver.H.


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