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
mixedSmagorinsky
mixedSmagorinsky.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::mixedSmagorinsky
26
27
Description
28
The mixed Isochoric Smagorinsky Model for incompressible flows.
29
30
The mixed model is a linear combination of an eddy viscosity model
31
(Smagorinsky) with a scale similarity model. Hence
32
@verbatim
33
B = (L + C) + R = (F(v*v) - F(v)*F(v)) + R
34
@endverbatim
35
36
The algebraic eddy viscosity SGS model is founded on the assumption
37
that local equilibrium prevails, hence
38
@verbatim
39
R = 2/3*k*I - 2*nuEff*dev(D)
40
where
41
k = cI*delta^2*||D||^2
42
nuEff = ck*sqrt(k)*delta + nu
43
@endverbatim
44
45
The Leonard and cross contributions are incorporated
46
by adding,
47
@verbatim
48
+ div(((filter(U*U) - filter(U)*filter(U)) -
49
0.333*I*tr(filter(U*U) - filter(U)*filter(U))))
50
+ div((filter(U*epsilon) - filter(U)*filter(epsilon)))
51
@endverbatim
52
to the rhs. of the equations.
53
54
SourceFiles
55
mixedSmagorinsky.C
56
57
\*---------------------------------------------------------------------------*/
58
59
#ifndef mixedSmagorinsky_H
60
#define mixedSmagorinsky_H
61
62
#include <
incompressibleLESModels/scaleSimilarity.H
>
63
#include <
incompressibleLESModels/Smagorinsky.H
>
64
65
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
66
67
namespace
Foam
68
{
69
namespace
incompressible
70
{
71
namespace
LESModels
72
{
73
74
/*---------------------------------------------------------------------------*\
75
Class mixedSmagorinsky Declaration
76
\*---------------------------------------------------------------------------*/
77
78
class
mixedSmagorinsky
79
:
80
public
scaleSimilarity
,
81
public
Smagorinsky
82
{
83
// Private Member Functions
84
85
// Disallow default bitwise copy construct and assignment
86
mixedSmagorinsky
(
const
mixedSmagorinsky
&);
87
mixedSmagorinsky
& operator=(
const
mixedSmagorinsky
&);
88
89
public
:
90
91
//- Runtime type information
92
TypeName
(
"mixedSmagorinsky"
);
93
94
95
// Constructors
96
97
//- Construct from components
98
mixedSmagorinsky
99
(
100
const
volVectorField
&
U
,
101
const
surfaceScalarField
&
phi
,
102
transportModel
&
transport
103
);
104
105
106
//- Destructor
107
virtual
~mixedSmagorinsky
()
108
{}
109
110
111
// Member Functions
112
113
//- Return the SGS turbulent kinetic energy.
114
virtual
tmp<volScalarField>
k
()
const
;
115
116
//- Return the SGS turbulent disipation rate.
117
virtual
tmp<volScalarField>
epsilon
()
const
;
118
119
//- Return the SGS viscosity.
120
virtual
tmp<volScalarField>
nuSgs
()
const
121
{
122
return
nuSgs_
;
123
}
124
125
//- Return the sub-grid stress tensor.
126
virtual
tmp<volSymmTensorField>
B
()
const
;
127
128
//- Return the effective sub-grid turbulence stress tensor
129
// including the laminar stress
130
virtual
tmp<volSymmTensorField>
devBeff
()
const
;
131
132
//- Implementation of div(B). This is necessary to override
133
// (and include) the div(B) terms from both the parent classes.
134
virtual
tmp<fvVectorMatrix>
divDevBeff
(
volVectorField
& U)
const
;
135
136
//- Correct Eddy-Viscosity and related properties
137
virtual
void
correct
(
const
tmp<volTensorField>
& gradU);
138
139
//- Read LESProperties dictionary
140
virtual
bool
read
();
141
};
142
143
144
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
145
146
}
// End namespace LESModels
147
}
// End namespace incompressible
148
}
// End namespace Foam
149
150
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
151
152
#endif
153
154
// ************************ vim: set sw=4 sts=4 et: ************************ //