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
EulerDdtScheme
EulerDdtScheme.C
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
\*---------------------------------------------------------------------------*/
25
26
#include "
EulerDdtScheme.H
"
27
#include <
finiteVolume/surfaceInterpolate.H
>
28
#include <
finiteVolume/fvcDiv.H
>
29
#include <
finiteVolume/fvMatrices.H
>
30
31
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32
33
namespace
Foam
34
{
35
36
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37
38
namespace
fv
39
{
40
41
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42
43
template
<
class
Type>
44
tmp<GeometricField<Type, fvPatchField, volMesh> >
45
EulerDdtScheme<Type>::fvcDdt
46
(
47
const
dimensioned<Type>
& dt
48
)
49
{
50
dimensionedScalar
rDeltaT = 1.0/
mesh
().time().deltaT();
51
52
IOobject
ddtIOobject
53
(
54
"ddt("
+dt.
name
()+
')'
,
55
mesh
().time().timeName(),
56
mesh
()
57
);
58
59
if
(
mesh
().moving())
60
{
61
tmp<GeometricField<Type, fvPatchField, volMesh>
> tdtdt
62
(
63
new
GeometricField<Type, fvPatchField, volMesh>
64
(
65
ddtIOobject,
66
mesh
(),
67
dimensioned<Type>
68
(
69
"0"
,
70
dt.
dimensions
()/
dimTime
,
71
pTraits<Type>::zero
72
)
73
)
74
);
75
76
tdtdt().internalField() =
77
rDeltaT.
value
()*dt.
value
()*(1.0 -
mesh
().V0()/
mesh
().V());
78
79
return
tdtdt;
80
}
81
else
82
{
83
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
84
(
85
new
GeometricField<Type, fvPatchField, volMesh>
86
(
87
ddtIOobject,
88
mesh
(),
89
dimensioned<Type>
90
(
91
"0"
,
92
dt.
dimensions
()/
dimTime
,
93
pTraits<Type>::zero
94
),
95
calculatedFvPatchField<Type>::typeName
96
)
97
);
98
}
99
}
100
101
102
template
<
class
Type>
103
tmp<GeometricField<Type, fvPatchField, volMesh>
>
104
EulerDdtScheme<Type>::fvcDdt
105
(
106
const
GeometricField<Type, fvPatchField, volMesh>
& vf
107
)
108
{
109
dimensionedScalar
rDeltaT = 1.0/
mesh
().time().deltaT();
110
111
IOobject
ddtIOobject
112
(
113
"ddt("
+vf.
name
()+
')'
,
114
mesh
().time().timeName(),
115
mesh
()
116
);
117
118
if
(
mesh
().moving())
119
{
120
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
121
(
122
new
GeometricField<Type, fvPatchField, volMesh>
123
(
124
ddtIOobject,
125
mesh
(),
126
rDeltaT.
dimensions
()*vf.
dimensions
(),
127
rDeltaT.
value
()*
128
(
129
vf.
internalField
()
130
- vf.
oldTime
().
internalField
()*
mesh
().V0()/
mesh
().V()
131
),
132
rDeltaT.
value
()*
133
(
134
vf.
boundaryField
() - vf.
oldTime
().
boundaryField
()
135
)
136
)
137
);
138
}
139
else
140
{
141
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
142
(
143
new
GeometricField<Type, fvPatchField, volMesh>
144
(
145
ddtIOobject,
146
rDeltaT*(vf - vf.
oldTime
())
147
)
148
);
149
}
150
}
151
152
153
template
<
class
Type>
154
tmp<GeometricField<Type, fvPatchField, volMesh>
>
155
EulerDdtScheme<Type>::fvcDdt
156
(
157
const
dimensionedScalar
&
rho
,
158
const
GeometricField<Type, fvPatchField, volMesh>
& vf
159
)
160
{
161
dimensionedScalar
rDeltaT = 1.0/
mesh
().time().deltaT();
162
163
IOobject
ddtIOobject
164
(
165
"ddt("
+rho.
name
()+
','
+vf.
name
()+
')'
,
166
mesh
().time().timeName(),
167
mesh
()
168
);
169
170
if
(
mesh
().moving())
171
{
172
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
173
(
174
new
GeometricField<Type, fvPatchField, volMesh>
175
(
176
ddtIOobject,
177
mesh
(),
178
rDeltaT.
dimensions
()*rho.
dimensions
()*vf.
dimensions
(),
179
rDeltaT.
value
()*rho.
value
()*
180
(
181
vf.
internalField
()
182
- vf.
oldTime
().
internalField
()*
mesh
().V0()/
mesh
().V()
183
),
184
rDeltaT.
value
()*rho.
value
()*
185
(
186
vf.
boundaryField
() - vf.
oldTime
().
boundaryField
()
187
)
188
)
189
);
190
}
191
else
192
{
193
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
194
(
195
new
GeometricField<Type, fvPatchField, volMesh>
196
(
197
ddtIOobject,
198
rDeltaT*rho*(vf - vf.
oldTime
())
199
)
200
);
201
}
202
}
203
204
205
template
<
class
Type>
206
tmp<GeometricField<Type, fvPatchField, volMesh>
>
207
EulerDdtScheme<Type>::fvcDdt
208
(
209
const
volScalarField
& rho,
210
const
GeometricField<Type, fvPatchField, volMesh>
& vf
211
)
212
{
213
dimensionedScalar
rDeltaT = 1.0/
mesh
().time().deltaT();
214
215
IOobject
ddtIOobject
216
(
217
"ddt("
+rho.
name
()+
','
+vf.
name
()+
')'
,
218
mesh
().time().timeName(),
219
mesh
()
220
);
221
222
if
(
mesh
().moving())
223
{
224
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
225
(
226
new
GeometricField<Type, fvPatchField, volMesh>
227
(
228
ddtIOobject,
229
mesh
(),
230
rDeltaT.
dimensions
()*rho.
dimensions
()*vf.
dimensions
(),
231
rDeltaT.
value
()*
232
(
233
rho.
internalField
()*vf.
internalField
()
234
- rho.
oldTime
().
internalField
()
235
*vf.
oldTime
().
internalField
()*
mesh
().V0()/
mesh
().V()
236
),
237
rDeltaT.
value
()*
238
(
239
rho.
boundaryField
()*vf.
boundaryField
()
240
- rho.
oldTime
().
boundaryField
()
241
*vf.
oldTime
().
boundaryField
()
242
)
243
)
244
);
245
}
246
else
247
{
248
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
249
(
250
new
GeometricField<Type, fvPatchField, volMesh>
251
(
252
ddtIOobject,
253
rDeltaT*(rho*vf - rho.
oldTime
()*vf.
oldTime
())
254
)
255
);
256
}
257
}
258
259
260
template
<
class
Type>
261
tmp<fvMatrix<Type>
>
262
EulerDdtScheme<Type>::fvmDdt
263
(
264
GeometricField<Type, fvPatchField, volMesh>
& vf
265
)
266
{
267
tmp<fvMatrix<Type>
> tfvm
268
(
269
new
fvMatrix<Type>
270
(
271
vf,
272
vf.
dimensions
()*
dimVol
/
dimTime
273
)
274
);
275
276
fvMatrix<Type>
& fvm = tfvm();
277
278
scalar rDeltaT = 1.0/
mesh
().time().deltaT().value();
279
280
fvm.
diag
() = rDeltaT*
mesh
().V();
281
282
if
(
mesh
().moving())
283
{
284
fvm.
source
() = rDeltaT*vf.
oldTime
().
internalField
()*
mesh
().V0();
285
}
286
else
287
{
288
fvm.
source
() = rDeltaT*vf.
oldTime
().
internalField
()*
mesh
().V();
289
}
290
291
return
tfvm;
292
}
293
294
295
template
<
class
Type>
296
tmp<fvMatrix<Type>
>
297
EulerDdtScheme<Type>::fvmDdt
298
(
299
const
dimensionedScalar
& rho,
300
GeometricField<Type, fvPatchField, volMesh>
& vf
301
)
302
{
303
tmp<fvMatrix<Type>
> tfvm
304
(
305
new
fvMatrix<Type>
306
(
307
vf,
308
rho.
dimensions
()*vf.
dimensions
()*
dimVol
/
dimTime
309
)
310
);
311
fvMatrix<Type>
& fvm = tfvm();
312
313
scalar rDeltaT = 1.0/
mesh
().time().deltaT().value();
314
315
fvm.
diag
() = rDeltaT*rho.
value
()*
mesh
().V();
316
317
if
(
mesh
().moving())
318
{
319
fvm.
source
() = rDeltaT
320
*rho.
value
()*vf.
oldTime
().
internalField
()*
mesh
().V0();
321
}
322
else
323
{
324
fvm.
source
() = rDeltaT
325
*rho.
value
()*vf.
oldTime
().
internalField
()*
mesh
().V();
326
}
327
328
return
tfvm;
329
}
330
331
332
template
<
class
Type>
333
tmp<fvMatrix<Type>
>
334
EulerDdtScheme<Type>::fvmDdt
335
(
336
const
volScalarField
& rho,
337
GeometricField<Type, fvPatchField, volMesh>
& vf
338
)
339
{
340
tmp<fvMatrix<Type>
> tfvm
341
(
342
new
fvMatrix<Type>
343
(
344
vf,
345
rho.
dimensions
()*vf.
dimensions
()*
dimVol
/
dimTime
346
)
347
);
348
fvMatrix<Type>
& fvm = tfvm();
349
350
scalar rDeltaT = 1.0/
mesh
().time().deltaT().value();
351
352
fvm.
diag
() = rDeltaT*rho.
internalField
()*
mesh
().V();
353
354
if
(
mesh
().moving())
355
{
356
fvm.
source
() = rDeltaT
357
*rho.
oldTime
().
internalField
()
358
*vf.
oldTime
().
internalField
()*
mesh
().V0();
359
}
360
else
361
{
362
fvm.
source
() = rDeltaT
363
*rho.
oldTime
().
internalField
()
364
*vf.
oldTime
().
internalField
()*
mesh
().V();
365
}
366
367
return
tfvm;
368
}
369
370
371
template
<
class
Type>
372
tmp<typename EulerDdtScheme<Type>::fluxFieldType
>
373
EulerDdtScheme<Type>::fvcDdtPhiCorr
374
(
375
const
volScalarField
& rA,
376
const
GeometricField<Type, fvPatchField, volMesh>
&
U
,
377
const
fluxFieldType
& phiAbs
378
)
379
{
380
dimensionedScalar
rDeltaT = 1.0/
mesh
().time().deltaT();
381
382
IOobject
ddtIOobject
383
(
384
"ddtPhiCorr("
+ rA.
name
() +
','
+ U.
name
() +
','
+ phiAbs.
name
() +
')'
,
385
mesh
().time().timeName(),
386
mesh
()
387
);
388
389
tmp<fluxFieldType>
phiCorr =
390
phiAbs.
oldTime
() - (
fvc::interpolate
(U.
oldTime
()) &
mesh
().Sf());
391
392
return
tmp<fluxFieldType>
393
(
394
new
fluxFieldType
395
(
396
ddtIOobject,
397
fvcDdtPhiCoeff(U.
oldTime
(), phiAbs.
oldTime
(), phiCorr())
398
*
fvc::interpolate
(rDeltaT*rA)*phiCorr
399
)
400
);
401
}
402
403
404
template
<
class
Type>
405
tmp<typename EulerDdtScheme<Type>::fluxFieldType
>
406
EulerDdtScheme<Type>::fvcDdtPhiCorr
407
(
408
const
volScalarField
& rA,
409
const
volScalarField
& rho,
410
const
GeometricField<Type, fvPatchField, volMesh>
& U,
411
const
fluxFieldType
& phiAbs
412
)
413
{
414
dimensionedScalar
rDeltaT = 1.0/
mesh
().time().deltaT();
415
416
IOobject
ddtIOobject
417
(
418
"ddtPhiCorr("
419
+ rA.
name
() +
','
420
+ rho.
name
() +
','
421
+ U.
name
() +
','
422
+ phiAbs.
name
() +
')'
,
423
mesh
().time().timeName(),
424
mesh
()
425
);
426
427
if
428
(
429
U.
dimensions
() ==
dimVelocity
430
&& phiAbs.
dimensions
() ==
dimVelocity
*
dimArea
431
)
432
{
433
return
tmp<fluxFieldType>
434
(
435
new
fluxFieldType
436
(
437
ddtIOobject,
438
rDeltaT
439
*fvcDdtPhiCoeff(U.
oldTime
(), phiAbs.
oldTime
())
440
*(
441
fvc::interpolate
(rA*rho.
oldTime
())*phiAbs.
oldTime
()
442
- (
fvc::interpolate
(rA*rho.
oldTime
()*U.
oldTime
())
443
&
mesh
().Sf())
444
)
445
)
446
);
447
}
448
else
if
449
(
450
U.
dimensions
() ==
dimVelocity
451
&& phiAbs.
dimensions
() == rho.
dimensions
()*
dimVelocity
*
dimArea
452
)
453
{
454
return
tmp<fluxFieldType>
455
(
456
new
fluxFieldType
457
(
458
ddtIOobject,
459
rDeltaT
460
*fvcDdtPhiCoeff
461
(
462
U.
oldTime
(),
463
phiAbs.
oldTime
()/
fvc::interpolate
(rho.
oldTime
())
464
)
465
*(
466
fvc::interpolate
(rA*rho.
oldTime
())
467
*phiAbs.
oldTime
()/
fvc::interpolate
(rho.
oldTime
())
468
- (
469
fvc::interpolate
470
(
471
rA*rho.
oldTime
()*U.
oldTime
()
472
) &
mesh
().Sf()
473
)
474
)
475
)
476
);
477
}
478
else
if
479
(
480
U.
dimensions
() == rho.
dimensions
()*
dimVelocity
481
&& phiAbs.
dimensions
() == rho.
dimensions
()*
dimVelocity
*
dimArea
482
)
483
{
484
return
tmp<fluxFieldType>
485
(
486
new
fluxFieldType
487
(
488
ddtIOobject,
489
rDeltaT
490
*fvcDdtPhiCoeff(rho.
oldTime
(), U.
oldTime
(), phiAbs.
oldTime
())
491
*(
492
fvc::interpolate
(rA)*phiAbs.
oldTime
()
493
- (
fvc::interpolate
(rA*U.
oldTime
()) &
mesh
().Sf())
494
)
495
)
496
);
497
}
498
else
499
{
500
FatalErrorIn
501
(
502
"EulerDdtScheme<Type>::fvcDdtPhiCorr"
503
) <<
"dimensions of phiAbs are not correct"
504
<<
abort
(
FatalError
);
505
506
return
fluxFieldType::null();
507
}
508
}
509
510
511
template
<
class
Type>
512
tmp<surfaceScalarField>
EulerDdtScheme<Type>::meshPhi
513
(
514
const
GeometricField<Type, fvPatchField, volMesh>
&
515
)
516
{
517
return
mesh
().phi();
518
}
519
520
521
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
522
523
}
// End namespace fv
524
525
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
526
527
}
// End namespace Foam
528
529
// ************************ vim: set sw=4 sts=4 et: ************************ //