FreeFOAM The Cross-Platform CFD Toolkit
epsilonWallFunctionFvPatchScalarField.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 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
43 
44 void epsilonWallFunctionFvPatchScalarField::checkType()
45 {
46  if (!isA<wallFvPatch>(patch()))
47  {
48  FatalErrorIn("epsilonWallFunctionFvPatchScalarField::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 
58 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
59 
61 (
62  const fvPatch& p,
64 )
65 :
67  UName_("U"),
68  kName_("k"),
69  GName_("RASModel::G"),
70  nuName_("nu"),
71  nutName_("nut"),
72  Cmu_(0.09),
73  kappa_(0.41),
74  E_(9.8)
75 {
76  checkType();
77 }
78 
79 
81 (
83  const fvPatch& p,
85  const fvPatchFieldMapper& mapper
86 )
87 :
88  fixedInternalValueFvPatchField<scalar>(ptf, p, iF, mapper),
89  UName_(ptf.UName_),
90  kName_(ptf.kName_),
91  GName_(ptf.GName_),
92  nuName_(ptf.nuName_),
93  nutName_(ptf.nutName_),
94  Cmu_(ptf.Cmu_),
95  kappa_(ptf.kappa_),
96  E_(ptf.E_)
97 {
98  checkType();
99 }
100 
101 
103 (
104  const fvPatch& p,
106  const dictionary& dict
107 )
108 :
110  UName_(dict.lookupOrDefault<word>("U", "U")),
111  kName_(dict.lookupOrDefault<word>("k", "k")),
112  GName_(dict.lookupOrDefault<word>("G", "RASModel::G")),
113  nuName_(dict.lookupOrDefault<word>("nu", "nu")),
114  nutName_(dict.lookupOrDefault<word>("nut", "nut")),
115  Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
116  kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
117  E_(dict.lookupOrDefault<scalar>("E", 9.8))
118 {
119  checkType();
120 }
121 
122 
124 (
126 )
127 :
129  UName_(ewfpsf.UName_),
130  kName_(ewfpsf.kName_),
131  GName_(ewfpsf.GName_),
132  nuName_(ewfpsf.nuName_),
133  nutName_(ewfpsf.nutName_),
134  Cmu_(ewfpsf.Cmu_),
135  kappa_(ewfpsf.kappa_),
136  E_(ewfpsf.E_)
137 {
138  checkType();
139 }
140 
141 
143 (
146 )
147 :
149  UName_(ewfpsf.UName_),
150  kName_(ewfpsf.kName_),
151  GName_(ewfpsf.GName_),
152  nuName_(ewfpsf.nuName_),
153  nutName_(ewfpsf.nutName_),
154  Cmu_(ewfpsf.Cmu_),
155  kappa_(ewfpsf.kappa_),
156  E_(ewfpsf.E_)
157 {
158  checkType();
159 }
160 
161 
162 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
163 
165 {
166  if (updated())
167  {
168  return;
169  }
170 
171  const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
172  const scalarField& y = rasModel.y()[patch().index()];
173 
174  const scalar Cmu25 = pow(Cmu_, 0.25);
175  const scalar Cmu75 = pow(Cmu_, 0.75);
176 
177  volScalarField& G = const_cast<volScalarField&>
178  (db().lookupObject<volScalarField>(GName_));
179 
180  volScalarField& epsilon = const_cast<volScalarField&>
182 
183  const volScalarField& k = db().lookupObject<volScalarField>(kName_);
184 
185  const scalarField& nuw =
186  patch().lookupPatchField<volScalarField, scalar>(nuName_);
187 
188  const scalarField& nutw =
189  patch().lookupPatchField<volScalarField, scalar>(nutName_);
190 
191  const fvPatchVectorField& Uw =
193 
194  const scalarField magGradUw = mag(Uw.snGrad());
195 
196  // Set epsilon and G
197  forAll(nutw, faceI)
198  {
199  label faceCellI = patch().faceCells()[faceI];
200 
201  epsilon[faceCellI] = Cmu75*pow(k[faceCellI], 1.5)/(kappa_*y[faceI]);
202 
203  G[faceCellI] =
204  (nutw[faceI] + nuw[faceI])
205  *magGradUw[faceI]
206  *Cmu25*sqrt(k[faceCellI])
207  /(kappa_*y[faceI]);
208  }
209 
210  // TODO: perform averaging for cells sharing more than one boundary face
211 
213 }
214 
215 
217 (
218  const Pstream::commsTypes commsType
219 )
220 {
222 }
223 
224 
226 {
228  writeEntryIfDifferent<word>(os, "U", "U", UName_);
229  writeEntryIfDifferent<word>(os, "k", "k", kName_);
230  writeEntryIfDifferent<word>(os, "G", "RASModel::G", GName_);
231  writeEntryIfDifferent<word>(os, "nu", "nu", nuName_);
232  writeEntryIfDifferent<word>(os, "nut", "nut", nutName_);
233  os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
234  os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
235  os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
236  writeEntry("value", os);
237 }
238 
239 
240 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
241 
243 (
246 );
247 
248 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
249 
250 } // End namespace RASModels
251 } // End namespace incompressible
252 } // End namespace Foam
253 
254 // ************************ vim: set sw=4 sts=4 et: ************************ //