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