SpalartAllmarasDES.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2015-2019 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Class
28 Foam::LESModels::SpalartAllmarasDES
29
30Group
31 grpDESTurbulence
32
33Description
34 SpalartAllmarasDES DES turbulence model for incompressible and
35 compressible flows
36
37 Reference:
38 \verbatim
39 Spalart, P. R., Jou, W. H., Strelets, M., & Allmaras, S. R. (1997).
40 Comments on the feasibility of LES for wings, and on a hybrid
41 RANS/LES approach.
42 Advances in DNS/LES, 1, 4-8.
43 \endverbatim
44
45 Including the low Reynolds number correction documented in
46 \verbatim
47 Spalart, P. R., Deck, S., Shur, M.L., Squires, K.D., Strelets, M.Kh,
48 Travin, A. (2006).
49 A new version of detached-eddy simulation, resistant to ambiguous grid
50 densities.
51 Theor. Comput. Fluid Dyn., 20, 181-195.
52 \endverbatim
53
54Note
55 The default value of the DES constant implemented was calibrated for
56 OpenFOAM using decaying isotropic turbulence and agrees with the value
57 suggested in the reference publication.
58
59SourceFiles
60 SpalartAllmarasDES.C
61
62\*---------------------------------------------------------------------------*/
63
64#ifndef SpalartAllmarasDES_H
65#define SpalartAllmarasDES_H
66
67#include "DESModel.H"
68
69// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
70
71namespace Foam
72{
73namespace LESModels
74{
75
76/*---------------------------------------------------------------------------*\
77 Class SpalartAllmarasDES Declaration
78\*---------------------------------------------------------------------------*/
79
80template<class BasicTurbulenceModel>
82:
83 public DESModel<BasicTurbulenceModel>
84{
85 // Private Member Functions
86
87 //- No copy construct
89
90 //- No copy assignment
91 void operator=(const SpalartAllmarasDES&) = delete;
92
93
94protected:
95
96 // Protected data
97
98 // Model constants
112
113
114 // Low Reynolds number correction
120
121
122 // Fields
125
126 //- Wall distance
127 // Note: different to wall distance in parent RASModel
128 // which is for near-wall cells only
129 const volScalarField& y_;
130
131
132 // Protected Member Functions
133
134 tmp<volScalarField> chi() const;
135
137
139 (
140 const volScalarField& chi,
141 const volScalarField& fv1
142 ) const;
143
145
146 tmp<volScalarField> Omega(const volTensorField& gradU) const;
147
149 (
150 const volScalarField& chi,
151 const volScalarField& fv1,
152 const volScalarField& Omega,
154 ) const;
155
157 (
158 const volScalarField& nur,
159 const volScalarField& Omega,
161 ) const;
162
164 (
165 const volScalarField& Omega,
167 ) const;
168
170 (
171 const volScalarField& chi,
172 const volScalarField& fv1
173 ) const;
174
175 //- Length scale
177 (
178 const volScalarField& chi,
179 const volScalarField& fv1,
180 const volTensorField& gradU
181 ) const;
182
183 void correctNut(const volScalarField& fv1);
184 virtual void correctNut();
185
186
187public:
189 typedef typename BasicTurbulenceModel::alphaField alphaField;
190 typedef typename BasicTurbulenceModel::rhoField rhoField;
191 typedef typename BasicTurbulenceModel::transportModel transportModel;
192
193
194 //- Runtime type information
195 TypeName("SpalartAllmarasDES");
196
197
198 // Constructors
199
200 //- Construct from components
202 (
203 const alphaField& alpha,
204 const rhoField& rho,
205 const volVectorField& U,
206 const surfaceScalarField& alphaRhoPhi,
207 const surfaceScalarField& phi,
208 const transportModel& transport,
209 const word& propertiesName = turbulenceModel::propertiesName,
210 const word& type = typeName
211 );
212
213
214 //- Destructor
215 virtual ~SpalartAllmarasDES() = default;
216
217
218 // Member Functions
219
220 //- Re-read model coefficients if they have changed
221 virtual bool read();
222
223 //- Return the effective diffusivity for nuTilda
225
226 //- Return SGS kinetic energy
227 virtual tmp<volScalarField> k() const;
230 {
231 return nuTilda_;
232 }
233
234 //- Return the LES field indicator
235 virtual tmp<volScalarField> LESRegion() const;
236
237 //- Correct nuTilda and related properties
238 virtual void correct();
239};
240
241
242// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
243
244} // End namespace LESModels
245} // End namespace Foam
246
247// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
248
249#ifdef NoRepository
250 #include "SpalartAllmarasDES.C"
251#endif
252
253// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
254
255#endif
256
257// ************************************************************************* //
surfaceScalarField & phi
Templated abstract base class for DES models.
Definition: DESModel.H:58
SpalartAllmarasDES DES turbulence model for incompressible and compressible flows.
BasicTurbulenceModel::alphaField alphaField
tmp< volScalarField > chi() const
BasicTurbulenceModel::rhoField rhoField
tmp< volScalarField > fw(const volScalarField &Omega, const volScalarField &dTilda) const
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.
virtual void correct()
Correct nuTilda and related properties.
tmp< volScalarField > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
tmp< volScalarField > Stilda(const volScalarField &chi, const volScalarField &fv1, const volScalarField &Omega, const volScalarField &dTilda) const
tmp< volScalarField > Omega(const volTensorField &gradU) const
tmp< volScalarField > nuTilda() const
virtual ~SpalartAllmarasDES()=default
Destructor.
tmp< volScalarField > fv1(const volScalarField &chi) const
virtual tmp< volScalarField > LESRegion() const
Return the LES field indicator.
const volScalarField & y_
Wall distance.
tmp< volScalarField > r(const volScalarField &nur, const volScalarField &Omega, const volScalarField &dTilda) const
TypeName("SpalartAllmarasDES")
Runtime type information.
BasicTurbulenceModel::transportModel transportModel
tmp< volScalarField > ft2(const volScalarField &chi) const
tmp< volScalarField > fv2(const volScalarField &chi, const volScalarField &fv1) const
virtual tmp< volScalarField > dTilda(const volScalarField &chi, const volScalarField &fv1, const volTensorField &gradU) const
Length scale.
virtual bool read()
Re-read model coefficients if they have changed.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:78
A class for managing temporary objects.
Definition: tmp.H:65
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from Foam::string.
Definition: word.H:68
U
Definition: pEqn.H:72
const volScalarField & psi
Namespace for OpenFOAM.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
volScalarField & alpha
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73