kLowReWallFunctionFvPatchScalarField.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) 2012-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 
27 #include "turbulenceModel.H"
28 #include "fvPatchFieldMapper.H"
29 #include "volFields.H"
31 #include "wallFvPatch.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 
38 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
39 
41 {
42  if (!isA<wallFvPatch>(patch()))
43  {
45  << "Invalid wall function specification" << nl
46  << " Patch type for patch " << patch().name()
47  << " must be wall" << nl
48  << " Current patch type is " << patch().type() << nl << endl
49  << abort(FatalError);
50  }
51 }
52 
53 
55 (
56  const scalar kappa,
57  const scalar E
58 )
59 {
60  scalar ypl = 11.0;
61 
62  for (int i=0; i<10; i++)
63  {
64  ypl = log(max(E*ypl, 1))/kappa;
65  }
66 
67  return ypl;
68 }
69 
70 
71 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
72 
74 (
75  const fvPatch& p,
77 )
78 :
80  Cmu_(0.09),
81  kappa_(0.41),
82  E_(9.8),
83  Ceps2_(1.9),
85 {
86  checkType();
87 }
88 
89 
91 (
93  const fvPatch& p,
95  const fvPatchFieldMapper& mapper
96 )
97 :
98  fixedValueFvPatchField<scalar>(ptf, p, iF, mapper),
99  Cmu_(ptf.Cmu_),
100  kappa_(ptf.kappa_),
101  E_(ptf.E_),
102  Ceps2_(ptf.Ceps2_),
103  yPlusLam_(ptf.yPlusLam_)
104 {
105  checkType();
106 }
107 
108 
110 (
111  const fvPatch& p,
113  const dictionary& dict
114 )
115 :
117  Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
118  kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
119  E_(dict.lookupOrDefault<scalar>("E", 9.8)),
120  Ceps2_(dict.lookupOrDefault<scalar>("Ceps2", 1.9)),
122 {
123  checkType();
124 }
125 
126 
128 (
130 )
131 :
133  Cmu_(kwfpsf.Cmu_),
134  kappa_(kwfpsf.kappa_),
135  E_(kwfpsf.E_),
136  Ceps2_(kwfpsf.Ceps2_),
137  yPlusLam_(kwfpsf.yPlusLam_)
138 {
139  checkType();
140 }
141 
142 
144 (
147 )
148 :
149  fixedValueFvPatchField<scalar>(kwfpsf, iF),
150  Cmu_(kwfpsf.Cmu_),
151  kappa_(kwfpsf.kappa_),
152  E_(kwfpsf.E_),
153  Ceps2_(kwfpsf.Ceps2_),
154  yPlusLam_(kwfpsf.yPlusLam_)
155 {
156  checkType();
157 }
158 
159 
160 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
161 
163 {
164  if (updated())
165  {
166  return;
167  }
168 
169  const label patchi = patch().index();
170 
171  const turbulenceModel& turbModel = db().lookupObject<turbulenceModel>
172  (
174  (
176  internalField().group()
177  )
178  );
179  const scalarField& y = turbModel.y()[patchi];
180 
181  const tmp<volScalarField> tk = turbModel.k();
182  const volScalarField& k = tk();
183 
184  const tmp<scalarField> tnuw = turbModel.nu(patchi);
185  const scalarField& nuw = tnuw();
186 
187  const scalar Cmu25 = pow025(Cmu_);
188 
189  scalarField& kw = *this;
190 
191  // Set k wall values
192  forAll(kw, facei)
193  {
194  label celli = patch().faceCells()[facei];
195 
196  scalar uTau = Cmu25*sqrt(k[celli]);
197 
198  scalar yPlus = uTau*y[facei]/nuw[facei];
199 
200  if (yPlus > yPlusLam_)
201  {
202  scalar Ck = -0.416;
203  scalar Bk = 8.366;
204  kw[facei] = Ck/kappa_*log(yPlus) + Bk;
205  }
206  else
207  {
208  scalar C = 11.0;
209  scalar Cf = (1.0/sqr(yPlus + C) + 2.0*yPlus/pow3(C) - 1.0/sqr(C));
210  kw[facei] = 2400.0/sqr(Ceps2_)*Cf;
211  }
212 
213  kw[facei] *= sqr(uTau);
214  }
215 
216  // Limit kw to avoid failure of the turbulence model due to division by kw
217  kw = max(kw, small);
218 
220 
221  // TODO: perform averaging for cells sharing more than one boundary face
222 }
223 
224 
226 (
227  const Pstream::commsTypes commsType
228 )
229 {
231 }
232 
233 
235 {
236  os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
237  os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
238  os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
239  os.writeKeyword("Ceps2") << Ceps2_ << token::END_STATEMENT << nl;
241 }
242 
243 
244 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
245 
247 (
250 );
251 
252 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
253 
254 } // End namespace Foam
255 
256 // ************************************************************************* //
const char *const group
Group name for atomic constants.
dictionary dict
Graphite solid properties.
Definition: C.H:48
virtual void checkType()
Check the type of the patch.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
virtual void write(Ostream &) const
Write.
dimensionedScalar log(const dimensionedScalar &ds)
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:179
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
bool updated() const
Return true if the boundary condition has already been updated.
Definition: fvPatchField.H:371
commsTypes
Types of communications.
Definition: UPstream.H:64
dimensionedSymmTensor sqr(const dimensionedVector &dv)
This boundary condition provides a turbulence kinetic energy wall function condition for low- and hig...
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
dimensionedScalar pow025(const dimensionedScalar &ds)
virtual void evaluate(const Pstream::commsTypes)
Evaluate the patchField.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
label k
Boltzmann constant.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Abstract base class for turbulence models (RAS, LES and laminar).
const word & name() const
Return name.
Definition: fvPatch.H:149
Macros for easy insertion into run-time selection tables.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
scalar uTau
scalar y
static const word propertiesName
Default name of the turbulence properties dictionary.
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:93
static word groupName(Name name, const word &group)
Foam::fvPatchFieldMapper.
const nearWallDist & y() const
Return the near wall distances.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:340
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:265
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field, sets Updated to false.
Definition: fvPatchField.C:331
dimensionedScalar pow3(const dimensionedScalar &ds)
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
scalar yPlusLam_
Y+ at the edge of the laminar sublayer.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
scalar yPlus
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:311
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:197
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
kLowReWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Namespace for OpenFOAM.
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:347
scalar yPlusLam(const scalar kappa, const scalar E)
Calculate the Y+ at the edge of the laminar sublayer.