reacting*EulerFoam/.../ThermalPhaseChange: Update to new definition of continuity error

The ThermalPhaseChangePhaseSystem stores the thermal phase change dmdt
used in the previous continuity error update and uses that to stabilize
the interfacial heat transfer calculations when phase fractions approach
zero.

Patch contributed by Juho Peltola, VTT.
This commit is contained in:
Will Bainbridge
2020-02-27 10:25:02 +00:00
parent 57e6515379
commit 047e9640b6
2 changed files with 45 additions and 13 deletions

View File

@ -68,7 +68,8 @@ ThermalPhaseChangePhaseSystem
)
:
BasePhaseSystem(mesh),
volatile_(this->template lookupOrDefault<word>("volatile", "none"))
volatile_(this->template lookupOrDefault<word>("volatile", "none")),
dmdt0s_(this->phases().size())
{
this->generatePairsAndSubModels
(
@ -325,9 +326,15 @@ Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::heatTransfer() const
const volScalarField dmdtf21(posPart(dmdtf));
const volScalarField dmdtf12(negPart(dmdtf));
*eqns[phase1.name()] += - fvm::Sp(dmdtf21, he1) + dmdtf21*(K2 - K1);
*eqns[phase1.name()] +=
fvm::Sp(dmdt0s_[phase1.index()],he1)
- fvm::Sp(dmdtf21, he1)
+ dmdtf21*K2 + dmdtf12*K1;
*eqns[phase2.name()] -= - fvm::Sp(dmdtf12, he2) + dmdtf12*(K1 - K2);
*eqns[phase2.name()] -=
- fvm::Sp(dmdt0s_[phase2.index()],he2)
- fvm::Sp(dmdtf12, he2)
+ dmdtf12*K1 + dmdtf21*K2;
if (this->saturationModels_.found(phasePairIter.key()))
{
@ -376,14 +383,12 @@ Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::heatTransfer() const
else
{
*eqns[phase1.name()] += dmdtf21*he2;
*eqns[phase2.name()] -= dmdtf12*he1;
}
if (this->nDmdtLfs_.found(phasePairIter.key()))
{
*eqns[phase1.name()] += negPart(*this->nDmdtLfs_[pair]);
*eqns[phase2.name()] -= posPart(*this->nDmdtLfs_[pair]);
}
@ -398,7 +403,6 @@ Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::heatTransfer() const
*eqns[phase2.name()] -=
phase2.thermo().p()*dmdtf/phase2.thermo().rho();
}
}
return eqnsPtr;
@ -457,6 +461,34 @@ Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::specieTransfer() const
}
template<class BasePhaseSystem>
void Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::
correctContinuityError()
{
dmdt0s_ = PtrList<volScalarField>(this->phases().size());
forAllConstIter(phaseSystem::dmdtfTable, dmdtfs_, dmdtfIter)
{
const phasePair& pair = this->phasePairs_[dmdtfIter.key()];
const volScalarField& dmdtf = *dmdtfIter();
addField(pair.phase1(), "dmdt", dmdtf, dmdt0s_);
addField(pair.phase2(), "dmdt", - dmdtf, dmdt0s_);
}
forAllConstIter(phaseSystem::dmdtfTable, nDmdtfs_, nDmdtfIter)
{
const phasePair& pair = this->phasePairs_[nDmdtfIter.key()];
const volScalarField& nDmdtf = *nDmdtfIter();
addField(pair.phase1(), "dmdt", nDmdtf, dmdt0s_);
addField(pair.phase2(), "dmdt", - nDmdtf, dmdt0s_);
}
BasePhaseSystem::correctContinuityError();
}
template<class BasePhaseSystem>
void
Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::correctInterfaceThermo()

View File

@ -91,16 +91,9 @@ class ThermalPhaseChangePhaseSystem
//- Name of the volatile specie
word volatile_;
//- The saturation model used to evaluate Tsat = Tf
//autoPtr<saturationModel> saturationModel_;
//- Saturation models used to evaluate Tsat = Tf
saturationModelTable saturationModels_;
//- PhasePairs between which phaseChange occurs, e.g.,
// "((gasI and liquid) (gasII and liquid))"
//List<phasePairKey> phaseChangePairKeys_;
//- Mass transfer rates
phaseSystem::dmdtfTable dmdtfs_;
@ -110,6 +103,10 @@ class ThermalPhaseChangePhaseSystem
//- Nucleate thermal energy transfer rates
phaseSystem::dmdtfTable nDmdtLfs_;
//- Previous continuity error update phase dmdts for the heat transfer
// function
PtrList<volScalarField> dmdt0s_;
// Private member functions
@ -153,6 +150,9 @@ public:
virtual autoPtr<phaseSystem::specieTransferTable>
specieTransfer() const;
//- Store phase dmdts at the during the continuity error update
virtual void correctContinuityError();
//- Correct the interface thermodynamics
virtual void correctInterfaceThermo();