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-2022 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<GeometricField<Type, fvPatchField, volMesh>>
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  {
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  {
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 (
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  {
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  {
124  (
125  ddtName,
126  rDeltaT*(vf - vf.oldTime())
127  );
128  }
129 }
130 
131 
132 template<class Type>
135 (
136  const dimensionedScalar& rho,
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  {
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  {
163  (
164  ddtName,
165  rDeltaT*rho*(vf - vf.oldTime())
166  );
167  }
168 }
169 
170 
171 template<class Type>
174 (
175  const volScalarField& rho,
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  {
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  {
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,
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  {
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  {
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 (
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 (
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,
328 )
329 {
330  tmp<fvMatrix<Type>> tfvm
331  (
332  new fvMatrix<Type>
333  (
334  vf,
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,
365 )
366 {
367  tmp<fvMatrix<Type>> tfvm
368  (
369  new fvMatrix<Type>
370  (
371  vf,
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,
405 )
406 {
407  tmp<fvMatrix<Type>> tfvm
408  (
409  new fvMatrix<Type>
410  (
411  vf,
412  alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
413  )
414  );
415  fvMatrix<Type>& fvm = tfvm.ref();
416 
417  const scalar rDeltaT = 1.0/mesh().time().deltaTValue();
418 
419  fvm.diag() =
420  rDeltaT*alpha.primitiveField()*rho.primitiveField()*mesh().Vsc();
421 
422  if (mesh().moving())
423  {
424  fvm.source() = rDeltaT
425  *alpha.oldTime().primitiveField()
426  *rho.oldTime().primitiveField()
427  *vf.oldTime().primitiveField()*mesh().Vsc0();
428  }
429  else
430  {
431  fvm.source() = rDeltaT
432  *alpha.oldTime().primitiveField()
433  *rho.oldTime().primitiveField()
434  *vf.oldTime().primitiveField()*mesh().Vsc();
435  }
436 
437  return tfvm;
438 }
439 
440 
441 template<class Type>
444 (
447 )
448 {
449  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
450 
451  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
452  fluxFieldType phiCorr
453  (
454  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
455  );
456 
457  return fluxFieldType::New
458  (
459  "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
460  this->fvcDdtPhiCoeff(U.oldTime(), phiUf0, phiCorr) *rDeltaT*phiCorr
461  );
462 }
463 
464 
465 template<class Type>
468 (
470  const fluxFieldType& phi
471 )
472 {
473  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
474 
475  fluxFieldType phiCorr
476  (
477  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
478  );
479 
480  return fluxFieldType::New
481  (
482  "ddtCorr(" + U.name() + ',' + phi.name() + ')',
483  this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime(), phiCorr)
484  *rDeltaT*phiCorr
485  );
486 }
487 
488 
489 template<class Type>
492 (
493  const volScalarField& rho,
496 )
497 {
498  const dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
499 
500  if
501  (
502  U.dimensions() == dimVelocity
503  && Uf.dimensions() == rho.dimensions()*dimVelocity
504  )
505  {
507  (
508  rho.oldTime()*U.oldTime()
509  );
510 
511  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
512  fluxFieldType phiCorr(phiUf0 - fvc::dotInterpolate(mesh().Sf(), rhoU0));
513 
514  return fluxFieldType::New
515  (
516  "ddtCorr(" + rho.name() + ',' + U.name() + ',' + Uf.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  && Uf.dimensions() == rho.dimensions()*dimVelocity
525  )
526  {
527  fluxFieldType phiUf0(mesh().Sf() & Uf.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() + ',' + Uf.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 Uf are not correct"
549  << abort(FatalError);
550 
551  return fluxFieldType::null();
552  }
553 }
554 
555 
556 template<class Type>
559 (
560  const volScalarField& rho,
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()*dimFlux
571  )
572  {
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()*dimFlux
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>
631 (
633 )
634 {
635  return mesh().phi();
636 }
637 
638 
639 template<class Type>
641 (
643  const label patchi
644 )
645 {
646  return mesh().phi().boundaryField()[patchi];
647 }
648 
649 
650 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
651 
652 } // End namespace fv
653 
654 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
655 
656 } // End namespace Foam
657 
658 // ************************************************************************* //
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
const word & name() const
Return name.
Definition: IOobject.H:315
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
static tmp< GeometricField< typename innerProduct< vector, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const surfaceVectorField &Sf, const GeometricField< Type, fvPatchField, volMesh > &tvf)
Interpolate field onto faces.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &)
Return a temporary field constructed from name,.
virtual tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
Generic dimensioned Type class.
fvMesh & mesh
const dimensionSet dimTime
const dimensionSet & dimensions() const
Return dimensions.
virtual tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
virtual tmp< surfaceScalarField > meshPhi(const GeometricField< Type, fvPatchField, volMesh > &)
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
const dimensionSet dimVol
const dimensionSet dimFlux
This boundary condition is not designed to be evaluated; it is assumed that the value is assigned via...
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
static const zero Zero
Definition: zero.H:97
errorManip< error > abort(error &err)
Definition: errorManip.H:131
virtual tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Calculate the divergence of the given field.
const dimensionSet dimVelocity
Field< Type > & source()
Definition: fvMatrix.H:300
const word & name() const
Return const reference to name.
label patchi
A special matrix type and solver, designed for finite volume solutions of scalar equations.
const dimensionSet & dimensions() const
Return const reference to dimensions.
A class for managing temporary objects.
Definition: PtrList.H:53
scalarField & diag()
Definition: lduMatrix.C:186
virtual tmp< fvMatrix< Type > > fvmDdt(const GeometricField< Type, fvPatchField, volMesh > &)
Namespace for OpenFOAM.