nutURoughWallFunctionFvPatchScalarField.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) 2019-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
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 "turbulenceModel.H"
31 #include "fvPatchFieldMapper.H"
32 #include "volFields.H"
34 
35 
36 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
37 
39 Foam::nutURoughWallFunctionFvPatchScalarField::calcNut() const
40 {
41  const label patchi = patch().index();
42 
43  const turbulenceModel& turbModel = db().lookupObject<turbulenceModel>
44  (
46  (
48  internalField().group()
49  )
50  );
51  const scalarField& y = turbModel.y()[patchi];
52  const fvPatchVectorField& Uw = U(turbModel).boundaryField()[patchi];
53  const tmp<scalarField> tnuw = turbModel.nu(patchi);
54  const scalarField& nuw = tnuw();
55 
56  // The flow velocity at the adjacent cell centre
57  const scalarField magUp(mag(Uw.patchInternalField() - Uw));
58 
59  tmp<scalarField> tyPlus = calcYPlus(magUp);
60  scalarField& yPlus = tyPlus.ref();
61 
62  tmp<scalarField> tnutw(new scalarField(patch().size(), Zero));
63  scalarField& nutw = tnutw.ref();
64 
65  forAll(yPlus, facei)
66  {
67  if (yPlusLam_ < yPlus[facei])
68  {
69  const scalar Re = magUp[facei]*y[facei]/nuw[facei] + ROOTVSMALL;
70  nutw[facei] = nuw[facei]*(sqr(yPlus[facei])/Re - 1);
71  }
72  }
73 
74  return tnutw;
75 }
76 
77 
79 Foam::nutURoughWallFunctionFvPatchScalarField::calcYPlus
80 (
81  const scalarField& magUp
82 ) const
83 {
84  const label patchi = patch().index();
85 
86  const turbulenceModel& turbModel = db().lookupObject<turbulenceModel>
87  (
89  (
91  internalField().group()
92  )
93  );
94  const scalarField& y = turbModel.y()[patchi];
95  const tmp<scalarField> tnuw = turbModel.nu(patchi);
96  const scalarField& nuw = tnuw();
97 
98  tmp<scalarField> tyPlus(new scalarField(patch().size(), Zero));
99  scalarField& yPlus = tyPlus.ref();
100 
101  if (0.0 < roughnessHeight_)
102  {
103  // Rough Walls
104  const scalar c_1 = 1/(90 - 2.25) + roughnessConstant_;
105  static const scalar c_2 = 2.25/(90 - 2.25);
106  static const scalar c_3 = 2.0*atan(1.0)/log(90/2.25);
107  static const scalar c_4 = c_3*log(2.25);
108 
109  //if (KsPlusBasedOnYPlus_)
110  {
111  // If KsPlus is based on YPlus the extra term added to the law
112  // of the wall will depend on yPlus
113  forAll(yPlus, facei)
114  {
115  const scalar magUpara = magUp[facei];
116  const scalar Re = magUpara*y[facei]/nuw[facei];
117  const scalar kappaRe = kappa_*Re;
118 
119  scalar yp = yPlusLam_;
120  const scalar ryPlusLam = 1.0/yp;
121 
122  int iter = 0;
123  scalar yPlusLast = 0.0;
124  scalar dKsPlusdYPlus = roughnessHeight_/y[facei];
125 
126  // Additional tuning parameter - nominally = 1
127  dKsPlusdYPlus *= roughnessFactor_;
128 
129  do
130  {
131  yPlusLast = yp;
132 
133  // The non-dimensional roughness height
134  scalar KsPlus = yp*dKsPlusdYPlus;
135 
136  // The extra term in the law-of-the-wall
137  scalar G = 0.0;
138 
139  scalar yPlusGPrime = 0.0;
140 
141  if (KsPlus >= 90)
142  {
143  const scalar t_1 = 1 + roughnessConstant_*KsPlus;
144  G = log(t_1);
145  yPlusGPrime = roughnessConstant_*KsPlus/t_1;
146  }
147  else if (KsPlus > 2.25)
148  {
149  const scalar t_1 = c_1*KsPlus - c_2;
150  const scalar t_2 = c_3*log(KsPlus) - c_4;
151  const scalar sint_2 = sin(t_2);
152  const scalar logt_1 = log(t_1);
153  G = logt_1*sint_2;
154  yPlusGPrime =
155  (c_1*sint_2*KsPlus/t_1) + (c_3*logt_1*cos(t_2));
156  }
157 
158  scalar denom = 1.0 + log(E_*yp) - G - yPlusGPrime;
159  if (mag(denom) > VSMALL)
160  {
161  yp = (kappaRe + yp*(1 - yPlusGPrime))/denom;
162  }
163  } while
164  (
165  mag(ryPlusLam*(yp - yPlusLast)) > tolerance_
166  && ++iter < maxIter_
167  && yp > VSMALL
168  );
169 
170  yPlus[facei] = max(0.0, yp);
171  }
172  }
173  }
174  else
175  {
176  // Smooth Walls
177  forAll(yPlus, facei)
178  {
179  const scalar magUpara = magUp[facei];
180  const scalar Re = magUpara*y[facei]/nuw[facei];
181  const scalar kappaRe = kappa_*Re;
182 
183  scalar yp = yPlusLam_;
184  const scalar ryPlusLam = 1.0/yp;
185 
186  int iter = 0;
187  scalar yPlusLast = 0.0;
188 
189  do
190  {
191  yPlusLast = yp;
192  yp = (kappaRe + yp)/(1.0 + log(E_*yp));
193 
194  }
195  while
196  (
197  mag(ryPlusLam*(yp - yPlusLast)) > tolerance_
198  && ++iter < maxIter_
199  );
200 
201  yPlus[facei] = max(0.0, yp);
202  }
203  }
204 
205  return tyPlus;
206 }
207 
208 
209 void Foam::nutURoughWallFunctionFvPatchScalarField::writeLocalEntries
210 (
211  Ostream& os
212 ) const
213 {
215  os.writeEntry("roughnessHeight", roughnessHeight_);
216  os.writeEntry("roughnessConstant", roughnessConstant_);
217  os.writeEntry("roughnessFactor", roughnessFactor_);
218  os.writeEntryIfDifferent<label>("maxIter", 10, maxIter_);
219  os.writeEntryIfDifferent<scalar>("tolerance", 0.0001, tolerance_);
220 }
221 
222 
223 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
224 
227 (
228  const fvPatch& p,
230 )
231 :
233  roughnessHeight_(Zero),
234  roughnessConstant_(Zero),
235  roughnessFactor_(Zero),
236  maxIter_(10),
237  tolerance_(0.0001)
238 {}
239 
240 
243 (
245  const fvPatch& p,
247  const fvPatchFieldMapper& mapper
248 )
249 :
250  nutWallFunctionFvPatchScalarField(ptf, p, iF, mapper),
251  roughnessHeight_(ptf.roughnessHeight_),
252  roughnessConstant_(ptf.roughnessConstant_),
253  roughnessFactor_(ptf.roughnessFactor_),
254  maxIter_(ptf.maxIter_),
255  tolerance_(ptf.tolerance_)
256 {}
257 
258 
261 (
262  const fvPatch& p,
264  const dictionary& dict
265 )
266 :
268  roughnessHeight_(dict.get<scalar>("roughnessHeight")),
269  roughnessConstant_(dict.get<scalar>("roughnessConstant")),
270  roughnessFactor_(dict.get<scalar>("roughnessFactor")),
271  maxIter_(dict.getOrDefault<label>("maxIter", 10)),
272  tolerance_(dict.getOrDefault<scalar>("tolerance", 0.0001))
273 {}
274 
275 
278 (
280 )
281 :
283  roughnessHeight_(rwfpsf.roughnessHeight_),
284  roughnessConstant_(rwfpsf.roughnessConstant_),
285  roughnessFactor_(rwfpsf.roughnessFactor_),
286  maxIter_(rwfpsf.maxIter_),
287  tolerance_(rwfpsf.tolerance_)
288 {}
289 
290 
293 (
296 )
297 :
299  roughnessHeight_(rwfpsf.roughnessHeight_),
300  roughnessConstant_(rwfpsf.roughnessConstant_),
301  roughnessFactor_(rwfpsf.roughnessFactor_),
302  maxIter_(rwfpsf.maxIter_),
303  tolerance_(rwfpsf.tolerance_)
304 {}
305 
306 
307 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
308 
311 {
312  const label patchi = patch().index();
313 
314  const turbulenceModel& turbModel = db().lookupObject<turbulenceModel>
315  (
317  (
319  internalField().group()
320  )
321  );
322  const fvPatchVectorField& Uw = U(turbModel).boundaryField()[patchi];
323  tmp<scalarField> magUp = mag(Uw.patchInternalField() - Uw);
324 
325  return calcYPlus(magUp());
326 }
327 
328 
330 (
331  Ostream& os
332 ) const
333 {
335  writeLocalEntries(os);
336  writeEntry("value", os);
337 }
338 
339 
340 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
341 
342 namespace Foam
343 {
345  (
348  );
349 }
350 
351 // ************************************************************************* //
Foam::fvPatchField< vector >
volFields.H
Foam::fvPatchField::write
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:384
Foam::Ostream::writeEntryIfDifferent
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:248
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::constant::atomic::group
constexpr const char *const group
Group name for atomic constants.
Definition: atomicConstants.H:52
Foam::constant::universal::G
const dimensionedScalar G
Newtonian constant of gravitation.
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::nutWallFunctionFvPatchScalarField::yPlusLam_
scalar yPlusLam_
Definition: nutWallFunctionFvPatchScalarField.H:297
Foam::sin
dimensionedScalar sin(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:264
Foam::nutURoughWallFunctionFvPatchScalarField
This boundary condition provides a wall constraint on the turbulent viscosity, i.e....
Definition: nutURoughWallFunctionFvPatchScalarField.H:138
Foam::turbulenceModel::propertiesName
static const word propertiesName
Default name of the turbulence properties dictionary.
Definition: turbulenceModel.H:100
Foam::nutWallFunctionFvPatchScalarField::U
virtual const volVectorField & U(const turbulenceModel &turb) const
Definition: nutWallFunctionFvPatchScalarField.C:75
fvPatchFieldMapper.H
Foam::dictionary::get
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:107
magUp
scalar magUp
Definition: evaluateNearWall.H:10
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::compressible::turbulenceModel
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
Definition: turbulentFluidThermoModel.H:63
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::nutURoughWallFunctionFvPatchScalarField::write
virtual void write(Ostream &os) const
Write.
Definition: nutURoughWallFunctionFvPatchScalarField.C:330
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::turbulenceModel
Abstract base class for turbulence models (RAS, LES and laminar).
Definition: turbulenceModel.H:63
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::nutURoughWallFunctionFvPatchScalarField::nutURoughWallFunctionFvPatchScalarField
nutURoughWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Definition: nutURoughWallFunctionFvPatchScalarField.C:227
Foam::log
dimensionedScalar log(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:262
os
OBJstream os(runTime.globalPath()/outputName)
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::nutURoughWallFunctionFvPatchScalarField::yPlus
virtual tmp< scalarField > yPlus() const
Calculate and return the yPlus at the boundary.
Definition: nutURoughWallFunctionFvPatchScalarField.C:310
Foam::fvPatchVectorField
fvPatchField< vector > fvPatchVectorField
Definition: fvPatchFieldsFwd.H:43
U
U
Definition: pEqn.H:72
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::nutWallFunctionFvPatchScalarField::nutw
static const nutWallFunctionFvPatchScalarField & nutw(const turbulenceModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
Definition: nutWallFunctionFvPatchScalarField.C:235
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
nutURoughWallFunctionFvPatchScalarField.H
Foam::Re
scalarField Re(const UList< complex > &cf)
Extract real component.
Definition: complexField.C:159
Foam::atan
dimensionedScalar atan(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:269
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
yPlus
scalar yPlus
Definition: evaluateNearWall.H:16
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
Foam::IOobject::groupName
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::makePatchTypeField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Foam::dictionary::getOrDefault
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:148
Foam::nutWallFunctionFvPatchScalarField
The class nutWallFunction is a base class that parents the derived boundary conditions which provide ...
Definition: nutWallFunctionFvPatchScalarField.H:251
turbulenceModel.H
Foam::nutWallFunctionFvPatchScalarField::writeLocalEntries
virtual void writeLocalEntries(Ostream &) const
Write local wall function variables.
Definition: nutWallFunctionFvPatchScalarField.C:91
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
y
scalar y
Definition: LISASMDCalcMethod1.H:14
Foam::cos
dimensionedScalar cos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:265