FreeFOAM The Cross-Platform CFD Toolkit
SpalartAllmaras.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) 1991-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 
26 #include "SpalartAllmaras.H"
28 #include <finiteVolume/wallDist.H>
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace compressible
35 {
36 namespace LESModels
37 {
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 defineTypeNameAndDebug(SpalartAllmaras, 0);
42 addToRunTimeSelectionTable(LESModel, SpalartAllmaras, dictionary);
43 
44 
45 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
46 
47 void SpalartAllmaras::updateSubGridScaleFields()
48 {
49  muSgs_.internalField() = rho()*fv1()*nuTilda_.internalField();
50  muSgs_.correctBoundaryConditions();
51 
52  alphaSgs_ = muSgs_/Prt_;
53  alphaSgs_.correctBoundaryConditions();
54 }
55 
56 
57 tmp<volScalarField> SpalartAllmaras::fv1() const
58 {
59  volScalarField chi3 = pow3(rho()*nuTilda_/mu());
60  return chi3/(chi3 + pow3(Cv1_));
61 }
62 
63 
64 tmp<volScalarField> SpalartAllmaras::fv2() const
65 {
66  return 1.0/pow3(scalar(1) + rho()*nuTilda_/(mu()*Cv2_));
67 }
68 
69 
70 tmp<volScalarField> SpalartAllmaras::fv3() const
71 {
72  volScalarField chi = rho()*nuTilda_/mu();
73  volScalarField chiByCv2 = (1/Cv2_)*chi;
74 
75  return
76  (scalar(1) + chi*fv1())
77  *(1/Cv2_)
78  *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
79  /pow3(scalar(1) + chiByCv2);
80 }
81 
82 
83 tmp<volScalarField> SpalartAllmaras::fw(const volScalarField& Stilda) const
84 {
86  (
87  nuTilda_
88  /(
89  max(Stilda, dimensionedScalar("SMALL", Stilda.dimensions(), SMALL))
90  *sqr(kappa_*dTilda_)
91  ),
92  scalar(10.0)
93  );
94  r.boundaryField() == 0.0;
95 
96  volScalarField g = r + Cw2_*(pow6(r) - r);
97 
98  return g*pow((1.0 + pow6(Cw3_))/(pow6(g) + pow6(Cw3_)), 1.0/6.0);
99 }
100 
101 
102 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
103 
104 SpalartAllmaras::SpalartAllmaras
105 (
106  const volScalarField& rho,
107  const volVectorField& U,
108  const surfaceScalarField& phi,
109  const basicThermo& thermoPhysicalModel
110 )
111 :
112  LESModel(typeName, rho, U, phi, thermoPhysicalModel),
113 
114  sigmaNut_
115  (
117  (
118  "sigmaNut",
119  coeffDict_,
120  0.66666
121  )
122  ),
123  Prt_
124  (
126  (
127  "Prt",
128  coeffDict_,
129  1.0
130  )
131  ),
132 
133  Cb1_
134  (
136  (
137  "Cb1",
138  coeffDict_,
139  0.1355
140  )
141  ),
142  Cb2_
143  (
145  (
146  "Cb2",
147  coeffDict_,
148  0.622
149  )
150  ),
151  Cv1_
152  (
154  (
155  "Cv1",
156  coeffDict_,
157  7.1
158  )
159  ),
160  Cv2_
161  (
163  (
164  "Cv2",
165  coeffDict_,
166  5.0
167  )
168  ),
169  CDES_
170  (
172  (
173  "CDES",
174  coeffDict_,
175  0.65
176  )
177  ),
178  ck_
179  (
181  (
182  "ck",
183  coeffDict_,
184  0.07
185  )
186  ),
187  kappa_
188  (
190  (
191  "kappa",
192  *this,
193  0.41
194  )
195  ),
196  Cw1_(Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_),
197  Cw2_
198  (
200  (
201  "Cw2",
202  coeffDict_,
203  0.3
204  )
205  ),
206  Cw3_
207  (
209  (
210  "Cw3",
211  coeffDict_,
212  2.0
213  )
214  ),
215 
216  nuTilda_
217  (
218  IOobject
219  (
220  "nuTilda",
221  runTime_.timeName(),
222  mesh_,
225  ),
226  mesh_
227  ),
228 
229  dTilda_(min(CDES_*delta(), wallDist(mesh_).y())),
230  muSgs_
231  (
232  IOobject
233  (
234  "muSgs",
235  runTime_.timeName(),
236  mesh_,
239  ),
240  mesh_
241  ),
242 
243  alphaSgs_
244  (
245  IOobject
246  (
247  "alphaSgs",
248  runTime_.timeName(),
249  mesh_,
252  ),
253  mesh_
254  )
255 {
256  updateSubGridScaleFields();
257 
258  printCoeffs();
259 }
260 
261 
262 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
263 
264 tmp<volSymmTensorField> SpalartAllmaras::B() const
265 {
266  return ((2.0/3.0)*I)*k() - (muSgs_/rho())*dev(twoSymm(fvc::grad(U())));
267 }
268 
269 
270 tmp<volSymmTensorField> SpalartAllmaras::devRhoBeff() const
271 {
272  return -muEff()*dev(twoSymm(fvc::grad(U())));
273 }
274 
275 
276 tmp<volScalarField> SpalartAllmaras::epsilon() const
277 {
278  return 2*muEff()/rho()*magSqr(symm(fvc::grad(U())));
279 }
280 
281 
282 tmp<fvVectorMatrix> SpalartAllmaras::divDevRhoBeff(volVectorField& U) const
283 {
284  return
285  (
286  - fvm::laplacian(muEff(), U) - fvc::div(muEff()*dev2(fvc::grad(U)().T()))
287  );
288 }
289 
290 
292 {
293  const volTensorField& gradU = tgradU();
294  LESModel::correct(gradU);
295 
296  if (mesh_.changing())
297  {
298  dTilda_ = min(CDES_*delta(), wallDist(mesh_).y());
299  }
300 
301  volScalarField Stilda =
302  fv3()*::sqrt(2.0)*mag(skew(gradU)) + fv2()*nuTilda_/sqr(kappa_*dTilda_);
303 
304  solve
305  (
306  fvm::ddt(rho(), nuTilda_)
307  + fvm::div(phi(), nuTilda_)
309  (
310  (nuTilda_*rho() + mu())/sigmaNut_,
311  nuTilda_,
312  "laplacian(DnuTildaEff,nuTilda)"
313  )
314  - rho()*Cb2_/sigmaNut_*magSqr(fvc::grad(nuTilda_))
315  ==
316  rho()*Cb1_*Stilda*nuTilda_
317  - fvm::Sp(rho()*Cw1_*fw(Stilda)*nuTilda_/sqr(dTilda_), nuTilda_)
318  );
319 
320  bound(nuTilda_, dimensionedScalar("zero", nuTilda_.dimensions(), 0.0));
321  nuTilda_.correctBoundaryConditions();
322 
323  updateSubGridScaleFields();
324 }
325 
326 
327 bool SpalartAllmaras::read()
328 {
329  if (LESModel::read())
330  {
331  sigmaNut_.readIfPresent(coeffDict());
332  Prt_.readIfPresent(coeffDict());
333  Cb1_.readIfPresent(coeffDict());
334  Cb2_.readIfPresent(coeffDict());
335  Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
336  Cw2_.readIfPresent(coeffDict());
337  Cw3_.readIfPresent(coeffDict());
338  Cv1_.readIfPresent(coeffDict());
339  Cv2_.readIfPresent(coeffDict());
340  CDES_.readIfPresent(coeffDict());
341  ck_.readIfPresent(coeffDict());
342  kappa_.readIfPresent(*this);
343 
344  return true;
345  }
346  else
347  {
348  return false;
349  }
350 }
351 
352 
353 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
354 
355 } // End namespace LESModels
356 } // End namespace compressible
357 } // End namespace Foam
358 
359 // ************************ vim: set sw=4 sts=4 et: ************************ //