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-2018 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 template<class Type>
353 (
356 )
357 {
358  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
359 
360  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
361  fluxFieldType phiCorr
362  (
363  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
364  );
365 
366  return tmp<fluxFieldType>
367  (
368  new fluxFieldType
369  (
370  IOobject
371  (
372  "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
373  mesh().time().timeName(),
374  mesh()
375  ),
376  this->fvcDdtPhiCoeff(U.oldTime(), phiUf0, phiCorr)
377  *rDeltaT*phiCorr
378  )
379  );
380 }
381 
382 
383 template<class Type>
386 (
388  const fluxFieldType& phi
389 )
390 {
391  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
392 
393  fluxFieldType phiCorr
394  (
395  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
396  );
397 
398  return tmp<fluxFieldType>
399  (
400  new fluxFieldType
401  (
402  IOobject
403  (
404  "ddtCorr(" + U.name() + ',' + phi.name() + ')',
405  mesh().time().timeName(),
406  mesh()
407  ),
408  this->fvcDdtPhiCoeff(U.oldTime(), phi.oldTime(), phiCorr)
409  *rDeltaT*phiCorr
410  )
411  );
412 }
413 
414 
415 template<class Type>
418 (
419  const volScalarField& rho,
422 )
423 {
424  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
425 
426  if
427  (
428  U.dimensions() == dimVelocity
430  )
431  {
433  (
434  rho.oldTime()*U.oldTime()
435  );
436 
437  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
438  fluxFieldType phiCorr(phiUf0 - fvc::dotInterpolate(mesh().Sf(), rhoU0));
439 
440  return tmp<fluxFieldType>
441  (
442  new fluxFieldType
443  (
444  IOobject
445  (
446  "ddtCorr("
447  + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
448  mesh().time().timeName(),
449  mesh()
450  ),
451  this->fvcDdtPhiCoeff(rhoU0, phiUf0, phiCorr, rho.oldTime())
452  *rDeltaT*phiCorr
453  )
454  );
455  }
456  else if
457  (
460  )
461  {
462  fluxFieldType phiUf0(mesh().Sf() & Uf.oldTime());
463  fluxFieldType phiCorr
464  (
465  phiUf0 - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
466  );
467 
468  return tmp<fluxFieldType>
469  (
470  new fluxFieldType
471  (
472  IOobject
473  (
474  "ddtCorr("
475  + rho.name() + ',' + U.name() + ',' + Uf.name() + ')',
476  mesh().time().timeName(),
477  mesh()
478  ),
479  this->fvcDdtPhiCoeff
480  (
481  U.oldTime(),
482  phiUf0,
483  phiCorr,
484  rho.oldTime()
485  )*rDeltaT*phiCorr
486  )
487  );
488  }
489  else
490  {
492  << "dimensions of Uf are not correct"
493  << abort(FatalError);
494 
495  return fluxFieldType::null();
496  }
497 }
498 
499 
500 template<class Type>
503 (
504  const volScalarField& rho,
506  const fluxFieldType& phi
507 )
508 {
509  const surfaceScalarField rDeltaT(fvc::interpolate(localRDeltaT()));
510 
511  if
512  (
513  U.dimensions() == dimVelocity
514  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
515  )
516  {
518  (
519  rho.oldTime()*U.oldTime()
520  );
521 
522  fluxFieldType phiCorr
523  (
524  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), rhoU0)
525  );
526 
527  return tmp<fluxFieldType>
528  (
529  new fluxFieldType
530  (
531  IOobject
532  (
533  "ddtCorr("
534  + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
535  mesh().time().timeName(),
536  mesh()
537  ),
538  this->fvcDdtPhiCoeff
539  (
540  rhoU0,
541  phi.oldTime(),
542  phiCorr,
543  rho.oldTime()
544  )*rDeltaT*phiCorr
545  )
546  );
547  }
548  else if
549  (
550  U.dimensions() == rho.dimensions()*dimVelocity
551  && phi.dimensions() == rho.dimensions()*dimVelocity*dimArea
552  )
553  {
554  fluxFieldType phiCorr
555  (
556  phi.oldTime() - fvc::dotInterpolate(mesh().Sf(), U.oldTime())
557  );
558 
559  return tmp<fluxFieldType>
560  (
561  new fluxFieldType
562  (
563  IOobject
564  (
565  "ddtCorr("
566  + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
567  mesh().time().timeName(),
568  mesh()
569  ),
570  this->fvcDdtPhiCoeff
571  (
572  U.oldTime(),
573  phi.oldTime(),
574  phiCorr,
575  rho.oldTime()
576  )*rDeltaT*phiCorr
577  )
578  );
579  }
580  else
581  {
583  << "dimensions of phi are not correct"
584  << abort(FatalError);
585 
586  return fluxFieldType::null();
587  }
588 }
589 
590 
591 template<class Type>
593 (
595 )
596 {
598  (
600  (
601  IOobject
602  (
603  "meshPhi",
604  mesh().time().timeName(),
605  mesh(),
608  false
609  ),
610  mesh(),
612  )
613  );
614 }
615 
616 
617 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
618 
619 } // End namespace fv
620 
621 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
622 
623 } // End namespace Foam
624 
625 // ************************************************************************* //
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
const word & name() const
Return name.
Definition: IOobject.H:297
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.
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.
Generic GeometricField class.
Generic dimensioned Type class.
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)
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
word timeName
Definition: getTimeIndex.H:3
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:294
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
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.
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