From d47a42458f4bbaacd4ec41c68aadab3f7e98ea11 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Thu, 16 Jul 2015 14:12:03 +0100 Subject: [PATCH] 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. --- .../HeatAndMassTransferPhaseSystem.C | 264 -------------- .../HeatAndMassTransferPhaseSystem.H | 23 +- ...terfaceCompositionPhaseChangePhaseSystem.C | 336 ++++++++++++++++++ ...terfaceCompositionPhaseChangePhaseSystem.H | 120 +++++++ .../ThermalPhaseChangePhaseSystem.C | 199 +++++++++++ .../ThermalPhaseChangePhaseSystem.H | 112 ++++++ .../twoPhaseSystem/twoPhaseSystems.C | 27 +- 7 files changed, 793 insertions(+), 288 deletions(-) create mode 100644 applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C create mode 100644 applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.H create mode 100644 applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C create mode 100644 applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.H diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C index 71c5a06724..ff4fb10e9d 100644 --- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C +++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C @@ -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::heatTransfer() const } -template -Foam::autoPtr -Foam::HeatAndMassTransferPhaseSystem::massTransfer() const -{ - // Create a mass transfer matrix for each species of each phase - autoPtr eqnsPtr - ( - new phaseSystem::massTransferTable() - ); - - phaseSystem::massTransferTable& eqns = eqnsPtr(); - - forAllConstIter - ( - phaseSystem::phaseModelTable, - this->phaseModels_, - phaseModelIter - ) - { - const phaseModel& phase(phaseModelIter()); - - const PtrList& 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::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 -void Foam::HeatAndMassTransferPhaseSystem::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 bool Foam::HeatAndMassTransferPhaseSystem::read() { diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.H b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.H index 27e4f32a77..3f20c575f8 100644 --- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.H +++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.H @@ -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, - 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 - massTransfer() const; + massTransfer() const = 0; //- Correct the thermodynamics - virtual void correctThermo(); + virtual void correctThermo() = 0; //- Read base phaseProperties dictionary virtual bool read(); diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C new file mode 100644 index 0000000000..47c994334d --- /dev/null +++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C @@ -0,0 +1,336 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / 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 . + +\*---------------------------------------------------------------------------*/ + +#include "InterfaceCompositionPhaseChangePhaseSystem.H" +#include "interfaceCompositionModel.H" + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template +Foam::InterfaceCompositionPhaseChangePhaseSystem:: +InterfaceCompositionPhaseChangePhaseSystem +( + const fvMesh& mesh +) +: + HeatAndMassTransferPhaseSystem(mesh) +{ + this->generatePairsAndSubModels + ( + "interfaceComposition", + interfaceCompositionModels_ + ); +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +template +Foam::InterfaceCompositionPhaseChangePhaseSystem:: +~InterfaceCompositionPhaseChangePhaseSystem() +{} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + +template +Foam::autoPtr +Foam::InterfaceCompositionPhaseChangePhaseSystem:: +massTransfer() const +{ + // Create a mass transfer matrix for each species of each phase + autoPtr eqnsPtr + ( + new phaseSystem::massTransferTable() + ); + + phaseSystem::massTransferTable& eqns = eqnsPtr(); + + forAllConstIter + ( + phaseSystem::phaseModelTable, + this->phaseModels_, + phaseModelIter + ) + { + const phaseModel& phase(phaseModelIter()); + + const PtrList& 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::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; +} + + +template +void Foam::InterfaceCompositionPhaseChangePhaseSystem:: +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(this->heatTransferModels_[pair][pair.first()]->K()); + volScalarField H2(this->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 = *this->Tf_[pair]; + + // Add latent heats from forward and backward models + if (this->interfaceCompositionModels_.found(key12)) + { + this->interfaceCompositionModels_[key12]->addMDotL + ( + this->massTransferModels_[pair][pair.first()]->K(), + Tf, + mDotL, + mDotLPrime + ); + } + if (this->interfaceCompositionModels_.found(key21)) + { + this->interfaceCompositionModels_[key21]->addMDotL + ( + this->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) + ); + + Tf.correctBoundaryConditions(); + + Info<< "Tf." << pair.name() + << ": min = " << min(Tf.internalField()) + << ", mean = " << average(Tf.internalField()) + << ", max = " << max(Tf.internalField()) + << endl; + + // Update the interface compositions + if (this->interfaceCompositionModels_.found(key12)) + { + this->interfaceCompositionModels_[key12]->update(Tf); + } + if (this->interfaceCompositionModels_.found(key21)) + { + this->interfaceCompositionModels_[key21]->update(Tf); + } + } +} + + +template +bool Foam::InterfaceCompositionPhaseChangePhaseSystem::read() +{ + if (BasePhaseSystem::read()) + { + bool readOK = true; + + // Models ... + + return readOK; + } + else + { + return false; + } +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.H b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.H new file mode 100644 index 0000000000..b87c5a709d --- /dev/null +++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.H @@ -0,0 +1,120 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / 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 . + +Class + Foam::InterfaceCompositionPhaseChangePhaseSystem + +Description + Class to provide interfacial heat and mass transfer between a number of + phases according to a interface composition model. + + 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. + +SourceFiles + InterfaceCompositionPhaseChangePhaseSystem.C + +\*---------------------------------------------------------------------------*/ + +#ifndef InterfaceCompositionPhaseChangePhaseSystem_H +#define InterfaceCompositionPhaseChangePhaseSystem_H + +#include "HeatAndMassTransferPhaseSystem.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +class interfaceCompositionModel; + +/*---------------------------------------------------------------------------*\ + Class InterfaceCompositionPhaseChangePhaseSystem Declaration +\*---------------------------------------------------------------------------*/ + +template +class InterfaceCompositionPhaseChangePhaseSystem +: + public HeatAndMassTransferPhaseSystem +{ +protected: + + // Protected typedefs + + typedef HashTable + < + autoPtr, + phasePairKey, + phasePairKey::hash + > interfaceCompositionModelTable; + + + // Protected data + + // Sub Models + + //- Interface composition models + interfaceCompositionModelTable interfaceCompositionModels_; + + +public: + + // Constructors + + //- Construct from fvMesh + InterfaceCompositionPhaseChangePhaseSystem(const fvMesh&); + + + //- Destructor + virtual ~InterfaceCompositionPhaseChangePhaseSystem(); + + + // Member Functions + + //- Return the mass transfer matrices + virtual autoPtr massTransfer() const; + + //- Correct the thermodynamics + virtual void correctThermo(); + + //- Read base phaseProperties dictionary + virtual bool read(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "InterfaceCompositionPhaseChangePhaseSystem.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C new file mode 100644 index 0000000000..bbce1962a1 --- /dev/null +++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C @@ -0,0 +1,199 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / 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 . + +\*---------------------------------------------------------------------------*/ + +#include "ThermalPhaseChangePhaseSystem.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template +Foam::ThermalPhaseChangePhaseSystem:: +ThermalPhaseChangePhaseSystem +( + const fvMesh& mesh +) +: + HeatAndMassTransferPhaseSystem(mesh), + volatile_(this->lookup("volatile")), + saturationModel_(saturationModel::New(this->subDict("saturationModel"))) +{} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +template +Foam::ThermalPhaseChangePhaseSystem:: +~ThermalPhaseChangePhaseSystem() +{} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + +template +Foam::autoPtr +Foam::ThermalPhaseChangePhaseSystem::massTransfer() const +{ + // Create a mass transfer matrix for each species of each phase + autoPtr eqnsPtr + ( + new phaseSystem::massTransferTable() + ); + + phaseSystem::massTransferTable& eqns = eqnsPtr(); + + forAllConstIter + ( + phaseSystem::phaseModelTable, + this->phaseModels_, + phaseModelIter + ) + { + const phaseModel& phase(phaseModelIter()); + + const PtrList& Yi = phase.Y(); + + forAll(Yi, i) + { + eqns.insert + ( + Yi[i].name(), + new fvScalarMatrix(Yi[i], dimMass/dimTime) + ); + } + } + + forAllConstIter + ( + phaseSystem::phasePairTable, + this->phasePairs_, + phasePairIter + ) + { + const phasePair& pair(phasePairIter()); + + if (pair.ordered()) + { + continue; + } + const phaseModel& phase = pair.phase1(); + const phaseModel& otherPhase = pair.phase2(); + + const word name + ( + IOobject::groupName(volatile_, phase.name()) + ); + + const word otherName + ( + IOobject::groupName(volatile_, otherPhase.name()) + ); + + const volScalarField dmdt(this->dmdt(pair)); + const volScalarField dmdt12(posPart(dmdt)); + const volScalarField dmdt21(negPart(dmdt)); + + *eqns[name] += fvm::Sp(dmdt21, eqns[name]->psi()) - dmdt21; + *eqns[otherName] += dmdt12 - fvm::Sp(dmdt12, eqns[otherName]->psi()); + } + + return eqnsPtr; +} + + +template +void Foam::ThermalPhaseChangePhaseSystem::correctThermo() +{ + BasePhaseSystem::correctThermo(); + + forAllConstIter + ( + phaseSystem::phasePairTable, + this->phasePairs_, + phasePairIter + ) + { + const phasePair& pair(phasePairIter()); + + if (pair.ordered()) + { + continue; + } + + const phaseModel& phase1 = pair.phase1(); + const phaseModel& phase2 = pair.phase2(); + + volScalarField& Tf = *this->Tf_[pair]; + Tf = saturationModel_->Tsat(phase1.thermo().p()); + + Info<< "Tf." << pair.name() + << ": min = " << min(Tf.internalField()) + << ", mean = " << average(Tf.internalField()) + << ", max = " << max(Tf.internalField()) + << endl; + + volScalarField& dmdt(*this->dmdt_[pair]); + + volScalarField H1(this->heatTransferModels_[pair][pair.first()]->K()); + volScalarField H2(this->heatTransferModels_[pair][pair.second()]->K()); + + const volScalarField& T1(phase1.thermo().T()); + const volScalarField& T2(phase2.thermo().T()); + + const volScalarField& he1(phase1.thermo().he()); + const volScalarField& he2(phase2.thermo().he()); + + volScalarField hef2(phase2.thermo().he(phase2.thermo().p(), Tf)); + volScalarField hef1(phase1.thermo().he(phase1.thermo().p(), Tf)); + + dmdt = + (H2*(Tf - T1) + H1*(Tf - T2)) + /min + ( + (pos(dmdt)*he2 + neg(dmdt)*hef2) + - (neg(dmdt)*he1 + pos(dmdt)*hef1), + 0.3*(hef1 - hef2) + ); + } +} + + +template +bool Foam::ThermalPhaseChangePhaseSystem::read() +{ + if (BasePhaseSystem::read()) + { + bool readOK = true; + + // Models ... + + return readOK; + } + else + { + return false; + } +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.H b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.H new file mode 100644 index 0000000000..9a0a0d109a --- /dev/null +++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.H @@ -0,0 +1,112 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / 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 . + +Class + Foam::ThermalPhaseChangePhaseSystem + +Description + Class to provide interfacial heat and mass transfer between a number of + phases according the interfacial temperature approximated by the saturation + temperature. + + Currently only a single specified specie is considered volatile and changes + phase, all other species are considered nonvolatile and do not + affect the mass-transfer. + +SourceFiles + ThermalPhaseChangePhaseSystem.C + +\*---------------------------------------------------------------------------*/ + +#ifndef ThermalPhaseChangePhaseSystem_H +#define ThermalPhaseChangePhaseSystem_H + +#include "HeatAndMassTransferPhaseSystem.H" +#include "saturationModel.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class ThermalPhaseChangePhaseSystem Declaration +\*---------------------------------------------------------------------------*/ + +template +class ThermalPhaseChangePhaseSystem +: + public HeatAndMassTransferPhaseSystem +{ + +protected: + + // Protected data + + //- Name of the volatile specie + word volatile_; + + //- The saturation model used to evaluate Tsat = Tf + autoPtr saturationModel_; + + +public: + + // Constructors + + //- Construct from fvMesh + ThermalPhaseChangePhaseSystem(const fvMesh&); + + + //- Destructor + virtual ~ThermalPhaseChangePhaseSystem(); + + + // Member Functions + + //- Return the mass transfer matrices + virtual autoPtr massTransfer() const; + + //- Correct the thermodynamics + virtual void correctThermo(); + + //- Read base phaseProperties dictionary + virtual bool read(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "ThermalPhaseChangePhaseSystem.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/twoPhaseSystem/twoPhaseSystems.C b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/twoPhaseSystem/twoPhaseSystems.C index 6b649afb64..075a578c23 100644 --- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/twoPhaseSystem/twoPhaseSystems.C +++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/twoPhaseSystem/twoPhaseSystems.C @@ -29,7 +29,8 @@ License #include "twoPhaseSystem.H" #include "MomentumTransferPhaseSystem.H" #include "HeatTransferPhaseSystem.H" -#include "HeatAndMassTransferPhaseSystem.H" +#include "InterfaceCompositionPhaseChangePhaseSystem.H" +#include "ThermalPhaseChangePhaseSystem.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -51,19 +52,35 @@ namespace Foam ); typedef - HeatAndMassTransferPhaseSystem + InterfaceCompositionPhaseChangePhaseSystem < MomentumTransferPhaseSystem > - heatMassAndMomentumTransferTwoPhaseSystem; + interfaceCompositionPhaseChangeTwoPhaseSystem; addNamedToRunTimeSelectionTable ( twoPhaseSystem, - heatMassAndMomentumTransferTwoPhaseSystem, + interfaceCompositionPhaseChangeTwoPhaseSystem, dictionary, - heatMassAndMomentumTransferTwoPhaseSystem + interfaceCompositionPhaseChangeTwoPhaseSystem + ); + + typedef + ThermalPhaseChangePhaseSystem + < + MomentumTransferPhaseSystem + > + thermalPhaseChangeTwoPhaseSystem; + + addNamedToRunTimeSelectionTable + ( + twoPhaseSystem, + thermalPhaseChangeTwoPhaseSystem, + dictionary, + thermalPhaseChangeTwoPhaseSystem ); } + // ************************************************************************* //