MomentumTransferPhaseSystem.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) 2015 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 
27 
28 #include "BlendedInterfacialModel.H"
29 #include "dragModel.H"
30 #include "virtualMassModel.H"
31 #include "liftModel.H"
32 #include "wallLubricationModel.H"
33 #include "turbulentDispersionModel.H"
34 
35 #include "HashPtrTable.H"
36 
37 #include "fvmDdt.H"
38 #include "fvmDiv.H"
39 #include "fvmSup.H"
40 #include "fvcDiv.H"
41 #include "fvcSnGrad.H"
42 #include "fvMatrix.H"
43 
44 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45 
46 template<class BasePhaseSystem>
49 (
50  const fvMesh& mesh
51 )
52 :
53  BasePhaseSystem(mesh)
54 {
55  this->generatePairsAndSubModels
56  (
57  "drag",
58  dragModels_
59  );
60 
61  this->generatePairsAndSubModels
62  (
63  "virtualMass",
64  virtualMassModels_
65  );
66 
67  this->generatePairsAndSubModels
68  (
69  "lift",
70  liftModels_
71  );
72 
73  this->generatePairsAndSubModels
74  (
75  "wallLubrication",
76  wallLubricationModels_
77  );
78 
79  this->generatePairsAndSubModels
80  (
81  "turbulentDispersion",
82  turbulentDispersionModels_
83  );
84 
86  (
87  dragModelTable,
88  dragModels_,
89  dragModelIter
90  )
91  {
92  const phasePair& pair(this->phasePairs_[dragModelIter.key()]);
93 
94  Kds_.insert
95  (
96  pair,
97  new volScalarField
98  (
99  IOobject::groupName("Kd", pair.name()),
100  dragModelIter()->K()
101  )
102  );
103  }
104 
106  (
107  virtualMassModelTable,
108  virtualMassModels_,
109  virtualMassModelIter
110  )
111  {
112  const phasePair& pair(this->phasePairs_[virtualMassModelIter.key()]);
113 
114  Vms_.insert
115  (
116  pair,
117  new volScalarField
118  (
119  IOobject::groupName("Vm", pair.name()),
120  virtualMassModelIter()->K()
121  )
122  );
123  }
124 }
125 
126 
127 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
128 
129 template<class BasePhaseSystem>
132 {}
133 
134 
135 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
136 
137 template<class BasePhaseSystem>
140 (
141  const phasePairKey& key
142 ) const
143 {
144  return dragModels_[key]->K();
145 }
146 
147 
148 template<class BasePhaseSystem>
151 (
152  const phasePairKey& key
153 ) const
154 {
155  return dragModels_[key]->Kf();
156 }
157 
158 
159 template<class BasePhaseSystem>
162 (
163  const Foam::phaseModel& phase
164 ) const
165 {
166  tmp<volScalarField> tKd
167  (
168  new volScalarField
169  (
170  IOobject
171  (
172  IOobject::groupName("Kd", phase.name()),
173  this->mesh_.time().timeName(),
174  this->mesh_
175  ),
176  this->mesh_,
178  (
179  IOobject::groupName("Kd", phase.name()),
180  dimensionSet(1, -3, -1, 0, 0),
181  0
182  )
183  )
184  );
185 
187  (
189  Kds_,
190  KdIter
191  )
192  {
193  const volScalarField& K(*KdIter());
194 
195  const phasePair& pair(this->phasePairs_[KdIter.key()]);
196 
197  const phaseModel* phase1 = &pair.phase1();
198  const phaseModel* phase2 = &pair.phase2();
199 
200  forAllConstIter(phasePair, pair, iter)
201  {
202  if (phase1 == &phase)
203  {
204  tKd() += K;
205  }
206 
207  Swap(phase1, phase2);
208  }
209  }
210 
211  return tKd;
212 }
213 
214 
215 template<class BasePhaseSystem>
218 (
219  const phasePairKey& key
220 ) const
221 {
222  if (virtualMassModels_.found(key))
223  {
224  return virtualMassModels_[key]->K();
225  }
226  else
227  {
228  return tmp<volScalarField>
229  (
230  new volScalarField
231  (
232  IOobject
233  (
234  virtualMassModel::typeName + ":K",
235  this->mesh_.time().timeName(),
236  this->mesh_,
239  false
240  ),
241  this->mesh_,
243  )
244  );
245  }
246 }
247 
248 
249 template<class BasePhaseSystem>
252 (
253  const phasePairKey& key
254 ) const
255 {
256  if (virtualMassModels_.found(key))
257  {
258  return virtualMassModels_[key]->Kf();
259  }
260  else
261  {
262  return tmp<surfaceScalarField>
263  (
265  (
266  IOobject
267  (
268  virtualMassModel::typeName + ":Kf",
269  this->mesh_.time().timeName(),
270  this->mesh_,
273  false
274  ),
275  this->mesh_,
277  )
278  );
279  }
280 }
281 
282 
283 template<class BasePhaseSystem>
286 (
287  const phasePairKey& key
288 ) const
289 {
290  if (liftModels_.found(key) && wallLubricationModels_.found(key))
291  {
292  return
293  liftModels_[key]->template F<vector>()
294  + wallLubricationModels_[key]->template F<vector>();
295  }
296  else if (liftModels_.found(key))
297  {
298  return liftModels_[key]->template F<vector>();
299  }
300  else if (wallLubricationModels_.found(key))
301  {
302  return wallLubricationModels_[key]->template F<vector>();
303  }
304  else
305  {
306  return tmp<volVectorField>
307  (
308  new volVectorField
309  (
310  IOobject
311  (
312  liftModel::typeName + ":F",
313  this->mesh_.time().timeName(),
314  this->mesh_,
317  false
318  ),
319  this->mesh_,
321  )
322  );
323  }
324 }
325 
326 
327 template<class BasePhaseSystem>
330 (
331  const phasePairKey& key
332 ) const
333 {
334  if (liftModels_.found(key) && wallLubricationModels_.found(key))
335  {
336  return
337  liftModels_[key]->Ff()
338  + wallLubricationModels_[key]->Ff();
339  }
340  else if (liftModels_.found(key))
341  {
342  return liftModels_[key]->Ff();
343  }
344  else if (wallLubricationModels_.found(key))
345  {
346  return wallLubricationModels_[key]->Ff();
347  }
348  else
349  {
350  return tmp<surfaceScalarField>
351  (
353  (
354  IOobject
355  (
356  liftModel::typeName + ":Ff",
357  this->mesh_.time().timeName(),
358  this->mesh_,
361  false
362  ),
363  this->mesh_,
365  )
366  );
367  }
368 }
369 
370 
371 template<class BasePhaseSystem>
374 (
375  const phasePairKey& key
376 ) const
377 {
378  if (turbulentDispersionModels_.found(key))
379  {
380  return turbulentDispersionModels_[key]->D();
381  }
382  else
383  {
384  return tmp<volScalarField>
385  (
386  new volScalarField
387  (
388  IOobject
389  (
390  turbulentDispersionModel::typeName + ":D",
391  this->mesh_.time().timeName(),
392  this->mesh_,
395  false
396  ),
397  this->mesh_,
399  )
400  );
401  }
402 }
403 
404 
405 template<class BasePhaseSystem>
408 {
409  // Create a momentum transfer matrix for each phase
410  autoPtr<phaseSystem::momentumTransferTable> eqnsPtr
411  (
413  );
414 
415  phaseSystem::momentumTransferTable& eqns = eqnsPtr();
416 
417  forAll(this->phaseModels_, phasei)
418  {
419  const phaseModel& phase = this->phaseModels_[phasei];
420 
421  eqns.insert
422  (
423  phase.name(),
424  new fvVectorMatrix(phase.U(), dimMass*dimVelocity/dimTime)
425  );
426  }
427 
428  // Update the drag coefficients
430  (
431  dragModelTable,
432  dragModels_,
433  dragModelIter
434  )
435  {
436  *Kds_[dragModelIter.key()] = dragModelIter()->K();
437  }
438 
439  // Add the implicit part of the drag force
441  (
443  Kds_,
444  KdIter
445  )
446  {
447  const volScalarField& K(*KdIter());
448 
449  const phasePair& pair(this->phasePairs_[KdIter.key()]);
450 
451  const phaseModel* phase = &pair.phase1();
452  const phaseModel* otherPhase = &pair.phase2();
453 
454  forAllConstIter(phasePair, pair, iter)
455  {
456  const volVectorField& U = phase->U();
457 
458  *eqns[phase->name()] -= fvm::Sp(K, U);
459 
460  Swap(phase, otherPhase);
461  }
462  }
463 
464  // Update the virtual mass coefficients
466  (
467  virtualMassModelTable,
468  virtualMassModels_,
469  virtualMassModelIter
470  )
471  {
472  *Vms_[virtualMassModelIter.key()] = virtualMassModelIter()->K();
473  }
474 
475  // Add the virtual mass force
477  (
479  Vms_,
480  VmIter
481  )
482  {
483  const volScalarField& Vm(*VmIter());
484 
485  const phasePair& pair(this->phasePairs_[VmIter.key()]);
486 
487  const phaseModel* phase = &pair.phase1();
488  const phaseModel* otherPhase = &pair.phase2();
489 
490  forAllConstIter(phasePair, pair, iter)
491  {
492  const volVectorField& U = phase->U();
493  const surfaceScalarField& phi = phase->phi();
494 
495  *eqns[phase->name()] -=
496  Vm
497  *(
498  fvm::ddt(U)
499  + fvm::div(phi, U)
500  - fvm::Sp(fvc::div(phi), U)
501  - otherPhase->DUDt()
502  )
503  + this->MRF_.DDt(Vm, U - otherPhase->U());
504 
505  Swap(phase, otherPhase);
506  }
507  }
508 
509  return eqnsPtr;
510 }
511 
512 
513 template<class BasePhaseSystem>
515 (
516  PtrList<volVectorField>& Fs, const label phasei
517 ) const
518 {
519  if (!Fs.set(phasei))
520  {
521  Fs.set
522  (
523  phasei,
524  new volVectorField
525  (
526  IOobject
527  (
528  liftModel::typeName + ":F",
529  this->mesh_.time().timeName(),
530  this->mesh_,
533  false
534  ),
535  this->mesh_,
537  )
538  );
539  }
540 
541  return Fs[phasei];
542 }
543 
544 
545 template<class BasePhaseSystem>
548 {
549  autoPtr<PtrList<volVectorField> > tFs
550  (
551  new PtrList<volVectorField>(this->phases().size())
552  );
553  PtrList<volVectorField>& Fs = tFs();
554 
555  // Add the lift force
557  (
558  liftModelTable,
559  liftModels_,
560  liftModelIter
561  )
562  {
563  const volVectorField F(liftModelIter()->F<vector>());
564 
565  const phasePair& pair(this->phasePairs_[liftModelIter.key()]);
566 
567  setF(Fs, pair.phase1().index()) += F;
568  setF(Fs, pair.phase2().index()) -= F;
569  }
570 
571  // Add the wall lubrication force
573  (
574  wallLubricationModelTable,
575  wallLubricationModels_,
576  wallLubricationModelIter
577  )
578  {
579  const volVectorField F(wallLubricationModelIter()->F<vector>());
580 
581  const phasePair&
582  pair(this->phasePairs_[wallLubricationModelIter.key()]);
583 
584  setF(Fs, pair.phase1().index()) += F;
585  setF(Fs, pair.phase2().index()) -= F;
586  }
587 
588  return tFs;
589 }
590 
591 
592 template<class BasePhaseSystem>
595 (
596  PtrList<surfaceScalarField>& phiDs, const label phasei
597 ) const
598 {
599  if (!phiDs.set(phasei))
600  {
601  phiDs.set
602  (
603  phasei,
605  (
606  IOobject
607  (
608  turbulentDispersionModel::typeName + ":phiD",
609  this->mesh_.time().timeName(),
610  this->mesh_,
613  false
614  ),
615  this->mesh_,
617  (
618  "zero",
620  0
621  )
622  )
623  );
624  }
625 
626  return phiDs[phasei];
627 }
628 
629 
630 template<class BasePhaseSystem>
633 (
634  const PtrList<volScalarField>& rAUs
635 ) const
636 {
637  autoPtr<PtrList<surfaceScalarField> > tphiDs
638  (
639  new PtrList<surfaceScalarField>(this->phases().size())
640  );
641  PtrList<surfaceScalarField>& phiDs = tphiDs();
642 
643  // Add the turbulent dispersion force
645  (
646  turbulentDispersionModelTable,
647  turbulentDispersionModels_,
648  turbulentDispersionModelIter
649  )
650  {
651  const phasePair&
652  pair(this->phasePairs_[turbulentDispersionModelIter.key()]);
653 
654  const volScalarField D(turbulentDispersionModelIter()->D());
656  (
657  fvc::snGrad(pair.phase1())*this->mesh_.magSf()
658  );
659 
660  setPhiD(phiDs, pair.phase1().index()) +=
661  fvc::interpolate(rAUs[pair.phase1().index()]*D)*snGradAlpha1;
662  setPhiD(phiDs, pair.phase2().index()) -=
663  fvc::interpolate(rAUs[pair.phase2().index()]*D)*snGradAlpha1;
664  }
665 
666  return tphiDs;
667 }
668 
669 
670 template<class BasePhaseSystem>
672 {
673  if (BasePhaseSystem::read())
674  {
675  bool readOK = true;
676 
677  // Read models ...
678 
679  return readOK;
680  }
681  else
682  {
683  return false;
684  }
685 }
686 
687 
688 // ************************************************************************* //
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &vf, const surfaceScalarField &faceFlux, Istream &schemeData)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
static const dimensionSet dimF
Force dimensions.
const word & name() const
Definition: phaseModel.H:109
Calculate the snGrad of the given volField.
HashPtrTable< volScalarField, phasePairKey, phasePairKey::hash > KdTable
Definition: phaseSystem.H:82
Calculate the matrix for implicit and explicit sources.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
label phasei
Definition: pEqn.H:37
CGAL::Exact_predicates_exact_constructions_kernel K
virtual tmp< volScalarField > D(const phasePairKey &key) const
Return the turbulent diffusivity.
virtual tmp< volScalarField > Vm(const phasePairKey &key) const
Return the virtual mass coefficient.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
virtual tmp< volScalarField > Kd(const phasePairKey &key) const
Return the drag coefficient.
virtual tmp< surfaceScalarField > Ff(const phasePairKey &key) const
Return the combined face-force (lift + wall-lubrication)
static const dimensionSet dimK
Coefficient dimensions.
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Class which models interfacial momenum transfer between a number of phases. Drag, virtual mass...
Calculate the matrix for the divergence of the given field and flux.
Calulate the matrix for the first temporal derivative.
virtual bool read()
Read base phaseProperties dictionary.
MomentumTransferPhaseSystem(const fvMesh &)
Construct from fvMesh.
virtual autoPtr< PtrList< volVectorField > > Fs() const
Return the combined force (lift + wall-lubrication)
const dimensionSet dimDensity
virtual autoPtr< PtrList< surfaceScalarField > > phiDs(const PtrList< volScalarField > &rAUs) const
Return the turbulent dispersion force on faces for phase pair.
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
fvMatrix< vector > fvVectorMatrix
Definition: fvMatricesFwd.H:45
static word groupName(Name name, const word &group)
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
surfaceScalarField snGradAlpha1(fvc::snGrad(alpha1)*mesh.magSf())
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
#define forAll(list, i)
Definition: UList.H:421
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
virtual tmp< surfaceScalarField > Vmf(const phasePairKey &key) const
Return the face virtual mass coefficient.
HashPtrTable< volScalarField, phasePairKey, phasePairKey::hash > VmTable
Definition: phaseSystem.H:91
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
static const dimensionSet dimF
Force dimensions.
Definition: liftModel.H:86
HashPtrTable< fvVectorMatrix, word, string::hash > momentumTransferTable
Definition: phaseSystem.H:100
void Swap(T &a, T &b)
Definition: Swap.H:43
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
static const dimensionSet dimD
Diffusivity dimensions.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
Calculate the divergence of the given field.
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransfer() const
Return the momentum transfer matrices.
virtual tmp< surfaceScalarField > Kdf(const phasePairKey &key) const
Return the face drag coefficient.
static const Vector zero
Definition: Vector.H:80
bool read(const char *, int32_t &)
Definition: int32IO.C:87
multiphaseSystem::phaseModelList & phases
Definition: createFields.H:11
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
const dimensionedScalar F
Faraday constant: default SI units: [C/mol].
const dimensionSet dimVelocity
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
A class for managing temporary objects.
Definition: PtrList.H:118
virtual tmp< volVectorField > F(const phasePairKey &key) const
Return the combined force (lift + wall-lubrication)
virtual ~MomentumTransferPhaseSystem()
Destructor.