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
LienCubicKE
LienCubicKE.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::incompressible::RASModels::LienCubicKE
26
27
Description
28
Lien cubic non-linear k-epsilon turbulence model for incompressible flows.
29
30
SourceFiles
31
LienCubicKE.C
32
33
\*---------------------------------------------------------------------------*/
34
35
#ifndef LienCubicKE_H
36
#define LienCubicKE_H
37
38
#include <
incompressibleRASModels/RASModel.H
>
39
40
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41
42
namespace
Foam
43
{
44
namespace
incompressible
45
{
46
namespace
RASModels
47
{
48
49
/*---------------------------------------------------------------------------*\
50
Class LienCubicKE Declaration
51
\*---------------------------------------------------------------------------*/
52
53
class
LienCubicKE
54
:
55
public
RASModel
56
{
57
// Private data
58
59
// Model coefficients
60
61
dimensionedScalar
C1_;
62
dimensionedScalar
C2_;
63
dimensionedScalar
sigmak_;
64
dimensionedScalar
sigmaEps_;
65
dimensionedScalar
A1_;
66
dimensionedScalar
A2_;
67
dimensionedScalar
Ctau1_;
68
dimensionedScalar
Ctau2_;
69
dimensionedScalar
Ctau3_;
70
dimensionedScalar
alphaKsi_;
71
72
73
// Fields
74
75
volScalarField
k_;
76
volScalarField
epsilon_;
77
78
volTensorField
gradU_;
79
volScalarField
eta_;
80
volScalarField
ksi_;
81
volScalarField
Cmu_;
82
volScalarField
fEta_;
83
volScalarField
C5viscosity_;
84
85
volScalarField
nut_;
86
87
volSymmTensorField
nonlinearStress_;
88
89
90
public
:
91
92
//- Runtime type information
93
TypeName
(
"LienCubicKE"
);
94
95
// Constructors
96
97
//- Construct from components
98
LienCubicKE
99
(
100
const
volVectorField
&
U
,
101
const
surfaceScalarField
&
phi
,
102
transportModel
&
transport
103
);
104
105
106
// Destructor
107
virtual
~LienCubicKE
()
108
{}
109
110
111
// Member Functions
112
113
//- Return the turbulence viscosity
114
virtual
tmp<volScalarField>
nut
()
const
115
{
116
return
nut_;
117
}
118
119
//- Return the effective diffusivity for k
120
tmp<volScalarField>
DkEff
()
const
121
{
122
return
tmp<volScalarField>
123
(
124
new
volScalarField
(
"DkEff"
, nut_/sigmak_ +
nu
())
125
);
126
}
127
128
//- Return the effective diffusivity for epsilon
129
tmp<volScalarField>
DepsilonEff
()
const
130
{
131
return
tmp<volScalarField>
132
(
133
new
volScalarField
(
"DepsilonEff"
, nut_/sigmaEps_ +
nu
())
134
);
135
}
136
137
//- Return the turbulence kinetic energy
138
virtual
tmp<volScalarField>
k
()
const
139
{
140
return
k_;
141
}
142
143
//- Return the turbulence kinetic energy dissipation rate
144
virtual
tmp<volScalarField>
epsilon
()
const
145
{
146
return
epsilon_;
147
}
148
149
//- Return the Reynolds stress tensor
150
virtual
tmp<volSymmTensorField>
R
()
const
;
151
152
//- Return the effective stress tensor including the laminar stress
153
virtual
tmp<volSymmTensorField>
devReff
()
const
;
154
155
//- Return the source term for the momentum equation
156
virtual
tmp<fvVectorMatrix>
divDevReff
(
volVectorField
& U)
const
;
157
158
//- Solve the turbulence equations and correct the turbulence viscosity
159
virtual
void
correct
();
160
161
//- Read RASProperties dictionary
162
virtual
bool
read
();
163
};
164
165
166
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
167
168
}
// End namespace RASModels
169
}
// End namespace incompressible
170
}
// End namespace Foam
171
172
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
173
174
#endif
175
176
// ************************ vim: set sw=4 sts=4 et: ************************ //