CrankNicolsonDdtScheme.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 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 "CrankNicolsonDdtScheme.H"
27 #include "surfaceInterpolate.H"
28 #include "fvcDiv.H"
29 #include "fvMatrices.H"
30 #include "Constant.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 namespace fv
40 {
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 template<class Type>
45 template<class GeoField>
46 CrankNicolsonDdtScheme<Type>::DDt0Field<GeoField>::DDt0Field
47 (
48  const IOobject& io,
49  const fvMesh& mesh
50 )
51 :
52  GeoField(io, mesh),
53  startTimeIndex_(-2) // This field is for a restart and thus correct so set
54  // the start time-index to correspond to a previous run
55 {
56  // Set the time-index to the beginning of the run to ensure the field
57  // is updated during the first time-step
58  this->timeIndex() = mesh.time().startTimeIndex();
59 }
60 
61 
62 template<class Type>
63 template<class GeoField>
64 CrankNicolsonDdtScheme<Type>::DDt0Field<GeoField>::DDt0Field
65 (
66  const IOobject& io,
67  const fvMesh& mesh,
68  const dimensioned<typename GeoField::value_type>& dimType
69 )
70 :
71  GeoField(io, mesh, dimType),
72  startTimeIndex_(mesh.time().timeIndex())
73 {}
74 
75 
76 template<class Type>
77 template<class GeoField>
78 label CrankNicolsonDdtScheme<Type>::DDt0Field<GeoField>::
79 startTimeIndex() const
80 {
81  return startTimeIndex_;
82 }
83 
84 
85 template<class Type>
86 template<class GeoField>
87 GeoField& CrankNicolsonDdtScheme<Type>::DDt0Field<GeoField>::
88 operator()()
89 {
90  return *this;
91 }
92 
93 
94 template<class Type>
95 template<class GeoField>
96 void CrankNicolsonDdtScheme<Type>::DDt0Field<GeoField>::
97 operator=(const GeoField& gf)
98 {
99  GeoField::operator=(gf);
100 }
101 
102 
103 template<class Type>
104 template<class GeoField>
105 typename CrankNicolsonDdtScheme<Type>::template DDt0Field<GeoField>&
106 CrankNicolsonDdtScheme<Type>::ddt0_
107 (
108  const word& name,
109  const dimensionSet& dims
110 )
111 {
112  if (!mesh().objectRegistry::template foundObject<GeoField>(name))
113  {
114  const Time& runTime = mesh().time();
115  word startTimeName = runTime.timeName(runTime.startTime().value());
116 
117  if
118  (
119  (
120  runTime.timeIndex() == runTime.startTimeIndex()
121  || runTime.timeIndex() == runTime.startTimeIndex() + 1
122  )
123  && IOobject
124  (
125  name,
126  startTimeName,
127  mesh()
128  ).typeHeaderOk<DDt0Field<GeoField>>()
129  )
130  {
132  (
133  new DDt0Field<GeoField>
134  (
135  IOobject
136  (
137  name,
138  startTimeName,
139  mesh(),
142  ),
143  mesh()
144  )
145  );
146  }
147  else
148  {
150  (
151  new DDt0Field<GeoField>
152  (
153  IOobject
154  (
155  name,
156  mesh().time().timeName(),
157  mesh(),
160  ),
161  mesh(),
163  (
164  "0",
165  dims/dimTime,
166  Zero
167  )
168  )
169  );
170  }
171  }
172 
173  DDt0Field<GeoField>& ddt0 = static_cast<DDt0Field<GeoField>&>
174  (
175  mesh().objectRegistry::template lookupObjectRef<GeoField>(name)
176  );
177 
178  return ddt0;
179 }
180 
181 
182 template<class Type>
183 template<class GeoField>
185 (
186  const DDt0Field<GeoField>& ddt0
187 ) const
188 {
189  return ddt0.timeIndex() != mesh().time().timeIndex();
190 }
191 
192 template<class Type>
193 template<class GeoField>
195 (
196  const DDt0Field<GeoField>& ddt0
197 ) const
198 {
199  if (mesh().time().timeIndex() > ddt0.startTimeIndex())
200  {
201  return 1 + ocCoeff();
202  }
203  else
204  {
205  return 1;
206  }
207 }
208 
209 
210 template<class Type>
211 template<class GeoField>
213 (
214  const DDt0Field<GeoField>& ddt0
215 ) const
216 {
217  if (mesh().time().timeIndex() > ddt0.startTimeIndex() + 1)
218  {
219  return 1 + ocCoeff();
220  }
221  else
222  {
223  return 1;
224  }
225 }
226 
227 
228 template<class Type>
229 template<class GeoField>
231 (
232  const DDt0Field<GeoField>& ddt0
233 ) const
234 {
235  return coef_(ddt0)/mesh().time().deltaT();
236 }
237 
238 
239 template<class Type>
240 template<class GeoField>
242 (
243  const DDt0Field<GeoField>& ddt0
244 ) const
245 {
246  return coef0_(ddt0)/mesh().time().deltaT0();
247 }
248 
249 
250 template<class Type>
251 template<class GeoField>
253 (
254  const GeoField& ddt0
255 ) const
256 {
257  if (ocCoeff() < 1)
258  {
259  return ocCoeff()*ddt0;
260  }
261  else
262  {
263  return ddt0;
264  }
265 }
266 
267 
268 template<class Type>
270 (
272 )
273 {
274  return bf;
275 }
276 
277 
278 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
279 
280 template<class Type>
282 :
283  ddtScheme<Type>(mesh),
284  ocCoeff_(new Function1Types::Constant<scalar>("ocCoeff", 1))
285 {
286  // Ensure the old-old-time cell volumes are available
287  // for moving meshes
288  if (mesh.moving())
289  {
290  mesh.V00();
291  }
292 }
293 
294 
295 template<class Type>
297 (
298  const fvMesh& mesh,
299  Istream& is
300 )
301 :
302  ddtScheme<Type>(mesh, is)
303 {
304  token firstToken(is);
305 
306  if (firstToken.isNumber())
307  {
308  const scalar ocCoeff = firstToken.number();
309  if (ocCoeff < 0 || ocCoeff > 1)
310  {
312  (
313  is
314  ) << "Off-centreing coefficient = " << ocCoeff
315  << " should be >= 0 and <= 1"
316  << exit(FatalIOError);
317  }
318 
319  ocCoeff_ = new Function1Types::Constant<scalar>
320  (
321  "ocCoeff",
322  ocCoeff
323  );
324  }
325  else
326  {
327  is.putBack(firstToken);
328  dictionary dict(is);
329  ocCoeff_ = Function1<scalar>::New("ocCoeff", dict);
330  }
331 
332  // Ensure the old-old-time cell volumes are available
333  // for moving meshes
334  if (mesh.moving())
335  {
336  mesh.V00();
337  }
338 }
339 
340 
341 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
342 
343 template<class Type>
346 (
347  const dimensioned<Type>& dt
348 )
349 {
350  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
351  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
352  (
353  "ddt0(" + dt.name() + ')',
354  dt.dimensions()
355  );
356 
357  IOobject ddtIOobject
358  (
359  "ddt(" + dt.name() + ')',
360  mesh().time().timeName(),
361  mesh()
362  );
363 
365  (
367  (
368  ddtIOobject,
369  mesh(),
371  (
372  "0",
373  dt.dimensions()/dimTime,
374  Zero
375  )
376  )
377  );
378 
379  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
380 
381  if (mesh().moving())
382  {
383  if (evaluate(ddt0))
384  {
385  dimensionedScalar rDtCoef0 = rDtCoef0_(ddt0);
386 
387  ddt0.ref() =
388  (
389  (rDtCoef0*dt)*(mesh().V0() - mesh().V00())
390  - mesh().V00()*offCentre_(ddt0.internalField())
391  )/mesh().V0();
392  }
393 
394  tdtdt.ref().ref() =
395  (
396  (rDtCoef*dt)*(mesh().V() - mesh().V0())
397  - mesh().V0()*offCentre_(ddt0.internalField())
398  )/mesh().V();
399  }
400 
401  return tdtdt;
402 }
403 
404 
405 template<class Type>
408 (
410 )
411 {
412  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
413  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
414  (
415  "ddt0(" + vf.name() + ')',
416  vf.dimensions()
417  );
418 
419  IOobject ddtIOobject
420  (
421  "ddt(" + vf.name() + ')',
422  mesh().time().timeName(),
423  mesh()
424  );
425 
426  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
427 
428  if (mesh().moving())
429  {
430  if (evaluate(ddt0))
431  {
432  const scalar rDtCoef0 = rDtCoef0_(ddt0).value();
433 
434  ddt0.primitiveFieldRef() =
435  (
436  rDtCoef0*
437  (
438  mesh().V0()*vf.oldTime().primitiveField()
439  - mesh().V00()*vf.oldTime().oldTime().primitiveField()
440  ) - mesh().V00()*offCentre_(ddt0.primitiveField())
441  )/mesh().V0();
442 
443  ddt0.boundaryFieldRef() =
444  (
445  rDtCoef0*
446  (
447  vf.oldTime().boundaryField()
448  - vf.oldTime().oldTime().boundaryField()
449  ) - offCentre_(ff(ddt0.boundaryField()))
450  );
451  }
452 
454  (
456  (
457  ddtIOobject,
458  (
459  rDtCoef*
460  (
461  mesh().V()*vf
462  - mesh().V0()*vf.oldTime()
463  ) - mesh().V0()*offCentre_(ddt0()())
464  )/mesh().V(),
465  rDtCoef.value()*
466  (
467  vf.boundaryField() - vf.oldTime().boundaryField()
468  ) - offCentre_(ff(ddt0.boundaryField()))
469  )
470  );
471  }
472  else
473  {
474  if (evaluate(ddt0))
475  {
476  ddt0 = rDtCoef0_(ddt0)*(vf.oldTime() - vf.oldTime().oldTime())
477  - offCentre_(ddt0());
478  }
479 
481  (
483  (
484  ddtIOobject,
485  rDtCoef*(vf - vf.oldTime()) - offCentre_(ddt0())
486  )
487  );
488  }
489 }
490 
491 
492 template<class Type>
495 (
496  const dimensionedScalar& rho,
498 )
499 {
500  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
501  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
502  (
503  "ddt0(" + rho.name() + ',' + vf.name() + ')',
504  rho.dimensions()*vf.dimensions()
505  );
506 
507  IOobject ddtIOobject
508  (
509  "ddt(" + rho.name() + ',' + vf.name() + ')',
510  mesh().time().timeName(),
511  mesh()
512  );
513 
514  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
515 
516  if (mesh().moving())
517  {
518  if (evaluate(ddt0))
519  {
520  const scalar rDtCoef0 = rDtCoef0_(ddt0).value();
521 
522  ddt0.primitiveFieldRef() =
523  (
524  rDtCoef0*rho.value()*
525  (
526  mesh().V0()*vf.oldTime().primitiveField()
527  - mesh().V00()*vf.oldTime().oldTime().primitiveField()
528  ) - mesh().V00()*offCentre_(ddt0.primitiveField())
529  )/mesh().V0();
530 
531  ddt0.boundaryFieldRef() =
532  (
533  rDtCoef0*rho.value()*
534  (
535  vf.oldTime().boundaryField()
536  - vf.oldTime().oldTime().boundaryField()
537  ) - offCentre_(ff(ddt0.boundaryField()))
538  );
539  }
540 
542  (
544  (
545  ddtIOobject,
546  mesh(),
547  rDtCoef.dimensions()*rho.dimensions()*vf.dimensions(),
548  (
549  rDtCoef.value()*rho.value()*
550  (
551  mesh().V()*vf.primitiveField()
552  - mesh().V0()*vf.oldTime().primitiveField()
553  ) - mesh().V0()*offCentre_(ddt0.primitiveField())
554  )/mesh().V(),
555  rDtCoef.value()*rho.value()*
556  (
557  vf.boundaryField() - vf.oldTime().boundaryField()
558  ) - offCentre_(ff(ddt0.boundaryField()))
559  )
560  );
561  }
562  else
563  {
564  if (evaluate(ddt0))
565  {
566  ddt0 = rDtCoef0_(ddt0)*rho*(vf.oldTime() - vf.oldTime().oldTime())
567  - offCentre_(ddt0());
568  }
569 
571  (
573  (
574  ddtIOobject,
575  rDtCoef*rho*(vf - vf.oldTime()) - offCentre_(ddt0())
576  )
577  );
578  }
579 }
580 
581 
582 template<class Type>
585 (
586  const volScalarField& rho,
588 )
589 {
590  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
591  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
592  (
593  "ddt0(" + rho.name() + ',' + vf.name() + ')',
594  rho.dimensions()*vf.dimensions()
595  );
596 
597  IOobject ddtIOobject
598  (
599  "ddt(" + rho.name() + ',' + vf.name() + ')',
600  mesh().time().timeName(),
601  mesh()
602  );
603 
604  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
605 
606  if (mesh().moving())
607  {
608  if (evaluate(ddt0))
609  {
610  const scalar rDtCoef0 = rDtCoef0_(ddt0).value();
611 
612  ddt0.primitiveFieldRef() =
613  (
614  rDtCoef0*
615  (
616  mesh().V0()*rho.oldTime().primitiveField()
617  *vf.oldTime().primitiveField()
618  - mesh().V00()*rho.oldTime().oldTime().primitiveField()
619  *vf.oldTime().oldTime().primitiveField()
620  ) - mesh().V00()*offCentre_(ddt0.primitiveField())
621  )/mesh().V0();
622 
623  ddt0.boundaryFieldRef() =
624  (
625  rDtCoef0*
626  (
627  rho.oldTime().boundaryField()
628  *vf.oldTime().boundaryField()
629  - rho.oldTime().oldTime().boundaryField()
630  *vf.oldTime().oldTime().boundaryField()
631  ) - offCentre_(ff(ddt0.boundaryField()))
632  );
633  }
634 
636  (
638  (
639  ddtIOobject,
640  mesh(),
641  rDtCoef.dimensions()*rho.dimensions()*vf.dimensions(),
642  (
643  rDtCoef.value()*
644  (
645  mesh().V()*rho.primitiveField()*vf.primitiveField()
646  - mesh().V0()*rho.oldTime().primitiveField()
647  *vf.oldTime().primitiveField()
648  ) - mesh().V00()*offCentre_(ddt0.primitiveField())
649  )/mesh().V(),
650  rDtCoef.value()*
651  (
652  rho.boundaryField()*vf.boundaryField()
653  - rho.oldTime().boundaryField()*vf.oldTime().boundaryField()
654  ) - offCentre_(ff(ddt0.boundaryField()))
655  )
656  );
657  }
658  else
659  {
660  if (evaluate(ddt0))
661  {
662  ddt0 = rDtCoef0_(ddt0)*
663  (
664  rho.oldTime()*vf.oldTime()
665  - rho.oldTime().oldTime()*vf.oldTime().oldTime()
666  ) - offCentre_(ddt0());
667  }
668 
670  (
672  (
673  ddtIOobject,
674  rDtCoef*(rho*vf - rho.oldTime()*vf.oldTime())
675  - offCentre_(ddt0())
676  )
677  );
678  }
679 }
680 
681 
682 template<class Type>
685 (
686  const volScalarField& alpha,
687  const volScalarField& rho,
689 )
690 {
691  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
692  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
693  (
694  "ddt0(" + alpha.name() + ',' + rho.name() + ',' + vf.name() + ')',
695  alpha.dimensions()*rho.dimensions()*vf.dimensions()
696  );
697 
698  IOobject ddtIOobject
699  (
700  "ddt(" + alpha.name() + ',' + rho.name() + ',' + vf.name() + ')',
701  mesh().time().timeName(),
702  mesh()
703  );
704 
705  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
706 
707  if (mesh().moving())
708  {
709  if (evaluate(ddt0))
710  {
711  const scalar rDtCoef0 = rDtCoef0_(ddt0).value();
712 
713  ddt0.primitiveFieldRef() =
714  (
715  rDtCoef0*
716  (
717  mesh().V0()
718  *alpha.oldTime().primitiveField()
719  *rho.oldTime().primitiveField()
720  *vf.oldTime().primitiveField()
721 
722  - mesh().V00()
723  *alpha.oldTime().oldTime().primitiveField()
724  *rho.oldTime().oldTime().primitiveField()
725  *vf.oldTime().oldTime().primitiveField()
726  ) - mesh().V00()*offCentre_(ddt0.primitiveField())
727  )/mesh().V0();
728 
729  ddt0.boundaryFieldRef() =
730  (
731  rDtCoef0*
732  (
733  alpha.oldTime().boundaryField()
734  *rho.oldTime().boundaryField()
735  *vf.oldTime().boundaryField()
736 
737  - alpha.oldTime().oldTime().boundaryField()
738  *rho.oldTime().oldTime().boundaryField()
739  *vf.oldTime().oldTime().boundaryField()
740  ) - offCentre_(ff(ddt0.boundaryField()))
741  );
742  }
743 
745  (
747  (
748  ddtIOobject,
749  mesh(),
750  rDtCoef.dimensions()
751  *alpha.dimensions()*rho.dimensions()*vf.dimensions(),
752  (
753  rDtCoef.value()*
754  (
755  mesh().V()
756  *alpha.primitiveField()
757  *rho.primitiveField()
758  *vf.primitiveField()
759 
760  - mesh().V0()
761  *alpha.oldTime().primitiveField()
762  *rho.oldTime().primitiveField()
763  *vf.oldTime().primitiveField()
764  ) - mesh().V00()*offCentre_(ddt0.primitiveField())
765  )/mesh().V(),
766  rDtCoef.value()*
767  (
768  alpha.boundaryField()
769  *rho.boundaryField()
770  *vf.boundaryField()
771 
772  - alpha.oldTime().boundaryField()
773  *rho.oldTime().boundaryField()
774  *vf.oldTime().boundaryField()
775  ) - offCentre_(ff(ddt0.boundaryField()))
776  )
777  );
778  }
779  else
780  {
781  if (evaluate(ddt0))
782  {
783  ddt0 = rDtCoef0_(ddt0)*
784  (
785  alpha.oldTime()
786  *rho.oldTime()
787  *vf.oldTime()
788 
789  - alpha.oldTime().oldTime()
790  *rho.oldTime().oldTime()
791  *vf.oldTime().oldTime()
792  ) - offCentre_(ddt0());
793  }
794 
796  (
798  (
799  ddtIOobject,
800  rDtCoef
801  *(
802  alpha*rho*vf
803  - alpha.oldTime()*rho.oldTime()*vf.oldTime()
804  )
805  - offCentre_(ddt0())
806  )
807  );
808  }
809 }
810 
811 
812 template<class Type>
815 (
817 )
818 {
819  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
820  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
821  (
822  "ddt0(" + vf.name() + ')',
823  vf.dimensions()
824  );
825 
826  tmp<fvMatrix<Type>> tfvm
827  (
828  new fvMatrix<Type>
829  (
830  vf,
832  )
833  );
834 
835  fvMatrix<Type>& fvm = tfvm.ref();
836 
837  const scalar rDtCoef = rDtCoef_(ddt0).value();
838  fvm.diag() = rDtCoef*mesh().V();
839 
840  vf.oldTime().oldTime();
841 
842  if (mesh().moving())
843  {
844  if (evaluate(ddt0))
845  {
846  const scalar rDtCoef0 = rDtCoef0_(ddt0).value();
847 
848  ddt0.primitiveFieldRef() =
849  (
850  rDtCoef0*
851  (
852  mesh().V0()*vf.oldTime().primitiveField()
853  - mesh().V00()*vf.oldTime().oldTime().primitiveField()
854  )
855  - mesh().V00()*offCentre_(ddt0.primitiveField())
856  )/mesh().V0();
857 
858  ddt0.boundaryFieldRef() =
859  (
860  rDtCoef0*
861  (
862  vf.oldTime().boundaryField()
863  - vf.oldTime().oldTime().boundaryField()
864  )
865  - offCentre_(ff(ddt0.boundaryField()))
866  );
867  }
868 
869  fvm.source() =
870  (
871  rDtCoef*vf.oldTime().primitiveField()
872  + offCentre_(ddt0.primitiveField())
873  )*mesh().V0();
874  }
875  else
876  {
877  if (evaluate(ddt0))
878  {
879  ddt0 = rDtCoef0_(ddt0)*(vf.oldTime() - vf.oldTime().oldTime())
880  - offCentre_(ddt0());
881  }
882 
883  fvm.source() =
884  (
885  rDtCoef*vf.oldTime().primitiveField()
886  + offCentre_(ddt0.primitiveField())
887  )*mesh().V();
888  }
889 
890  return tfvm;
891 }
892 
893 
894 template<class Type>
897 (
898  const dimensionedScalar& rho,
900 )
901 {
902  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
903  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
904  (
905  "ddt0(" + rho.name() + ',' + vf.name() + ')',
906  rho.dimensions()*vf.dimensions()
907  );
908 
909  tmp<fvMatrix<Type>> tfvm
910  (
911  new fvMatrix<Type>
912  (
913  vf,
915  )
916  );
917  fvMatrix<Type>& fvm = tfvm.ref();
918 
919  const scalar rDtCoef = rDtCoef_(ddt0).value();
920  fvm.diag() = rDtCoef*rho.value()*mesh().V();
921 
922  vf.oldTime().oldTime();
923 
924  if (mesh().moving())
925  {
926  if (evaluate(ddt0))
927  {
928  const scalar rDtCoef0 = rDtCoef0_(ddt0).value();
929 
930  ddt0.primitiveFieldRef() =
931  (
932  rDtCoef0*rho.value()*
933  (
934  mesh().V0()*vf.oldTime().primitiveField()
935  - mesh().V00()*vf.oldTime().oldTime().primitiveField()
936  )
937  - mesh().V00()*offCentre_(ddt0.primitiveField())
938  )/mesh().V0();
939 
940  ddt0.boundaryFieldRef() =
941  (
942  rDtCoef0*rho.value()*
943  (
944  vf.oldTime().boundaryField()
945  - vf.oldTime().oldTime().boundaryField()
946  )
947  - offCentre_(ff(ddt0.boundaryField()))
948  );
949  }
950 
951  fvm.source() =
952  (
953  rDtCoef*rho.value()*vf.oldTime().primitiveField()
954  + offCentre_(ddt0.primitiveField())
955  )*mesh().V0();
956  }
957  else
958  {
959  if (evaluate(ddt0))
960  {
961  ddt0 = rDtCoef0_(ddt0)*rho*(vf.oldTime() - vf.oldTime().oldTime())
962  - offCentre_(ddt0());
963  }
964 
965  fvm.source() =
966  (
967  rDtCoef*rho.value()*vf.oldTime().primitiveField()
968  + offCentre_(ddt0.primitiveField())
969  )*mesh().V();
970  }
971 
972  return tfvm;
973 }
974 
975 
976 template<class Type>
979 (
980  const volScalarField& rho,
982 )
983 {
984  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
985  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
986  (
987  "ddt0(" + rho.name() + ',' + vf.name() + ')',
988  rho.dimensions()*vf.dimensions()
989  );
990 
991  tmp<fvMatrix<Type>> tfvm
992  (
993  new fvMatrix<Type>
994  (
995  vf,
997  )
998  );
999  fvMatrix<Type>& fvm = tfvm.ref();
1000 
1001  const scalar rDtCoef = rDtCoef_(ddt0).value();
1002  fvm.diag() = rDtCoef*rho.primitiveField()*mesh().V();
1003 
1004  vf.oldTime().oldTime();
1005  rho.oldTime().oldTime();
1006 
1007  if (mesh().moving())
1008  {
1009  if (evaluate(ddt0))
1010  {
1011  const scalar rDtCoef0 = rDtCoef0_(ddt0).value();
1012 
1013  ddt0.primitiveFieldRef() =
1014  (
1015  rDtCoef0*
1016  (
1017  mesh().V0()*rho.oldTime().primitiveField()
1018  *vf.oldTime().primitiveField()
1019  - mesh().V00()*rho.oldTime().oldTime().primitiveField()
1020  *vf.oldTime().oldTime().primitiveField()
1021  )
1022  - mesh().V00()*offCentre_(ddt0.primitiveField())
1023  )/mesh().V0();
1024 
1025  ddt0.boundaryFieldRef() =
1026  (
1027  rDtCoef0*
1028  (
1029  rho.oldTime().boundaryField()
1030  *vf.oldTime().boundaryField()
1031  - rho.oldTime().oldTime().boundaryField()
1032  *vf.oldTime().oldTime().boundaryField()
1033  )
1034  - offCentre_(ff(ddt0.boundaryField()))
1035  );
1036  }
1037 
1038  fvm.source() =
1039  (
1040  rDtCoef*rho.oldTime().primitiveField()*vf.oldTime().primitiveField()
1041  + offCentre_(ddt0.primitiveField())
1042  )*mesh().V0();
1043  }
1044  else
1045  {
1046  if (evaluate(ddt0))
1047  {
1048  ddt0 = rDtCoef0_(ddt0)*
1049  (
1050  rho.oldTime()*vf.oldTime()
1051  - rho.oldTime().oldTime()*vf.oldTime().oldTime()
1052  ) - offCentre_(ddt0());
1053  }
1054 
1055  fvm.source() =
1056  (
1057  rDtCoef*rho.oldTime().primitiveField()*vf.oldTime().primitiveField()
1058  + offCentre_(ddt0.primitiveField())
1059  )*mesh().V();
1060  }
1061 
1062  return tfvm;
1063 }
1064 
1065 
1066 template<class Type>
1070  const volScalarField& alpha,
1071  const volScalarField& rho,
1073 )
1074 {
1075  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
1076  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
1077  (
1078  "ddt0(" + alpha.name() + ',' + rho.name() + ',' + vf.name() + ')',
1079  alpha.dimensions()*rho.dimensions()*vf.dimensions()
1080  );
1081 
1082  tmp<fvMatrix<Type>> tfvm
1083  (
1084  new fvMatrix<Type>
1085  (
1086  vf,
1087  alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
1088  )
1089  );
1090  fvMatrix<Type>& fvm = tfvm.ref();
1091 
1092  const scalar rDtCoef = rDtCoef_(ddt0).value();
1093  fvm.diag() = rDtCoef*alpha.primitiveField()*rho.primitiveField()*mesh().V();
1094 
1095  vf.oldTime().oldTime();
1096  alpha.oldTime().oldTime();
1097  rho.oldTime().oldTime();
1098 
1099  if (mesh().moving())
1100  {
1101  if (evaluate(ddt0))
1102  {
1103  const scalar rDtCoef0 = rDtCoef0_(ddt0).value();
1104 
1105  ddt0.primitiveFieldRef() =
1106  (
1107  rDtCoef0*
1108  (
1109  mesh().V0()
1110  *alpha.oldTime().primitiveField()
1111  *rho.oldTime().primitiveField()
1112  *vf.oldTime().primitiveField()
1113 
1114  - mesh().V00()
1115  *alpha.oldTime().oldTime().primitiveField()
1116  *rho.oldTime().oldTime().primitiveField()
1117  *vf.oldTime().oldTime().primitiveField()
1118  )
1119  - mesh().V00()*offCentre_(ddt0.primitiveField())
1120  )/mesh().V0();
1121 
1122  ddt0.boundaryFieldRef() =
1123  (
1124  rDtCoef0*
1125  (
1126  alpha.oldTime().boundaryField()
1127  *rho.oldTime().boundaryField()
1128  *vf.oldTime().boundaryField()
1129 
1130  - alpha.oldTime().oldTime().boundaryField()
1131  *rho.oldTime().oldTime().boundaryField()
1132  *vf.oldTime().oldTime().boundaryField()
1133  )
1134  - offCentre_(ff(ddt0.boundaryField()))
1135  );
1136  }
1137 
1138  fvm.source() =
1139  (
1140  rDtCoef
1141  *alpha.oldTime().primitiveField()
1142  *rho.oldTime().primitiveField()
1143  *vf.oldTime().primitiveField()
1144  + offCentre_(ddt0.primitiveField())
1145  )*mesh().V0();
1146  }
1147  else
1148  {
1149  if (evaluate(ddt0))
1150  {
1151  ddt0 = rDtCoef0_(ddt0)*
1152  (
1153  alpha.oldTime()
1154  *rho.oldTime()
1155  *vf.oldTime()
1156 
1157  - alpha.oldTime().oldTime()
1158  *rho.oldTime().oldTime()
1159  *vf.oldTime().oldTime()
1160  ) - offCentre_(ddt0());
1161  }
1162 
1163  fvm.source() =
1164  (
1165  rDtCoef
1166  *alpha.oldTime().primitiveField()
1167  *rho.oldTime().primitiveField()
1168  *vf.oldTime().primitiveField()
1169  + offCentre_(ddt0.primitiveField())
1170  )*mesh().V();
1171  }
1172 
1173  return tfvm;
1174 }
1175 
1176 
1177 template<class Type>
1183 )
1184 {
1185  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
1186  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
1187  (
1188  "ddtCorrDdt0(" + U.name() + ')',
1189  U.dimensions()
1190  );
1191 
1192  DDt0Field<GeometricField<Type, fvsPatchField, surfaceMesh>>& dUfdt0 =
1193  ddt0_<GeometricField<Type, fvsPatchField, surfaceMesh>>
1194  (
1195  "ddtCorrDdt0(" + Uf.name() + ')',
1196  Uf.dimensions()
1197  );
1198 
1199  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
1200 
1201  if (evaluate(ddt0))
1202  {
1203  ddt0 =
1204  rDtCoef0_(ddt0)*(U.oldTime() - U.oldTime().oldTime())
1205  - offCentre_(ddt0());
1206  }
1207 
1208  if (evaluate(dUfdt0))
1209  {
1210  dUfdt0 =
1211  rDtCoef0_(dUfdt0)*(Uf.oldTime() - Uf.oldTime().oldTime())
1212  - offCentre_(dUfdt0());
1213  }
1214 
1215  return tmp<fluxFieldType>
1216  (
1217  new fluxFieldType
1218  (
1219  IOobject
1220  (
1221  "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
1222  mesh().time().timeName(),
1223  mesh()
1224  ),
1225  this->fvcDdtPhiCoeff(U.oldTime(), mesh().Sf() & Uf.oldTime())
1226  *(
1227  mesh().Sf()
1228  & (
1229  (rDtCoef*Uf.oldTime() + offCentre_(dUfdt0()))
1230  - fvc::interpolate(rDtCoef*U.oldTime() + offCentre_(ddt0()))
1231  )
1232  )
1233  )
1234  );
1235 }
1236 
1237 
1238 template<class Type>
1243  const fluxFieldType& phi
1244 )
1245 {
1246  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
1247  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
1248  (
1249  "ddtCorrDdt0(" + U.name() + ')',
1250  U.dimensions()
1251  );
1252 
1253  DDt0Field<fluxFieldType>& dphidt0 =
1254  ddt0_<fluxFieldType>
1255  (
1256  "ddtCorrDdt0(" + phi.name() + ')',
1257  phi.dimensions()
1258  );
1259 
1260  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
1261 
1262  if (evaluate(ddt0))
1263  {
1264  ddt0 =
1265  rDtCoef0_(ddt0)*(U.oldTime() - U.oldTime().oldTime())
1266  - offCentre_(ddt0());
1267  }
1268 
1269  if (evaluate(dphidt0))
1270  {
1271  dphidt0 =
1272  rDtCoef0_(dphidt0)*(phi.oldTime() - phi.oldTime().oldTime())
1273  - offCentre_(dphidt0());
1274  }
1275 
1276  return tmp<fluxFieldType>
1277  (
1278  new fluxFieldType
1279  (
1280  IOobject
1281  (
1282  "ddtCorr(" + U.name() + ',' + phi.name() + ')',
1283  mesh().time().timeName(),
1284  mesh()
1285  ),
1286  this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())
1287  *(
1288  (rDtCoef*phi.oldTime() + offCentre_(dphidt0()))
1290  (
1291  mesh().Sf(),
1292  rDtCoef*U.oldTime() + offCentre_(ddt0())
1293  )
1294  )
1295  )
1296  );
1297 }
1298 
1299 
1300 template<class Type>
1304  const volScalarField& rho,
1307 )
1308 {
1309  if
1310  (
1311  U.dimensions() == dimVelocity
1312  && Uf.dimensions() == rho.dimensions()*dimVelocity
1313  )
1314  {
1315  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
1316  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
1317  (
1318  "ddtCorrDdt0(" + rho.name() + ',' + U.name() + ')',
1319  rho.dimensions()*U.dimensions()
1320  );
1321 
1322  DDt0Field<GeometricField<Type, fvsPatchField, surfaceMesh>>& dUfdt0 =
1323  ddt0_<GeometricField<Type, fvsPatchField, surfaceMesh>>
1324  (
1325  "ddtCorrDdt0(" + Uf.name() + ')',
1326  Uf.dimensions()
1327  );
1328 
1329  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
1330 
1332  (
1333  rho.oldTime()*U.oldTime()
1334  );
1335 
1336  if (evaluate(ddt0))
1337  {
1338  ddt0 =
1339  rDtCoef0_(ddt0)
1340  *(rhoU0 - rho.oldTime().oldTime()*U.oldTime().oldTime())
1341  - offCentre_(ddt0());
1342  }
1343 
1344  if (evaluate(dUfdt0))
1345  {
1346  dUfdt0 =
1347  rDtCoef0_(dUfdt0)
1348  *(Uf.oldTime() - Uf.oldTime().oldTime())
1349  - offCentre_(dUfdt0());
1350  }
1351 
1353  (
1354  new fluxFieldType
1355  (
1356  IOobject
1357  (
1358  "ddtCorr("
1359  + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
1360  mesh().time().timeName(),
1361  mesh()
1362  ),
1363  this->fvcDdtPhiCoeff
1364  (
1365  rhoU0,
1366  mesh().Sf() & Uf.oldTime(),
1367  rho.oldTime()
1368  )
1369  *(
1370  mesh().Sf()
1371  & (
1372  (rDtCoef*Uf.oldTime() + offCentre_(dUfdt0()))
1373  - fvc::interpolate(rDtCoef*rhoU0 + offCentre_(ddt0()))
1374  )
1375  )
1376  )
1377  );
1378 
1379  return ddtCorr;
1380  }
1381  else if
1382  (
1383  U.dimensions() == rho.dimensions()*dimVelocity
1384  && Uf.dimensions() == rho.dimensions()*dimVelocity
1385  )
1386  {
1387  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
1388  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
1389  (
1390  "ddtCorrDdt0(" + U.name() + ')',
1391  U.dimensions()
1392  );
1393 
1394  DDt0Field<GeometricField<Type, fvsPatchField, surfaceMesh>>& dUfdt0 =
1395  ddt0_<GeometricField<Type, fvsPatchField, surfaceMesh>>
1396  (
1397  "ddtCorrDdt0(" + Uf.name() + ')',
1398  Uf.dimensions()
1399  );
1400 
1401  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
1402 
1403  if (evaluate(ddt0))
1404  {
1405  ddt0 =
1406  rDtCoef0_(ddt0)*(U.oldTime() - U.oldTime().oldTime())
1407  - offCentre_(ddt0());
1408  }
1409 
1410  if (evaluate(dUfdt0))
1411  {
1412  dUfdt0 =
1413  rDtCoef0_(dUfdt0)*(Uf.oldTime() - Uf.oldTime().oldTime())
1414  - offCentre_(dUfdt0());
1415  }
1416 
1417  return tmp<fluxFieldType>
1418  (
1419  new fluxFieldType
1420  (
1421  IOobject
1422  (
1423  "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
1424  mesh().time().timeName(),
1425  mesh()
1426  ),
1427  this->fvcDdtPhiCoeff
1428  (
1429  U.oldTime(),
1430  mesh().Sf() & Uf.oldTime(),
1431  rho.oldTime()
1432  )
1433  *(
1434  mesh().Sf()
1435  & (
1436  (rDtCoef*Uf.oldTime() + offCentre_(dUfdt0()))
1438  (
1439  rDtCoef*U.oldTime() + offCentre_(ddt0())
1440  )
1441  )
1442  )
1443  )
1444  );
1445  }
1446  else
1447  {
1449  << "dimensions of Uf are not correct"
1450  << abort(FatalError);
1451 
1452  return fluxFieldType::null();
1453  }
1454 }
1455 
1456 
1457 template<class Type>
1461  const volScalarField& rho,
1463  const fluxFieldType& phi
1464 )
1465 {
1466  if
1467  (
1468  U.dimensions() == dimVelocity
1469  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
1470  )
1471  {
1472  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
1473  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
1474  (
1475  "ddtCorrDdt0(" + rho.name() + ',' + U.name() + ')',
1476  rho.dimensions()*U.dimensions()
1477  );
1478 
1479  DDt0Field<fluxFieldType>& dphidt0 =
1480  ddt0_<fluxFieldType>
1481  (
1482  "ddtCorrDdt0(" + phi.name() + ')',
1483  phi.dimensions()
1484  );
1485 
1486  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
1487 
1489  (
1490  rho.oldTime()*U.oldTime()
1491  );
1492 
1493  if (evaluate(ddt0))
1494  {
1495  ddt0 =
1496  rDtCoef0_(ddt0)
1497  *(rhoU0 - rho.oldTime().oldTime()*U.oldTime().oldTime())
1498  - offCentre_(ddt0());
1499  }
1500 
1501  if (evaluate(dphidt0))
1502  {
1503  dphidt0 =
1504  rDtCoef0_(dphidt0)
1505  *(phi.oldTime() - phi.oldTime().oldTime())
1506  - offCentre_(dphidt0());
1507  }
1508 
1510  (
1511  new fluxFieldType
1512  (
1513  IOobject
1514  (
1515  "ddtCorr("
1516  + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
1517  mesh().time().timeName(),
1518  mesh()
1519  ),
1520  this->fvcDdtPhiCoeff(rhoU0, phi.oldTime(), rho.oldTime())
1521  *(
1522  (rDtCoef*phi.oldTime() + offCentre_(dphidt0()))
1524  (
1525  mesh().Sf(),
1526  rDtCoef*rhoU0 + offCentre_(ddt0())
1527  )
1528  )
1529  )
1530  );
1531 
1532  return ddtCorr;
1533  }
1534  else if
1535  (
1536  U.dimensions() == rho.dimensions()*dimVelocity
1537  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
1538  )
1539  {
1540  DDt0Field<GeometricField<Type, fvPatchField, volMesh>>& ddt0 =
1541  ddt0_<GeometricField<Type, fvPatchField, volMesh>>
1542  (
1543  "ddtCorrDdt0(" + U.name() + ')',
1544  U.dimensions()
1545  );
1546 
1547  DDt0Field<fluxFieldType>& dphidt0 =
1548  ddt0_<fluxFieldType>
1549  (
1550  "ddtCorrDdt0(" + phi.name() + ')',
1551  phi.dimensions()
1552  );
1553 
1554  dimensionedScalar rDtCoef = rDtCoef_(ddt0);
1555 
1556  if (evaluate(ddt0))
1557  {
1558  ddt0 =
1559  rDtCoef0_(ddt0)*(U.oldTime() - U.oldTime().oldTime())
1560  - offCentre_(ddt0());
1561  }
1562 
1563  if (evaluate(dphidt0))
1564  {
1565  dphidt0 =
1566  rDtCoef0_(dphidt0)*(phi.oldTime() - phi.oldTime().oldTime())
1567  - offCentre_(dphidt0());
1568  }
1569 
1570  return tmp<fluxFieldType>
1571  (
1572  new fluxFieldType
1573  (
1574  IOobject
1575  (
1576  "ddtCorr(" + U.name() + ',' + phi.name() + ')',
1577  mesh().time().timeName(),
1578  mesh()
1579  ),
1580  this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime(), rho.oldTime())
1581  *(
1582  (rDtCoef*phi.oldTime() + offCentre_(dphidt0()))
1584  (
1585  mesh().Sf(),
1586  rDtCoef*U.oldTime() + offCentre_(ddt0())
1587  )
1588  )
1589  )
1590  );
1591  }
1592  else
1593  {
1595  << "dimensions of phi are not correct"
1596  << abort(FatalError);
1597 
1598  return fluxFieldType::null();
1599  }
1600 }
1601 
1602 
1603 template<class Type>
1607 )
1608 {
1609  DDt0Field<surfaceScalarField>& meshPhi0 = ddt0_<surfaceScalarField>
1610  (
1611  "meshPhiCN_0",
1612  dimVolume
1613  );
1614 
1615  if (evaluate(meshPhi0))
1616  {
1617  meshPhi0 =
1618  coef0_(meshPhi0)*mesh().phi().oldTime() - offCentre_(meshPhi0());
1619  }
1620 
1622  (
1623  new surfaceScalarField
1624  (
1625  IOobject
1626  (
1627  mesh().phi().name(),
1628  mesh().time().timeName(),
1629  mesh(),
1632  false
1633  ),
1634  coef_(meshPhi0)*mesh().phi() - offCentre_(meshPhi0())
1635  )
1636  );
1637 }
1638 
1639 
1640 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1641 
1642 } // End namespace fv
1643 
1644 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1645 
1646 } // End namespace Foam
1647 
1648 // ************************************************************************* //
tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
const fvMesh & mesh() const
Return mesh reference.
dictionary dict
Second-oder Crank-Nicolson implicit ddt using the current and previous time-step fields as well as th...
const surfaceVectorField & Sf() const
Return cell face area vectors.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const word & name() const
Return name.
Definition: IOobject.H:297
dimensionedScalar deltaT0() const
Return old time step.
Definition: TimeStateI.H:59
surfaceScalarField & phi
bool moving() const
Is mesh moving.
Definition: polyMesh.H:502
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const Boundary & boundaryField() const
Return const-reference to the boundary field.
scalar number() const
Definition: tokenI.H:382
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
scalar ocCoeff() const
Return the current off-centreing coefficient.
A token holds items read from Istream.
Definition: token.H:69
static tmp< GeometricField< typename innerProduct< vector, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const surfaceVectorField &Sf, const GeometricField< Type, fvPatchField, volMesh > &tvf)
Interpolate field onto faces.
tmp< GeometricField< typename flux< Type >::type, fvsPatchField, surfaceMesh > > ddtCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Definition: fvcDdt.C:170
virtual dimensionedScalar startTime() const
Return start time.
Definition: Time.C:781
void putBack(const token &)
Put back token.
Definition: Istream.C:30
const surfaceScalarField & phi() const
Return cell face motion fluxes.
engineTime & runTime
const dimensionSet dimVol(dimVolume)
Definition: dimensionSets.H:59
const DimensionedField< scalar, volMesh > & V00() const
Return old-old-time cell volumes.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
Generic dimensioned Type class.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
Generic field type.
Definition: FieldField.H:51
Abstract base class for ddt schemes.
Definition: ddtScheme.H:64
bool isNumber() const
Definition: tokenI.H:377
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
const dimensionSet & dimensions() const
Return dimensions.
Dimension set for the base types.
Definition: dimensionSet.H:120
ddtScheme< Type >::fluxFieldType fluxFieldType
const DimensionedField< scalar, volMesh > & V0() const
Return old-time cell volumes.
A class for handling words, derived from string.
Definition: word.H:59
tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
labelList fv(nPoints)
Templated function that returns a constant value.
Definition: Constant.H:56
const Type & value() const
Return const reference to value.
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:72
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:97
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Calculate the divergence of the given field.
void store()
Transfer ownership of this object to its registry.
Definition: regIOobjectI.H:34
tmp< fvMatrix< Type > > fvmDdt(const GeometricField< Type, fvPatchField, volMesh > &)
Field< Type > & source()
Definition: fvMatrix.H:294
const word & name() const
Return const reference to name.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
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.
static const GeometricField< Type, PatchField, GeoMesh > & null()
Return a null geometric field.
tmp< surfaceScalarField > fvcDdtPhiCoeff(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi, const fluxFieldType &phiCorr)
Definition: ddtScheme.C:153
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label timeIndex
Definition: getTimeIndex.H:4
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
A special matrix type and solver, designed for finite volume solutions of scalar equations.
const dimensionSet & dimensions() const
Return const reference to dimensions.
tmp< surfaceScalarField > meshPhi(const GeometricField< Type, fvPatchField, volMesh > &)
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
scalarField & diag()
Definition: lduMatrix.C:186
virtual label startTimeIndex() const
Return start time index.
Definition: Time.C:775
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:53
static autoPtr< Function1< Type > > New(const word &entryName, const dictionary &dict)
Selector.
Definition: Function1New.C:32
const FieldField< fvPatchField, Type > & ff(const FieldField< fvPatchField, Type > &bf)
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
Namespace for OpenFOAM.
IOerror FatalIOError
const dimensionSet dimVelocity