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