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
finiteVolume
finiteVolume
ddtSchemes
SLTSDdtScheme
SLTSDdtScheme.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::fv::SLTSDdtScheme
26
27
Description
28
Stabilised local time-step first-order Euler implicit/explicit ddt.
29
The time-step is adjusted locally so that an advective equations remains
30
diagonally dominant.
31
32
This scheme should only be used for steady-state computations
33
using transient codes where local time-stepping is preferably to
34
under-relaxation for transport consistency reasons.
35
36
See also CoEulerDdtScheme.
37
38
SourceFiles
39
SLTSDdtScheme.C
40
41
\*---------------------------------------------------------------------------*/
42
43
#ifndef SLTSDdtScheme_H
44
#define SLTSDdtScheme_H
45
46
#include <
finiteVolume/ddtScheme.H
>
47
48
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
49
50
namespace
Foam
51
{
52
53
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
54
55
namespace
fv
56
{
57
58
/*---------------------------------------------------------------------------*\
59
Class SLTSDdtScheme Declaration
60
\*---------------------------------------------------------------------------*/
61
62
template
<
class
Type>
63
class
SLTSDdtScheme
64
:
65
public
fv::ddtScheme
<Type>
66
{
67
// Private Data
68
69
//- Name of the flux field used to calculate the local time-step
70
word
phiName_;
71
72
//- Name of the density field used to obtain the volumetric flux
73
// from the mass flux if required
74
word
rhoName_;
75
76
//- Under-relaxation factor
77
scalar alpha_;
78
79
80
// Private Member Functions
81
82
//- Disallow default bitwise copy construct
83
SLTSDdtScheme
(
const
SLTSDdtScheme
&);
84
85
//- Disallow default bitwise assignment
86
void
operator=(
const
SLTSDdtScheme
&);
87
88
//- Calculate a relaxed diagonal from the given flux field
89
void
relaxedDiag(
scalarField
& rD,
const
surfaceScalarField
&
phi
)
const
;
90
91
//- Return the reciprocal of the stabilised local time-step
92
tmp<volScalarField>
SLrDeltaT()
const
;
93
94
95
public
:
96
97
//- Runtime type information
98
TypeName
(
"SLTS"
);
99
100
101
// Constructors
102
103
//- Construct from mesh and Istream
104
SLTSDdtScheme
(
const
fvMesh
&
mesh
,
Istream
& is)
105
:
106
ddtScheme
<Type>(mesh, is),
107
phiName_(is),
108
rhoName_(is),
109
alpha_(
readScalar
(is))
110
{}
111
112
113
// Member Functions
114
115
//- Return mesh reference
116
const
fvMesh
&
mesh
()
const
117
{
118
return
fv::ddtScheme<Type>::mesh
();
119
}
120
121
tmp<GeometricField<Type, fvPatchField, volMesh>
>
fvcDdt
122
(
123
const
dimensioned<Type>
&
124
);
125
126
tmp<GeometricField<Type, fvPatchField, volMesh>
>
fvcDdt
127
(
128
const
GeometricField<Type, fvPatchField, volMesh>
&
129
);
130
131
tmp<GeometricField<Type, fvPatchField, volMesh>
>
fvcDdt
132
(
133
const
dimensionedScalar
&,
134
const
GeometricField<Type, fvPatchField, volMesh>
&
135
);
136
137
tmp<GeometricField<Type, fvPatchField, volMesh>
>
fvcDdt
138
(
139
const
volScalarField
&,
140
const
GeometricField<Type, fvPatchField, volMesh>
&
141
);
142
143
tmp<fvMatrix<Type>
>
fvmDdt
144
(
145
GeometricField<Type, fvPatchField, volMesh>
&
146
);
147
148
tmp<fvMatrix<Type>
>
fvmDdt
149
(
150
const
dimensionedScalar
&,
151
GeometricField<Type, fvPatchField, volMesh>
&
152
);
153
154
tmp<fvMatrix<Type>
>
fvmDdt
155
(
156
const
volScalarField
&,
157
GeometricField<Type, fvPatchField, volMesh>
&
158
);
159
160
typedef
typename
ddtScheme<Type>::fluxFieldType
fluxFieldType
;
161
162
tmp<fluxFieldType>
fvcDdtPhiCorr
163
(
164
const
volScalarField
& rA,
165
const
GeometricField<Type, fvPatchField, volMesh>
&
U
,
166
const
fluxFieldType
& phi
167
);
168
169
tmp<fluxFieldType>
fvcDdtPhiCorr
170
(
171
const
volScalarField
& rA,
172
const
volScalarField
&
rho
,
173
const
GeometricField<Type, fvPatchField, volMesh>
&
U
,
174
const
fluxFieldType
& phi
175
);
176
177
tmp<surfaceScalarField>
meshPhi
178
(
179
const
GeometricField<Type, fvPatchField, volMesh>
&
180
);
181
};
182
183
184
template
<>
185
tmp<surfaceScalarField>
SLTSDdtScheme<scalar>::fvcDdtPhiCorr
186
(
187
const
volScalarField
& rA,
188
const
volScalarField
&
U
,
189
const
surfaceScalarField
&
phi
190
);
191
192
193
template
<>
194
tmp<surfaceScalarField>
SLTSDdtScheme<scalar>::fvcDdtPhiCorr
195
(
196
const
volScalarField
& rA,
197
const
volScalarField
&
rho
,
198
const
volScalarField
&
U
,
199
const
surfaceScalarField
&
phi
200
);
201
202
203
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
204
205
}
// End namespace fv
206
207
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
208
209
}
// End namespace Foam
210
211
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
212
213
#ifdef NoRepository
214
# include <
finiteVolume/SLTSDdtScheme.C
>
215
#endif
216
217
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
218
219
#endif
220
221
// ************************ vim: set sw=4 sts=4 et: ************************ //