FreeFOAM The Cross-Platform CFD Toolkit
kappatJayatillekeWallFunctionFvPatchScalarField.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) 2010-2010 OpenCFD Ltd.
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 
29 #include <finiteVolume/volFields.H>
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace incompressible
38 {
39 namespace RASModels
40 {
41 
42 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
43 
46 
47 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
48 
50 {
51  if (!isA<wallFvPatch>(patch()))
52  {
54  (
55  "kappatJayatillekeWallFunctionFvPatchScalarField::checkType()"
56  ) << "Invalid wall function specification" << nl
57  << " Patch type for patch " << patch().name()
58  << " must be wall" << nl
59  << " Current patch type is " << patch().type() << nl << endl
60  << abort(FatalError);
61  }
62 }
63 
64 
66 (
67  const scalar Prat
68 ) const
69 {
70  return 9.24*(pow(Prat, 0.75) - 1.0)*(1.0 + 0.28*exp(-0.007*Prat));
71 }
72 
73 
75 (
76  const scalar P,
77  const scalar Prat
78 ) const
79 {
80  scalar ypt = 11.0;
81 
82  for (int i=0; i<maxIters_; i++)
83  {
84  scalar f = ypt - (log(E_*ypt)/kappa_ + P)/Prat;
85  scalar df = 1.0 - 1.0/(ypt*kappa_*Prat);
86  scalar yptNew = ypt - f/df;
87 
88  if (yptNew < VSMALL)
89  {
90  return 0;
91  }
92  else if (mag(yptNew - ypt) < tolerance_)
93  {
94  return yptNew;
95  }
96  else
97  {
98  ypt = yptNew;
99  }
100  }
101 
102  return ypt;
103 }
104 
105 
106 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
107 
110 (
111  const fvPatch& p,
113 )
114 :
115  fixedValueFvPatchScalarField(p, iF),
116  Prt_(0.85),
117  Cmu_(0.09),
118  kappa_(0.41),
119  E_(9.8)
120 {
121  checkType();
122 }
123 
124 
127 (
129  const fvPatch& p,
131  const fvPatchFieldMapper& mapper
132 )
133 :
134  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
135  Prt_(ptf.Prt_),
136  Cmu_(ptf.Cmu_),
137  kappa_(ptf.kappa_),
138  E_(ptf.E_)
139 {
140  checkType();
141 }
142 
143 
146 (
147  const fvPatch& p,
149  const dictionary& dict
150 )
151 :
152  fixedValueFvPatchScalarField(p, iF, dict),
153  Prt_(readScalar(dict.lookup("Prt"))), // force read to avoid ambiguity
154  Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
155  kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
156  E_(dict.lookupOrDefault<scalar>("E", 9.8))
157 {
158  checkType();
159 }
160 
161 
164 (
166 )
167 :
168  fixedValueFvPatchScalarField(wfpsf),
169  Prt_(wfpsf.Prt_),
170  Cmu_(wfpsf.Cmu_),
171  kappa_(wfpsf.kappa_),
172  E_(wfpsf.E_)
173 {
174  checkType();
175 }
176 
177 
180 (
183 )
184 :
185  fixedValueFvPatchScalarField(wfpsf, iF),
186  Prt_(wfpsf.Prt_),
187  Cmu_(wfpsf.Cmu_),
188  kappa_(wfpsf.kappa_),
189  E_(wfpsf.E_)
190 {
191  checkType();
192 }
193 
194 
195 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
196 
198 {
199  if (updated())
200  {
201  return;
202  }
203 
204  const label patchI = patch().index();
205 
206  // Retrieve turbulence properties from model
207  const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
208  const scalar Cmu25 = pow(Cmu_, 0.25);
209  const scalarField& y = rasModel.y()[patchI];
210  const scalarField& nuw = rasModel.nu().boundaryField()[patchI];
211  const tmp<volScalarField> tk = rasModel.k();
212  const volScalarField& k = tk();
213 
214  // Molecular Prandtl number
215  const scalar
216  Pr(dimensionedScalar(rasModel.transport().lookup("Pr")).value());
217 
218  // Populate boundary values
219  scalarField& kappatw = *this;
220  forAll(kappatw, faceI)
221  {
222  label faceCellI = patch().faceCells()[faceI];
223 
224  // y+
225  scalar yPlus = Cmu25*sqrt(k[faceCellI])*y[faceI]/nuw[faceI];
226 
227  // Molecular-to-turbulent Prandtl number ratio
228  scalar Prat = Pr/Prt_;
229 
230  // Thermal sublayer thickness
231  scalar P = Psmooth(Prat);
232  scalar yPlusTherm = this->yPlusTherm(P, Prat);
233 
234  // Update turbulent thermal conductivity
235  if (yPlus > yPlusTherm)
236  {
237  scalar nu = nuw[faceI];
238  scalar kt = nu*(yPlus/(Prt_*(log(E_*yPlus)/kappa_ + P)) - 1/Pr);
239  kappatw[faceI] = max(0.0, kt);
240  }
241  else
242  {
243  kappatw[faceI] = 0.0;
244  }
245  }
246 
248 }
249 
250 
252 {
254  os.writeKeyword("Prt") << Prt_ << token::END_STATEMENT << nl;
255  os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
256  os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
257  os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
258  writeEntry("value", os);
259 }
260 
261 
262 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
263 
265 
266 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
267 
268 } // End namespace RASModels
269 } // End namespace incompressible
270 } // End namespace Foam
271 
272 // ************************************************************************* //