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