kLowReWallFunctionFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2012-2016 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 // * * * * * * * * * * * * * Private 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 
54 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
55 
57 (
58  const scalar kappa,
59  const scalar E
60 )
61 {
62  scalar ypl = 11.0;
63 
64  for (int i=0; i<10; i++)
65  {
66  ypl = log(max(E*ypl, 1))/kappa;
67  }
68 
69  return ypl;
70 }
71 
72 
73 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
74 
76 (
77  const fvPatch& p,
79 )
80 :
82  Cmu_(0.09),
83  kappa_(0.41),
84  E_(9.8),
85  Ceps2_(1.9),
87 {
88  checkType();
89 }
90 
91 
93 (
95  const fvPatch& p,
97  const fvPatchFieldMapper& mapper
98 )
99 :
100  fixedValueFvPatchField<scalar>(ptf, p, iF, mapper),
101  Cmu_(ptf.Cmu_),
102  kappa_(ptf.kappa_),
103  E_(ptf.E_),
104  Ceps2_(ptf.Ceps2_),
105  yPlusLam_(ptf.yPlusLam_)
106 {
107  checkType();
108 }
109 
110 
112 (
113  const fvPatch& p,
115  const dictionary& dict
116 )
117 :
119  Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
120  kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
121  E_(dict.lookupOrDefault<scalar>("E", 9.8)),
122  Ceps2_(dict.lookupOrDefault<scalar>("Ceps2", 1.9)),
124 {
125  checkType();
126 }
127 
128 
130 (
132 )
133 :
135  Cmu_(kwfpsf.Cmu_),
136  kappa_(kwfpsf.kappa_),
137  E_(kwfpsf.E_),
138  Ceps2_(kwfpsf.Ceps2_),
139  yPlusLam_(kwfpsf.yPlusLam_)
140 {
141  checkType();
142 }
143 
144 
146 (
149 )
150 :
151  fixedValueFvPatchField<scalar>(kwfpsf, iF),
152  Cmu_(kwfpsf.Cmu_),
153  kappa_(kwfpsf.kappa_),
154  E_(kwfpsf.E_),
155  Ceps2_(kwfpsf.Ceps2_),
156  yPlusLam_(kwfpsf.yPlusLam_)
157 {
158  checkType();
159 }
160 
161 
162 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
163 
165 {
166  if (updated())
167  {
168  return;
169  }
170 
171  const label patchi = patch().index();
172 
173  const turbulenceModel& turbModel = db().lookupObject<turbulenceModel>
174  (
176  (
178  internalField().group()
179  )
180  );
181  const scalarField& y = turbModel.y()[patchi];
182 
183  const tmp<volScalarField> tk = turbModel.k();
184  const volScalarField& k = tk();
185 
186  const tmp<scalarField> tnuw = turbModel.nu(patchi);
187  const scalarField& nuw = tnuw();
188 
189  const scalar Cmu25 = pow025(Cmu_);
190 
191  scalarField& kw = *this;
192 
193  // Set k wall values
194  forAll(kw, facei)
195  {
196  label faceCelli = patch().faceCells()[facei];
197 
198  scalar uTau = Cmu25*sqrt(k[faceCelli]);
199 
200  scalar yPlus = uTau*y[facei]/nuw[facei];
201 
202  if (yPlus > yPlusLam_)
203  {
204  scalar Ck = -0.416;
205  scalar Bk = 8.366;
206  kw[facei] = Ck/kappa_*log(yPlus) + Bk;
207  }
208  else
209  {
210  scalar C = 11.0;
211  scalar Cf = (1.0/sqr(yPlus + C) + 2.0*yPlus/pow3(C) - 1.0/sqr(C));
212  kw[facei] = 2400.0/sqr(Ceps2_)*Cf;
213  }
214 
215  kw[facei] *= sqr(uTau);
216  }
217 
218  // Limit kw to avoid failure of the turbulence model due to division by kw
219  kw = max(kw, SMALL);
220 
222 
223  // TODO: perform averaging for cells sharing more than one boundary face
224 }
225 
226 
228 (
229  const Pstream::commsTypes commsType
230 )
231 {
233 }
234 
235 
237 {
238  os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
239  os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
240  os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
241  os.writeKeyword("Ceps2") << Ceps2_ << token::END_STATEMENT << nl;
243 }
244 
245 
246 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
247 
249 (
252 );
253 
254 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
255 
256 } // End namespace Foam
257 
258 // ************************************************************************* //
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:339
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::blocking)
Evaluate the patch field, sets Updated to false.
Definition: fvPatchField.C:332
const char *const group
Group name for atomic constants.
dictionary dict
Graphite solid properties.
Definition: C.H:57
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
dimensionedScalar log(const dimensionedScalar &ds)
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
virtual void write(Ostream &) const
Write.
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.
const word & name() const
Return name.
Definition: fvPatch.H:149
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
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:65
label k
Boltzmann constant.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Abstract base class for turbulence models (RAS, LES and laminar).
Macros for easy insertion into run-time selection tables.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
scalar uTau
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
scalar y
static const word propertiesName
Default name of the turbulence properties dictionary.
static word groupName(Name name, const word &group)
Foam::fvPatchFieldMapper.
bool updated() const
Return true if the boundary condition has already been updated.
Definition: fvPatchField.H:370
errorManip< error > abort(error &err)
Definition: errorManip.H:131
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:93
static const char nl
Definition: Ostream.H:262
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:198
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:179
dimensionedScalar pow3(const dimensionedScalar &ds)
label patchi
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:312
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:54
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
kLowReWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
const nearWallDist & y() const
Return the near wall distances.
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:346
Namespace for OpenFOAM.
scalar yPlusLam(const scalar kappa, const scalar E)
Calculate the Y+ at the edge of the laminar sublayer.