FreeFOAM The Cross-Platform CFD Toolkit
GAMGSolver.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::GAMGSolver
26 
27 Description
28  Geometric agglomerated algebraic multigrid solver.
29 
30  Characteristics:
31  - Requires positive definite, diagonally dominant matrix.
32  - Agglomeration algorithm: selectable and optionally cached.
33  - Restriction operator: summation.
34  - Prolongation operator: injection.
35  - Smoother: Gauss-Seidel.
36  - Coarse matrix creation: central coefficient: summation of fine grid
37  central coefficients with the removal of intra-cluster face;
38  off-diagonal coefficient: summation of off-diagonal faces.
39  - Coarse matrix scaling: performed by correction scaling, using steepest
40  descent optimisation.
41  - Type of cycle: V-cycle with optional pre-smoothing.
42  - Coarsest-level matrix solved using ICCG or BICCG.
43 
44 SourceFiles
45  GAMGSolver.C
46  GAMGSolverCalcAgglomeration.C
47  GAMGSolverMakeCoarseMatrix.C
48  GAMGSolverOperations.C
49  GAMGSolverSolve.C
50 
51 \*---------------------------------------------------------------------------*/
52 
53 #ifndef GAMGSolver_H
54 #define GAMGSolver_H
55 
57 #include <OpenFOAM/lduMatrix.H>
58 #include <OpenFOAM/labelField.H>
61 #include <OpenFOAM/Switch.H>
62 
63 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
64 
65 namespace Foam
66 {
67 
68 /*---------------------------------------------------------------------------*\
69  Class GAMGSolver Declaration
70 \*---------------------------------------------------------------------------*/
71 
73 :
74  public lduMatrix::solver
75 {
76  // Private data
77 
78  Switch cacheAgglomeration_;
79 
80  //- Number of pre-smoothing sweeps
81  label nPreSweeps_;
82 
83  //- Number of post-smoothing sweeps
84  label nPostSweeps_;
85 
86  //- Number of smoothing sweeps on finest mesh
87  label nFinestSweeps_;
88 
89  //- Choose if the corrections should be scaled.
90  // By default corrections for symmetric matrices are scaled
91  // but not for asymmetric matrices.
92  bool scaleCorrection_;
93 
94  //- Direct or iteratively solve the coarsest level
95  bool directSolveCoarsest_;
96 
97  //- The agglomeration
98  const GAMGAgglomeration& agglomeration_;
99 
100  //- Hierarchy of matrix levels
101  PtrList<lduMatrix> matrixLevels_;
102 
103  //- Hierarchy of interfaces.
104  // Warning: Needs to be deleted explicitly.
105  PtrList<lduInterfaceFieldPtrsList> interfaceLevels_;
106 
107  //- Hierarchy of interface boundary coefficients
108  PtrList<FieldField<Field, scalar> > interfaceLevelsBouCoeffs_;
109 
110  //- Hierarchy of interface internal coefficients
111  PtrList<FieldField<Field, scalar> > interfaceLevelsIntCoeffs_;
112 
113  //- LU decompsed coarsest matrix
114  autoPtr<LUscalarMatrix> coarsestLUMatrixPtr_;
115 
116 
117  // Private Member Functions
118 
119  //- Read control parameters from the control dictionary
120  virtual void readControls();
121 
122  //- Simplified access to interface level
123  const lduInterfaceFieldPtrsList& interfaceLevel
124  (
125  const label i
126  ) const;
127 
128  //- Simplified access to matrix level
129  const lduMatrix& matrixLevel(const label i) const;
130 
131  //- Simplified access to interface boundary coeffs level
132  const FieldField<Field, scalar>& interfaceBouCoeffsLevel
133  (
134  const label i
135  ) const;
136 
137  //- Simplified access to interface internal coeffs level
138  const FieldField<Field, scalar>& interfaceIntCoeffsLevel
139  (
140  const label i
141  ) const;
142 
143  //- Agglomerate coarse matrix
144  void agglomerateMatrix(const label fineLevelIndex);
145 
146  //- Calculate and return the scaling factor from Acf, coarseSource
147  // and coarseField.
148  // At the same time do a Jacobi iteration on the coarseField using
149  // the Acf provided after the coarseField values are used for the
150  // scaling factor.
151  scalar scalingFactor
152  (
153  scalarField& field,
154  const scalarField& source,
155  const scalarField& Acf,
156  const scalarField& D
157  ) const;
158 
159  //- Calculate Acf and calculate and return the scaling factor.
160  scalar scalingFactor
161  (
162  scalarField& Acf,
163  const lduMatrix& A,
164  scalarField& field,
165  const FieldField<Field, scalar>& interfaceLevelBouCoeffs,
166  const lduInterfaceFieldPtrsList& interfaceLevel,
167  const scalarField& source,
168  const direction cmpt
169  ) const;
170 
171 
172  //- Initialise the data structures for the V-cycle
173  void initVcycle
174  (
175  PtrList<scalarField>& coarseCorrFields,
176  PtrList<scalarField>& coarseSources,
178  ) const;
179 
180 
181  //- Perform a single GAMG V-cycle with pre, post and finest smoothing.
182  void Vcycle
183  (
184  const PtrList<lduMatrix::smoother>& smoothers,
185  scalarField& psi,
186  const scalarField& source,
187  scalarField& Apsi,
188  scalarField& finestCorrection,
189  scalarField& finestResidual,
190  PtrList<scalarField>& coarseCorrFields,
191  PtrList<scalarField>& coarseSources,
192  const direction cmpt=0
193  ) const;
194 
195 
196  //- Solve the coarsest level with either an iterative or direct solver
197  void solveCoarsestLevel
198  (
199  scalarField& coarsestCorrField,
200  const scalarField& coarsestSource
201  ) const;
202 
203 
204 public:
205 
206  friend class GAMGPreconditioner;
207 
208  //- Runtime type information
209  TypeName("GAMG");
210 
211 
212  // Constructors
213 
214  //- Construct from lduMatrix and solver controls
215  GAMGSolver
216  (
217  const word& fieldName,
218  const lduMatrix& matrix,
222  const dictionary& solverControls
223  );
224 
225 
226  // Destructor
227 
228  virtual ~GAMGSolver();
229 
230 
231  // Member Functions
232 
233  //- Solve
235  (
236  scalarField& psi,
237  const scalarField& source,
238  const direction cmpt=0
239  ) const;
240 };
241 
242 
243 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
244 
245 } // End namespace Foam
246 
247 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
248 
249 #endif
250 
251 // ************************ vim: set sw=4 sts=4 et: ************************ //