phaseSystem.H
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) 2015-2019 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 Class
25  Foam::phaseSystem
26 
27 Description
28  Class to represent a system of phases and model interfacial transfers
29  between them.
30 
31 SourceFiles
32  phaseSystem.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef phaseSystem_H
37 #define phaseSystem_H
38 
39 #include "IOdictionary.H"
40 
41 #include "phaseModel.H"
42 #include "phasePair.H"
43 #include "orderedPhasePair.H"
44 #include "HashPtrTable.H"
45 #include "PtrListDictionary.H"
46 
47 #include "IOMRFZoneList.H"
48 #include "fvOptions.H"
49 
50 #include "volFields.H"
51 #include "surfaceFields.H"
52 #include "fvMatricesFwd.H"
53 
54 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55 
56 namespace Foam
57 {
58 
59 class blendingMethod;
60 template<class modelType> class BlendedInterfacialModel;
61 class surfaceTensionModel;
62 class aspectRatioModel;
63 
64 /*---------------------------------------------------------------------------*\
65  Class phaseSystem Declaration
66 \*---------------------------------------------------------------------------*/
67 
68 class phaseSystem
69 :
70  public IOdictionary
71 {
72 public:
73 
74  // Public Typedefs
75 
77 
79 
81 
83 
85 
86  typedef
89 
90 
91 protected:
92 
93  // Protected typedefs
94 
95  typedef
97  dictTable;
98 
99  typedef
102 
103  typedef
104  HashTable
105  <
107  phasePairKey,
109  >
111 
112  typedef
113  HashTable
114  <
116  phasePairKey,
118  >
120 
121 
122  // Protected data
123 
124  //- Reference to the mesh
125  const fvMesh& mesh_;
126 
127  //- Phase models
128  phaseModelList phaseModels_;
129 
130  //- Moving phase models
131  phaseModelPartialList movingPhaseModels_;
132 
133  //- Stationary phase models
134  phaseModelPartialList stationaryPhaseModels_;
135 
136  //- Anisothermal phase models
137  phaseModelPartialList anisothermalPhaseModels_;
138 
139  //- Multi-component phase models
140  phaseModelPartialList multiComponentPhaseModels_;
141 
142  //- Phase pairs
143  phasePairTable phasePairs_;
144 
145  //- Total volumetric flux
147 
148  //- Rate of change of pressure
150 
151  //- Optional MRF zones
153 
154  //- Blending methods
155  blendingMethodTable blendingMethods_;
156 
157 
158  // Sub Models
159 
160  //- Surface tension models
162 
163  //- Aspect ratio models
165 
166 
167  // Protected member functions
168 
169  //- Calculate and return the mixture flux
171  (
172  const phaseModelList& phaseModels
173  ) const;
174 
175  //- Generate pairs
176  void generatePairs
177  (
178  const dictTable& modelDicts
179  );
180 
181  //- Generate pairs and sub-model tables
182  template<class modelType>
183  void createSubModels
184  (
185  const dictTable& modelDicts,
186  HashTable
187  <
189  phasePairKey,
191  >& models
192  );
193 
194  //- Generate pairs and sub-model tables
195  template<class modelType>
197  (
198  const word& modelName,
199  HashTable
200  <
202  phasePairKey,
204  >& models
205  );
206 
207  //- Generate pairs and blended sub-model tables
208  template<class modelType>
210  (
211  const word& modelName,
212  HashTable
213  <
215  phasePairKey,
217  >& models,
218  const bool correctFixedFluxBCs = true
219  );
220 
221  //- Generate pairs and two-sided sub-model tables
222  template<class modelType>
224  (
225  const word& modelName,
226  HashTable
227  <
229  phasePairKey,
231  >& models,
232  const bool correctFixedFluxBCs = true
233  );
234 
235  //- Add the field to a phase-indexed list, with the given name,
236  // constructing if necessary
237  template<class GeoField>
238  void addField
239  (
240  const phaseModel& phase,
241  const word& fieldName,
242  tmp<GeoField> field,
243  PtrList<GeoField>& fieldList
244  ) const;
245 
246  //- Add the field to a phase-indexed list, with the given name,
247  // constructing if necessary
248  template<class GeoField>
249  void addField
250  (
251  const phaseModel& phase,
252  const word& fieldName,
253  const GeoField& field,
254  PtrList<GeoField>& fieldList
255  ) const;
256 
257  //- Add the field to a phase-indexed table, with the given name,
258  // constructing if necessary
259  template<class GeoField>
260  void addField
261  (
262  const phaseModel& phase,
263  const word& fieldName,
264  tmp<GeoField> field,
265  HashPtrTable<GeoField>& fieldTable
266  ) const;
267 
268  //- Add the field to a phase-indexed table, with the given name,
269  // constructing if necessary
270  template<class GeoField>
271  void addField
272  (
273  const phaseModel& phase,
274  const word& fieldName,
275  const GeoField& field,
276  HashPtrTable<GeoField>& fieldTable
277  ) const;
278 
279 
280 public:
281 
282  //- Runtime type information
283  TypeName("phaseSystem");
284 
285  //- Default name of the phase properties dictionary
286  static const word propertiesName;
287 
288 
289  // Constructors
290 
291  //- Construct from fvMesh
292  phaseSystem(const fvMesh& mesh);
293 
294 
295  //- Destructor
296  virtual ~phaseSystem();
297 
298 
299  // Member Functions
300 
301  // Access
302 
303  //- Return the mesh
304  inline const fvMesh& mesh() const;
305 
306  //- Return the phase models
307  inline const phaseModelList& phases() const;
308 
309  //- Access the phase models
310  inline phaseModelList& phases();
311 
312  //- Return the models for phases that are moving
313  inline const phaseModelPartialList& movingPhases() const;
314 
315  //- Access the models for phases that are moving
316  inline phaseModelPartialList& movingPhases();
317 
318  //- Return the models for phases that are stationary
319  inline const phaseModelPartialList& stationaryPhases() const;
320 
321  //- Access the models for phases that are stationary
322  inline phaseModelPartialList& stationaryPhases();
323 
324  //- Return the models for phases that have variable temperature
325  inline const phaseModelPartialList& anisothermalPhases() const;
326 
327  //- Access the models for phases that have variable temperature
328  inline phaseModelPartialList& anisothermalPhases();
329 
330  //- Return the models for phases that have multiple species
331  inline const phaseModelPartialList& multiComponentPhases() const;
332 
333  //- Access the models for phases that have multiple species
334  inline phaseModelPartialList& multiComponentPhases();
335 
336  //- Return the phase pairs
337  inline const phasePairTable& phasePairs() const;
338 
339  //- Return the mixture flux
340  inline const surfaceScalarField& phi() const;
341 
342  //- Access the mixture flux
343  inline surfaceScalarField& phi();
344 
345  //- Return the rate of change of the pressure
346  inline const volScalarField& dpdt() const;
347 
348  //- Access the rate of change of the pressure
349  inline volScalarField& dpdt();
350 
351  //- Return MRF zones
352  inline const IOMRFZoneList& MRF() const;
353 
354  //- Access the fvOptions
355  inline fv::options& fvOptions() const;
356 
357 
358  // Sub-model lookup
359 
360  //- Check availability of a sub model for a given phase pair
361  template<class modelType>
362  bool foundSubModel(const phasePair& key) const;
363 
364  //- Return a sub model between a phase pair
365  template<class modelType>
366  const modelType& lookupSubModel(const phasePair& key) const;
367 
368  //- Check availability of a sub model between two phases
369  template<class modelType>
370  bool foundSubModel
371  (
372  const phaseModel& dispersed,
373  const phaseModel& continuous
374  ) const;
375 
376  //- Return a sub model between two phases
377  template<class modelType>
378  const modelType& lookupSubModel
379  (
380  const phaseModel& dispersed,
381  const phaseModel& continuous
382  ) const;
383 
384  //- Check availability of a blended sub model for a given phase pair
385  template<class modelType>
386  bool foundBlendedSubModel(const phasePair& key) const;
387 
388  //- Return a blended sub model between a phase pair
389  template<class modelType>
391  lookupBlendedSubModel(const phasePair& key) const;
392 
393 
394  // Field construction
395 
396  //- Fill up gaps in a phase-indexed list of fields with zeros
397  template
398  <
399  class Type,
400  template<class> class PatchField,
401  class GeoMesh
402  >
403  void fillFields
404  (
405  const word& name,
406  const dimensionSet& dims,
408  ) const;
409 
410  //- Fill up gaps in a phase-indexed table of fields with zeros
411  template
412  <
413  class Type,
414  template<class> class PatchField,
415  class GeoMesh
416  >
417  void fillFields
418  (
419  const word& name,
420  const dimensionSet& dims,
422  fieldTable
423  ) const;
424 
425 
426  // Properties
427 
428  //- Return the mixture density
429  tmp<volScalarField> rho() const;
430 
431  //- Return the mixture velocity
432  tmp<volVectorField> U() const;
433 
434  //- Return the aspect-ratio for a pair
435  tmp<volScalarField> E(const phasePairKey& key) const;
436 
437  //- Return the surface tension coefficient for a pair
438  tmp<volScalarField> sigma(const phasePairKey& key) const;
439 
440  //- Return the mass transfer rate for a pair
441  virtual tmp<volScalarField> dmdt(const phasePairKey& key) const;
442 
443  //- Return the mass transfer rates for each phase
444  virtual PtrList<volScalarField> dmdts() const;
445 
446 
447  // Transfers
448 
449  //- Return the momentum transfer matrices for the cell-based
450  // algorithm
452 
453  //- Return the momentum transfer matrices for the face-based
454  // algorithm
456 
457  //- Return the implicit force coefficients for the face-based
458  // algorithm
459  virtual PtrList<surfaceScalarField> AFfs() const = 0;
460 
461  //- Return the force fluxes for the cell-based algorithm
463  (
465  ) = 0;
466 
467  //- Return the force fluxes for the face-based algorithm
469  (
471  ) = 0;
472 
473  //- Return the force fluxes for the cell-based algorithm
475  (
477  ) const = 0;
478 
479  //- Return the force fluxes for the face-based algorithm
481  (
483  ) const = 0;
484 
485  //- Return the explicit part of the drag force
487  (
489  ) const = 0;
490 
491  //- Solve the drag system for the new velocities and fluxes
492  virtual void partialElimination
493  (
495  ) = 0;
496 
497  //- Solve the drag system for the new fluxes
498  virtual void partialEliminationf
499  (
501  ) = 0;
502 
503  //- Return the flux corrections for the cell-based algorithm
505  (
507  const bool includeVirtualMass = false
508  ) const = 0;
509 
510  //- Return the phase diffusivities divided by the momentum
511  // coefficients
512  virtual const HashPtrTable<surfaceScalarField>& DByAfs() const = 0;
513 
514  //- Return the heat transfer matrices
515  virtual autoPtr<heatTransferTable> heatTransfer() const = 0;
516 
517  //- Return the mass transfer matrices
518  virtual autoPtr<massTransferTable> massTransfer() const = 0;
519 
520 
521  // Evolution
522 
523  //- Solve for the phase fractions
524  virtual void solve();
525 
526  //- Correct the fluid properties other than those listed below
527  virtual void correct();
528 
529  //- Correct the kinematics
530  virtual void correctKinematics();
531 
532  //- Correct the thermodynamics
533  virtual void correctThermo();
534 
535  //- Correct the turbulence
536  virtual void correctTurbulence();
537 
538  //- Correct the energy transport e.g. alphat
539  virtual void correctEnergyTransport();
540 
541 
542  // IO
543 
544  //- Read base phaseProperties dictionary
545  virtual bool read();
546 };
547 
548 
549 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
550 
553 
554 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
555 
556 } // End namespace Foam
557 
558 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
559 
560 #include "phaseSystemI.H"
561 
562 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
563 
564 #ifdef NoRepository
565  #include "phaseSystemTemplates.C"
566 #endif
567 
568 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
569 
570 #endif
571 
572 // ************************************************************************* //
Foam::surfaceFields.
const BlendedInterfacialModel< modelType > & lookupBlendedSubModel(const phasePair &key) const
Return a blended sub model between a phase pair.
static const word propertiesName
Default name of the phase properties dictionary.
Definition: phaseSystem.H:285
Hashing function class, shared by all the derived classes.
Definition: string.H:90
virtual PtrList< surfaceScalarField > AFfs() const =0
Return the implicit force coefficients for the face-based.
void generatePairs(const dictTable &modelDicts)
Generate pairs.
fv::options & fvOptions() const
Access the fvOptions.
Definition: phaseSystemI.H:141
const phaseModelPartialList & movingPhases() const
Return the models for phases that are moving.
Definition: phaseSystemI.H:49
surfaceTensionModelTable surfaceTensionModels_
Surface tension models.
Definition: phaseSystem.H:160
tmp< volScalarField > E(const phasePairKey &key) const
Return the aspect-ratio for a pair.
virtual void correctTurbulence()
Correct the turbulence.
virtual PtrList< surfaceScalarField > phiFfs(const PtrList< surfaceScalarField > &rAUfs)=0
Return the force fluxes for the face-based algorithm.
PtrList< volScalarField > rAUs
Definition: pEqn.H:4
HashPtrTable< fvScalarMatrix > massTransferTable
Definition: phaseSystem.H:79
HashTable< autoPtr< phasePair >, phasePairKey, phasePairKey::hash > phasePairTable
Definition: phaseSystem.H:87
HashPtrTable< fvVectorMatrix > momentumTransferTable
Definition: phaseSystem.H:75
virtual void solve()
Solve for the phase fractions.
phasePairTable phasePairs_
Phase pairs.
Definition: phaseSystem.H:142
const modelType & lookupSubModel(const phasePair &key) const
Return a sub model between a phase pair.
virtual void partialElimination(const PtrList< volScalarField > &rAUs)=0
Solve the drag system for the new velocities and fluxes.
bool foundBlendedSubModel(const phasePair &key) const
Check availability of a blended sub model for a given phase pair.
virtual autoPtr< heatTransferTable > heatTransfer() const =0
Return the heat transfer matrices.
tmp< volScalarField > sigma(const phasePairKey &key) const
Return the surface tension coefficient for a pair.
volScalarField dpdt_
Rate of change of pressure.
Definition: phaseSystem.H:148
TypeName("phaseSystem")
Runtime type information.
virtual autoPtr< momentumTransferTable > momentumTransferf()=0
Return the momentum transfer matrices for the face-based.
HashPtrTable< fvScalarMatrix > heatTransferTable
Definition: phaseSystem.H:77
phaseModelPartialList stationaryPhaseModels_
Stationary phase models.
Definition: phaseSystem.H:133
virtual PtrList< surfaceScalarField > ddtCorrByAs(const PtrList< volScalarField > &rAUs, const bool includeVirtualMass=false) const =0
Return the flux corrections for the cell-based algorithm.
A HashTable specialization for hashing pointers.
Definition: HashPtrTable.H:50
void createSubModels(const dictTable &modelDicts, HashTable< autoPtr< modelType >, phasePairKey, phasePairKey::hash > &models)
Generate pairs and sub-model tables.
virtual tmp< volScalarField > dmdt(const phasePairKey &key) const
Return the mass transfer rate for a pair.
HashTable< autoPtr< surfaceTensionModel >, phasePairKey, phasePairKey::hash > surfaceTensionModelTable
Definition: phaseSystem.H:109
Finite-volume options.
Definition: fvOptions.H:52
virtual autoPtr< momentumTransferTable > momentumTransfer()=0
Return the momentum transfer matrices for the cell-based.
phaseModelList phaseModels_
Phase models.
Definition: phaseSystem.H:127
phaseModelPartialList multiComponentPhaseModels_
Multi-component phase models.
Definition: phaseSystem.H:139
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
const phaseModelPartialList & stationaryPhases() const
Return the models for phases that are stationary.
Definition: phaseSystemI.H:63
const phasePairTable & phasePairs() const
Return the phase pairs.
Definition: phaseSystemI.H:105
const phaseModelList & phases() const
Return the phase models.
Definition: phaseSystemI.H:35
virtual PtrList< volVectorField > KdUByAs(const PtrList< volScalarField > &rAUs) const =0
Return the explicit part of the drag force.
virtual void correctKinematics()
Correct the kinematics.
Dimension set for the base types.
Definition: dimensionSet.H:120
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
phaseSystem(const fvMesh &mesh)
Construct from fvMesh.
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:67
virtual void partialEliminationf(const PtrList< surfaceScalarField > &rAUfs)=0
Solve the drag system for the new fluxes.
A class for handling words, derived from string.
Definition: word.H:59
tmp< volVectorField > U() const
Return the mixture velocity.
const IOMRFZoneList & MRF() const
Return MRF zones.
Definition: phaseSystemI.H:135
const phaseModelPartialList & anisothermalPhases() const
Return the models for phases that have variable temperature.
Definition: phaseSystemI.H:77
virtual void correctThermo()
Correct the thermodynamics.
PtrListDictionary< phaseModel > phaseModelList
Definition: phaseSystem.H:81
UPtrList< phaseModel > phaseModelPartialList
Definition: phaseSystem.H:83
const fvMesh & mesh_
Reference to the mesh.
Definition: phaseSystem.H:124
virtual PtrList< surfaceScalarField > phiFs(const PtrList< volScalarField > &rAUs)=0
Return the force fluxes for the cell-based algorithm.
bool foundSubModel(const phasePair &key) const
Check availability of a sub model for a given phase pair.
const word & name() const
Name function is needed to disambiguate those inherited.
An STL-conforming hash table.
Definition: HashTable.H:61
PtrList< surfaceScalarField > rAUfs
Definition: pEqn.H:30
phaseModelPartialList anisothermalPhaseModels_
Anisothermal phase models.
Definition: phaseSystem.H:136
phaseModelPartialList movingPhaseModels_
Moving phase models.
Definition: phaseSystem.H:130
virtual void correct()
Correct the fluid properties other than those listed below.
tmp< volScalarField > rho() const
Return the mixture density.
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
aspectRatioModelTable aspectRatioModels_
Aspect ratio models.
Definition: phaseSystem.H:163
IOMRFZoneList MRF_
Optional MRF zones.
Definition: phaseSystem.H:151
volScalarField sf(fieldObject, mesh)
const volScalarField & dpdt() const
Return the rate of change of the pressure.
Definition: phaseSystemI.H:123
tmp< volScalarField > byDt(const volScalarField &vf)
virtual void correctEnergyTransport()
Correct the energy transport e.g. alphat.
void fillFields(const word &name, const dimensionSet &dims, PtrList< GeometricField< Type, PatchField, GeoMesh >> &fieldList) const
Fill up gaps in a phase-indexed list of fields with zeros.
void generatePairsAndSubModels(const word &modelName, HashTable< autoPtr< modelType >, phasePairKey, phasePairKey::hash > &models)
Generate pairs and sub-model tables.
Forward declarations of fvMatrix specializations.
tmp< surfaceScalarField > calcPhi(const phaseModelList &phaseModels) const
Calculate and return the mixture flux.
const surfaceScalarField & phi() const
Return the mixture flux.
Definition: phaseSystemI.H:111
virtual ~phaseSystem()
Destructor.
void addField(const phaseModel &phase, const word &fieldName, tmp< GeoField > field, PtrList< GeoField > &fieldList) const
Add the field to a phase-indexed list, with the given name,.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
virtual bool read()
Read base phaseProperties dictionary.
blendingMethodTable blendingMethods_
Blending methods.
Definition: phaseSystem.H:154
HashTable< autoPtr< blendingMethod >, word, word::hash > blendingMethodTable
Definition: phaseSystem.H:100
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:52
virtual autoPtr< massTransferTable > massTransfer() const =0
Return the mass transfer matrices.
HashTable< autoPtr< aspectRatioModel >, phasePairKey, phasePairKey::hash > aspectRatioModelTable
Definition: phaseSystem.H:118
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
surfaceScalarField phi_
Total volumetric flux.
Definition: phaseSystem.H:145
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
virtual PtrList< surfaceScalarField > phiKdPhifs(const PtrList< surfaceScalarField > &rAUfs) const =0
Return the force fluxes for the face-based algorithm.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual PtrList< surfaceScalarField > phiKdPhis(const PtrList< volScalarField > &rAUs) const =0
Return the force fluxes for the cell-based algorithm.
const phaseModelPartialList & multiComponentPhases() const
Return the models for phases that have multiple species.
Definition: phaseSystemI.H:91
const fvMesh & mesh() const
Return the mesh.
Definition: phaseSystemI.H:28
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
List of MRF zones with IO functionality. MRF zones are specified by a list of dictionary entries...
Definition: IOMRFZoneList.H:66
HashTable< dictionary, phasePairKey, phasePairKey::hash > dictTable
Definition: phaseSystem.H:96
virtual const HashPtrTable< surfaceScalarField > & DByAfs() const =0
Return the phase diffusivities divided by the momentum.
Namespace for OpenFOAM.