localEulerDdtScheme.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 | Copyright (C) 2011-2016 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 tmp<GeometricField<Type, fvPatchField, volMesh>>
52 (
53  const dimensioned<Type>& dt
54 )
55 {
56  const volScalarField& rDeltaT = localRDeltaT();
57 
58  IOobject ddtIOobject
59  (
60  "ddt(" + dt.name() + ')',
61  mesh().time().timeName(),
62  mesh()
63  );
64 
65  if (mesh().moving())
66  {
68  (
70  (
71  ddtIOobject,
72  mesh(),
74  (
75  "0",
76  dt.dimensions()/dimTime,
77  Zero
78  )
79  )
80  );
81 
82  tdtdt.ref().primitiveFieldRef() =
83  rDeltaT.primitiveField()*dt.value()
84  *(1.0 - mesh().Vsc0()/mesh().Vsc());
85 
86  return tdtdt;
87  }
88  else
89  {
91  (
93  (
94  ddtIOobject,
95  mesh(),
97  (
98  "0",
99  dt.dimensions()/dimTime,
100  Zero
101  ),
103  )
104  );
105  }
106 }
107 
108 
109 template<class Type>
112 (
114 )
115 {
116  const volScalarField& rDeltaT = localRDeltaT();
117 
118  IOobject ddtIOobject
119  (
120  "ddt(" + vf.name() + ')',
121  mesh().time().timeName(),
122  mesh()
123  );
124 
125  if (mesh().moving())
126  {
128  (
130  (
131  ddtIOobject,
132  mesh(),
133  rDeltaT.dimensions()*vf.dimensions(),
134  rDeltaT.primitiveField()*
135  (
136  vf.primitiveField()
137  - vf.oldTime().primitiveField()*mesh().Vsc0()/mesh().Vsc()
138  ),
139  rDeltaT.boundaryField()*
140  (
141  vf.boundaryField() - vf.oldTime().boundaryField()
142  )
143  )
144  );
145  }
146  else
147  {
149  (
151  (
152  ddtIOobject,
153  rDeltaT*(vf - vf.oldTime())
154  )
155  );
156  }
157 }
158 
159 
160 template<class Type>
163 (
164  const dimensionedScalar& rho,
166 )
167 {
168  const volScalarField& rDeltaT = localRDeltaT();
169 
170  IOobject ddtIOobject
171  (
172  "ddt(" + rho.name() + ',' + vf.name() + ')',
173  mesh().time().timeName(),
174  mesh()
175  );
176 
177  if (mesh().moving())
178  {
180  (
182  (
183  ddtIOobject,
184  mesh(),
185  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
186  rDeltaT.primitiveField()*rho.value()*
187  (
188  vf.primitiveField()
189  - vf.oldTime().primitiveField()*mesh().Vsc0()/mesh().Vsc()
190  ),
191  rDeltaT.boundaryField()*rho.value()*
192  (
193  vf.boundaryField() - vf.oldTime().boundaryField()
194  )
195  )
196  );
197  }
198  else
199  {
201  (
203  (
204  ddtIOobject,
205  rDeltaT*rho*(vf - vf.oldTime())
206  )
207  );
208  }
209 }
210 
211 
212 template<class Type>
215 (
216  const volScalarField& rho,
218 )
219 {
220  const volScalarField& rDeltaT = localRDeltaT();
221 
222  IOobject ddtIOobject
223  (
224  "ddt(" + rho.name() + ',' + vf.name() + ')',
225  mesh().time().timeName(),
226  mesh()
227  );
228 
229  if (mesh().moving())
230  {
232  (
234  (
235  ddtIOobject,
236  mesh(),
237  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
238  rDeltaT.primitiveField()*
239  (
240  rho.primitiveField()*vf.primitiveField()
241  - rho.oldTime().primitiveField()
242  *vf.oldTime().primitiveField()*mesh().Vsc0()/mesh().Vsc()
243  ),
244  rDeltaT.boundaryField()*
245  (
246  rho.boundaryField()*vf.boundaryField()
247  - rho.oldTime().boundaryField()
248  *vf.oldTime().boundaryField()
249  )
250  )
251  );
252  }
253  else
254  {
256  (
258  (
259  ddtIOobject,
260  rDeltaT*(rho*vf - rho.oldTime()*vf.oldTime())
261  )
262  );
263  }
264 }
265 
266 
267 template<class Type>
270 (
271  const volScalarField& alpha,
272  const volScalarField& rho,
274 )
275 {
276  const volScalarField& rDeltaT = localRDeltaT();
277 
278  IOobject ddtIOobject
279  (
280  "ddt("+alpha.name()+','+rho.name()+','+vf.name()+')',
281  mesh().time().timeName(),
282  mesh()
283  );
284 
285  if (mesh().moving())
286  {
288  (
290  (
291  ddtIOobject,
292  mesh(),
293  rDeltaT.dimensions()
294  *alpha.dimensions()*rho.dimensions()*vf.dimensions(),
295  rDeltaT.primitiveField()*
296  (
297  alpha.primitiveField()
298  *rho.primitiveField()
299  *vf.primitiveField()
300 
301  - alpha.oldTime().primitiveField()
302  *rho.oldTime().primitiveField()
303  *vf.oldTime().primitiveField()*mesh().Vsc0()/mesh().Vsc()
304  ),
305  rDeltaT.boundaryField()*
306  (
307  alpha.boundaryField()
308  *rho.boundaryField()
309  *vf.boundaryField()
310 
311  - alpha.oldTime().boundaryField()
312  *rho.oldTime().boundaryField()
313  *vf.oldTime().boundaryField()
314  )
315  )
316  );
317  }
318  else
319  {
321  (
323  (
324  ddtIOobject,
325  rDeltaT
326  *(
327  alpha*rho*vf
328  - alpha.oldTime()*rho.oldTime()*vf.oldTime()
329  )
330  )
331  );
332  }
333 }
334 
335 
336 template<class Type>
339 (
341 )
342 {
343  tmp<fvMatrix<Type>> tfvm
344  (
345  new fvMatrix<Type>
346  (
347  vf,
349  )
350  );
351 
352  fvMatrix<Type>& fvm = tfvm.ref();
353 
354  const scalarField& rDeltaT = localRDeltaT();
355 
356  fvm.diag() = rDeltaT*mesh().Vsc();
357 
358  if (mesh().moving())
359  {
360  fvm.source() = rDeltaT*vf.oldTime().primitiveField()*mesh().Vsc0();
361  }
362  else
363  {
364  fvm.source() = rDeltaT*vf.oldTime().primitiveField()*mesh().Vsc();
365  }
366 
367  return tfvm;
368 }
369 
370 
371 template<class Type>
374 (
375  const dimensionedScalar& rho,
377 )
378 {
379  tmp<fvMatrix<Type>> tfvm
380  (
381  new fvMatrix<Type>
382  (
383  vf,
385  )
386  );
387  fvMatrix<Type>& fvm = tfvm.ref();
388 
389  const scalarField& rDeltaT = localRDeltaT();
390 
391  fvm.diag() = rDeltaT*rho.value()*mesh().Vsc();
392 
393  if (mesh().moving())
394  {
395  fvm.source() = rDeltaT
396  *rho.value()*vf.oldTime().primitiveField()*mesh().Vsc0();
397  }
398  else
399  {
400  fvm.source() = rDeltaT
401  *rho.value()*vf.oldTime().primitiveField()*mesh().Vsc();
402  }
403 
404  return tfvm;
405 }
406 
407 
408 template<class Type>
411 (
412  const volScalarField& rho,
414 )
415 {
416  tmp<fvMatrix<Type>> tfvm
417  (
418  new fvMatrix<Type>
419  (
420  vf,
422  )
423  );
424  fvMatrix<Type>& fvm = tfvm.ref();
425 
426  const scalarField& rDeltaT = localRDeltaT();
427 
428  fvm.diag() = rDeltaT*rho.primitiveField()*mesh().Vsc();
429 
430  if (mesh().moving())
431  {
432  fvm.source() = rDeltaT
433  *rho.oldTime().primitiveField()
434  *vf.oldTime().primitiveField()*mesh().Vsc0();
435  }
436  else
437  {
438  fvm.source() = rDeltaT
439  *rho.oldTime().primitiveField()
440  *vf.oldTime().primitiveField()*mesh().Vsc();
441  }
442 
443  return tfvm;
444 }
445 
446 
447 template<class Type>
450 (
451  const volScalarField& alpha,
452  const volScalarField& rho,
454 )
455 {
456  tmp<fvMatrix<Type>> tfvm
457  (
458  new fvMatrix<Type>
459  (
460  vf,
461  alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
462  )
463  );
464  fvMatrix<Type>& fvm = tfvm.ref();
465 
466  const scalarField& rDeltaT = localRDeltaT();
467 
468  fvm.diag() =
469  rDeltaT*alpha.primitiveField()*rho.primitiveField()*mesh().Vsc();
470 
471  if (mesh().moving())
472  {
473  fvm.source() = rDeltaT
474  *alpha.oldTime().primitiveField()
475  *rho.oldTime().primitiveField()
476  *vf.oldTime().primitiveField()*mesh().Vsc0();
477  }
478  else
479  {
480  fvm.source() = rDeltaT
481  *alpha.oldTime().primitiveField()
482  *rho.oldTime().primitiveField()
483  *vf.oldTime().primitiveField()*mesh().Vsc();
484  }
485 
486  return tfvm;
487 }
488 
489 
490 template<class Type>
493 (
496 )
497 {
498  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
499 
500  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
501  fluxFieldType phiCorr
502  (
503  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
504  );
505 
506  return tmp<fluxFieldType>
507  (
508  new fluxFieldType
509  (
510  IOobject
511  (
512  "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
513  mesh().time().timeName(),
514  mesh()
515  ),
516  this->fvcDdtPhiCoeff(U.oldTime(), phiUf0, phiCorr)
517  *rDeltaT*phiCorr
518  )
519  );
520 }
521 
522 
523 template<class Type>
526 (
528  const fluxFieldType& phi
529 )
530 {
531  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
532 
533  fluxFieldType phiCorr
534  (
535  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
536  );
537 
538  return tmp<fluxFieldType>
539  (
540  new fluxFieldType
541  (
542  IOobject
543  (
544  "ddtCorr(" + U.name() + ',' + phi.name() + ')',
545  mesh().time().timeName(),
546  mesh()
547  ),
548  this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime(), phiCorr)
549  *rDeltaT*phiCorr
550  )
551  );
552 }
553 
554 
555 template<class Type>
558 (
559  const volScalarField& rho,
562 )
563 {
564  if
565  (
566  U.dimensions() == dimVelocity
568  )
569  {
570  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
571 
573  (
574  rho.oldTime()*U.oldTime()
575  );
576 
577  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
578  fluxFieldType phiCorr(phiUf0 - fvc::dotInterpolate(mesh().Sf(), rhoU0));
579 
580  return tmp<fluxFieldType>
581  (
582  new fluxFieldType
583  (
584  IOobject
585  (
586  "ddtCorr("
587  + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
588  mesh().time().timeName(),
589  mesh()
590  ),
591  this->fvcDdtPhiCoeff(rhoU0, phiUf0, phiCorr)*rDeltaT*phiCorr
592  )
593  );
594  }
595  else if
596  (
599  )
600  {
601  return fvcDdtUfCorr(U, Uf);
602  }
603  else
604  {
606  << "dimensions of Uf are not correct"
607  << abort(FatalError);
608 
609  return fluxFieldType::null();
610  }
611 }
612 
613 
614 template<class Type>
617 (
618  const volScalarField& rho,
620  const fluxFieldType& phi
621 )
622 {
623  if
624  (
625  U.dimensions() == dimVelocity
626  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
627  )
628  {
629  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
630 
632  (
633  rho.oldTime()*U.oldTime()
634  );
635 
636  fluxFieldType phiCorr
637  (
638  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), rhoU0)
639  );
640 
641  return tmp<fluxFieldType>
642  (
643  new fluxFieldType
644  (
645  IOobject
646  (
647  "ddtCorr("
648  + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
649  mesh().time().timeName(),
650  mesh()
651  ),
652  this->fvcDdtPhiCoeff(rhoU0, phi.oldTime(), phiCorr)
653  *rDeltaT*phiCorr
654  )
655  );
656  }
657  else if
658  (
659  U.dimensions() == rho.dimensions()*dimVelocity
660  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
661  )
662  {
663  return fvcDdtPhiCorr(U, phi);
664  }
665  else
666  {
668  << "dimensions of phi are not correct"
669  << abort(FatalError);
670 
671  return fluxFieldType::null();
672  }
673 }
674 
675 
676 template<class Type>
678 (
680 )
681 {
683  (
685  (
686  IOobject
687  (
688  "meshPhi",
689  mesh().time().timeName(),
690  mesh(),
693  false
694  ),
695  mesh(),
697  )
698  );
699 }
700 
701 
702 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
703 
704 } // End namespace fv
705 
706 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
707 
708 } // End namespace Foam
709 
710 // ************************************************************************* //
const dimensionSet & dimensions() const
Return const reference to dimensions.
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static tmp< GeometricField< typename innerProduct< vector, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const surfaceVectorField &Sf, const GeometricField< Type, fvPatchField, volMesh > &tvf)
Interpolate field onto faces.
tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
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.
const Type & value() const
Return const reference to value.
Generic dimensioned Type class.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const word & name() const
Return const reference to name.
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
dynamicFvMesh & mesh
const Boundary & boundaryField() const
Return const-reference to the boundary field.
labelList fv(nPoints)
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...
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:71
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:91
static const volScalarField & localRDeltaT(const fvMesh &mesh)
Return the reciprocal of the local time-step.
Definition: localEulerDdt.C:45
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const dimensionSet & dimensions() const
Return dimensions.
Field< Type > & source()
Definition: fvMatrix.H:291
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)
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
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.
A class for managing temporary objects.
Definition: PtrList.H:54
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
scalarField & diag()
Definition: lduMatrix.C:183
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
const word & name() const
Return name.
Definition: IOobject.H:260
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
Namespace for OpenFOAM.
const dimensionSet dimVelocity