FreeFOAM The Cross-Platform CFD Toolkit
EulerDdtScheme.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) 1991-2010 OpenCFD Ltd.
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"
28 #include <finiteVolume/fvcDiv.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  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
51 
52  IOobject ddtIOobject
53  (
54  "ddt("+dt.name()+')',
55  mesh().time().timeName(),
56  mesh()
57  );
58 
59  if (mesh().moving())
60  {
62  (
64  (
65  ddtIOobject,
66  mesh(),
68  (
69  "0",
70  dt.dimensions()/dimTime,
72  )
73  )
74  );
75 
76  tdtdt().internalField() =
77  rDeltaT.value()*dt.value()*(1.0 - mesh().V0()/mesh().V());
78 
79  return tdtdt;
80  }
81  else
82  {
84  (
86  (
87  ddtIOobject,
88  mesh(),
90  (
91  "0",
92  dt.dimensions()/dimTime,
94  ),
96  )
97  );
98  }
99 }
100 
101 
102 template<class Type>
105 (
107 )
108 {
109  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
110 
111  IOobject ddtIOobject
112  (
113  "ddt("+vf.name()+')',
114  mesh().time().timeName(),
115  mesh()
116  );
117 
118  if (mesh().moving())
119  {
121  (
123  (
124  ddtIOobject,
125  mesh(),
126  rDeltaT.dimensions()*vf.dimensions(),
127  rDeltaT.value()*
128  (
129  vf.internalField()
130  - vf.oldTime().internalField()*mesh().V0()/mesh().V()
131  ),
132  rDeltaT.value()*
133  (
134  vf.boundaryField() - vf.oldTime().boundaryField()
135  )
136  )
137  );
138  }
139  else
140  {
142  (
144  (
145  ddtIOobject,
146  rDeltaT*(vf - vf.oldTime())
147  )
148  );
149  }
150 }
151 
152 
153 template<class Type>
156 (
157  const dimensionedScalar& rho,
159 )
160 {
161  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
162 
163  IOobject ddtIOobject
164  (
165  "ddt("+rho.name()+','+vf.name()+')',
166  mesh().time().timeName(),
167  mesh()
168  );
169 
170  if (mesh().moving())
171  {
173  (
175  (
176  ddtIOobject,
177  mesh(),
178  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
179  rDeltaT.value()*rho.value()*
180  (
181  vf.internalField()
182  - vf.oldTime().internalField()*mesh().V0()/mesh().V()
183  ),
184  rDeltaT.value()*rho.value()*
185  (
186  vf.boundaryField() - vf.oldTime().boundaryField()
187  )
188  )
189  );
190  }
191  else
192  {
194  (
196  (
197  ddtIOobject,
198  rDeltaT*rho*(vf - vf.oldTime())
199  )
200  );
201  }
202 }
203 
204 
205 template<class Type>
208 (
209  const volScalarField& rho,
211 )
212 {
213  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
214 
215  IOobject ddtIOobject
216  (
217  "ddt("+rho.name()+','+vf.name()+')',
218  mesh().time().timeName(),
219  mesh()
220  );
221 
222  if (mesh().moving())
223  {
225  (
227  (
228  ddtIOobject,
229  mesh(),
230  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
231  rDeltaT.value()*
232  (
233  rho.internalField()*vf.internalField()
234  - rho.oldTime().internalField()
235  *vf.oldTime().internalField()*mesh().V0()/mesh().V()
236  ),
237  rDeltaT.value()*
238  (
239  rho.boundaryField()*vf.boundaryField()
240  - rho.oldTime().boundaryField()
241  *vf.oldTime().boundaryField()
242  )
243  )
244  );
245  }
246  else
247  {
249  (
251  (
252  ddtIOobject,
253  rDeltaT*(rho*vf - rho.oldTime()*vf.oldTime())
254  )
255  );
256  }
257 }
258 
259 
260 template<class Type>
263 (
265 )
266 {
267  tmp<fvMatrix<Type> > tfvm
268  (
269  new fvMatrix<Type>
270  (
271  vf,
273  )
274  );
275 
276  fvMatrix<Type>& fvm = tfvm();
277 
278  scalar rDeltaT = 1.0/mesh().time().deltaT().value();
279 
280  fvm.diag() = rDeltaT*mesh().V();
281 
282  if (mesh().moving())
283  {
284  fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V0();
285  }
286  else
287  {
288  fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V();
289  }
290 
291  return tfvm;
292 }
293 
294 
295 template<class Type>
298 (
299  const dimensionedScalar& rho,
301 )
302 {
303  tmp<fvMatrix<Type> > tfvm
304  (
305  new fvMatrix<Type>
306  (
307  vf,
309  )
310  );
311  fvMatrix<Type>& fvm = tfvm();
312 
313  scalar rDeltaT = 1.0/mesh().time().deltaT().value();
314 
315  fvm.diag() = rDeltaT*rho.value()*mesh().V();
316 
317  if (mesh().moving())
318  {
319  fvm.source() = rDeltaT
320  *rho.value()*vf.oldTime().internalField()*mesh().V0();
321  }
322  else
323  {
324  fvm.source() = rDeltaT
325  *rho.value()*vf.oldTime().internalField()*mesh().V();
326  }
327 
328  return tfvm;
329 }
330 
331 
332 template<class Type>
335 (
336  const volScalarField& rho,
338 )
339 {
340  tmp<fvMatrix<Type> > tfvm
341  (
342  new fvMatrix<Type>
343  (
344  vf,
346  )
347  );
348  fvMatrix<Type>& fvm = tfvm();
349 
350  scalar rDeltaT = 1.0/mesh().time().deltaT().value();
351 
352  fvm.diag() = rDeltaT*rho.internalField()*mesh().V();
353 
354  if (mesh().moving())
355  {
356  fvm.source() = rDeltaT
357  *rho.oldTime().internalField()
358  *vf.oldTime().internalField()*mesh().V0();
359  }
360  else
361  {
362  fvm.source() = rDeltaT
363  *rho.oldTime().internalField()
364  *vf.oldTime().internalField()*mesh().V();
365  }
366 
367  return tfvm;
368 }
369 
370 
371 template<class Type>
374 (
375  const volScalarField& rA,
377  const fluxFieldType& phiAbs
378 )
379 {
380  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
381 
382  IOobject ddtIOobject
383  (
384  "ddtPhiCorr(" + rA.name() + ',' + U.name() + ',' + phiAbs.name() + ')',
385  mesh().time().timeName(),
386  mesh()
387  );
388 
389  tmp<fluxFieldType> phiCorr =
390  phiAbs.oldTime() - (fvc::interpolate(U.oldTime()) & mesh().Sf());
391 
392  return tmp<fluxFieldType>
393  (
394  new fluxFieldType
395  (
396  ddtIOobject,
397  fvcDdtPhiCoeff(U.oldTime(), phiAbs.oldTime(), phiCorr())
398  *fvc::interpolate(rDeltaT*rA)*phiCorr
399  )
400  );
401 }
402 
403 
404 template<class Type>
407 (
408  const volScalarField& rA,
409  const volScalarField& rho,
411  const fluxFieldType& phiAbs
412 )
413 {
414  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
415 
416  IOobject ddtIOobject
417  (
418  "ddtPhiCorr("
419  + rA.name() + ','
420  + rho.name() + ','
421  + U.name() + ','
422  + phiAbs.name() + ')',
423  mesh().time().timeName(),
424  mesh()
425  );
426 
427  if
428  (
429  U.dimensions() == dimVelocity
430  && phiAbs.dimensions() == dimVelocity*dimArea
431  )
432  {
433  return tmp<fluxFieldType>
434  (
435  new fluxFieldType
436  (
437  ddtIOobject,
438  rDeltaT
439  *fvcDdtPhiCoeff(U.oldTime(), phiAbs.oldTime())
440  *(
441  fvc::interpolate(rA*rho.oldTime())*phiAbs.oldTime()
442  - (fvc::interpolate(rA*rho.oldTime()*U.oldTime())
443  & mesh().Sf())
444  )
445  )
446  );
447  }
448  else if
449  (
450  U.dimensions() == dimVelocity
451  && phiAbs.dimensions() == rho.dimensions()*dimVelocity*dimArea
452  )
453  {
454  return tmp<fluxFieldType>
455  (
456  new fluxFieldType
457  (
458  ddtIOobject,
459  rDeltaT
460  *fvcDdtPhiCoeff
461  (
462  U.oldTime(),
463  phiAbs.oldTime()/fvc::interpolate(rho.oldTime())
464  )
465  *(
466  fvc::interpolate(rA*rho.oldTime())
467  *phiAbs.oldTime()/fvc::interpolate(rho.oldTime())
468  - (
470  (
471  rA*rho.oldTime()*U.oldTime()
472  ) & mesh().Sf()
473  )
474  )
475  )
476  );
477  }
478  else if
479  (
480  U.dimensions() == rho.dimensions()*dimVelocity
481  && phiAbs.dimensions() == rho.dimensions()*dimVelocity*dimArea
482  )
483  {
484  return tmp<fluxFieldType>
485  (
486  new fluxFieldType
487  (
488  ddtIOobject,
489  rDeltaT
490  *fvcDdtPhiCoeff(rho.oldTime(), U.oldTime(), phiAbs.oldTime())
491  *(
492  fvc::interpolate(rA)*phiAbs.oldTime()
493  - (fvc::interpolate(rA*U.oldTime()) & mesh().Sf())
494  )
495  )
496  );
497  }
498  else
499  {
501  (
502  "EulerDdtScheme<Type>::fvcDdtPhiCorr"
503  ) << "dimensions of phiAbs are not correct"
504  << abort(FatalError);
505 
506  return fluxFieldType::null();
507  }
508 }
509 
510 
511 template<class Type>
513 (
515 )
516 {
517  return mesh().phi();
518 }
519 
520 
521 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
522 
523 } // End namespace fv
524 
525 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
526 
527 } // End namespace Foam
528 
529 // ************************ vim: set sw=4 sts=4 et: ************************ //