SpalartAllmarasDES.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) 2011-2015 OpenFOAM Foundation
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 "SpalartAllmarasDES.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 namespace LESModels
33 {
34 
35 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
36 
37 template<class BasicTurbulenceModel>
39 {
40  return nuTilda_/this->nu();
41 }
42 
43 
44 template<class BasicTurbulenceModel>
46 (
47  const volScalarField& chi
48 ) const
49 {
50  const volScalarField chi3("chi3", pow3(chi));
51  return chi3/(chi3 + pow3(Cv1_));
52 }
53 
54 
55 template<class BasicTurbulenceModel>
57 (
58  const volScalarField& chi,
59  const volScalarField& fv1
60 ) const
61 {
62  return 1.0 - chi/(1.0 + chi*fv1);
63 }
64 
65 
66 template<class BasicTurbulenceModel>
68 (
69  const volTensorField& gradU
70 ) const
71 {
72  return sqrt(2.0)*mag(symm(gradU));
73 }
74 
75 
76 template<class BasicTurbulenceModel>
78 (
79  const volTensorField& gradU
80 ) const
81 {
82  return sqrt(2.0)*mag(skew(gradU));
83 }
84 
85 
86 template<class BasicTurbulenceModel>
88 (
89  const volScalarField& chi,
90  const volScalarField& fv1,
91  const volScalarField& Omega,
92  const volScalarField& dTilda
93 ) const
94 {
95  return
96  (
97  max
98  (
99  Omega
100  + fv2(chi, fv1)*nuTilda_/sqr(kappa_*dTilda),
101  Cs_*Omega
102  )
103  );
104 }
105 
106 
107 template<class BasicTurbulenceModel>
109 (
110  const volScalarField& nur,
111  const volScalarField& Omega,
112  const volScalarField& dTilda
113 ) const
114 {
116  (
117  min
118  (
119  nur
120  /(
121  max
122  (
123  Omega,
124  dimensionedScalar("SMALL", Omega.dimensions(), SMALL)
125  )
126  *sqr(kappa_*dTilda)
127  ),
128  scalar(10)
129  )
130  );
131  tr().boundaryField() == 0.0;
132 
133  return tr;
134 }
135 
136 
137 template<class BasicTurbulenceModel>
139 (
140  const volScalarField& Omega,
141  const volScalarField& dTilda
142 ) const
143 {
144  const volScalarField r(this->r(nuTilda_, Omega, dTilda));
145  const volScalarField g(r + Cw2_*(pow6(r) - r));
146 
147  return g*pow((1 + pow6(Cw3_))/(pow6(g) + pow6(Cw3_)), 1.0/6.0);
148 }
149 
150 
151 template<class BasicTurbulenceModel>
153 (
154  const volScalarField& chi,
155  const volScalarField& fv1,
156  const volTensorField& gradU
157 ) const
158 {
159  tmp<volScalarField> tdTilda(CDES_*this->delta());
160  min(tdTilda().dimensionedInternalField(), tdTilda(), y_);
161  return tdTilda;
162 }
163 
164 
165 template<class BasicTurbulenceModel>
167 (
168  const volScalarField& fv1
169 )
170 {
171  this->nut_ = nuTilda_*fv1;
172  this->nut_.correctBoundaryConditions();
173 
174  BasicTurbulenceModel::correctNut();
175 }
176 
177 
178 template<class BasicTurbulenceModel>
180 {
181  correctNut(fv1(this->chi()));
182 }
183 
184 
185 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
186 
187 template<class BasicTurbulenceModel>
189 (
190  const alphaField& alpha,
191  const rhoField& rho,
192  const volVectorField& U,
193  const surfaceScalarField& alphaRhoPhi,
194  const surfaceScalarField& phi,
195  const transportModel& transport,
196  const word& propertiesName,
197  const word& type
198 )
199 :
201  (
202  type,
203  alpha,
204  rho,
205  U,
206  alphaRhoPhi,
207  phi,
208  transport,
209  propertiesName
210  ),
211 
212  sigmaNut_
213  (
215  (
216  "sigmaNut",
217  this->coeffDict_,
218  0.66666
219  )
220  ),
221  kappa_
222  (
224  (
225  "kappa",
226  this->coeffDict_,
227  0.41
228  )
229  ),
230  Cb1_
231  (
233  (
234  "Cb1",
235  this->coeffDict_,
236  0.1355
237  )
238  ),
239  Cb2_
240  (
242  (
243  "Cb2",
244  this->coeffDict_,
245  0.622
246  )
247  ),
248  Cw1_(Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_),
249  Cw2_
250  (
252  (
253  "Cw2",
254  this->coeffDict_,
255  0.3
256  )
257  ),
258  Cw3_
259  (
261  (
262  "Cw3",
263  this->coeffDict_,
264  2.0
265  )
266  ),
267  Cv1_
268  (
270  (
271  "Cv1",
272  this->coeffDict_,
273  7.1
274  )
275  ),
276  Cs_
277  (
279  (
280  "Cs",
281  this->coeffDict_,
282  0.3
283  )
284  ),
285  CDES_
286  (
288  (
289  "CDES",
290  this->coeffDict_,
291  0.65
292  )
293  ),
294  ck_
295  (
297  (
298  "ck",
299  this->coeffDict_,
300  0.07
301  )
302  ),
303 
304  nuTilda_
305  (
306  IOobject
307  (
308  "nuTilda",
309  this->runTime_.timeName(),
310  this->mesh_,
313  ),
314  this->mesh_
315  ),
316 
317  y_(wallDist::New(this->mesh_).y())
318 {
319  if (type == typeName)
320  {
321  correctNut();
322  this->printCoeffs(type);
323  }
324 }
325 
326 
327 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
328 
329 template<class BasicTurbulenceModel>
331 {
333  {
334  sigmaNut_.readIfPresent(this->coeffDict());
335  kappa_.readIfPresent(*this);
336 
337  Cb1_.readIfPresent(this->coeffDict());
338  Cb2_.readIfPresent(this->coeffDict());
339  Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
340  Cw2_.readIfPresent(this->coeffDict());
341  Cw3_.readIfPresent(this->coeffDict());
342  Cv1_.readIfPresent(this->coeffDict());
343  Cs_.readIfPresent(this->coeffDict());
344 
345  CDES_.readIfPresent(this->coeffDict());
346  ck_.readIfPresent(this->coeffDict());
347 
348  return true;
349  }
350  else
351  {
352  return false;
353  }
354 }
355 
356 
357 template<class BasicTurbulenceModel>
359 DnuTildaEff() const
360 {
361  return tmp<volScalarField>
362  (
363  new volScalarField("DnuTildaEff", (nuTilda_ + this->nu())/sigmaNut_)
364  );
365 }
366 
367 
368 template<class BasicTurbulenceModel>
370 {
371  const volScalarField chi(this->chi());
372  const volScalarField fv1(this->fv1(chi));
373  return sqr(this->nut()/ck_/dTilda(chi, fv1, fvc::grad(this->U_)));
374 }
375 
376 
377 template<class BasicTurbulenceModel>
379 {
380  const volScalarField chi(this->chi());
381  const volScalarField fv1(this->fv1(chi));
382 
383  tmp<volScalarField> tLESRegion
384  (
385  new volScalarField
386  (
387  IOobject
388  (
389  "DES::LESRegion",
390  this->mesh_.time().timeName(),
391  this->mesh_,
394  ),
395  neg(dTilda(chi, fv1, fvc::grad(this->U_)) - y_)
396  )
397  );
398 
399  return tLESRegion;
400 }
401 
402 
403 template<class BasicTurbulenceModel>
405 {
406  if (!this->turbulence_)
407  {
408  return;
409  }
410 
411  // Local references
412  const alphaField& alpha = this->alpha_;
413  const rhoField& rho = this->rho_;
414  const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_;
415  const volVectorField& U = this->U_;
416 
418 
419  const volScalarField chi(this->chi());
420  const volScalarField fv1(this->fv1(chi));
421 
422  tmp<volTensorField> tgradU = fvc::grad(U);
423  const volScalarField Omega(this->Omega(tgradU()));
424  const volScalarField dTilda(this->dTilda(chi, fv1, tgradU()));
425  const volScalarField Stilda(this->Stilda(chi, fv1, Omega, dTilda));
426 
427  tmp<fvScalarMatrix> nuTildaEqn
428  (
429  fvm::ddt(alpha, rho, nuTilda_)
430  + fvm::div(alphaRhoPhi, nuTilda_)
431  - fvm::laplacian(alpha*rho*DnuTildaEff(), nuTilda_)
432  - Cb2_/sigmaNut_*alpha*rho*magSqr(fvc::grad(nuTilda_))
433  ==
434  Cb1_*alpha*rho*Stilda*nuTilda_
435  - fvm::Sp
436  (
437  Cw1_*alpha*rho*fw(Stilda, dTilda)*nuTilda_/sqr(dTilda),
438  nuTilda_
439  )
440  );
441 
442  nuTildaEqn().relax();
443  solve(nuTildaEqn);
444  bound(nuTilda_, dimensionedScalar("zero", nuTilda_.dimensions(), 0.0));
445  nuTilda_.correctBoundaryConditions();
446 
447  correctNut(fv1);
448 }
449 
450 
451 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
452 
453 } // End namespace LESModels
454 } // End namespace Foam
455 
456 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
dimensionedScalar pow3(const dimensionedScalar &ds)
tmp< volScalarField > S(const volTensorField &gradU) const
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar neg(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
dimensionedTensor skew(const dimensionedTensor &dt)
tmp< volScalarField > fv1(const volScalarField &chi) const
dimensioned< scalar > mag(const dimensioned< Type > &)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
static const wallDist & New(const fvMesh &mesh)
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:33
A class for handling words, derived from string.
Definition: word.H:59
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
virtual void correct()
Correct nuTilda and related properties.
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField > chi() const
tmp< volScalarField > Stilda(const volScalarField &chi, const volScalarField &fv1, const volScalarField &Omega, const volScalarField &dTilda) const
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
Namespace for OpenFOAM.
tmp< volScalarField > r(const volScalarField &nur, const volScalarField &Omega, const volScalarField &dTilda) const
dimensionedScalar pow6(const dimensionedScalar &ds)
BasicTurbulenceModel::alphaField alphaField
Generic dimensioned Type class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
scalar nut
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
BasicTurbulenceModel::rhoField rhoField
virtual tmp< volScalarField > dTilda(const volScalarField &chi, const volScalarField &fv1, const volTensorField &gradU) const
Length scale.
bool readIfPresent(const dictionary &)
Update the value of dimensioned<Type> if found in the dictionary.
scalar delta
tmp< volScalarField > fw(const volScalarField &Omega, const volScalarField &dTilda) const
const dimensionSet & dimensions() const
Return dimensions.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
surfaceScalarField & phi
rDeltaT dimensionedInternalField()
tmp< volScalarField > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
BasicTurbulenceModel::transportModel transportModel
solverPerformance solve(fvMatrix< Type > &, const dictionary &)
Solve returning the solution statistics given convergence tolerance.
const dimensionedVector & g
virtual tmp< volScalarField > LESRegion() const
Return the LES field indicator.
scalar y
void correctBoundaryConditions()
Correct boundary field.
volScalarField & nu
virtual bool read()
Read LESProperties dictionary.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
Eddy viscosity LES SGS model base class.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< volScalarField > fv2(const volScalarField &chi, const volScalarField &fv1) const
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
SpalartAllmarasDES DES turbulence model for incompressible and compressible flows.
tmp< volScalarField > Omega(const volTensorField &gradU) const
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
U
Definition: pEqn.H:82