Commit d47a4245 authored by Henry Weller's avatar Henry Weller
Browse files

reactingTwoPhaseEulerFoam: Added support for thermally-driven phase-change (boiling)

The interfacial temperature is assumed equal to the saturation
temperature.  Only a single species is considered volatile and the other
species to not affect the mass-transfer.
parent 82241e64
......@@ -28,7 +28,6 @@ License
#include "BlendedInterfacialModel.H"
#include "heatTransferModel.H"
#include "massTransferModel.H"
#include "interfaceCompositionModel.H"
#include "HashPtrTable.H"
......@@ -59,12 +58,6 @@ HeatAndMassTransferPhaseSystem
massTransferModels_
);
this->generatePairsAndSubModels
(
"interfaceComposition",
interfaceCompositionModels_
);
forAllConstIter
(
phaseSystem::phasePairTable,
......@@ -332,263 +325,6 @@ Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::heatTransfer() const
}
template<class BasePhaseSystem>
Foam::autoPtr<Foam::phaseSystem::massTransferTable>
Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::massTransfer() const
{
// Create a mass transfer matrix for each species of each phase
autoPtr<phaseSystem::massTransferTable> eqnsPtr
(
new phaseSystem::massTransferTable()
);
phaseSystem::massTransferTable& eqns = eqnsPtr();
forAllConstIter
(
phaseSystem::phaseModelTable,
this->phaseModels_,
phaseModelIter
)
{
const phaseModel& phase(phaseModelIter());
const PtrList<volScalarField>& Yi = phase.Y();
forAll(Yi, i)
{
eqns.insert
(
Yi[i].name(),
new fvScalarMatrix(Yi[i], dimMass/dimTime)
);
}
}
// Reset the interfacial mass flow rates
forAllConstIter
(
phaseSystem::phasePairTable,
this->phasePairs_,
phasePairIter
)
{
const phasePair& pair(phasePairIter());
if (pair.ordered())
{
continue;
}
*dmdt_[pair] =
*dmdtExplicit_[pair];
*dmdtExplicit_[pair] =
dimensionedScalar("zero", dimDensity/dimTime, 0);
}
// Sum up the contribution from each interface composition model
forAllConstIter
(
interfaceCompositionModelTable,
interfaceCompositionModels_,
interfaceCompositionModelIter
)
{
const interfaceCompositionModel& compositionModel
(
interfaceCompositionModelIter()
);
const phasePair& pair
(
this->phasePairs_[interfaceCompositionModelIter.key()]
);
const phaseModel& phase = pair.phase1();
const phaseModel& otherPhase = pair.phase2();
const phasePairKey key(phase.name(), otherPhase.name());
const volScalarField& Tf(*Tf_[key]);
volScalarField& dmdtExplicit(*dmdtExplicit_[key]);
volScalarField& dmdt(*dmdt_[key]);
scalar dmdtSign(Pair<word>::compare(dmdt_.find(key).key(), key));
const volScalarField K
(
massTransferModels_[key][phase.name()]->K()
);
forAllConstIter
(
hashedWordList,
compositionModel.species(),
memberIter
)
{
const word& member = *memberIter;
const word name
(
IOobject::groupName(member, phase.name())
);
const word otherName
(
IOobject::groupName(member, otherPhase.name())
);
const volScalarField KD
(
K*compositionModel.D(member)
);
const volScalarField Yf
(
compositionModel.Yf(member, Tf)
);
// Implicit transport through the phase
*eqns[name] +=
phase.rho()*KD*Yf
- fvm::Sp(phase.rho()*KD, eqns[name]->psi());
// Sum the mass transfer rate
dmdtExplicit += dmdtSign*phase.rho()*KD*Yf;
dmdt -= dmdtSign*phase.rho()*KD*eqns[name]->psi();
// Explicit transport out of the other phase
if (eqns.found(otherName))
{
*eqns[otherName] -=
otherPhase.rho()*KD*compositionModel.dY(member, Tf);
}
}
}
return eqnsPtr;
}
template<class BasePhaseSystem>
void Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::correctThermo()
{
BasePhaseSystem::correctThermo();
// This loop solves for the interface temperatures, Tf, and updates the
// interface composition models.
//
// The rate of heat transfer to the interface must equal the latent heat
// consumed at the interface, i.e.:
//
// H1*(T1 - Tf) + H2*(T2 - Tf) == mDotL
// == K*rho*(Yfi - Yi)*Li
//
// Yfi is likely to be a strong non-linear (typically exponential) function
// of Tf, so the solution for the temperature is newton-accelerated
forAllConstIter
(
phaseSystem::phasePairTable,
this->phasePairs_,
phasePairIter
)
{
const phasePair& pair(phasePairIter());
if (pair.ordered())
{
continue;
}
const phasePairKey key12(pair.first(), pair.second(), true);
const phasePairKey key21(pair.second(), pair.first(), true);
volScalarField H1(heatTransferModels_[pair][pair.first()]->K());
volScalarField H2(heatTransferModels_[pair][pair.second()]->K());
dimensionedScalar HSmall("small", heatTransferModel::dimK, SMALL);
volScalarField mDotL
(
IOobject
(
"mDotL",
this->mesh().time().timeName(),
this->mesh()
),
this->mesh(),
dimensionedScalar("zero", dimEnergy/dimVolume/dimTime, 0)
);
volScalarField mDotLPrime
(
IOobject
(
"mDotLPrime",
this->mesh().time().timeName(),
this->mesh()
),
this->mesh(),
dimensionedScalar("zero", mDotL.dimensions()/dimTemperature, 0)
);
volScalarField& Tf = *Tf_[pair];
// Add latent heats from forward and backward models
if (interfaceCompositionModels_.found(key12))
{
interfaceCompositionModels_[key12]->addMDotL
(
massTransferModels_[pair][pair.first()]->K(),
Tf,
mDotL,
mDotLPrime
);
}
if (interfaceCompositionModels_.found(key21))
{
interfaceCompositionModels_[key21]->addMDotL
(
massTransferModels_[pair][pair.second()]->K(),
Tf,
mDotL,
mDotLPrime
);
}
// Update the interface temperature by applying one step of newton's
// method to the interface relation
Tf -=
(
H1*(Tf - pair.phase1().thermo().T())
+ H2*(Tf - pair.phase2().thermo().T())
+ mDotL
)
/(
max(H1 + H2 + mDotLPrime, HSmall)
);
// Update the interface compositions
if (interfaceCompositionModels_.found(key12))
{
interfaceCompositionModels_[key12]->update(Tf);
}
if (interfaceCompositionModels_.found(key21))
{
interfaceCompositionModels_[key21]->update(Tf);
}
Tf.correctBoundaryConditions();
Info<< "Tf." << pair.name()
<< ": min = " << min(Tf.internalField())
<< ", mean = " << average(Tf.internalField())
<< ", max = " << max(Tf.internalField())
<< endl;
}
}
template<class BasePhaseSystem>
bool Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::read()
{
......
......@@ -25,12 +25,8 @@ Class
Foam::HeatAndMassTransferPhaseSystem
Description
Class which models interfacial heat and mass transfer between a number of
phases. Mass is transferred to or from a phase according to a composition
model at the surface of another phase. Heat is transferred from a phase to
an interfacial temperature. The interface temperature is calculated such
that the net rate at which the heat is transferred to the interface is
equal to the latent heat consumed by the mass transfer.
Base class to support interfacial heat and mass transfer between a number
of phases.
SourceFiles
HeatAndMassTransferPhaseSystem.C
......@@ -53,7 +49,6 @@ class BlendedInterfacialModel;
class blendingMethod;
class heatTransferModel;
class massTransferModel;
class interfaceCompositionModel;
/*---------------------------------------------------------------------------*\
Class HeatAndMassTransferPhaseSystem Declaration
......@@ -88,13 +83,6 @@ protected:
phasePairKey::hash
> massTransferModelTable;
typedef HashTable
<
autoPtr<interfaceCompositionModel>,
phasePairKey,
phasePairKey::hash
> interfaceCompositionModelTable;
// Protected data
......@@ -117,9 +105,6 @@ protected:
//- Mass transfer models
massTransferModelTable massTransferModels_;
//- Interface composition models
interfaceCompositionModelTable interfaceCompositionModels_;
public:
......@@ -151,10 +136,10 @@ public:
//- Return the mass transfer matrices
virtual autoPtr<phaseSystem::massTransferTable>
massTransfer() const;
massTransfer() const = 0;
//- Correct the thermodynamics
virtual void correctThermo();
virtual void correctThermo() = 0;
//- Read base phaseProperties dictionary
virtual bool read();
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "InterfaceCompositionPhaseChangePhaseSystem.H"
#include "interfaceCompositionModel.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class BasePhaseSystem>
Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::
InterfaceCompositionPhaseChangePhaseSystem
(
const fvMesh& mesh
)
:
HeatAndMassTransferPhaseSystem<BasePhaseSystem>(mesh)
{
this->generatePairsAndSubModels
(
"interfaceComposition",
interfaceCompositionModels_
);
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
template<class BasePhaseSystem>
Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::
~InterfaceCompositionPhaseChangePhaseSystem()
{}
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
template<class BasePhaseSystem>
Foam::autoPtr<Foam::phaseSystem::massTransferTable>
Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::
massTransfer() const
{
// Create a mass transfer matrix for each species of each phase
autoPtr<phaseSystem::massTransferTable> eqnsPtr
(
new phaseSystem::massTransferTable()
);
phaseSystem::massTransferTable& eqns = eqnsPtr();
forAllConstIter
(
phaseSystem::phaseModelTable,
this->phaseModels_,
phaseModelIter
)
{
const phaseModel& phase(phaseModelIter());
const PtrList<volScalarField>& Yi = phase.Y();
forAll(Yi, i)
{
eqns.insert
(
Yi[i].name(),
new fvScalarMatrix(Yi[i], dimMass/dimTime)
);
}
}
// Reset the interfacial mass flow rates
forAllConstIter
(
phaseSystem::phasePairTable,
this->phasePairs_,
phasePairIter
)
{
const phasePair& pair(phasePairIter());
if (pair.ordered())
{
continue;
}
*this->dmdt_[pair] =
*this->dmdtExplicit_[pair];
*this->dmdtExplicit_[pair] =
dimensionedScalar("zero", dimDensity/dimTime, 0);
}
// Sum up the contribution from each interface composition model
forAllConstIter
(
interfaceCompositionModelTable,
interfaceCompositionModels_,
interfaceCompositionModelIter
)
{
const interfaceCompositionModel& compositionModel
(
interfaceCompositionModelIter()
);
const phasePair& pair
(
this->phasePairs_[interfaceCompositionModelIter.key()]
);
const phaseModel& phase = pair.phase1();
const phaseModel& otherPhase = pair.phase2();
const phasePairKey key(phase.name(), otherPhase.name());
const volScalarField& Tf(*this->Tf_[key]);
volScalarField& dmdtExplicit(*this->dmdtExplicit_[key]);
volScalarField& dmdt(*this->dmdt_[key]);
scalar dmdtSign(Pair<word>::compare(this->dmdt_.find(key).key(), key));
const volScalarField K
(
this->massTransferModels_[key][phase.name()]->K()
);
forAllConstIter
(
hashedWordList,
compositionModel.species(),
memberIter
)
{
const word& member = *memberIter;
const word name
(
IOobject::groupName(member, phase.name())
);
const word otherName
(
IOobject::groupName(member, otherPhase.name())
);
const volScalarField KD
(
K*compositionModel.D(member)
);
const volScalarField Yf
(
compositionModel.Yf(member, Tf)
);
// Implicit transport through the phase
*eqns[name] +=
phase.rho()*KD*Yf
- fvm::Sp(phase.rho()*KD, eqns[name]->psi());
// Sum the mass transfer rate
dmdtExplicit += dmdtSign*phase.rho()*KD*Yf;
dmdt -= dmdtSign*phase.rho()*KD*eqns[name]->psi();
// Explicit transport out of the other phase
if (eqns.found(otherName))
{
*eqns[otherName] -=
otherPhase.rho()*KD*compositionModel.dY(member, Tf);
}
}
}
return eqnsPtr;
}