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
incompressible
LES
kOmegaSSTSAS
kOmegaSSTSAS.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) 2008-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::LESmodels::kOmegaSSTSAS
26
27
Description
28
kOmegaSSTSAS LES turbulence model for incompressible flows
29
30
Reference:
31
32
DESider A European Effort on Hybrid RANS-LES Modelling:
33
Results of the European-Union Funded Project, 2004 - 2007
34
(Notes on Numerical Fluid Mechanics and Multidisciplinary Design).
35
Chapter 8 Formulation of the Scale-Adaptive Simulation (SAS) Model during
36
the DESIDER Project. Published in Springer-Verlag Berlin Heidelberg 2009.
37
F. R. Menter and Y. Egorov.
38
39
SourceFiles
40
kOmegaSSTSAS.C
41
42
\*---------------------------------------------------------------------------*/
43
44
#ifndef kOmegaSSTSAS_H
45
#define kOmegaSSTSAS_H
46
47
#include <
incompressibleLESModels/LESModel.H
>
48
#include <
finiteVolume/volFields.H
>
49
#include <
finiteVolume/wallDist.H
>
50
51
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
52
53
namespace
Foam
54
{
55
namespace
incompressible
56
{
57
namespace
LESModels
58
{
59
60
/*---------------------------------------------------------------------------*\
61
Class kOmegaSSTSAS Declaration
62
\*---------------------------------------------------------------------------*/
63
64
class
kOmegaSSTSAS
65
:
66
public
LESModel
67
{
68
// Private member functions
69
70
//- Update sub-grid scale fields
71
void
updateSubGridScaleFields(
const
volScalarField
& D);
72
73
// Disallow default bitwise copy construct and assignment
74
kOmegaSSTSAS
(
const
kOmegaSSTSAS
&);
75
kOmegaSSTSAS
& operator=(
const
kOmegaSSTSAS
&);
76
77
78
protected
:
79
80
// Protected data
81
82
// Model constants
83
84
dimensionedScalar
alphaK1_
;
85
dimensionedScalar
alphaK2_
;
86
87
dimensionedScalar
alphaOmega1_
;
88
dimensionedScalar
alphaOmega2_
;
89
90
dimensionedScalar
gamma1_
;
91
dimensionedScalar
gamma2_
;
92
93
dimensionedScalar
beta1_
;
94
dimensionedScalar
beta2_
;
95
96
dimensionedScalar
betaStar_
;
97
98
dimensionedScalar
a1_
;
99
dimensionedScalar
c1_
;
100
dimensionedScalar
Cs_
;
101
102
dimensionedScalar
alphaPhi_
;
103
dimensionedScalar
zetaTilda2_
;
104
dimensionedScalar
FSAS_
;
105
106
dimensionedScalar
omega0_
;
107
dimensionedScalar
omegaSmall_
;
108
109
wallDist
y_
;
110
dimensionedScalar
Cmu_
;
111
dimensionedScalar
kappa_
;
112
113
114
// Fields
115
116
volScalarField
k_
;
117
volScalarField
omega_
;
118
volScalarField
nuSgs_
;
119
120
121
// Protected member functions
122
123
tmp<volScalarField>
Lvk2
124
(
125
const
volScalarField
& S2
126
)
const
;
127
128
tmp<volScalarField>
F1
(
const
volScalarField
& CDkOmega)
const
;
129
tmp<volScalarField>
F2
()
const
;
130
131
tmp<volScalarField>
blend
132
(
133
const
volScalarField
&
F1
,
134
const
dimensionedScalar
& psi1,
135
const
dimensionedScalar
& psi2
136
)
const
137
{
138
return
F1*(psi1 - psi2) + psi2;
139
}
140
141
tmp<volScalarField>
alphaK
142
(
143
const
volScalarField
& F1
144
)
const
145
{
146
return
blend
(F1,
alphaK1_
,
alphaK2_
);
147
}
148
149
tmp<volScalarField>
alphaOmega
150
(
151
const
volScalarField
& F1
152
)
const
153
{
154
return
blend
(F1,
alphaOmega1_
,
alphaOmega2_
);
155
}
156
157
tmp<volScalarField>
beta
158
(
159
const
volScalarField
& F1
160
)
const
161
{
162
return
blend
(F1,
beta1_
,
beta2_
);
163
}
164
165
tmp<volScalarField>
gamma
166
(
167
const
volScalarField
& F1
168
)
const
169
{
170
return
blend
(F1,
gamma1_
,
gamma2_
);
171
}
172
173
174
public
:
175
176
//- Runtime type information
177
TypeName
(
"kOmegaSSTSAS"
);
178
179
180
// Constructors
181
182
//- Construct from components
183
kOmegaSSTSAS
184
(
185
const
volVectorField
&
U
,
186
const
surfaceScalarField
&
phi
,
187
transportModel
&
transport
,
188
const
word
& modelName = typeName
189
);
190
191
192
//- Destructor
193
virtual
~kOmegaSSTSAS
()
194
{}
195
196
197
// Member Functions
198
199
//- Return SGS kinetic energy
200
virtual
tmp<volScalarField>
k
()
const
201
{
202
return
k_
;
203
}
204
205
//- Return omega
206
virtual
tmp<volScalarField>
omega
()
const
207
{
208
return
omega_
;
209
}
210
211
//- Return the effective diffusivity for k
212
tmp<volScalarField>
DkEff
(
const
volScalarField
& F1)
const
213
{
214
return
tmp<volScalarField>
215
(
216
new
volScalarField
(
"DkEff"
,
alphaK
(F1)*
nuSgs_
+
nu
())
217
);
218
}
219
220
//- Return the effective diffusivity for omega
221
tmp<volScalarField>
DomegaEff
(
const
volScalarField
& F1)
const
222
{
223
return
tmp<volScalarField>
224
(
225
new
volScalarField
(
"DomegaEff"
,
alphaOmega
(F1)*
nuSgs_
+
nu
())
226
);
227
}
228
229
//- Return sub-grid disipation rate
230
virtual
tmp<volScalarField>
epsilon
()
const
;
231
232
//- Return SGS viscosity
233
virtual
tmp<volScalarField>
nuSgs
()
const
234
{
235
return
nuSgs_
;
236
}
237
238
//- Return the sub-grid stress tensor.
239
virtual
tmp<volSymmTensorField>
B
()
const
;
240
241
//- Return the effective sub-grid turbulence stress tensor
242
// including the laminar stress
243
virtual
tmp<volSymmTensorField>
devBeff
()
const
;
244
245
//- Return the deviatoric part of the divergence of Beff
246
// i.e. the additional term in the filtered NSE.
247
virtual
tmp<fvVectorMatrix>
divDevBeff
(
volVectorField
&
U
)
const
;
248
249
//- Solve the turbulence equations (k-w) and correct the turbulence
250
// viscosity
251
virtual
void
correct
(
const
tmp<volTensorField>
& gradU);
252
253
//- Read LESProperties dictionary
254
virtual
bool
read
();
255
};
256
257
258
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
259
260
}
// End namespace LESModels
261
}
// End namespace incompressible
262
}
// End namespace Foam
263
264
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
265
266
#endif
267
268
// ************************ vim: set sw=4 sts=4 et: ************************ //