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
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::compressible::LESModels
26
27
Description
28
Namespace for compressible LES models.
29
30
31
Class
32
Foam::compressible::LESModel
33
34
Description
35
Base class for all compressible flow LES SGS models.
36
37
This class defines the basic interface for a compressible flow SGS
38
model, and encapsulates data of value to all possible models.
39
In particular this includes references to all the dependent fields
40
(rho, U, phi), the physical viscosity mu, and the LESProperties
41
dictionary, which contains the model selection and model coefficients.
42
43
SourceFiles
44
LESModel.C
45
46
\*---------------------------------------------------------------------------*/
47
48
#ifndef compressibleLESModel_H
49
#define compressibleLESModel_H
50
51
#include <
compressibleTurbulenceModel/turbulenceModel.H
>
52
#include <
LESdeltas/LESdelta.H
>
53
#include <
finiteVolume/fvm.H
>
54
#include <
finiteVolume/fvc.H
>
55
#include <
finiteVolume/fvMatrices.H
>
56
#include <
basicThermophysicalModels/basicThermo.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
compressible
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
volScalarField
&
rho
,
122
const
volVectorField
&
U
,
123
const
surfaceScalarField
&
phi
,
124
const
basicThermo
& thermoPhysicalModel
125
),
126
(rho, U, phi, thermoPhysicalModel)
127
);
128
129
130
// Constructors
131
132
//- Construct from components
133
LESModel
134
(
135
const
word
&
type
,
136
const
volScalarField
& rho,
137
const
volVectorField
& U,
138
const
surfaceScalarField
& phi,
139
const
basicThermo
& thermoPhysicalModel
140
);
141
142
143
// Selectors
144
145
//- Return a reference to the selected LES model
146
static
autoPtr<LESModel>
New
147
(
148
const
volScalarField
& rho,
149
const
volVectorField
& U,
150
const
surfaceScalarField
& phi,
151
const
basicThermo
& thermoPhysicalModel
152
);
153
154
155
//- Destructor
156
virtual
~LESModel
()
157
{}
158
159
160
// Member Functions
161
162
// Access
163
164
//- Const access to the coefficients dictionary,
165
// which provides info. about choice of models,
166
// and all related data (particularly model coefficients).
167
inline
const
dictionary
&
coeffDict
()
const
168
{
169
return
coeffDict_
;
170
}
171
172
//- Return the value of k0 which k is not allowed to be less than
173
const
dimensionedScalar
&
k0
()
const
174
{
175
return
k0_
;
176
}
177
178
//- Allow k0 to be changed
179
dimensionedScalar
&
k0
()
180
{
181
return
k0_
;
182
}
183
184
//- Access function to filter width
185
inline
const
volScalarField
&
delta
()
const
186
{
187
return
delta_
();
188
}
189
190
191
//- Return the SGS turbulent kinetic energy.
192
virtual
tmp<volScalarField>
k
()
const
= 0;
193
194
//- Return the SGS turbulent dissipation.
195
virtual
tmp<volScalarField>
epsilon()
const
= 0;
196
197
//- Return the SGS turbulent viscosity
198
virtual
tmp<volScalarField>
muSgs()
const
= 0;
199
200
//- Return the effective viscosity
201
virtual
tmp<volScalarField>
muEff
()
const
202
{
203
return
tmp<volScalarField>
204
(
205
new
volScalarField
(
"muEff"
,
muSgs
() +
mu
())
206
);
207
}
208
209
//- Return the SGS turbulent thermal diffusivity
210
virtual
tmp<volScalarField>
alphaSgs()
const
= 0;
211
212
//- Return the SGS thermal conductivity.
213
virtual
tmp<volScalarField>
alphaEff()
const
= 0;
214
215
//- Return the sub-grid stress tensor.
216
virtual
tmp<volSymmTensorField>
B()
const
= 0;
217
218
//- Return the deviatoric part of the effective sub-grid
219
// turbulence stress tensor including the laminar stress
220
virtual
tmp<volSymmTensorField>
devRhoBeff()
const
= 0;
221
222
//- Returns div(rho*dev(B)).
223
// This is the additional term due to the filtering of the NSE.
224
virtual
tmp<fvVectorMatrix>
divDevRhoBeff(
volVectorField
& U)
const
= 0;
225
226
227
// RAS compatibility functions for the turbulenceModel base class
228
229
//- Return the turbulence viscosity
230
virtual
tmp<volScalarField>
mut
()
const
231
{
232
return
muSgs
();
233
}
234
235
//- Return the turbulence thermal diffusivity
236
virtual
tmp<volScalarField>
alphat
()
const
237
{
238
return
alphaSgs
();
239
}
240
241
//- Return the Reynolds stress tensor
242
virtual
tmp<volSymmTensorField>
R
()
const
243
{
244
return
B
();
245
}
246
247
//- Return the effective stress tensor including the laminar stress
248
virtual
tmp<volSymmTensorField>
devRhoReff
()
const
249
{
250
return
devRhoBeff
();
251
}
252
253
//- Return the source term for the momentum equation
254
virtual
tmp<fvVectorMatrix>
divDevRhoReff
(
volVectorField
& U)
const
255
{
256
return
divDevRhoBeff
(U);
257
}
258
259
260
//- Correct Eddy-Viscosity and related properties.
261
// This calls correct(const tmp<volTensorField>& gradU) by supplying
262
// gradU calculated locally.
263
void
correct
();
264
265
//- Correct Eddy-Viscosity and related properties
266
virtual
void
correct
(
const
tmp<volTensorField>
& gradU);
267
268
//- Read LESProperties dictionary
269
virtual
bool
read() = 0;
270
};
271
272
273
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
274
275
}
// End namespace compressible
276
}
// End namespace Foam
277
278
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
279
280
#endif
281
282
// ************************ vim: set sw=4 sts=4 et: ************************ //