EulerDdtScheme.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 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "EulerDdtScheme.H"
29 #include "surfaceInterpolate.H"
30 #include "fvcDiv.H"
31 #include "fvMatrices.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 
38 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
39 
40 namespace fv
41 {
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 template<class Type>
46 tmp<GeometricField<Type, fvPatchField, volMesh>>
48 (
49  const dimensioned<Type>& dt
50 )
51 {
52  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
53 
54  IOobject ddtIOobject
55  (
56  "ddt("+dt.name()+')',
57  mesh().time().timeName(),
58  mesh()
59  );
60 
61  if (mesh().moving())
62  {
64  (
66  (
67  ddtIOobject,
68  mesh(),
70  )
71  );
72 
73  tdtdt.ref().primitiveFieldRef() =
74  rDeltaT.value()*dt.value()*(1.0 - mesh().Vsc0()/mesh().Vsc());
75 
76  return tdtdt;
77  }
78 
80  (
81  ddtIOobject,
82  mesh(),
85  );
86 }
87 
88 
89 template<class Type>
92 (
94 )
95 {
96  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
97 
98  IOobject ddtIOobject
99  (
100  "ddt("+vf.name()+')',
101  mesh().time().timeName(),
102  mesh()
103  );
104 
105  if (mesh().moving())
106  {
108  (
110  (
111  ddtIOobject,
112  rDeltaT*
113  (
114  vf()
115  - vf.oldTime()()*mesh().Vsc0()/mesh().Vsc()
116  ),
117  rDeltaT.value()*
118  (
119  vf.boundaryField() - vf.oldTime().boundaryField()
120  )
121  )
122  );
123  }
124  else
125  {
127  (
129  (
130  ddtIOobject,
131  rDeltaT*(vf - vf.oldTime())
132  )
133  );
134  }
135 }
136 
137 
138 template<class Type>
141 (
142  const dimensionedScalar& rho,
144 )
145 {
146  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
147 
148  IOobject ddtIOobject
149  (
150  "ddt("+rho.name()+','+vf.name()+')',
151  mesh().time().timeName(),
152  mesh()
153  );
154 
155  if (mesh().moving())
156  {
158  (
160  (
161  ddtIOobject,
162  rDeltaT*rho*
163  (
164  vf()
165  - vf.oldTime()()*mesh().Vsc0()/mesh().Vsc()
166  ),
167  rDeltaT.value()*rho.value()*
168  (
169  vf.boundaryField() - vf.oldTime().boundaryField()
170  )
171  )
172  );
173  }
174  else
175  {
177  (
179  (
180  ddtIOobject,
181  rDeltaT*rho*(vf - vf.oldTime())
182  )
183  );
184  }
185 }
186 
187 
188 template<class Type>
191 (
192  const volScalarField& rho,
194 )
195 {
196  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
197 
198  IOobject ddtIOobject
199  (
200  "ddt("+rho.name()+','+vf.name()+')',
201  mesh().time().timeName(),
202  mesh()
203  );
204 
205  if (mesh().moving())
206  {
208  (
210  (
211  ddtIOobject,
212  rDeltaT*
213  (
214  rho()*vf()
215  - rho.oldTime()()
216  *vf.oldTime()()*mesh().Vsc0()/mesh().Vsc()
217  ),
218  rDeltaT.value()*
219  (
220  rho.boundaryField()*vf.boundaryField()
221  - rho.oldTime().boundaryField()
222  *vf.oldTime().boundaryField()
223  )
224  )
225  );
226  }
227  else
228  {
230  (
232  (
233  ddtIOobject,
234  rDeltaT*(rho*vf - rho.oldTime()*vf.oldTime())
235  )
236  );
237  }
238 }
239 
240 
241 template<class Type>
244 (
245  const volScalarField& alpha,
246  const volScalarField& rho,
248 )
249 {
250  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
251 
252  IOobject ddtIOobject
253  (
254  "ddt("+alpha.name()+','+rho.name()+','+vf.name()+')',
255  mesh().time().timeName(),
256  mesh()
257  );
258 
259  if (mesh().moving())
260  {
262  (
264  (
265  ddtIOobject,
266  rDeltaT*
267  (
268  alpha()
269  *rho()
270  *vf()
271 
272  - alpha.oldTime()()
273  *rho.oldTime()()
274  *vf.oldTime()()*mesh().Vsc0()/mesh().Vsc()
275  ),
276  rDeltaT.value()*
277  (
278  alpha.boundaryField()
279  *rho.boundaryField()
280  *vf.boundaryField()
281 
282  - alpha.oldTime().boundaryField()
283  *rho.oldTime().boundaryField()
284  *vf.oldTime().boundaryField()
285  )
286  )
287  );
288  }
289  else
290  {
292  (
294  (
295  ddtIOobject,
296  rDeltaT
297  *(
298  alpha*rho*vf
299  - alpha.oldTime()*rho.oldTime()*vf.oldTime()
300  )
301  )
302  );
303  }
304 }
305 
306 
307 template<class Type>
310 (
312 )
313 {
314  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
315 
316  IOobject ddtIOobject
317  (
318  "ddt("+sf.name()+')',
319  mesh().time().timeName(),
320  mesh()
321  );
322 
324  (
326  (
327  ddtIOobject,
328  rDeltaT*(sf - sf.oldTime())
329  )
330  );
331 }
332 
333 
334 template<class Type>
337 (
339 )
340 {
341  tmp<fvMatrix<Type>> tfvm
342  (
343  new fvMatrix<Type>
344  (
345  vf,
346  vf.dimensions()*dimVol/dimTime
347  )
348  );
349 
350  fvMatrix<Type>& fvm = tfvm.ref();
351 
352  scalar rDeltaT = 1.0/mesh().time().deltaTValue();
353 
354  fvm.diag() = rDeltaT*mesh().Vsc();
355 
356  if (mesh().moving())
357  {
358  fvm.source() = rDeltaT*vf.oldTime().primitiveField()*mesh().Vsc0();
359  }
360  else
361  {
362  fvm.source() = rDeltaT*vf.oldTime().primitiveField()*mesh().Vsc();
363  }
364 
365  return tfvm;
366 }
367 
368 
369 template<class Type>
372 (
373  const dimensionedScalar& rho,
375 )
376 {
377  tmp<fvMatrix<Type>> tfvm
378  (
379  new fvMatrix<Type>
380  (
381  vf,
382  rho.dimensions()*vf.dimensions()*dimVol/dimTime
383  )
384  );
385  fvMatrix<Type>& fvm = tfvm.ref();
386 
387  scalar rDeltaT = 1.0/mesh().time().deltaTValue();
388 
389  fvm.diag() = rDeltaT*rho.value()*mesh().Vsc();
390 
391  if (mesh().moving())
392  {
393  fvm.source() = rDeltaT
394  *rho.value()*vf.oldTime().primitiveField()*mesh().Vsc0();
395  }
396  else
397  {
398  fvm.source() = rDeltaT
399  *rho.value()*vf.oldTime().primitiveField()*mesh().Vsc();
400  }
401 
402  return tfvm;
403 }
404 
405 
406 template<class Type>
409 (
410  const volScalarField& rho,
412 )
413 {
414  tmp<fvMatrix<Type>> tfvm
415  (
416  new fvMatrix<Type>
417  (
418  vf,
419  rho.dimensions()*vf.dimensions()*dimVol/dimTime
420  )
421  );
422  fvMatrix<Type>& fvm = tfvm.ref();
423 
424  scalar rDeltaT = 1.0/mesh().time().deltaTValue();
425 
426  fvm.diag() = rDeltaT*rho.primitiveField()*mesh().Vsc();
427 
428  if (mesh().moving())
429  {
430  fvm.source() = rDeltaT
431  *rho.oldTime().primitiveField()
432  *vf.oldTime().primitiveField()*mesh().Vsc0();
433  }
434  else
435  {
436  fvm.source() = rDeltaT
437  *rho.oldTime().primitiveField()
438  *vf.oldTime().primitiveField()*mesh().Vsc();
439  }
440 
441  return tfvm;
442 }
443 
444 
445 template<class Type>
448 (
449  const volScalarField& alpha,
450  const volScalarField& rho,
452 )
453 {
454  tmp<fvMatrix<Type>> tfvm
455  (
456  new fvMatrix<Type>
457  (
458  vf,
459  alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
460  )
461  );
462  fvMatrix<Type>& fvm = tfvm.ref();
463 
464  scalar rDeltaT = 1.0/mesh().time().deltaTValue();
465 
466  fvm.diag() =
467  rDeltaT*alpha.primitiveField()*rho.primitiveField()*mesh().Vsc();
468 
469  if (mesh().moving())
470  {
471  fvm.source() = rDeltaT
472  *alpha.oldTime().primitiveField()
473  *rho.oldTime().primitiveField()
474  *vf.oldTime().primitiveField()*mesh().Vsc0();
475  }
476  else
477  {
478  fvm.source() = rDeltaT
479  *alpha.oldTime().primitiveField()
480  *rho.oldTime().primitiveField()
481  *vf.oldTime().primitiveField()*mesh().Vsc();
482  }
483 
484  return tfvm;
485 }
486 
487 
488 template<class Type>
491 (
494 )
495 {
496  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
497 
498  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
499  fluxFieldType phiCorr
500  (
501  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
502  );
503 
504  return tmp<fluxFieldType>
505  (
506  new fluxFieldType
507  (
508  IOobject
509  (
510  "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
511  mesh().time().timeName(),
512  mesh()
513  ),
514  this->fvcDdtPhiCoeff(U.oldTime(), phiUf0, phiCorr)
515  *rDeltaT*phiCorr
516  )
517  );
518 }
519 
520 
521 template<class Type>
524 (
526  const fluxFieldType& phi
527 )
528 {
529  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
530 
531  fluxFieldType phiCorr
532  (
533  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
534  );
535 
536  return tmp<fluxFieldType>
537  (
538  new fluxFieldType
539  (
540  IOobject
541  (
542  "ddtCorr(" + U.name() + ',' + phi.name() + ')',
543  mesh().time().timeName(),
544  mesh()
545  ),
546  this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime(), phiCorr)
547  *rDeltaT*phiCorr
548  )
549  );
550 }
551 
552 
553 template<class Type>
556 (
557  const volScalarField& rho,
560 )
561 {
562  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
563 
564  if
565  (
566  U.dimensions() == dimVelocity
567  && Uf.dimensions() == rho.dimensions()*dimVelocity
568  )
569  {
571  (
572  rho.oldTime()*U.oldTime()
573  );
574 
575  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
576  fluxFieldType phiCorr(phiUf0 - fvc::dotInterpolate(mesh().Sf(), rhoU0));
577 
578  return tmp<fluxFieldType>
579  (
580  new fluxFieldType
581  (
582  IOobject
583  (
584  "ddtCorr("
585  + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
586  mesh().time().timeName(),
587  mesh()
588  ),
589  this->fvcDdtPhiCoeff(rhoU0, phiUf0, phiCorr, rho.oldTime())
590  *rDeltaT*phiCorr
591  )
592  );
593  }
594  else if
595  (
596  U.dimensions() == rho.dimensions()*dimVelocity
597  && Uf.dimensions() == rho.dimensions()*dimVelocity
598  )
599  {
600  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
601  fluxFieldType phiCorr
602  (
603  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
604  );
605 
606  return tmp<fluxFieldType>
607  (
608  new fluxFieldType
609  (
610  IOobject
611  (
612  "ddtCorr("
613  + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
614  mesh().time().timeName(),
615  mesh()
616  ),
617  this->fvcDdtPhiCoeff
618  (
619  U.oldTime(),
620  phiUf0,
621  phiCorr,
622  rho.oldTime()
623  )*rDeltaT*phiCorr
624  )
625  );
626  }
627  else
628  {
630  << "dimensions of Uf are not correct"
631  << abort(FatalError);
632 
633  return fluxFieldType::null();
634  }
635 }
636 
637 
638 template<class Type>
641 (
642  const volScalarField& rho,
644  const fluxFieldType& phi
645 )
646 {
647  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
648 
649  if
650  (
651  U.dimensions() == dimVelocity
652  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
653  )
654  {
656  (
657  rho.oldTime()*U.oldTime()
658  );
659 
660  fluxFieldType phiCorr
661  (
662  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), rhoU0)
663  );
664 
665  return tmp<fluxFieldType>
666  (
667  new fluxFieldType
668  (
669  IOobject
670  (
671  "ddtCorr("
672  + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
673  mesh().time().timeName(),
674  mesh()
675  ),
676  this->fvcDdtPhiCoeff
677  (
678  rhoU0,
679  phi.oldTime(),
680  phiCorr,
681  rho.oldTime()
682  )*rDeltaT*phiCorr
683  )
684  );
685  }
686  else if
687  (
688  U.dimensions() == rho.dimensions()*dimVelocity
689  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
690  )
691  {
692  fluxFieldType phiCorr
693  (
694  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
695  );
696 
697  return tmp<fluxFieldType>
698  (
699  new fluxFieldType
700  (
701  IOobject
702  (
703  "ddtCorr("
704  + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
705  mesh().time().timeName(),
706  mesh()
707  ),
708  this->fvcDdtPhiCoeff
709  (
710  U.oldTime(),
711  phi.oldTime(),
712  phiCorr,
713  rho.oldTime()
714  )*rDeltaT*phiCorr
715  )
716  );
717  }
718  else
719  {
721  << "dimensions of phi are not correct"
722  << abort(FatalError);
723 
724  return fluxFieldType::null();
725  }
726 }
727 
728 
729 template<class Type>
731 (
733 )
734 {
735  return mesh().phi();
736 }
737 
738 
739 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
740 
741 } // End namespace fv
742 
743 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
744 
745 } // End namespace Foam
746 
747 // ************************************************************************* //
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::fv::EulerDdtScheme::fvcDdtPhiCorr
tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
Definition: EulerDdtScheme.C:524
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::fv::EulerDdtScheme::meshPhi
tmp< surfaceScalarField > meshPhi(const GeometricField< Type, fvPatchField, volMesh > &)
Definition: EulerDdtScheme.C:731
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
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::fv::EulerDdtScheme::fvcDdtUfCorr
tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Definition: EulerDdtScheme.C:491
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
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
EulerDdtScheme.H
phi
surfaceScalarField & phi
Definition: setRegionFluidFields.H:8
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::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
Foam::fv::EulerDdtScheme::fvmDdt
tmp< fvMatrix< Type > > fvmDdt(const GeometricField< Type, fvPatchField, volMesh > &)
Definition: EulerDdtScheme.C:337
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::dimVol
const dimensionSet dimVol(dimVolume)
Older spelling for dimVolume.
Definition: dimensionSets.H:64
Foam::dimensioned::dimensions
const dimensionSet & dimensions() const
Return const reference to dimensions.
Definition: dimensionedType.C:420
Foam::GeometricField< Type, fvPatchField, volMesh >
Foam::fv::EulerDdtScheme::fvcDdt
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)
Definition: EulerDdtScheme.C:48
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62