EulerDdtScheme.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-2024 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 "EulerDdtScheme.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 tmp<VolField<Type>>
46 (
47  const dimensioned<Type>& dt
48 )
49 {
50  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
51 
52  const word ddtName("ddt("+dt.name()+')');
53 
54  if (mesh().moving())
55  {
56  tmp<VolField<Type>> tdtdt
57  (
59  (
60  ddtName,
61  mesh(),
63  (
64  "0",
65  dt.dimensions()/dimTime,
66  Zero
67  )
68  )
69  );
70 
71  tdtdt.ref().primitiveFieldRef() =
72  rDeltaT.value()*dt.value()*(1.0 - mesh().Vsc0()/mesh().Vsc());
73 
74  return tdtdt;
75  }
76  else
77  {
78  return VolField<Type>::New
79  (
80  ddtName,
81  mesh(),
83  (
84  "0",
85  dt.dimensions()/dimTime,
86  Zero
87  ),
89  );
90  }
91 }
92 
93 
94 template<class Type>
97 (
98  const VolField<Type>& vf
99 )
100 {
101  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
102 
103  const word ddtName("ddt("+vf.name()+')');
104 
105  if (mesh().moving())
106  {
107  return VolField<Type>::New
108  (
109  ddtName,
110  rDeltaT*
111  (
112  vf()
113  - vf.oldTime()()*mesh().Vsc0()/mesh().Vsc()
114  ),
115  rDeltaT.value()*
116  (
117  vf.boundaryField() - vf.oldTime().boundaryField()
118  )
119  );
120  }
121  else
122  {
123  return VolField<Type>::New
124  (
125  ddtName,
126  rDeltaT*(vf - vf.oldTime())
127  );
128  }
129 }
130 
131 
132 template<class Type>
135 (
136  const dimensionedScalar& rho,
137  const VolField<Type>& vf
138 )
139 {
140  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
141 
142  const word ddtName("ddt("+rho.name()+','+vf.name()+')');
143 
144  if (mesh().moving())
145  {
146  return VolField<Type>::New
147  (
148  ddtName,
149  rDeltaT*rho*
150  (
151  vf()
152  - vf.oldTime()()*mesh().Vsc0()/mesh().Vsc()
153  ),
154  rDeltaT.value()*rho.value()*
155  (
156  vf.boundaryField() - vf.oldTime().boundaryField()
157  )
158  );
159  }
160  else
161  {
162  return VolField<Type>::New
163  (
164  ddtName,
165  rDeltaT*rho*(vf - vf.oldTime())
166  );
167  }
168 }
169 
170 
171 template<class Type>
174 (
175  const volScalarField& rho,
176  const VolField<Type>& vf
177 )
178 {
179  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
180 
181  const word ddtName("ddt("+rho.name()+','+vf.name()+')');
182 
183  if (mesh().moving())
184  {
185  return VolField<Type>::New
186  (
187  ddtName,
188  rDeltaT*
189  (
190  rho()*vf()
191  - rho.oldTime()()
192  *vf.oldTime()()*mesh().Vsc0()/mesh().Vsc()
193  ),
194  rDeltaT.value()*
195  (
196  rho.boundaryField()*vf.boundaryField()
197  - rho.oldTime().boundaryField()
198  *vf.oldTime().boundaryField()
199  )
200  );
201  }
202  else
203  {
204  return VolField<Type>::New
205  (
206  ddtName,
207  rDeltaT*(rho*vf - rho.oldTime()*vf.oldTime())
208  );
209  }
210 }
211 
212 
213 template<class Type>
216 (
217  const volScalarField& alpha,
218  const volScalarField& rho,
219  const VolField<Type>& vf
220 )
221 {
222  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
223 
224  const word ddtName("ddt("+alpha.name()+','+rho.name()+','+vf.name()+')');
225 
226  if (mesh().moving())
227  {
228  return VolField<Type>::New
229  (
230  ddtName,
231  rDeltaT*
232  (
233  alpha()
234  *rho()
235  *vf()
236 
237  - alpha.oldTime()()
238  *rho.oldTime()()
239  *vf.oldTime()()*mesh().Vsc0()/mesh().Vsc()
240  ),
241  rDeltaT.value()*
242  (
243  alpha.boundaryField()
244  *rho.boundaryField()
245  *vf.boundaryField()
246 
247  - alpha.oldTime().boundaryField()
248  *rho.oldTime().boundaryField()
249  *vf.oldTime().boundaryField()
250  )
251  );
252  }
253  else
254  {
255  return VolField<Type>::New
256  (
257  ddtName,
258  rDeltaT
259  *(
260  alpha*rho*vf
261  - alpha.oldTime()*rho.oldTime()*vf.oldTime()
262  )
263  );
264  }
265 }
266 
267 
268 template<class Type>
271 (
272  const SurfaceField<Type>& sf
273 )
274 {
275  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
276 
277  const word ddtName("ddt("+sf.name()+')');
278 
280  (
281  ddtName,
282  rDeltaT*(sf - sf.oldTime())
283  );
284 }
285 
286 
287 template<class Type>
290 (
291  const VolField<Type>& vf
292 )
293 {
294  tmp<fvMatrix<Type>> tfvm
295  (
296  new fvMatrix<Type>
297  (
298  vf,
300  )
301  );
302 
303  fvMatrix<Type>& fvm = tfvm.ref();
304 
305  const scalar rDeltaT = 1.0/mesh().time().deltaTValue();
306 
307  fvm.diag() = rDeltaT*mesh().Vsc();
308 
309  if (mesh().moving())
310  {
311  fvm.source() = rDeltaT*vf.oldTime().primitiveField()*mesh().Vsc0();
312  }
313  else
314  {
315  fvm.source() = rDeltaT*vf.oldTime().primitiveField()*mesh().Vsc();
316  }
317 
318  return tfvm;
319 }
320 
321 
322 template<class Type>
325 (
326  const dimensionedScalar& rho,
327  const VolField<Type>& vf
328 )
329 {
330  tmp<fvMatrix<Type>> tfvm
331  (
332  new fvMatrix<Type>
333  (
334  vf,
335  rho.dimensions()*vf.dimensions()*dimVolume/dimTime
336  )
337  );
338  fvMatrix<Type>& fvm = tfvm.ref();
339 
340  const scalar rDeltaT = 1.0/mesh().time().deltaTValue();
341 
342  fvm.diag() = rDeltaT*rho.value()*mesh().Vsc();
343 
344  if (mesh().moving())
345  {
346  fvm.source() = rDeltaT
347  *rho.value()*vf.oldTime().primitiveField()*mesh().Vsc0();
348  }
349  else
350  {
351  fvm.source() = rDeltaT
352  *rho.value()*vf.oldTime().primitiveField()*mesh().Vsc();
353  }
354 
355  return tfvm;
356 }
357 
358 
359 template<class Type>
362 (
363  const volScalarField& rho,
364  const VolField<Type>& vf
365 )
366 {
367  tmp<fvMatrix<Type>> tfvm
368  (
369  new fvMatrix<Type>
370  (
371  vf,
372  rho.dimensions()*vf.dimensions()*dimVolume/dimTime
373  )
374  );
375  fvMatrix<Type>& fvm = tfvm.ref();
376 
377  const scalar rDeltaT = 1.0/mesh().time().deltaTValue();
378 
379  fvm.diag() = rDeltaT*rho.primitiveField()*mesh().Vsc();
380 
381  if (mesh().moving())
382  {
383  fvm.source() = rDeltaT
384  *rho.oldTime().primitiveField()
385  *vf.oldTime().primitiveField()*mesh().Vsc0();
386  }
387  else
388  {
389  fvm.source() = rDeltaT
390  *rho.oldTime().primitiveField()
391  *vf.oldTime().primitiveField()*mesh().Vsc();
392  }
393 
394  return tfvm;
395 }
396 
397 
398 template<class Type>
401 (
402  const volScalarField& alpha,
403  const volScalarField& rho,
404  const VolField<Type>& vf
405 )
406 {
407  tmp<fvMatrix<Type>> tfvm
408  (
409  new fvMatrix<Type>
410  (
411  vf,
412  alpha.dimensions()
413  *rho.dimensions()
414  *vf.dimensions()
415  *dimVolume
416  /dimTime
417  )
418  );
419  fvMatrix<Type>& fvm = tfvm.ref();
420 
421  const scalar rDeltaT = 1.0/mesh().time().deltaTValue();
422 
423  fvm.diag() =
424  rDeltaT*alpha.primitiveField()*rho.primitiveField()*mesh().Vsc();
425 
426  if (mesh().moving())
427  {
428  fvm.source() = rDeltaT
429  *alpha.oldTime().primitiveField()
430  *rho.oldTime().primitiveField()
431  *vf.oldTime().primitiveField()*mesh().Vsc0();
432  }
433  else
434  {
435  fvm.source() = rDeltaT
436  *alpha.oldTime().primitiveField()
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 VolField<Type>& U,
450  const SurfaceField<Type>& Uf
451 )
452 {
453  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
454 
455  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
456  fluxFieldType phiCorr
457  (
458  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
459  );
460 
461  return fluxFieldType::New
462  (
463  "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
464  this->fvcDdtPhiCoeff(U.oldTime(), phiUf0, phiCorr)*rDeltaT*phiCorr
465  );
466 }
467 
468 
469 template<class Type>
472 (
473  const VolField<Type>& U,
474  const fluxFieldType& phi
475 )
476 {
477  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
478 
479  fluxFieldType phiCorr
480  (
481  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
482  );
483 
484  return fluxFieldType::New
485  (
486  "ddtCorr(" + U.name() + ',' + phi.name() + ')',
487  this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime(), phiCorr)
488  *rDeltaT*phiCorr
489  );
490 }
491 
492 
493 template<class Type>
496 (
497  const volScalarField& rho,
498  const VolField<Type>& U,
499  const SurfaceField<Type>& rhoUf
500 )
501 {
502  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
503 
504  if
505  (
506  U.dimensions() == dimVelocity
507  && rhoUf.dimensions() == rho.dimensions()*dimVelocity
508  )
509  {
510  VolField<Type> rhoU0(rho.oldTime()*U.oldTime());
511  fluxFieldType phiUf0(mesh().Sf() & rhoUf.oldTime());
512  fluxFieldType phiCorr(phiUf0 - fvc::dotInterpolate(mesh().Sf(), rhoU0));
513 
514  return fluxFieldType::New
515  (
516  "ddtCorr(" + rho.name() + ',' + U.name() + ',' + rhoUf.name() + ')',
517  this->fvcDdtPhiCoeff(rhoU0, phiUf0, phiCorr, rho.oldTime())
518  *rDeltaT*phiCorr
519  );
520  }
521  else if
522  (
523  U.dimensions() == rho.dimensions()*dimVelocity
524  && rhoUf.dimensions() == rho.dimensions()*dimVelocity
525  )
526  {
527  fluxFieldType phiUf0(mesh().Sf() & rhoUf.oldTime());
528  fluxFieldType phiCorr
529  (
530  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
531  );
532 
533  return fluxFieldType::New
534  (
535  "ddtCorr(" + rho.name() + ',' + U.name() + ',' + rhoUf.name() + ')',
536  this->fvcDdtPhiCoeff
537  (
538  U.oldTime(),
539  phiUf0,
540  phiCorr,
541  rho.oldTime()
542  )*rDeltaT*phiCorr
543  );
544  }
545  else
546  {
548  << "dimensions of rhoUf are not correct"
549  << abort(FatalError);
550 
551  return fluxFieldType::null();
552  }
553 }
554 
555 
556 template<class Type>
559 (
560  const volScalarField& rho,
561  const VolField<Type>& U,
562  const fluxFieldType& phi
563 )
564 {
565  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
566 
567  if
568  (
569  U.dimensions() == dimVelocity
570  && phi.dimensions() == rho.dimensions()*dimVolumetricFlux
571  )
572  {
573  VolField<Type> rhoU0
574  (
575  rho.oldTime()*U.oldTime()
576  );
577 
578  fluxFieldType phiCorr
579  (
580  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), rhoU0)
581  );
582 
583  return fluxFieldType::New
584  (
585  "ddtCorr(" + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
586  this->fvcDdtPhiCoeff
587  (
588  rhoU0,
589  phi.oldTime(),
590  phiCorr,
591  rho.oldTime()
592  )*rDeltaT*phiCorr
593  );
594  }
595  else if
596  (
597  U.dimensions() == rho.dimensions()*dimVelocity
598  && phi.dimensions() == rho.dimensions()*dimVolumetricFlux
599  )
600  {
601  fluxFieldType phiCorr
602  (
603  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
604  );
605 
606  return fluxFieldType::New
607  (
608  "ddtCorr(" + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
609  this->fvcDdtPhiCoeff
610  (
611  U.oldTime(),
612  phi.oldTime(),
613  phiCorr,
614  rho.oldTime()
615  )*rDeltaT*phiCorr
616  );
617  }
618  else
619  {
621  << "dimensions of phi are not correct"
622  << abort(FatalError);
623 
624  return fluxFieldType::null();
625  }
626 }
627 
628 
629 template<class Type>
632 (
633  const volScalarField& alpha,
634  const volScalarField& rho,
635  const VolField<Type>& U,
636  const SurfaceField<Type>& Uf
637 )
638 {
639  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
640 
641  if (U.dimensions() == dimVelocity && Uf.dimensions() == dimVelocity)
642  {
643  const volScalarField alphaRho0(alpha.oldTime()*rho.oldTime());
644 
645  return fluxFieldType::New
646  (
647  "ddtCorr("
648  + alpha.name() + rho.name() + ',' + U.name() + ',' + Uf.name()
649  + ')',
650  this->fvcDdtPhiCoeff(U.oldTime(), mesh().Sf() & Uf.oldTime())
651  *rDeltaT
652  *(
653  mesh().Sf()
654  & (
655  fvc::interpolate(alphaRho0)*Uf.oldTime()
656  - fvc::interpolate(alphaRho0*U.oldTime())
657  )
658  )
659  );
660  }
661  else
662  {
664  << "dimensions of Uf are not correct"
665  << abort(FatalError);
666 
667  return fluxFieldType::null();
668  }
669 }
670 
671 
672 template<class Type>
675 (
676  const volScalarField& alpha,
677  const volScalarField& rho,
678  const VolField<Type>& U,
679  const fluxFieldType& phi
680 )
681 {
682  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
683 
684  if (U.dimensions() == dimVelocity && phi.dimensions() == dimVolumetricFlux)
685  {
686  const volScalarField alphaRho0(alpha.oldTime()*rho.oldTime());
687 
688  return fluxFieldType::New
689  (
690  "ddtCorr("
691  + alpha.name() + rho.name() + ',' + U.name() + ',' + phi.name()
692  + ')',
693  this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())
694  *rDeltaT
695  *(
696  fvc::interpolate(alphaRho0)*phi.oldTime()
697  -fvc::dotInterpolate(mesh().Sf(), alphaRho0*U.oldTime())
698  )
699  );
700  }
701  else
702  {
704  << "dimensions of phi are not correct"
705  << abort(FatalError);
706 
707  return fluxFieldType::null();
708  }
709 }
710 
711 
712 template<class Type>
714 (
715  const VolField<Type>&
716 )
717 {
718  return mesh().phi();
719 }
720 
721 
722 template<class Type>
724 (
725  const VolField<Type>&,
726  const label patchi
727 )
728 {
729  return mesh().phi().boundaryField()[patchi];
730 }
731 
732 
733 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
734 
735 } // End namespace fv
736 
737 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
738 
739 } // End namespace Foam
740 
741 // ************************************************************************* //
const dimensionSet & dimensions() const
Return dimensions.
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
const word & name() const
Return name.
Definition: IOobject.H:310
const FieldType & oldTime() const
Return the old time field.
Definition: OldTimeField.C:315
This boundary condition is not designed to be evaluated; it is assumed that the value is assigned via...
Generic dimensioned Type class.
const dimensionSet & dimensions() const
Return const reference to dimensions.
const Type & value() const
Return const reference to value.
const word & name() const
Return const reference to name.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
Field< Type > & source()
Definition: fvMatrix.H:307
virtual tmp< fluxFieldType > fvcDdtPhiCorr(const VolField< Type > &U, const fluxFieldType &phi)
virtual tmp< fvMatrix< Type > > fvmDdt(const VolField< Type > &)
virtual tmp< surfaceScalarField > meshPhi(const VolField< Type > &)
virtual tmp< fluxFieldType > fvcDdtUfCorr(const VolField< Type > &U, const SurfaceField< Type > &Uf)
virtual tmp< VolField< Type > > fvcDdt(const dimensioned< Type > &)
scalarField & diag()
Definition: lduMatrix.C:186
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
A special matrix type and solver, designed for finite volume solutions of scalar equations.
Calculate the divergence of the given field.
label patchi
volScalarField sf(fieldObject, mesh)
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
static tmp< SurfaceField< typename innerProduct< vector, Type >::type > > dotInterpolate(const surfaceVectorField &Sf, const VolField< Type > &tvf)
Interpolate field onto faces.
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
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 dimensionSet dimVolumetricFlux
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const dimensionSet dimTime
const dimensionSet dimVolume
const dimensionSet dimVelocity
error FatalError
labelList fv(nPoints)