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
LESModel
LESModel.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::LESModels
26
27
Description
28
Namespace for incompressible LES models.
29
30
Class
31
Foam::incompressible::LESModel
32
33
Description
34
Base class for all incompressible flow LES SGS models.
35
36
This class defines the basic interface for an incompressible flow SGS
37
model, and encapsulates data of value to all possible models.
38
In particular this includes references to all the dependent fields
39
(U, phi), the physical viscosity nu, and the LESProperties
40
dictionary, which contains the model selection and model coefficients.
41
42
SourceFiles
43
LESModel.C
44
45
\*---------------------------------------------------------------------------*/
46
47
#ifndef LESModel_H
48
#define LESModel_H
49
50
#include <
incompressibleTurbulenceModel/turbulenceModel.H
>
51
#include <
LESdeltas/LESdelta.H
>
52
#include <
finiteVolume/fvm.H
>
53
#include <
finiteVolume/fvc.H
>
54
#include <
finiteVolume/fvMatrices.H
>
55
#include <
incompressibleTransportModels/transportModel.H
>
56
#include <
finiteVolume/wallFvPatch.H
>
57
#include <
finiteVolume/bound.H
>
58
#include <
OpenFOAM/autoPtr.H
>
59
#include <
OpenFOAM/runTimeSelectionTables.H
>
60
61
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
62
63
namespace
Foam
64
{
65
namespace
incompressible
66
{
67
68
/*---------------------------------------------------------------------------*\
69
Class LESModel Declaration
70
\*---------------------------------------------------------------------------*/
71
72
class
LESModel
73
:
74
public
turbulenceModel
,
75
public
IOdictionary
76
{
77
78
protected
:
79
80
// Protected data
81
82
Switch
printCoeffs_
;
83
dictionary
coeffDict_
;
84
85
dimensionedScalar
k0_
;
86
87
autoPtr<LESdelta>
delta_
;
88
89
90
// Protected member functions
91
92
//- Print model coefficients
93
virtual
void
printCoeffs
();
94
95
96
private
:
97
98
// Private Member Functions
99
100
//- Disallow default bitwise copy construct
101
LESModel
(
const
LESModel
&);
102
103
//- Disallow default bitwise assignment
104
LESModel
&
operator=
(
const
LESModel
&);
105
106
107
public
:
108
109
//- Runtime type information
110
TypeName
(
"LESModel"
);
111
112
113
// Declare run-time constructor selection table
114
115
declareRunTimeSelectionTable
116
(
117
autoPtr
,
118
LESModel
,
119
dictionary
,
120
(
121
const
volVectorField
&
U
,
122
const
surfaceScalarField
&
phi
,
123
transportModel
& lamTransportModel
124
),
125
(U, phi, lamTransportModel)
126
);
127
128
129
// Constructors
130
131
//- Construct from components
132
LESModel
133
(
134
const
word
&
type
,
135
const
volVectorField
& U,
136
const
surfaceScalarField
& phi,
137
transportModel
& lamTransportModel
138
);
139
140
141
// Selectors
142
143
//- Return a reference to the selected LES model
144
static
autoPtr<LESModel>
New
145
(
146
const
volVectorField
& U,
147
const
surfaceScalarField
& phi,
148
transportModel
& lamTransportModel
149
);
150
151
152
//- Destructor
153
virtual
~LESModel
()
154
{}
155
156
157
// Member Functions
158
159
//- Const access to the coefficients dictionary,
160
// which provides info. about choice of models,
161
// and all related data (particularly model coefficients).
162
inline
const
dictionary
&
coeffDict
()
const
163
{
164
return
coeffDict_
;
165
}
166
167
//- Access function to filter width
168
virtual
const
volScalarField
&
delta
()
const
169
{
170
return
delta_
();
171
}
172
173
//- Return the value of k0 which k is not allowed to be less than
174
const
dimensionedScalar
&
k0
()
const
175
{
176
return
k0_
;
177
}
178
179
//- Allow k0 to be changed
180
dimensionedScalar
&
k0
()
181
{
182
return
k0_
;
183
}
184
185
186
//- Return the SGS turbulent kinetic energy.
187
virtual
tmp<volScalarField>
k
()
const
= 0;
188
189
//- Return the SGS turbulent dissipation.
190
virtual
tmp<volScalarField>
epsilon()
const
= 0;
191
192
//- Return the SGS viscosity.
193
virtual
tmp<volScalarField>
nuSgs()
const
= 0;
194
195
//- Return the effective viscosity
196
virtual
tmp<volScalarField>
nuEff
()
const
197
{
198
return
tmp<volScalarField>
199
(
200
new
volScalarField
(
"nuEff"
,
nuSgs
() +
nu
())
201
);
202
}
203
204
//- Return the sub-grid stress tensor.
205
virtual
tmp<volSymmTensorField>
B()
const
= 0;
206
207
//- Return the deviatoric part of the effective sub-grid
208
// turbulence stress tensor including the laminar stress
209
virtual
tmp<volSymmTensorField>
devBeff()
const
= 0;
210
211
//- Returns div(dev(Beff)).
212
// This is the additional term due to the filtering of the NSE.
213
virtual
tmp<fvVectorMatrix>
divDevBeff(
volVectorField
& U)
const
= 0;
214
215
216
// RAS compatibility functions for the turbulenceModel base class
217
218
//- Return the turbulence viscosity
219
virtual
tmp<volScalarField>
nut
()
const
220
{
221
return
nuSgs
();
222
}
223
224
//- Return the Reynolds stress tensor
225
virtual
tmp<volSymmTensorField>
R
()
const
226
{
227
return
B
();
228
}
229
230
//- Return the effective stress tensor including the laminar stress
231
virtual
tmp<volSymmTensorField>
devReff
()
const
232
{
233
return
devBeff
();
234
}
235
236
//- Return the source term for the momentum equation
237
virtual
tmp<fvVectorMatrix>
divDevReff
(
volVectorField
& U)
const
238
{
239
return
divDevBeff
(U);
240
}
241
242
243
//- Correct Eddy-Viscosity and related properties.
244
// This calls correct(const tmp<volTensorField>& gradU) by supplying
245
// gradU calculated locally.
246
void
correct
();
247
248
//- Correct Eddy-Viscosity and related properties
249
virtual
void
correct
(
const
tmp<volTensorField>
& gradU);
250
251
//- Read LESProperties dictionary
252
virtual
bool
read() = 0;
253
};
254
255
256
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
257
258
}
// End namespace incompressible
259
}
// End namespace Foam
260
261
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
262
263
#endif
264
265
// ************************ vim: set sw=4 sts=4 et: ************************ //