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
RAS
RASModel
RASModel.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
Namespace
25
Foam::incompressible::RASModels
26
27
Description
28
Namespace for incompressible RAS turbulence models.
29
30
Class
31
Foam::incompressible::RASModel
32
33
Description
34
Abstract base class for incompressible turbulence models.
35
36
SourceFiles
37
RASModel.C
38
39
\*---------------------------------------------------------------------------*/
40
41
#ifndef RASModel_H
42
#define RASModel_H
43
44
#include <
incompressibleTurbulenceModel/turbulenceModel.H
>
45
#include <
finiteVolume/volFields.H
>
46
#include <
finiteVolume/surfaceFields.H
>
47
#include <
finiteVolume/nearWallDist.H
>
48
#include <
finiteVolume/fvm.H
>
49
#include <
finiteVolume/fvc.H
>
50
#include <
finiteVolume/fvMatrices.H
>
51
#include <
incompressibleTransportModels/transportModel.H
>
52
#include <
OpenFOAM/IOdictionary.H
>
53
#include <
OpenFOAM/Switch.H
>
54
#include <
finiteVolume/bound.H
>
55
#include <
OpenFOAM/autoPtr.H
>
56
#include <
OpenFOAM/runTimeSelectionTables.H
>
57
58
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
59
60
namespace
Foam
61
{
62
namespace
incompressible
63
{
64
65
/*---------------------------------------------------------------------------*\
66
Class RASModel Declaration
67
\*---------------------------------------------------------------------------*/
68
69
class
RASModel
70
:
71
public
turbulenceModel
,
72
public
IOdictionary
73
{
74
75
protected
:
76
77
// Protected data
78
79
//- Turbulence on/off flag
80
Switch
turbulence_
;
81
82
//- Flag to print the model coeffs at run-time
83
Switch
printCoeffs_
;
84
85
//- Model coefficients dictionary
86
dictionary
coeffDict_
;
87
88
//- Lower limit of k
89
dimensionedScalar
k0_
;
90
91
//- Lower limit of epsilon
92
dimensionedScalar
epsilon0_
;
93
94
//- Small epsilon value used to avoid divide by zero
95
dimensionedScalar
epsilonSmall_
;
96
97
//- Lower limit for omega
98
dimensionedScalar
omega0_
;
99
100
//- Small omega value used to avoid divide by zero
101
dimensionedScalar
omegaSmall_
;
102
103
//- Near wall distance boundary field
104
nearWallDist
y_
;
105
106
107
// Protected member functions
108
109
//- Print model coefficients
110
virtual
void
printCoeffs
();
111
112
113
private
:
114
115
// Private Member Functions
116
117
//- Disallow default bitwise copy construct
118
RASModel
(
const
RASModel
&);
119
120
//- Disallow default bitwise assignment
121
void
operator=
(
const
RASModel
&);
122
123
124
public
:
125
126
//- Runtime type information
127
TypeName
(
"RASModel"
);
128
129
130
// Declare run-time constructor selection table
131
132
declareRunTimeSelectionTable
133
(
134
autoPtr
,
135
RASModel
,
136
dictionary
,
137
(
138
const
volVectorField
&
U
,
139
const
surfaceScalarField
&
phi
,
140
transportModel
& lamTransportModel
141
),
142
(U, phi, lamTransportModel)
143
);
144
145
146
// Constructors
147
148
//- Construct from components
149
RASModel
150
(
151
const
word
&
type
,
152
const
volVectorField
& U,
153
const
surfaceScalarField
& phi,
154
transportModel
& lamTransportModel
155
);
156
157
158
// Selectors
159
160
//- Return a reference to the selected RAS model
161
static
autoPtr<RASModel>
New
162
(
163
const
volVectorField
& U,
164
const
surfaceScalarField
& phi,
165
transportModel
& lamTransportModel
166
);
167
168
169
//- Destructor
170
virtual
~RASModel
()
171
{}
172
173
174
// Member Functions
175
176
// Access
177
178
//- Return the value of k0 which k is not allowed to be less than
179
const
dimensionedScalar
&
k0
()
const
180
{
181
return
k0_
;
182
}
183
184
//- Return the value of epsilon0 which epsilon is not allowed to be
185
// less than
186
const
dimensionedScalar
&
epsilon0
()
const
187
{
188
return
epsilon0_
;
189
}
190
191
//- Return the value of epsilonSmall which is added to epsilon when
192
// calculating nut
193
const
dimensionedScalar
&
epsilonSmall
()
const
194
{
195
return
epsilonSmall_
;
196
}
197
198
//- Return the value of omega0 which epsilon is not allowed to be
199
// less than
200
const
dimensionedScalar
&
omega0
()
const
201
{
202
return
omega0_
;
203
}
204
205
//- Return the value of omegaSmall which is added to epsilon when
206
// calculating nut
207
const
dimensionedScalar
&
omegaSmall
()
const
208
{
209
return
omegaSmall_
;
210
}
211
212
//- Allow k0 to be changed
213
dimensionedScalar
&
k0
()
214
{
215
return
k0_
;
216
}
217
218
//- Allow epsilon0 to be changed
219
dimensionedScalar
&
epsilon0
()
220
{
221
return
epsilon0_
;
222
}
223
224
//- Allow epsilonSmall to be changed
225
dimensionedScalar
&
epsilonSmall
()
226
{
227
return
epsilonSmall_
;
228
}
229
230
//- Allow omega0 to be changed
231
dimensionedScalar
&
omega0
()
232
{
233
return
omega0_
;
234
}
235
236
//- Allow omegaSmall to be changed
237
dimensionedScalar
&
omegaSmall
()
238
{
239
return
omegaSmall_
;
240
}
241
242
//- Return the near wall distances
243
const
nearWallDist
&
y
()
const
244
{
245
return
y_
;
246
}
247
248
//- Calculate y+ at the edge of the laminar sublayer
249
scalar yPlusLam(
const
scalar
kappa
,
const
scalar
E
)
const
;
250
251
//- Const access to the coefficients dictionary
252
const
dictionary
&
coeffDict
()
const
253
{
254
return
coeffDict_
;
255
}
256
257
258
//- Return the turbulence viscosity
259
virtual
tmp<volScalarField>
nut()
const
= 0;
260
261
//- Return the effective viscosity
262
virtual
tmp<volScalarField>
nuEff
()
const
263
{
264
return
tmp<volScalarField>
265
(
266
new
volScalarField
(
"nuEff"
,
nut
() +
nu
())
267
);
268
}
269
270
//- Return the turbulence kinetic energy
271
virtual
tmp<volScalarField>
k
()
const
= 0;
272
273
//- Return the turbulence kinetic energy dissipation rate
274
virtual
tmp<volScalarField>
epsilon()
const
= 0;
275
276
//- Return the Reynolds stress tensor
277
virtual
tmp<volSymmTensorField>
R
()
const
= 0;
278
279
//- Return the effective stress tensor including the laminar stress
280
virtual
tmp<volSymmTensorField>
devReff()
const
= 0;
281
282
//- Return the source term for the momentum equation
283
virtual
tmp<fvVectorMatrix>
divDevReff(
volVectorField
& U)
const
= 0;
284
285
//- Return yPlus for the given patch
286
virtual
tmp<scalarField>
yPlus
287
(
288
const
label patchI,
289
const
scalar
Cmu
290
)
const
;
291
292
//- Solve the turbulence equations and correct the turbulence viscosity
293
virtual
void
correct
() = 0;
294
295
//- Read RASProperties dictionary
296
virtual
bool
read() = 0;
297
};
298
299
300
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
301
302
}
// End namespace incompressible
303
}
// End namespace Foam
304
305
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
306
307
#endif
308
309
// ************************ vim: set sw=4 sts=4 et: ************************ //