FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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::incompressible::RASModels::kOmegaSST
26 
27 Description
28  Implementation of the k-omega-SST turbulence model for incompressible
29  flows.
30 
31  Turbulence model described in:
32  @verbatim
33  Menter, F., Esch, T.
34  "Elements of Industrial Heat Transfer Prediction"
35  16th Brazilian Congress of Mechanical Engineering (COBEM),
36  Nov. 2001
37  @endverbatim
38 
39  Note that this implementation is written in terms of alpha diffusion
40  coefficients rather than the more traditional sigma (alpha = 1/sigma) so
41  that the blending can be applied to all coefficuients in a consistent
42  manner. The paper suggests that sigma is blended but this would not be
43  consistent with the blending of the k-epsilon and k-omega models.
44 
45  Also note that the error in the last term of equation (2) relating to
46  sigma has been corrected.
47 
48  Wall-functions are applied in this implementation by using equations (14)
49  to specify the near-wall omega as appropriate.
50 
51  The blending functions (15) and (16) are not currently used because of the
52  uncertainty in their origin, range of applicability and that is y+ becomes
53  sufficiently small blending u_tau in this manner clearly becomes nonsense.
54 
55  The default model coefficients correspond to the following:
56  @verbatim
57  kOmegaSSTCoeffs
58  {
59  alphaK1 0.85034;
60  alphaK2 1.0;
61  alphaOmega1 0.5;
62  alphaOmega2 0.85616;
63  beta1 0.075;
64  beta2 0.0828;
65  betaStar 0.09;
66  gamma1 0.5532;
67  gamma2 0.4403;
68  a1 0.31;
69  c1 10.0;
70  }
71  @endverbatim
72 
73 SourceFiles
74  kOmegaSST.C
75 
76 \*---------------------------------------------------------------------------*/
77 
78 #ifndef kOmegaSST_H
79 #define kOmegaSST_H
80 
82 #include <finiteVolume/wallDist.H>
83 
84 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
85 
86 namespace Foam
87 {
88 namespace incompressible
89 {
90 namespace RASModels
91 {
92 
93 /*---------------------------------------------------------------------------*\
94  Class kOmega Declaration
95 \*---------------------------------------------------------------------------*/
96 
97 class kOmegaSST
98 :
99  public RASModel
100 {
101  // Private data
102 
103  // Model coefficients
104  dimensionedScalar alphaK1_;
105  dimensionedScalar alphaK2_;
106 
107  dimensionedScalar alphaOmega1_;
108  dimensionedScalar alphaOmega2_;
109 
110  dimensionedScalar gamma1_;
111  dimensionedScalar gamma2_;
112 
113  dimensionedScalar beta1_;
114  dimensionedScalar beta2_;
115 
116  dimensionedScalar betaStar_;
117 
118  dimensionedScalar a1_;
119  dimensionedScalar c1_;
120 
121  //- Wall distance field
122  // Note: different to wall distance in parent RASModel
123  wallDist y_;
124 
125  // Fields
126 
127  volScalarField k_;
128  volScalarField omega_;
129  volScalarField nut_;
130 
131 
132  // Private member functions
133 
134  tmp<volScalarField> F1(const volScalarField& CDkOmega) const;
135  tmp<volScalarField> F2() const;
136 
137  tmp<volScalarField> blend
138  (
139  const volScalarField& F1,
140  const dimensionedScalar& psi1,
141  const dimensionedScalar& psi2
142  ) const
143  {
144  return F1*(psi1 - psi2) + psi2;
145  }
146 
147  tmp<volScalarField> alphaK
148  (
149  const volScalarField& F1
150  ) const
151  {
152  return blend(F1, alphaK1_, alphaK2_);
153  }
154 
155  tmp<volScalarField> alphaOmega
156  (
157  const volScalarField& F1
158  ) const
159  {
160  return blend(F1, alphaOmega1_, alphaOmega2_);
161  }
162 
164  (
165  const volScalarField& F1
166  ) const
167  {
168  return blend(F1, beta1_, beta2_);
169  }
170 
171  tmp<volScalarField> gamma
172  (
173  const volScalarField& F1
174  ) const
175  {
176  return blend(F1, gamma1_, gamma2_);
177  }
178 
179 
180 public:
181 
182  //- Runtime type information
183  TypeName("kOmegaSST");
184 
185 
186  // Constructors
187 
188  //- Construct from components
189  kOmegaSST
190  (
191  const volVectorField& U,
192  const surfaceScalarField& phi,
194  );
195 
196 
197  //- Destructor
198  virtual ~kOmegaSST()
199  {}
200 
201 
202  // Member Functions
203 
204  //- Return the turbulence viscosity
205  virtual tmp<volScalarField> nut() const
206  {
207  return nut_;
208  }
209 
210  //- Return the effective diffusivity for k
212  {
213  return tmp<volScalarField>
214  (
215  new volScalarField("DkEff", alphaK(F1)*nut_ + nu())
216  );
217  }
218 
219  //- Return the effective diffusivity for omega
221  {
222  return tmp<volScalarField>
223  (
224  new volScalarField("DomegaEff", alphaOmega(F1)*nut_ + nu())
225  );
226  }
227 
228  //- Return the turbulence kinetic energy
229  virtual tmp<volScalarField> k() const
230  {
231  return k_;
232  }
233 
234  //- Return the turbulence specific dissipation rate
235  virtual tmp<volScalarField> omega() const
236  {
237  return omega_;
238  }
239 
240  //- Return the turbulence kinetic energy dissipation rate
241  virtual tmp<volScalarField> epsilon() const
242  {
243  return tmp<volScalarField>
244  (
245  new volScalarField
246  (
247  IOobject
248  (
249  "epsilon",
250  mesh_.time().timeName(),
251  mesh_
252  ),
253  betaStar_*k_*omega_,
254  omega_.boundaryField().types()
255  )
256  );
257  }
258 
259  //- Return the Reynolds stress tensor
260  virtual tmp<volSymmTensorField> R() const;
261 
262  //- Return the effective stress tensor including the laminar stress
263  virtual tmp<volSymmTensorField> devReff() const;
264 
265  //- Return the source term for the momentum equation
267 
268  //- Solve the turbulence equations and correct the turbulence viscosity
269  virtual void correct();
270 
271  //- Read RASProperties dictionary
272  virtual bool read();
273 };
274 
275 
276 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
277 
278 } // End namespace RASModels
279 } // namespace incompressible
280 } // End namespace Foam
281 
282 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
283 
284 #endif
285 
286 // ************************ vim: set sw=4 sts=4 et: ************************ //