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