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