backwardDdtScheme.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2018 OpenFOAM Foundation
9 Copyright (C) 2017 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "backwardDdtScheme.H"
30#include "surfaceInterpolate.H"
31#include "fvcDiv.H"
32#include "fvMatrices.H"
33
34// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35
36namespace Foam
37{
38
39// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40
41namespace fv
42{
43
44// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45
46template<class Type>
47scalar backwardDdtScheme<Type>::deltaT_() const
48{
49 return mesh().time().deltaTValue();
50}
51
52
53template<class Type>
54scalar backwardDdtScheme<Type>::deltaT0_() const
55{
56 return mesh().time().deltaT0Value();
57}
58
59
60template<class Type>
61template<class GeoField>
62scalar backwardDdtScheme<Type>::deltaT0_(const GeoField& vf) const
63{
64 if (mesh().time().timeIndex() < 2)
65 {
66 return GREAT;
67 }
68 else
69 {
70 return deltaT0_();
71 }
72}
73
74
75// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
76
77template<class Type>
78tmp<GeometricField<Type, fvPatchField, volMesh>>
79backwardDdtScheme<Type>::fvcDdt
80(
81 const dimensioned<Type>& dt
82)
83{
84 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
85
86 IOobject ddtIOobject
87 (
88 "ddt("+dt.name()+')',
89 mesh().time().timeName(),
90 mesh()
91 );
92
93 scalar deltaT = deltaT_();
94 scalar deltaT0 = deltaT0_();
95
96 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
97 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
98 scalar coefft0 = coefft + coefft00;
99
100 if (mesh().moving())
101 {
102 tmp<GeometricField<Type, fvPatchField, volMesh>> tdtdt
103 (
104 new GeometricField<Type, fvPatchField, volMesh>
105 (
106 ddtIOobject,
107 mesh(),
108 dimensioned<Type>(dt.dimensions()/dimTime, Zero)
109 )
110 );
111
112 tdtdt.ref().primitiveFieldRef() = rDeltaT.value()*dt.value()*
113 (
114 coefft - (coefft0*mesh().V0() - coefft00*mesh().V00())/mesh().V()
115 );
116
117 return tdtdt;
118 }
119
120 return tmp<GeometricField<Type, fvPatchField, volMesh>>::New
121 (
122 ddtIOobject,
123 mesh(),
124 dimensioned<Type>(dt.dimensions()/dimTime, Zero),
125 calculatedFvPatchField<Type>::typeName
126 );
127}
128
129
130template<class Type>
131tmp<GeometricField<Type, fvPatchField, volMesh>>
132backwardDdtScheme<Type>::fvcDdt
133(
134 const GeometricField<Type, fvPatchField, volMesh>& vf
135)
136{
137 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
138
139 IOobject ddtIOobject
140 (
141 "ddt("+vf.name()+')',
142 mesh().time().timeName(),
143 mesh()
144 );
145
146 scalar deltaT = deltaT_();
147 scalar deltaT0 = deltaT0_(vf);
148
149 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
150 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
151 scalar coefft0 = coefft + coefft00;
152
153 if (mesh().moving())
154 {
155 return tmp<GeometricField<Type, fvPatchField, volMesh>>
156 (
157 new GeometricField<Type, fvPatchField, volMesh>
158 (
159 ddtIOobject,
160 mesh(),
161 rDeltaT.dimensions()*vf.dimensions(),
162 rDeltaT.value()*
163 (
164 coefft*vf.primitiveField() -
165 (
166 coefft0*vf.oldTime().primitiveField()*mesh().V0()
167 - coefft00*vf.oldTime().oldTime().primitiveField()
168 *mesh().V00()
169 )/mesh().V()
170 ),
171 rDeltaT.value()*
172 (
173 coefft*vf.boundaryField() -
174 (
175 coefft0*vf.oldTime().boundaryField()
176 - coefft00*vf.oldTime().oldTime().boundaryField()
177 )
178 )
179 )
180 );
181 }
182 else
183 {
184 return tmp<GeometricField<Type, fvPatchField, volMesh>>
185 (
186 new GeometricField<Type, fvPatchField, volMesh>
187 (
188 ddtIOobject,
189 rDeltaT*
190 (
191 coefft*vf
192 - coefft0*vf.oldTime()
193 + coefft00*vf.oldTime().oldTime()
194 )
195 )
196 );
197 }
198}
199
200
201template<class Type>
202tmp<GeometricField<Type, fvPatchField, volMesh>>
203backwardDdtScheme<Type>::fvcDdt
204(
205 const dimensionedScalar& rho,
206 const GeometricField<Type, fvPatchField, volMesh>& vf
207)
208{
209 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
210
211 IOobject ddtIOobject
212 (
213 "ddt("+rho.name()+','+vf.name()+')',
214 mesh().time().timeName(),
215 mesh()
216 );
217
218 scalar deltaT = deltaT_();
219 scalar deltaT0 = deltaT0_(vf);
220
221 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
222 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
223 scalar coefft0 = coefft + coefft00;
224
225 if (mesh().moving())
226 {
227 return tmp<GeometricField<Type, fvPatchField, volMesh>>
228 (
229 new GeometricField<Type, fvPatchField, volMesh>
230 (
231 ddtIOobject,
232 mesh(),
233 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
234 rDeltaT.value()*rho.value()*
235 (
236 coefft*vf.primitiveField() -
237 (
238 coefft0*vf.oldTime().primitiveField()*mesh().V0()
239 - coefft00*vf.oldTime().oldTime().primitiveField()
240 *mesh().V00()
241 )/mesh().V()
242 ),
243 rDeltaT.value()*rho.value()*
244 (
245 coefft*vf.boundaryField() -
246 (
247 coefft0*vf.oldTime().boundaryField()
248 - coefft00*vf.oldTime().oldTime().boundaryField()
249 )
250 )
251 )
252 );
253 }
254 else
255 {
256 return tmp<GeometricField<Type, fvPatchField, volMesh>>
257 (
258 new GeometricField<Type, fvPatchField, volMesh>
259 (
260 ddtIOobject,
261 rDeltaT*rho*
262 (
263 coefft*vf
264 - coefft0*vf.oldTime()
265 + coefft00*vf.oldTime().oldTime()
266 )
267 )
268 );
269 }
270}
271
272
273template<class Type>
274tmp<GeometricField<Type, fvPatchField, volMesh>>
275backwardDdtScheme<Type>::fvcDdt
276(
277 const volScalarField& rho,
278 const GeometricField<Type, fvPatchField, volMesh>& vf
279)
280{
281 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
282
283 IOobject ddtIOobject
284 (
285 "ddt("+rho.name()+','+vf.name()+')',
286 mesh().time().timeName(),
287 mesh()
288 );
289
290 scalar deltaT = deltaT_();
291 scalar deltaT0 = deltaT0_(vf);
292
293 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
294 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
295 scalar coefft0 = coefft + coefft00;
296
297 if (mesh().moving())
298 {
299 return tmp<GeometricField<Type, fvPatchField, volMesh>>
300 (
301 new GeometricField<Type, fvPatchField, volMesh>
302 (
303 ddtIOobject,
304 mesh(),
305 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
306 rDeltaT.value()*
307 (
308 coefft*rho.primitiveField()*vf.primitiveField() -
309 (
310 coefft0*rho.oldTime().primitiveField()
311 *vf.oldTime().primitiveField()*mesh().V0()
312 - coefft00*rho.oldTime().oldTime().primitiveField()
313 *vf.oldTime().oldTime().primitiveField()*mesh().V00()
314 )/mesh().V()
315 ),
316 rDeltaT.value()*
317 (
318 coefft*rho.boundaryField()*vf.boundaryField() -
319 (
320 coefft0*rho.oldTime().boundaryField()
321 *vf.oldTime().boundaryField()
322 - coefft00*rho.oldTime().oldTime().boundaryField()
323 *vf.oldTime().oldTime().boundaryField()
324 )
325 )
326 )
327 );
328 }
329 else
330 {
331 return tmp<GeometricField<Type, fvPatchField, volMesh>>
332 (
333 new GeometricField<Type, fvPatchField, volMesh>
334 (
335 ddtIOobject,
336 rDeltaT*
337 (
338 coefft*rho*vf
339 - coefft0*rho.oldTime()*vf.oldTime()
340 + coefft00*rho.oldTime().oldTime()*vf.oldTime().oldTime()
341 )
342 )
343 );
344 }
345}
346
347
348template<class Type>
349tmp<GeometricField<Type, fvPatchField, volMesh>>
350backwardDdtScheme<Type>::fvcDdt
351(
352 const volScalarField& alpha,
353 const volScalarField& rho,
354 const GeometricField<Type, fvPatchField, volMesh>& vf
355)
356{
357 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
358
359 IOobject ddtIOobject
360 (
361 "ddt("+alpha.name()+','+rho.name()+','+vf.name()+')',
362 mesh().time().timeName(),
363 mesh()
364 );
365
366 scalar deltaT = deltaT_();
367 scalar deltaT0 = deltaT0_(vf);
368
369 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
370 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
371 scalar coefft0 = coefft + coefft00;
372
373 if (mesh().moving())
374 {
375 return tmp<GeometricField<Type, fvPatchField, volMesh>>
376 (
377 new GeometricField<Type, fvPatchField, volMesh>
378 (
379 ddtIOobject,
380 mesh(),
381 rDeltaT.dimensions()
382 *alpha.dimensions()*rho.dimensions()*vf.dimensions(),
383 rDeltaT.value()*
384 (
385 coefft
386 *alpha.primitiveField()
387 *rho.primitiveField()
388 *vf.primitiveField() -
389 (
390 coefft0
391 *alpha.oldTime().primitiveField()
392 *rho.oldTime().primitiveField()
393 *vf.oldTime().primitiveField()*mesh().V0()
394
395 - coefft00
396 *alpha.oldTime().oldTime().primitiveField()
397 *rho.oldTime().oldTime().primitiveField()
398 *vf.oldTime().oldTime().primitiveField()*mesh().V00()
399 )/mesh().V()
400 ),
401 rDeltaT.value()*
402 (
403 coefft
404 *alpha.boundaryField()
405 *rho.boundaryField()
406 *vf.boundaryField() -
407 (
408 coefft0
409 *alpha.oldTime().boundaryField()
410 *rho.oldTime().boundaryField()
411 *vf.oldTime().boundaryField()
412
413 - coefft00
414 *alpha.oldTime().oldTime().boundaryField()
415 *rho.oldTime().oldTime().boundaryField()
416 *vf.oldTime().oldTime().boundaryField()
417 )
418 )
419 )
420 );
421 }
422 else
423 {
424 return tmp<GeometricField<Type, fvPatchField, volMesh>>
425 (
426 new GeometricField<Type, fvPatchField, volMesh>
427 (
428 ddtIOobject,
429 rDeltaT*
430 (
431 coefft*alpha*rho*vf
432 - coefft0*alpha.oldTime()*rho.oldTime()*vf.oldTime()
433 + coefft00*alpha.oldTime().oldTime()
434 *rho.oldTime().oldTime()*vf.oldTime().oldTime()
435 )
436 )
437 );
438 }
439}
440
441
442template<class Type>
443tmp<fvMatrix<Type>>
444backwardDdtScheme<Type>::fvmDdt
445(
446 const GeometricField<Type, fvPatchField, volMesh>& vf
447)
448{
449 tmp<fvMatrix<Type>> tfvm
450 (
451 new fvMatrix<Type>
452 (
453 vf,
454 vf.dimensions()*dimVol/dimTime
455 )
456 );
457
458 fvMatrix<Type>& fvm = tfvm.ref();
459
460 scalar rDeltaT = 1.0/deltaT_();
461
462 scalar deltaT = deltaT_();
463 scalar deltaT0 = deltaT0_(vf);
464
465 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
466 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
467 scalar coefft0 = coefft + coefft00;
468
469 fvm.diag() = (coefft*rDeltaT)*mesh().V();
470
471 if (mesh().moving())
472 {
473 fvm.source() = rDeltaT*
474 (
475 coefft0*vf.oldTime().primitiveField()*mesh().V0()
476 - coefft00*vf.oldTime().oldTime().primitiveField()
477 *mesh().V00()
478 );
479 }
480 else
481 {
482 fvm.source() = rDeltaT*mesh().V()*
483 (
484 coefft0*vf.oldTime().primitiveField()
485 - coefft00*vf.oldTime().oldTime().primitiveField()
486 );
487 }
488
489 return tfvm;
490}
491
492
493template<class Type>
494tmp<fvMatrix<Type>>
495backwardDdtScheme<Type>::fvmDdt
496(
497 const dimensionedScalar& rho,
498 const GeometricField<Type, fvPatchField, volMesh>& vf
499)
500{
501 tmp<fvMatrix<Type>> tfvm
502 (
503 new fvMatrix<Type>
504 (
505 vf,
506 rho.dimensions()*vf.dimensions()*dimVol/dimTime
507 )
508 );
509 fvMatrix<Type>& fvm = tfvm.ref();
510
511 scalar rDeltaT = 1.0/deltaT_();
512
513 scalar deltaT = deltaT_();
514 scalar deltaT0 = deltaT0_(vf);
515
516 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
517 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
518 scalar coefft0 = coefft + coefft00;
519
520 fvm.diag() = (coefft*rDeltaT*rho.value())*mesh().V();
521
522 if (mesh().moving())
523 {
524 fvm.source() = rDeltaT*rho.value()*
525 (
526 coefft0*vf.oldTime().primitiveField()*mesh().V0()
527 - coefft00*vf.oldTime().oldTime().primitiveField()
528 *mesh().V00()
529 );
530 }
531 else
532 {
533 fvm.source() = rDeltaT*mesh().V()*rho.value()*
534 (
535 coefft0*vf.oldTime().primitiveField()
536 - coefft00*vf.oldTime().oldTime().primitiveField()
537 );
538 }
539
540 return tfvm;
541}
542
543
544template<class Type>
545tmp<fvMatrix<Type>>
546backwardDdtScheme<Type>::fvmDdt
547(
548 const volScalarField& rho,
549 const GeometricField<Type, fvPatchField, volMesh>& vf
550)
551{
552 tmp<fvMatrix<Type>> tfvm
553 (
554 new fvMatrix<Type>
555 (
556 vf,
557 rho.dimensions()*vf.dimensions()*dimVol/dimTime
558 )
559 );
560 fvMatrix<Type>& fvm = tfvm.ref();
561
562 scalar rDeltaT = 1.0/deltaT_();
563
564 scalar deltaT = deltaT_();
565 scalar deltaT0 = deltaT0_(vf);
566
567 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
568 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
569 scalar coefft0 = coefft + coefft00;
570
571 fvm.diag() = (coefft*rDeltaT)*rho.primitiveField()*mesh().V();
572
573 if (mesh().moving())
574 {
575 fvm.source() = rDeltaT*
576 (
577 coefft0*rho.oldTime().primitiveField()
578 *vf.oldTime().primitiveField()*mesh().V0()
579 - coefft00*rho.oldTime().oldTime().primitiveField()
580 *vf.oldTime().oldTime().primitiveField()*mesh().V00()
581 );
582 }
583 else
584 {
585 fvm.source() = rDeltaT*mesh().V()*
586 (
587 coefft0*rho.oldTime().primitiveField()
588 *vf.oldTime().primitiveField()
589 - coefft00*rho.oldTime().oldTime().primitiveField()
590 *vf.oldTime().oldTime().primitiveField()
591 );
592 }
593
594 return tfvm;
595}
596
597
598template<class Type>
599tmp<fvMatrix<Type>>
600backwardDdtScheme<Type>::fvmDdt
601(
602 const volScalarField& alpha,
603 const volScalarField& rho,
604 const GeometricField<Type, fvPatchField, volMesh>& vf
605)
606{
607 tmp<fvMatrix<Type>> tfvm
608 (
609 new fvMatrix<Type>
610 (
611 vf,
612 alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
613 )
614 );
615 fvMatrix<Type>& fvm = tfvm.ref();
616
617 scalar rDeltaT = 1.0/deltaT_();
618
619 scalar deltaT = deltaT_();
620 scalar deltaT0 = deltaT0_(vf);
621
622 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
623 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
624 scalar coefft0 = coefft + coefft00;
625
626 fvm.diag() =
627 (coefft*rDeltaT)*alpha.primitiveField()*rho.primitiveField()*mesh().V();
628
629 if (mesh().moving())
630 {
631 fvm.source() = rDeltaT*
632 (
633 coefft0
634 *alpha.oldTime().primitiveField()
635 *rho.oldTime().primitiveField()
636 *vf.oldTime().primitiveField()*mesh().V0()
637
638 - coefft00
639 *alpha.oldTime().oldTime().primitiveField()
640 *rho.oldTime().oldTime().primitiveField()
641 *vf.oldTime().oldTime().primitiveField()*mesh().V00()
642 );
643 }
644 else
645 {
646 fvm.source() = rDeltaT*mesh().V()*
647 (
648 coefft0
649 *alpha.oldTime().primitiveField()
650 *rho.oldTime().primitiveField()
651 *vf.oldTime().primitiveField()
652
653 - coefft00
654 *alpha.oldTime().oldTime().primitiveField()
655 *rho.oldTime().oldTime().primitiveField()
656 *vf.oldTime().oldTime().primitiveField()
657 );
658 }
659
660 return tfvm;
661}
662
663
664template<class Type>
665tmp<typename backwardDdtScheme<Type>::fluxFieldType>
666backwardDdtScheme<Type>::fvcDdtUfCorr
667(
668 const GeometricField<Type, fvPatchField, volMesh>& U,
669 const GeometricField<Type, fvsPatchField, surfaceMesh>& Uf
670)
671{
672 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
673
674 scalar deltaT = deltaT_();
675 scalar deltaT0 = deltaT0_(U);
676
677 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
678 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
679 scalar coefft0 = coefft + coefft00;
680
681 return tmp<fluxFieldType>
682 (
683 new fluxFieldType
684 (
685 IOobject
686 (
687 "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
688 mesh().time().timeName(),
689 mesh()
690 ),
691 this->fvcDdtPhiCoeff(U.oldTime(), (mesh().Sf() & Uf.oldTime()))
692 *rDeltaT
693 *(
694 mesh().Sf()
695 & (
696 (coefft0*Uf.oldTime() - coefft00*Uf.oldTime().oldTime())
697 - fvc::interpolate
698 (
699 coefft0*U.oldTime() - coefft00*U.oldTime().oldTime()
700 )
701 )
702 )
703 )
704 );
705}
706
707
708template<class Type>
709tmp<typename backwardDdtScheme<Type>::fluxFieldType>
710backwardDdtScheme<Type>::fvcDdtPhiCorr
711(
712 const GeometricField<Type, fvPatchField, volMesh>& U,
713 const fluxFieldType& phi
714)
715{
716 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
717
718 scalar deltaT = deltaT_();
719 scalar deltaT0 = deltaT0_(U);
720
721 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
722 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
723 scalar coefft0 = coefft + coefft00;
724
725 return tmp<fluxFieldType>
726 (
727 new fluxFieldType
728 (
729 IOobject
730 (
731 "ddtCorr(" + U.name() + ',' + phi.name() + ')',
732 mesh().time().timeName(),
733 mesh()
734 ),
735 this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())
736 *rDeltaT
737 *(
738 (coefft0*phi.oldTime() - coefft00*phi.oldTime().oldTime())
739 - fvc::dotInterpolate
740 (
741 mesh().Sf(),
742 coefft0*U.oldTime() - coefft00*U.oldTime().oldTime()
743 )
744 )
745 )
746 );
747}
748
749
750template<class Type>
751tmp<typename backwardDdtScheme<Type>::fluxFieldType>
752backwardDdtScheme<Type>::fvcDdtUfCorr
753(
754 const volScalarField& rho,
755 const GeometricField<Type, fvPatchField, volMesh>& U,
756 const GeometricField<Type, fvsPatchField, surfaceMesh>& Uf
757)
758{
759 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
760
761 scalar deltaT = deltaT_();
762 scalar deltaT0 = deltaT0_(U);
763
764 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
765 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
766 scalar coefft0 = coefft + coefft00;
767
768 if
769 (
770 U.dimensions() == dimVelocity
771 && Uf.dimensions() == rho.dimensions()*dimVelocity
772 )
773 {
774 GeometricField<Type, fvPatchField, volMesh> rhoU0
775 (
776 rho.oldTime()*U.oldTime()
777 );
778
779 GeometricField<Type, fvPatchField, volMesh> rhoU00
780 (
781 rho.oldTime().oldTime()*U.oldTime().oldTime()
782 );
783
784 return tmp<fluxFieldType>
785 (
786 new fluxFieldType
787 (
788 IOobject
789 (
790 "ddtCorr("
791 + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
792 mesh().time().timeName(),
793 mesh()
794 ),
795 this->fvcDdtPhiCoeff
796 (
797 rhoU0,
798 mesh().Sf() & Uf.oldTime(),
799 rho.oldTime()
800 )
801 *rDeltaT
802 *(
803 mesh().Sf()
804 & (
805 (coefft0*Uf.oldTime() - coefft00*Uf.oldTime().oldTime())
806 - fvc::interpolate(coefft0*rhoU0 - coefft00*rhoU00)
807 )
808 )
809 )
810 );
811 }
812 else if
813 (
814 U.dimensions() == rho.dimensions()*dimVelocity
815 && Uf.dimensions() == rho.dimensions()*dimVelocity
816 )
817 {
818 return tmp<fluxFieldType>
819 (
820 new fluxFieldType
821 (
822 IOobject
823 (
824 "ddtCorr("
825 + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
826 mesh().time().timeName(),
827 mesh()
828 ),
829 this->fvcDdtPhiCoeff
830 (
831 U.oldTime(),
832 mesh().Sf() & Uf.oldTime(),
833 rho.oldTime()
834 )
835 *rDeltaT
836 *(
837 mesh().Sf()
838 & (
839 (coefft0*Uf.oldTime() - coefft00*Uf.oldTime().oldTime())
840 - fvc::interpolate
841 (
842 coefft0*U.oldTime()
843 - coefft00*U.oldTime().oldTime()
844 )
845 )
846 )
847 )
848 );
849 }
850 else
851 {
852 FatalErrorInFunction
853 << "dimensions of phi are not correct"
854 << abort(FatalError);
855
856 return fluxFieldType::null();
857 }
858}
859
860
861template<class Type>
862tmp<typename backwardDdtScheme<Type>::fluxFieldType>
863backwardDdtScheme<Type>::fvcDdtPhiCorr
864(
865 const volScalarField& rho,
866 const GeometricField<Type, fvPatchField, volMesh>& U,
867 const fluxFieldType& phi
868)
869{
870 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
871
872 scalar deltaT = deltaT_();
873 scalar deltaT0 = deltaT0_(U);
874
875 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
876 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
877 scalar coefft0 = coefft + coefft00;
878
879 if
880 (
881 U.dimensions() == dimVelocity
882 && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
883 )
884 {
885 GeometricField<Type, fvPatchField, volMesh> rhoU0
886 (
887 rho.oldTime()*U.oldTime()
888 );
889
890 GeometricField<Type, fvPatchField, volMesh> rhoU00
891 (
892 rho.oldTime().oldTime()*U.oldTime().oldTime()
893 );
894
895 return tmp<fluxFieldType>
896 (
897 new fluxFieldType
898 (
899 IOobject
900 (
901 "ddtCorr("
902 + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
903 mesh().time().timeName(),
904 mesh()
905 ),
906 this->fvcDdtPhiCoeff(rhoU0, phi.oldTime(), rho.oldTime())
907 *rDeltaT
908 *(
909 (coefft0*phi.oldTime() - coefft00*phi.oldTime().oldTime())
910 - fvc::dotInterpolate
911 (
912 mesh().Sf(),
913 coefft0*rhoU0 - coefft00*rhoU00
914 )
915 )
916 )
917 );
918 }
919 else if
920 (
921 U.dimensions() == rho.dimensions()*dimVelocity
922 && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
923 )
924 {
925 return tmp<fluxFieldType>
926 (
927 new fluxFieldType
928 (
929 IOobject
930 (
931 "ddtCorr("
932 + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
933 mesh().time().timeName(),
934 mesh()
935 ),
936 this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime(), rho.oldTime())
937 *rDeltaT
938 *(
939 (coefft0*phi.oldTime() - coefft00*phi.oldTime().oldTime())
940 - fvc::dotInterpolate
941 (
942 mesh().Sf(),
943 coefft0*U.oldTime() - coefft00*U.oldTime().oldTime()
944 )
945 )
946 )
947 );
948 }
949 else
950 {
951 FatalErrorInFunction
952 << "dimensions of phi are not correct"
953 << abort(FatalError);
954
955 return fluxFieldType::null();
956 }
957}
958
959
960template<class Type>
961tmp<surfaceScalarField> backwardDdtScheme<Type>::meshPhi
962(
963 const GeometricField<Type, fvPatchField, volMesh>& vf
964)
965{
966 scalar deltaT = deltaT_();
967 scalar deltaT0 = deltaT0_(vf);
968
969 // Coefficient for t-3/2 (between times 0 and 00)
970 scalar coefft0_00 = deltaT/(deltaT + deltaT0);
971
972 // Coefficient for t-1/2 (between times n and 0)
973 scalar coefftn_0 = 1 + coefft0_00;
974
975 return tmp<surfaceScalarField>
976 (
977 new surfaceScalarField
978 (
979 IOobject
980 (
981 mesh().phi().name(),
982 mesh().time().timeName(),
983 mesh(),
984 IOobject::NO_READ,
985 IOobject::NO_WRITE,
986 false
987 ),
988 coefftn_0*mesh().phi() - coefft0_00*mesh().phi().oldTime()
989 )
990 );
991}
992
993
994// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
995
996} // End namespace fv
997
998// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
999
1000} // End namespace Foam
1001
1002// ************************************************************************* //
phi
surfaceScalarField & phi
Definition: setRegionFluidFields.H:8
Foam::DimensionedField::dimensions
const dimensionSet & dimensions() const
Return dimensions.
Definition: DimensionedFieldI.H:49
Foam::GeometricField
Generic GeometricField class.
Definition: GeometricField.H:80
Foam::GeometricField::oldTime
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
Definition: GeometricField.C:978
Foam::GeometricField::primitiveField
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Definition: GeometricFieldI.H:53
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
Foam::IOobject::name
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
Foam::calculatedFvPatchField
This boundary condition is not designed to be evaluated; it is assmued that the value is assigned via...
Definition: calculatedFvPatchField.H:69
Foam::dimensioned
Generic dimensioned Type class.
Definition: dimensionedType.H:67
Foam::dimensioned::dimensions
const dimensionSet & dimensions() const
Return const reference to dimensions.
Definition: dimensionedType.C:420
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
Foam::dimensioned::name
const word & name() const
Return const reference to name.
Definition: dimensionedType.C:406
Foam::fvMatrix
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
Foam::fvMatrix::source
Field< Type > & source() noexcept
Definition: fvMatrix.H:458
Foam::fv::backwardDdtScheme::fvcDdtUfCorr
tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Definition: backwardDdtScheme.C:667
Foam::fv::backwardDdtScheme::fvcDdt
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)
Definition: backwardDdtScheme.C:80
Foam::fv::backwardDdtScheme::fvcDdtPhiCorr
tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
Definition: backwardDdtScheme.C:711
Foam::fv::backwardDdtScheme::fvmDdt
tmp< fvMatrix< Type > > fvmDdt(const GeometricField< Type, fvPatchField, volMesh > &)
Definition: backwardDdtScheme.C:445
Foam::fv::backwardDdtScheme::meshPhi
tmp< surfaceScalarField > meshPhi(const GeometricField< Type, fvPatchField, volMesh > &)
Definition: backwardDdtScheme.C:962
Foam::lduMatrix::diag
scalarField & diag()
Definition: lduMatrix.C:192
Foam::tmp
A class for managing temporary objects.
Definition: tmp.H:65
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
U
U
Definition: pEqn.H:72
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
oldTime
Info<< "Creating field kinetic energy K\n"<< endl;volScalarField K("K", 0.5 *magSqr(U));if(U.nOldTimes()){ volVectorField *Uold=&U.oldTime();volScalarField *Kold=&K.oldTime();*Kold==0.5 *magSqr(*Uold);while(Uold->nOldTimes()) { Uold=&Uold-> oldTime()
Uf
autoPtr< surfaceVectorField > Uf
Definition: createUfIfPresent.H:33
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
fvMatrices.H
A special matrix type and solver, designed for finite volume solutions of scalar equations.
fvcDiv.H
Calculate the divergence of the given field.
timeName
word timeName
Definition: getTimeIndex.H:3
Foam::fvc::interpolate
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Foam::fvc::dotInterpolate
static tmp< GeometricField< typename innerProduct< vector, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const surfaceVectorField &Sf, const GeometricField< Type, fvPatchField, volMesh > &tvf)
Interpolate field onto faces.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:34
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::dimArea
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
Foam::dimVelocity
const dimensionSet dimVelocity
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
Foam::FatalError
error FatalError
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::dimVol
const dimensionSet dimVol(dimVolume)
Older spelling for dimVolume.
Definition: dimensionSets.H:64
timeIndex
label timeIndex
Definition: getTimeIndex.H:30
fv
labelList fv(nPoints)
alpha
volScalarField & alpha
Definition: readThermalProperties.H:212