localEulerDdtScheme.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-2019 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 "localEulerDdtScheme.H"
27 #include "surfaceInterpolate.H"
28 #include "fvMatrices.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace fv
38 {
39 
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 
42 template<class Type>
43 const volScalarField& localEulerDdtScheme<Type>::localRDeltaT() const
44 {
46 }
47 
48 
49 template<class Type>
50 const surfaceScalarField& localEulerDdtScheme<Type>::localRDeltaTf() const
51 {
53 }
54 
55 
56 template<class Type>
57 tmp<GeometricField<Type, fvPatchField, volMesh>>
59 (
60  const dimensioned<Type>& dt
61 )
62 {
63  IOobject ddtIOobject
64  (
65  "ddt(" + dt.name() + ')',
66  mesh().time().timeName(),
67  mesh()
68  );
69 
71  (
73  (
74  ddtIOobject,
75  mesh(),
77  (
78  "0",
79  dt.dimensions()/dimTime,
80  Zero
81  ),
83  )
84  );
85 }
86 
87 
88 template<class Type>
91 (
93 )
94 {
95  const volScalarField& rDeltaT = localRDeltaT();
96 
97  IOobject ddtIOobject
98  (
99  "ddt(" + vf.name() + ')',
100  mesh().time().timeName(),
101  mesh()
102  );
103 
105  (
107  (
108  ddtIOobject,
109  rDeltaT*(vf - vf.oldTime())
110  )
111  );
112 }
113 
114 
115 template<class Type>
118 (
119  const dimensionedScalar& rho,
121 )
122 {
123  const volScalarField& rDeltaT = localRDeltaT();
124 
125  IOobject ddtIOobject
126  (
127  "ddt(" + rho.name() + ',' + vf.name() + ')',
128  mesh().time().timeName(),
129  mesh()
130  );
131 
133  (
135  (
136  ddtIOobject,
137  rDeltaT*rho*(vf - vf.oldTime())
138  )
139  );
140 }
141 
142 
143 template<class Type>
146 (
147  const volScalarField& rho,
149 )
150 {
151  const volScalarField& rDeltaT = localRDeltaT();
152 
153  IOobject ddtIOobject
154  (
155  "ddt(" + rho.name() + ',' + vf.name() + ')',
156  mesh().time().timeName(),
157  mesh()
158  );
159 
161  (
163  (
164  ddtIOobject,
165  rDeltaT*(rho*vf - rho.oldTime()*vf.oldTime())
166  )
167  );
168 }
169 
170 
171 template<class Type>
174 (
175  const volScalarField& alpha,
176  const volScalarField& rho,
178 )
179 {
180  const volScalarField& rDeltaT = localRDeltaT();
181 
182  IOobject ddtIOobject
183  (
184  "ddt("+alpha.name()+','+rho.name()+','+vf.name()+')',
185  mesh().time().timeName(),
186  mesh()
187  );
188 
190  (
192  (
193  ddtIOobject,
194  rDeltaT
195  *(
196  alpha*rho*vf
197  - alpha.oldTime()*rho.oldTime()*vf.oldTime()
198  )
199  )
200  );
201 }
202 
203 
204 template<class Type>
207 (
209 )
210 {
211  const surfaceScalarField& rDeltaT = localRDeltaTf();
212 
213  IOobject ddtIOobject
214  (
215  "ddt("+sf.name()+')',
216  mesh().time().timeName(),
217  mesh()
218  );
219 
221  (
223  (
224  ddtIOobject,
225  rDeltaT*(sf - sf.oldTime())
226  )
227  );
228 }
229 
230 
231 template<class Type>
234 (
236 )
237 {
238  tmp<fvMatrix<Type>> tfvm
239  (
240  new fvMatrix<Type>
241  (
242  vf,
244  )
245  );
246 
247  fvMatrix<Type>& fvm = tfvm.ref();
248 
249  const scalarField& rDeltaT = localRDeltaT();
250 
251  fvm.diag() = rDeltaT*mesh().Vsc();
252  fvm.source() = rDeltaT*vf.oldTime().primitiveField()*mesh().Vsc();
253 
254  return tfvm;
255 }
256 
257 
258 template<class Type>
261 (
262  const dimensionedScalar& rho,
264 )
265 {
266  tmp<fvMatrix<Type>> tfvm
267  (
268  new fvMatrix<Type>
269  (
270  vf,
272  )
273  );
274  fvMatrix<Type>& fvm = tfvm.ref();
275 
276  const scalarField& rDeltaT = localRDeltaT();
277 
278  fvm.diag() = rDeltaT*rho.value()*mesh().Vsc();
279 
280  fvm.source() =
281  rDeltaT*rho.value()*vf.oldTime().primitiveField()*mesh().Vsc();
282 
283  return tfvm;
284 }
285 
286 
287 template<class Type>
290 (
291  const volScalarField& rho,
293 )
294 {
295  tmp<fvMatrix<Type>> tfvm
296  (
297  new fvMatrix<Type>
298  (
299  vf,
301  )
302  );
303  fvMatrix<Type>& fvm = tfvm.ref();
304 
305  const scalarField& rDeltaT = localRDeltaT();
306 
307  fvm.diag() = rDeltaT*rho.primitiveField()*mesh().Vsc();
308 
309  fvm.source() = rDeltaT
310  *rho.oldTime().primitiveField()
311  *vf.oldTime().primitiveField()*mesh().Vsc();
312 
313  return tfvm;
314 }
315 
316 
317 template<class Type>
320 (
321  const volScalarField& alpha,
322  const volScalarField& rho,
324 )
325 {
326  tmp<fvMatrix<Type>> tfvm
327  (
328  new fvMatrix<Type>
329  (
330  vf,
331  alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
332  )
333  );
334  fvMatrix<Type>& fvm = tfvm.ref();
335 
336  const scalarField& rDeltaT = localRDeltaT();
337 
338  fvm.diag() =
339  rDeltaT*alpha.primitiveField()*rho.primitiveField()*mesh().Vsc();
340 
341  fvm.source() = rDeltaT
342  *alpha.oldTime().primitiveField()
343  *rho.oldTime().primitiveField()
344  *vf.oldTime().primitiveField()*mesh().Vsc();
345 
346  return tfvm;
347 }
348 
349 
350 /*
351 // Courant number limited formulation
352 template<class Type>
353 tmp<surfaceScalarField> localEulerDdtScheme<Type>::fvcDdtPhiCoeff
354 (
355  const GeometricField<Type, fvPatchField, volMesh>& U,
356  const fluxFieldType& phi,
357  const fluxFieldType& phiCorr
358 )
359 {
360  // Courant number limited formulation
361  tmp<surfaceScalarField> tddtCouplingCoeff = scalar(1)
362  - min
363  (
364  mag(phiCorr)*mesh().deltaCoeffs()
365  /(fvc::interpolate(localRDeltaT())*mesh().magSf()),
366  scalar(1)
367  );
368 
369  surfaceScalarField& ddtCouplingCoeff = tddtCouplingCoeff.ref();
370 
371  surfaceScalarField::Boundary& ccbf = ddtCouplingCoeff.boundaryFieldRef();
372 
373  forAll(U.boundaryField(), patchi)
374  {
375  if
376  (
377  U.boundaryField()[patchi].fixesValue()
378  || isA<cyclicAMIFvPatch>(mesh().boundary()[patchi])
379  )
380  {
381  ccbf[patchi] = 0.0;
382  }
383  }
384 
385  if (debug > 1)
386  {
387  InfoInFunction
388  << "ddtCouplingCoeff mean max min = "
389  << gAverage(ddtCouplingCoeff.primitiveField())
390  << " " << gMax(ddtCouplingCoeff.primitiveField())
391  << " " << gMin(ddtCouplingCoeff.primitiveField())
392  << endl;
393  }
394 
395  return tddtCouplingCoeff;
396 }
397 */
398 
399 
400 template<class Type>
403 (
406 )
407 {
408  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
409 
410  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
411  fluxFieldType phiCorr
412  (
413  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
414  );
415 
416  return tmp<fluxFieldType>
417  (
418  new fluxFieldType
419  (
420  IOobject
421  (
422  "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
423  mesh().time().timeName(),
424  mesh()
425  ),
426  this->fvcDdtPhiCoeff(U.oldTime(), phiUf0, phiCorr)
427  *rDeltaT*phiCorr
428  )
429  );
430 }
431 
432 
433 template<class Type>
436 (
438  const fluxFieldType& phi
439 )
440 {
441  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
442 
443  fluxFieldType phiCorr
444  (
445  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
446  );
447 
448  return tmp<fluxFieldType>
449  (
450  new fluxFieldType
451  (
452  IOobject
453  (
454  "ddtCorr(" + U.name() + ',' + phi.name() + ')',
455  mesh().time().timeName(),
456  mesh()
457  ),
458  this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime(), phiCorr)
459  *rDeltaT*phiCorr
460  )
461  );
462 }
463 
464 
465 template<class Type>
468 (
469  const volScalarField& rho,
472 )
473 {
474  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
475 
476  if
477  (
478  U.dimensions() == dimVelocity
480  )
481  {
483  (
484  rho.oldTime()*U.oldTime()
485  );
486 
487  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
488  fluxFieldType phiCorr(phiUf0 - fvc::dotInterpolate(mesh().Sf(), rhoU0));
489 
490  return tmp<fluxFieldType>
491  (
492  new fluxFieldType
493  (
494  IOobject
495  (
496  "ddtCorr("
497  + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
498  mesh().time().timeName(),
499  mesh()
500  ),
501  this->fvcDdtPhiCoeff(rhoU0, phiUf0, phiCorr, rho.oldTime())
502  *rDeltaT*phiCorr
503  )
504  );
505  }
506  else if
507  (
510  )
511  {
512  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
513  fluxFieldType phiCorr
514  (
515  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
516  );
517 
518  return tmp<fluxFieldType>
519  (
520  new fluxFieldType
521  (
522  IOobject
523  (
524  "ddtCorr("
525  + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
526  mesh().time().timeName(),
527  mesh()
528  ),
529  this->fvcDdtPhiCoeff
530  (
531  U.oldTime(),
532  phiUf0,
533  phiCorr,
534  rho.oldTime()
535  )*rDeltaT*phiCorr
536  )
537  );
538  }
539  else
540  {
542  << "dimensions of Uf are not correct"
543  << abort(FatalError);
544 
545  return fluxFieldType::null();
546  }
547 }
548 
549 
550 template<class Type>
553 (
554  const volScalarField& rho,
556  const fluxFieldType& phi
557 )
558 {
559  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
560 
561  if
562  (
563  U.dimensions() == dimVelocity
564  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
565  )
566  {
568  (
569  rho.oldTime()*U.oldTime()
570  );
571 
572  fluxFieldType phiCorr
573  (
574  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), rhoU0)
575  );
576 
577  return tmp<fluxFieldType>
578  (
579  new fluxFieldType
580  (
581  IOobject
582  (
583  "ddtCorr("
584  + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
585  mesh().time().timeName(),
586  mesh()
587  ),
588  this->fvcDdtPhiCoeff
589  (
590  rhoU0,
591  phi.oldTime(),
592  phiCorr,
593  rho.oldTime()
594  )*rDeltaT*phiCorr
595  )
596  );
597  }
598  else if
599  (
600  U.dimensions() == rho.dimensions()*dimVelocity
601  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
602  )
603  {
604  fluxFieldType phiCorr
605  (
606  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
607  );
608 
609  return tmp<fluxFieldType>
610  (
611  new fluxFieldType
612  (
613  IOobject
614  (
615  "ddtCorr("
616  + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
617  mesh().time().timeName(),
618  mesh()
619  ),
620  this->fvcDdtPhiCoeff
621  (
622  U.oldTime(),
623  phi.oldTime(),
624  phiCorr,
625  rho.oldTime()
626  )*rDeltaT*phiCorr
627  )
628  );
629  }
630  else
631  {
633  << "dimensions of phi are not correct"
634  << abort(FatalError);
635 
636  return fluxFieldType::null();
637  }
638 }
639 
640 
641 template<class Type>
643 (
645 )
646 {
648  (
649  "meshPhi",
650  mesh(),
652  );
653 }
654 
655 
656 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
657 
658 } // End namespace fv
659 
660 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
661 
662 } // End namespace Foam
663 
664 // ************************************************************************* //
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
const word & name() const
Return name.
Definition: IOobject.H:295
static const surfaceScalarField & localRDeltaTf(const fvMesh &mesh)
Return the reciprocal of the local face time-step.
Definition: localEulerDdt.C:58
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
static tmp< GeometricField< typename innerProduct< vector, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const surfaceVectorField &Sf, const GeometricField< Type, fvPatchField, volMesh > &tvf)
Interpolate field onto faces.
virtual tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
virtual tmp< fvMatrix< Type > > fvmDdt(const GeometricField< Type, fvPatchField, volMesh > &)
const dimensionSet dimVol(dimVolume)
Definition: dimensionSets.H:59
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
Generic dimensioned Type class.
static tmp< GeometricField< scalar, fvsPatchField, surfaceMesh > > New(const word &name, const Mesh &, const dimensionSet &, const word &patchFieldType=fvsPatchField< scalar >::calculatedType())
Return a temporary field constructed from name, mesh, dimensionSet.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
const dimensionSet & dimensions() const
Return dimensions.
dynamicFvMesh & mesh
labelList fv(nPoints)
virtual tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
This boundary condition is not designed to be evaluated; it is assmued 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
static const volScalarField & localRDeltaT(const fvMesh &mesh)
Return the reciprocal of the local time-step.
Definition: localEulerDdt.C:46
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Field< Type > & source()
Definition: fvMatrix.H:292
virtual tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)
const word & name() const
Return const reference to name.
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.
const dimensionSet dimDensity
virtual tmp< surfaceScalarField > meshPhi(const GeometricField< Type, fvPatchField, volMesh > &)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
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.
mesh Sf()
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
scalarField & diag()
Definition: lduMatrix.C:186
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
Namespace for OpenFOAM.
const dimensionSet dimVelocity