nutUWallFunctionFvPatchScalarField.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-2021 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 
40 {
41  const label patchi = patch().index();
42 
43  const turbulenceModel& turbModel = db().lookupObject<turbulenceModel>
44  (
46  (
48  internalField().group()
49  )
50  );
51  const fvPatchVectorField& Uw = U(turbModel).boundaryField()[patchi];
52  const scalarField magUp(mag(Uw.patchInternalField() - Uw));
53  const tmp<scalarField> tnuw = turbModel.nu(patchi);
54  const scalarField& nuw = tnuw();
55 
57  const scalarField& yPlus = tyPlus();
58 
59  tmp<scalarField> tnutw(new scalarField(patch().size(), Zero));
60  scalarField& nutw = tnutw.ref();
61 
62  forAll(yPlus, facei)
63  {
64  // Viscous sublayer contribution
65  const scalar nutVis = 0.0;
66 
67  // Inertial sublayer contribution
68  const scalar nutLog =
69  nuw[facei]
70  *(yPlus[facei]*kappa_/log(max(E_*yPlus[facei], 1 + 1e-4)) - 1.0);
71 
72  nutw[facei] = blend(nutVis, nutLog, yPlus[facei]);
73  }
74 
75  return tnutw;
76 }
77 
78 
81 (
82  const scalarField& magUp
83 ) const
84 {
85  const label patchi = patch().index();
86 
87  const turbulenceModel& turbModel = db().lookupObject<turbulenceModel>
88  (
90  (
92  internalField().group()
93  )
94  );
95  const scalarField& y = turbModel.y()[patchi];
96  const tmp<scalarField> tnuw = turbModel.nu(patchi);
97  const scalarField& nuw = tnuw();
98 
99  tmp<scalarField> tyPlus(new scalarField(patch().size(), Zero));
100  scalarField& yPlus = tyPlus.ref();
101 
102  forAll(yPlus, facei)
103  {
104  const scalar kappaRe = kappa_*magUp[facei]*y[facei]/nuw[facei];
105 
106  scalar yp = yPlusLam_;
107  const scalar ryPlusLam = 1.0/yp;
108 
109  int iter = 0;
110  scalar yPlusLast = 0.0;
111 
112  do
113  {
114  yPlusLast = yp;
115  yp = (kappaRe + yp)/(1.0 + log(E_*yp));
116 
117  } while (mag(ryPlusLam*(yp - yPlusLast)) > 0.01 && ++iter < 10 );
118 
119  yPlus[facei] = max(0.0, yp);
120  }
121 
122  return tyPlus;
123 }
124 
125 
126 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
127 
129 (
130  const fvPatch& p,
132 )
133 :
135 {}
136 
137 
139 (
141  const fvPatch& p,
143  const fvPatchFieldMapper& mapper
144 )
145 :
146  nutWallFunctionFvPatchScalarField(ptf, p, iF, mapper)
147 {}
148 
149 
151 (
152  const fvPatch& p,
154  const dictionary& dict
155 )
156 :
158 {}
159 
160 
162 (
164 )
165 :
167 {}
168 
169 
171 (
172  const nutUWallFunctionFvPatchScalarField& sawfpsf,
174 )
175 :
177 {}
178 
179 
180 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
181 
184 {
185  const label patchi = patch().index();
186 
187  const auto& turbModel = db().lookupObject<turbulenceModel>
188  (
190  (
192  internalField().group()
193  )
194  );
195 
196  const scalarField& y = turbModel.y()[patchi];
197 
198  tmp<scalarField> tnuw = turbModel.nu(patchi);
199  const scalarField& nuw = tnuw();
200 
201  tmp<scalarField> tnuEff = turbModel.nuEff(patchi);
202  const scalarField& nuEff = tnuEff();
203 
204  const fvPatchVectorField& Uw = U(turbModel).boundaryField()[patchi];
205  const scalarField magUp(mag(Uw.patchInternalField() - Uw));
206  const scalarField magGradUw(mag(Uw.snGrad()));
207 
208  tmp<scalarField> tyPlus = calcYPlus(magUp);
209  scalarField& yPlus = tyPlus.ref();
210 
211  forAll(yPlus, facei)
212  {
213  if (yPlusLam_ > yPlus[facei])
214  {
215  // viscous sublayer
216  yPlus[facei] =
217  y[facei]*sqrt(nuEff[facei]*magGradUw[facei])/nuw[facei];
218  }
219  }
220 
221  return tyPlus;
222 }
223 
224 
226 (
227  Ostream& os
228 ) const
229 {
231  writeLocalEntries(os);
232  writeEntry("value", os);
233 }
234 
235 
236 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
237 
238 namespace Foam
239 {
241  (
244  );
245 }
246 
247 
248 // ************************************************************************* //
Foam::nutUWallFunctionFvPatchScalarField::yPlus
virtual tmp< scalarField > yPlus() const
Calculate and return the yPlus at the boundary.
Definition: nutUWallFunctionFvPatchScalarField.C:183
Foam::fvPatchField< vector >
volFields.H
Foam::fvPatchField::write
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:384
Foam::fvPatchField::snGrad
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:225
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
Foam::turbulenceModel::nu
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
p
volScalarField & p
Definition: createFieldRefs.H:8
nutUWallFunctionFvPatchScalarField.H
Foam::constant::atomic::group
constexpr const char *const group
Group name for atomic constants.
Definition: atomicConstants.H:52
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::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::turbulenceModel::y
const nearWallDist & y() const
Return the near wall distances.
Definition: turbulenceModel.H:159
magUp
scalar magUp
Definition: evaluateNearWall.H:10
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::nutWallFunctionFvPatchScalarField::kappa_
scalar kappa_
von Kármán constant
Definition: nutWallFunctionFvPatchScalarField.H:290
Foam::nutUWallFunctionFvPatchScalarField::nutUWallFunctionFvPatchScalarField
nutUWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Definition: nutUWallFunctionFvPatchScalarField.C:129
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
Foam::Field< scalar >
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::fvPatchField::patchInternalField
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:233
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
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::nutWallFunctionFvPatchScalarField::E_
scalar E_
Wall roughness parameter.
Definition: nutWallFunctionFvPatchScalarField.H:293
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::nutUWallFunctionFvPatchScalarField::calcYPlus
virtual tmp< scalarField > calcYPlus(const scalarField &magUp) const
Calculate yPlus.
Definition: nutUWallFunctionFvPatchScalarField.C:81
U
U
Definition: pEqn.H:72
Foam::nutUWallFunctionFvPatchScalarField
This boundary condition provides a wall constraint on the turbulent viscosity, i.e....
Definition: nutUWallFunctionFvPatchScalarField.H:98
Foam::nutWallFunctionFvPatchScalarField::blend
scalar blend(const scalar nutVis, const scalar nutLog, const scalar yPlus) const
Return the blended nut according to the chosen blending treatment.
Definition: nutWallFunctionFvPatchScalarField.C:267
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::nutUWallFunctionFvPatchScalarField::calcNut
virtual tmp< scalarField > calcNut() const
Calculate the turbulent viscosity.
Definition: nutUWallFunctionFvPatchScalarField.C:39
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::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
Foam::nutUWallFunctionFvPatchScalarField::write
virtual void write(Ostream &os) const
Write.
Definition: nutUWallFunctionFvPatchScalarField.C:226
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::nutWallFunctionFvPatchScalarField
The class nutWallFunction is a base class that parents the derived boundary conditions which provide ...
Definition: nutWallFunctionFvPatchScalarField.H:251
turbulenceModel.H
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