Go to the documentation of this file.
30 #include "BlendedInterfacialModel.H"
31 #include "dragModel.H"
32 #include "virtualMassModel.H"
33 #include "liftModel.H"
34 #include "wallLubricationModel.H"
35 #include "turbulentDispersionModel.H"
52 template<
class BasePhaseSystem>
56 const phasePairKey& key
59 if (dragModels_.found(key))
61 return dragModels_[key]->K();
67 dragModel::typeName +
":K",
75 template<
class BasePhaseSystem>
79 const phasePairKey& key
82 if (dragModels_.found(key))
84 return dragModels_[key]->Kf();
90 dragModel::typeName +
":K",
98 template<
class BasePhaseSystem>
102 const phasePairKey& key
105 if (virtualMassModels_.found(key))
107 return virtualMassModels_[key]->K();
113 virtualMassModel::typeName +
":K",
121 template<
class BasePhaseSystem>
127 phaseSystem::phasePairTable,
132 const phasePair& pair(phasePairIter());
145 if (!pair.phase1().stationary())
150 eqn += dmdt21*pair.phase2().U() -
fvm::Sp(dmdt21, eqn.psi());
153 if (!pair.phase2().stationary())
158 eqn -= dmdt12*pair.phase1().U() -
fvm::Sp(dmdt12, eqn.psi());
166 template<
class BasePhaseSystem>
173 BasePhaseSystem(
mesh)
175 this->generatePairsAndSubModels
181 this->generatePairsAndSubModels
187 this->generatePairsAndSubModels
193 this->generatePairsAndSubModels
196 wallLubricationModels_
199 this->generatePairsAndSubModels
201 "turbulentDispersion",
202 turbulentDispersionModels_
212 const phasePair& pair(this->phasePairs_[dragModelIter.key()]);
219 IOobject::groupName(
"Kd", pair.name()),
229 IOobject::groupName(
"Kdf", pair.name()),
230 dragModelIter()->Kf()
242 const phasePair& pair(this->phasePairs_[virtualMassModelIter.key()]);
249 IOobject::groupName(
"Vm", pair.name()),
250 virtualMassModelIter()->K()
259 IOobject::groupName(
"Vmf", pair.name()),
260 virtualMassModelIter()->Kf()
269 template<
class BasePhaseSystem>
277 template<
class BasePhaseSystem>
289 forAll(this->movingPhases(), movingPhasei)
308 *Kds_[dragModelIter.key()] = dragModelIter()->K();
309 *Kdfs_[dragModelIter.key()] = dragModelIter()->Kf();
316 const phasePair& pair(this->phasePairs_[KdIter.key()]);
320 if (!iter().stationary())
337 *Vms_[virtualMassModelIter.key()] = virtualMassModelIter()->K();
338 *Vmfs_[virtualMassModelIter.key()] = virtualMassModelIter()->Kf();
345 const phasePair& pair(this->phasePairs_[VmIter.key()]);
352 if (!
phase.stationary())
367 + this->MRF_.DDt(Vm,
U - otherPhase.
U());
373 addMassTransferMomentumTransfer(eqns);
379 template<
class BasePhaseSystem>
391 forAll(this->movingPhases(), movingPhasei)
408 if (!
phase.stationary())
429 const phasePair& pair(this->phasePairs_[VmIter.key()]);
436 if (!
phase.stationary())
441 UgradUs[
phase.index()]
442 - (UgradUs[otherPhase.
index()] & otherPhase.
U())
449 addMassTransferMomentumTransfer(eqns);
455 template<
class BasePhaseSystem>
465 const phasePair& pair(this->phasePairs_[KdfIter.key()]);
469 this->addField(iter(),
"AFf", Kf,
AFfs);
477 const phasePair& pair(this->phasePairs_[VmfIter.key()]);
485 if (this->fillFields_)
494 template<
class BasePhaseSystem>
512 const phasePair& pair(this->phasePairs_[liftModelIter.key()]);
534 wallLubricationModels_,
535 wallLubricationModelIter
540 pair(this->phasePairs_[wallLubricationModelIter.key()]);
576 const bool implicitPhasePressure =
577 this->mesh_.solverDict(
phase.volScalarField::name()).
578 template lookupOrDefault<Switch>
580 "implicitPhasePressure",
584 if (implicitPhasePressure)
586 this->addField(
phase,
"DByAf", pPrimeByAf, DByAfs_);
594 turbulentDispersionModels_,
595 turbulentDispersionModelIter
599 pair(this->phasePairs_[turbulentDispersionModelIter.key()]);
620 if (DByAfs_.found(pair.phase1().name()))
622 this->addField(pair.phase1(),
"DByAf", DByA1f, DByAfs_);
626 if (this->fillFields_)
635 template<
class BasePhaseSystem>
648 const phasePair& pair(this->phasePairs_[VmfIter.key()]);
659 + iter.otherPhase().DUDtf()
675 const phasePair& pair(this->phasePairs_[liftModelIter.key()]);
681 rAUfs[pair.phase1().index()]*
Ff,
688 -
rAUfs[pair.phase2().index()]*
Ff,
697 wallLubricationModels_,
698 wallLubricationModelIter
703 pair(this->phasePairs_[wallLubricationModelIter.key()]);
709 rAUfs[pair.phase1().index()]*
Ff,
716 -
rAUfs[pair.phase2().index()]*
Ff,
739 const bool implicitPhasePressure =
740 this->mesh_.solverDict(
phase.volScalarField::name()).
741 template lookupOrDefault<Switch>
743 "implicitPhasePressure",
747 if (implicitPhasePressure)
749 this->addField(
phase,
"DByAf", pPrimeByAf, DByAfs_);
757 turbulentDispersionModels_,
758 turbulentDispersionModelIter
762 pair(this->phasePairs_[turbulentDispersionModelIter.key()]);
783 if (DByAfs_.found(pair.phase1().name()))
785 this->addField(pair.phase1(),
"DByAf", DByAf1, DByAfs_);
789 if (this->fillFields_)
798 template<
class BasePhaseSystem>
811 const phasePair& pair(this->phasePairs_[KdIter.key()]);
826 if (this->fillFields_)
840 template<
class BasePhaseSystem>
853 const phasePair& pair(this->phasePairs_[KdfIter.key()]);
861 -
rAUfs[iter().index()]*Kf
868 if (this->fillFields_)
882 template<
class BasePhaseSystem>
895 const phasePair& pair(this->phasePairs_[KdIter.key()]);
903 -
rAUs[iter().index()]*
K*iter.otherPhase().U(),
909 if (this->fillFields_)
918 template<
class BasePhaseSystem>
923 const bool includeVirtualMass
956 surfaceScalarField::Boundary& phiCorrCoeffBf =
957 phiCorrCoeff.
ref().boundaryFieldRef();
959 forAll(this->mesh_.boundary(), patchi)
964 !this->mesh_.boundary()[patchi].coupled()
965 || isA<cyclicAMIFvPatch>(this->mesh_.boundary()[patchi])
968 phiCorrCoeffBf[patchi] = 0;
985 if (includeVirtualMass)
990 const phasePair& pair(this->phasePairs_[VmIter.key()]);
1003 phiCorrs[
phase.index()]
1004 + this->
MRF().absolute(otherPhase.
phi())
1006 - phiCorrs[otherPhase.
index()]
1018 template<
class BasePhaseSystem>
1024 Info<<
"Inverting drag systems: ";
1050 const phasePair& pair(this->phasePairs_[KdIter.key()]);
1053 const label phase2i = pair.phase2().index();
1100 KdByAs[i][j] /= KdByAs[i][i];
1101 phiKds[i][j] /= phiKds[i][i];
1104 KdByAs[j][
k] -= KdByAs[j][i]*KdByAs[i][
k];
1105 phiKds[j][
k] -= phiKds[j][i]*phiKds[i][
k];
1114 detKdByAs *= KdByAs[i][i];
1115 detPhiKdfs *= phiKds[i][i];
1124 if (!
phases[i].stationary())
1126 for (
label j = 0; j < i; j ++)
1135 if (!
phases[i].stationary())
1142 phases[i].URef() /= KdByAs[i][i];
1143 phases[i].phiRef() /= phiKds[i][i];
1149 template<
class BasePhaseSystem>
1155 Info<<
"Inverting drag system: ";
1174 const phasePair& pair(this->phasePairs_[KdfIter.key()]);
1177 const label phase2i = pair.phase2().index();
1207 phiKdfs[i][j] /= phiKdfs[i][i];
1210 phiKdfs[j][
k] -= phiKdfs[j][i]*phiKdfs[i][
k];
1218 detPhiKdfs *= phiKdfs[i][i];
1226 if (!
phases[i].stationary())
1228 for (
label j = 0; j < i; j ++)
1236 if (!
phases[i].stationary())
1242 phases[i].phiRef() /= phiKdfs[i][i];
1248 template<
class BasePhaseSystem>
1256 template<
class BasePhaseSystem>
virtual PtrList< surfaceScalarField > phiFfs(const PtrList< surfaceScalarField > &rAUfs)
As phiFs, but for the face-based algorithm.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
label index() const
Return the index of the phase.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
virtual ~MomentumTransferPhaseSystem()
Destructor.
PtrList< volScalarField > rAUs
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
A class for managing temporary objects.
const dimensionSet dimVelocity
const dimensionSet dimDensity
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
dimensionedScalar posPart(const dimensionedScalar &ds)
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Calculate the snGrad of the given volField.
virtual PtrList< volVectorField > KdUByAs(const PtrList< volScalarField > &rAUs) const
Return the explicit part of the drag force for the cell-based.
Info<< "Creating field kinetic energy K\n"<< endl;volScalarField K("K", 0.5 *magSqr(U));if(U.nOldTimes()){ volVectorField *Uold=&U.oldTime();volScalarField *Kold=&K.oldTime();*Kold==0.5 *magSqr(*Uold);while(Uold->nOldTimes()) { Uold=&Uold-> oldTime()
Calculate the divergence of the given field.
const phaseModel & otherPhase() const
Return the other phase in this two-phase system.
surfaceScalarField Ff(fluid.Ff())
volVectorField F(fluid.F())
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar pos0(const dimensionedScalar &ds)
const dimensionSet dimForce
tmp< GeometricField< Type, faPatchField, areaMesh > > div(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
const GeometricField< Type, fvPatchField, volMesh > & psi() const
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object.
Calculate the matrix for the divergence of the given field and flux.
#define forAll(list, i)
Loop across all elements in list.
CGAL::Exact_predicates_exact_constructions_kernel K
surfaceScalarField Vmf("Vmf", fluid.Vmf())
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
fvMatrix< vector > fvVectorMatrix
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
bool set(const Key &key, T *ptr)
Assign a new entry, overwriting existing entries.
tmp< GeometricField< Type, faPatchField, areaMesh > > average(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
Area-weighted average a edgeField creating a areaField.
PtrList< surfaceScalarField > AFfs(fluid.AFfs())
messageStream Info
Information stream (uses stdout - output is on the master only)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
GeometricField< scalar, fvPatchField, volMesh > volScalarField
virtual void partialElimination(const PtrList< volScalarField > &rAUs)
Solve the drag system for the velocities and fluxes.
virtual PtrList< surfaceScalarField > phiKdPhis(const PtrList< volScalarField > &rAUs) const
Return the explicit drag force fluxes for the cell-based algorithm.
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransferf()
As momentumTransfer, but for the face-based algorithm.
MomentumTransferPhaseSystem(const fvMesh &)
Construct from fvMesh.
PtrList< surfaceScalarField > phiFfs(fluid.phiFfs(rAUfs))
Mesh data needed to do the Finite Volume discretisation.
virtual PtrList< surfaceScalarField > phiKdPhifs(const PtrList< surfaceScalarField > &rAUfs) const
As phiKdPhis, but for the face-based algorithm.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
PtrList< surfaceScalarField > phiFs(fluid.phiFs(rAUs))
surfaceScalarField snGradAlpha1(fvc::snGrad(alpha1) *mesh.magSf())
Calculate the matrix for implicit and explicit sources.
const word & name() const
virtual tmp< surfaceScalarField > phi() const =0
Return the volumetric flux.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
fluid fillFields("rAU", dimTime/dimDensity, rAUs)
virtual bool read()
Read base phaseProperties dictionary.
virtual PtrList< surfaceScalarField > phiFs(const PtrList< volScalarField > &rAUs)
Return the explicit force fluxes for the cell-based algorithm, that.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
dimensionedScalar negPart(const dimensionedScalar &ds)
PtrList< surfaceScalarField > rAUfs
A HashTable of pointers to objects of type <T>.
tmp< GeometricField< Type, faPatchField, areaMesh > > ddt(const dimensioned< Type > dt, const faMesh &mesh)
Calculate the face-flux of the given field.
multiphaseSystem::phaseModelList & phases
label k
Boltzmann constant.
const dimensionedScalar & D
A special matrix type and solver, designed for finite volume solutions of scalar equations....
virtual tmp< volVectorField > U() const =0
Return the velocity.
virtual void partialEliminationf(const PtrList< surfaceScalarField > &rAUfs)
As partialElimination, but for the face-based algorithm. Only solves.
Calculate the first temporal derivative.
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransfer()
Return the momentum transfer matrices for the cell-based algorithm.
virtual PtrList< surfaceScalarField > AFfs() const
Return implicit force coefficients on the faces, for the face-based.
Type gMin(const FieldField< Field, Type > &f)
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
tmp< volScalarField > byDt(const volScalarField &vf)
const phaseModel & phase1() const
Return phase 1.
Area-weighted average a surfaceField creating a volField.
Class which models interfacial momenum transfer between a number of phases. Drag, virtual mass,...
static tmp< GeometricField< Type, faePatchField, edgeMesh > > interpolate(const GeometricField< Type, faPatchField, areaMesh > &tvf, const edgeScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
virtual PtrList< surfaceScalarField > ddtCorrByAs(const PtrList< volScalarField > &rAUs, const bool includeVirtualMass=false) const
Return the flux corrections for the cell-based algorithm. These.
Calulate the matrix for the first temporal derivative.
virtual const HashPtrTable< surfaceScalarField > & DByAfs() const
Return the phase diffusivities divided by the momentum coefficients.
virtual tmp< volVectorField > DUDt() const =0
Return the substantive acceleration.
const dimensionedScalar & rho() const
Return const-access to phase1 density.