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