28 #include "BlendedInterfacialModel.H" 29 #include "dragModel.H" 30 #include "virtualMassModel.H" 31 #include "liftModel.H" 32 #include "wallLubricationModel.H" 33 #include "turbulentDispersionModel.H" 50 template<
class BasePhaseSystem>
54 const phasePairKey& key
57 if (dragModels_.found(key))
59 return dragModels_[key]->K();
65 dragModel::typeName +
":K",
73 template<
class BasePhaseSystem>
77 const phasePairKey& key
80 if (dragModels_.found(key))
82 return dragModels_[key]->Kf();
88 dragModel::typeName +
":K",
96 template<
class BasePhaseSystem>
100 const phasePairKey& key
103 if (virtualMassModels_.found(key))
105 return virtualMassModels_[key]->K();
111 virtualMassModel::typeName +
":K",
119 template<
class BasePhaseSystem>
125 phaseSystem::phasePairTable,
130 const phasePair& pair(phasePairIter());
143 if (!pair.phase1().stationary())
148 eqn += dmdt21*pair.phase2().U() -
fvm::Sp(dmdt21, eqn.psi());
151 if (!pair.phase2().stationary())
156 eqn -= dmdt12*pair.phase1().U() -
fvm::Sp(dmdt12, eqn.psi());
164 template<
class BasePhaseSystem>
171 BasePhaseSystem(mesh)
173 this->generatePairsAndSubModels
179 this->generatePairsAndSubModels
185 this->generatePairsAndSubModels
191 this->generatePairsAndSubModels
194 wallLubricationModels_
197 this->generatePairsAndSubModels
199 "turbulentDispersion",
200 turbulentDispersionModels_
210 const phasePair& pair(this->phasePairs_[dragModelIter.key()]);
228 dragModelIter()->Kf()
235 virtualMassModelTable,
240 const phasePair& pair(this->phasePairs_[virtualMassModelIter.key()]);
248 virtualMassModelIter()->K()
258 virtualMassModelIter()->Kf()
267 template<
class BasePhaseSystem>
275 template<
class BasePhaseSystem>
280 autoPtr<phaseSystem::momentumTransferTable> eqnsPtr
287 forAll(this->movingPhases(), movingPhasei)
289 const phaseModel& phase = this->movingPhases()[movingPhasei];
306 *Kds_[dragModelIter.key()] = dragModelIter()->K();
307 *Kdfs_[dragModelIter.key()] = dragModelIter()->Kf();
314 const phasePair& pair(this->phasePairs_[KdIter.key()]);
318 if (!iter().stationary())
330 virtualMassModelTable,
335 *Vms_[virtualMassModelIter.key()] = virtualMassModelIter()->K();
336 *Vmfs_[virtualMassModelIter.key()] = virtualMassModelIter()->Kf();
343 const phasePair& pair(this->phasePairs_[VmIter.key()]);
347 const phaseModel& phase = iter();
348 const phaseModel& otherPhase = iter.otherPhase();
350 if (!phase.stationary())
365 + this->MRF_.DDt(Vm, U - otherPhase.U());
371 addMassTransferMomentumTransfer(eqns);
377 template<
class BasePhaseSystem>
382 autoPtr<phaseSystem::momentumTransferTable> eqnsPtr
389 forAll(this->movingPhases(), movingPhasei)
391 const phaseModel& phase = this->movingPhases()[movingPhasei];
401 PtrList<fvVectorMatrix> UgradUs(this->phaseModels_.size());
404 const phaseModel& phase = this->phaseModels_[
phasei];
406 if (!phase.stationary())
427 const phasePair& pair(this->phasePairs_[VmIter.key()]);
431 const phaseModel& phase = iter();
432 const phaseModel& otherPhase = iter.otherPhase();
434 if (!phase.stationary())
436 *eqns[phase.name()] -=
439 UgradUs[phase.index()]
440 - (UgradUs[otherPhase.index()] & otherPhase.U())
447 addMassTransferMomentumTransfer(eqns);
453 template<
class BasePhaseSystem>
457 PtrList<surfaceScalarField>
AFfs(this->phaseModels_.size());
463 const phasePair& pair(this->phasePairs_[KdfIter.key()]);
467 this->addField(iter(),
"AFf", Kf,
AFfs);
475 const phasePair& pair(this->phasePairs_[VmfIter.key()]);
483 if (this->fillFields_)
492 template<
class BasePhaseSystem>
496 const PtrList<volScalarField>& rAUs
499 PtrList<surfaceScalarField>
phiFs(this->phaseModels_.size());
510 const phasePair& pair(this->phasePairs_[liftModelIter.key()]);
531 wallLubricationModelTable,
532 wallLubricationModels_,
533 wallLubricationModelIter
538 pair(this->phasePairs_[wallLubricationModelIter.key()]);
560 const phaseModel& phase = this->phaseModels_[
phasei];
572 this->addField(phase,
"phiF", pPrimeByAf*snGradAlpha1,
phiFs);
574 const bool implicitPhasePressure =
575 this->mesh_.solverDict(phase.volScalarField::name()).
576 template lookupOrDefault<Switch>
578 "implicitPhasePressure",
582 if (implicitPhasePressure)
584 this->addField(phase,
"DByAf", pPrimeByAf, DByAfs_);
591 turbulentDispersionModelTable,
592 turbulentDispersionModels_,
593 turbulentDispersionModelIter
597 pair(this->phasePairs_[turbulentDispersionModelIter.key()]);
618 if (DByAfs_.found(pair.phase1().name()))
620 this->addField(pair.phase1(),
"DByAf", DByA1f, DByAfs_);
624 if (this->fillFields_)
633 template<
class BasePhaseSystem>
637 const PtrList<surfaceScalarField>& rAUfs
640 PtrList<surfaceScalarField>
phiFfs(this->phaseModels_.size());
646 const phasePair& pair(this->phasePairs_[VmfIter.key()]);
654 - rAUfs[iter().index()]*
Vmf 657 + iter.otherPhase().DUDtf()
673 const phasePair& pair(this->phasePairs_[liftModelIter.key()]);
679 rAUfs[pair.phase1().index()]*
Ff,
686 - rAUfs[pair.phase2().index()]*
Ff,
694 wallLubricationModelTable,
695 wallLubricationModels_,
696 wallLubricationModelIter
701 pair(this->phasePairs_[wallLubricationModelIter.key()]);
707 rAUfs[pair.phase1().index()]*
Ff,
714 - rAUfs[pair.phase2().index()]*
Ff,
723 const phaseModel& phase = this->phaseModels_[
phasei];
735 this->addField(phase,
"phiFf", pPrimeByAf*snGradAlpha1,
phiFfs);
737 const bool implicitPhasePressure =
738 this->mesh_.solverDict(phase.volScalarField::name()).
739 template lookupOrDefault<Switch>
741 "implicitPhasePressure",
745 if (implicitPhasePressure)
747 this->addField(phase,
"DByAf", pPrimeByAf, DByAfs_);
754 turbulentDispersionModelTable,
755 turbulentDispersionModels_,
756 turbulentDispersionModelIter
760 pair(this->phasePairs_[turbulentDispersionModelIter.key()]);
781 if (DByAfs_.found(pair.phase1().name()))
783 this->addField(pair.phase1(),
"DByAf", DByAf1, DByAfs_);
787 if (this->fillFields_)
796 template<
class BasePhaseSystem>
800 const PtrList<volScalarField>& rAUs
803 PtrList<surfaceScalarField> phiKdPhis(this->phaseModels_.size());
809 const phasePair& pair(this->phasePairs_[KdIter.key()]);
824 if (this->fillFields_)
838 template<
class BasePhaseSystem>
842 const PtrList<surfaceScalarField>& rAUfs
845 PtrList<surfaceScalarField> phiKdPhifs(this->phaseModels_.size());
851 const phasePair& pair(this->phasePairs_[KdfIter.key()]);
859 - rAUfs[iter().index()]*Kf
866 if (this->fillFields_)
880 template<
class BasePhaseSystem>
884 const PtrList<volScalarField>& rAUs
887 PtrList<volVectorField> KdUByAs(this->phaseModels_.size());
893 const phasePair& pair(this->phasePairs_[KdIter.key()]);
901 - rAUs[iter().index()]*K*iter.otherPhase().U(),
907 if (this->fillFields_)
916 template<
class BasePhaseSystem>
920 const PtrList<volScalarField>& rAUs,
921 const bool includeVirtualMass
924 PtrList<surfaceScalarField> ddtCorrByAs(this->phaseModels_.size());
927 PtrList<surfaceScalarField> phiCorrs(this->phaseModels_.size());
930 const phaseModel& phase = this->phaseModels_[
phasei];
943 const phaseModel& phase = this->phaseModels_[
phasei];
952 tmp<surfaceScalarField> phiCorrCoeff =
pos0(alphafBar - 0.99);
954 surfaceScalarField::Boundary& phiCorrCoeffBf =
955 phiCorrCoeff.ref().boundaryFieldRef();
962 !this->mesh_.boundary()[
patchi].coupled()
963 || isA<cyclicAMIFvPatch>(this->mesh_.boundary()[
patchi])
966 phiCorrCoeffBf[
patchi] = 0;
983 if (includeVirtualMass)
988 const phasePair& pair(this->phasePairs_[VmIter.key()]);
992 const phaseModel& phase = iter();
993 const phaseModel& otherPhase = iter.otherPhase();
1001 phiCorrs[phase.index()]
1002 + this->
MRF().absolute(otherPhase.phi())
1004 - phiCorrs[otherPhase.index()]
1016 template<
class BasePhaseSystem>
1019 const PtrList<volScalarField>& rAUs
1022 Info<<
"Inverting drag systems: ";
1027 PtrList<PtrList<volScalarField>> KdByAs(phases.size());
1028 PtrList<PtrList<surfaceScalarField>> phiKds(phases.size());
1035 new PtrList<volScalarField>(phases.size())
1041 new PtrList<surfaceScalarField>(phases.size())
1048 const phasePair& pair(this->phasePairs_[KdIter.key()]);
1050 const label phase1i = pair.phase1().index();
1051 const label phase2i = pair.phase2().index();
1094 for (
label i = 0; i < phases.size(); ++ i)
1096 for (
label j = i + 1; j < phases.size(); ++ j)
1098 KdByAs[i][j] /= KdByAs[i][i];
1099 phiKds[i][j] /= phiKds[i][i];
1100 for (
label k = i + 1;
k < phases.size(); ++
k)
1102 KdByAs[j][
k] -= KdByAs[j][i]*KdByAs[i][
k];
1103 phiKds[j][
k] -= phiKds[j][i]*phiKds[i][
k];
1110 for (
label i = 1; i < phases.size(); ++ i)
1112 detKdByAs *= KdByAs[i][i];
1113 detPhiKdfs *= phiKds[i][i];
1115 Info<<
"Min cell/face det = " <<
gMin(detKdByAs.primitiveField())
1116 <<
"/" <<
gMin(detPhiKdfs.primitiveField()) <<
endl;
1120 for (
label i = 1; i < phases.size(); ++ i)
1122 if (!phases[i].stationary())
1124 for (
label j = 0; j < i; j ++)
1126 phases[i].URef() -= KdByAs[i][j]*phases[j].U();
1127 phases[i].phiRef() -= phiKds[i][j]*phases[j].phi();
1131 for (
label i = phases.size() - 1; i >= 0; i --)
1133 if (!phases[i].stationary())
1135 for (
label j = phases.size() - 1; j > i; j --)
1137 phases[i].URef() -= KdByAs[i][j]*phases[j].U();
1138 phases[i].phiRef() -= phiKds[i][j]*phases[j].phi();
1140 phases[i].URef() /= KdByAs[i][i];
1141 phases[i].phiRef() /= phiKds[i][i];
1147 template<
class BasePhaseSystem>
1150 const PtrList<surfaceScalarField>& rAUfs
1153 Info<<
"Inverting drag system: ";
1158 PtrList<PtrList<surfaceScalarField>> phiKdfs(phases.size());
1165 new PtrList<surfaceScalarField>(phases.size())
1172 const phasePair& pair(this->phasePairs_[KdfIter.key()]);
1174 const label phase1i = pair.phase1().index();
1175 const label phase2i = pair.phase2().index();
1201 for (
label i = 0; i < phases.size(); ++ i)
1203 for (
label j = i + 1; j < phases.size(); ++ j)
1205 phiKdfs[i][j] /= phiKdfs[i][i];
1206 for (
label k = i + 1;
k < phases.size(); ++
k)
1208 phiKdfs[j][
k] -= phiKdfs[j][i]*phiKdfs[i][
k];
1214 for (
label i = 1; i < phases.size(); ++ i)
1216 detPhiKdfs *= phiKdfs[i][i];
1218 Info<<
"Min face det = " <<
gMin(detPhiKdfs.primitiveField()) <<
endl;
1222 for (
label i = 1; i < phases.size(); ++ i)
1224 if (!phases[i].stationary())
1226 for (
label j = 0; j < i; j ++)
1228 phases[i].phiRef() -= phiKdfs[i][j]*phases[j].phi();
1232 for (
label i = phases.size() - 1; i >= 0; i --)
1234 if (!phases[i].stationary())
1236 for (
label j = phases.size() - 1; j > i; j --)
1238 phases[i].phiRef() -= phiKdfs[i][j]*phases[j].phi();
1240 phases[i].phiRef() /= phiKdfs[i][i];
1246 template<
class BasePhaseSystem>
1254 template<
class BasePhaseSystem>
Class which models interfacial momenum transfer between a number of phases. Drag, virtual mass...
surfaceScalarField Vmf("Vmf", fluid.Vmf())
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
virtual PtrList< surfaceScalarField > phiFs(const PtrList< volScalarField > &rAUs)
Return the explicit force fluxes for the cell-based algorithm, that.
Type gMin(const FieldField< Field, Type > &f)
HashPtrTable< fvVectorMatrix > momentumTransferTable
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransferf()
As momentumTransfer, but for the face-based algorithm.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
virtual void partialEliminationf(const PtrList< surfaceScalarField > &rAUfs)
As partialElimination, but for the face-based algorithm. Only solves.
virtual bool read()
Read base phaseProperties dictionary.
Ostream & endl(Ostream &os)
Add newline and flush stream.
MomentumTransferPhaseSystem(const fvMesh &)
Construct from fvMesh.
Calculate the snGrad of the given volField.
label k
Boltzmann constant.
A HashTable specialization for hashing pointers.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
PtrList< surfaceScalarField > AFfs(fluid.AFfs())
GeometricField< vector, fvPatchField, volMesh > volVectorField
CGAL::Exact_predicates_exact_constructions_kernel K
dimensionedScalar posPart(const dimensionedScalar &ds)
Area-weighted average a surfaceField creating a volField.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Calculate the first temporal derivative.
PtrList< surfaceScalarField > phiFs(fluid.phiFs(rAUs))
virtual const HashPtrTable< surfaceScalarField > & DByAfs() const
Return the phase diffusivities divided by the momentum coefficients.
Calculate the face-flux of the given field.
static word groupName(Name name, const word &group)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Calculate the matrix for the first temporal derivative.
PtrListDictionary< phaseModel > phaseModelList
virtual PtrList< volVectorField > KdUByAs(const PtrList< volScalarField > &rAUs) const
Return the explicit part of the drag force for the cell-based.
volVectorField F(fluid.F())
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
virtual PtrList< surfaceScalarField > phiKdPhis(const PtrList< volScalarField > &rAUs) const
Return the explicit drag force fluxes for the cell-based algorithm.
Calculate the divergence of the given field.
dimensionedScalar pos0(const dimensionedScalar &ds)
const dimensionSet dimForce
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
virtual PtrList< surfaceScalarField > ddtCorrByAs(const PtrList< volScalarField > &rAUs, const bool includeVirtualMass=false) const
Return the flux corrections for the cell-based algorithm. These.
tmp< volScalarField > byDt(const volScalarField &vf)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
const dimensionSet dimDensity
Calculate the matrix for the divergence of the given field and flux.
virtual ~MomentumTransferPhaseSystem()
Destructor.
virtual PtrList< surfaceScalarField > phiFfs(const PtrList< surfaceScalarField > &rAUfs)
As phiFs, but for the face-based algorithm.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
surfaceScalarField snGradAlpha1(fvc::snGrad(alpha1) *mesh.magSf())
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
fvMatrix< vector > fvVectorMatrix
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
tmp< GeometricField< Type, fvPatchField, volMesh > > average(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Area-weighted average a surfaceField creating a volField.
fluid fillFields("rAU", dimTime/dimDensity, rAUs)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
A class for managing temporary objects.
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
surfaceScalarField Ff(fluid.Ff())
virtual PtrList< surfaceScalarField > AFfs() const
Return implicit force coefficients on the faces, for the face-based.
PtrList< surfaceScalarField > phiFfs(fluid.phiFfs(rAUfs))
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Calculate the matrix for implicit and explicit sources.
dimensionedScalar negPart(const dimensionedScalar &ds)
virtual void partialElimination(const PtrList< volScalarField > &rAUs)
Solve the drag system for the velocities and fluxes.
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransfer()
Return the momentum transfer matrices for the cell-based algorithm.
virtual PtrList< surfaceScalarField > phiKdPhifs(const PtrList< surfaceScalarField > &rAUfs) const
As phiKdPhis, but for the face-based algorithm.
const dimensionSet dimVelocity