All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
SpalartAllmaras.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::RASModels::SpalartAllmaras
26 
27 Description
28  Spalart-Allmaras one-eqn mixing-length model for incompressible and
29  compressible external flows.
30 
31  Reference:
32  \verbatim
33  Spalart, P.R., & Allmaras, S.R. (1994).
34  A one-equation turbulence model for aerodynamic flows.
35  La Recherche Aerospatiale, 1, 5-21.
36  \endverbatim
37 
38  The model is implemented without the trip-term and hence the ft2 term is
39  not needed.
40 
41  It is necessary to limit the Stilda generation term as the model generates
42  unphysical results if this term becomes negative which occurs for complex
43  flow. Several approaches have been proposed to limit Stilda but it is not
44  clear which is the most appropriate. Here the limiter proposed by Spalart
45  is implemented in which Stilda is clipped at Cs*Omega with the default value
46  of Cs = 0.3.
47 
48  The default model coefficients are
49  \verbatim
50  SpalartAllmarasCoeffs
51  {
52  Cb1 0.1355;
53  Cb2 0.622;
54  Cw2 0.3;
55  Cw3 2.0;
56  Cv1 7.1;
57  Cs 0.3;
58  sigmaNut 0.66666;
59  kappa 0.41;
60  }
61  \endverbatim
62 
63 SourceFiles
64  SpalartAllmaras.C
65 
66 \*---------------------------------------------------------------------------*/
67 
68 #ifndef SpalartAllmaras_H
69 #define SpalartAllmaras_H
70 
71 #include "RASModel.H"
72 #include "eddyViscosity.H"
73 
74 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
75 
76 namespace Foam
77 {
78 namespace RASModels
79 {
80 
81 /*---------------------------------------------------------------------------*\
82  Class SpalartAllmaras Declaration
83 \*---------------------------------------------------------------------------*/
84 
85 template<class BasicTurbulenceModel>
86 class SpalartAllmaras
87 :
88  public eddyViscosity<RASModel<BasicTurbulenceModel>>
89 {
90 protected:
91 
92  // Protected data
93 
94  // Model coefficients
95 
98 
106 
107 
108  // Fields
111 
112  //- Wall distance
113  // Note: different to wall distance in parent RASModel
114  // which is for near-wall cells only
115  const volScalarField& y_;
116 
117 
118  // Protected Member Functions
119 
120  tmp<volScalarField> chi() const;
121 
123 
125  (
126  const volScalarField& chi,
127  const volScalarField& fv1
128  ) const;
129 
131  (
132  const volScalarField& chi,
133  const volScalarField& fv1
134  ) const;
135 
137 
138  void correctNut(const volScalarField& fv1);
139  virtual void correctNut();
140 
141 
142 public:
144  typedef typename BasicTurbulenceModel::alphaField alphaField;
145  typedef typename BasicTurbulenceModel::rhoField rhoField;
146  typedef typename BasicTurbulenceModel::transportModel transportModel;
147 
148 
149  //- Runtime type information
150  TypeName("SpalartAllmaras");
151 
152 
153  // Constructors
154 
155  //- Construct from components
157  (
158  const alphaField& alpha,
159  const rhoField& rho,
160  const volVectorField& U,
161  const surfaceScalarField& alphaRhoPhi,
162  const surfaceScalarField& phi,
163  const transportModel& transport,
164  const word& propertiesName = turbulenceModel::propertiesName,
165  const word& type = typeName
166  );
167 
168  //- Disallow default bitwise copy construction
169  SpalartAllmaras(const SpalartAllmaras&) = delete;
170 
171 
172  //- Destructor
173  virtual ~SpalartAllmaras()
174  {}
175 
176 
177  // Member Functions
178 
179  //- Read RASProperties dictionary
180  virtual bool read();
181 
182  //- Return the effective diffusivity for nuTilda
184 
185  //- Return the turbulence kinetic energy
186  virtual tmp<volScalarField> k() const;
187 
188  //- Return the turbulence kinetic energy dissipation rate
189  virtual tmp<volScalarField> epsilon() const;
190 
191  //- Solve the turbulence equations and correct the turbulence viscosity
192  virtual void correct();
193 
194 
195  // Member Operators
196 
197  //- Disallow default bitwise assignment
198  void operator=(const SpalartAllmaras&) = delete;
199 };
200 
201 
202 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
203 
204 } // End namespace RASModels
205 } // End namespace Foam
206 
207 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
208 
209 #ifdef NoRepository
210  #include "SpalartAllmaras.C"
211 #endif
212 
213 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
214 
215 #endif
216 
217 // ************************************************************************* //
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
surfaceScalarField & phi
TypeName("SpalartAllmaras")
Runtime type information.
tmp< volScalarField > Stilda(const volScalarField &chi, const volScalarField &fv1) const
BasicTurbulenceModel::alphaField alphaField
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:49
tmp< volScalarField > chi() const
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
const volScalarField & y_
Wall distance.
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 > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
virtual bool read()
Read RASProperties dictionary.
BasicTurbulenceModel::rhoField rhoField
BasicTurbulenceModel::transportModel transportModel
virtual ~SpalartAllmaras()
Destructor.
Spalart-Allmaras one-eqn mixing-length model for incompressible and compressible external flows...
U
Definition: pEqn.H:72
tmp< volScalarField > fv1(const volScalarField &chi) const
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
void operator=(const SpalartAllmaras &)=delete
Disallow default bitwise assignment.
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 > fv2(const volScalarField &chi, const volScalarField &fv1) const
tmp< volScalarField > fw(const volScalarField &Stilda) const
Namespace for OpenFOAM.
SpalartAllmaras(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.