From 68ea75a37cd3306769c45303ba040b9201487826 Mon Sep 17 00:00:00 2001
From: Henry Weller <http://cfd.direct>
Date: Thu, 2 Jul 2015 22:51:06 +0100
Subject: [PATCH] reactingTwoPhaseEulerFoam: Corrected handling of
heat-transfer caused by mass-transfer
---
.../HeatAndMassTransferPhaseSystem.C | 32 +++++++++++++++++++
.../AnisothermalPhaseModel.C | 8 +++++
.../AnisothermalPhaseModel.H | 7 ++--
.../MovingPhaseModel/MovingPhaseModel.C | 8 -----
.../MovingPhaseModel/MovingPhaseModel.H | 3 --
.../phaseModel/phaseModel/phaseModel.C | 7 ++++
.../phaseModel/phaseModel/phaseModel.H | 6 ++--
7 files changed, 55 insertions(+), 16 deletions(-)
diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C
index ba153f166a8..ab62b151362 100644
--- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C
+++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C
@@ -288,6 +288,38 @@ Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::heatTransfer() const
}
}
+ // Source term due to mass trasfer
+ forAllConstIter
+ (
+ phaseSystem::phasePairTable,
+ this->phasePairs_,
+ phasePairIter
+ )
+ {
+ const phasePair& pair(phasePairIter());
+
+ if (pair.ordered())
+ {
+ continue;
+ }
+
+ const volScalarField& he1(pair.phase1().thermo().he());
+ const volScalarField& he2(pair.phase2().thermo().he());
+
+ const volScalarField& K1(pair.phase1().K());
+ const volScalarField& K2(pair.phase2().K());
+
+ const volScalarField dmdt(this->dmdt(pair));
+ const volScalarField dmdt12(dmdt*pos(dmdt));
+ const volScalarField dmdt21(dmdt*neg(dmdt));
+
+ *eqns[pair.phase1().name()] +=
+ fvm::Sp(dmdt21, he1) + dmdt21*K1 - dmdt21*(he2 + K2);
+
+ *eqns[pair.phase2().name()] +=
+ dmdt12*(he1 + K1) - fvm::Sp(dmdt12, he2) - dmdt12*K2;
+ }
+
return eqnsPtr;
}
diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.C b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.C
index cbefea2852e..3a2786d75c1 100644
--- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.C
+++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.C
@@ -154,4 +154,12 @@ Foam::AnisothermalPhaseModel<BasePhaseModel>::divU(const volScalarField& divU)
}
+template<class BasePhaseModel>
+const Foam::volScalarField&
+Foam::AnisothermalPhaseModel<BasePhaseModel>::K() const
+{
+ return K_;
+}
+
+
// ************************************************************************* //
diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.H b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.H
index 48b0cbd4991..e7a76bf496e 100644
--- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.H
+++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.H
@@ -89,11 +89,14 @@ public:
//- Return true if the phase is compressible otherwise false
virtual bool compressible() const;
- //- Phase dilatation rate (d(alpha)/dt + div(alpha*phi))
+ //- Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
virtual const volScalarField& divU() const;
- //- Set phase dilatation rate (d(alpha)/dt + div(alpha*phi))
+ //- Set the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
virtual void divU(const volScalarField& divU);
+
+ //- Return the phase kinetic energy
+ virtual const volScalarField& K() const;
};
diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.C b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.C
index aeeef98149c..b1869b9fbfe 100644
--- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.C
+++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.C
@@ -314,14 +314,6 @@ Foam::MovingPhaseModel<BasePhaseModel>::continuityError() const
}
-template<class BasePhaseModel>
-Foam::volScalarField&
-Foam::MovingPhaseModel<BasePhaseModel>::continuityError()
-{
- return continuityError_;
-}
-
-
template<class BasePhaseModel>
Foam::tmp<Foam::surfaceScalarField>
Foam::MovingPhaseModel<BasePhaseModel>::phi() const
diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.H b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.H
index 96e8942f6fc..14e9c057b7e 100644
--- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.H
+++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.H
@@ -152,9 +152,6 @@ public:
//- Constant access the continuity error
virtual tmp<volScalarField> continuityError() const;
- //- Access the continuity error
- virtual volScalarField& continuityError();
-
//- Constant access the volumetric flux
virtual tmp<surfaceScalarField> phi() const;
diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/phaseModel/phaseModel.C b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/phaseModel/phaseModel.C
index 2697587b6a5..61ea75e586c 100644
--- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/phaseModel/phaseModel.C
+++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/phaseModel/phaseModel.C
@@ -155,6 +155,13 @@ void Foam::phaseModel::divU(const volScalarField& divU)
}
+const Foam::volScalarField& Foam::phaseModel::K() const
+{
+ notImplemented("Foam::phaseModel::K()");
+ return volScalarField::null();
+}
+
+
const Foam::surfaceScalarField& Foam::phaseModel::DbyA() const
{
return surfaceScalarField::null();
diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/phaseModel/phaseModel.H b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/phaseModel/phaseModel.H
index de85779d550..082bb6c0191 100644
--- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/phaseModel/phaseModel.H
+++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/phaseSystems/phaseModel/phaseModel/phaseModel.H
@@ -175,6 +175,9 @@ public:
//- Set the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
virtual void divU(const volScalarField& divU);
+ //- Return the phase kinetic energy
+ virtual const volScalarField& K() const;
+
// Implicit phase pressure and dispersion support
@@ -218,9 +221,6 @@ public:
//- Constant access the continuity error
virtual tmp<volScalarField> continuityError() const = 0;
- //- Access the continuity error
- virtual volScalarField& continuityError() = 0;
-
//- Constant access the volumetric flux
virtual tmp<surfaceScalarField> phi() const = 0;
--
GitLab