SmoothSolver.C
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-2018 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 \*---------------------------------------------------------------------------*/
25 
26 #include "SmoothSolver.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class Type, class DType, class LUType>
32 (
33  const word& fieldName,
34  const LduMatrix<Type, DType, LUType>& matrix,
35  const dictionary& solverDict
36 )
37 :
38  LduMatrix<Type, DType, LUType>::solver
39  (
40  fieldName,
41  matrix,
42  solverDict
43  ),
44  nSweeps_(1)
45 {
46  readControls();
47 }
48 
49 
50 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
51 
52 template<class Type, class DType, class LUType>
54 {
56  this->readControl(this->controlDict_, nSweeps_, "nSweeps");
57 }
58 
59 
60 template<class Type, class DType, class LUType>
63 {
64  // --- Setup class containing solver performance data
65  SolverPerformance<Type> solverPerf
66  (
67  typeName,
68  this->fieldName_
69  );
70 
71  label nIter = 0;
72 
73  // If the nSweeps_ is negative do a fixed number of sweeps
74  if (nSweeps_ < 0)
75  {
78  (
79  this->fieldName_,
80  this->matrix_,
81  this->controlDict_
82  );
83 
84  smootherPtr->smooth(psi, -nSweeps_);
85 
86  nIter -= nSweeps_;
87  }
88  else
89  {
90  Type normFactor = Zero;
91 
92  {
93  Field<Type> Apsi(psi.size());
94  Field<Type> temp(psi.size());
95 
96  // Calculate A.psi
97  this->matrix_.Amul(Apsi, psi);
98 
99  // Calculate normalisation factor
100  normFactor = this->normFactor(psi, Apsi, temp);
101 
102  // Calculate residual magnitude
103  solverPerf.initialResidual() = cmptDivide
104  (
105  gSumCmptMag(this->matrix_.source() - Apsi),
106  normFactor
107  );
108  solverPerf.finalResidual() = solverPerf.initialResidual();
109  }
110 
112  {
113  Info<< " Normalisation factor = " << normFactor << endl;
114  }
115 
116 
117  // Check convergence, solve if not converged
118  if
119  (
120  this->minIter_ > 0
121  || !solverPerf.checkConvergence(this->tolerance_, this->relTol_)
122  )
123  {
126  (
127  this->fieldName_,
128  this->matrix_,
129  this->controlDict_
130  );
131 
132  // Smoothing loop
133  do
134  {
135  smootherPtr->smooth
136  (
137  psi,
138  nSweeps_
139  );
140 
141  // Calculate the residual to check convergence
142  solverPerf.finalResidual() = cmptDivide
143  (
144  gSumCmptMag(this->matrix_.residual(psi)),
145  normFactor
146  );
147  } while
148  (
149  (
150  (nIter += nSweeps_) < this->maxIter_
151  && !solverPerf.checkConvergence(this->tolerance_, this->relTol_)
152  )
153  || nIter < this->minIter_
154  );
155  }
156  }
157 
158  solverPerf.nIterations() =
160 
161  return solverPerf;
162 }
163 
164 
165 // ************************************************************************* //
Pre-declare SubField and related Field type.
Definition: Field.H:83
static autoPtr< smoother > New(const word &fieldName, const LduMatrix< Type, DType, LUType > &matrix, const dictionary &smootherDict)
Return a new smoother.
virtual void readControls()
Read the control parameters from the controlDict_.
LduMatrix is a general matrix class in which the coefficients are stored as three arrays,...
Definition: LduMatrix.H:85
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
SmoothSolver(const word &fieldName, const LduMatrix< Type, DType, LUType > &matrix, const dictionary &solverDict)
Construct from matrix components and solver data dictionary.
Definition: SmoothSolver.C:32
virtual void readControls()
Read the control parameters from the controlDict_.
Definition: SmoothSolver.C:53
virtual SolverPerformance< Type > solve(Field< Type > &psi) const
Solve the matrix with this solver.
Definition: SmoothSolver.C:62
SolverPerformance is the class returned by the LduMatrix solver containing performance statistics.
const labelType & nIterations() const
Return number of iterations.
const Type & initialResidual() const
Return initial residual.
const Type & finalResidual() const
Return final residual.
bool checkConvergence(const Type &tolerance, const Type &relTolerance)
Check, store and return convergence.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
A class representing the concept of 1 (scalar(1)) used to avoid unnecessary manipulations for objects...
Definition: one.H:51
Traits class for primitives.
Definition: pTraits.H:53
Abstract base class for run-time selectable region solvers.
Definition: solver.H:56
A class for handling words, derived from string.
Definition: word.H:62
const volScalarField & psi
static const zero Zero
Definition: zero.H:97
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
messageStream Info
dimensioned< Type > cmptDivide(const dimensioned< Type > &, const dimensioned< Type > &)
Type gSumCmptMag(const UList< Type > &f, const label comm)