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
locDynOneEqEddy
locDynOneEqEddy.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::locDynOneEqEddy
26
27
Description
28
Localised Dynamic One Equation Eddy Viscosity Model for incompressible
29
flows
30
31
Eddy viscosity SGS model using a modeled balance equation to simulate
32
the behaviour of k, hence
33
@verbatim
34
d/dt(k) + div(U*k) - div(nuSgs*grad(k))
35
=
36
-B*L - ce*rho*k^3/2/delta
37
and
38
B = 2/3*k*I - 2*nuSgs*dev(D)
39
Beff = 2/3*k*I - 2*nuEff*dev(D)
40
where
41
nuSgs = cD*delta^2*||D||
42
nuEff = nuSgs + nu
43
@endverbatim
44
45
A dynamic procedure is here applied to evaluate ck and ce
46
@verbatim
47
ck=<L.M>/<M.M>
48
and
49
ce=<e*m>/<m*m>
50
where
51
K = 0.5*(F(U.U) - F(U).F(U))
52
L = (F(U*U) - F(U)*F(U) - 0.33*K*I)
53
M = delta*(F(sqrt(k)*D) - 2*sqrt(K + filter(k))*F(D))
54
m = pow(K + F(k), 3.0/2.0)/(2*delta) - F(pow(k, 3.0/2.0))/delta
55
e = 2*delta*ck*(F(sqrt(k)*(D && D)) - 2*sqrt(K + F(k))*(F(D) && F(D)))/
56
@endverbatim
57
58
SourceFiles
59
locDynOneEqEddy.C
60
61
\*---------------------------------------------------------------------------*/
62
63
#ifndef locDynOneEqEddy_H
64
#define locDynOneEqEddy_H
65
66
#include <
incompressibleLESModels/GenEddyVisc.H
>
67
#include <
LESfilters/simpleFilter.H
>
68
#include <
LESfilters/LESfilter.H
>
69
70
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
71
72
namespace
Foam
73
{
74
namespace
incompressible
75
{
76
namespace
LESModels
77
{
78
79
/*---------------------------------------------------------------------------*\
80
Class locDynOneEqEddy Declaration
81
\*---------------------------------------------------------------------------*/
82
83
class
locDynOneEqEddy
84
:
85
public
GenEddyVisc
86
{
87
// Private data
88
89
volScalarField
k_;
90
91
simpleFilter
simpleFilter_;
92
autoPtr<LESfilter>
filterPtr_;
93
LESfilter
& filter_;
94
95
96
// Private Member Functions
97
98
//- Update sub-grid scale fields
99
void
updateSubGridScaleFields
100
(
101
const
volSymmTensorField
& D,
102
const
volScalarField
& KK
103
);
104
105
//- Calculate ck, ce by filtering the velocity field U.
106
volScalarField
ck
107
(
108
const
volSymmTensorField
&,
109
const
volScalarField
&
110
)
const
;
111
112
volScalarField
ce
113
(
114
const
volSymmTensorField
&,
115
const
volScalarField
&
116
)
const
;
117
118
// Disallow default bitwise copy construct and assignment
119
locDynOneEqEddy
(
const
locDynOneEqEddy
&);
120
locDynOneEqEddy
& operator=(
const
locDynOneEqEddy
&);
121
122
123
public
:
124
125
//- Runtime type information
126
TypeName
(
"locDynOneEqEddy"
);
127
128
// Constructors
129
130
//- Construct from components
131
locDynOneEqEddy
132
(
133
const
volVectorField
&
U
,
134
const
surfaceScalarField
&
phi
,
135
transportModel
&
transport
136
);
137
138
139
//- Destructor
140
virtual
~locDynOneEqEddy
()
141
{}
142
143
144
// Member Functions
145
146
//- Return SGS kinetic energy
147
virtual
tmp<volScalarField>
k
()
const
148
{
149
return
k_;
150
}
151
152
//- Return the effective diffusivity for k
153
tmp<volScalarField>
DkEff
()
const
154
{
155
return
tmp<volScalarField>
156
(
157
new
volScalarField
(
"DkEff"
,
nuSgs_
+
nu
())
158
);
159
}
160
161
//- Correct Eddy-Viscosity and related properties
162
virtual
void
correct
(
const
tmp<volTensorField>
& gradU);
163
164
//- Read LESProperties dictionary
165
virtual
bool
read
();
166
};
167
168
169
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
170
171
}
// End namespace LESModels
172
}
// End namespace incompressible
173
}
// End namespace Foam
174
175
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
176
177
#endif
178
179
// ************************ vim: set sw=4 sts=4 et: ************************ //