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 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace incompressible
34 {
35 namespace RASModels
36 {
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 defineTypeNameAndDebug(SpalartAllmaras, 0);
41 addToRunTimeSelectionTable(RASModel, SpalartAllmaras, dictionary);
42 
43 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
44 
45 tmp<volScalarField> SpalartAllmaras::chi() const
46 {
47  return nuTilda_/nu();
48 }
49 
50 
51 tmp<volScalarField> SpalartAllmaras::fv1(const volScalarField& chi) const
52 {
53  volScalarField chi3 = pow3(chi);
54  return chi3/(chi3 + pow3(Cv1_));
55 }
56 
57 
58 tmp<volScalarField> SpalartAllmaras::fv2
59 (
60  const volScalarField& chi,
61  const volScalarField& fv1
62 ) const
63 {
64  return 1.0/pow3(scalar(1) + chi/Cv2_);
65 }
66 
67 
68 tmp<volScalarField> SpalartAllmaras::fv3
69 (
70  const volScalarField& chi,
71  const volScalarField& fv1
72 ) const
73 {
74  volScalarField chiByCv2 = (1/Cv2_)*chi;
75 
76  return
77  (scalar(1) + chi*fv1)
78  *(1/Cv2_)
79  *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
80  /pow3(scalar(1) + chiByCv2);
81 }
82 
83 
84 tmp<volScalarField> SpalartAllmaras::fw(const volScalarField& Stilda) const
85 {
87  (
88  nuTilda_
89  /(
90  max(Stilda, dimensionedScalar("SMALL", Stilda.dimensions(), SMALL))
91  *sqr(kappa_*d_)
92  ),
93  scalar(10.0)
94  );
95  r.boundaryField() == 0.0;
96 
97  volScalarField g = r + Cw2_*(pow6(r) - r);
98 
99  return g*pow((1.0 + pow6(Cw3_))/(pow6(g) + pow6(Cw3_)), 1.0/6.0);
100 }
101 
102 
103 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
104 
106 (
107  const volVectorField& U,
108  const surfaceScalarField& phi,
109  transportModel& lamTransportModel
110 )
111 :
112  RASModel(typeName, U, phi, lamTransportModel),
113 
114  sigmaNut_
115  (
117  (
118  "sigmaNut",
119  coeffDict_,
120  0.66666
121  )
122  ),
123  kappa_
124  (
126  (
127  "kappa",
128  coeffDict_,
129  0.41
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  Cw1_(Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_),
152  Cw2_
153  (
155  (
156  "Cw2",
157  coeffDict_,
158  0.3
159  )
160  ),
161  Cw3_
162  (
164  (
165  "Cw3",
166  coeffDict_,
167  2.0
168  )
169  ),
170  Cv1_
171  (
173  (
174  "Cv1",
175  coeffDict_,
176  7.1
177  )
178  ),
179  Cv2_
180  (
182  (
183  "Cv2",
184  coeffDict_,
185  5.0
186  )
187  ),
188 
189  nuTilda_
190  (
191  IOobject
192  (
193  "nuTilda",
194  runTime_.timeName(),
195  mesh_,
198  ),
199  mesh_
200  ),
201 
202  nut_
203  (
204  IOobject
205  (
206  "nut",
207  runTime_.timeName(),
208  mesh_,
211  ),
212  mesh_
213  ),
214 
215  d_(mesh_)
216 {
217  printCoeffs();
218 }
219 
220 
221 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
222 
224 {
225  return tmp<volScalarField>
226  (
227  new volScalarField("DnuTildaEff", nuTilda_/sigmaNut_ + nu())
228  );
229 }
230 
231 
233 {
234  return tmp<volScalarField>
235  (
236  new volScalarField
237  (
238  IOobject
239  (
240  "k",
241  runTime_.timeName(),
242  mesh_
243  ),
244  mesh_,
245  dimensionedScalar("0", dimensionSet(0, 2, -2, 0, 0), 0)
246  )
247  );
248 }
249 
250 
252 {
253  return tmp<volScalarField>
254  (
255  new volScalarField
256  (
257  IOobject
258  (
259  "epsilon",
260  runTime_.timeName(),
261  mesh_
262  ),
263  mesh_,
264  dimensionedScalar("0", dimensionSet(0, 2, -3, 0, 0), 0)
265  )
266  );
267 }
268 
269 
271 {
273  (
275  (
276  IOobject
277  (
278  "R",
279  runTime_.timeName(),
280  mesh_,
283  ),
284  ((2.0/3.0)*I)*k() - nut()*twoSymm(fvc::grad(U_))
285  )
286  );
287 }
288 
289 
291 {
293  (
295  (
296  IOobject
297  (
298  "devRhoReff",
299  runTime_.timeName(),
300  mesh_,
303  ),
305  )
306  );
307 }
308 
309 
311 {
312  volScalarField nuEff_ = nuEff();
313 
314  return
315  (
316  - fvm::laplacian(nuEff_, U)
317  - fvc::div(nuEff_*dev(fvc::grad(U)().T()))
318  );
319 }
320 
321 
323 {
324  if (RASModel::read())
325  {
326  sigmaNut_.readIfPresent(coeffDict());
327 
328  Cb1_.readIfPresent(coeffDict());
329  Cb2_.readIfPresent(coeffDict());
330  Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
331  Cw2_.readIfPresent(coeffDict());
332  Cw3_.readIfPresent(coeffDict());
333  Cv1_.readIfPresent(coeffDict());
334  Cv2_.readIfPresent(coeffDict());
335 
336  return true;
337  }
338  else
339  {
340  return false;
341  }
342 }
343 
344 
346 {
348 
349  if (!turbulence_)
350  {
351  return;
352  }
353 
354  if (mesh_.changing())
355  {
356  d_.correct();
357  }
358 
359  volScalarField chi = this->chi();
360  volScalarField fv1 = this->fv1(chi);
361 
362  volScalarField Stilda =
363  fv3(chi, fv1)*::sqrt(2.0)*mag(skew(fvc::grad(U_)))
364  + fv2(chi, fv1)*nuTilda_/sqr(kappa_*d_);
365 
366  tmp<fvScalarMatrix> nuTildaEqn
367  (
368  fvm::ddt(nuTilda_)
369  + fvm::div(phi_, nuTilda_)
370  - fvm::Sp(fvc::div(phi_), nuTilda_)
371  - fvm::laplacian(DnuTildaEff(), nuTilda_)
372  - Cb2_/sigmaNut_*magSqr(fvc::grad(nuTilda_))
373  ==
374  Cb1_*Stilda*nuTilda_
375  - fvm::Sp(Cw1_*fw(Stilda)*nuTilda_/sqr(d_), nuTilda_)
376  );
377 
378  nuTildaEqn().relax();
379  solve(nuTildaEqn);
380  bound(nuTilda_, dimensionedScalar("0", nuTilda_.dimensions(), 0.0));
381  nuTilda_.correctBoundaryConditions();
382 
383  nut_.internalField() = fv1*nuTilda_.internalField();
385 }
386 
387 
388 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
389 
390 } // End namespace RASModels
391 } // End namespace incompressible
392 } // End namespace Foam
393 
394 // ************************ vim: set sw=4 sts=4 et: ************************ //