SpalartAllmarasDES.H
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-2019 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 Class
25  Foam::LESModels::SpalartAllmarasDES
26 
27 Description
28  SpalartAllmarasDES DES turbulence model for incompressible and
29  compressible flows
30 
31  Reference:
32  \verbatim
33  Spalart, P. R., Jou, W. H., Strelets, M., & Allmaras, S. R. (1997).
34  Comments on the feasibility of LES for wings, and on a hybrid
35  RANS/LES approach.
36  Advances in DNS/LES, 1, 4-8.
37  \endverbatim
38 
39 SourceFiles
40  SpalartAllmarasDES.C
41 
42 \*---------------------------------------------------------------------------*/
43 
44 #ifndef SpalartAllmarasDES_H
45 #define SpalartAllmarasDES_H
46 
47 #include "LESeddyViscosity.H"
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 namespace Foam
52 {
53 namespace LESModels
54 {
55 
56 /*---------------------------------------------------------------------------*\
57  Class SpalartAllmarasDES Declaration
58 \*---------------------------------------------------------------------------*/
59 
60 template<class BasicTurbulenceModel>
62 :
63  public LESeddyViscosity<BasicTurbulenceModel>
64 {
65 protected:
66 
67  // Protected data
68 
69  // Model constants
70 
73 
83 
84  // Fields
85 
87 
88  //- Wall distance
89  // Note: different to wall distance in parent RASModel
90  // which is for near-wall cells only
91  const volScalarField& y_;
92 
93 
94  // Protected Member Functions
95 
96  tmp<volScalarField> chi() const;
97 
99 
101  (
102  const volScalarField& chi,
103  const volScalarField& fv1
104  ) const;
105 
106  tmp<volScalarField> S(const volTensorField& gradU) const;
107 
108  tmp<volScalarField> Omega(const volTensorField& gradU) const;
109 
111  (
112  const volScalarField& chi,
113  const volScalarField& fv1,
114  const volScalarField& Omega,
115  const volScalarField& dTilda
116  ) const;
117 
119  (
120  const volScalarField& nur,
121  const volScalarField& Omega,
122  const volScalarField& dTilda
123  ) const;
124 
126  (
127  const volScalarField& Omega,
128  const volScalarField& dTilda
129  ) const;
130 
131  //- Length scale
133  (
134  const volScalarField& chi,
135  const volScalarField& fv1,
136  const volTensorField& gradU
137  ) const;
138 
139  void correctNut(const volScalarField& fv1);
140  virtual void correctNut();
141 
142 
143 public:
145  typedef typename BasicTurbulenceModel::alphaField alphaField;
146  typedef typename BasicTurbulenceModel::rhoField rhoField;
147  typedef typename BasicTurbulenceModel::transportModel transportModel;
148 
149 
150  //- Runtime type information
151  TypeName("SpalartAllmarasDES");
152 
153 
154  // Constructors
155 
156  //- Construct from components
158  (
159  const alphaField& alpha,
160  const rhoField& rho,
161  const volVectorField& U,
162  const surfaceScalarField& alphaRhoPhi,
163  const surfaceScalarField& phi,
164  const transportModel& transport,
165  const word& propertiesName = turbulenceModel::propertiesName,
166  const word& type = typeName
167  );
168 
169  //- Disallow default bitwise copy construction
170  SpalartAllmarasDES(const SpalartAllmarasDES&) = delete;
171 
172 
173  //- Destructor
174  virtual ~SpalartAllmarasDES()
175  {}
176 
177 
178  // Member Functions
179 
180  //- Read model coefficients if they have changed
181  virtual bool read();
182 
183  //- Return the effective diffusivity for nuTilda
185 
186  //- Return SGS kinetic energy
187  virtual tmp<volScalarField> k() const;
190  {
191  return nuTilda_;
192  }
193 
194  //- Return the LES field indicator
195  virtual tmp<volScalarField> LESRegion() const;
196 
197  //- Correct nuTilda and related properties
198  virtual void correct();
199 
200 
201  // Member Operators
202 
203  //- Disallow default bitwise assignment
204  void operator=(const SpalartAllmarasDES&) = delete;
205 };
206 
207 
208 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
209 
210 } // End namespace LESModels
211 } // End namespace Foam
212 
213 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
214 
215 #ifdef NoRepository
216  #include "SpalartAllmarasDES.C"
217 #endif
218 
219 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
220 
221 #endif
222 
223 // ************************************************************************* //
Eddy viscosity LES SGS model base class.
TypeName("SpalartAllmarasDES")
Runtime type information.
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.
void operator=(const SpalartAllmarasDES &)=delete
Disallow default bitwise assignment.
surfaceScalarField & phi
tmp< volScalarField > Stilda(const volScalarField &chi, const volScalarField &fv1, const volScalarField &Omega, const volScalarField &dTilda) const
const volScalarField & y_
Wall distance.
tmp< volScalarField > fv2(const volScalarField &chi, const volScalarField &fv1) const
tmp< volScalarField > r(const volScalarField &nur, const volScalarField &Omega, const volScalarField &dTilda) const
virtual bool read()
Read model coefficients if they have changed.
SpalartAllmarasDES DES turbulence model for incompressible and compressible flows.
BasicTurbulenceModel::transportModel transportModel
virtual tmp< volScalarField > LESRegion() const
Return the LES field indicator.
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
tmp< volScalarField > fv1(const volScalarField &chi) const
tmp< volScalarField > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
tmp< volScalarField > S(const volTensorField &gradU) const
BasicTurbulenceModel::rhoField rhoField
U
Definition: pEqn.H:72
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
tmp< volScalarField > nuTilda() const
BasicTurbulenceModel::alphaField alphaField
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
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.
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.