44 #include "phaseModel.H" 45 #include "phasePair.H" 46 #include "orderedPhasePair.H" 63 template<
class modelType>
class BlendedInterfacialModel;
64 class aspectRatioModel;
66 namespace reactingMultiphaseEuler
68 class surfaceTensionModel;
187 template<
class modelType>
200 template<
class modelType>
203 const word& modelName,
213 template<
class modelType>
216 const word& modelName,
223 const bool correctFixedFluxBCs =
true 227 template<
class modelType>
230 const word& modelName,
237 const bool correctFixedFluxBCs =
true 242 template<
class GeoField>
246 const word& fieldName,
253 template<
class GeoField>
257 const word& fieldName,
258 const GeoField&
field,
264 template<
class GeoField>
268 const word& fieldName,
275 template<
class GeoField>
279 const word& fieldName,
280 const GeoField&
field,
366 template<
class modelType>
370 template<
class modelType>
374 template<
class modelType>
382 template<
class modelType>
390 template<
class modelType>
394 template<
class modelType>
405 template<
class>
class PatchField,
419 template<
class>
class PatchField,
512 const bool includeVirtualMass =
false 529 virtual void solve();
virtual void correctKinematics()
Correct the kinematics.
const BlendedInterfacialModel< modelType > & lookupBlendedSubModel(const phasePair &key) const
Return a blended sub model between a phase pair.
virtual PtrList< surfaceScalarField > AFfs() const =0
Return the implicit force coefficients for the face-based.
virtual void correctEnergyTransport()
Correct the energy transport e.g. alphat.
void generatePairs(const dictTable &modelDicts)
Generate pairs.
fv::options & fvOptions() const
Access the fvOptions.
const phaseModelPartialList & movingPhases() const
Return the models for phases that are moving.
surfaceTensionModelTable surfaceTensionModels_
Surface tension models.
virtual PtrList< surfaceScalarField > phiFfs(const PtrList< surfaceScalarField > &rAUfs)=0
Return the force fluxes for the face-based algorithm.
HashPtrTable< fvScalarMatrix > massTransferTable
HashTable< autoPtr< phasePair >, phasePairKey, phasePairKey::hash > phasePairTable
virtual tmp< volScalarField > dmdt(const phasePairKey &key) const
Return the mass transfer rate for a pair.
HashPtrTable< fvVectorMatrix > momentumTransferTable
phasePairTable phasePairs_
Phase pairs.
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.
volScalarField dpdt_
Rate of change of pressure.
TypeName("phaseSystem")
Runtime type information.
virtual autoPtr< momentumTransferTable > momentumTransferf()=0
Return the momentum transfer matrices for the face-based.
HashPtrTable< fvScalarMatrix > heatTransferTable
phaseModelPartialList stationaryPhaseModels_
Stationary phase models.
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 of pointers to objects of type <T>, with deallocation management of the pointers...
void createSubModels(const dictTable &modelDicts, HashTable< autoPtr< modelType >, phasePairKey, phasePairKey::hash > &models)
Generate pairs and sub-model tables.
virtual autoPtr< momentumTransferTable > momentumTransfer()=0
Return the momentum transfer matrices for the cell-based.
phaseModelList phaseModels_
Phase models.
Hashing functor for phasePairKey.
virtual void correctTurbulence()
Correct the turbulence.
phaseModelPartialList multiComponentPhaseModels_
Multi-component phase models.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
const phaseModelPartialList & stationaryPhases() const
Return the models for phases that are stationary.
tmp< volScalarField > rho() const
Return the mixture density.
virtual void correctThermo()
Correct the thermodynamics.
const phasePairTable & phasePairs() const
Return the phase pairs.
tmp< volScalarField > sigma(const phasePairKey &key) const
Return the surface tension coefficient for a pair.
const phaseModelList & phases() const
Return the phase models.
virtual PtrList< volVectorField > KdUByAs(const PtrList< volScalarField > &rAUs) const =0
Return the explicit part of the drag force.
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
An ordered pair of two objects of type <T> with first() and second() elements.
phaseSystem(const fvMesh &mesh)
Construct from fvMesh.
virtual bool read()
Read base phaseProperties dictionary.
Class to represent a system of phases and model interfacial transfers between them.
virtual void partialEliminationf(const PtrList< surfaceScalarField > &rAUfs)=0
Solve the drag system for the new fluxes.
A class for handling words, derived from Foam::string.
const IOMRFZoneList & MRF() const
Return MRF zones.
PtrList< surfaceScalarField > rAUfs
const phaseModelPartialList & anisothermalPhases() const
Return the models for phases that have variable temperature.
PtrListDictionary< phaseModel > phaseModelList
tmp< volScalarField > E(const phasePairKey &key) const
Return the aspect-ratio for a pair.
tmp< volScalarField > byDt(const volScalarField &vf)
UPtrList< phaseModel > phaseModelPartialList
const fvMesh & mesh_
Reference to the mesh.
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 from regIOobject and dictionary.
A HashTable similar to std::unordered_map.
phaseModelPartialList anisothermalPhaseModels_
Anisothermal phase models.
virtual void solve()
Solve for the phase fractions.
const modelType & lookupSubModel(const phasePair &key) const
Return a sub model between a phase pair.
An ordered or unorder pair of phase names. Typically specified as follows.
phaseModelPartialList movingPhaseModels_
Moving phase models.
aspectRatioModelTable aspectRatioModels_
Aspect ratio models.
IOMRFZoneList MRF_
Optional MRF zones.
static const word propertiesName
Default name of the phase properties dictionary.
const volScalarField & dpdt() const
Return the rate of change of the pressure.
HashTable< autoPtr< reactingMultiphaseEuler::surfaceTensionModel >, phasePairKey, phasePairKey::hash > surfaceTensionModelTable
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.
PtrList< volScalarField > rAUs
const surfaceScalarField & phi() const
Return the mixture flux.
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,.
auto key(const Type &t) -> typename std::enable_if< std::is_enum< Type >::value, typename std::underlying_type< Type >::type >::type
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
virtual void correct()
Correct the fluid properties other than those listed below.
blendingMethodTable blendingMethods_
Blending methods.
Mesh data needed to do the Finite Volume discretisation.
HashTable< autoPtr< blendingMethod > > blendingMethodTable
tmp< volVectorField > U() const
Return the mixture velocity.
tmp< surfaceScalarField > calcPhi(const phaseModelList &phaseModels) const
Calculate and return the mixture flux.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
virtual autoPtr< massTransferTable > massTransfer() const =0
Return the mass transfer matrices.
virtual ~phaseSystem()
Destructor.
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
HashTable< autoPtr< aspectRatioModel >, phasePairKey, phasePairKey::hash > aspectRatioModelTable
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
surfaceScalarField phi_
Total volumetric flux.
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
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.
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.
const fvMesh & mesh() const
Return the mesh.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
List of MRF zones with IO functionality. MRF zones are specified by a list of dictionary entries...
HashTable< dictionary, phasePairKey, phasePairKey::hash > dictTable
virtual const HashPtrTable< surfaceScalarField > & DByAfs() const =0
Return the phase diffusivities divided by the momentum.