multiphaseEulerFoam: Prevent population balances from interfering

Population balance models now own their mass transfer rates, rather than
taking a non-constant reference to rates held by the phase system. This
means that they cannot reset or modify rates that relate to other
population balances.
This commit is contained in:
Will Bainbridge
2022-08-23 08:19:41 +01:00
parent c85a01a6e3
commit b6818dd901
7 changed files with 152 additions and 141 deletions

View File

@ -38,72 +38,9 @@ PopulationBalancePhaseSystem
populationBalances_
(
this->lookup("populationBalances"),
diameterModels::populationBalanceModel::iNew(*this, dmdtfs_)
diameterModels::populationBalanceModel::iNew(*this)
)
{
forAll(populationBalances_, popBali)
{
const diameterModels::populationBalanceModel& popBal =
populationBalances_[popBali];
forAllConstIter
(
HashTable<const diameterModels::velocityGroup*>,
popBal.velocityGroupPtrs(),
iter1
)
{
const diameterModels::velocityGroup& velGrp1 = *iter1();
forAllConstIter
(
HashTable<const diameterModels::velocityGroup*>,
popBal.velocityGroupPtrs(),
iter2
)
{
const diameterModels::velocityGroup& velGrp2 = *iter2();
const phaseInterface interface
(
velGrp1.phase(),
velGrp2.phase()
);
if
(
&velGrp1 != &velGrp2
&&
!dmdtfs_.found(interface)
)
{
this->template validateMassTransfer
<diameterModels::populationBalanceModel>(interface);
dmdtfs_.insert
(
interface,
new volScalarField
(
IOobject
(
IOobject::groupName
(
"populationBalance:dmdtf",
interface.name()
),
this->mesh().time().timeName(),
this->mesh()
),
this->mesh(),
dimensionedScalar(dimDensity/dimTime, 0)
)
);
}
}
}
}
}
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
@ -125,9 +62,12 @@ Foam::PopulationBalancePhaseSystem<BasePhaseSystem>::dmdtf
{
tmp<volScalarField> tDmdtf = BasePhaseSystem::dmdtf(key);
if (dmdtfs_.found(key))
forAll(populationBalances_, popBali)
{
tDmdtf.ref() += *dmdtfs_[key];
if (populationBalances_[popBali].dmdtfs().found(key))
{
tDmdtf.ref() += *populationBalances_[popBali].dmdtfs()[key];
}
}
return tDmdtf;
@ -140,12 +80,20 @@ Foam::PopulationBalancePhaseSystem<BasePhaseSystem>::dmdts() const
{
PtrList<volScalarField> dmdts(BasePhaseSystem::dmdts());
forAllConstIter(phaseSystem::dmdtfTable, dmdtfs_, dmdtfIter)
forAll(populationBalances_, popBali)
{
const phaseInterface interface(*this, dmdtfIter.key());
forAllConstIter
(
phaseSystem::dmdtfTable,
populationBalances_[popBali].dmdtfs(),
dmdtfIter
)
{
const phaseInterface interface(*this, dmdtfIter.key());
addField(interface.phase1(), "dmdt", *dmdtfIter(), dmdts);
addField(interface.phase2(), "dmdt", - *dmdtfIter(), dmdts);
addField(interface.phase1(), "dmdt", *dmdtfIter(), dmdts);
addField(interface.phase2(), "dmdt", - *dmdtfIter(), dmdts);
}
}
return dmdts;
@ -161,7 +109,10 @@ Foam::PopulationBalancePhaseSystem<BasePhaseSystem>::momentumTransfer()
phaseSystem::momentumTransferTable& eqns = eqnsPtr();
this->addDmdtUfs(dmdtfs_, eqns);
forAll(populationBalances_, popBali)
{
this->addDmdtUfs(populationBalances_[popBali].dmdtfs(), eqns);
}
return eqnsPtr;
}
@ -176,7 +127,10 @@ Foam::PopulationBalancePhaseSystem<BasePhaseSystem>::momentumTransferf()
phaseSystem::momentumTransferTable& eqns = eqnsPtr();
this->addDmdtUfs(dmdtfs_, eqns);
forAll(populationBalances_, popBali)
{
this->addDmdtUfs(populationBalances_[popBali].dmdtfs(), eqns);
}
return eqnsPtr;
}
@ -191,7 +145,10 @@ Foam::PopulationBalancePhaseSystem<BasePhaseSystem>::heatTransfer() const
phaseSystem::heatTransferTable& eqns = eqnsPtr();
this->addDmdtHefs(dmdtfs_, eqns);
forAll(populationBalances_, popBali)
{
this->addDmdtHefs(populationBalances_[popBali].dmdtfs(), eqns);
}
return eqnsPtr;
}
@ -206,7 +163,10 @@ Foam::PopulationBalancePhaseSystem<BasePhaseSystem>::specieTransfer() const
phaseSystem::specieTransferTable& eqns = eqnsPtr();
this->addDmdtYfs(dmdtfs_, eqns);
forAll(populationBalances_, popBali)
{
this->addDmdtYfs(populationBalances_[popBali].dmdtfs(), eqns);
}
return eqnsPtr;
}

View File

@ -62,9 +62,6 @@ class PopulationBalancePhaseSystem
//- Population balances
PtrList<diameterModels::populationBalanceModel> populationBalances_;
//- Mass transfer rates
phaseSystem::dmdtfTable dmdtfs_;
public:

View File

@ -181,10 +181,6 @@ protected:
const phaseModelList& phaseModels
) const;
//- Check that mass transfer is supported across the given interface
template<class ModelType>
void validateMassTransfer(const phaseInterface& interface) const;
//- Return the sum of the phase fractions of the moving phases
tmp<volScalarField> sumAlphaMoving() const;
@ -413,6 +409,10 @@ public:
template<class ModelType>
static const dictionary& modelSubDict(const dictionary& dict);
//- Check that mass transfer is supported across the given interface
template<class ModelType>
void validateMassTransfer(const phaseInterface& interface) const;
// Sub-model lookup

View File

@ -28,23 +28,6 @@ License
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
template<class ModelType>
void Foam::phaseSystem::validateMassTransfer
(
const phaseInterface& interface
) const
{
if (interface.phase1().stationary() || interface.phase2().stationary())
{
FatalErrorInFunction
<< "A " << ModelType::typeName << " was specified for pair "
<< interface.name() << ", but one of these phases is stationary. "
<< "Mass transfer is not supported on stationary phases"
<< exit(FatalError);
}
}
namespace Foam
{
@ -455,6 +438,23 @@ const Foam::dictionary& Foam::phaseSystem::modelSubDict
}
template<class ModelType>
void Foam::phaseSystem::validateMassTransfer
(
const phaseInterface& interface
) const
{
if (interface.phase1().stationary() || interface.phase2().stationary())
{
FatalErrorInFunction
<< "A " << ModelType::typeName << " was specified for pair "
<< interface.name() << ", but one of these phases is stationary. "
<< "Mass transfer is not supported on stationary phases"
<< exit(FatalError);
}
}
template<class ModelType>
bool Foam::phaseSystem::foundInterfacialModel
(

View File

@ -192,6 +192,63 @@ void Foam::diameterModels::populationBalanceModel::registerSizeGroups
}
void Foam::diameterModels::populationBalanceModel::initialiseDmdtfs()
{
forAllConstIter
(
HashTable<const diameterModels::velocityGroup*>,
velocityGroupPtrs(),
iter1
)
{
const diameterModels::velocityGroup& velGrp1 = *iter1();
forAllConstIter
(
HashTable<const diameterModels::velocityGroup*>,
velocityGroupPtrs(),
iter2
)
{
const diameterModels::velocityGroup& velGrp2 = *iter2();
const phaseInterface interface(velGrp1.phase(), velGrp2.phase());
if
(
&velGrp1 != &velGrp2
&&
!dmdtfs_.found(interface)
)
{
fluid_.validateMassTransfer
<diameterModels::populationBalanceModel>(interface);
dmdtfs_.insert
(
interface,
new volScalarField
(
IOobject
(
IOobject::groupName
(
"populationBalance:dmdtf",
interface.name()
),
mesh().time().timeName(),
mesh()
),
mesh(),
dimensionedScalar(dimDensity/dimTime, 0)
)
);
}
}
}
}
void Foam::diameterModels::populationBalanceModel::precompute()
{
forAll(coalescenceModels_, model)
@ -260,22 +317,22 @@ void Foam::diameterModels::populationBalanceModel::birthByCoalescence
const phaseInterface interfaceij(fi.phase(), fj.phase());
if (pDmdt_.found(interfaceij))
if (dmdtfs_.found(interfaceij))
{
const scalar dmdtSign =
interfaceij.index(fi.phase()) == 0 ? +1 : -1;
*pDmdt_[interfaceij] += dmdtSign*fj.x()/v*Sui_*fj.phase().rho();
*dmdtfs_[interfaceij] += dmdtSign*fj.x()/v*Sui_*fj.phase().rho();
}
const phaseInterface interfaceik(fi.phase(), fk.phase());
if (pDmdt_.found(interfaceik))
if (dmdtfs_.found(interfaceik))
{
const scalar dmdtSign =
interfaceik.index(fi.phase()) == 0 ? +1 : -1;
*pDmdt_[interfaceik] += dmdtSign*fk.x()/v*Sui_*fk.phase().rho();
*dmdtfs_[interfaceik] += dmdtSign*fk.x()/v*Sui_*fk.phase().rho();
}
sizeGroups_[i].shapeModelPtr()->addCoalescence(Sui_, fj, fk);
@ -321,12 +378,12 @@ void Foam::diameterModels::populationBalanceModel::birthByBreakup
const phaseInterface interface(fi.phase(), fk.phase());
if (pDmdt_.found(interface))
if (dmdtfs_.found(interface))
{
const scalar dmdtSign =
interface.index(fi.phase()) == 0 ? +1 : -1;
*pDmdt_[interface] += dmdtSign*Sui_*fk.phase().rho();
*dmdtfs_[interface] += dmdtSign*Sui_*fk.phase().rho();
}
sizeGroups_[i].shapeModelPtr()->addBreakup(Sui_, fk);
@ -398,12 +455,12 @@ void Foam::diameterModels::populationBalanceModel::birthByBinaryBreakup
const phaseInterface interfaceij(fi.phase(), fj.phase());
if (pDmdt_.found(interfaceij))
if (dmdtfs_.found(interfaceij))
{
const scalar dmdtSign =
interfaceij.index(fi.phase()) == 0 ? +1 : -1;
*pDmdt_[interfaceij] += dmdtSign*Sui_*fj.phase().rho();
*dmdtfs_[interfaceij] += dmdtSign*Sui_*fj.phase().rho();
}
dimensionedScalar Eta;
@ -425,12 +482,12 @@ void Foam::diameterModels::populationBalanceModel::birthByBinaryBreakup
const phaseInterface interfacekj(fk.phase(), fj.phase());
if (pDmdt_.found(interfacekj))
if (dmdtfs_.found(interfacekj))
{
const scalar dmdtSign =
interfacekj.index(fk.phase()) == 0 ? +1 : -1;
*pDmdt_[interfacekj] += dmdtSign*Suk*fj.phase().rho();
*dmdtfs_[interfacekj] += dmdtSign*Suk*fj.phase().rho();
}
sizeGroups_[k].shapeModelPtr()->addBreakup(Suk, fj);
@ -472,12 +529,12 @@ void Foam::diameterModels::populationBalanceModel::drift
const phaseInterface interfacepe(fp.phase(), fe.phase());
if (pDmdt_.found(interfacepe))
if (dmdtfs_.found(interfacepe))
{
const scalar dmdtSign =
interfacepe.index(fp.phase()) == 0 ? +1 : -1;
*pDmdt_[interfacepe] -= dmdtSign*Sue*fp.phase().rho();
*dmdtfs_[interfacepe] -= dmdtSign*Sue*fp.phase().rho();
}
sizeGroups_[i+1].shapeModelPtr()->addDrift(Sue, fp, model);
@ -502,12 +559,12 @@ void Foam::diameterModels::populationBalanceModel::drift
const phaseInterface interfacepw(fp.phase(), fw.phase());
if (pDmdt_.found(interfacepw))
if (dmdtfs_.found(interfacepw))
{
const scalar dmdtSign =
interfacepw.index(fp.phase()) == 0 ? +1 : -1;
*pDmdt_[interfacepw] -= dmdtSign*Suw*fp.phase().rho();
*dmdtfs_[interfacepw] -= dmdtSign*Suw*fp.phase().rho();
}
sizeGroups_[i-1].shapeModelPtr()->addDrift(Suw, fp, model);
@ -547,7 +604,7 @@ void Foam::diameterModels::populationBalanceModel::sources()
Sp_[i] = Zero;
}
forAllIter(phaseSystem::dmdtfTable, pDmdt_, pDmdtIter)
forAllIter(phaseSystem::dmdtfTable, dmdtfs_, pDmdtIter)
{
*pDmdtIter() = Zero;
}
@ -720,8 +777,7 @@ bool Foam::diameterModels::populationBalanceModel::updateSources()
Foam::diameterModels::populationBalanceModel::populationBalanceModel
(
const phaseSystem& fluid,
const word& name,
phaseSystem::dmdtfTable& pDmdt
const word& name
)
:
regIOobject
@ -734,7 +790,7 @@ Foam::diameterModels::populationBalanceModel::populationBalanceModel
)
),
fluid_(fluid),
pDmdt_(pDmdt),
dmdtfs_(),
mesh_(fluid_.mesh()),
name_(name),
dict_
@ -812,6 +868,7 @@ Foam::diameterModels::populationBalanceModel::populationBalanceModel
<< exit(FatalError);
}
this->initialiseDmdtfs();
if (coalescenceModels_.size() != 0)
{

View File

@ -240,7 +240,7 @@ private:
const phaseSystem& fluid_;
//- Interfacial mass transfer rates
phaseSystem::dmdtfTable& pDmdt_;
phaseSystem::dmdtfTable dmdtfs_;
//- Reference to the mesh
const fvMesh& mesh_;
@ -336,6 +336,8 @@ private:
void registerSizeGroups(sizeGroup& group);
void initialiseDmdtfs();
void createPhasePairs();
void precompute();
@ -380,19 +382,10 @@ public:
TypeName("populationBalanceModel");
// Constructor
// Constructors
populationBalanceModel
(
const phaseSystem& fluid,
const word& name,
HashPtrTable
<
volScalarField,
phaseInterfaceKey,
phaseInterfaceKey::hash
>& pDmdt
);
//- Construct for a fluid
populationBalanceModel(const phaseSystem& fluid, const word& name);
//- Return clone
autoPtr<populationBalanceModel> clone() const;
@ -403,29 +396,22 @@ public:
{
const phaseSystem& fluid_;
phaseSystem::dmdtfTable& pDmdt_;
public:
iNew
(
const phaseSystem& fluid,
phaseSystem::dmdtfTable& pDmdt
)
iNew(const phaseSystem& fluid)
:
fluid_(fluid),
pDmdt_(pDmdt)
fluid_(fluid)
{}
autoPtr<populationBalanceModel> operator()(Istream& is) const
{
word name(is);
const word name(is);
Info << "Setting up population balance: " << name << endl;
return autoPtr<populationBalanceModel>
(
new populationBalanceModel(fluid_, name, pDmdt_)
new populationBalanceModel(fluid_, name)
);
}
};
@ -434,6 +420,7 @@ public:
//- Destructor
virtual ~populationBalanceModel();
// Member Functions
//- Dummy write for regIOobject
@ -442,6 +429,9 @@ public:
//- Return reference to the phaseSystem
inline const phaseSystem& fluid() const;
//- Return reference to the interfacial mass transfer rates
inline const phaseSystem::dmdtfTable& dmdtfs() const;
//- Return reference to the mesh
inline const fvMesh& mesh() const;

View File

@ -43,6 +43,13 @@ Foam::diameterModels::populationBalanceModel::fluid() const
}
inline const Foam::phaseSystem::dmdtfTable&
Foam::diameterModels::populationBalanceModel::dmdtfs() const
{
return dmdtfs_;
}
inline const Foam::fvMesh&
Foam::diameterModels::populationBalanceModel::mesh() const
{