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 | Copyright (C) 2011-2016 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 Group
28  grpDESTurbulence
29 
30 Description
31  SpalartAllmarasDES DES turbulence model for incompressible and
32  compressible flows
33 
34  Reference:
35  \verbatim
36  Spalart, P. R., Jou, W. H., Strelets, M., & Allmaras, S. R. (1997).
37  Comments on the feasibility of LES for wings, and on a hybrid
38  RANS/LES approach.
39  Advances in DNS/LES, 1, 4-8.
40  \endverbatim
41 
42 SourceFiles
43  SpalartAllmarasDES.C
44 
45 \*---------------------------------------------------------------------------*/
46 
47 #ifndef SpalartAllmarasDES_H
48 #define SpalartAllmarasDES_H
49 
50 #include "LESeddyViscosity.H"
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 namespace LESModels
57 {
58 
59 /*---------------------------------------------------------------------------*\
60  Class SpalartAllmarasDES Declaration
61 \*---------------------------------------------------------------------------*/
62 
63 template<class BasicTurbulenceModel>
65 :
66  public LESeddyViscosity<BasicTurbulenceModel>
67 {
68  // Private Member Functions
69 
70  // Disallow default bitwise copy construct and assignment
72  void operator=(const SpalartAllmarasDES&);
73 
74 
75 protected:
76 
77  // Protected data
78 
79  // Model constants
80 
83 
93 
94  // Fields
95 
97 
98  //- Wall distance
99  // Note: different to wall distance in parent RASModel
100  // which is for near-wall cells only
101  const volScalarField& y_;
102 
103 
104  // Protected Member Functions
105 
106  tmp<volScalarField> chi() const;
107 
109 
111  (
112  const volScalarField& chi,
113  const volScalarField& fv1
114  ) const;
115 
116  tmp<volScalarField> S(const volTensorField& gradU) const;
117 
118  tmp<volScalarField> Omega(const volTensorField& gradU) const;
119 
121  (
122  const volScalarField& chi,
123  const volScalarField& fv1,
124  const volScalarField& Omega,
125  const volScalarField& dTilda
126  ) const;
127 
129  (
130  const volScalarField& nur,
131  const volScalarField& Omega,
132  const volScalarField& dTilda
133  ) const;
134 
136  (
137  const volScalarField& Omega,
138  const volScalarField& dTilda
139  ) const;
140 
141  //- Length scale
143  (
144  const volScalarField& chi,
145  const volScalarField& fv1,
146  const volTensorField& gradU
147  ) const;
148 
149  void correctNut(const volScalarField& fv1);
150  virtual void correctNut();
151 
152 
153 public:
155  typedef typename BasicTurbulenceModel::alphaField alphaField;
156  typedef typename BasicTurbulenceModel::rhoField rhoField;
157  typedef typename BasicTurbulenceModel::transportModel transportModel;
158 
159 
160  //- Runtime type information
161  TypeName("SpalartAllmarasDES");
162 
163 
164  // Constructors
165 
166  //- Construct from components
168  (
169  const alphaField& alpha,
170  const rhoField& rho,
171  const volVectorField& U,
172  const surfaceScalarField& alphaRhoPhi,
173  const surfaceScalarField& phi,
174  const transportModel& transport,
175  const word& propertiesName = turbulenceModel::propertiesName,
176  const word& type = typeName
177  );
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  //- Return the LES field indicator
202  virtual tmp<volScalarField> LESRegion() const;
203 
204  //- Correct nuTilda and related properties
205  virtual void correct();
206 };
207 
208 
209 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
210 
211 } // End namespace LESModels
212 } // End namespace Foam
213 
214 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
215 
216 #ifdef NoRepository
217  #include "SpalartAllmarasDES.C"
218 #endif
219 
220 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
221 
222 #endif
223 
224 // ************************************************************************* //
Eddy viscosity LES SGS model base class.
TypeName("SpalartAllmarasDES")
Runtime type information.
surfaceScalarField & phi
U
Definition: pEqn.H:83
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
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
BasicTurbulenceModel::rhoField rhoField
tmp< volScalarField > nuTilda() const
BasicTurbulenceModel::alphaField alphaField
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
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.