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-2021 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 BasicMomentumTransportModel>
62 :
63  public LESeddyViscosity<BasicMomentumTransportModel>
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  using IOobject::modelName;
97 
98  tmp<volScalarField> chi() const;
99 
101 
103  (
106  ) const;
107 
109  (
110  const volTensorField::Internal& gradU
111  ) const;
112 
114  (
119  ) const;
120 
122  (
123  const volScalarField::Internal& nur,
126  ) const;
127 
129  (
132  ) const;
133 
134  //- Length scale
136  (
139  const volTensorField::Internal& gradU
140  ) const;
141 
142  //- Cache the LES region indicator field
143  virtual void cacheLESRegion
144  (
146  ) const;
147 
148  void correctNut(const volScalarField& fv1);
149  virtual void correctNut();
150 
151 
152 public:
154  typedef typename BasicMomentumTransportModel::alphaField alphaField;
155  typedef typename BasicMomentumTransportModel::rhoField rhoField;
156 
157 
158  //- Runtime type information
159  TypeName("SpalartAllmarasDES");
160 
161 
162  // Constructors
163 
164  //- Construct from components
166  (
167  const alphaField& alpha,
168  const rhoField& rho,
169  const volVectorField& U,
170  const surfaceScalarField& alphaRhoPhi,
171  const surfaceScalarField& phi,
172  const viscosity& viscosity,
173  const word& type = typeName
174  );
175 
176  //- Disallow default bitwise copy construction
177  SpalartAllmarasDES(const SpalartAllmarasDES&) = delete;
178 
179 
180  //- Destructor
181  virtual ~SpalartAllmarasDES()
182  {}
183 
184 
185  // Member Functions
186 
187  //- Read model coefficients if they have changed
188  virtual bool read();
189 
190  //- Return the effective diffusivity for nuTilda
192 
193  //- Return SGS kinetic energy
194  virtual tmp<volScalarField> k() const;
197  {
198  return nuTilda_;
199  }
200 
201  //- Correct nuTilda and related properties
202  virtual void correct();
203 
204 
205  // Member Operators
206 
207  //- Disallow default bitwise assignment
208  void operator=(const SpalartAllmarasDES&) = delete;
209 };
210 
211 
212 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
213 
214 } // End namespace LESModels
215 } // End namespace Foam
216 
217 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
218 
219 #ifdef NoRepository
220  #include "SpalartAllmarasDES.C"
221 #endif
222 
223 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
224 
225 #endif
226 
227 // ************************************************************************* //
Eddy viscosity LES SGS model base class.
TypeName("SpalartAllmarasDES")
Runtime type information.
void operator=(const SpalartAllmarasDES &)=delete
Disallow default bitwise assignment.
U
Definition: pEqn.H:72
const volScalarField & y_
Wall distance.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
virtual bool read()
Read model coefficients if they have changed.
SpalartAllmarasDES DES turbulence model for incompressible and compressible flows.
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
A class for handling words, derived from string.
Definition: word.H:59
BasicMomentumTransportModel::alphaField alphaField
virtual tmp< volScalarField::Internal > dTilda(const volScalarField::Internal &chi, const volScalarField::Internal &fv1, const volTensorField::Internal &gradU) const
Length scale.
phi
Definition: correctPhi.H:3
tmp< volScalarField::Internal > r(const volScalarField::Internal &nur, const volScalarField::Internal &Omega, const volScalarField::Internal &dTilda) const
tmp< volScalarField > fv1(const volScalarField &chi) const
Abstract base class for all fluid physical properties.
Definition: viscosity.H:49
tmp< volScalarField > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
BasicMomentumTransportModel::rhoField rhoField
tmp< volScalarField::Internal > Omega(const volTensorField::Internal &gradU) const
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
virtual void cacheLESRegion(const volScalarField::Internal &dTilda) const
Cache the LES region indicator field.
tmp< volScalarField::Internal > Stilda(const volScalarField::Internal &chi, const volScalarField::Internal &fv1, const volScalarField::Internal &Omega, const volScalarField::Internal &dTilda) const
SpalartAllmarasDES(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity, const word &type=typeName)
Construct from components.
static word modelName(Name name, const word &model)
Return the name of the object within the given model.
A class for managing temporary objects.
Definition: PtrList.H:53
tmp< volScalarField > chi() const
virtual void correct()
Correct nuTilda and related properties.
tmp< volScalarField::Internal > fw(const volScalarField::Internal &Omega, const volScalarField::Internal &dTilda) const
tmp< volScalarField::Internal > fv2(const volScalarField::Internal &chi, const volScalarField::Internal &fv1) const
Namespace for OpenFOAM.
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.