OpenFOAM logo
Open Source CFD Toolkit

amgSymSolver.H

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2005 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software; you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by the
00013     Free Software Foundation; either version 2 of the License, or (at your
00014     option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM; if not, write to the Free Software Foundation,
00023     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
00024 
00025 Class
00026     amgSymSolver
00027 
00028 Description
00029     Agglomerated algebraic multigrid solver tuned for the FV elliptic matrices.
00030 
00031   Characteristics:
00032       - Specific for elliptic matrices (symmetric positive definite, diagonally
00033         equal or better).
00034       - Choice of smoother: Gauss-Seidel (assymetric).
00035       - Restriction operator: summation (zero-order).
00036       - Prolongation operator: injection (zero-order).
00037       - Agglomeration algorithm: single-pass pairwise agglomeration with change
00038         of direction between levels.
00039       - Coarse matrix creation: central coefficient: summation of fine grid
00040         central coefficients with the removal of intra-cluster face;
00041         off-diagonal coefficient: summation of off-diagonal faces.
00042       - Coarse matrix scaling: performed by correction scaling, using steepest
00043         descent optimisation.
00044       - Type of cycle: Saw-tooth (V-cycle with no pre-smoothing). Top-level
00045         matrices below 100 equations solved using CG.
00046 
00047 SourceFiles
00048     amgSymSolver.C
00049     amgSymSolverCalcAgglomeration.C
00050     amgSymSolverMakeCoarseMatrix.C
00051     amgSymSolverOperations.C
00052     amgSymSolverSolve.C
00053 
00054 \*---------------------------------------------------------------------------*/
00055 
00056 #ifndef amgSymSolver_H
00057 #define amgSymSolver_H
00058 
00059 #include "lduMatrix.H"
00060 #include "PtrList.H"
00061 #include "labelField.H"
00062 #include "lduAddressingStore.H"
00063 #include "primitiveFields.H"
00064 #include "lduCoupledInterfacePtrsList.H"
00065 #include "cpuTime.H"
00066 
00067 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00068 
00069 namespace Foam
00070 {
00071 
00072 /*---------------------------------------------------------------------------*\
00073                            Class amgSymSolver Declaration
00074 \*---------------------------------------------------------------------------*/
00075 
00076 class amgSymSolver
00077 :
00078     public lduMatrix::solver
00079 {
00080     // Static data members
00081 
00082         //- Max number of levels
00083         static label maxLevels_;
00084 
00085         //- Max number of cycles
00086         static label maxCycles_;
00087 
00088         //- Number of post-smoothing sweeps
00089         static label nPostSweeps_;
00090 
00091         //- Number of smoothing sweeps on finest mesh
00092         static label nBottomSweeps_;
00093 
00094 
00095     // Private data
00096 
00097         //- Final convergence tolerance
00098         scalar tolerance_;
00099 
00100         //- Convergence tolerance relative to the initial
00101         scalar relTol_;
00102 
00103         //- Min number of cells in top level
00104         label nCellsInTopLevel_;
00105 
00106 
00107         //- Restrict addressing array
00108         // Note: restrict addressing array n tells you how to grab the next
00109         // finer level and create the n-th level of the hierarchy. It has the
00110         // same size as meshLevels; the zeroth level tells you how to grab the
00111         // finest matrix and create the zeroth level of the hierarchy
00112         PtrList<labelField> restrictAddressing_;
00113 
00114         //- Hierarchy of matrix addressing
00115         PtrList<lduAddressingStore> addrLevels_;
00116 
00117         //- Hierarchy of matrix levels
00118         PtrList<lduMatrix> matrixLevels_;
00119 
00120         //- Hierarchy of coupled interfaces.
00121         //  Warning! Needs to be cleared by hand.  
00122         PtrList<lduCoupledInterfacePtrsList> interfaceLevels_;
00123 
00124         //- Hierarchy of coupled interface coefficients
00125         PtrList<FieldField<Field, scalar> > interfaceCoeffs_;
00126 
00127         //- Timing
00128         mutable cpuTime cpu_;
00129 
00130 
00131     // Private Member Functions
00132 
00133         //- Disallow default bitwise copy construct
00134         amgSymSolver(const amgSymSolver&);
00135 
00136         //- Disallow default bitwise assignment
00137         void operator=(const amgSymSolver&);
00138 
00139         //- Simplified access to matrix level
00140         const lduMatrix& matrixLevel(const label i) const;
00141 
00142         //- Simplified access to interface level
00143         const lduCoupledInterfacePtrsList& interfaceLevel
00144         (
00145             const label i
00146         ) const;
00147 
00148         //- Simplified access to interface coeffs level
00149         const FieldField<Field, scalar>& interfaceCoeffsLevel
00150         (
00151             const label i
00152         ) const;
00153 
00154         //- Calculate agglomeration and return success
00155         bool calcAgglomeration(const label fineLevelIndex);
00156 
00157         //- Assemble coarse matrix
00158         void makeCoarseMatrix(const label fineLevelIndex);
00159 
00160         //- Restrict cell field
00161         void restrictField
00162         (
00163             scalarField& result,
00164             const scalarField& f,
00165             const label fineLevelIndex
00166         ) const;
00167 
00168         //- Prolong cell field
00169         void prolongField
00170         (
00171             scalarField& result,
00172             const scalarField& f,
00173             const label coarseLevelIndex
00174         ) const;
00175 
00176         //- Return the scaling factor
00177         scalar scalingFactor
00178         (
00179             const scalarField& coarseField,
00180             const scalarField& coarseSource,
00181             const scalarField& ACf
00182         ) const;
00183 
00184 
00185 public:
00186 
00187     //- Runtime type information
00188     TypeName("AMG");
00189 
00190 
00191     // Constructors
00192 
00193         //- Construct from lduMatrix
00194         amgSymSolver
00195         (
00196             const word& fieldName,
00197             scalarField& psi,
00198             const lduMatrix& matrix,
00199             const scalarField& source,
00200             const FieldField<Field, scalar>& coupleBouCoeffs,
00201             const FieldField<Field, scalar>& coupleIntCoeffs,
00202             const lduCoupledInterfacePtrsList& interfaces,
00203             const direction cmpt,
00204             Istream& solverData
00205         );
00206 
00207 
00208     // Destructor
00209 
00210         ~amgSymSolver();
00211 
00212 
00213     // Member Functions
00214 
00215         //- Solve
00216         lduMatrix::solverPerformance solve();
00217 
00218         //- Set maximum number of level
00219         static label setMaxLevels(const label mLevels)
00220         {
00221             label oldMaxLevels = maxLevels_;
00222             maxLevels_ = mLevels;
00223 
00224             return oldMaxLevels;
00225         }
00226 
00227         //- Set maximum number of cycles
00228         static label setMaxCycles(const label mCycles)
00229         {
00230             label oldMaxCycles = maxCycles_;
00231             maxCycles_ = mCycles;
00232 
00233             return oldMaxCycles;
00234         }
00235 
00236         //- Set number of post-smoothing sweeps
00237         static label setPostSweeps(const label pSweeps)
00238         {
00239             label oldPostSweeps = nPostSweeps_;
00240             nPostSweeps_ = pSweeps;
00241 
00242             return oldPostSweeps;
00243         }
00244 
00245         //- Set number of bottom-smoothing sweeps
00246         static label setBottomSweeps(const label bSweeps)
00247         {
00248             label oldBottomSweeps = nBottomSweeps_;
00249             nBottomSweeps_ = bSweeps;
00250 
00251             return oldBottomSweeps;
00252         }
00253 };
00254 
00255 
00256 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00257 
00258 } // End namespace Foam
00259 
00260 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00261 
00262 #endif
00263 
00264 // ************************************************************************* //
For further information go to www.openfoam.org