GAMGSolver.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
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 PCG or PBiCGStab.
43 
44 SourceFiles
45  GAMGSolver.C
46  GAMGSolverAgglomerateMatrix.C
47  GAMGSolverInterpolate.C
48  GAMGSolverScale.C
49  GAMGSolverSolve.C
50 
51 \*---------------------------------------------------------------------------*/
52 
53 #ifndef GAMGSolver_H
54 #define GAMGSolver_H
55 
56 #include "GAMGAgglomeration.H"
57 #include "lduMatrix.H"
58 #include "labelField.H"
59 #include "primitiveFields.H"
60 #include "LUscalarMatrix.H"
61 
62 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
63 
64 namespace Foam
65 {
66 
67 /*---------------------------------------------------------------------------*\
68  Class GAMGSolver Declaration
69 \*---------------------------------------------------------------------------*/
70 
71 class GAMGSolver
72 :
73  public lduMatrix::solver
74 {
75  // Private Data
76 
77  bool cacheAgglomeration_;
78 
79  //- Number of pre-smoothing sweeps
80  label nPreSweeps_;
81 
82  //- Lever multiplier for the number of pre-smoothing sweeps
83  label preSweepsLevelMultiplier_;
84 
85  //- Maximum number of pre-smoothing sweeps
86  label maxPreSweeps_;
87 
88  //- Number of post-smoothing sweeps
89  label nPostSweeps_;
90 
91  //- Lever multiplier for the number of post-smoothing sweeps
92  label postSweepsLevelMultiplier_;
93 
94  //- Maximum number of post-smoothing sweeps
95  label maxPostSweeps_;
96 
97  //- Number of smoothing sweeps on finest mesh
98  label nFinestSweeps_;
99 
100  //- Choose if the corrections should be interpolated after injection.
101  // By default corrections are not interpolated.
102  bool interpolateCorrection_;
103 
104  //- Choose if the corrections should be scaled.
105  // By default corrections for symmetric matrices are scaled
106  // but not for asymmetric matrices.
107  bool scaleCorrection_;
108 
109  //- Direct or iteratively solve the coarsest level
110  bool directSolveCoarsest_;
111 
112  //- The agglomeration
113  const GAMGAgglomeration& agglomeration_;
114 
115  //- Hierarchy of matrix levels
116  PtrList<lduMatrix> matrixLevels_;
117 
118  //- Hierarchy of interfaces.
119  PtrList<PtrList<lduInterfaceField>> primitiveInterfaceLevels_;
120 
121  //- Hierarchy of interfaces in lduInterfaceFieldPtrs form
122  PtrList<lduInterfaceFieldPtrsList> interfaceLevels_;
123 
124  //- Hierarchy of interface boundary coefficients
125  PtrList<FieldField<Field, scalar>> interfaceLevelsBouCoeffs_;
126 
127  //- Hierarchy of interface internal coefficients
128  PtrList<FieldField<Field, scalar>> interfaceLevelsIntCoeffs_;
129 
130  //- LU decompsed coarsest matrix
131  autoPtr<LUscalarMatrix> coarsestLUMatrixPtr_;
132 
133 
134  // Private Member Functions
135 
136  //- Read control parameters from the control dictionary
137  virtual void readControls();
138 
139  //- Simplified access to interface level
140  const lduInterfaceFieldPtrsList& interfaceLevel
141  (
142  const label i
143  ) const;
144 
145  //- Simplified access to matrix level
146  const lduMatrix& matrixLevel(const label i) const;
147 
148  //- Simplified access to interface boundary coeffs level
149  const FieldField<Field, scalar>& interfaceBouCoeffsLevel
150  (
151  const label i
152  ) const;
153 
154  //- Simplified access to interface internal coeffs level
155  const FieldField<Field, scalar>& interfaceIntCoeffsLevel
156  (
157  const label i
158  ) const;
159 
160  //- Agglomerate coarse matrix. Supply mesh to use - so we can
161  // construct temporary matrix on the fine mesh (instead of the coarse
162  // mesh)
163  void agglomerateMatrix
164  (
165  const label fineLevelIndex,
166  const lduMesh& coarseMesh,
167  const lduInterfacePtrsList& coarseMeshInterfaces
168  );
169 
170  //- Agglomerate coarse interface coefficients
171  void agglomerateInterfaceCoefficients
172  (
173  const label fineLevelIndex,
174  const lduInterfacePtrsList& coarseMeshInterfaces,
175  PtrList<lduInterfaceField>& coarsePrimInterfaces,
176  lduInterfaceFieldPtrsList& coarseInterfaces,
177  FieldField<Field, scalar>& coarseInterfaceBouCoeffs,
178  FieldField<Field, scalar>& coarseInterfaceIntCoeffs
179  ) const;
180 
181  //- Collect matrices from other processors
182  void gatherMatrices
183  (
184  const labelList& procIDs,
185  const lduMesh& dummyMesh,
186  const label meshComm,
187 
188  const lduMatrix& mat,
192 
193  PtrList<lduMatrix>& otherMats,
194  PtrList<FieldField<Field, scalar>>& otherBouCoeffs,
195  PtrList<FieldField<Field, scalar>>& otherIntCoeffs,
196  List<boolList>& otherTransforms,
197  List<List<label>>& otherRanks
198  ) const;
199 
200  //- Agglomerate processor matrices
201  void procAgglomerateMatrix
202  (
203  // Agglomeration information
204  const labelList& procAgglomMap,
205  const List<label>& agglomProcIDs,
206 
207  const label levelI,
208 
209  // Resulting matrix
210  autoPtr<lduMatrix>& allMatrixPtr,
211  FieldField<Field, scalar>& allInterfaceBouCoeffs,
212  FieldField<Field, scalar>& allInterfaceIntCoeffs,
213  PtrList<lduInterfaceField>& allPrimitiveInterfaces,
214  lduInterfaceFieldPtrsList& allInterfaces
215  ) const;
216 
217  //- Agglomerate processor matrices
218  void procAgglomerateMatrix
219  (
220  const labelList& procAgglomMap,
221  const List<label>& agglomProcIDs,
222  const label levelI
223  );
224 
225  //- Interpolate the correction after injected prolongation
226  void interpolate
227  (
228  scalarField& psi,
229  scalarField& Apsi,
230  const lduMatrix& m,
231  const FieldField<Field, scalar>& interfaceBouCoeffs,
232  const lduInterfaceFieldPtrsList& interfaces,
233  const direction cmpt
234  ) const;
235 
236  //- Interpolate the correction after injected prolongation and
237  // re-normalise
238  void interpolate
239  (
240  scalarField& psi,
241  scalarField& Apsi,
242  const lduMatrix& m,
243  const FieldField<Field, scalar>& interfaceBouCoeffs,
244  const lduInterfaceFieldPtrsList& interfaces,
245  const labelList& restrictAddressing,
246  const scalarField& psiC,
247  const direction cmpt
248  ) const;
249 
250  //- Calculate and apply the scaling factor from Acf, coarseSource
251  // and coarseField.
252  // At the same time do a Jacobi iteration on the coarseField using
253  // the Acf provided after the coarseField values are used for the
254  // scaling factor.
255  void scale
256  (
258  scalarField& Acf,
259  const lduMatrix& A,
260  const FieldField<Field, scalar>& interfaceLevelBouCoeffs,
261  const lduInterfaceFieldPtrsList& interfaceLevel,
262  const scalarField& source,
263  const direction cmpt
264  ) const;
265 
266  //- Initialise the data structures for the V-cycle
267  void initVcycle
268  (
269  PtrList<scalarField>& coarseCorrFields,
270  PtrList<scalarField>& coarseSources,
271  PtrList<lduMatrix::smoother>& smoothers,
272  scalarField& scratch1,
273  scalarField& scratch2
274  ) const;
275 
276 
277  //- Perform a single GAMG V-cycle with pre, post and finest smoothing.
278  void Vcycle
279  (
280  const PtrList<lduMatrix::smoother>& smoothers,
281  scalarField& psi,
282  const scalarField& source,
283  scalarField& Apsi,
284  scalarField& finestCorrection,
285  scalarField& finestResidual,
286 
287  scalarField& scratch1,
288  scalarField& scratch2,
289 
290  PtrList<scalarField>& coarseCorrFields,
291  PtrList<scalarField>& coarseSources,
292  const direction cmpt=0
293  ) const;
294 
295  //- Create and return the dictionary to specify the PCG solver
296  // to solve the coarsest level
297  dictionary PCGsolverDict
298  (
299  const scalar tol,
300  const scalar relTol
301  ) const;
302 
303  //- Create and return the dictionary to specify the PBiCGStab solver
304  // to solve the coarsest level
305  dictionary PBiCGStabSolverDict
306  (
307  const scalar tol,
308  const scalar relTol
309  ) const;
310 
311  //- Solve the coarsest level with either an iterative or direct solver
312  void solveCoarsestLevel
313  (
314  scalarField& coarsestCorrField,
315  const scalarField& coarsestSource
316  ) const;
317 
318 
319 public:
321  friend class GAMGPreconditioner;
322 
323  //- Runtime type information
324  TypeName("GAMG");
325 
326 
327  // Constructors
328 
329  //- Construct from lduMatrix and solver controls
330  GAMGSolver
331  (
332  const word& fieldName,
333  const lduMatrix& matrix,
334  const FieldField<Field, scalar>& interfaceBouCoeffs,
335  const FieldField<Field, scalar>& interfaceIntCoeffs,
336  const lduInterfaceFieldPtrsList& interfaces,
337  const dictionary& solverControls
338  );
339 
340 
341  //- Destructor
342  virtual ~GAMGSolver();
343 
344 
345  // Member Functions
346 
347  //- Solve
348  virtual solverPerformance solve
349  (
350  scalarField& psi,
351  const scalarField& source,
352  const direction cmpt=0
353  ) const;
354 };
355 
356 
357 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
358 
359 } // End namespace Foam
360 
361 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
362 
363 #endif
364 
365 // ************************************************************************* //
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.
Definition: GAMGSolver.C:46
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
Geometric agglomerated algebraic multigrid preconditioner.
uint8_t direction
Definition: direction.H:45
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
Definition: lduMesh.H:59
const lduMatrix & matrix() const
Definition: lduMatrix.H:226
virtual solverPerformance solve(scalarField &psi, const scalarField &source, const direction cmpt=0) const
Solve.
const FieldField< Field, scalar > & interfaceBouCoeffs() const
Definition: lduMatrix.H:231
Generic field type.
Definition: FieldField.H:51
const FieldField< Field, scalar > & interfaceIntCoeffs() const
Definition: lduMatrix.H:236
A class for handling words, derived from string.
Definition: word.H:59
SolverPerformance is the class returned by the LduMatrix solver containing performance statistics...
Abstract base-class for lduMatrix solvers.
Definition: lduMatrix.H:93
const lduInterfaceFieldPtrsList & interfaces() const
Definition: lduMatrix.H:241
const volScalarField & psi
Geometric agglomerated algebraic multigrid solver.
Definition: GAMGSolver.H:70
virtual ~GAMGSolver()
Destructor.
Definition: GAMGSolver.C:282
Specialisations of Field<T> for scalar, vector and tensor.
lduMatrix is a general matrix class in which the coefficients are stored as three arrays...
Definition: lduMatrix.H:79
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
TypeName("GAMG")
Runtime type information.
const word & fieldName() const
Definition: lduMatrix.H:221
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
rDeltaTY field()
Geometric agglomerated algebraic multigrid agglomeration class.
Namespace for OpenFOAM.