boundedBackwardFaDdtScheme.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) 2016-2017 Wikki Ltd
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
29#include "facDiv.H"
30#include "faMatrices.H"
31
32// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33
34namespace Foam
35{
36
37// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38
39namespace fa
40{
41
42// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43
44scalar boundedBackwardFaDdtScheme::deltaT_() const
45{
46 return mesh().time().deltaT().value();
47}
48
49
50scalar boundedBackwardFaDdtScheme::deltaT0_() const
51{
52 return mesh().time().deltaT0().value();
53}
54
55
56// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
57
59(
60 const dimensionedScalar dt
61)
62{
63 // No change compared to backward
64
65 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
66
67 IOobject ddtIOobject
68 (
69 "ddt("+dt.name()+')',
70 mesh()().time().timeName(),
71 mesh()(),
74 );
75
76 scalar deltaT = deltaT_();
77 scalar deltaT0 = deltaT0_();
78
79 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
80 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
81 scalar coefft0 = coefft + coefft00;
82
83 if (mesh().moving())
84 {
86 (
88 (
89 ddtIOobject,
90 mesh(),
92 )
93 );
94
95 tdtdt.ref().primitiveFieldRef() = rDeltaT.value()*dt.value()*
96 (
97 coefft - (coefft0*mesh().S0() - coefft00*mesh().S00())/mesh().S()
98 );
99
100 return tdtdt;
101 }
102
104 (
105 ddtIOobject,
106 mesh(),
108 calculatedFaPatchScalarField::typeName
109 );
110}
111
112
114(
115 const dimensionedScalar dt
116)
117{
118 // No change compared to backward
119
120 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
121
122 IOobject ddtIOobject
123 (
124 "ddt("+dt.name()+')',
125 mesh()().time().timeName(),
126 mesh()(),
129 );
130
131 scalar deltaT = deltaT_();
132 scalar deltaT0 = deltaT0_();
133
134 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
135 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
136 scalar coefft0 = coefft + coefft00;
137
139 (
141 (
142 ddtIOobject,
143 mesh(),
144 -rDeltaT*(coefft0 - coefft00)*dt
145 )
146 );
147
148 if (mesh().moving())
149 {
150 tdtdt0.ref().primitiveFieldRef() = (-rDeltaT.value()*dt.value())*
151 (
152 (coefft0*mesh().S0() - coefft00*mesh().S00())/mesh().S()
153 );
154 }
155
156 return tdtdt0;
157}
158
159
161(
162 const areaScalarField& vf
163)
164{
165 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
166
167 IOobject ddtIOobject
168 (
169 "ddt("+vf.name()+')',
170 mesh()().time().timeName(),
171 mesh()(),
174 );
175
176 scalar deltaT = deltaT_();
177 scalar deltaT0 = deltaT0_(vf);
178
179 // Calculate unboundedness indicator
180 // Note: all times moved by one because access to internal field
181 // copies current field into the old-time level.
182 areaScalarField phict
183 (
184 mag
185 (
186 vf.oldTime().oldTime()
187 - vf.oldTime().oldTime().oldTime()
188 )/
189 (
190 mag
191 (
192 vf.oldTime()
193 - vf.oldTime().oldTime()
194 )
195 + dimensionedScalar("small", vf.dimensions(), SMALL)
196 )
197 );
198
199 areaScalarField limiter(pos(phict) - pos(phict - scalar(1)));
200
201 areaScalarField coefft(scalar(1) + limiter*deltaT/(deltaT + deltaT0));
202 areaScalarField coefft00(limiter*sqr(deltaT)/(deltaT0*(deltaT + deltaT0)));
203 areaScalarField coefft0(coefft + coefft00);
204
205 if (mesh().moving())
206 {
208 (
210 (
211 ddtIOobject,
212 mesh(),
213 rDeltaT.dimensions()*vf.dimensions(),
214 rDeltaT.value()*
215 (
216 coefft*vf.primitiveField() -
217 (
218 coefft0.primitiveField()
219 *vf.oldTime().primitiveField()*mesh().S0()
220 - coefft00.primitiveField()
222 *mesh().S00()
223 )/mesh().S()
224 ),
225 rDeltaT.value()*
226 (
227 coefft.boundaryField()*vf.boundaryField() -
228 (
229 coefft0.boundaryField()*
231 - coefft00.boundaryField()*
233 )
234 )
235 )
236 );
237 }
238 else
239 {
241 (
243 (
244 ddtIOobject,
245 rDeltaT*
246 (
247 coefft*vf
248 - coefft0*vf.oldTime()
249 + coefft00*vf.oldTime().oldTime()
250 )
251 )
252 );
253 }
254}
255
256
258(
259 const areaScalarField& vf
260)
261{
262 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
263
264 IOobject ddtIOobject
265 (
266 "ddt0("+vf.name()+')',
267 mesh()().time().timeName(),
268 mesh()(),
271 );
272
273 scalar deltaT = deltaT_();
274 scalar deltaT0 = deltaT0_(vf);
275
276 // Calculate unboundedness indicator
277 // Note: all times moved by one because access to internal field
278 // copies current field into the old-time level.
279 areaScalarField phict
280 (
281 mag
282 (
283 vf.oldTime().oldTime()
284 - vf.oldTime().oldTime().oldTime()
285 )/
286 (
287 mag
288 (
289 vf.oldTime()
290 - vf.oldTime().oldTime()
291 )
292 + dimensionedScalar("small", vf.dimensions(), SMALL)
293 )
294 );
295
296 areaScalarField limiter(pos(phict) - pos(phict - scalar(1)));
297
298 areaScalarField coefft(scalar(1) + limiter*deltaT/(deltaT + deltaT0));
299 areaScalarField coefft00(limiter*sqr(deltaT)/(deltaT0*(deltaT + deltaT0)));
300 areaScalarField coefft0(coefft + coefft00);
301
302 if (mesh().moving())
303 {
305 (
307 (
308 ddtIOobject,
309 mesh(),
310 rDeltaT.dimensions()*vf.dimensions(),
311 rDeltaT.value()*
312 (
313 - (
314 coefft0.primitiveField()*
315 vf.oldTime().primitiveField()*mesh().S0()
316 - coefft00.primitiveField()*
318 *mesh().S00()
319 )/mesh().S()
320 ),
321 rDeltaT.value()*
322 (
323 - (
324 coefft0.boundaryField()*
326 - coefft00.boundaryField()*
328 )
329 )
330 )
331 );
332 }
333 else
334 {
336 (
338 (
339 ddtIOobject,
340 rDeltaT*
341 (
342 - coefft0*vf.oldTime()
343 + coefft00*vf.oldTime().oldTime()
344 )
345 )
346 );
347 }
348}
349
350
352(
353 const dimensionedScalar& rho,
354 const areaScalarField& vf
355)
356{
357 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
358
359 IOobject ddtIOobject
360 (
361 "ddt("+rho.name()+','+vf.name()+')',
362 mesh()().time().timeName(),
363 mesh()(),
366 );
367
368 scalar deltaT = deltaT_();
369 scalar deltaT0 = deltaT0_(vf);
370
371 // Calculate unboundedness indicator
372 // Note: all times moved by one because access to internal field
373 // copies current field into the old-time level.
374 areaScalarField phict
375 (
376 mag
377 (
378 vf.oldTime().oldTime()
379 - vf.oldTime().oldTime().oldTime()
380 )/
381 (
382 mag
383 (
384 vf.oldTime()
385 - vf.oldTime().oldTime()
386 )
387 + dimensionedScalar("small", vf.dimensions(), SMALL)
388 )
389 );
390
391 areaScalarField limiter(pos(phict) - pos(phict - scalar(1)));
392
393 areaScalarField coefft(scalar(1) + limiter*deltaT/(deltaT + deltaT0));
394 areaScalarField coefft00(limiter*sqr(deltaT)/(deltaT0*(deltaT + deltaT0)));
395 areaScalarField coefft0(coefft + coefft00);
396
397 if (mesh().moving())
398 {
400 (
402 (
403 ddtIOobject,
404 mesh(),
405 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
406 rDeltaT.value()*rho.value()*
407 (
408 coefft*vf.primitiveField() -
409 (
410 coefft0.primitiveField()*
411 vf.oldTime().primitiveField()*mesh().S0()
412 - coefft00.primitiveField()*
414 *mesh().S00()
415 )/mesh().S()
416 ),
417 rDeltaT.value()*rho.value()*
418 (
419 coefft.boundaryField()*vf.boundaryField() -
420 (
421 coefft0.boundaryField()*
423 - coefft00.boundaryField()*
425 )
426 )
427 )
428 );
429 }
430 else
431 {
433 (
435 (
436 ddtIOobject,
437 rDeltaT*rho*
438 (
439 coefft*vf
440 - coefft0*vf.oldTime()
441 + coefft00*vf.oldTime().oldTime()
442 )
443 )
444 );
445 }
446}
447
449(
450 const dimensionedScalar& rho,
451 const areaScalarField& vf
452)
453{
454 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
455
456 IOobject ddtIOobject
457 (
458 "ddt0("+rho.name()+','+vf.name()+')',
459 mesh()().time().timeName(),
460 mesh()(),
463 );
464
465 scalar deltaT = deltaT_();
466 scalar deltaT0 = deltaT0_(vf);
467
468 // Calculate unboundedness indicator
469 // Note: all times moved by one because access to internal field
470 // copies current field into the old-time level.
471 areaScalarField phict
472 (
473 mag
474 (
475 vf.oldTime().oldTime()
476 - vf.oldTime().oldTime().oldTime()
477 )/
478 (
479 mag
480 (
481 vf.oldTime()
482 - vf.oldTime().oldTime()
483 )
484 + dimensionedScalar("small", vf.dimensions(), SMALL)
485 )
486 );
487
488 areaScalarField limiter(pos(phict) - pos(phict - scalar(1)));
489
490 areaScalarField coefft(scalar(1) + limiter*deltaT/(deltaT + deltaT0));
491 areaScalarField coefft00(limiter*sqr(deltaT)/(deltaT0*(deltaT + deltaT0)));
492 areaScalarField coefft0(coefft + coefft00);
493
494 if (mesh().moving())
495 {
497 (
499 (
500 ddtIOobject,
501 mesh(),
502 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
503 rDeltaT.value()*rho.value()*
504 (
505 -(
506 coefft0.primitiveField()*
507 vf.oldTime().primitiveField()*mesh().S0()
508 - coefft00.primitiveField()*
510 *mesh().S00()
511 )/mesh().S()
512 ),
513 rDeltaT.value()*rho.value()*
514 (
515 -(
516 coefft0.boundaryField()*
518 - coefft00.boundaryField()*
520 )
521 )
522 )
523 );
524 }
525 else
526 {
528 (
530 (
531 ddtIOobject,
532 rDeltaT*rho*
533 (
534 - coefft0*vf.oldTime()
535 + coefft00*vf.oldTime().oldTime()
536 )
537 )
538 );
539 }
540}
541
542
544(
545 const areaScalarField& rho,
546 const areaScalarField& vf
547)
548{
549 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
550
551 IOobject ddtIOobject
552 (
553 "ddt("+rho.name()+','+vf.name()+')',
554 mesh()().time().timeName(),
555 mesh()(),
558 );
559
560 scalar deltaT = deltaT_();
561 scalar deltaT0 = deltaT0_(vf);
562
563 // Calculate unboundedness indicator
564 // Note: all times moved by one because access to internal field
565 // copies current field into the old-time level.
566 areaScalarField phict
567 (
568 mag
569 (
570 rho.oldTime().oldTime()*vf.oldTime().oldTime()
571 - rho.oldTime().oldTime().oldTime()*vf.oldTime().oldTime().oldTime()
572 )/
573 (
574 mag
575 (
576 rho.oldTime()*vf.oldTime()
577 - rho.oldTime().oldTime()*vf.oldTime().oldTime()
578 )
579 + dimensionedScalar("small", rho.dimensions()*vf.dimensions(), SMALL)
580 )
581 );
582
583 areaScalarField limiter(pos(phict) - pos(phict - scalar(1)));
584
585 areaScalarField coefft(scalar(1) + limiter*deltaT/(deltaT + deltaT0));
586 areaScalarField coefft00(limiter*sqr(deltaT)/(deltaT0*(deltaT + deltaT0)));
587 areaScalarField coefft0(coefft + coefft00);
588
589 if (mesh().moving())
590 {
592 (
594 (
595 ddtIOobject,
596 mesh(),
597 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
598 rDeltaT.value()*
599 (
600 coefft*rho.primitiveField()*vf.primitiveField() -
601 (
602 coefft0.primitiveField()*
603 rho.oldTime().primitiveField()*
604 vf.oldTime().primitiveField()*mesh().S0()
605 - coefft00.primitiveField()*
606 rho.oldTime().oldTime().primitiveField()
607 *vf.oldTime().oldTime().primitiveField()*mesh().S00()
608 )/mesh().S()
609 ),
610 rDeltaT.value()*
611 (
612 coefft.boundaryField()*vf.boundaryField() -
613 (
614 coefft0.boundaryField()*
615 rho.oldTime().boundaryField()*
617 - coefft00.boundaryField()*
618 rho.oldTime().oldTime().boundaryField()*
620 )
621 )
622 )
623 );
624 }
625 else
626 {
628 (
630 (
631 ddtIOobject,
632 rDeltaT*
633 (
634 coefft*rho*vf
635 - coefft0*rho.oldTime()*vf.oldTime()
636 + coefft00*rho.oldTime().oldTime()*vf.oldTime().oldTime()
637 )
638 )
639 );
640 }
641}
642
643
645(
646 const areaScalarField& rho,
647 const areaScalarField& vf
648)
649{
650 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
651
652 IOobject ddtIOobject
653 (
654 "ddt0("+rho.name()+','+vf.name()+')',
655 mesh()().time().timeName(),
656 mesh()(),
659 );
660
661 scalar deltaT = deltaT_();
662 scalar deltaT0 = deltaT0_(vf);
663
664 // Calculate unboundedness indicator
665 // Note: all times moved by one because access to internal field
666 // copies current field into the old-time level.
667 areaScalarField phict
668 (
669 mag
670 (
671 rho.oldTime().oldTime()*vf.oldTime().oldTime()
672 - rho.oldTime().oldTime().oldTime()*vf.oldTime().oldTime().oldTime()
673 )/
674 (
675 mag
676 (
677 rho.oldTime()*vf.oldTime()
678 - rho.oldTime().oldTime()*vf.oldTime().oldTime()
679 )
680 + dimensionedScalar("small", rho.dimensions()*vf.dimensions(), SMALL)
681 )
682 );
683
684 areaScalarField limiter(pos(phict) - pos(phict - scalar(1)));
685
686 areaScalarField coefft(scalar(1) + limiter*deltaT/(deltaT + deltaT0));
687 areaScalarField coefft00(limiter*sqr(deltaT)/(deltaT0*(deltaT + deltaT0)));
688 areaScalarField coefft0(coefft + coefft00);
689
690 if (mesh().moving())
691 {
693 (
695 (
696 ddtIOobject,
697 mesh(),
698 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
699 rDeltaT.value()*
700 (
701 - (
702 coefft0.primitiveField()*
703 rho.oldTime().primitiveField()*
704 vf.oldTime().primitiveField()*mesh().S0()
705 - coefft00.primitiveField()*
706 rho.oldTime().oldTime().primitiveField()*
707 vf.oldTime().oldTime().primitiveField()*mesh().S00()
708 )/mesh().S()
709 ),
710 rDeltaT.value()*
711 (
712 - (
713 coefft0.boundaryField()*
714 rho.oldTime().boundaryField()*
716 - coefft00.boundaryField()*
717 rho.oldTime().oldTime().boundaryField()*
719 )
720 )
721 )
722 );
723 }
724 else
725 {
727 (
729 (
730 ddtIOobject,
731 rDeltaT*
732 (
733 - coefft0*rho.oldTime()*vf.oldTime()
734 + coefft00*rho.oldTime().oldTime()*vf.oldTime().oldTime()
735 )
736 )
737 );
738 }
739}
740
741
743(
744 const areaScalarField& vf
745)
746{
748 (
750 (
751 vf,
753 )
754 );
755
756 faScalarMatrix& fam = tfam.ref();
757
758 scalar rDeltaT = 1.0/deltaT_();
759
760 scalar deltaT = deltaT_();
761 scalar deltaT0 = deltaT0_(vf);
762
763 // Calculate unboundedness indicator
764 // Note: all times moved by one because access to internal field
765 // copies current field into the old-time level.
766 scalarField phict
767 (
768 mag
769 (
772 )/
773 (
774 mag
775 (
777 - vf.oldTime().oldTime().internalField()
778 )
779 + SMALL
780 )
781 );
782
783 scalarField limiter(pos(phict) - pos(phict - 1.0));
784
785 scalarField coefft(1.0 + limiter*deltaT/(deltaT + deltaT0));
786 scalarField coefft00(limiter*deltaT*deltaT/(deltaT0*(deltaT + deltaT0)));
787 scalarField coefft0(coefft + coefft00);
788
789 fam.diag() = (coefft*rDeltaT)*mesh().S();
790
791 if (mesh().moving())
792 {
793 fam.source() = rDeltaT*
794 (
795 coefft0*vf.oldTime().primitiveField()*mesh().S0()
796 - coefft00*vf.oldTime().oldTime().primitiveField()
797 *mesh().S00()
798 );
799 }
800 else
801 {
802 fam.source() = rDeltaT*mesh().S()*
803 (
804 coefft0*vf.oldTime().primitiveField()
805 - coefft00*vf.oldTime().oldTime().primitiveField()
806 );
807 }
808
809 return tfam;
810}
811
812
814(
815 const dimensionedScalar& rho,
816 const areaScalarField& vf
817)
818{
820 (
822 (
823 vf,
824 rho.dimensions()*vf.dimensions()*dimArea/dimTime
825 )
826 );
827 faScalarMatrix& fam = tfam.ref();
828
829 scalar rDeltaT = 1.0/deltaT_();
830
831 scalar deltaT = deltaT_();
832 scalar deltaT0 = deltaT0_(vf);
833
834 // Calculate unboundedness indicator
835 // Note: all times moved by one because access to internal field
836 // copies current field into the old-time level.
837 scalarField phict
838 (
839 mag
840 (
843 )/
844 (
845 mag
846 (
848 - vf.oldTime().oldTime().internalField()
849 )
850 + SMALL
851 )
852 );
853
854 scalarField limiter(pos(phict) - pos(phict - 1.0));
855
856 scalarField coefft(1.0 + limiter*deltaT/(deltaT + deltaT0));
857 scalarField coefft00(limiter*deltaT*deltaT/(deltaT0*(deltaT + deltaT0)));
858 scalarField coefft0(coefft + coefft00);
859
860 fam.diag() = (coefft*rDeltaT*rho.value())*mesh().S();
861
862 if (mesh().moving())
863 {
864 fam.source() = rDeltaT*rho.value()*
865 (
866 coefft0*vf.oldTime().primitiveField()*mesh().S0()
867 - coefft00*vf.oldTime().oldTime().primitiveField()
868 *mesh().S00()
869 );
870 }
871 else
872 {
873 fam.source() = rDeltaT*mesh().S()*rho.value()*
874 (
875 coefft0*vf.oldTime().primitiveField()
876 - coefft00*vf.oldTime().oldTime().primitiveField()
877 );
878 }
879
880 return tfam;
881}
882
883
885(
886 const areaScalarField& rho,
887 const areaScalarField& vf
888)
889{
891 (
893 (
894 vf,
895 rho.dimensions()*vf.dimensions()*dimArea/dimTime
896 )
897 );
898 faScalarMatrix& fam = tfam.ref();
899
900 scalar rDeltaT = 1.0/deltaT_();
901
902 scalar deltaT = deltaT_();
903 scalar deltaT0 = deltaT0_(vf);
904
905 // Calculate unboundedness indicator
906 // Note: all times moved by one because access to internal field
907 // copies current field into the old-time level.
908 scalarField phict
909 (
910 mag
911 (
912 rho.oldTime().oldTime().internalField()*
914 - rho.oldTime().oldTime().oldTime().internalField()*
916 )/
917 (
918 mag
919 (
920 rho.oldTime().internalField()*
922 - rho.oldTime().oldTime().internalField()*
924 )
925 + SMALL
926 )
927 );
928
929 scalarField limiter(pos(phict) - pos(phict - 1.0));
930
931 scalarField coefft(1.0 + limiter*deltaT/(deltaT + deltaT0));
932 scalarField coefft00(limiter*deltaT*deltaT/(deltaT0*(deltaT + deltaT0)));
933 scalarField coefft0(coefft + coefft00);
934
935 fam.diag() = (coefft*rDeltaT)*rho.primitiveField()*mesh().S();
936
937 if (mesh().moving())
938 {
939 fam.source() = rDeltaT*
940 (
941 coefft0*rho.oldTime().primitiveField()
942 *vf.oldTime().primitiveField()*mesh().S0()
943 - coefft00*rho.oldTime().oldTime().primitiveField()
945 );
946 }
947 else
948 {
949 fam.source() = rDeltaT*mesh().S()*
950 (
951 coefft0*rho.oldTime().primitiveField()
952 *vf.oldTime().primitiveField()
953 - coefft00*rho.oldTime().oldTime().primitiveField()
955 );
956 }
957
958 return tfam;
959}
960
961
962// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
963
965
968
969
970// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
971
972} // End namespace fa
973
974// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
975
976} // End namespace Foam
977
978// ************************************************************************* //
const dimensionSet & dimensions() const
Return dimensions.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const Internal & internalField() const
Return a const-reference to the dimensioned internal field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
dimensionedScalar deltaT0() const
Return old time step.
Definition: TimeStateI.H:61
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:55
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
const dimensionSet & dimensions() const
Return const reference to dimensions.
const Type & value() const
Return const reference to value.
const word & name() const
Return const reference to name.
A special matrix type and solver, designed for finite area solutions of scalar equations....
Definition: faMatrix.H:76
const Time & time() const
Return reference to time.
Definition: faMesh.C:673
const DimensionedField< scalar, areaMesh > & S00() const
Return old-old-time face areas.
Definition: faMesh.C:806
const DimensionedField< scalar, areaMesh > & S() const
Return face areas.
Definition: faMesh.C:780
const DimensionedField< scalar, areaMesh > & S0() const
Return old-time face areas.
Definition: faMesh.C:792
bool moving() const
Is mesh moving.
Definition: faMesh.H:776
tmp< areaScalarField > facDdt(const dimensionedScalar)
tmp< areaScalarField > facDdt0(const dimensionedScalar)
const faMesh & mesh() const
Return mesh reference.
tmp< faScalarMatrix > famDdt(const areaScalarField &)
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
Calculate the divergence of the given field.
word timeName
Definition: getTimeIndex.H:3
faDdtScheme< scalar >::addIstreamConstructorToTable< boundedBackwardFaDdtScheme > addboundedBackwardFaDdtSchemeIstreamConstructorToTable_
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pos(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
tmp< areaScalarField > limiter(const areaScalarField &phi)
Definition: faNVDscheme.C:38
Calculate the matrix for the second temporal derivative.