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