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
homogeneousDynSmagorinsky
homogeneousDynSmagorinsky.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-2011 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::homogeneousDynSmagorinsky
26
27
Description
28
The Isochoric homogeneous dynamic Smagorinsky Model for
29
incompressible flows.
30
31
Algebraic eddy viscosity SGS model founded on the assumption that
32
local equilibrium prevails.
33
Thus,
34
@verbatim
35
B = 2/3*k*I - 2*nuSgs*dev(D)
36
Beff = 2/3*k*I - 2*nuEff*dev(D)
37
38
where
39
40
k = cI*delta^2*||D||^2
41
nuSgs = cD*delta^2*||D||
42
nuEff = nuSgs + nu
43
44
In the dynamic version of the choric Smagorinsky model
45
the coefficients cI and cD are calculated during the simulation,
46
47
cI=<K*m>/<m*m>
48
49
and
50
51
cD=<L.M>/<M.M>,
52
53
where
54
55
K = 0.5*(F(U.U) - F(U).F(U))
56
m = delta^2*(4*||F(D)||^2 - F(||D||^2))
57
L = dev(F(U*U) - F(U)*F(U))
58
M = delta^2*(F(||D||*dev(D)) - 4*||F(D)||*F(dev(D)))
59
60
The averaging <...> is over the whole domain, i.e. homogeneous turbulence
61
is assumed.
62
@endverbatim
63
64
SourceFiles
65
homogeneousDynSmagorinsky.C
66
67
\*---------------------------------------------------------------------------*/
68
69
#ifndef homogeneousDynSmagorinsky_H
70
#define homogeneousDynSmagorinsky_H
71
72
#include <
incompressibleLESModels/Smagorinsky.H
>
73
#include <
LESfilters/LESfilter.H
>
74
75
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
76
77
namespace
Foam
78
{
79
namespace
incompressible
80
{
81
namespace
LESModels
82
{
83
84
/*---------------------------------------------------------------------------*\
85
Class homogeneousDynSmagorinsky Declaration
86
\*---------------------------------------------------------------------------*/
87
88
class
homogeneousDynSmagorinsky
89
:
90
public
GenEddyVisc
91
{
92
// Private data
93
94
volScalarField
k_;
95
96
autoPtr<LESfilter>
filterPtr_;
97
LESfilter
& filter_;
98
99
100
// Private Member Functions
101
102
//- Update sub-grid scale fields
103
void
updateSubGridScaleFields(
const
volSymmTensorField
& D);
104
105
//- Calculate coefficients cD, cI from filtering velocity field
106
dimensionedScalar
cD(
const
volSymmTensorField
& D)
const
;
107
dimensionedScalar
cI(
const
volSymmTensorField
& D)
const
;
108
109
// Disallow default bitwise copy construct and assignment
110
homogeneousDynSmagorinsky
(
const
homogeneousDynSmagorinsky
&);
111
homogeneousDynSmagorinsky
& operator=(
const
homogeneousDynSmagorinsky
&);
112
113
114
public
:
115
116
//- Runtime type information
117
TypeName
(
"homogeneousDynSmagorinsky"
);
118
119
// Constructors
120
121
//- Construct from components
122
homogeneousDynSmagorinsky
123
(
124
const
volVectorField
&
U
,
125
const
surfaceScalarField
&
phi
,
126
transportModel
&
transport
127
);
128
129
130
//- Destructor
131
virtual
~homogeneousDynSmagorinsky
()
132
{}
133
134
135
// Member Functions
136
137
//- Return SGS kinetic energy
138
virtual
tmp<volScalarField>
k
()
const
139
{
140
return
k_;
141
}
142
143
//- Correct Eddy-Viscosity and related properties
144
virtual
void
correct
(
const
tmp<volTensorField>
& gradU);
145
146
//- Read LESProperties dictionary
147
virtual
bool
read
();
148
};
149
150
151
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
152
153
}
// End namespace LESModels
154
}
// End namespace incompressible
155
}
// End namespace Foam
156
157
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
158
159
#endif
160
161
// ************************ vim: set sw=4 sts=4 et: ************************ //