limitedSnGrad.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-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 \*---------------------------------------------------------------------------*/
25 
26 #include "limitedSnGrad.H"
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace fv
38 {
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
42 template<class Type>
44 :
45  snGradScheme<Type>(mesh),
46  correctedScheme_(new correctedSnGrad<Type>(this->mesh())),
47  limitCoeff_(1)
48 {}
49 
50 
51 template<class Type>
53 :
54  snGradScheme<Type>(mesh),
55  correctedScheme_(lookupCorrectedScheme(schemeData))
56 {
57  if (limitCoeff_ < 0 || limitCoeff_ > 1)
58  {
60  (
61  schemeData
62  ) << "limitCoeff is specified as " << limitCoeff_
63  << " but should be >= 0 && <= 1"
64  << exit(FatalIOError);
65  }
66 }
67 
68 
69 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
70 
71 template<class Type>
73 {}
74 
75 
76 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
77 
78 template<class Type>
80 (
81  const VolField<Type>& vf
82 ) const
83 {
84  return correctedScheme_->deltaCoeffs(vf);
85 }
86 
87 
88 template<class Type>
91 (
92  const VolField<Type>& vf
93 ) const
94 {
95  const SurfaceField<Type> corr
96  (
97  correctedScheme_().correction(vf)
98  );
99 
101  (
102  min
103  (
104  limitCoeff_
105  *mag(snGradScheme<Type>::snGrad(vf, deltaCoeffs(vf), "SndGrad"))
106  /(
107  (1 - limitCoeff_)*mag(corr)
108  + dimensionedScalar(corr.dimensions(), small)
109  ),
111  )
112  );
113 
114  if (fv::debug)
115  {
117  << "limiter min: " << min(limiter.primitiveField())
118  << " max: "<< max(limiter.primitiveField())
119  << " avg: " << average(limiter.primitiveField()) << endl;
120  }
121 
122  return limiter*corr;
123 }
124 
125 
126 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
127 
128 } // End namespace fv
129 
130 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
131 
132 } // End namespace Foam
133 
134 // ************************************************************************* //
const dimensionSet & dimensions() const
Return dimensions.
Generic GeometricField class.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Simple central-difference snGrad scheme with non-orthogonal correction.
virtual ~limitedSnGrad()
Destructor.
Definition: limitedSnGrad.C:72
virtual tmp< SurfaceField< Type > > correction(const VolField< Type > &) const
Return the explicit correction to the limitedSnGrad.
Definition: limitedSnGrad.C:91
virtual tmp< surfaceScalarField > deltaCoeffs(const VolField< Type > &) const
Return the interpolation weighting factors for the given field.
Definition: limitedSnGrad.C:80
limitedSnGrad(const fvMesh &mesh)
Construct from mesh.
Definition: limitedSnGrad.C:43
Abstract base class for snGrad schemes.
Definition: snGradScheme.H:63
A class for managing temporary objects.
Definition: tmp.H:55
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:346
#define InfoInFunction
Report an information message using Foam::Info.
void limiter(surfaceScalarField &lambda, const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phiBD, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su, const PsiMaxType &psiMax, const PsiMinType &psiMin)
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
const dimensionSet dimless
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
IOerror FatalIOError
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
labelList fv(nPoints)
Foam::surfaceFields.