SpalartAllmarasDES.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 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 #include "fvOptions.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace LESModels
34 {
35 
36 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
37 
38 template<class BasicTurbulenceModel>
40 {
41  return nuTilda_/this->nu();
42 }
43 
44 
45 template<class BasicTurbulenceModel>
47 (
48  const volScalarField& chi
49 ) const
50 {
51  const volScalarField chi3("chi3", pow3(chi));
52  return chi3/(chi3 + pow3(Cv1_));
53 }
54 
55 
56 template<class BasicTurbulenceModel>
58 (
59  const volScalarField& chi,
60  const volScalarField& fv1
61 ) const
62 {
63  return 1.0 - chi/(1.0 + chi*fv1);
64 }
65 
66 
67 template<class BasicTurbulenceModel>
69 (
70  const volTensorField& gradU
71 ) const
72 {
73  return sqrt(2.0)*mag(symm(gradU));
74 }
75 
76 
77 template<class BasicTurbulenceModel>
79 (
80  const volTensorField& gradU
81 ) const
82 {
83  return sqrt(2.0)*mag(skew(gradU));
84 }
85 
86 
87 template<class BasicTurbulenceModel>
89 (
90  const volScalarField& chi,
91  const volScalarField& fv1,
92  const volScalarField& Omega,
93  const volScalarField& dTilda
94 ) const
95 {
96  return
97  (
98  max
99  (
100  Omega
101  + fv2(chi, fv1)*nuTilda_/sqr(kappa_*dTilda),
102  Cs_*Omega
103  )
104  );
105 }
106 
107 
108 template<class BasicTurbulenceModel>
110 (
111  const volScalarField& nur,
112  const volScalarField& Omega,
113  const volScalarField& dTilda
114 ) const
115 {
117  (
118  min
119  (
120  nur
121  /(
122  max
123  (
124  Omega,
125  dimensionedScalar(Omega.dimensions(), small)
126  )
127  *sqr(kappa_*dTilda)
128  ),
129  scalar(10)
130  )
131  );
132  tr.ref().boundaryFieldRef() == 0.0;
133 
134  return tr;
135 }
136 
137 
138 template<class BasicTurbulenceModel>
140 (
141  const volScalarField& Omega,
142  const volScalarField& dTilda
143 ) const
144 {
145  const volScalarField r(this->r(nuTilda_, Omega, dTilda));
146  const volScalarField g(r + Cw2_*(pow6(r) - r));
147 
148  return g*pow((1 + pow6(Cw3_))/(pow6(g) + pow6(Cw3_)), 1.0/6.0);
149 }
150 
151 
152 template<class BasicTurbulenceModel>
154 (
155  const volScalarField& chi,
156  const volScalarField& fv1,
157  const volTensorField& gradU
158 ) const
159 {
160  tmp<volScalarField> tdTilda(CDES_*this->delta());
161  min(tdTilda.ref().ref(), tdTilda(), y_);
162  return tdTilda;
163 }
164 
165 
166 template<class BasicTurbulenceModel>
168 (
169  const volScalarField& fv1
170 )
171 {
172  this->nut_ = nuTilda_*fv1;
173  this->nut_.correctBoundaryConditions();
174  fv::options::New(this->mesh_).correct(this->nut_);
175 
176  BasicTurbulenceModel::correctNut();
177 }
178 
179 
180 template<class BasicTurbulenceModel>
182 {
183  correctNut(fv1(this->chi()));
184 }
185 
186 
187 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
188 
189 template<class BasicTurbulenceModel>
191 (
192  const alphaField& alpha,
193  const rhoField& rho,
194  const volVectorField& U,
195  const surfaceScalarField& alphaRhoPhi,
196  const surfaceScalarField& phi,
197  const transportModel& transport,
198  const word& propertiesName,
199  const word& type
200 )
201 :
203  (
204  type,
205  alpha,
206  rho,
207  U,
208  alphaRhoPhi,
209  phi,
210  transport,
211  propertiesName
212  ),
213 
214  sigmaNut_
215  (
217  (
218  "sigmaNut",
219  this->coeffDict_,
220  0.66666
221  )
222  ),
223  kappa_
224  (
226  (
227  "kappa",
228  this->coeffDict_,
229  0.41
230  )
231  ),
232  Cb1_
233  (
235  (
236  "Cb1",
237  this->coeffDict_,
238  0.1355
239  )
240  ),
241  Cb2_
242  (
244  (
245  "Cb2",
246  this->coeffDict_,
247  0.622
248  )
249  ),
250  Cw1_(Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_),
251  Cw2_
252  (
254  (
255  "Cw2",
256  this->coeffDict_,
257  0.3
258  )
259  ),
260  Cw3_
261  (
263  (
264  "Cw3",
265  this->coeffDict_,
266  2.0
267  )
268  ),
269  Cv1_
270  (
272  (
273  "Cv1",
274  this->coeffDict_,
275  7.1
276  )
277  ),
278  Cs_
279  (
281  (
282  "Cs",
283  this->coeffDict_,
284  0.3
285  )
286  ),
287  CDES_
288  (
290  (
291  "CDES",
292  this->coeffDict_,
293  0.65
294  )
295  ),
296  ck_
297  (
299  (
300  "ck",
301  this->coeffDict_,
302  0.07
303  )
304  ),
305 
306  nuTilda_
307  (
308  IOobject
309  (
310  "nuTilda",
311  this->runTime_.timeName(),
312  this->mesh_,
315  ),
316  this->mesh_
317  ),
318 
319  y_(wallDist::New(this->mesh_).y())
320 {
321  if (type == typeName)
322  {
323  this->printCoeffs(type);
324  }
325 }
326 
327 
328 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
329 
330 template<class BasicTurbulenceModel>
332 {
334  {
335  sigmaNut_.readIfPresent(this->coeffDict());
336  kappa_.readIfPresent(*this);
337 
338  Cb1_.readIfPresent(this->coeffDict());
339  Cb2_.readIfPresent(this->coeffDict());
340  Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
341  Cw2_.readIfPresent(this->coeffDict());
342  Cw3_.readIfPresent(this->coeffDict());
343  Cv1_.readIfPresent(this->coeffDict());
344  Cs_.readIfPresent(this->coeffDict());
345 
346  CDES_.readIfPresent(this->coeffDict());
347  ck_.readIfPresent(this->coeffDict());
348 
349  return true;
350  }
351  else
352  {
353  return false;
354  }
355 }
356 
357 
358 template<class BasicTurbulenceModel>
360 DnuTildaEff() const
361 {
362  return volScalarField::New
363  (
364  "DnuTildaEff",
365  (nuTilda_ + this->nu())/sigmaNut_
366  );
367 }
368 
369 
370 template<class BasicTurbulenceModel>
372 {
373  const volScalarField chi(this->chi());
374  const volScalarField fv1(this->fv1(chi));
375  return sqr(this->nut()/ck_/dTilda(chi, fv1, fvc::grad(this->U_)));
376 }
377 
378 
379 template<class BasicTurbulenceModel>
381 {
382  const volScalarField chi(this->chi());
383  const volScalarField fv1(this->fv1(chi));
384 
385  tmp<volScalarField> tLESRegion
386  (
388  (
389  "DES::LESRegion",
390  neg(dTilda(chi, fv1, fvc::grad(this->U_)) - y_)
391  )
392  );
393 
394  return tLESRegion;
395 }
396 
397 
398 template<class BasicTurbulenceModel>
400 {
401  if (!this->turbulence_)
402  {
403  return;
404  }
405 
406  // Local references
407  const alphaField& alpha = this->alpha_;
408  const rhoField& rho = this->rho_;
409  const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_;
410  const volVectorField& U = this->U_;
411  fv::options& fvOptions(fv::options::New(this->mesh_));
412 
414 
415  const volScalarField chi(this->chi());
416  const volScalarField fv1(this->fv1(chi));
417 
418  tmp<volTensorField> tgradU = fvc::grad(U);
419  const volScalarField Omega(this->Omega(tgradU()));
420  const volScalarField dTilda(this->dTilda(chi, fv1, tgradU()));
421  const volScalarField Stilda(this->Stilda(chi, fv1, Omega, dTilda));
422 
423  tmp<fvScalarMatrix> nuTildaEqn
424  (
425  fvm::ddt(alpha, rho, nuTilda_)
426  + fvm::div(alphaRhoPhi, nuTilda_)
427  - fvm::laplacian(alpha*rho*DnuTildaEff(), nuTilda_)
428  - Cb2_/sigmaNut_*alpha*rho*magSqr(fvc::grad(nuTilda_))
429  ==
430  Cb1_*alpha*rho*Stilda*nuTilda_
431  - fvm::Sp
432  (
433  Cw1_*alpha*rho*fw(Stilda, dTilda)*nuTilda_/sqr(dTilda),
434  nuTilda_
435  )
436  + fvOptions(alpha, rho, nuTilda_)
437  );
438 
439  nuTildaEqn.ref().relax();
440  fvOptions.constrain(nuTildaEqn.ref());
441  solve(nuTildaEqn);
442  fvOptions.correct(nuTilda_);
443  bound(nuTilda_, dimensionedScalar(nuTilda_.dimensions(), 0));
444  nuTilda_.correctBoundaryConditions();
445 
446  correctNut();
447 }
448 
449 
450 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
451 
452 } // End namespace LESModels
453 } // End namespace Foam
454 
455 // ************************************************************************* //
Eddy viscosity LES SGS model base class.
scalar delta
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
SpalartAllmarasDES(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
fv::options & fvOptions
surfaceScalarField & phi
BasicTurbulenceModel::transportModel transportModel
tmp< volScalarField > Stilda(const volScalarField &chi, const volScalarField &fv1, const volScalarField &Omega, const volScalarField &dTilda) const
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedTensor skew(const dimensionedTensor &dt)
tmp< volScalarField > fv2(const volScalarField &chi, const volScalarField &fv1) const
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
tmp< volScalarField > r(const volScalarField &nur, const volScalarField &Omega, const volScalarField &dTilda) const
bool readIfPresent(const dictionary &)
Update the value of dimensioned<Type> if found in the dictionary.
virtual bool read()
Read model coefficients if they have changed.
Generic dimensioned Type class.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Finite-volume options.
Definition: fvOptions.H:52
dimensionedScalar neg(const dimensionedScalar &ds)
BasicTurbulenceModel::alphaField alphaField
static const wallDist & New(const fvMesh &mesh)
Definition: MeshObject.C:44
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Mesh &, const dimensionSet &, const word &patchFieldType=fvPatchField< scalar >::calculatedType())
Return a temporary field constructed from name, mesh, dimensionSet.
virtual tmp< volScalarField > LESRegion() const
Return the LES field indicator.
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
scalar y
const dimensionSet & dimensions() const
Return dimensions.
A class for handling words, derived from string.
Definition: word.H:59
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
tmp< volScalarField > fv1(const volScalarField &chi) const
tmp< volScalarField > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
tmp< volScalarField > S(const volTensorField &gradU) const
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
U
Definition: pEqn.H:72
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:33
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
void correctBoundaryConditions()
Correct boundary field.
dimensionedScalar pow6(const dimensionedScalar &ds)
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
BasicTurbulenceModel::rhoField rhoField
dimensioned< scalar > mag(const dimensioned< Type > &)
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
Definition: fvOptions.C:101
const dimensionedVector & g
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
tmp< volScalarField > chi() const
virtual tmp< volScalarField > dTilda(const volScalarField &chi, const volScalarField &fv1, const volTensorField &gradU) const
Length scale.
virtual void correct()
Correct nuTilda and related properties.
volScalarField & nu
scalar nut
tmp< volScalarField > fw(const volScalarField &Omega, const volScalarField &dTilda) const
tmp< volScalarField > Omega(const volTensorField &gradU) const
Namespace for OpenFOAM.
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.