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 LESModels
36 {
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 defineTypeNameAndDebug(SpalartAllmaras, 0);
41 addToRunTimeSelectionTable(LESModel, SpalartAllmaras, dictionary);
42 
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
44 
45 void SpalartAllmaras::updateSubGridScaleFields()
46 {
49 }
50 
51 
52 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
53 
55 {
56  volScalarField chi3 = pow3(nuTilda_/nu());
57  return chi3/(chi3 + pow3(Cv1_));
58 }
59 
60 
62 {
63  return 1/pow3(scalar(1) + nuTilda_/(Cv2_*nu()));
64 }
65 
66 
68 {
69  volScalarField chi = nuTilda_/nu();
70  volScalarField chiByCv2 = (1/Cv2_)*chi;
71 
72  return
73  (scalar(1) + chi*fv1())
74  *(1/Cv2_)
75  *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
76  /pow3(scalar(1) + chiByCv2);
77 }
78 
79 
81 {
82  return sqrt(2.0)*mag(skew(gradU));
83 }
84 
85 
87 (
88  const volScalarField& S,
89  const volScalarField& dTilda
90 ) const
91 {
92  return fv3()*S + fv2()*nuTilda_/sqr(kappa_*dTilda);
93 }
94 
95 
97 (
98  const volScalarField& visc,
99  const volScalarField& S,
100  const volScalarField& dTilda
101 ) const
102 {
103  return min
104  (
105  visc
106  /(
107  max
108  (
109  S,
110  dimensionedScalar("SMALL", S.dimensions(), SMALL)
111  )
112  *sqr(kappa_*dTilda)
114  (
115  "ROOTVSMALL",
116  dimensionSet(0, 2 , -1, 0, 0),
117  ROOTVSMALL
118  )
119  ),
120  scalar(10)
121  );
122 }
123 
124 
126 (
127  const volScalarField& S,
128  const volScalarField& dTilda
129 ) const
130 {
131  volScalarField r = this->r(nuTilda_, S, dTilda);
132 
133  volScalarField g = r + Cw2_*(pow6(r) - r);
134 
135  return g*pow((1 + pow6(Cw3_))/(pow6(g) + pow6(Cw3_)), 1.0/6.0);
136 }
137 
138 
140 {
141  return min(CDES_*delta(), y_);
142 }
143 
144 
145 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
146 
147 SpalartAllmaras::SpalartAllmaras
148 (
149  const volVectorField& U,
150  const surfaceScalarField& phi,
151  transportModel& transport,
152  const word& modelName
153 )
154 :
155  LESModel(modelName, U, phi, transport),
156 
157  sigmaNut_
158  (
160  (
161  "sigmaNut",
162  coeffDict_,
163  0.66666
164  )
165  ),
166  kappa_
167  (
169  (
170  "kappa",
171  coeffDict_,
172  0.41
173  )
174  ),
175  Cb1_
176  (
178  (
179  "Cb1",
180  coeffDict_,
181  0.1355
182  )
183  ),
184  Cb2_
185  (
187  (
188  "Cb2",
189  coeffDict_,
190  0.622
191  )
192  ),
193  Cv1_
194  (
196  (
197  "Cv1",
198  coeffDict_,
199  7.1
200  )
201  ),
202  Cv2_
203  (
205  (
206  "Cv2",
207  coeffDict_,
208  5.0
209  )
210  ),
211  CDES_
212  (
214  (
215  "CDES",
216  coeffDict_,
217  0.65
218  )
219  ),
220  ck_
221  (
223  (
224  "ck",
225  coeffDict_,
226  0.07
227  )
228  ),
229  Cw1_(Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_),
230  Cw2_
231  (
233  (
234  "Cw2",
235  coeffDict_,
236  0.3
237  )
238  ),
239  Cw3_
240  (
242  (
243  "Cw3",
244  coeffDict_,
245  2.0
246  )
247  ),
248 
249  y_(mesh_),
250 
251  nuTilda_
252  (
253  IOobject
254  (
255  "nuTilda",
256  runTime_.timeName(),
257  mesh_,
260  ),
261  mesh_
262  ),
263 
264  nuSgs_
265  (
266  IOobject
267  (
268  "nuSgs",
269  runTime_.timeName(),
270  mesh_,
273  ),
274  mesh_
275  )
276 {
277  updateSubGridScaleFields();
278 }
279 
280 
281 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
282 
284 {
285  LESModel::correct(gradU);
286 
287  if (mesh_.changing())
288  {
289  y_.correct();
290  y_.boundaryField() = max(y_.boundaryField(), VSMALL);
291  }
292 
293  const volScalarField S = this->S(gradU);
294  const volScalarField dTilda = this->dTilda(S);
295  const volScalarField STilda = this->STilda(S, dTilda);
296 
297  fvScalarMatrix nuTildaEqn
298  (
300  + fvm::div(phi(), nuTilda_)
302  (
303  (nuTilda_ + nu())/sigmaNut_,
304  nuTilda_,
305  "laplacian(DnuTildaEff,nuTilda)"
306  )
308  ==
309  Cb1_*STilda*nuTilda_
310  - fvm::Sp(Cw1_*fw(STilda, dTilda)*nuTilda_/sqr(dTilda), nuTilda_)
311  );
312 
313  nuTildaEqn.relax();
314  nuTildaEqn.solve();
315 
318 
319  updateSubGridScaleFields();
320 }
321 
322 
324 {
325  return sqr(nuSgs()/ck_/dTilda(S(fvc::grad(U()))));
326 }
327 
328 
330 {
331  return 2*nuEff()*magSqr(symm(fvc::grad(U())));
332 }
333 
334 
336 {
337  return ((2.0/3.0)*I)*k() - nuSgs()*twoSymm(fvc::grad(U()));
338 }
339 
340 
342 {
343  return -nuEff()*dev(twoSymm(fvc::grad(U())));
344 }
345 
346 
348 {
349  return
350  (
351  - fvm::laplacian(nuEff(), U) - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
352  );
353 }
354 
355 
357 {
358  if (LESModel::read())
359  {
361  kappa_.readIfPresent(*this);
368  Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
371 
372  return true;
373  }
374  else
375  {
376  return false;
377  }
378 }
379 
380 
381 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
382 
383 } // End namespace LESModels
384 } // End namespace incompressible
385 } // End namespace Foam
386 
387 // ************************ vim: set sw=4 sts=4 et: ************************ //