kOmegaSST.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-2015 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::kOmegaSST
26 
27 Group
28  grpRASTurbulence
29 
30 Description
31  Implementation of the k-omega-SST turbulence model for
32  incompressible and compressible flows.
33 
34  Turbulence model described in:
35  \verbatim
36  Menter, F. R. & Esch, T. (2001).
37  Elements of Industrial Heat Transfer Prediction.
38  16th Brazilian Congress of Mechanical Engineering (COBEM).
39  \endverbatim
40 
41  with updated coefficients from
42  \verbatim
43  Menter, F. R., Kuntz, M., and Langtry, R. (2003).
44  Ten Years of Industrial Experience with the SST Turbulence Model.
45  Turbulence, Heat and Mass Transfer 4, ed: K. Hanjalic, Y. Nagano,
46  & M. Tummers, Begell House, Inc., 625 - 632.
47  \endverbatim
48 
49  but with the consistent production terms from the 2001 paper as form in the
50  2003 paper is a typo, see
51  \verbatim
52  http://turbmodels.larc.nasa.gov/sst.html
53  \endverbatim
54 
55  and the addition of the optional F3 term for rough walls from
56  \verbatim
57  Hellsten, A. (1998).
58  "Some Improvements in Menter’s k-omega-SST turbulence model"
59  29th AIAA Fluid Dynamics Conference, AIAA-98-2554.
60  \endverbatim
61 
62  Note that this implementation is written in terms of alpha diffusion
63  coefficients rather than the more traditional sigma (alpha = 1/sigma) so
64  that the blending can be applied to all coefficuients in a consistent
65  manner. The paper suggests that sigma is blended but this would not be
66  consistent with the blending of the k-epsilon and k-omega models.
67 
68  Also note that the error in the last term of equation (2) relating to
69  sigma has been corrected.
70 
71  Wall-functions are applied in this implementation by using equations (14)
72  to specify the near-wall omega as appropriate.
73 
74  The blending functions (15) and (16) are not currently used because of the
75  uncertainty in their origin, range of applicability and that if y+ becomes
76  sufficiently small blending u_tau in this manner clearly becomes nonsense.
77 
78  The default model coefficients are
79  \verbatim
80  kOmegaSSTCoeffs
81  {
82  alphaK1 0.85;
83  alphaK2 1.0;
84  alphaOmega1 0.5;
85  alphaOmega2 0.856;
86  beta1 0.075;
87  beta2 0.0828;
88  betaStar 0.09;
89  gamma1 5/9;
90  gamma2 0.44;
91  a1 0.31;
92  b1 1.0;
93  c1 10.0;
94  F3 no;
95  }
96  \endverbatim
97 
98 SourceFiles
99  kOmegaSST.C
100 
101 \*---------------------------------------------------------------------------*/
102 
103 #ifndef kOmegaSST_H
104 #define kOmegaSST_H
105 
106 #include "RASModel.H"
107 #include "eddyViscosity.H"
108 
109 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
110 
111 namespace Foam
112 {
113 namespace RASModels
114 {
115 
116 /*---------------------------------------------------------------------------*\
117  Class kOmegaSST Declaration
118 \*---------------------------------------------------------------------------*/
119 
120 template<class BasicTurbulenceModel>
121 class kOmegaSST
122 :
123  public eddyViscosity<RASModel<BasicTurbulenceModel> >
124 {
125  // Private Member Functions
126 
127  // Disallow default bitwise copy construct and assignment
128  kOmegaSST(const kOmegaSST&);
129  kOmegaSST& operator=(const kOmegaSST&);
130 
131 
132 protected:
133 
134  // Protected data
135 
136  // Model coefficients
156  Switch F3_;
157 
158 
159  // Fields
160 
161  //- Wall distance
162  // Note: different to wall distance in parent RASModel
163  // which is for near-wall cells only
164  const volScalarField& y_;
168 
169 
170  // Private Member Functions
171 
172  tmp<volScalarField> F1(const volScalarField& CDkOmega) const;
173  tmp<volScalarField> F2() const;
174  tmp<volScalarField> F3() const;
175  tmp<volScalarField> F23() const;
176 
178  (
179  const volScalarField& F1,
180  const dimensionedScalar& psi1,
181  const dimensionedScalar& psi2
182  ) const
183  {
184  return F1*(psi1 - psi2) + psi2;
185  }
188  {
189  return blend(F1, alphaK1_, alphaK2_);
190  }
193  {
194  return blend(F1, alphaOmega1_, alphaOmega2_);
195  }
198  {
199  return blend(F1, beta1_, beta2_);
200  }
203  {
204  return blend(F1, gamma1_, gamma2_);
205  }
206 
207  void correctNut(const volScalarField& S2);
208 
209 
210  // Protected Member Functions
211 
212  virtual void correctNut();
213  virtual tmp<fvScalarMatrix> kSource() const;
214  virtual tmp<fvScalarMatrix> omegaSource() const;
215  virtual tmp<fvScalarMatrix> Qsas
216  (
217  const volScalarField& S2,
218  const volScalarField& gamma,
219  const volScalarField& beta
220  ) const;
221 
222 
223 public:
225  typedef typename BasicTurbulenceModel::alphaField alphaField;
226  typedef typename BasicTurbulenceModel::rhoField rhoField;
227  typedef typename BasicTurbulenceModel::transportModel transportModel;
228 
229 
230  //- Runtime type information
231  TypeName("kOmegaSST");
232 
233 
234  // Constructors
235 
236  //- Construct from components
237  kOmegaSST
238  (
239  const alphaField& alpha,
240  const rhoField& rho,
241  const volVectorField& U,
242  const surfaceScalarField& alphaRhoPhi,
243  const surfaceScalarField& phi,
244  const transportModel& transport,
245  const word& propertiesName = turbulenceModel::propertiesName,
246  const word& type = typeName
247  );
248 
249 
250  //- Destructor
251  virtual ~kOmegaSST()
252  {}
253 
254 
255  // Member Functions
256 
257  //- Re-read model coefficients if they have changed
258  virtual bool read();
259 
260  //- Return the effective diffusivity for k
262  {
263  return tmp<volScalarField>
264  (
265  new volScalarField("DkEff", alphaK(F1)*this->nut_ + this->nu())
266  );
267  }
268 
269  //- Return the effective diffusivity for omega
271  {
272  return tmp<volScalarField>
273  (
274  new volScalarField
275  (
276  "DomegaEff",
277  alphaOmega(F1)*this->nut_ + this->nu()
278  )
279  );
280  }
281 
282  //- Return the turbulence kinetic energy
283  virtual tmp<volScalarField> k() const
284  {
285  return k_;
286  }
287 
288  //- Return the turbulence kinetic energy dissipation rate
289  virtual tmp<volScalarField> epsilon() const
290  {
291  return tmp<volScalarField>
292  (
293  new volScalarField
294  (
295  IOobject
296  (
297  "epsilon",
298  this->mesh_.time().timeName(),
299  this->mesh_
300  ),
301  betaStar_*k_*omega_,
302  omega_.boundaryField().types()
303  )
304  );
305  }
306 
307  //- Return the turbulence kinetic energy dissipation rate
308  virtual tmp<volScalarField> omega() const
309  {
310  return omega_;
311  }
312 
313  //- Solve the turbulence equations and correct the turbulence viscosity
314  virtual void correct();
315 };
316 
317 
318 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
319 
320 } // End namespace RASModels
321 } // End namespace Foam
322 
323 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
324 #ifdef NoRepository
325 # include "kOmegaSST.C"
326 #endif
327 
328 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
329 #endif
330 
331 // ************************************************************************* //
dimensionedScalar alphaK2_
Definition: kOmegaSST.H:138
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
Definition: kOmegaSST.H:288
virtual tmp< fvScalarMatrix > omegaSource() const
Definition: kOmegaSST.C:144
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kOmegaSST.C:390
dimensionedScalar beta1_
Definition: kOmegaSST.H:146
virtual void correctNut()
Definition: kOmegaSST.C:123
Implementation of the k-omega-SST turbulence model for incompressible and compressible flows...
Definition: kOmegaSST.H:120
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
tmp< volScalarField > F2() const
tmp< volScalarField > F23() const
BasicTurbulenceModel::rhoField rhoField
Definition: kOmegaSST.H:225
tmp< volScalarField > gamma(const volScalarField &F1) const
Definition: kOmegaSST.H:201
A class for handling words, derived from string.
Definition: word.H:59
BasicTurbulenceModel::transportModel transportModel
Definition: kOmegaSST.H:226
psi2
Definition: TEqns.H:35
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
tmp< volScalarField > DkEff(const volScalarField &F1) const
Return the effective diffusivity for k.
Definition: kOmegaSST.H:260
Namespace for OpenFOAM.
virtual tmp< fvScalarMatrix > Qsas(const volScalarField &S2, const volScalarField &gamma, const volScalarField &beta) const
Definition: kOmegaSST.C:159
dimensionedScalar gamma2_
Definition: kOmegaSST.H:144
dimensionedScalar b1_
Definition: kOmegaSST.H:152
dimensionedScalar c1_
Definition: kOmegaSST.H:153
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
dimensionedScalar a1_
Definition: kOmegaSST.H:151
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
Definition: kOmegaSST.H:282
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:52
const volScalarField & y_
Wall distance.
Definition: kOmegaSST.H:163
virtual tmp< volScalarField > omega() const
Return the turbulence kinetic energy dissipation rate.
Definition: kOmegaSST.H:307
dimensionedScalar alphaOmega1_
Definition: kOmegaSST.H:140
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
surfaceScalarField & phi
tmp< volScalarField > F3() const
virtual tmp< fvScalarMatrix > kSource() const
Definition: kOmegaSST.C:130
dimensionedScalar beta2_
Definition: kOmegaSST.H:147
tmp< volScalarField > beta(const volScalarField &F1) const
Definition: kOmegaSST.H:196
dimensionedScalar alphaOmega2_
Definition: kOmegaSST.H:141
TypeName("kOmegaSST")
Runtime type information.
tmp< volScalarField > blend(const volScalarField &F1, const dimensionedScalar &psi1, const dimensionedScalar &psi2) const
Definition: kOmegaSST.H:177
psi1
Definition: TEqns.H:34
static const word propertiesName
Default name of the turbulence properties dictionary.
volScalarField & nu
dimensionedScalar gamma1_
Definition: kOmegaSST.H:143
virtual ~kOmegaSST()
Destructor.
Definition: kOmegaSST.H:250
tmp< volScalarField > alphaOmega(const volScalarField &F1) const
Definition: kOmegaSST.H:191
BasicTurbulenceModel::alphaField alphaField
Definition: kOmegaSST.H:224
volScalarField omega_
Definition: kOmegaSST.H:166
tmp< volScalarField > alphaK(const volScalarField &F1) const
Definition: kOmegaSST.H:186
A class for managing temporary objects.
Definition: PtrList.H:118
U
Definition: pEqn.H:82
tmp< volScalarField > DomegaEff(const volScalarField &F1) const
Return the effective diffusivity for omega.
Definition: kOmegaSST.H:269
virtual bool read()
Re-read model coefficients if they have changed.
Definition: kOmegaSST.C:362
dimensionedScalar alphaK1_
Definition: kOmegaSST.H:137
dimensionedScalar betaStar_
Definition: kOmegaSST.H:149
tmp< volScalarField > F1(const volScalarField &CDkOmega) const