nutUSpaldingWallFunctionFvPatchScalarField.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) 2011-2020 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 "momentumTransportModel.H"
28 #include "fvPatchFieldMapper.H"
29 #include "volFields.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
38 
40 {
41  const label patchi = patch().index();
42 
43  const momentumTransportModel& turbModel =
44  db().lookupObject<momentumTransportModel>
45  (
47  (
48  momentumTransportModel::typeName,
49  internalField().group()
50  )
51  );
52  const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];
53  const scalarField magGradU(mag(Uw.snGrad()));
54  const tmp<scalarField> tnuw = turbModel.nu(patchi);
55  const scalarField& nuw = tnuw();
56 
57  return max
58  (
59  scalar(0),
60  sqr(calcUTau(magGradU))/(magGradU + rootVSmall) - nuw
61  );
62 }
63 
64 
66 (
67  const scalarField& magGradU
68 ) const
69 {
70  const label patchi = patch().index();
71 
72  const momentumTransportModel& turbModel =
73  db().lookupObject<momentumTransportModel>
74  (
76  (
77  momentumTransportModel::typeName,
78  internalField().group()
79  )
80  );
81  const scalarField& y = turbModel.y()[patchi];
82 
83  const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];
84  const scalarField magUp(mag(Uw.patchInternalField() - Uw));
85 
86  const tmp<scalarField> tnuw = turbModel.nu(patchi);
87  const scalarField& nuw = tnuw();
88 
89  const scalarField& nutw = *this;
90 
91  tmp<scalarField> tuTau(new scalarField(patch().size(), 0.0));
92  scalarField& uTau = tuTau.ref();
93 
94  forAll(uTau, facei)
95  {
96  scalar ut = sqrt((nutw[facei] + nuw[facei])*magGradU[facei]);
97 
98  if (ut > rootVSmall)
99  {
100  int iter = 0;
101  scalar err = great;
102 
103  do
104  {
105  scalar kUu = min(kappa_*magUp[facei]/ut, 50);
106  scalar fkUu = exp(kUu) - 1 - kUu*(1 + 0.5*kUu);
107 
108  scalar f =
109  - ut*y[facei]/nuw[facei]
110  + magUp[facei]/ut
111  + 1/E_*(fkUu - 1.0/6.0*kUu*sqr(kUu));
112 
113  scalar df =
114  y[facei]/nuw[facei]
115  + magUp[facei]/sqr(ut)
116  + 1/E_*kUu*fkUu/ut;
117 
118  scalar uTauNew = ut + f/df;
119  err = mag((ut - uTauNew)/ut);
120  ut = uTauNew;
121 
122  } while (ut > rootVSmall && err > 0.01 && ++iter < 10);
123 
124  uTau[facei] = max(0.0, ut);
125  }
126  }
127 
128  return tuTau;
129 }
130 
131 
132 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
133 
136 (
137  const fvPatch& p,
139 )
140 :
142 {}
143 
144 
147 (
149  const fvPatch& p,
151  const fvPatchFieldMapper& mapper
152 )
153 :
154  nutWallFunctionFvPatchScalarField(ptf, p, iF, mapper)
155 {}
156 
157 
160 (
161  const fvPatch& p,
163  const dictionary& dict
164 )
165 :
167 {}
168 
169 
172 (
174 )
175 :
177 {}
178 
179 
182 (
185 )
186 :
188 {}
189 
190 
191 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
192 
194 {
195  const label patchi = patch().index();
196 
197  const momentumTransportModel& turbModel =
198  db().lookupObject<momentumTransportModel>
199  (
201  (
202  momentumTransportModel::typeName,
203  internalField().group()
204  )
205  );
206  const scalarField& y = turbModel.y()[patchi];
207  const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];
208  const tmp<scalarField> tnuw = turbModel.nu(patchi);
209  const scalarField& nuw = tnuw();
210 
211  return y*calcUTau(mag(Uw.snGrad()))/nuw;
212 }
213 
214 
216 {
218  writeLocalEntries(os);
219  writeEntry(os, "value", *this);
220 }
221 
222 
223 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
224 
226 (
229 );
230 
231 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
232 
233 } // End namespace Foam
234 
235 // ************************************************************************* //
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
const char *const group
Group name for atomic constants.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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 tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:186
virtual void writeLocalEntries(Ostream &) const
Write local wall function variables.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
const nearWallDist & y() const
Return the near wall distances.
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:280
virtual tmp< scalarField > calcUTau(const scalarField &magGradU) const
Calculate the friction velocity.
Macros for easy insertion into run-time selection tables.
scalar magUp
scalar uTau
static const nutWallFunctionFvPatchScalarField & nutw(const momentumTransportModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
scalar y
virtual tmp< scalarField > nut() const
Calculate the turbulence viscosity.
This boundary condition provides a turbulent kinematic viscosity condition when using wall functions ...
dimensionedScalar exp(const dimensionedScalar &ds)
nutUSpaldingWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static word groupName(Name name, const word &group)
Foam::fvPatchFieldMapper.
virtual tmp< scalarField > yPlus() const
Calculate and return the yPlus at the boundary.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:194
labelList f(nPoints)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Abstract base class for turbulence models (RAS, LES and laminar).
const volVectorField & U() const
Access function to velocity field.
label patchi
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensioned< scalar > mag(const dimensioned< Type > &)
nutWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
A class for managing temporary objects.
Definition: PtrList.H:53
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Namespace for OpenFOAM.