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
oneEqEddy
oneEqEddy.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::LESModels::oneEqEddy
26
27
Description
28
One Equation Eddy Viscosity Model for incompressible flows
29
30
Eddy viscosity SGS model using a modeled balance equation to simulate the
31
behaviour of k, hence,
32
@verbatim
33
d/dt(k) + div(U*k) - div(nuEff*grad(k))
34
=
35
-B*L - ce*k^3/2/delta
36
37
and
38
39
B = 2/3*k*I - 2*nuSgs*dev(D)
40
Beff = 2/3*k*I - 2*nuEff*dev(D)
41
42
where
43
44
D = symm(grad(U));
45
nuSgs = ck*sqrt(k)*delta
46
nuEff = nuSgs + nu
47
@endverbatim
48
49
SourceFiles
50
oneEqEddy.C
51
52
\*---------------------------------------------------------------------------*/
53
54
#ifndef oneEqEddy_H
55
#define oneEqEddy_H
56
57
#include <
incompressibleLESModels/GenEddyVisc.H
>
58
59
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
60
61
namespace
Foam
62
{
63
namespace
incompressible
64
{
65
namespace
LESModels
66
{
67
68
/*---------------------------------------------------------------------------*\
69
Class oneEqEddy Declaration
70
\*---------------------------------------------------------------------------*/
71
72
class
oneEqEddy
73
:
74
public
GenEddyVisc
75
{
76
// Private data
77
78
volScalarField
k_;
79
80
dimensionedScalar
ck_;
81
82
83
// Private Member Functions
84
85
//- Update sub-grid scale fields
86
void
updateSubGridScaleFields();
87
88
// Disallow default bitwise copy construct and assignment
89
oneEqEddy
(
const
oneEqEddy
&);
90
oneEqEddy
& operator=(
const
oneEqEddy
&);
91
92
93
public
:
94
95
//- Runtime type information
96
TypeName
(
"oneEqEddy"
);
97
98
// Constructors
99
100
//- Construct from components
101
oneEqEddy
102
(
103
const
volVectorField
&
U
,
104
const
surfaceScalarField
&
phi
,
105
transportModel
&
transport
106
);
107
108
109
//- Destructor
110
virtual
~oneEqEddy
()
111
{}
112
113
114
// Member Functions
115
116
//- Return SGS kinetic energy
117
virtual
tmp<volScalarField>
k
()
const
118
{
119
return
k_;
120
}
121
122
//- Return the effective diffusivity for k
123
tmp<volScalarField>
DkEff
()
const
124
{
125
return
tmp<volScalarField>
126
(
127
new
volScalarField
(
"DkEff"
,
nuSgs_
+
nu
())
128
);
129
}
130
131
//- Correct Eddy-Viscosity and related properties
132
virtual
void
correct
(
const
tmp<volTensorField>
& gradU);
133
134
//- Read LESProperties dictionary
135
virtual
bool
read
();
136
};
137
138
139
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
140
141
}
// End namespace LESModels
142
}
// End namespace incompressible
143
}
// End namespace Foam
144
145
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
146
147
#endif
148
149
// ************************ vim: set sw=4 sts=4 et: ************************ //