EulerFaDdtScheme.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2016-2017 Wikki Ltd
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "EulerFaDdtScheme.H"
29 #include "facDiv.H"
30 #include "faMatrices.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 namespace fa
40 {
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 template<class Type>
45 tmp<GeometricField<Type, faPatchField, areaMesh>>
47 (
48  const dimensioned<Type> dt
49 )
50 {
51  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
52 
53  IOobject ddtIOobject
54  (
55  "ddt("+dt.name()+')',
56  mesh()().time().timeName(),
57  mesh()(),
60  );
61 
62  if (mesh().moving())
63  {
65  (
67  (
68  ddtIOobject,
69  mesh(),
71  )
72  );
73 
74  tdtdt.ref().primitiveFieldRef() =
75  rDeltaT.value()*dt.value()*(1.0 - mesh().S0()/mesh().S());
76 
77  return tdtdt;
78  }
79 
80 
82  (
83  ddtIOobject,
84  mesh(),
87  );
88 }
89 
90 
91 template<class Type>
94 (
95  const dimensioned<Type> dt
96 )
97 {
98  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
99 
100  IOobject ddtIOobject
101  (
102  "ddt("+dt.name()+')',
103  mesh()().time().timeName(),
104  mesh()(),
107  );
108 
110  (
112  (
113  ddtIOobject,
114  mesh(),
115  -rDeltaT*dt
116  )
117  );
118 
119  if (mesh().moving())
120  {
121  tdtdt0.ref().primitiveFieldRef() =
122  (-rDeltaT.value()*dt.value())*mesh().S0()/mesh().S();
123  }
124 
125  return tdtdt0;
126 }
127 
128 
129 template<class Type>
132 (
134 )
135 {
136  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
137 
138  IOobject ddtIOobject
139  (
140  "ddt("+vf.name()+')',
141  mesh()().time().timeName(),
142  mesh()(),
145  );
146 
147  if (mesh().moving())
148  {
150  (
152  (
153  ddtIOobject,
154  mesh(),
155  rDeltaT.dimensions()*vf.dimensions(),
156  rDeltaT.value()*
157  (
158  vf()
159  - vf.oldTime()()*mesh().S0()/mesh().S()
160  ),
161  rDeltaT.value()*
162  (
163  vf.boundaryField() - vf.oldTime().boundaryField()
164  )
165  )
166  );
167  }
168  else
169  {
171  (
173  (
174  ddtIOobject,
175  rDeltaT*(vf - vf.oldTime())
176  )
177  );
178  }
179 }
180 
181 
182 template<class Type>
185 (
187 )
188 {
189  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
190 
191  IOobject ddtIOobject
192  (
193  "ddt0("+vf.name()+')',
194  mesh()().time().timeName(),
195  mesh()(),
198  );
199 
200  if (mesh().moving())
201  {
203  (
205  (
206  ddtIOobject,
207  mesh(),
208  rDeltaT.dimensions()*vf.dimensions(),
209  (-rDeltaT.value())*vf.oldTime().internalField(),
210  (-rDeltaT.value())*vf.oldTime().boundaryField()
211  )
212  );
213  }
214  else
215  {
217  (
219  (
220  ddtIOobject,
221  (-rDeltaT)*vf.oldTime()
222  )
223  );
224  }
225 }
226 
227 
228 template<class Type>
231 (
232  const dimensionedScalar& rho,
234 )
235 {
236  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
237 
238  IOobject ddtIOobject
239  (
240  "ddt("+rho.name()+','+vf.name()+')',
241  mesh()().time().timeName(),
242  mesh()(),
245  );
246 
247  if (mesh().moving())
248  {
250  (
252  (
253  ddtIOobject,
254  mesh(),
255  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
256  rDeltaT.value()*rho.value()*
257  (
258  vf()
259  - vf.oldTime()()*mesh().S0()/mesh().S()
260  ),
261  rDeltaT.value()*rho.value()*
262  (
263  vf.boundaryField() - vf.oldTime().boundaryField()
264  )
265  )
266  );
267  }
268  else
269  {
271  (
273  (
274  ddtIOobject,
275  rDeltaT*rho*(vf - vf.oldTime())
276  )
277  );
278  }
279 }
280 
281 template<class Type>
284 (
285  const dimensionedScalar& rho,
287 )
288 {
289  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
290 
291  IOobject ddtIOobject
292  (
293  "ddt0("+rho.name()+','+vf.name()+')',
294  mesh()().time().timeName(),
295  mesh()(),
298  );
299 
300  if (mesh().moving())
301  {
303  (
305  (
306  ddtIOobject,
307  mesh(),
308  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
309  (-rDeltaT.value())*rho.value()*
310  vf.oldTime()()*mesh().S0()/mesh().S(),
311  (-rDeltaT.value())*rho.value()*
312  vf.oldTime().boundaryField()
313  )
314  );
315  }
316  else
317  {
319  (
321  (
322  ddtIOobject,
323  (-rDeltaT)*rho*vf.oldTime()
324  )
325  );
326  }
327 }
328 
329 
330 template<class Type>
333 (
334  const areaScalarField& rho,
336 )
337 {
338  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
339 
340  IOobject ddtIOobject
341  (
342  "ddt("+rho.name()+','+vf.name()+')',
343  mesh()().time().timeName(),
344  mesh()(),
347  );
348 
349  if (mesh().moving())
350  {
352  (
354  (
355  ddtIOobject,
356  mesh(),
357  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
358  rDeltaT.value()*
359  (
360  rho()*vf()
361  - rho.oldTime()()
362  *vf.oldTime()()*mesh().S0()/mesh().S()
363  ),
364  rDeltaT.value()*
365  (
366  rho.boundaryField()*vf.boundaryField()
367  - rho.oldTime().boundaryField()
368  *vf.oldTime().boundaryField()
369  )
370  )
371  );
372  }
373  else
374  {
376  (
378  (
379  ddtIOobject,
380  rDeltaT*(rho*vf - rho.oldTime()*vf.oldTime())
381  )
382  );
383  }
384 }
385 
386 
387 template<class Type>
390 (
391  const areaScalarField& rho,
393 )
394 {
395  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
396 
397  IOobject ddtIOobject
398  (
399  "ddt0("+rho.name()+','+vf.name()+')',
400  mesh()().time().timeName(),
401  mesh()(),
404  );
405 
406  if (mesh().moving())
407  {
409  (
411  (
412  ddtIOobject,
413  mesh(),
414  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
415  rDeltaT.value()*
416  (
417  - rho.oldTime()()
418  *vf.oldTime()()*mesh().S0()/mesh().S()
419  ),
420  rDeltaT.value()*
421  (
422  - rho.oldTime().boundaryField()
423  *vf.oldTime().boundaryField()
424  )
425  )
426  );
427  }
428  else
429  {
431  (
433  (
434  ddtIOobject,
435  (-rDeltaT)*rho.oldTime()*vf.oldTime()
436  )
437  );
438  }
439 }
440 
441 template<class Type>
444 (
446 )
447 {
448  tmp<faMatrix<Type>> tfam
449  (
450  new faMatrix<Type>
451  (
452  vf,
453  vf.dimensions()*dimArea/dimTime
454  )
455  );
456 
457  faMatrix<Type>& fam = tfam.ref();
458 
459  scalar rDeltaT = 1.0/mesh().time().deltaT().value();
460 
461  fam.diag() = rDeltaT*mesh().S();
462 
463  if (mesh().moving())
464  {
465  fam.source() = rDeltaT*vf.oldTime().primitiveField()*mesh().S0();
466  }
467  else
468  {
469  fam.source() = rDeltaT*vf.oldTime().primitiveField()*mesh().S();
470  }
471 
472  return tfam;
473 }
474 
475 
476 template<class Type>
479 (
480  const dimensionedScalar& rho,
482 )
483 {
484  tmp<faMatrix<Type>> tfam
485  (
486  new faMatrix<Type>
487  (
488  vf,
489  rho.dimensions()*vf.dimensions()*dimArea/dimTime
490  )
491  );
492  faMatrix<Type>& fam = tfam.ref();
493 
494  scalar rDeltaT = 1.0/mesh().time().deltaT().value();
495 
496  fam.diag() = rDeltaT*rho.value()*mesh().S();
497 
498  if (mesh().moving())
499  {
500  fam.source() = rDeltaT
501  *rho.value()*vf.oldTime().primitiveField()*mesh().S0();
502  }
503  else
504  {
505  fam.source() = rDeltaT
506  *rho.value()*vf.oldTime().primitiveField()*mesh().S();
507  }
508 
509  return tfam;
510 }
511 
512 
513 template<class Type>
516 (
517  const areaScalarField& rho,
519 )
520 {
521  tmp<faMatrix<Type>> tfam
522  (
523  new faMatrix<Type>
524  (
525  vf,
526  rho.dimensions()*vf.dimensions()*dimArea/dimTime
527  )
528  );
529  faMatrix<Type>& fam = tfam.ref();
530 
531  scalar rDeltaT = 1.0/mesh().time().deltaT().value();
532 
533  fam.diag() = rDeltaT*rho.primitiveField()*mesh().S();
534 
535  if (mesh().moving())
536  {
537  fam.source() = rDeltaT
538  *rho.oldTime().primitiveField()
539  *vf.oldTime().primitiveField()*mesh().S0();
540  }
541  else
542  {
543  fam.source() = rDeltaT
544  *rho.oldTime().primitiveField()
545  *vf.oldTime().primitiveField()*mesh().S();
546  }
547 
548  return tfam;
549 }
550 
551 
552 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
553 
554 } // End namespace fa
555 
556 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
557 
558 } // End namespace Foam
559 
560 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::fa::EulerFaDdtScheme::facDdt
tmp< GeometricField< Type, faPatchField, areaMesh > > facDdt(const dimensioned< Type >)
Definition: EulerFaDdtScheme.C:47
Foam::faMatrix
A special matrix type and solver, designed for finite area solutions of scalar equations....
Definition: faMatricesFwd.H:43
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
EulerFaDdtScheme.H
Foam::GeometricField::oldTime
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
Definition: GeometricField.C:850
Foam::dimensioned::name
const word & name() const
Return const reference to name.
Definition: dimensionedType.C:406
faMatrices.H
facDiv.H
Calculate the divergence of the given field.
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
rho
rho
Definition: readInitialConditions.H:88
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
Foam::dimArea
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
Foam::calculatedFaPatchField
Author Zeljko Tukovic, FMENA Hrvoje Jasak, Wikki Ltd.
Definition: calculatedFaPatchField.H:54
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::dimensioned
Generic dimensioned Type class.
Definition: dimensionedScalarFwd.H:42
Foam::fa::EulerFaDdtScheme::facDdt0
tmp< GeometricField< Type, faPatchField, areaMesh > > facDdt0(const dimensioned< Type >)
Definition: EulerFaDdtScheme.C:94
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
Foam::fa::EulerFaDdtScheme::famDdt
tmp< faMatrix< Type > > famDdt(const GeometricField< Type, faPatchField, areaMesh > &)
Definition: EulerFaDdtScheme.C:444
fam
Calculate the matrix for the second temporal derivative.
Foam::dimensioned::dimensions
const dimensionSet & dimensions() const
Return const reference to dimensions.
Definition: dimensionedType.C:420
Foam::GeometricField
Generic GeometricField class.
Definition: areaFieldsFwd.H:53
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62