limitedSnGrad.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::fv::limitedSnGrad
26 
27 Description
28  Run-time selected snGrad scheme with limited non-orthogonal correction.
29 
30  The limiter is controlled by a coefficient with a value between 0 and 1
31  which when 0 switches the correction off and the scheme behaves as
32  uncorrectedSnGrad, when set to 1 the full correction of the selected scheme
33  is used and when set to 0.5 the limiter is calculated such that the
34  non-orthogonal contribution does not exceed the orthogonal part.
35 
36  Format:
37  limited <corrected scheme> <coefficient>;
38 
39  or
40 
41  limited <coefficient>; // Backward compatibility
42 
43 SourceFiles
44  limitedSnGrad.C
45 
46 \*---------------------------------------------------------------------------*/
47 
48 #ifndef limitedSnGrad_H
49 #define limitedSnGrad_H
50 
51 #include "correctedSnGrad.H"
52 
53 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
54 
55 namespace Foam
56 {
57 
58 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
59 
60 namespace fv
61 {
62 
63 /*---------------------------------------------------------------------------*\
64  Class limitedSnGrad Declaration
65 \*---------------------------------------------------------------------------*/
66 
67 template<class Type>
68 class limitedSnGrad
69 :
70  public snGradScheme<Type>
71 {
72  // Private Data
73 
74  tmp<snGradScheme<Type>> correctedScheme_;
75 
76  scalar limitCoeff_;
77 
78 
79  // Private Member Functions
80 
81  //- Lookup function for the corrected to support backward compatibility
82  // of dictionary specification
83  tmp<snGradScheme<Type>> lookupCorrectedScheme(Istream& schemeData)
84  {
85  token nextToken(schemeData);
86 
87  if (nextToken.isNumber())
88  {
89  limitCoeff_ = nextToken.number();
91  (
92  new correctedSnGrad<Type>(this->mesh())
93  );
94  }
95  else
96  {
97  schemeData.putBack(nextToken);
98  tmp<snGradScheme<Type>> tcorrectedScheme
99  (
100  fv::snGradScheme<Type>::New(this->mesh(), schemeData)
101  );
102 
103  schemeData >> limitCoeff_;
104 
105  return tcorrectedScheme;
106  }
107  }
108 
109 
110 public:
111 
112  //- Runtime type information
113  TypeName("limited");
114 
115 
116  // Constructors
117 
118  //- Construct from mesh
119  limitedSnGrad(const fvMesh& mesh)
120  :
121  snGradScheme<Type>(mesh),
122  correctedScheme_(new correctedSnGrad<Type>(this->mesh())),
123  limitCoeff_(1)
124  {}
125 
126 
127  //- Construct from mesh and data stream
128  limitedSnGrad(const fvMesh& mesh, Istream& schemeData)
129  :
130  snGradScheme<Type>(mesh),
131  correctedScheme_(lookupCorrectedScheme(schemeData))
132  {
133  if (limitCoeff_ < 0 || limitCoeff_ > 1)
134  {
136  (
137  schemeData
138  ) << "limitCoeff is specified as " << limitCoeff_
139  << " but should be >= 0 && <= 1"
140  << exit(FatalIOError);
141  }
142  }
143 
144 
145  //- Destructor
146  virtual ~limitedSnGrad();
147 
148 
149  // Member Functions
150 
151  //- Return the interpolation weighting factors for the given field
153  (
155  ) const
156  {
157  return this->mesh().nonOrthDeltaCoeffs();
158  }
159 
160  //- Return true if this scheme uses an explicit correction
161  virtual bool corrected() const
162  {
163  return true;
164  }
165 
166  //- Return the explicit correction to the limitedSnGrad
167  // for the given field
170 
171 
172  // Member Operators
173 
174  //- Disallow default bitwise assignment
175  void operator=(const limitedSnGrad&) = delete;
176 };
177 
178 
179 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
180 
181 } // End namespace fv
182 
183 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
184 
185 } // End namespace Foam
186 
187 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
188 
189 #ifdef NoRepository
190  #include "limitedSnGrad.C"
191 #endif
192 
193 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
194 
195 #endif
196 
197 // ************************************************************************* //
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
scalar number() const
Definition: tokenI.H:503
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
A token holds items read from Istream.
Definition: token.H:72
void putBack(const token &)
Put back token.
Definition: Istream.C:30
Generic GeometricField class.
TypeName("limited")
Runtime type information.
limitedSnGrad(const fvMesh &mesh)
Construct from mesh.
bool isNumber() const
Definition: tokenI.H:498
const surfaceScalarField & nonOrthDeltaCoeffs() const
Return reference to non-orthogonal cell-centre difference.
virtual tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > correction(const GeometricField< Type, fvPatchField, volMesh > &) const
Return the explicit correction to the limitedSnGrad.
Definition: limitedSnGrad.C:54
Run-time selected snGrad scheme with limited non-orthogonal correction.
Definition: limitedSnGrad.H:67
labelList fv(nPoints)
Simple central-difference snGrad scheme with non-orthogonal correction.
void operator=(const limitedSnGrad &)=delete
Disallow default bitwise assignment.
virtual bool corrected() const
Return true if this scheme uses an explicit correction.
Abstract base class for snGrad schemes.
Definition: snGradScheme.H:60
const fvMesh & mesh() const
Return mesh reference.
Definition: snGradScheme.H:117
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:335
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual tmp< surfaceScalarField > deltaCoeffs(const GeometricField< Type, fvPatchField, volMesh > &) const
Return the interpolation weighting factors for the given field.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual ~limitedSnGrad()
Destructor.
Definition: limitedSnGrad.C:45
Namespace for OpenFOAM.
IOerror FatalIOError