alphatJayatillekeWallFunctionFvPatchScalarField.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2017-2022 OpenCFD Ltd
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
30#include "fvPatchFieldMapper.H"
31#include "volFields.H"
32#include "wallFvPatch.H"
35
36// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37
38namespace Foam
39{
40namespace incompressible
41{
42
43// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
44
47
48// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
49
51{
52 if (!isA<wallFvPatch>(patch()))
53 {
55 << "Invalid wall function specification" << nl
56 << " Patch type for patch " << patch().name()
57 << " must be wall" << nl
58 << " Current patch type is " << patch().type() << nl << endl
59 << abort(FatalError);
60 }
61}
62
63
65(
66 const turbulenceModel& turbModel
67) const
68{
69 const label patchi = patch().index();
70 const tmp<volScalarField> tnut = turbModel.nut();
71 const volScalarField& nut = tnut();
72
73 if (isA<nutWallFunctionFvPatchScalarField>(nut.boundaryField()[patchi]))
74 {
75 const auto& nutPf =
76 dynamic_cast<const nutWallFunctionFvPatchScalarField&>
77 (
78 nut.boundaryField()[patchi]
79 );
80
81 return nutPf.yPlus();
82 }
83 else
84 {
85 const scalarField& y = turbModel.y()[patchi];
86 const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];
87
88 return
89 y*sqrt(turbModel.nuEff(patchi)*mag(Uw.snGrad()))
90 /turbModel.nu(patchi);
91 }
92}
93
94
96(
97 const scalar Prat
98) const
99{
100 return 9.24*(pow(Prat, 0.75) - 1.0)*(1.0 + 0.28*exp(-0.007*Prat));
101}
102
103
105(
106 const scalar P,
107 const scalar Prat
108) const
109{
110 scalar ypt = 11;
111
112 for (int iter = 0; iter < maxIters_; ++iter)
113 {
114 const scalar f = ypt - (log(E_*ypt)/kappa_ + P)/Prat;
115 const scalar df = 1.0 - 1.0/(ypt*kappa_*Prat);
116 const scalar yptNew = ypt - f/df;
117
118 if (yptNew < VSMALL)
119 {
120 return 0;
121 }
122 else if (mag(yptNew - ypt) < tolerance_)
123 {
124 return yptNew;
125 }
126 else
127 {
128 ypt = yptNew;
129 }
130 }
131
132 return ypt;
133}
134
135
136// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
137
138alphatJayatillekeWallFunctionFvPatchScalarField::
139alphatJayatillekeWallFunctionFvPatchScalarField
140(
141 const fvPatch& p,
143)
144:
145 fixedValueFvPatchScalarField(p, iF),
146 Prt_(0.85),
147 kappa_(0.41),
148 E_(9.8)
149{
150 checkType();
151}
152
153
154alphatJayatillekeWallFunctionFvPatchScalarField::
155alphatJayatillekeWallFunctionFvPatchScalarField
156(
158 const fvPatch& p,
160 const fvPatchFieldMapper& mapper
161)
162:
163 fixedValueFvPatchScalarField(ptf, p, iF, mapper),
164 Prt_(ptf.Prt_),
165 kappa_(ptf.kappa_),
166 E_(ptf.E_)
167{
168 checkType();
169}
170
171
172alphatJayatillekeWallFunctionFvPatchScalarField::
173alphatJayatillekeWallFunctionFvPatchScalarField
174(
175 const fvPatch& p,
177 const dictionary& dict
178)
179:
180 fixedValueFvPatchScalarField(p, iF, dict),
181 Prt_(dict.get<scalar>("Prt")), // force read to avoid ambiguity
182 kappa_(dict.getOrDefault<scalar>("kappa", 0.41)),
183 E_(dict.getOrDefault<scalar>("E", 9.8))
184{
185 checkType();
186}
187
188
189alphatJayatillekeWallFunctionFvPatchScalarField::
190alphatJayatillekeWallFunctionFvPatchScalarField
191(
193)
194:
195 fixedValueFvPatchScalarField(wfpsf),
196 Prt_(wfpsf.Prt_),
197 kappa_(wfpsf.kappa_),
198 E_(wfpsf.E_)
199{
200 checkType();
201}
202
203
204alphatJayatillekeWallFunctionFvPatchScalarField::
205alphatJayatillekeWallFunctionFvPatchScalarField
206(
209)
210:
211 fixedValueFvPatchScalarField(wfpsf, iF),
212 Prt_(wfpsf.Prt_),
213 kappa_(wfpsf.kappa_),
214 E_(wfpsf.E_)
215{
216 checkType();
217}
218
219
220// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
221
223{
224 if (updated())
225 {
226 return;
227 }
228
229 const label patchi = patch().index();
230
231 // Retrieve turbulence properties from model
232
233 const auto& turbModel = db().lookupObject<turbulenceModel>
234 (
236 (
238 internalField().group()
239 )
240 );
241
242 const scalarField yPlusp(yPlus(turbModel));
243
244 const tmp<volScalarField> tnu = turbModel.nu();
245 const volScalarField& nu = tnu();
246 const scalarField& nuw = nu.boundaryField()[patchi];
247
248 const auto& transportProperties =
249 db().lookupObject<IOdictionary>("transportProperties");
250
251 // Molecular Prandtl number
252 const scalar Pr
253 (
255 );
256
257 // Populate boundary values
258 scalarField& alphatw = *this;
259 forAll(alphatw, facei)
260 {
261 const scalar yPlus = yPlusp[facei];
262
263 // Molecular-to-turbulent Prandtl number ratio
264 const scalar Prat = Pr/Prt_;
265
266 // Thermal sublayer thickness
267 const scalar P = Psmooth(Prat);
268 const scalar yPlusTherm = this->yPlusTherm(P, Prat);
269
270 // Update turbulent thermal conductivity
271 if (yPlus > yPlusTherm)
272 {
273 const scalar nu = nuw[facei];
274 const scalar kt =
275 nu*(yPlus/(Prt_*(log(E_*yPlus)/kappa_ + P)) - 1.0/Pr);
276
277 alphatw[facei] = max(scalar(0), kt);
278 }
279 else
280 {
281 alphatw[facei] = 0.0;
282 }
283 }
284
286}
287
288
290{
292 os.writeEntry("Prt", Prt_);
293 os.writeEntryIfDifferent<scalar>("kappa", 0.41, kappa_);
294 os.writeEntryIfDifferent<scalar>("E", 9.8, E_);
295 writeEntry("value", os);
296}
297
298
299// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
300
302(
305);
306
307// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
308
309} // End namespace incompressible
310} // End namespace Foam
311
312// ************************************************************************* //
scalar y
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Boundary & boundaryField() const
Return const-reference to the boundary field.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Templated abstract base class for single-phase incompressible turbulence models.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:239
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Definition: Ostream.H:251
virtual tmp< volScalarField > nu() const
Return the laminar viscosity.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
virtual bool write()
Write the output fields.
Computes the near wall for turbulence models.
Definition: yPlus.H:160
A FieldMapper for finite-volume patch fields.
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:229
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
This boundary condition provides a kinematic turbulent thermal conductivity for using wall functions,...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
scalar yPlusTherm(const scalar P, const scalar Prat) const
Calculate y+ at the edge of the thermal laminar sublayer.
The class nutWallFunction is an abstract base class that hosts calculation methods and common functi...
virtual tmp< scalarField > yPlus() const =0
Calculate and return the yPlus at the boundary.
static const word propertiesName
Default name of the phase properties dictionary.
Definition: phaseSystem.H:290
A class for managing temporary objects.
Definition: tmp.H:65
const volVectorField & U() const
Access function to velocity field.
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
const nearWallDist & y() const
Return the near wall distances.
virtual tmp< volScalarField > nuEff() const =0
Return the effective viscosity.
volScalarField & p
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
scalar yPlus
scalar nut
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
Definition: fvPatchField.H:676
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensionedScalar exp(const dimensionedScalar &ds)
const dimensionSet dimless
Dimensionless.
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManip< error > abort(error &err)
Definition: errorManip.H:144
error FatalError
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
labelList f(nPoints)
volScalarField & nu
dictionary dict
IOdictionary transportProperties(IOobject("transportProperties", runTime.constant(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
dimensionedScalar Pr("Pr", dimless, laminarTransport)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333