FreeFOAM The Cross-Platform CFD Toolkit
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) 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 "localEulerDdtScheme.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 {
45  return mesh().objectRegistry::lookupObject<volScalarField>(rDeltaTName_);
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,
78  )
79  )
80  );
81 
82  tdtdt().internalField() =
83  rDeltaT.internalField()*dt.value()*(1.0 - mesh().V0()/mesh().V());
84 
85  return tdtdt;
86  }
87  else
88  {
90  (
92  (
93  ddtIOobject,
94  mesh(),
96  (
97  "0",
98  dt.dimensions()/dimTime,
100  ),
102  )
103  );
104  }
105 }
106 
107 
108 template<class Type>
111 (
113 )
114 {
115  const volScalarField& rDeltaT = localRDeltaT();
116 
117  IOobject ddtIOobject
118  (
119  "ddt(" + vf.name() + ')',
120  mesh().time().timeName(),
121  mesh()
122  );
123 
124  if (mesh().moving())
125  {
127  (
129  (
130  ddtIOobject,
131  mesh(),
132  rDeltaT.dimensions()*vf.dimensions(),
133  rDeltaT.internalField()*
134  (
135  vf.internalField()
136  - vf.oldTime().internalField()*mesh().V0()/mesh().V()
137  ),
138  rDeltaT.boundaryField()*
139  (
140  vf.boundaryField() - vf.oldTime().boundaryField()
141  )
142  )
143  );
144  }
145  else
146  {
148  (
150  (
151  ddtIOobject,
152  rDeltaT*(vf - vf.oldTime())
153  )
154  );
155  }
156 }
157 
158 
159 template<class Type>
162 (
163  const dimensionedScalar& rho,
165 )
166 {
167  const volScalarField& rDeltaT = localRDeltaT();
168 
169  IOobject ddtIOobject
170  (
171  "ddt(" + rho.name() + ',' + vf.name() + ')',
172  mesh().time().timeName(),
173  mesh()
174  );
175 
176  if (mesh().moving())
177  {
179  (
181  (
182  ddtIOobject,
183  mesh(),
184  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
185  rDeltaT.internalField()*rho.value()*
186  (
187  vf.internalField()
188  - vf.oldTime().internalField()*mesh().V0()/mesh().V()
189  ),
190  rDeltaT.boundaryField()*rho.value()*
191  (
192  vf.boundaryField() - vf.oldTime().boundaryField()
193  )
194  )
195  );
196  }
197  else
198  {
200  (
202  (
203  ddtIOobject,
204  rDeltaT*rho*(vf - vf.oldTime())
205  )
206  );
207  }
208 }
209 
210 
211 template<class Type>
214 (
215  const volScalarField& rho,
217 )
218 {
219  const volScalarField& rDeltaT = localRDeltaT();
220 
221  IOobject ddtIOobject
222  (
223  "ddt(" + rho.name() + ',' + vf.name() + ')',
224  mesh().time().timeName(),
225  mesh()
226  );
227 
228  if (mesh().moving())
229  {
231  (
233  (
234  ddtIOobject,
235  mesh(),
236  rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
237  rDeltaT.internalField()*
238  (
239  rho.internalField()*vf.internalField()
240  - rho.oldTime().internalField()
241  *vf.oldTime().internalField()*mesh().V0()/mesh().V()
242  ),
243  rDeltaT.boundaryField()*
244  (
245  rho.boundaryField()*vf.boundaryField()
246  - rho.oldTime().boundaryField()
247  *vf.oldTime().boundaryField()
248  )
249  )
250  );
251  }
252  else
253  {
255  (
257  (
258  ddtIOobject,
259  rDeltaT*(rho*vf - rho.oldTime()*vf.oldTime())
260  )
261  );
262  }
263 }
264 
265 
266 template<class Type>
269 (
271 )
272 {
273  tmp<fvMatrix<Type> > tfvm
274  (
275  new fvMatrix<Type>
276  (
277  vf,
279  )
280  );
281 
282  fvMatrix<Type>& fvm = tfvm();
283 
284  const scalarField& rDeltaT = localRDeltaT().internalField();
285 
286  fvm.diag() = rDeltaT*mesh().V();
287 
288  if (mesh().moving())
289  {
290  fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V0();
291  }
292  else
293  {
294  fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V();
295  }
296 
297  return tfvm;
298 }
299 
300 
301 template<class Type>
304 (
305  const dimensionedScalar& rho,
307 )
308 {
309  tmp<fvMatrix<Type> > tfvm
310  (
311  new fvMatrix<Type>
312  (
313  vf,
315  )
316  );
317  fvMatrix<Type>& fvm = tfvm();
318 
319  const scalarField& rDeltaT = localRDeltaT().internalField();
320 
321  fvm.diag() = rDeltaT*rho.value()*mesh().V();
322 
323  if (mesh().moving())
324  {
325  fvm.source() = rDeltaT
326  *rho.value()*vf.oldTime().internalField()*mesh().V0();
327  }
328  else
329  {
330  fvm.source() = rDeltaT
331  *rho.value()*vf.oldTime().internalField()*mesh().V();
332  }
333 
334  return tfvm;
335 }
336 
337 
338 template<class Type>
341 (
342  const volScalarField& rho,
344 )
345 {
346  tmp<fvMatrix<Type> > tfvm
347  (
348  new fvMatrix<Type>
349  (
350  vf,
352  )
353  );
354  fvMatrix<Type>& fvm = tfvm();
355 
356  const scalarField& rDeltaT = localRDeltaT().internalField();
357 
358  fvm.diag() = rDeltaT*rho.internalField()*mesh().V();
359 
360  if (mesh().moving())
361  {
362  fvm.source() = rDeltaT
363  *rho.oldTime().internalField()
364  *vf.oldTime().internalField()*mesh().V0();
365  }
366  else
367  {
368  fvm.source() = rDeltaT
369  *rho.oldTime().internalField()
370  *vf.oldTime().internalField()*mesh().V();
371  }
372 
373  return tfvm;
374 }
375 
376 
377 template<class Type>
380 (
381  const volScalarField& rA,
383  const fluxFieldType& phi
384 )
385 {
386  IOobject ddtIOobject
387  (
388  "ddtPhiCorr(" + rA.name() + ',' + U.name() + ',' + phi.name() + ')',
389  mesh().time().timeName(),
390  mesh()
391  );
392 
393  if (mesh().moving())
394  {
395  return tmp<fluxFieldType>
396  (
397  new fluxFieldType
398  (
399  ddtIOobject,
400  mesh(),
401  dimensioned<typename flux<Type>::type>
402  (
403  "0",
404  rA.dimensions()*phi.dimensions()/dimTime,
406  )
407  )
408  );
409  }
410  else
411  {
412  const volScalarField& rDeltaT = localRDeltaT();
413 
414  return tmp<fluxFieldType>
415  (
416  new fluxFieldType
417  (
418  ddtIOobject,
419  fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())*
420  (
421  fvc::interpolate(rDeltaT*rA)*phi.oldTime()
422  - (fvc::interpolate(rDeltaT*rA*U.oldTime()) & mesh().Sf())
423  )
424  )
425  );
426  }
427 }
428 
429 
430 template<class Type>
433 (
434  const volScalarField& rA,
435  const volScalarField& rho,
437  const fluxFieldType& phi
438 )
439 {
440  IOobject ddtIOobject
441  (
442  "ddtPhiCorr("
443  + rA.name() + ',' + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
444  mesh().time().timeName(),
445  mesh()
446  );
447 
448  if (mesh().moving())
449  {
450  return tmp<fluxFieldType>
451  (
452  new fluxFieldType
453  (
454  ddtIOobject,
455  mesh(),
456  dimensioned<typename flux<Type>::type>
457  (
458  "0",
459  rA.dimensions()*rho.dimensions()*phi.dimensions()/dimTime,
461  )
462  )
463  );
464  }
465  else
466  {
467  const volScalarField& rDeltaT = localRDeltaT();
468 
469  if
470  (
471  U.dimensions() == dimVelocity
472  && phi.dimensions() == dimVelocity*dimArea
473  )
474  {
475  return tmp<fluxFieldType>
476  (
477  new fluxFieldType
478  (
479  ddtIOobject,
480 // (
481 // scalar(1)
482 // - (min(rDeltaT)/fvc::interpolate(rDeltaT))
483 // *(scalar(1) - fvcDdtPhiCoeff(U.oldTime(), phi.oldTime()))
484 // )
485  fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())
486  //0.95
487  *(
488  fvc::interpolate(rDeltaT*rA*rho.oldTime())*phi.oldTime()
489  - (fvc::interpolate(rDeltaT*rA*rho.oldTime()*U.oldTime())
490  & mesh().Sf())
491  )
492  )
493  );
494  }
495  else if
496  (
497  U.dimensions() == dimVelocity
499  )
500  {
501  return tmp<fluxFieldType>
502  (
503  new fluxFieldType
504  (
505  ddtIOobject,
506  fvcDdtPhiCoeff
507  (
508  U.oldTime(),
509  phi.oldTime()/fvc::interpolate(rho.oldTime())
510  )
511  *(
512  fvc::interpolate(rDeltaT*rA*rho.oldTime())
513  *phi.oldTime()/fvc::interpolate(rho.oldTime())
514  - (
516  (
517  rDeltaT*rA*rho.oldTime()*U.oldTime()
518  ) & mesh().Sf()
519  )
520  )
521  )
522  );
523  }
524  else if
525  (
528  )
529  {
530  return tmp<fluxFieldType>
531  (
532  new fluxFieldType
533  (
534  ddtIOobject,
535  fvcDdtPhiCoeff(rho.oldTime(), U.oldTime(), phi.oldTime())
536  *(
537  fvc::interpolate(rDeltaT*rA)*phi.oldTime()
538  - (
539  fvc::interpolate(rDeltaT*rA*U.oldTime())&mesh().Sf()
540  )
541  )
542  )
543  );
544  }
545  else
546  {
548  (
549  "localEulerDdtScheme<Type>::fvcDdtPhiCorr"
550  ) << "dimensions of phi are not correct"
551  << abort(FatalError);
552 
553  return fluxFieldType::null();
554  }
555  }
556 }
557 
558 
559 template<class Type>
561 (
563 )
564 {
566  (
568  (
569  IOobject
570  (
571  "meshPhi",
572  mesh().time().timeName(),
573  mesh()
574  ),
575  mesh(),
577  )
578  );
579 }
580 
581 
582 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
583 
584 } // End namespace fv
585 
586 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
587 
588 } // End namespace Foam
589 
590 // ************************************************************************* //