CoBlended.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) 2012-2022 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::CoBlended
26 
27 Description
28  Two-scheme Courant number based blending interpolation scheme.
29 
30  Similar to localBlended but uses a blending factor computed from the
31  face-based Courant number and the lower and upper Courant number limits
32  supplied:
33  \f[
34  weight = 1 - max(min((Co - Co1)/(Co2 - Co1), 1), 0)
35  \f]
36  where
37  \vartable
38  Co1 | Courant number below which scheme1 is used
39  Co2 | Courant number above which scheme2 is used
40  \endvartable
41 
42  The weight applies to the first scheme and 1-weight to the second scheme.
43 
44  Example of the CoBlended scheme specification using LUST for Courant numbers
45  less than 1 and linearUpwind for Courant numbers greater than 10:
46  \verbatim
47  divSchemes
48  {
49  .
50  .
51  div(phi,U) Gauss CoBlended 1 LUST grad(U) 10 linearUpwind grad(U);
52  .
53  .
54  }
55  \endverbatim
56 
57 SourceFiles
58  CoBlended.C
59 
60 \*---------------------------------------------------------------------------*/
61 
62 #ifndef CoBlended_H
63 #define CoBlended_H
64 
66 #include "blendedSchemeBase.H"
67 #include "surfaceInterpolate.H"
68 
69 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
70 
71 namespace Foam
72 {
73 
74 /*---------------------------------------------------------------------------*\
75  Class CoBlended Declaration
76 \*---------------------------------------------------------------------------*/
77 
78 template<class Type>
79 class CoBlended
80 :
81  public surfaceInterpolationScheme<Type>,
82  public blendedSchemeBase<Type>
83 {
84  // Private Data
85 
86  //- Courant number below which scheme1 is used
87  const scalar Co1_;
88 
89  //- Scheme 1
91 
92  //- Courant number above which scheme2 is used
93  const scalar Co2_;
94 
95  //- Scheme 2
97 
98  //- The face-flux used to compute the face Courant number
99  const surfaceScalarField& faceFlux_;
100 
101 
102 public:
103 
104  //- Runtime type information
105  TypeName("CoBlended");
106 
107 
108  // Constructors
109 
110  //- Construct from mesh and Istream.
111  // The name of the flux field is read from the Istream and looked-up
112  // from the mesh objectRegistry
113  CoBlended
114  (
115  const fvMesh& mesh,
116  Istream& is
117  )
118  :
119  surfaceInterpolationScheme<Type>(mesh),
120  Co1_(readScalar(is)),
121  tScheme1_
122  (
123  surfaceInterpolationScheme<Type>::New(mesh, is)
124  ),
125  Co2_(readScalar(is)),
126  tScheme2_
127  (
128  surfaceInterpolationScheme<Type>::New(mesh, is)
129  ),
130  faceFlux_
131  (
132  mesh.lookupObject<surfaceScalarField>(word(is))
133  )
134  {
135  if (Co1_ < 0 || Co2_ < 0 || Co1_ >= Co2_)
136  {
138  << "coefficients = " << Co1_ << " and " << Co2_
139  << " should be > 0 and Co2 > Co1"
140  << exit(FatalIOError);
141  }
142  }
143 
144 
145  //- Construct from mesh, faceFlux and Istream
146  CoBlended
147  (
148  const fvMesh& mesh,
149  const surfaceScalarField& faceFlux,
150  Istream& is
151  )
152  :
154  Co1_(readScalar(is)),
155  tScheme1_
156  (
157  surfaceInterpolationScheme<Type>::New(mesh, faceFlux, is)
158  ),
159  Co2_(readScalar(is)),
160  tScheme2_
161  (
162  surfaceInterpolationScheme<Type>::New(mesh, faceFlux, is)
163  ),
164  faceFlux_(faceFlux)
165  {
166  if (Co1_ < 0 || Co2_ < 0 || Co1_ >= Co2_)
167  {
169  << "coefficients = " << Co1_ << " and " << Co2_
170  << " should be > 0 and Co2 > Co1"
171  << exit(FatalIOError);
172  }
173  }
174 
175  //- Disallow default bitwise copy construction
176  CoBlended(const CoBlended&) = delete;
177 
178 
179  // Member Functions
180 
181  //- Return the face-based blending factor
183  (
185  ) const
186  {
187  const fvMesh& mesh = this->mesh();
188  tmp<surfaceScalarField> tUflux = faceFlux_;
189 
190  if (faceFlux_.dimensions() == dimMassFlux)
191  {
192  // Currently assume that the density field
193  // corresponding to the mass-flux is named "rho"
194  const volScalarField& rho =
195  mesh.objectRegistry::template lookupObject<volScalarField>
196  ("rho");
197 
198  tUflux = faceFlux_/fvc::interpolate(rho);
199  }
200  else if (faceFlux_.dimensions() != dimFlux)
201  {
203  << "dimensions of faceFlux are not correct"
204  << exit(FatalError);
205  }
206 
208  (
209  vf.name() + "BlendingFactor",
210  scalar(1)
211  - max
212  (
213  min
214  (
215  (
216  mesh.time().deltaT()*mesh.deltaCoeffs()
217  *mag(tUflux)/mesh.magSf()
218  - Co1_
219  )/(Co2_ - Co1_),
220  scalar(1)
221  ),
222  scalar(0)
223  )
224  );
225  }
226 
227 
228  //- Return the interpolation weighting factors
230  weights
231  (
233  ) const
234  {
236 
237  return
238  bf*tScheme1_().weights(vf)
239  + (scalar(1) - bf)*tScheme2_().weights(vf);
240  }
241 
242 
243  //- Return the face-interpolate of the given cell field
244  // with explicit correction
247  (
249  ) const
250  {
252 
253  return
254  bf*tScheme1_().interpolate(vf)
255  + (scalar(1) - bf)*tScheme2_().interpolate(vf);
256  }
257 
258 
259  //- Return true if this scheme uses an explicit correction
260  virtual bool corrected() const
261  {
262  return tScheme1_().corrected() || tScheme2_().corrected();
263  }
264 
265 
266  //- Return the explicit correction to the face-interpolate
267  // for the given field
269  correction
270  (
272  ) const
273  {
275 
276  if (tScheme1_().corrected())
277  {
278  if (tScheme2_().corrected())
279  {
280  return
281  (
282  bf
283  * tScheme1_().correction(vf)
284  + (scalar(1) - bf)
285  * tScheme2_().correction(vf)
286  );
287  }
288  else
289  {
290  return
291  (
292  bf
293  * tScheme1_().correction(vf)
294  );
295  }
296  }
297  else if (tScheme2_().corrected())
298  {
299  return
300  (
301  (scalar(1) - bf)
302  * tScheme2_().correction(vf)
303  );
304  }
305  else
306  {
308  (
309  nullptr
310  );
311  }
312  }
313 
314 
315  // Member Operators
316 
317  //- Disallow default bitwise assignment
318  void operator=(const CoBlended&) = delete;
319 };
320 
321 
322 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
323 
324 } // End namespace Foam
325 
326 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
327 
328 #endif
329 
330 // ************************************************************************* //
void operator=(const CoBlended &)=delete
Disallow default bitwise assignment.
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
const word & name() const
Return name.
Definition: IOobject.H:315
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
virtual tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > correction(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the explicit correction to the face-interpolate.
Definition: CoBlended.H:277
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
static tmp< GeometricField< scalar, fvsPatchField, surfaceMesh > > New(const word &name, const Internal &, const PtrList< fvsPatchField< scalar >> &)
Return a temporary field constructed from name,.
const fvMesh & mesh() const
Return mesh reference.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:372
const dimensionSet & dimensions() const
Return dimensions.
A class for handling words, derived from string.
Definition: word.H:59
static tmp< surfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
const dimensionSet dimFlux
Two-scheme Courant number based blending interpolation scheme.
Definition: CoBlended.H:86
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:75
TypeName("CoBlended")
Runtime type information.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
const dimensionSet dimMassFlux
virtual tmp< surfaceScalarField > blendingFactor(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the face-based blending factor.
Definition: CoBlended.H:190
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:318
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
virtual bool corrected() const
Return true if this scheme uses an explicit correction.
Definition: CoBlended.H:267
dimensioned< scalar > mag(const dimensioned< Type > &)
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the face-interpolate of the given cell field.
Definition: CoBlended.H:254
A class for managing temporary objects.
Definition: PtrList.H:53
CoBlended(const fvMesh &mesh, Istream &is)
Construct from mesh and Istream.
Definition: CoBlended.H:121
Abstract base class for surface interpolation schemes.
const surfaceScalarField & deltaCoeffs() const
Return reference to cell-centre difference coefficients.
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:46
tmp< surfaceScalarField > weights(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the interpolation weighting factors.
Definition: CoBlended.H:238
Namespace for OpenFOAM.
IOerror FatalIOError