FreeFOAM The Cross-Platform CFD Toolkit
nutWallFunctionFvPatchScalarField.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) 2008-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 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
43 
45 {
46  if (!isA<wallFvPatch>(patch()))
47  {
48  FatalErrorIn("nutWallFunctionFvPatchScalarField::checkType()")
49  << "Invalid wall function specification" << nl
50  << " Patch type for patch " << patch().name()
51  << " must be wall" << nl
52  << " Current patch type is " << patch().type() << nl << endl
53  << abort(FatalError);
54  }
55 }
56 
57 
59 (
60  const scalar kappa,
61  const scalar E
62 ) const
63 {
64  scalar ypl = 11.0;
65 
66  for (int i=0; i<10; i++)
67  {
68  ypl = log(E*ypl)/kappa;
69  }
70 
71  return ypl;
72 }
73 
74 
76 {
77  const label patchI = patch().index();
78 
79  const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
80  const scalarField& y = rasModel.y()[patchI];
81  const tmp<volScalarField> tk = rasModel.k();
82  const volScalarField& k = tk();
83  const scalarField& nuw = rasModel.nu().boundaryField()[patchI];
84 
85  const scalar Cmu25 = pow(Cmu_, 0.25);
86 
87  tmp<scalarField> tnutw(new scalarField(patch().size(), 0.0));
88  scalarField& nutw = tnutw();
89 
90  forAll(nutw, faceI)
91  {
92  label faceCellI = patch().faceCells()[faceI];
93 
94  scalar yPlus = Cmu25*y[faceI]*sqrt(k[faceCellI])/nuw[faceI];
95 
96  if (yPlus > yPlusLam_)
97  {
98  nutw[faceI] = nuw[faceI]*(yPlus*kappa_/log(E_*yPlus) - 1.0);
99  }
100  }
101 
102  return tnutw;
103 }
104 
105 
107 {
108  os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
109  os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
110  os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
111 }
112 
113 
114 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
115 
117 (
118  const fvPatch& p,
120 )
121 :
122  fixedValueFvPatchScalarField(p, iF),
123  Cmu_(0.09),
124  kappa_(0.41),
125  E_(9.8),
126  yPlusLam_(calcYPlusLam(kappa_, E_))
127 {
128  checkType();
129 }
130 
131 
133 (
135  const fvPatch& p,
137  const fvPatchFieldMapper& mapper
138 )
139 :
140  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
141  Cmu_(ptf.Cmu_),
142  kappa_(ptf.kappa_),
143  E_(ptf.E_),
144  yPlusLam_(ptf.yPlusLam_)
145 {
146  checkType();
147 }
148 
149 
151 (
152  const fvPatch& p,
154  const dictionary& dict
155 )
156 :
157  fixedValueFvPatchScalarField(p, iF, dict),
158  Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
159  kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
160  E_(dict.lookupOrDefault<scalar>("E", 9.8)),
161  yPlusLam_(calcYPlusLam(kappa_, E_))
162 {
163  checkType();
164 }
165 
166 
168 (
170 )
171 :
172  fixedValueFvPatchScalarField(wfpsf),
173  Cmu_(wfpsf.Cmu_),
174  kappa_(wfpsf.kappa_),
175  E_(wfpsf.E_),
176  yPlusLam_(wfpsf.yPlusLam_)
177 {
178  checkType();
179 }
180 
181 
183 (
186 )
187 :
188  fixedValueFvPatchScalarField(wfpsf, iF),
189  Cmu_(wfpsf.Cmu_),
190  kappa_(wfpsf.kappa_),
191  E_(wfpsf.E_),
192  yPlusLam_(wfpsf.yPlusLam_)
193 {
194  checkType();
195 }
196 
197 
198 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
199 
201 {
202  operator==(calcNut());
203 
205 }
206 
207 
209 {
210  const label patchI = patch().index();
211 
212  const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
213  const scalarField& y = rasModel.y()[patchI];
214 
215  const tmp<volScalarField> tk = rasModel.k();
216  const volScalarField& k = tk();
217  const scalarField kwc = k.boundaryField()[patchI].patchInternalField();
218  const scalarField& nuw = rasModel.nu().boundaryField()[patchI];
219 
220  return pow(Cmu_, 0.25)*y*sqrt(kwc)/nuw;
221 }
222 
223 
225 {
227  writeLocalEntries(os);
228  writeEntry("value", os);
229 }
230 
231 
232 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
233 
235 
236 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
237 
238 } // End namespace RASModels
239 } // End namespace incompressible
240 } // End namespace Foam
241 
242 // ************************ vim: set sw=4 sts=4 et: ************************ //