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
CrankNicholsonDdtScheme
CrankNicholsonDdtScheme.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::CrankNicholsonDdtScheme
26
27
Description
28
Second-oder CrankNicholson implicit ddt using the current and
29
previous time-step fields as well as the previous time-step ddt.
30
31
SourceFiles
32
CrankNicholsonDdtScheme.C
33
34
\*---------------------------------------------------------------------------*/
35
36
#ifndef CrankNicholsonDdtScheme_H
37
#define CrankNicholsonDdtScheme_H
38
39
#include <
finiteVolume/ddtScheme.H
>
40
41
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42
43
namespace
Foam
44
{
45
46
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47
48
namespace
fv
49
{
50
51
/*---------------------------------------------------------------------------*\
52
Class CrankNicholsonDdtScheme Declaration
53
\*---------------------------------------------------------------------------*/
54
55
template
<
class
Type>
56
class
CrankNicholsonDdtScheme
57
:
58
public
fv::ddtScheme
<Type>
59
{
60
// Private Data
61
62
//- Class to store the ddt0 fields on the objectRegistry for use in the
63
// next time-step. The start-time index of the CN scheme is also
64
// stored to help handle the transition from Euler to CN
65
template
<
class
GeoField>
66
class
DDt0Field
67
:
68
public
GeoField
69
{
70
label startTimeIndex_;
71
72
public
:
73
74
//- Constructor from file for restart.
75
DDt0Field
76
(
77
const
IOobject
& io,
78
const
fvMesh
&
mesh
79
);
80
81
//- Constructor from components, initisalised to zero with given
82
// dimensions.
83
DDt0Field
84
(
85
const
IOobject
& io,
86
const
fvMesh
& mesh,
87
const
dimensioned<typename GeoField::value_type>
& dimType
88
);
89
90
//- Return the start-time index
91
label startTimeIndex()
const
;
92
93
//- Cast to the underlying GeoField
94
GeoField& operator()();
95
96
//- Assignment to a GeoField
97
void
operator=(
const
GeoField& gf);
98
};
99
100
101
//- Off-centering coefficient, 1 -> CN, less than one blends with EI
102
scalar ocCoeff_;
103
104
105
// Private Member Functions
106
107
//- Disallow default bitwise copy construct
108
CrankNicholsonDdtScheme
(
const
CrankNicholsonDdtScheme
&);
109
110
//- Disallow default bitwise assignment
111
void
operator=(
const
CrankNicholsonDdtScheme
&);
112
113
template
<
class
GeoField>
114
DDt0Field<GeoField>& ddt0_
115
(
116
const
word
&
name
,
117
const
dimensionSet
& dims
118
);
119
120
//- Check if the ddt0 needs to be evaluated for this time-step
121
template
<
class
GeoField>
122
bool
evaluate(
const
DDt0Field<GeoField>& ddt0)
const
;
123
124
//- Return the coefficient for Euler scheme for the first time-step
125
// for and CN thereafter
126
template
<
class
GeoField>
127
scalar coef_(
const
DDt0Field<GeoField>&)
const
;
128
129
//- Return the old time-step coefficient for Euler scheme for the
130
// second time-step and for CN thereafter
131
template
<
class
GeoField>
132
scalar coef0_(
const
DDt0Field<GeoField>&)
const
;
133
134
//- Return the reciprocal time-step coefficient for Euler for the
135
// first time-step and CN thereafter
136
template
<
class
GeoField>
137
dimensionedScalar
rDtCoef_(
const
DDt0Field<GeoField>&)
const
;
138
139
//- Return the reciprocal old time-step coefficient for Euler for the
140
// second time-step and CN thereafter
141
template
<
class
GeoField>
142
dimensionedScalar
rDtCoef0_(
const
DDt0Field<GeoField>&)
const
;
143
144
//- Return ddt0 multiplied by the off-centreing coefficient
145
template
<
class
GeoField>
146
tmp<GeoField>
offCentre_(
const
GeoField& ddt0)
const
;
147
148
149
public
:
150
151
//- Runtime type information
152
TypeName
(
"CrankNicholson"
);
153
154
155
// Constructors
156
157
//- Construct from mesh
158
CrankNicholsonDdtScheme
(
const
fvMesh
&
mesh
)
159
:
160
ddtScheme
<Type>(mesh),
161
ocCoeff_(1.0)
162
{}
163
164
//- Construct from mesh and Istream
165
CrankNicholsonDdtScheme
(
const
fvMesh
&
mesh
,
Istream
& is)
166
:
167
ddtScheme
<Type>(mesh, is),
168
ocCoeff_(
readScalar
(is))
169
{
170
if
(ocCoeff_ < 0 || ocCoeff_ > 1)
171
{
172
FatalIOErrorIn
173
(
174
"CrankNicholsonDdtScheme(const fvMesh& mesh, Istream& is)"
,
175
is
176
) <<
"coefficient = "
<< ocCoeff_
177
<<
" should be >= 0 and <= 1"
178
<<
exit
(
FatalIOError
);
179
}
180
}
181
182
183
// Member Functions
184
185
//- Return mesh reference
186
const
fvMesh
&
mesh
()
const
187
{
188
return
fv::ddtScheme<Type>::mesh
();
189
}
190
191
tmp<GeometricField<Type, fvPatchField, volMesh>
>
fvcDdt
192
(
193
const
dimensioned<Type>
&
194
);
195
196
tmp<GeometricField<Type, fvPatchField, volMesh>
>
fvcDdt
197
(
198
const
GeometricField<Type, fvPatchField, volMesh>
&
199
);
200
201
tmp<GeometricField<Type, fvPatchField, volMesh>
>
fvcDdt
202
(
203
const
dimensionedScalar
&,
204
const
GeometricField<Type, fvPatchField, volMesh>
&
205
);
206
207
tmp<GeometricField<Type, fvPatchField, volMesh>
>
fvcDdt
208
(
209
const
volScalarField
&,
210
const
GeometricField<Type, fvPatchField, volMesh>
&
211
);
212
213
tmp<fvMatrix<Type>
>
fvmDdt
214
(
215
GeometricField<Type, fvPatchField, volMesh>
&
216
);
217
218
tmp<fvMatrix<Type>
>
fvmDdt
219
(
220
const
dimensionedScalar
&,
221
GeometricField<Type, fvPatchField, volMesh>
&
222
);
223
224
tmp<fvMatrix<Type>
>
fvmDdt
225
(
226
const
volScalarField
&,
227
GeometricField<Type, fvPatchField, volMesh>
&
228
);
229
230
typedef
typename
ddtScheme<Type>::fluxFieldType
fluxFieldType
;
231
232
tmp<fluxFieldType>
fvcDdtPhiCorr
233
(
234
const
volScalarField
& rA,
235
const
GeometricField<Type, fvPatchField, volMesh>
&
U
,
236
const
fluxFieldType
&
phi
237
);
238
239
tmp<fluxFieldType>
fvcDdtPhiCorr
240
(
241
const
volScalarField
& rA,
242
const
volScalarField
&
rho
,
243
const
GeometricField<Type, fvPatchField, volMesh>
&
U
,
244
const
fluxFieldType
&
phi
245
);
246
247
248
tmp<surfaceScalarField>
meshPhi
249
(
250
const
GeometricField<Type, fvPatchField, volMesh>
&
251
);
252
};
253
254
255
template
<>
256
tmp<surfaceScalarField>
CrankNicholsonDdtScheme<scalar>::fvcDdtPhiCorr
257
(
258
const
volScalarField
& rA,
259
const
volScalarField
&
U
,
260
const
surfaceScalarField
&
phi
261
);
262
263
264
template
<>
265
tmp<surfaceScalarField>
CrankNicholsonDdtScheme<scalar>::fvcDdtPhiCorr
266
(
267
const
volScalarField
& rA,
268
const
volScalarField
&
rho
,
269
const
volScalarField
&
U
,
270
const
surfaceScalarField
&
phi
271
);
272
273
274
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
275
276
}
// End namespace fv
277
278
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
279
280
}
// End namespace Foam
281
282
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
283
284
#ifdef NoRepository
285
# include <
finiteVolume/CrankNicholsonDdtScheme.C
>
286
#endif
287
288
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
289
290
#endif
291
292
// ************************ vim: set sw=4 sts=4 et: ************************ //