Loading...
Searching...
No Matches
backwardFaDdtScheme.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 "backwardFaDdtScheme.H"
29#include "faMatrices.H"
30
31// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32
33namespace Foam
34{
35
36// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37
38namespace fa
39{
40
41// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42
43template<class Type>
44scalar backwardFaDdtScheme<Type>::deltaT_() const
45{
46 return mesh().time().deltaTValue();
47}
48
49
50template<class Type>
51scalar backwardFaDdtScheme<Type>::deltaT0_() const
52{
53 return mesh().time().deltaT0Value();
54}
55
56
57template<class Type>
58template<class GeoField>
59scalar backwardFaDdtScheme<Type>::deltaT0_(const GeoField& vf) const
60{
61 if (vf.oldTime().timeIndex() == vf.oldTime().oldTime().timeIndex())
62 {
63 return GREAT;
64 }
65 else
66 {
67 return mesh().time().deltaT0Value();
68 }
70
71
72// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
73
74template<class Type>
77(
78 const dimensioned<Type> dt
79)
80{
81 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
82
83 const IOobject ddtIOobject
84 (
85 mesh().thisDb().newIOobject
86 (
87 "ddt("+dt.name()+')',
88 { IOobject::REGISTER }
89 )
90 );
91
92 scalar deltaT = deltaT_();
93 scalar deltaT0 = deltaT0_();
94
95 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
96 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
97 scalar coefft0 = coefft + coefft00;
98
99 if (mesh().moving())
100 {
102 (
104 (
105 ddtIOobject,
106 mesh(),
108 )
109 );
110
111 tdtdt.ref().primitiveFieldRef() = rDeltaT.value()*dt.value()*
112 (
113 coefft - (coefft0*mesh().S0() - coefft00*mesh().S00())/mesh().S()
114 );
115
116 return tdtdt;
117 }
118
120 (
121 ddtIOobject,
122 mesh(),
125 );
126}
127
128
129template<class Type>
132(
133 const dimensioned<Type> dt
134)
135{
136 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
137
138 const IOobject ddtIOobject
139 (
140 mesh().thisDb().newIOobject
141 (
142 "ddt("+dt.name()+')',
143 { IOobject::REGISTER }
144 )
145 );
146
147 scalar deltaT = deltaT_();
148 scalar deltaT0 = deltaT0_();
149
150 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
151 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
152 scalar coefft0 = coefft + coefft00;
153
155 (
157 (
158 ddtIOobject,
159 mesh(),
160 -rDeltaT*(coefft0 - coefft00)*dt
161 )
162 );
163
164 if (mesh().moving())
165 {
166 tdtdt0.ref().primitiveFieldRef() = (-rDeltaT.value()*dt.value())*
167 (
168 (coefft0*mesh().S0() - coefft00*mesh().S00())/mesh().S()
169 );
170 }
172 return tdtdt0;
173}
174
175
176template<class Type>
179(
181)
182{
183 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
184
185 const IOobject ddtIOobject
186 (
187 mesh().thisDb().newIOobject
188 (
189 "ddt("+vf.name()+')',
190 { IOobject::REGISTER }
191 )
192 );
193
194 scalar deltaT = deltaT_();
195 scalar deltaT0 = deltaT0_(vf);
196
197 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
198 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
199 scalar coefft0 = coefft + coefft00;
200
201 if (mesh().moving())
202 {
204 (
206 (
207 ddtIOobject,
208 mesh(),
209 rDeltaT.dimensions()*vf.dimensions(),
210 rDeltaT.value()*
211 (
212 coefft*vf() -
213 (
214 coefft0*vf.oldTime()()*mesh().S0()
215 - coefft00*vf.oldTime().oldTime()()
216 *mesh().S00()
217 )/mesh().S()
218 ),
219 rDeltaT.value()*
220 (
221 coefft*vf.boundaryField() -
222 (
223 coefft0*vf.oldTime().boundaryField()
224 - coefft00*vf.oldTime().oldTime().boundaryField()
225 )
226 )
227 )
228 );
229 }
230 else
231 {
233 (
235 (
236 ddtIOobject,
237 rDeltaT*
238 (
239 coefft*vf
240 - coefft0*vf.oldTime()
241 + coefft00*vf.oldTime().oldTime()
242 )
243 )
244 );
245 }
246}
247
248
249template<class Type>
252(
254)
255{
256 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
257
258 const IOobject ddtIOobject
259 (
260 mesh().thisDb().newIOobject
261 (
262 "ddt0("+vf.name()+')',
263 { IOobject::REGISTER }
264 )
265 );
266
267 scalar deltaT = deltaT_();
268 scalar deltaT0 = deltaT0_(vf);
269
270 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
271 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
272 scalar coefft0 = coefft + coefft00;
273
274 if (mesh().moving())
275 {
277 (
279 (
280 ddtIOobject,
281 mesh(),
282 rDeltaT.dimensions()*vf.dimensions(),
283 rDeltaT.value()*
284 (
285 - (
286 coefft0*vf.oldTime()()*mesh().S0()
287 - coefft00*vf.oldTime().oldTime()()
288 *mesh().S00()
289 )/mesh().S()
290 ),
291 rDeltaT.value()*
292 (
293 - (
294 coefft0*vf.oldTime().boundaryField()
295 - coefft00*vf.oldTime().oldTime().boundaryField()
296 )
297 )
298 )
299 );
300 }
301 else
302 {
303 return tmp<GeometricField<Type, faPatchField, areaMesh>>
304 (
305 new GeometricField<Type, faPatchField, areaMesh>
306 (
307 ddtIOobject,
308 rDeltaT*
309 (
310 - coefft0*vf.oldTime()
311 + coefft00*vf.oldTime().oldTime()
312 )
313 )
314 );
315 }
316}
317
318
319template<class Type>
322(
323 const dimensionedScalar& rho,
325)
326{
327 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
328
329 const IOobject ddtIOobject
330 (
331 mesh().thisDb().newIOobject
332 (
333 "ddt("+rho.name()+','+vf.name()+')',
334 { IOobject::REGISTER }
335 )
336 );
337
338 scalar deltaT = deltaT_();
339 scalar deltaT0 = deltaT0_(vf);
340
341 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
342 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
343 scalar coefft0 = coefft + coefft00;
344
345 if (mesh().moving())
346 {
348 (
350 (
351 ddtIOobject,
352 mesh(),
353 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
354 rDeltaT.value()*rho.value()*
355 (
356 coefft*vf.internalField() -
357 (
358 coefft0*vf.oldTime()()*mesh().S0()
359 - coefft00*vf.oldTime().oldTime()()
360 *mesh().S00()
361 )/mesh().S()
362 ),
363 rDeltaT.value()*rho.value()*
364 (
365 coefft*vf.boundaryField() -
366 (
367 coefft0*vf.oldTime().boundaryField()
368 - coefft00*vf.oldTime().oldTime().boundaryField()
369 )
370 )
371 )
372 );
373 }
374 else
375 {
376 return tmp<GeometricField<Type, faPatchField, areaMesh>>
377 (
378 new GeometricField<Type, faPatchField, areaMesh>
379 (
380 ddtIOobject,
381 rDeltaT*rho*
382 (
383 coefft*vf
384 - coefft0*vf.oldTime()
385 + coefft00*vf.oldTime().oldTime()
386 )
388 );
389 }
390}
391
392template<class Type>
395(
396 const dimensionedScalar& rho,
398)
399{
400 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
401
402 const IOobject ddtIOobject
403 (
404 mesh().thisDb().newIOobject
405 (
406 "ddt0("+rho.name()+','+vf.name()+')',
407 { IOobject::REGISTER }
408 )
409 );
410
411 scalar deltaT = deltaT_();
412 scalar deltaT0 = deltaT0_(vf);
413
414 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
415 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
416 scalar coefft0 = coefft + coefft00;
417
418 if (mesh().moving())
419 {
421 (
423 (
424 ddtIOobject,
425 mesh(),
426 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
427 rDeltaT.value()*rho.value()*
428 (
429 -(
430 coefft0*vf.oldTime()()*mesh().S0()
431 - coefft00*vf.oldTime().oldTime()()
432 *mesh().S00()
433 )/mesh().S()
434 ),
435 rDeltaT.value()*rho.value()*
436 (
437 -(
438 coefft0*vf.oldTime().boundaryField()
439 - coefft00*vf.oldTime().oldTime().boundaryField()
440 )
441 )
442 )
443 );
444 }
445 else
446 {
447 return tmp<GeometricField<Type, faPatchField, areaMesh>>
448 (
449 new GeometricField<Type, faPatchField, areaMesh>
450 (
451 ddtIOobject,
452 rDeltaT*rho*
453 (
454 - coefft0*vf.oldTime()
455 + coefft00*vf.oldTime().oldTime()
456 )
457 )
458 );
459 }
460}
461
462
463template<class Type>
466(
467 const areaScalarField& rho,
469)
470{
471 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
472
473 const IOobject ddtIOobject
474 (
475 mesh().thisDb().newIOobject
476 (
477 "ddt("+rho.name()+','+vf.name()+')',
478 { IOobject::REGISTER }
479 )
480 );
481
482 scalar deltaT = deltaT_();
483 scalar deltaT0 = deltaT0_(vf);
484
485 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
486 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
487 scalar coefft0 = coefft + coefft00;
488
489 if (mesh().moving())
490 {
492 (
494 (
495 ddtIOobject,
496 mesh(),
497 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
498 rDeltaT.value()*
499 (
500 coefft*rho.internalField()*vf.internalField() -
501 (
502 coefft0*rho.oldTime()()
503 *vf.oldTime()()*mesh().S0()
504 - coefft00*rho.oldTime().oldTime()()
505 *vf.oldTime().oldTime()()*mesh().S00()
506 )/mesh().S()
507 ),
508 rDeltaT.value()*
509 (
510 coefft*vf.boundaryField() -
511 (
512 coefft0*rho.oldTime().boundaryField()
513 *vf.oldTime().boundaryField()
514 - coefft00*rho.oldTime().oldTime().boundaryField()
515 *vf.oldTime().oldTime().boundaryField()
516 )
517 )
518 )
519 );
520 }
521 else
522 {
523 return tmp<GeometricField<Type, faPatchField, areaMesh>>
524 (
525 new GeometricField<Type, faPatchField, areaMesh>
526 (
527 ddtIOobject,
528 rDeltaT*
529 (
530 coefft*rho*vf
531 - coefft0*rho.oldTime()*vf.oldTime()
532 + coefft00*rho.oldTime().oldTime()*vf.oldTime().oldTime()
533 )
534 )
535 );
536 }
537}
538
539
540template<class Type>
543(
544 const areaScalarField& rho,
546)
547{
548 dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
549
550 const IOobject ddtIOobject
551 (
552 mesh().thisDb().newIOobject
553 (
554 "ddt0("+rho.name()+','+vf.name()+')',
555 { IOobject::REGISTER }
556 )
557 );
558
559 scalar deltaT = deltaT_();
560 scalar deltaT0 = deltaT0_(vf);
561
562 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
563 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
564 scalar coefft0 = coefft + coefft00;
565
566 if (mesh().moving())
567 {
569 (
571 (
572 ddtIOobject,
573 mesh(),
574 rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
575 rDeltaT.value()*
576 (
577 - (
578 coefft0*rho.oldTime()()
579 *vf.oldTime()()*mesh().S0()
580 - coefft00*rho.oldTime().oldTime()()
581 *vf.oldTime().oldTime()()*mesh().S00()
582 )/mesh().S()
583 ),
584 rDeltaT.value()*
585 (
586 - (
587 coefft0*rho.oldTime().boundaryField()
588 *vf.oldTime().boundaryField()
589 - coefft00*rho.oldTime().oldTime().boundaryField()
590 *vf.oldTime().oldTime().boundaryField()
591 )
592 )
593 )
594 );
595 }
596 else
597 {
598 return tmp<GeometricField<Type, faPatchField, areaMesh>>
599 (
600 new GeometricField<Type, faPatchField, areaMesh>
601 (
602 ddtIOobject,
603 rDeltaT*
604 (
605 - coefft0*rho.oldTime()*vf.oldTime()
606 + coefft00*rho.oldTime().oldTime()*vf.oldTime().oldTime()
607 )
608 )
609 );
610 }
611}
612
613
614template<class Type>
617(
619)
620{
622 (
624 (
625 vf,
627 )
628 );
629
630 faMatrix<Type>& fam = tfam.ref();
631
632 scalar rDeltaT = 1.0/deltaT_();
633
634 scalar deltaT = deltaT_();
635 scalar deltaT0 = deltaT0_(vf);
636
637 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
638 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
639 scalar coefft0 = coefft + coefft00;
640
641 fam.diag() = (coefft*rDeltaT)*mesh().S();
642
643 if (mesh().moving())
644 {
645 fam.source() = rDeltaT*
646 (
647 coefft0*vf.oldTime().primitiveField()*mesh().S0()
648 - coefft00*vf.oldTime().oldTime().primitiveField()
649 *mesh().S00()
650 );
651 }
652 else
653 {
654 fam.source() = rDeltaT*mesh().S()*
655 (
656 coefft0*vf.oldTime().primitiveField()
657 - coefft00*vf.oldTime().oldTime().primitiveField()
658 );
659 }
661 return tfam;
662}
663
664
665template<class Type>
668(
669 const dimensionedScalar& rho,
671)
672{
674 (
676 (
677 vf,
678 rho.dimensions()*vf.dimensions()*dimArea/dimTime
679 )
680 );
681 faMatrix<Type>& fam = tfam.ref();
682
683 scalar rDeltaT = 1.0/deltaT_();
684
685 scalar deltaT = deltaT_();
686 scalar deltaT0 = deltaT0_(vf);
687
688 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
689 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
690 scalar coefft0 = coefft + coefft00;
691
692 fam.diag() = (coefft*rDeltaT*rho.value())*mesh().S();
693
694 if (mesh().moving())
695 {
696 fam.source() = rDeltaT*rho.value()*
697 (
698 coefft0*vf.oldTime().primitiveField()*mesh().S0()
699 - coefft00*vf.oldTime().oldTime().primitiveField()
700 *mesh().S00()
701 );
702 }
703 else
704 {
705 fam.source() = rDeltaT*mesh().S()*rho.value()*
706 (
707 coefft0*vf.oldTime().primitiveField()
708 - coefft00*vf.oldTime().oldTime().primitiveField()
709 );
710 }
712 return tfam;
713}
714
715
716template<class Type>
719(
720 const areaScalarField& rho,
722)
723{
725 (
727 (
728 vf,
729 rho.dimensions()*vf.dimensions()*dimArea/dimTime
730 )
731 );
732 faMatrix<Type>& fam = tfam.ref();
733
734 scalar rDeltaT = 1.0/deltaT_();
735
736 scalar deltaT = deltaT_();
737 scalar deltaT0 = deltaT0_(vf);
738
739 scalar coefft = 1 + deltaT/(deltaT + deltaT0);
740 scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
741 scalar coefft0 = coefft + coefft00;
742
743 fam.diag() = (coefft*rDeltaT)*rho.primitiveField()*mesh().S();
744
745 if (mesh().moving())
746 {
747 fam.source() = rDeltaT*
748 (
749 coefft0*rho.oldTime().primitiveField()
750 *vf.oldTime().primitiveField()*mesh().S0()
751 - coefft00*rho.oldTime().oldTime().primitiveField()
752 *vf.oldTime().oldTime().primitiveField()*mesh().S00()
753 );
754 }
755 else
756 {
757 fam.source() = rDeltaT*mesh().S()*
758 (
759 coefft0*rho.oldTime().primitiveField()
760 *vf.oldTime().primitiveField()
761 - coefft00*rho.oldTime().oldTime().primitiveField()
762 *vf.oldTime().oldTime().primitiveField()
763 );
764 }
765
766 return tfam;
767}
768
769
770// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
771
772} // End namespace fa
773
774// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
775
776} // End namespace Foam
777
778// ************************************************************************* //
const dimensionSet & dimensions() const noexcept
Return dimensions.
Generic GeometricField class.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
const Internal & internalField() const noexcept
Return a const-reference to the dimensioned internal field.
const Boundary & boundaryField() const noexcept
Return const-reference to the boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition IOobject.H:191
const word & name() const noexcept
Return the object name.
Definition IOobjectI.H:205
Generic dimensioned Type class.
const dimensionSet & dimensions() const noexcept
Return const reference to dimensions.
const word & name() const noexcept
Return const reference to name.
const Type & value() const noexcept
Return const reference to value.
A special matrix type and solver, designed for finite area solutions of scalar equations....
Definition faMatrix.H:108
static const word & calculatedType() noexcept
The type name for calculated patch fields.
tmp< GeometricField< Type, faPatchField, areaMesh > > facDdt(const dimensioned< Type >)
tmp< faMatrix< Type > > famDdt(const GeometricField< Type, faPatchField, areaMesh > &)
const faMesh & mesh() const
Return mesh reference.
tmp< GeometricField< Type, faPatchField, areaMesh > > facDdt0(const dimensioned< Type >)
A class for managing temporary objects.
Definition tmp.H:75
T & ref() const
Return non-const reference to the contents of a non-null managed pointer.
Definition tmpI.H:235
dynamicFvMesh & mesh
Namespace for finite-area.
Definition limitHeight.C:30
Namespace for OpenFOAM.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
const dimensionSet dimArea(sqr(dimLength))
GeometricField< scalar, faPatchField, areaMesh > areaScalarField
static constexpr const zero Zero
Global zero (0).
Definition zero.H:127
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Calculate the matrix for the second temporal derivative.