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
applications
solvers
combustion
PDRFoam
PDRModels
turbulence
PDRkEpsilon
PDRkEpsilon.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::PDRkEpsilon
26
27
Description
28
Standard k-epsilon turbulence model with additional source terms
29
corresponding to PDR basic drag model (\link basic.H \endlink)
30
31
The default model coefficients correspond to the following:
32
@verbatim
33
kEpsilonCoeffs
34
{
35
Cmu 0.09;
36
C1 1.44;
37
C2 1.92;
38
C3 -0.33; // only for compressible
39
sigmak 1.0; // only for compressible
40
sigmaEps 1.3;
41
Prt 1.0; // only for compressible
42
}
43
@endverbatim
44
45
The turbulence source term \f$ G_{R} \f$ appears in the
46
\f$ \kappa-\epsilon \f$ equation for the generation of turbulence due to
47
interaction with unresolved obstacles.
48
49
In the \f$ \epsilon \f$ equation \f$ C_{1} G_{R} \f$ is added as a source
50
term.
51
52
In the \f$ \kappa \f$ equation \f$ G_{R} \f$ is added as a source term.
53
54
SourceFiles
55
PDRkEpsilon.C
56
PDRkEpsilonCorrect.C
57
58
\*---------------------------------------------------------------------------*/
59
60
#ifndef compressiblePDRkEpsilon_H
61
#define compressiblePDRkEpsilon_H
62
63
#include <
compressibleRASModels/RASModel.H
>
64
65
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
66
67
namespace
Foam
68
{
69
namespace
compressible
70
{
71
namespace
RASModels
72
{
73
74
/*---------------------------------------------------------------------------*\
75
Class PDRkEpsilon Declaration
76
\*---------------------------------------------------------------------------*/
77
78
class
PDRkEpsilon
79
:
80
public
RASModel
81
{
82
// Private data
83
84
// Model coefficients
85
86
dimensionedScalar
Cmu_;
87
dimensionedScalar
C1_;
88
dimensionedScalar
C2_;
89
dimensionedScalar
sigmak_;
90
dimensionedScalar
sigmaEps_;
91
dimensionedScalar
Prt_;
92
93
// Fields
94
95
volScalarField
k_;
96
volScalarField
epsilon_;
97
volScalarField
mut_;
98
volScalarField
alphat_;
99
100
101
public
:
102
103
//- Runtime type information
104
TypeName
(
"PDRkEpsilon"
);
105
106
107
// Constructors
108
109
//- Construct from components
110
PDRkEpsilon
111
(
112
const
volScalarField
&
rho
,
113
const
volVectorField
&
U
,
114
const
surfaceScalarField
&
phi
,
115
const
basicThermo
& thermophysicalModel
116
);
117
118
119
//- Destructor
120
virtual
~
PDRkEpsilon
()
121
{}
122
123
124
// Member Functions
125
126
tmp<volScalarField>
mut
()
const
127
{
128
return
mut_;
129
}
130
131
//- Return the effective diffusivity for k
132
tmp<volScalarField>
DkEff
()
const
133
{
134
return
tmp<volScalarField>
135
(
136
new
volScalarField
(
"DkEff"
, mut_/sigmak_ +
mu
())
137
);
138
}
139
140
//- Return the effective diffusivity for epsilon
141
tmp<volScalarField>
DepsilonEff
()
const
142
{
143
return
tmp<volScalarField>
144
(
145
new
volScalarField
(
"DepsilonEff"
, mut_/sigmaEps_ +
mu
())
146
);
147
}
148
149
//- Return the effective turbulent thermal diffusivity
150
tmp<volScalarField>
alphaEff
()
const
151
{
152
return
tmp<volScalarField>
153
(
154
new
volScalarField
(
"alphaEff"
, alphat_ +
alpha
())
155
);
156
}
157
158
//- Return the turbulence kinetic energy
159
tmp<volScalarField>
k
()
const
160
{
161
return
k_;
162
}
163
164
//- Return the turbulence kinetic energy dissipation rate
165
tmp<volScalarField>
epsilon
()
const
166
{
167
return
epsilon_;
168
}
169
170
//- Return the Reynolds stress tensor
171
tmp<volSymmTensorField>
R()
const
;
172
173
//- Return the effective stress tensor including the laminar stress
174
tmp<volSymmTensorField>
devRhoReff()
const
;
175
176
//- Return the source term for the momentum equation
177
tmp<fvVectorMatrix>
divDevRhoReff(
volVectorField
& U)
const
;
178
179
//- Solve the turbulence equations and correct the turbulence viscosity
180
void
correct
();
181
182
//- Read turbulenceProperties dictionary
183
bool
read();
184
};
185
186
187
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
188
189
}
// End namespace RASModels
190
}
// End namespace compressible
191
}
// End namespace Foam
192
193
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
194
195
#endif
196
197
// ************************ vim: set sw=4 sts=4 et: ************************ //