From a7a6ac5104044383bb650b7e8c745764e094013c Mon Sep 17 00:00:00 2001
From: sergio <sergio>
Date: Fri, 27 Oct 2017 15:36:46 -0700
Subject: [PATCH] ENH: Final stages for mass transfer models and laser
modification for interface reflection
---
.../icoReactingMultiPhaseInterFoam/TEqn.H | 7 +-
.../icoReactingMultiphaseInterFoam.C | 4 +-
.../laserDTRM/DTRMParticle/DTRMParticle.C | 100 ++--
.../laserDTRM/DTRMParticle/DTRMParticle.H | 30 +-
.../laserDTRM/DTRMParticle/DTRMParticleI.H | 10 +-
.../laserDTRM/Make/options | 2 +-
.../laserDTRM/laserDTRM.C | 213 ++++---
.../laserDTRM/laserDTRM.H | 43 +-
.../localDensityAbsorptionEmission.C | 18 +-
.../FresnelLaser/FresnelLaser.C | 3 +-
.../FresnelLaser/FresnelLaser.H | 4 -
.../reflectionModel/reflectionModelNew.C | 2 +-
.../InterfaceCompositionModels.C | 144 +++--
.../massTransferModels/Lee/Lee.C | 46 +-
.../massTransferModels/Lee/Lee.H | 13 +-
.../interfaceCompositionModel.C | 27 +-
.../interfaceCompositionModel.H | 30 +-
.../kineticGasEvaporation.C | 301 +++-------
.../kineticGasEvaporation.H | 42 +-
.../MassTransferPhaseSystem.C | 521 +++---------------
.../MassTransferPhaseSystem.H | 7 +-
.../MultiComponentPhaseModel.C | 12 +-
.../MultiComponentPhaseModel.H | 3 -
.../phaseModel/phaseModel/phaseModel.C | 9 +-
.../phaseModel/phaseModel/phaseModel.H | 9 -
.../multiphaseSystem/multiphaseSystem.C | 58 +-
.../phasesSystem/phaseSystem/phaseSystem.C | 20 +-
.../phasesSystem/phaseSystem/phaseSystem.H | 2 +-
.../setDeltaT.H | 3 +-
.../radiationModel/radiationModel.C | 14 +
.../radiationModel/radiationModel.H | 7 +
31 files changed, 545 insertions(+), 1159 deletions(-)
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/TEqn.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/TEqn.H
index 805d4ee5a6..6bad14f11e 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/TEqn.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/TEqn.H
@@ -10,8 +10,6 @@
fluid.kappa() + fluid.Cp()*turbulence->mut()/fluid.Prt()
);
- //dimensionedScalar S("S", dimEnergy/dimVolume/dimTime, 1.225e8);
-
fvScalarMatrix TEqn
(
fvm::ddt(rhoCp, T)
@@ -24,6 +22,11 @@
+ fvOptions(rhoCp, T)
);
+ if (T.mesh().time().outputTime())
+ {
+ kappaEff.write();
+ }
+
TEqn.relax();
fvOptions.constrain(TEqn);
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/icoReactingMultiphaseInterFoam.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/icoReactingMultiphaseInterFoam.C
index 52aef44822..5e1e8cbe42 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/icoReactingMultiphaseInterFoam.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/icoReactingMultiphaseInterFoam.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2015 OpenCFD Ltd.
+ \\ / A nd | Copyright (C) 2017 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -106,8 +106,6 @@ int main(int argc, char *argv[])
}
}
-
-
rho = fluid.rho();
runTime.write();
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticle.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticle.C
index 54df98a54e..9f399b4ead 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticle.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticle.C
@@ -38,8 +38,7 @@ Foam::DTRMParticle::DTRMParticle
const scalar I,
const label cellI,
const scalar dA,
- const label reflectedId,
- const scalar Imin,
+ const label transmissiveId,
bool doCellFacePt
)
:
@@ -49,8 +48,7 @@ Foam::DTRMParticle::DTRMParticle
I0_(I),
I_(I),
dA_(dA),
- reflectedId_(reflectedId),
- Imin_(Imin)
+ transmissiveId_(transmissiveId)
{}
@@ -62,8 +60,7 @@ Foam::DTRMParticle::DTRMParticle(const DTRMParticle& p)
I0_(p.I0_),
I_(p.I_),
dA_(p.dA_),
- reflectedId_(p.reflectedId_),
- Imin_(p.Imin_)
+ transmissiveId_(p.transmissiveId_)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@@ -94,43 +91,28 @@ bool Foam::DTRMParticle::move
// Boltzman constant
const scalar sigma = physicoChemical::sigma.value();
+
+ label reflectedZoneId = td.relfectedCells()[celli];
+
if
(
- (!td.relfectedCells()[celli] > 0 && reflectedId_ == 0)
- || reflectedId_ > 0
+ (reflectedZoneId > -1)
+ && (
+ (transmissiveId_ == -1)
+ || (transmissiveId_ != reflectedZoneId)
+ )
)
- {
- scalar a = td.aInterp().interpolate(position(), cell0);
- scalar e = td.eInterp().interpolate(position(), cell0);
- scalar E = td.EInterp().interpolate(position(), cell0);
- scalar T = td.TInterp().interpolate(position(), cell0);
-
- const scalar I1 =
- (
- I_
- + ds*(e*sigma*pow4(T)/mathematical::pi + E)
- ) / (1 + ds*a);
-
- td.Q(cell0) += (I_ - I1)*dA_;
-
- I_ = I1;
-
- if ((I_ <= 0.01*I0_) || (I_ < Imin_))
- {
- break;
- }
- }
- else
{
scalar rho(0);
+
// Create a new reflected particle when the particles is not
// transmissive and larger than an absolute I
- if (reflectedId_ == 0 && I_ > Imin_)
+ if (I_ > 0.01*I0_)
{
vector pDir = dsv/ds;
cellPointWeight cpw(mesh_, position(), celli, face());
- //vector nHat = td.nHatCells()[celli];
+
vector nHat = td.nHatInterp().interpolate(cpw);
nHat /= mag(nHat);
@@ -138,12 +120,26 @@ bool Foam::DTRMParticle::move
// Only new incoming rays
if (cosTheta > SMALL)
{
- vector newDir = td.reflection().R(pDir, nHat);
-
- //scalar theta = acos(-pDir & nHat);
+ vector newDir =
+ td.reflection()
+ [
+ td.relfectedCells()[celli]
+ ].R(pDir, nHat);
// reflectivity
- rho = min(max(td.reflection().rho(cosTheta), 0.0), 0.98);
+ rho =
+ min
+ (
+ max
+ (
+ td.reflection()
+ [
+ td.relfectedCells()[celli]
+ ].rho(cosTheta)
+ , 0.0
+ )
+ , 0.98
+ );
scalar delaM = sqrt(mesh_.cellVolumes()[cell0]);
@@ -155,8 +151,7 @@ bool Foam::DTRMParticle::move
I_*rho,
cell0,
dA_,
- reflectedId_,
- Imin_,
+ transmissiveId_,
true
);
// Add to cloud
@@ -164,7 +159,8 @@ bool Foam::DTRMParticle::move
}
}
- reflectedId_++;
+ // Change transmissiveId of the particle
+ transmissiveId_ = reflectedZoneId;
const point p0 = position();
@@ -186,11 +182,33 @@ bool Foam::DTRMParticle::move
+ ds*(e*sigma*pow4(T)/mathematical::pi + E)
) / (1 + ds*a);
- td.Q(celli) += (Itran - I1)*dA_;
+ td.Q(celli) += (Itran - max(I1, 0.0))*dA_;
+
+ I_ = I1;
+
+ if (I_ <= 0.01*I0_)
+ {
+ break;
+ }
+ }
+ else
+ {
+ scalar a = td.aInterp().interpolate(position(), cell0);
+ scalar e = td.eInterp().interpolate(position(), cell0);
+ scalar E = td.EInterp().interpolate(position(), cell0);
+ scalar T = td.TInterp().interpolate(position(), cell0);
+
+ const scalar I1 =
+ (
+ I_
+ + ds*(e*sigma*pow4(T)/mathematical::pi + E)
+ ) / (1 + ds*a);
+
+ td.Q(cell0) += (I_ - max(I1, 0.0))*dA_;
I_ = I1;
- if (I_ <= 0.01*I0_ || I_ < Imin_)
+ if ((I_ <= 0.01*I0_))
{
break;
}
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticle.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticle.H
index e99925b2cf..c8af27360a 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticle.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticle.H
@@ -43,7 +43,6 @@ SourceFiles
#include "reflectionModel.H"
#include "interpolationCellPoint.H"
-
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
@@ -65,6 +64,7 @@ class DTRMParticle
{
private:
+
// Private data
//- Size in bytes of the fields
@@ -85,11 +85,8 @@ private:
//- Area of radiation
scalar dA_;
- //- Reflected index
- label reflectedId_;
-
- //- Minimum radiation intensity to which the particle is tracked [W/m2]
- scalar Imin_;
+ //- Trasnmissive index
+ label transmissiveId_;
public:
@@ -101,6 +98,8 @@ public:
:
public particle::TrackingData<Cloud<DTRMParticle>>
{
+
+
// Interpolators for continuous phase fields
const interpolationCell<scalar>& aInterp_;
@@ -111,7 +110,8 @@ public:
const interpolationCellPoint<vector>& nHatInterp_;
const labelField& relfectedCells_;
- const reflectionModel& reflection_;
+ //const reflectionModelTable& reflection_;
+ UPtrList<reflectionModel> reflection_;
//- Heat source term
volScalarField& Q_;
@@ -130,7 +130,7 @@ public:
const interpolationCell<scalar>& TInterp,
const interpolationCellPoint<vector>& nHatInterp,
const labelField&,
- const reflectionModel&,
+ const UPtrList<reflectionModel>&,
volScalarField& Q
);
@@ -147,7 +147,7 @@ public:
inline const interpolationCell<scalar>& TInterp() const;
inline const interpolationCellPoint<vector>& nHatInterp() const;
inline const labelField& relfectedCells() const;
- inline const reflectionModel& reflection() const;
+ inline const UPtrList<reflectionModel>& reflection() const;
inline scalar& Q(label celli);
@@ -190,8 +190,7 @@ public:
const scalar I,
const label cellI,
const scalar dA,
- const label reflectedId,
- const scalar Imin,
+ const label transmissiveId,
bool doCellFacePt = true
);
@@ -274,15 +273,6 @@ public:
//- Return access to reflectedId
inline label& reflectedId();
- //- Return access to initial tet face
- //inline label& tetFace0();
-
- //- Return access to initial tet point
- //inline label& tetPt0();
-
- //- Return access to initial proc Id
- //inline label& origProc0();
-
// Tracking
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticleI.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticleI.H
index 5aea5a45a0..52ca2f31cc 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticleI.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/DTRMParticle/DTRMParticleI.H
@@ -34,7 +34,7 @@ inline Foam::DTRMParticle::trackingData::trackingData
const interpolationCell<scalar>& TInterp,
const interpolationCellPoint<vector>& nHatInterp,
const labelField& relfectedCell,
- const reflectionModel& reflection,
+ const UPtrList<reflectionModel>& reflection,
volScalarField& Q
)
:
@@ -92,7 +92,7 @@ Foam::DTRMParticle::trackingData::relfectedCells() const
}
-inline const Foam::reflectionModel&
+inline const Foam::UPtrList<Foam::radiation::reflectionModel>&
Foam::DTRMParticle::trackingData::reflection() const
{
return reflection_;
@@ -147,12 +147,6 @@ inline Foam::point& Foam::DTRMParticle::p0()
}
-inline Foam::label& Foam::DTRMParticle::reflectedId()
-{
- return reflectedId_;
-}
-
-
inline Foam::point& Foam::DTRMParticle::p1()
{
return p1_;
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/Make/options b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/Make/options
index 1a31e07a96..165a136380 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/Make/options
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/Make/options
@@ -1,5 +1,5 @@
EXE_INC = \
- -DFULLDEBUG -g -O0 \
+ -I../phasesSystem/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/laserDTRM.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/laserDTRM.C
index c89e0ad3e6..c4ecf1d56c 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/laserDTRM.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/laserDTRM.C
@@ -124,19 +124,13 @@ Foam::tmp<Foam::volVectorField> Foam::radiation::laserDTRM::nHatfv
const dimensionedScalar deltaN
(
"deltaN",
- 1e-8/pow(average(mesh_.V()), 1.0/3.0)
- );
-
- const dimensionedScalar minAlpha
- (
- "minAlpha", dimless, 1e-3
+ 1e-7/pow(average(mesh_.V()), 1.0/3.0)
);
volVectorField gradAlphaf
(
- "gradAlphaf",
- (alpha2 + minAlpha)*fvc::grad(alpha1)
- - (alpha1 + minAlpha)*fvc::grad(alpha2)
+ alpha2*fvc::grad(alpha1)
+ - alpha1*fvc::grad(alpha2)
);
// Face unit interface normal
@@ -156,6 +150,48 @@ Foam::tmp<Foam::volScalarField>Foam::radiation::laserDTRM::nearInterface
}
+Foam::tmp<Foam::volScalarField> Foam::radiation::laserDTRM::sPhase
+(
+ const volScalarField& alpha1,
+ const volScalarField& alpha2
+) const
+{
+ // positive : the phases have similar orientation other
+ // negative: the phases oppose each
+ return fvc::grad(alpha1) & fvc::grad(alpha2);
+}
+
+
+
+void Foam::radiation::laserDTRM::initialiseReflection()
+{
+ if (found("reflectionModel"))
+ {
+ dictTable modelDicts(lookup("reflectionModel"));
+
+ forAllConstIter(dictTable, modelDicts, iter)
+ {
+ const phasePairKey& key = iter.key();
+
+ reflections_.insert
+ (
+ key,
+ reflectionModel::New
+ (
+ *iter,
+ mesh_
+ )
+ );
+ }
+
+ if (reflections_.size() > 0)
+ {
+ reflectionSwitch_ = true;
+ }
+ }
+}
+
+
void Foam::radiation::laserDTRM::initialise()
{
// Initialise the DTRM particles
@@ -203,19 +239,6 @@ void Foam::radiation::laserDTRM::initialise()
new interpolation2DTable<scalar>(*this)
);
- // Check dimensions ndr and ndTheta
-// if
-// (
-// (powerDistribution_->size() != ndTheta_)
-// || (powerDistribution_().first().second().size() != ndr_)
-// )
-// {
-// FatalErrorInFunction
-// << " The table dimensions should correspond with ndTheta "
-// << " and ndr "
-// << exit(FatalError);
-// }
-
break;
}
case pdUniform:
@@ -235,9 +258,6 @@ void Foam::radiation::laserDTRM::initialise()
p1 = p0;
if (mesh_.nGeometricD() == 3)
{
- //scalar r0 = dr/2.0;
- //scalar r1Max0 = drMax/2.0;
-
for (label ri = 0; ri < ndr_; ri++)
{
scalar r1 = SMALL + dr*ri;
@@ -275,8 +295,6 @@ void Foam::radiation::laserDTRM::initialise()
// calculate target point using new deviation rl
p1 = lPosition + finalPos + (0.5*maxTrackLength_*lDir);
- //scalar p = magSqr(p0 - lPosition);
-
scalar Ip = calculateIp(rP, thetaP);
scalar dAi = (sqr(r2) - sqr(r1))*(theta2 - theta1)/2.0;
@@ -290,7 +308,7 @@ void Foam::radiation::laserDTRM::initialise()
{
// Create a new particle
DTRMParticle* pPtr = new DTRMParticle
- (mesh_, p0, p1, Ip, cellI, dAi, 0 , 0.01*Ip, true);
+ (mesh_, p0, p1, Ip, cellI, dAi, -1, true);
// Add to cloud
DTRMCloud_.addParticle(pPtr);
@@ -338,12 +356,12 @@ Foam::radiation::laserDTRM::laserDTRM(const volScalarField& T)
(
Function1<point>::New("focalLaserPosition", *this)
),
+
laserDirection_
(
Function1<vector>::New("laserDirection", *this)
),
-
focalLaserRadius_(readScalar(lookup("focalLaserRadius"))),
qualityBeamLaser_
(
@@ -354,11 +372,10 @@ Foam::radiation::laserDTRM::laserDTRM(const volScalarField& T)
laserPower_(Function1<scalar>::New("laserPower", *this)),
powerDistribution_(),
- reflection_(reflectionModel::New(*this, mesh_).ptr()),
- reflectionSwitch_(lookupOrDefault("reflection", false)),
- initialPhase_(lookupOrDefault("initialPhase", word::null)),
- alpha1_(lookupOrDefault("alpha1", word::null)),
- alpha2_(lookupOrDefault("alpha2", word::null)),
+ reflectionSwitch_(false),
+
+ alphaCut_( lookupOrDefault<scalar>("alphaCut", 0.5)),
+
Qin_
(
IOobject
@@ -413,7 +430,7 @@ Foam::radiation::laserDTRM::laserDTRM(const volScalarField& T)
),
Q_
(
- IOobject
+ IOobject
(
"Q",
mesh_.time().timeName(),
@@ -425,6 +442,8 @@ Foam::radiation::laserDTRM::laserDTRM(const volScalarField& T)
dimensionedScalar("Q", dimPower/dimVolume, 0.0)
)
{
+ initialiseReflection();
+
initialise();
}
@@ -462,11 +481,10 @@ Foam::radiation::laserDTRM::laserDTRM
laserPower_(Function1<scalar>::New("laserPower", *this)),
powerDistribution_(),
- reflection_(reflectionModel::New(*this, mesh_).ptr()),
- reflectionSwitch_(dict.lookupOrDefault("reflection", false)),
- initialPhase_(lookupOrDefault("initialPhase", word::null)),
- alpha1_(lookupOrDefault("alpha1", word::null)),
- alpha2_(lookupOrDefault("alpha2", word::null)),
+ reflectionSwitch_(false),
+
+ alphaCut_( lookupOrDefault<scalar>("alphaCut", 0.5)),
+
Qin_
(
IOobject
@@ -533,6 +551,7 @@ Foam::radiation::laserDTRM::laserDTRM
dimensionedScalar("Q", dimPower/pow3(dimLength), 0.0)
)
{
+ initialiseReflection();
initialise();
}
@@ -573,25 +592,11 @@ void Foam::radiation::laserDTRM::calculate()
IOobject::NO_WRITE
),
mesh_,
- dimensionedScalar("zero", dimless, 0)
+ dimensionedScalar("zero", dimless, -1)
)
);
volScalarField& reflectingCellsVol = treflectingCells.ref();
- // Reset the fields
- Qin_ == dimensionedScalar("zero", Qin_.dimensions(), 0);
- Q_ == dimensionedScalar("zero", Q_.dimensions(), 0);
-
- a_ = absorptionEmission_->a();
- e_ = absorptionEmission_->e();
- E_ = absorptionEmission_->E();
-
- const interpolationCell<scalar> aInterp(a_);
- const interpolationCell<scalar> eInterp(e_);
- const interpolationCell<scalar> EInterp(E_);
- const interpolationCell<scalar> TInterp(T_);
-
- labelField reflectingCells(mesh_.nCells(), 0);
tmp<volVectorField> tnHat
(
@@ -611,54 +616,91 @@ void Foam::radiation::laserDTRM::calculate()
);
volVectorField& nHat = tnHat.ref();
+ // Reset the fields
+ Qin_ == dimensionedScalar("zero", Qin_.dimensions(), 0);
+ Q_ == dimensionedScalar("zero", Q_.dimensions(), 0);
+
+ a_ = absorptionEmission_->a();
+ e_ = absorptionEmission_->e();
+ E_ = absorptionEmission_->E();
+
+ const interpolationCell<scalar> aInterp(a_);
+ const interpolationCell<scalar> eInterp(e_);
+ const interpolationCell<scalar> EInterp(E_);
+ const interpolationCell<scalar> TInterp(T_);
+
+ labelField reflectingCells(mesh_.nCells(), -1);
+
autoPtr<interpolationCellPoint<vector>> nHatIntrPtr;
+ UPtrList<reflectionModel> reflactionUPtr;
+
if (reflectionSwitch_)
{
- const volScalarField& initialPhase =
- mesh_.lookupObject<volScalarField>(initialPhase_);
+ reflactionUPtr.resize(reflections_.size());
- if (alpha1_ != word::null)
+ label reflectionModelId(0);
+ forAllIter(reflectionModelTable, reflections_, iter1)
{
- const volScalarField& alpha1 =
- mesh_.lookupObject<volScalarField>(alpha1_);
+ reflectionModel& model = iter1()();
- nHat = nHatfv(initialPhase, alpha1);
+ reflactionUPtr.set(reflectionModelId, &model);
- forAll(alpha1, cellI)
- {
- if ((alpha1[cellI] > 0.9) && (mag(nHat[cellI]) > 0.9))
- {
- reflectingCells[cellI] = 1;
- reflectingCellsVol[cellI] = 1.0;
- }
- }
- }
+ const word alpha1Name = "alpha." + iter1.key().first();
+ const word alpha2Name = "alpha." + iter1.key().second();
- if (alpha2_ != word::null)
- {
- const volScalarField& alpha2 =
- mesh_.lookupObject<volScalarField>(alpha2_);
+ const volScalarField& alphaFrom =
+ mesh_.lookupObject<volScalarField>(alpha1Name);
+
+ const volScalarField& alphaTo =
+ mesh_.lookupObject<volScalarField>(alpha2Name);
+
+ const volVectorField nHatPhase(nHatfv(alphaFrom, alphaTo));
- nHat += nHatfv(initialPhase, alpha2);
+ volScalarField gradAlphaf
+ (
+ fvc::grad(alphaFrom)
+ & fvc::grad(alphaTo)
+ );
- forAll(alpha2, cellI)
+
+ const volScalarField nearInterface(pos(alphaTo - alphaCut_));
+
+ const volScalarField mask(nearInterface*gradAlphaf);
+
+ forAll(alphaFrom, cellI)
{
- if ((alpha2[cellI] > 0.9) && (mag(nHat[cellI]) > 0.9))
+ if
+ (
+ nearInterface[cellI]
+ && mag(nHatPhase[cellI]) > 0.99
+ && mask[cellI] < 0
+ )
{
- reflectingCells[cellI] = 1;
- reflectingCellsVol[cellI] = 1.0;
+ reflectingCells[cellI] = reflectionModelId;
+ reflectingCellsVol[cellI] = reflectionModelId;
+ if (!mag(nHat[cellI]) > 0.0)
+ {
+ nHat[cellI] += nHatPhase[cellI];
+ }
}
}
+ reflectionModelId++;
}
+ }
- nHatIntrPtr.reset
- (
- new interpolationCellPoint<vector>(nHat)
- );
+ nHatIntrPtr.reset
+ (
+ new interpolationCellPoint<vector>(nHat)
+ );
+ if (mesh_.time().outputTime())
+ {
+ reflectingCellsVol.write();
+ nHat.write();
}
+
DTRMParticle::trackingData td
(
DTRMCloud_,
@@ -668,11 +710,11 @@ void Foam::radiation::laserDTRM::calculate()
TInterp,
nHatIntrPtr,
reflectingCells,
- reflection_,
+ reflactionUPtr,
Q_
);
- Info << "Move particles..."
+ Info << "Move initial particles..."
<< returnReduce(DTRMCloud_.size(), sumOp<label>()) << endl;
DTRMCloud_.move(td, mesh_.time().deltaTValue());
@@ -706,6 +748,7 @@ void Foam::radiation::laserDTRM::calculate()
if (mesh_.time().outputTime())
{
reflectingCellsVol.write();
+ nHat.write();
}
}
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/laserDTRM.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/laserDTRM.H
index 1d93d7eb04..1f5255b4d2 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/laserDTRM.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/laserDTRM.H
@@ -50,6 +50,7 @@ SourceFiles
#include "Function1.H"
#include "interpolation2DTable.H"
#include "labelField.H"
+#include "phasePairKey.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -78,9 +79,28 @@ public:
pdUniform
};
+
private:
+ // Private types
+
+ typedef
+ HashTable<dictionary, phasePairKey, phasePairKey::hash> dictTable;
+
+
+
+ // Private types
+
+ typedef
+ HashTable
+ <
+ autoPtr<reflectionModel>,
+ phasePairKey,
+ phasePairKey::hash
+ > reflectionModelTable;
+
+
// Private data
@@ -104,6 +124,7 @@ private:
//- Maximum tracking length for particles
scalar maxTrackLength_;
+
// Laser parameters
//- Focal laser position
@@ -138,19 +159,13 @@ private:
// Reflection sub-model
//- Reflection model
- autoPtr<reflectionModel> reflection_;
+ reflectionModelTable reflections_;
//- Reflection switch
bool reflectionSwitch_;
- //- Phase in which the particles are inserted
- word initialPhase_;
-
- //- Phase name for absorbing medium 1
- word alpha1_;
-
- //- Phase name for absorbing medium 2
- word alpha2_;
+ //- Alpha value at which reflection is set
+ scalar alphaCut_;
// Fields
@@ -177,6 +192,9 @@ private:
//- Initialise
void initialise();
+ //- Initialise reflection model
+ void initialiseReflection();
+
//- Calculate Intensity of the laser at p(t, theta) [W/m2]
scalar calculateIp(scalar r, scalar theta);
@@ -196,6 +214,13 @@ private:
const volScalarField& alpha2
) const;
+ //- Indicator if the normal direction for each phase
+ tmp<volScalarField> sPhase
+ (
+ const volScalarField& alpha1,
+ const volScalarField& alpha2
+ ) const;
+
//- Disallow default bitwise copy construct
laserDTRM(const laserDTRM&);
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/localDensityAbsorptionEmission/localDensityAbsorptionEmission.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/localDensityAbsorptionEmission/localDensityAbsorptionEmission.C
index f049c29805..a0c5e4ef00 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/localDensityAbsorptionEmission/localDensityAbsorptionEmission.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/localDensityAbsorptionEmission/localDensityAbsorptionEmission.C
@@ -109,16 +109,19 @@ Foam::radiation::localDensityAbsorptionEmission::aCont(const label bandI) const
false
),
mesh_,
- dimensionedScalar("zero", dimless/dimLength, 0)
+ dimensionedScalar("zero", inv(dimLength), 0)
)
);
- scalarField& a = ta.ref().primitiveFieldRef();
+ volScalarField& a = ta.ref();
forAll(alphaNames_, i)
{
dimensionedScalar aPhase("a", dimless/dimLength, aCoeff_[i]);
- a += max(alpha(alphaNames_[i]), 0.0)*aPhase;
+ a +=
+ max(alpha(alphaNames_[i]), 0.0)
+ *pos(alpha(alphaNames_[i]) - 0.1)
+ *aPhase;
}
return ta;
@@ -142,16 +145,19 @@ Foam::radiation::localDensityAbsorptionEmission::eCont(const label bandI) const
false
),
mesh_,
- dimensionedScalar("zero", dimless/dimLength, 0)
+ dimensionedScalar("zero", inv(dimLength), 0)
)
);
- scalarField& e = te.ref().primitiveFieldRef();
+ volScalarField& e = te.ref();
forAll(alphaNames_, i)
{
dimensionedScalar ePhase("e", dimless/dimLength, eCoeff_[i]);
- e += max(alpha(alphaNames_[i]), 0.0)*ePhase;
+ e +=
+ max(alpha(alphaNames_[i]), 0.0)
+ *pos(alpha(alphaNames_[i]) - 0.1)
+ *ePhase;
}
return te;
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/FresnelLaser/FresnelLaser.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/FresnelLaser/FresnelLaser.C
index 19c9e02249..2d3b77c846 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/FresnelLaser/FresnelLaser.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/FresnelLaser/FresnelLaser.C
@@ -53,8 +53,7 @@ Foam::radiation::FresnelLaser::FresnelLaser
)
:
reflectionModel(dict, mesh),
- coeffsDict_(dict.subDict(typeName + "Coeffs")),
- epsilon_(readScalar(coeffsDict_.lookup("epsilon")))
+ epsilon_(readScalar(dict.lookup("epsilon")))
{}
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/FresnelLaser/FresnelLaser.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/FresnelLaser/FresnelLaser.H
index ad7e18af66..67104ecc31 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/FresnelLaser/FresnelLaser.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/FresnelLaser/FresnelLaser.H
@@ -58,10 +58,6 @@ class FresnelLaser
// Private data
- //- Coefficients dictionary
- dictionary coeffsDict_;
-
-
//- Model constant
scalar epsilon_;
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/reflectionModel/reflectionModelNew.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/reflectionModel/reflectionModelNew.C
index 2f95cb89f4..f94f18d4e0 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/reflectionModel/reflectionModelNew.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/laserDTRM/reflectionModel/reflectionModel/reflectionModelNew.C
@@ -36,7 +36,7 @@ New
const fvMesh& mesh
)
{
- const word modelType(dict.lookup("reflectionModel"));
+ const word modelType(dict.lookup("type"));
Info<< "Selecting reflectionModel " << modelType << endl;
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/InterfaceCompositionModel/InterfaceCompositionModels.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/InterfaceCompositionModel/InterfaceCompositionModels.C
index fe69ba56f7..e653451074 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/InterfaceCompositionModel/InterfaceCompositionModels.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/InterfaceCompositionModel/InterfaceCompositionModels.C
@@ -91,7 +91,6 @@ namespace Foam
using namespace meltingEvaporationModels;
//NOTE: First thermo (from) and second otherThermo (to)
- // in the phaseProperties: (from to to)
// kineticGasEvaporation model definitions
/*
@@ -110,7 +109,7 @@ namespace Foam
);
*/
- // pure from phase to a multi-component to phase
+ // From pure liquid (rhoConst) to a multi-component gas incomp phase
makeInterfaceContSpecieMixtureType
(
kineticGasEvaporation,
@@ -124,22 +123,36 @@ namespace Foam
constIncompressibleGasHThermoPhysics
);
+ // From pure liquid (BoussinesqFluid) to a multi-component gas incomp phase
+ makeInterfaceContSpecieMixtureType
+ (
+ kineticGasEvaporation,
+ heRhoThermo,
+ rhoThermo,
+ pureMixture,
+ BoussinesqFluidEThermoPhysics,
+ heRhoThermo,
+ rhoReactionThermo,
+ multiComponentMixture,
+ constIncompressibleGasHThermoPhysics
+ );
- // pure from phase and pure to phase with incompressible gas
+
+ // From pure liquid (rhoConst) to pure gas (incompressible ideal gas)
makeInterfacePureType
(
kineticGasEvaporation,
heRhoThermo,
rhoThermo,
pureMixture,
- constIncompressibleGasHThermoPhysics,
+ constFluidHThermoPhysics,
heRhoThermo,
rhoThermo,
pureMixture,
- constFluidHThermoPhysics
+ constIncompressibleGasHThermoPhysics
);
- // pure from phase and pure to phase with rhoConst gas
+ // From pure liquid (const rho) to pure gas (rhoConst) gas
makeInterfacePureType
(
kineticGasEvaporation,
@@ -154,9 +167,38 @@ namespace Foam
);
+ // From pure liquid (Boussinesq) to pure gas (incompressible ideal gas)
+ makeInterfacePureType
+ (
+ kineticGasEvaporation,
+ heRhoThermo,
+ rhoThermo,
+ pureMixture,
+ BoussinesqFluidEThermoPhysics,
+ heRhoThermo,
+ rhoThermo,
+ pureMixture,
+ constIncompressibleGasHThermoPhysics
+ );
+
+ // From pure liquid (Boussinesq) to pure gas (rho const)
+ makeInterfacePureType
+ (
+ kineticGasEvaporation,
+ heRhoThermo,
+ rhoThermo,
+ pureMixture,
+ BoussinesqFluidEThermoPhysics,
+ heRhoThermo,
+ rhoThermo,
+ pureMixture,
+ constFluidHThermoPhysics
+ );
+
+
// Lee model definitions
- // pure from phase and a pure to phase
+ // From pure phase (rho const) to phase (rho const)
makeInterfacePureType
(
Lee,
@@ -170,6 +212,7 @@ namespace Foam
constFluidHThermoPhysics
);
+ // From pure phase (rho const) to phase (Boussinesq)
makeInterfacePureType
(
Lee,
@@ -184,6 +227,7 @@ namespace Foam
);
+ // From pure phase (solidThermo) to phase (Boussinesq)
makeInterfacePureType
(
Lee,
@@ -197,19 +241,22 @@ namespace Foam
BoussinesqFluidEThermoPhysics
);
+ // From pure phase (solidThermo) to phase (rho const)
makeInterfacePureType
(
Lee,
- heRhoThermo,
- rhoThermo,
+ heSolidThermo,
+ solidThermo,
pureMixture,
- constFluidHThermoPhysics,
+ hConstSolidThermoPhysics,
heRhoThermo,
rhoThermo,
pureMixture,
- constIncompressibleGasHThermoPhysics
+ constFluidHThermoPhysics
);
+
+ // From pure phase (const rho) to multi phase (incomp ideal gas)
makeInterfaceContSpecieMixtureType
(
Lee,
@@ -223,22 +270,6 @@ namespace Foam
constIncompressibleGasHThermoPhysics
);
-/*
- makeInterfaceDispSpecieMixtureType
- (
- Lee,
- heRhoThermo,
- rhoReactionThermo,
- multiComponentMixture,
- constIncompressibleGasHThermoPhysics,
- heRhoThermo,
- rhoThermo,
- pureMixture,
- constFluidHThermoPhysics
- );
-*/
- // Lee model definitions
-
// pure from phase and a pure to phase
/*
makeInterfacePureType
@@ -256,6 +287,7 @@ namespace Foam
*/
+ // From pure phase (Boussinesq) to phase (solidThermo)
makeInterfacePureType
(
Lee,
@@ -269,65 +301,19 @@ namespace Foam
hConstSolidThermoPhysics
);
-/*
- // saturatedEvaporation model definitions
-
- // multi-component from phase and a pure to phase
- makeInterfaceDispSpecieMixtureType
- (
- saturatedEvaporation,
- heRhoThermo,
- rhoReactionThermo,
- multiComponentMixture,
- constIncompressibleGasHThermoPhysics,
- heRhoThermo,
- rhoThermo,
- pureMixture,
- constFluidHThermoPhysics
- );
-
- // pure from phase and a multi-component to phase
- makeInterfaceContSpecieMixtureType
- (
- saturatedEvaporation,
- heRhoThermo,
- rhoThermo,
- pureMixture,
- constFluidHThermoPhysics,
- heRhoThermo,
- rhoReactionThermo,
- multiComponentMixture,
- constIncompressibleGasHThermoPhysics
- );
-
- // multi-component from phase and a multi-componen to phase
- makeSpecieInterfaceSpecieMixtures
- (
- saturatedEvaporation,
- heRhoThermo,
- rhoReactionThermo,
- multiComponentMixture,
- constIncompressibleGasHThermoPhysics,
- heRhoThermo,
- rhoReactionThermo,
- multiComponentMixture,
- constIncompressibleGasHThermoPhysics
- );
-
- // pure from phase and pure to phase
+ // From pure phase (rho const) to phase (solidThermo)
makeInterfacePureType
(
- saturatedEvaporation,
+ Lee,
heRhoThermo,
rhoThermo,
pureMixture,
- constIncompressibleGasHThermoPhysics,
- heRhoThermo,
- rhoThermo,
+ constFluidHThermoPhysics,
+ heSolidThermo,
+ solidThermo,
pureMixture,
- constIncompressibleGasHThermoPhysics
+ hConstSolidThermoPhysics
);
-*/
}
// ************************************************************************* //
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/Lee/Lee.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/Lee/Lee.C
index 0c2e3f93f4..01dc7cf31a 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/Lee/Lee.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/Lee/Lee.C
@@ -49,7 +49,7 @@ Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>
template<class Thermo, class OtherThermo>
Foam::tmp<Foam::volScalarField>
Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>
-::Kexp(label variable, const volScalarField& refValue) const
+::Kexp(label variable, const volScalarField& refValue)
{
if (this->modelVariable_ == variable)
{
@@ -65,9 +65,9 @@ Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>
C_
* from
* this->pair().from().rho()
- * (refValue - Tactivate_)
+ * (refValue.oldTime() - Tactivate_)
* pos(from - alphaMin_)
- * pos(refValue - Tactivate_)/Tactivate_
+ * pos(refValue.oldTime() - Tactivate_)/Tactivate_
);
}
else
@@ -78,8 +78,8 @@ Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>
* from
* this->pair().from().rho()
* pos(from - alphaMin_)
- * (Tactivate_ - refValue)
- * pos(Tactivate_ - refValue)/Tactivate_
+ * (Tactivate_ - refValue.oldTime())
+ * pos(Tactivate_ - refValue.oldTime())/Tactivate_
);
}
@@ -91,33 +91,6 @@ Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>
}
-template<class Thermo, class OtherThermo>
-Foam::tmp<Foam::volScalarField>
-Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>
-::Kimp(label variable, const volScalarField& refValue) const
-{
- if (this->modelVariable_ == variable)
- {
- volScalarField limitedDispersed
- (
- min(max(this->pair().from(), scalar(0)), scalar(1))
- );
-
- return
- (
- C_
- *limitedDispersed
- *this->pair().from().rho()
- *pos(Tactivate_ - refValue.oldTime())
- );
- }
- else
- {
- return tmp<volScalarField> ();
- }
-}
-
-
template<class Thermo, class OtherThermo>
const Foam::dimensionedScalar&
Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>
@@ -127,13 +100,4 @@ Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>
}
-template<class Thermo, class OtherThermo>
-Foam::label
-Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>
-::dSdVariable()
-{
- return label(-1);
-}
-
-
// ************************************************************************* //
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/Lee/Lee.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/Lee/Lee.H
index 734470793b..5540c7f454 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/Lee/Lee.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/Lee/Lee.H
@@ -105,22 +105,11 @@ public:
(
label variable,
const volScalarField& field
- ) const;
-
- //- Semi implicit species mass transfer coefficient
- virtual tmp<volScalarField> Kimp
- (
- label variable,
- const volScalarField& field
- ) const;
+ );
//- Return T transition between phases
virtual const dimensionedScalar& Tactivate() const;
- //- Derivate sign of the source term:
- // d(Sk-i)/d(variable) > 0 or d(Si-k)/d(variable) < 0
- virtual label dSdVariable();
-
};
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/interfaceCompositionModel/interfaceCompositionModel.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/interfaceCompositionModel/interfaceCompositionModel.C
index 87cb330276..3851cc7f38 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/interfaceCompositionModel/interfaceCompositionModel.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/interfaceCompositionModel/interfaceCompositionModel.C
@@ -70,16 +70,10 @@ Foam::interfaceCompositionModel::interfaceCompositionModel
)
:
modelVariable_(modelVariableNames.read(dict.lookup("variable"))),
- semiImplicit_(dict.lookupOrDefault<bool>("semiImplicit", false)),
pair_(pair),
speciesName_(dict.lookupOrDefault<word>("species", "none")),
mesh_(pair_.from().mesh())
-{/*
- if (dict.found("species"))
- {
- speciesNames_ = ;
- }*/
-}
+{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
@@ -107,23 +101,4 @@ const Foam::word Foam::interfaceCompositionModel::variable() const
return modelVariableNames[modelVariable_];
}
-
-bool Foam::interfaceCompositionModel::semiImplicit() const
-{
- return semiImplicit_;
-}
-
-
-Foam::tmp<Foam::volScalarField> Foam::interfaceCompositionModel::KexpEnergy
-(
- label modelVariable,
- const volScalarField& field
-)
-const
-{
- NotImplemented;
- return tmp<volScalarField>();
-}
-
-
// ************************************************************************* //
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/interfaceCompositionModel/interfaceCompositionModel.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/interfaceCompositionModel/interfaceCompositionModel.H
index ee28fe69ac..303c6c0571 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/interfaceCompositionModel/interfaceCompositionModel.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/interfaceCompositionModel/interfaceCompositionModel.H
@@ -73,9 +73,6 @@ public:
//- Enumeration for model variables
modelVariable modelVariable_;
- //- Semi-Implicit or Explicit model
- bool semiImplicit_;
-
protected:
@@ -176,22 +173,7 @@ public:
(
label modelVariable,
const volScalarField& field
- ) const = 0;
-
- //- Implicit mass transfer coefficient of the shape
- // Kimp*(variable - refValue)
- virtual tmp<volScalarField> Kimp
- (
- label modelVariable,
- const volScalarField& field
- ) const = 0;
-
- //- Explicit mass transfer for energy
- virtual tmp<volScalarField> KexpEnergy
- (
- label modelVariable,
- const volScalarField& field
- ) const;
+ ) = 0;
//- Reference value
virtual const dimensionedScalar& Tactivate() const = 0;
@@ -199,16 +181,6 @@ public:
//- Returns the variable on which the model is based
const word variable() const;
- //- Returns type of model. false is fully Explicit
- // true is treated semiImplicit as :
- // Sk-i = Kimp()(variable - Tactivate) for positive k to i
- // Si-k = Kimp()(Tactivate - variable) for positive i to k
- bool semiImplicit() const;
-
- //- Derivate sign of the source term:
- // d(Sk-i)/d(variable) > 0 or d(Si-k)/d(variable) < 0
- virtual label dSdVariable() = 0;
-
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.C
index 76c127a393..cc558b0100 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.C
@@ -26,6 +26,7 @@ License
#include "kineticGasEvaporation.H"
#include "constants.H"
#include "fvcGrad.H"
+#include "fvcSnGrad.H"
#include "fvcDiv.H"
#include "surfaceInterpolate.H"
#include "fvcReconstruct.H"
@@ -59,13 +60,6 @@ Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
dimMass/dimMoles,
dict.lookupOrDefault<scalar>("Mv", -1)
),
- saturationPressureModel_
- (
- saturationPressureModel::New
- (
- dict.subDict("saturationPressure")
- )
- ),
phi0_
(
IOobject
@@ -73,35 +67,21 @@ Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
"phi0",
this->mesh_.time().timeName(),
this->mesh_,
- IOobject::READ_IF_PRESENT,
- IOobject::AUTO_WRITE
- ),
- this->mesh_,
- dimensionedScalar("zero", dimMass/dimTime/dimVolume, 0),
- zeroGradientFvPatchScalarField::typeName
- ),
- sPhi_
- (
- IOobject
- (
- "sPhi",
- this->mesh_.time().timeName(),
- this->mesh_,
- IOobject::READ_IF_PRESENT,
+ IOobject::NO_READ,
IOobject::AUTO_WRITE
),
this->mesh_,
dimensionedScalar("zero", dimMass/dimTime/dimVolume, 0),
zeroGradientFvPatchScalarField::typeName
),
- D_
- (
- "D",
- dimArea/dimTime,
- dict.lookupOrDefault<scalar>("D", 1)
- )
+ alphaMax_(dict.lookupOrDefault<scalar>("alphaMax", 0.03)),
+ alphaMin_(dict.lookupOrDefault<scalar>("alphaMin", 0.01))
{
- word speciesName = this->transferSpecie();
+ word fullSpeciesName = this->transferSpecie();
+
+ const label tempOpen(fullSpeciesName.find('.'));
+
+ const word speciesName(fullSpeciesName(0, tempOpen));
// Get the continuous thermo
const typename OtherThermo::thermoType& localThermo =
@@ -136,26 +116,18 @@ Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
template<class Thermo, class OtherThermo>
Foam::tmp<Foam::volScalarField>
Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
-::Kexp(label variable, const volScalarField& field) const
+::Kexp(label variable, const volScalarField& field)
{
if (this->modelVariable_ == variable)
{
- volScalarField limitedContinous
- (
- min(max(this->pair().to(), scalar(0)), scalar(1))
- );
+ const volScalarField& to = this->pair().to();
- volScalarField limitedDispersed
- (
- min(max(this->pair().from(), scalar(0)), scalar(1))
- );
+ const volScalarField& from = this->pair().from();
const fvMesh& mesh = this->mesh_;
-
- const volScalarField& T = mesh.lookupObject<volScalarField>("T");
- //const volScalarField& p = mesh.lookupObject<volScalarField>("p");
+ const volScalarField& T = mesh.lookupObject<volScalarField>("T").oldTime();
const dimensionedScalar HerztKnudsConst
(
@@ -168,68 +140,70 @@ Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
/pow3(Tactivate_)
)
);
-/*
- volScalarField pSat
- (
- "pSat",
- saturationPressureModel_->pSat(T)
- );
-*/
-/*
- volScalarField rhoMean
- (
- this->pair().to().rho()*this->pair().from().rho()
- / (this->pair().from().rho() - this->pair().to().rho())
- );
-*/
- word species(this->transferSpecie());
- tmp<volScalarField> L = mag(this->L(species, field));
+ word fullSpeciesName = this->transferSpecie();
+
+ const label tempOpen(fullSpeciesName.find('.'));
+
+ const word speciesName(fullSpeciesName(0, tempOpen));
+
+ tmp<volScalarField> L = this->L(speciesName, field);
+
+ const volVectorField gradFrom(fvc::grad(from));
+ const volVectorField gradTo(fvc::grad(to));
- DebugVar(Mv_);
- //dimensionedScalar psat("psat", dimPressure, 1e5);
- //scalarField maxAreaDen(sqrt(3.0)*pow(mesh.V(), -1.0/3.0));
- //scalarField avergageDelta(pow(mesh.V(), 1.0/3.0));
-/*
- volScalarField areaDensity
- (
- IOobject
- (
- "areaDensity",
- mesh.time().timeName(),
- mesh,
- IOobject::NO_READ,
- IOobject::NO_WRITE
- ),
- mesh,
- dimensionedScalar("zero", dimArea/dimVolume, 0.0),
- zeroGradientFvPatchScalarField::typeName
- );
-*/
volScalarField areaDensity
(
"areaDensity",
- mag(fvc::grad(limitedDispersed))
+ sqrt(mag(gradFrom)*mag(gradTo))
);
-/*
- volScalarField tempInterFace
+
+ volScalarField gradAlphaf
(
- "tempInterFace",
- pos(0.55 - limitedDispersed) * pos(limitedDispersed - 0.45)
- * pos(limitedContinous - 0.45) * pos(0.55 - limitedContinous)
- * T
+ gradFrom
+ & gradTo
);
-*/
- dimensionedScalar dMgc = max(areaDensity);
+// const volScalarField trigger
+// (
+// pos(alphaMax_ - from)
+// );
+
+ // const volScalarField mask(nearInterface*gradAlphaf);
- volScalarField Tave
+// const volScalarField Tmask
+// (
+// "Tmask",
+// pos(alphaMax_ - from)*pos(from - alphaMin_)
+// );
+
+ volScalarField Tmask
(
- "Tave", T*pos(areaDensity - 0.1*dMgc)
+ "Tmask",
+ neg(from + 1.0)
);
- const volScalarField massFluxEvap
+ Tmask *= 0.0;
+
+ forAll (Tmask, celli)
+ {
+ if (gradAlphaf[celli] < 0)
+ {
+ if (from[celli] > alphaMin_ && from[celli] < alphaMax_)
+ {
+ scalar alphaRes = 1.0 - from[celli] - to[celli];
+ if (alphaRes < 0.01)
+ {
+ Tmask[celli] = 1.0;
+ }
+ }
+ }
+ }
+
+ const volScalarField Tave(T*Tmask);
+
+ volScalarField massFluxEvap
(
"massFluxEvap",
2*C_/(2 - C_)
@@ -243,140 +217,29 @@ Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
)
);
- //scalar alphaC(0.9);
- //scalar alphaCMax(0.9);
-
-
- //.weightedAverage(mesh.V());
-/*
- volScalarField nearInterFace
- (
- "nearInterFace",
- pos(areaDensity - 0.1*dMgc)
- * pos(alphaC - limitedDispersed)
- * pos(limitedDispersed - (1 - alphaC))
- * pos(limitedContinous - (1 - alphaC))
- * pos(alphaC - limitedContinous)
- );
-
- volScalarField farDispersedInterFace
- (
- pos(limitedDispersed - alphaC)
- );
-
- volScalarField farContinousInterFace
- (
- pos(limitedContinous - alphaC)
- );
-
- volScalarField farInterFace
- (
- "farInterFace", farDispersedInterFace + farContinousInterFace
- );
-*/
- dimensionedScalar mIntDotPhi0
- (
- "mIntDotPhi0",
- gSum((pos(areaDensity - 0.1*dMgc)*areaDensity*mesh.V())())
- /gSum
- (
- (pos(areaDensity - 0.1*dMgc)*limitedDispersed*areaDensity*mesh.V())()
- )
- );
-
- DebugVar(mIntDotPhi0);
-
- // Local density rate (kg/m3/s)
- phi0_ = massFluxEvap*areaDensity; //mIntDotPhi0*limitedDispersed
-
-
- dimensionedScalar mIntDot("mIntDot", gSum((phi0_*mesh.V())()));
-
- DebugVar(mIntDot);
-
- volScalarField sPhi("sPhi", phi0_);
-
- phi0_ = 0.0*phi0_;
-
- fvScalarMatrix sPhiEq
- (
- fvm::ddt(sPhi)
- - fvm::laplacian(D_, sPhi)
- );
-
- sPhiEq.relax();
- sPhiEq.solve();
-
- // Liquid normalization
+ // Liquid normalization
const dimensionedScalar Nl
(
- gSum((sPhi*mesh.V())())
+ gSum((areaDensity*mesh.V())())
/(
gSum
(
- (sPhi*mesh.V()*limitedDispersed)()
+ ((areaDensity*from)*mesh.V())()
)
- + dimensionedScalar("SMALL", dimless, VSMALL)
+ + dimensionedScalar("SMALL", inv(dimless), VSMALL)
)
);
- // Vapour normalization
- const dimensionedScalar Nv
- (
- gSum((sPhi*mesh.V())())
- /(
- gSum
- (
- (sPhi*mesh.V()*limitedContinous)()
- )
- + dimensionedScalar("SMALL", dimless, VSMALL)
- )
- );
-
- // Spread density source
- sPhi_ =
- 0.5*
- (
- limitedContinous*Nv*sPhi
- + limitedDispersed*Nl*sPhi
- );
-
- dimensionedScalar msPhiIntDot("msPhiIntDot", gSum((sPhi_*mesh.V())()));
-
- DebugVar(msPhiIntDot);
+ // Local density rate (kg/m3/s)
+ phi0_ = massFluxEvap*Nl*areaDensity*from; //
if (this->pair().from().mesh().time().outputTime())
{
- massFluxEvap.write();
areaDensity.write();
- Tave.write();
+ Tmask.write();
}
- return
- (
- 0.5*
- (
- limitedContinous*Nv*sPhi
- + limitedDispersed*Nl*sPhi
- )
- );
- }
- else
- {
- return tmp<volScalarField> ();
- }
-}
-
-
-template<class Thermo, class OtherThermo>
-Foam::tmp<Foam::volScalarField>
-Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
-::Kimp(label variable, const volScalarField& field) const
-{
- if (this->modelVariable_ == variable)
- {
- return tmp<volScalarField> ();
-
+ return (1.0*phi0_);
}
else
{
@@ -394,26 +257,4 @@ Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
}
-template<class Thermo, class OtherThermo>
-Foam::tmp<Foam::volScalarField>
-Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
-::KexpEnergy
-(
- label modelVariable,
- const volScalarField& field
-)
-const
-{
- return dimensionedScalar("one", dimless, 1.0)*phi0_;
-}
-
-
-template<class Thermo, class OtherThermo>
-Foam::label
-Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
-::dSdVariable()
-{
- return label(1);
-}
-
// ************************************************************************* //
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.H
index 39ffaf1f3e..b44edbbf6b 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.H
@@ -37,7 +37,6 @@ SourceFiles
#include "InterfaceCompositionModel.H"
-#include "saturationPressureModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *//
@@ -69,23 +68,14 @@ class kineticGasEvaporation
//- Molar weight of the vapour in the continous phase
dimensionedScalar Mv_;
- //- Pointer to saturation pressure model
- autoPtr<saturationPressureModel> saturationPressureModel_;
-
//- Source
- mutable volScalarField phi0_;
-
- //- Spread source
- mutable volScalarField sPhi_;
+ volScalarField phi0_;
- //- Diffusivity for source term
- dimensionedScalar D_;
+ //- Alpha maximum value for the mass transfer
+ scalar alphaMax_;
- // Private Member Functions
-
- //- Constant of propotionality between partial pressure and mass
- // fraction
- tmp<volScalarField> wRatio() const;
+ //- Alpha minumum value for the mass transfer
+ scalar alphaMin_;
public:
@@ -116,31 +106,11 @@ public:
(
label variable,
const volScalarField& field
- ) const;
-
- //- Semi implicit species mass transfer coefficient
- virtual tmp<volScalarField> Kimp
- (
- label variable,
- const volScalarField& field
- ) const;
+ );
//- Return Tref
virtual const dimensionedScalar& Tactivate() const;
-
- //- Explicit mass transfer for energy
- virtual tmp<volScalarField> KexpEnergy
- (
- label modelVariable,
- const volScalarField& field
- ) const;
-
-
- //- Derivate sign of the source term:
- // d(Sk-i)/d(variable) > 0 or d(Si-k)/d(variable) < 0
- virtual label dSdVariable();
-
};
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/MassTransferPhaseSystem/MassTransferPhaseSystem.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/MassTransferPhaseSystem/MassTransferPhaseSystem.C
index 7d616bfec6..073f4b5e36 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/MassTransferPhaseSystem/MassTransferPhaseSystem.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/MassTransferPhaseSystem/MassTransferPhaseSystem.C
@@ -119,6 +119,7 @@ template<class BasePhaseSystem>
Foam::tmp<Foam::volScalarField>
Foam::MassTransferPhaseSystem<BasePhaseSystem>::calculateL
(
+ const volScalarField& dmdtNetki,
const phasePairKey& keyik,
const phasePairKey& keyki,
const volScalarField& T
@@ -153,11 +154,10 @@ Foam::MassTransferPhaseSystem<BasePhaseSystem>::calculateL
const word species(speciesName(0, tempOpen));
- L = -interfacePtr->L(species, T);
-
- return tL;
+ L -= neg(dmdtNetki)*interfacePtr->L(species, T);
}
- else if (massTransferModels_.found(keyki))
+
+ if (massTransferModels_.found(keyki))
{
const autoPtr<interfaceCompositionModel>& interfacePtr =
massTransferModels_[keyki];
@@ -168,9 +168,7 @@ Foam::MassTransferPhaseSystem<BasePhaseSystem>::calculateL
const word species(speciesName(0, tempOpen));
- L = interfacePtr->L(species, T);
-
- return tL;
+ L += pos(dmdtNetki)*interfacePtr->L(species, T);
}
return tL;
@@ -399,10 +397,6 @@ Foam::MassTransferPhaseSystem<BasePhaseSystem>::heatTransfer
const volScalarField& T
)
{
-// tmp<fvScalarMatrix> tEqnPtr
-// (
-// new fvScalarMatrix(T, dimVolume*dimTemperature/dimTime)
-// );
tmp<fvScalarMatrix> tEqnPtr
(
new fvScalarMatrix(T, dimEnergy/dimTime)
@@ -413,10 +407,6 @@ Foam::MassTransferPhaseSystem<BasePhaseSystem>::heatTransfer
forAllIter(phaseSystem::phaseModelTable, this->phaseModels_, iteri)
{
phaseModel& phasei = iteri()();
-// const volScalarField invRhoCpi
-// (
-// 1.0/phasei.rho()/phasei.Cp()
-// );
phaseSystem::phaseModelTable::iterator iterk = iteri;
iterk++;
@@ -431,63 +421,12 @@ Foam::MassTransferPhaseSystem<BasePhaseSystem>::heatTransfer
{
phaseModel& phasek = iterk()();
- //const volScalarField invRhoCpk(1.0/phasek.rho()/phasek.Cp());
-
// Phase i to phase k
const phasePairKey keyik(phasei.name(), phasek.name(), true);
// Phase k to phase i
const phasePairKey keyki(phasek.name(), phasei.name(), true);
-
- // mass transfer models coeffs*(T - Tref) for phase i
-// tmp<volScalarField> tcoeffsi
-// (
-// new volScalarField
-// (
-// IOobject
-// (
-// "tcoeffsi",
-// this->mesh().time().timeName(),
-// this->mesh(),
-// IOobject::NO_READ,
-// IOobject::NO_WRITE
-// ),
-// this->mesh(),
-// dimensionedScalar
-// (
-// "zero",
-// dimDensity/dimTemperature/dimTime,
-// 0
-// )
-// )
-// );
-// volScalarField& coeffsi = tcoeffsi.ref();
-
- // mass transfer models coeffs*(T - Tref) for phase k
-// tmp<volScalarField> tcoeffsk
-// (
-// new volScalarField
-// (
-// IOobject
-// (
-// "tcoeffsk",
-// this->mesh().time().timeName(),
-// this->mesh(),
-// IOobject::NO_READ,
-// IOobject::NO_WRITE
-// ),
-// this->mesh(),
-// dimensionedScalar
-// (
-// "zero",
-// dimDensity/dimTemperature/dimTime,
-// 0
-// )
-// )
-// );
-// volScalarField& coeffsk = tcoeffsk.ref();
-
// Net mass transfer from k to i phase
tmp<volScalarField> tdmdtNetki
(
@@ -518,284 +457,61 @@ Foam::MassTransferPhaseSystem<BasePhaseSystem>::heatTransfer
autoPtr<interfaceCompositionModel>& interfacePtr =
massTransferModels_[keyik];
- //word speciesName = interfacePtr->transferSpecie();
-
- //dSikdVar = interfacePtr->dSdVariable();
-
- // Look for mass transfer i to k specie driven by Y and p
-// if (dmdtYi_.found(keyik) || speciesName != "none")
-// {
-// Pair<word> YSpecie
-// (
-// interfaceCompositionModel::modelVariableNames
-// [
-// interfaceCompositionModel::Y
-// ],
-// speciesName
-// );
-//
-// if (dmdtYi_[keyik].found(YSpecie))
-// {
-// dmdtYki -= *dmdtYi_[keyik][YSpecie];
-// }
-//
-// Pair<word> PSpecie
-// (
-// interfaceCompositionModel::modelVariableNames
-// [
-// interfaceCompositionModel::P
-// ],
-// speciesName
-// );
-//
-// if (dmdtYi_[keyik].found(PSpecie))
-// {
-// dmdtYki -= *dmdtYi_[keyik][PSpecie];
-// }
-// }
-
-
-// volScalarField& dmdtYiT =
-// dmdtYi
-// (
-// keyik,
-// interfaceCompositionModel::modelVariableNames
-// [
-// interfaceCompositionModel::T
-// ],
-// speciesName
-// );
-
- //dimensionedScalar Tref(interfacePtr->Tactivate());
-
- // Fill dmdt for alpha's
- //volScalarField& dmdtik = *dmdt_[keyik];
-
- if (!interfacePtr->semiImplicit())
- {
- // Explicit temperature mass transfer rate
- tmp<volScalarField> Kexp =
- interfacePtr->Kexp
- (
- interfaceCompositionModel::T,
- T
- );
-
- //if (Kexp.valid())
- {
- Info << "Explicit temperature mass transfer.." << endl;
- Info << "keyik :" << keyik << endl;
- // Add explicit T based to all the other explixit terms
- dmdtNetki -= Kexp.ref();
- //dmdtYki -= *dmdt_[keyik];
- *dmdt_[keyik] = Kexp.ref();
-
-
-
- // Add explicit source to dmdt_
- //dmdtik = Kexp.ref();
- // Add explicit source to dmdtYi to be used in YEq's
- //dmdtYiT = Kexp.ref();
-
-// Info<< " max(dmdtik) "
-// << max(dmdtik.internalField()) << endl;
-// Info<< "temperature based Mass rate [kg/s]: " << keyik
-// << gSum((dmdtik*this->mesh().V())()) << endl;
- }
- }
- else if (interfacePtr->semiImplicit())
+ // Explicit temperature mass transfer rate
+ tmp<volScalarField> Kexp =
+ interfacePtr->Kexp
+ (
+ interfaceCompositionModel::T,
+ T
+ );
+
+ if (Kexp.valid())
{
- // Implicit temperature based mass transfer type :
- // Kimp*(T - Tref)
- // NOTE: dmdtYiT for species and dmdt for alpha's are
- // lagged in T
- tmp<volScalarField> Kimp =
- interfacePtr->Kimp
- (
- interfaceCompositionModel::T,
- T
- );
-
- if (Kimp.valid())
- {
- Info << "Implicit temperature mass transfer.." << endl;
- Info << "keyik :" << keyik << endl;
- // Fill alpha mass transfer
- //dmdtik = Kimp.ref()*mag(T.oldTime() - Tref);
- // Fill species mass transfer
- //dmdtYiT = Kimp.ref()*mag(T.oldTime() - Tref);
- //coeffsi = Kimp.ref();
- ///Trefi = Tref;
-
-// Info<< " max(dmdtik) "
-// << max(dmdtik.internalField()) << endl;
-// Info<< "temperature based Mass rate [kg/s]: " << keyik << " "
-// << gSum((dmdtik*this->mesh().V())()) << endl;
- }
+ Info << "Explicit temperature mass transfer.." << endl;
+ Info << "keyik :" << keyik << endl;
+
+ // Add explicit T based to all the other explixit terms
+ dmdtNetki -= Kexp.ref();
+ *dmdt_[keyik] = Kexp.ref();
+
}
}
- // derivate of source on variable (dSki/dVar) from k to i
- //label dSkidVar(0);
-
// Looking for mass transfer in the other direction (k to i)
if (massTransferModels_.found(keyki))
{
autoPtr<interfaceCompositionModel>& interfacePtr =
massTransferModels_[keyki];
- //word speciesName = interfacePtr->transferSpecie();
-
- //dSkidVar = interfacePtr->dSdVariable();
-
- // Look for mass transfer from k to i specie and pressure driven
-// if (dmdtYi_.found(keyki) && speciesName != "none")
-// {
-// Pair<word> YSpecie
-// (
-// interfaceCompositionModel::modelVariableNames
-// [
-// interfaceCompositionModel::Y
-// ],
-// speciesName
-// );
-//
-// if (dmdtYi_[keyki].found(YSpecie))
-// {
-// dmdtYki += *dmdtYi_[keyki][YSpecie];
-// }
-//
-// Pair<word> PSpecie
-// (
-// interfaceCompositionModel::modelVariableNames
-// [
-// interfaceCompositionModel::P
-// ],
-// speciesName
-// );
-//
-// if (dmdtYi_[keyki].found(PSpecie))
-// {
-// dmdtYki += *dmdtYi_[keyki][PSpecie];
-// }
-// }
-
- //dimensionedScalar Tref(interfacePtr->Tactivate());
- //volScalarField& dmdtki = *dmdt_(keyki);
-
-// volScalarField& dmdtYiT =
-// dmdtYi
-// (
-// keyki,
-// interfaceCompositionModel::modelVariableNames
-// [
-// interfaceCompositionModel::T
-// ],
-// speciesName
-// );
-
- if (!interfacePtr->semiImplicit())
- {
- // Explicit temperature mass transfer rate
- const tmp<volScalarField> Kexp =
- interfacePtr->Kexp
- (
- interfaceCompositionModel::T,
- T
- );
-
- if (Kexp.valid())
- {
- Info << "Explicit temperature mass transfer.." << endl;
- Info << "keyki :" << keyki << endl;
-
- dmdtNetki += Kexp.ref();
- //dmdtYki += *dmdt_[keyki];
- *dmdt_[keyki] = Kexp.ref();
-
-
- //dmdtki = Kexp.ref();
- //dmdtYiT = Kexp.ref();
-
-// Info<< " max(dmdtki) "
-// << max(dmdtki.internalField()) << endl;
- Info<< " max(dmdtYki) "
- << max(dmdtNetki.internalField()) << endl;
- //Info<< " max(dmdtYiT) "
- // << max(dmdtYiT.internalField()) << endl;
-
-// Info<< "temperature based Mass rate [kg/s]: "
-// << gSum((dmdtki*this->mesh().V())()) << endl;
- }
- }
- else if (interfacePtr->semiImplicit())
+ // Explicit temperature mass transfer rate
+ const tmp<volScalarField> Kexp =
+ interfacePtr->Kexp
+ (
+ interfaceCompositionModel::T,
+ T
+ );
+
+ if (Kexp.valid())
{
- // Implicit temperature based mass transfer
- // NOTE: dmdtYiT and dmdt are lagged in T
- tmp<volScalarField> Kimp =
- interfacePtr->Kimp
- (
- interfaceCompositionModel::T,
- T
- );
-
- if (Kimp.valid())
- {
- Info << "Implicit temperature mass transfer.." << endl;
-
- //dmdtki = Kimp.ref()*mag(T.oldTime() - Tref);
- //dmdtYiT = Kimp.ref()*mag(T.oldTime() - Tref);
- //coeffsk = Kimp.ref();
- //Trefk = Tref;
-
-// Info<< " max(dmdtki) "
-// << max(dmdtki.internalField()) << endl;
-// Info<< "temperature based Mass rate [kg/s]: "
-// << keyki << " "
-// << gSum((dmdtki*this->mesh().V())()) << endl;
- }
+ Info << "Explicit temperature mass transfer.." << endl;
+ Info << "keyki :" << keyki << endl;
+
+ dmdtNetki += Kexp.ref();
+ *dmdt_[keyki] = Kexp.ref();
}
+
}
-/*
- eqn +=
- (
- dmdtYki
-// + (-dSkidVar)*coeffsk*Trefk
-// + (dSikdVar)*coeffsi*Trefi
- )*tL.ref()
- * (invRhoCpk - invRhoCpi)
-
- - (
- dmdtYki
-// + (-dSkidVar)*coeffsk*Trefk
-// + (dSikdVar)*coeffsi*Trefi
- )*
- (
- invRhoCpk/iterk()->rho()
- - invRhoCpi/iteri()->rho()
- )*phasei.thermo().p();
-*/
-/*
- + fvm::SuSp
- (
- ((dSkidVar)*coeffsk + (-dSikdVar)*coeffsi)*L
- *(invRhoCpk - invRhoCpi)
- - ((dSkidVar)*coeffsk + (-dSikdVar)*coeffsi)*
- (
- invRhoCpk/iterk()->rho()
- - invRhoCpi/iteri()->rho()
- )*phasei.thermo().p(),
- T
- );
-*/
+
+ word keyikName(phasei.name() + phasek.name());
+ word keykiName(phasek.name() + phasei.name());
eqn -=
(
dmdtNetki
*(
- calculateL(keyik, keyki, T)
- - (phasek.Cp() - phasei.Cp())
- *dimensionedScalar("T0", dimTemperature, 298.0)
+ calculateL(dmdtNetki, keyik, keyki, T)
+ - (phasek.Cp() - phasei.Cp())
+ * dimensionedScalar("T0", dimTemperature, 298.0)
)
);
}
@@ -820,60 +536,29 @@ void Foam::MassTransferPhaseSystem<BasePhaseSystem>::massSpeciesTransfer
{
if (iter()->transferSpecie() == speciesName)
{
- //forAll (iter()->transferSpecie(), specieI)
- {
- //word speciesName = iter()->transferSpecie()[specieI];
-
- //const phasePair& pair(this->phasePairs_[iter.key()]);
-
- //const word from(pair.from().name());
- //const word to(pair.to().name());
-
- //const phasePairKey key(from, to, true);
-
- //const phaseModel& phaseFrom = this->phases()(from);
- //const volScalarField& phaseTo = phaseSystem::phases()(to);
+ // Extract alpha*div(u) added in alpha's
+ tmp<volScalarField> divU = fvc::div(this->phi());
-// volScalarField& dmdtYiY =
-// dmdtYi
-// (
-// key,
-// interfaceCompositionModel::modelVariableNames
-// [
-// interfaceCompositionModel::T
-// ],
-// speciesName
-// );
- DebugVar(max(this->Su()[phase.name()]));
- DebugVar(min(this->Su()[phase.name()]));
+ // Explicit source
+ Su =
+ this->Su()[phase.name()]
+ - divU.ref().internalField()*phase.oldTime()
+ + this->Sp()[phase.name()]*phase;
- DebugVar(max(this->Sp()[phase.name()]));
- DebugVar(min(this->Sp()[phase.name()]));
-
- // Extract alpha*div(u) added in alpha's
- tmp<volScalarField> divU = fvc::div(this->phi());
-
- // Explicit source
- Su =
- this->Su()[phase.name()]
- - divU.ref().internalField()*phase.oldTime()
- - this->Sp()[phase.name()]*phase;
-
- // Implicit source
- //Sp = this->Sp()[phase.name()];//*phase/phase.oldTime();
- }
+ // Implicit source
+ //Sp = this->Sp()[phase.name()];
}
}
}
-
+/*
template<class BasePhaseSystem>
void Foam::MassTransferPhaseSystem<BasePhaseSystem>::massTransfer
(
const volScalarField& T
)
{
- forAllConstIter(massTransferModelTable, massTransferModels_, iter)
+ forAllIter(massTransferModelTable, massTransferModels_, iter)
{
const phasePair& pair
(
@@ -894,35 +579,21 @@ void Foam::MassTransferPhaseSystem<BasePhaseSystem>::massTransfer
// Mass transfer with species
if (iter()->transferSpecie() != "none")
{
- //if (!iter()->semiImplicit())
+ if (Kexp.valid())
{
- //forAll (iter()->transferSpecie(), specieI)
- {
- //word speciesName = iter()->transferSpecie();
-
- //word nameDisp(IOobject::groupName(speciesName, from.name()));
- //word nameCont(IOobject::groupName(speciesName, to.name()));
-
- //Info<< "key : " << key << " nameDisp to nameCont " << " "
- // << nameDisp << " " << nameCont << endl;
+ volScalarField& dmdtYiY =
+ dmdtYi
+ (
+ iter.key(),
+ interfaceCompositionModel::modelVariableNames
+ [
+ interfaceCompositionModel::T
+ ],
+ speciesName
+ );
- if (Kexp.valid())
- {
-// volScalarField& dmdtYiY =
-// dmdtYi
-// (
-// iter.key(),
-// interfaceCompositionModel::modelVariableNames
-// [
-// interfaceCompositionModel::T
-// ],
-// speciesName
-// );
-
- *dmdt_[key] = Kexp.ref();
- //dmdtYiY = Kexp.ref();
- }
- }
+ *dmdt_[key] = Kexp.ref();
+ dmdtYiY = Kexp.ref();
}
}
else
@@ -932,65 +603,7 @@ void Foam::MassTransferPhaseSystem<BasePhaseSystem>::massTransfer
*dmdt_[key] = Kexp.ref();
}
}
- /*
- else if (iter()->semiImplicit())
- {
- // Implicit species based mass transfer
- tmp<volScalarField> Kimp = iter()->Kimp
- (
- interfaceCompositionModel::Y,
- T
- );
-
- if (Kimp.valid())
- {
- Info << "Semi-implicit species based mass transfer.." << endl;
-
- volScalarField& dmdtYiY =
- dmdtYi
- (
- iter.key(),
- interfaceCompositionModel::modelVariableNames
- [
- interfaceCompositionModel::Y
- ],
- speciesName
- );
-
- // Fill dm used in alpha's
- dmdt = Kimp.ref()*
- (
- max
- (
- eqns[nameCont]->psi()
- - iter()->Yf(speciesName, T),
- 0.0
- )
- );
- dmdtYiY = dmdt;
-
- *eqns[nameCont] +=
- - Kimp.ref()*iter()->Yf(speciesName, T)
- + fvm::Sp
- (
- Kimp.ref(),
- eqns[nameCont]->psi()
- );
-
- // Explicit in the from phase
- if (eqns.found(nameDisp))
- {
- *eqns[nameDisp] -= dmdt;
- }
- }
- }
- */
- //}
-
- Info<< "Species based Mass rate [kg/s]: "
- << key << " "
- << gSum((dmdt*this->mesh().V())()) << endl;
}
}
-
+*/
// ************************************************************************* //
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/MassTransferPhaseSystem/MassTransferPhaseSystem.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/MassTransferPhaseSystem/MassTransferPhaseSystem.H
index e86127d642..dcd0024cd7 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/MassTransferPhaseSystem/MassTransferPhaseSystem.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/MassTransferPhaseSystem/MassTransferPhaseSystem.H
@@ -116,6 +116,7 @@ protected:
//- Calculate L between phases
tmp<volScalarField> calculateL
(
+ const volScalarField& dmdtNetki,
const phasePairKey& keyik,
const phasePairKey& keyki,
const volScalarField& T
@@ -164,13 +165,13 @@ public:
);
- // Explicit and Implicit mass transfer matrices
+ // mass transfer functions
- //- Return the heat transfer matrix
+ //- Return the heat transfer matrix and fill dmdt for phases
virtual tmp<fvScalarMatrix> heatTransfer(const volScalarField& T);
//- Calculate mass transfer
- virtual void massTransfer(const volScalarField& T);
+ //virtual void massTransfer(const volScalarField& T);
//- Calculate mass transfer for species
virtual void massSpeciesTransfer
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.C
index d9b05114a2..791d2f4cf5 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.C
@@ -50,12 +50,6 @@ MultiComponentPhaseModel
)
:
BasePhaseModel(fluid, phaseName),
- residualAlpha_
- (
- "residualAlpha",
- dimless,
- fluid.mesh().solverDict("Yi").lookup("residualAlpha")
- ),
species_(),
inertIndex_(-1)
{
@@ -375,6 +369,12 @@ void Foam::MultiComponentPhaseModel<BasePhaseModel, phaseThermo>::solveYi
X_[inertIndex_].max(0.0);
calculateMassFractions();
+
+ if (mesh.time().outputTime())
+ {
+ X_[0].write();
+ X_[1].write();
+ }
}
template<class BasePhaseModel, class phaseThermo>
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.H
index 7001a8e934..068bf475da 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.H
@@ -60,9 +60,6 @@ protected:
// Protected data
- //- Residual phase fraction
- dimensionedScalar residualAlpha_;
-
//- Species table
hashedWordList species_;
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/phaseModel/phaseModel.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/phaseModel/phaseModel.C
index 8928564723..a4f7a2f09f 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/phaseModel/phaseModel.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/phaseModel/phaseModel.C
@@ -56,8 +56,7 @@ Foam::phaseModel::phaseModel
dimensionedScalar("alpha", dimless, 0)
),
fluid_(fluid),
- name_(phaseName),
- alphaMax_(fluid.subDict(phaseName).lookupOrDefault("alphaMax", 1.0))
+ name_(phaseName)
{}
@@ -81,12 +80,6 @@ const Foam::phaseSystem& Foam::phaseModel::fluid() const
}
-Foam::scalar Foam::phaseModel::alphaMax() const
-{
- return alphaMax_;
-}
-
-
void Foam::phaseModel::correct()
{
thermo().correct();
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/phaseModel/phaseModel.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/phaseModel/phaseModel.H
index ed4ae81fcf..bc7831859c 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/phaseModel/phaseModel.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseModel/phaseModel/phaseModel.H
@@ -64,9 +64,6 @@ class phaseModel
//- Name of phase
word name_;
- //- Optional maximum phase-fraction (e.g. packing limit)
- scalar alphaMax_;
-
public:
@@ -114,12 +111,6 @@ public:
//- Return the system to which this phase belongs
const phaseSystem& fluid() const;
- //- Return the other phase in a two-phase system
- //const phaseModel& otherPhase() const;
-
- //- Return the maximum phase-fraction (e.g. packing limit)
- scalar alphaMax() const;
-
//- Correct phase thermo
virtual void correct();
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/multiphaseSystem/multiphaseSystem.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/multiphaseSystem/multiphaseSystem.C
index 96d5ef4622..f2a9f9a0c1 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/multiphaseSystem/multiphaseSystem.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/multiphaseSystem/multiphaseSystem.C
@@ -283,7 +283,11 @@ void Foam::multiphaseSystem::calculateSuSp()
// Update ddtAlphaMax
ddtAlphaMax_ =
- max(gMax((dmdt21*coeffs1)()), gMax((dmdt12*coeffs2)()));
+ max
+ (
+ ddtAlphaMax_.value(),
+ max(gMax((dmdt21*coeffs1)()), gMax((dmdt12*coeffs2)()))
+ );
}
}
@@ -297,6 +301,9 @@ void Foam::multiphaseSystem::solve()
label nAlphaCorr(readLabel(alphaControls.lookup("nAlphaCorr")));
mesh.solverDict("alpha").lookup("cAlphas") >> cAlphas_;
+ // Reset ddtAlphaMax
+ ddtAlphaMax_ = dimensionedScalar("zero", dimless, 0.0);
+
PtrList<surfaceScalarField> phiAlphaCorrs(phases_.size());
const surfaceScalarField& phi = this->phi();
@@ -378,7 +385,7 @@ void Foam::multiphaseSystem::solve()
phirScheme
);
}
-/*
+
// Ensure that the flux at inflow BCs is preserved
forAll(phiAlphaCorr.boundaryField(), patchi)
{
@@ -400,8 +407,6 @@ void Foam::multiphaseSystem::solve()
}
}
}
-*/
-// limitedPhiAlphas_.set(phase1.name(), phiAlphaCorr);
phasei++;
}
@@ -488,8 +493,8 @@ void Foam::multiphaseSystem::solve()
phi,
upwind<scalar>(mesh, phi)
).fvmDiv(phi, alpha1)
- ==
- Su + fvm::Sp(Sp, alpha1)
+ ==
+ Su + fvm::Sp(Sp, alpha1)
);
alpha1Eqn.solve();
@@ -508,8 +513,6 @@ void Foam::multiphaseSystem::solve()
!(++alphaSubCycle).end();
)
{
- //surfaceScalarField phiAlphaCorrs0(phiAlphaCorrs[phasei]);
-
MULES::explicitSolve
(
geometricOneField(),
@@ -518,17 +521,17 @@ void Foam::multiphaseSystem::solve()
phiAlpha,
(alphaSubCycle.index()*Sp)(),
(Su - (alphaSubCycle.index() - 1)*Sp*alpha1)(),
- phase.alphaMax(),
+ 1,
0
);
if (alphaSubCycle.index() == 1)
{
- phase.alphaPhi() = phiAlpha;//phiAlphaCorrs0;
+ phase.alphaPhi() = phiAlpha;
}
else
{
- phase.alphaPhi() += phiAlpha;//phiAlphaCorrs0;
+ phase.alphaPhi() += phiAlpha;
}
}
@@ -540,8 +543,6 @@ void Foam::multiphaseSystem::solve()
phaseModel& phase = iter();
volScalarField& alpha1 = phase;
- //surfaceScalarField& phiAlpha = phiAlphaCorrs[phasei];
-
MULES::explicitSolve
(
geometricOneField(),
@@ -550,11 +551,10 @@ void Foam::multiphaseSystem::solve()
phiAlpha,
Sp,
Su,
- phase.alphaMax(),
+ 1,
0
);
-
phase.alphaPhi() = phiAlpha;
}
@@ -588,36 +588,38 @@ void Foam::multiphaseSystem::solve()
rhoPhi_ +=
fvc::interpolate(phase.rho())*phase.alphaPhi();
- Info<< alpha1.name() << " volume fraction = "
- << alpha1.weightedAverage(mesh.V()).value()
- << " Min(alpha) = " << min(alpha1).value()
- << " Max(alpha) = " << max(alpha1).value()
- << endl;
-
if (mesh.time().outputTime())
{
volScalarField dAlphadt("dAlphadt", fvc::ddt(alpha1));
volScalarField divrhoPhi("divrhoPhi", fvc::div(rhoPhi_));
+
divrhoPhi.write();
dAlphadt.write();
+
}
}
- // Normalize alphas
+ Info<< "Phase-sum volume fraction, min, max = "
+ << sumAlpha.weightedAverage(mesh_.V()).value()
+ << ' ' << min(sumAlpha).value()
+ << ' ' << max(sumAlpha).value()
+ << endl;
+
volScalarField sumCorr(1.0 - sumAlpha);
+
forAllIter(UPtrList<phaseModel>, phases_, iter)
{
phaseModel& phase = iter();
volScalarField& alpha = phase;
alpha += alpha*sumCorr;
- }
- Info<< "Phase-sum volume fraction, min, max = "
- << sumAlpha.weightedAverage(mesh_.V()).value()
- << ' ' << min(sumAlpha).value()
- << ' ' << max(sumAlpha).value()
- << endl;
+ Info<< alpha.name() << " volume fraction = "
+ << alpha.weightedAverage(mesh.V()).value()
+ << " Min(alpha) = " << min(alpha).value()
+ << " Max(alpha) = " << max(alpha).value()
+ << endl;
+ }
}
}
}
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/phaseSystem.C b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/phaseSystem.C
index 03a4cc32ff..8240133c03 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/phaseSystem.C
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/phaseSystem.C
@@ -217,11 +217,11 @@ Foam::phaseSystem::phaseSystem
"phi",
mesh_.time().timeName(),
mesh_,
- IOobject::NO_READ,
- IOobject::NO_WRITE
+ IOobject::READ_IF_PRESENT,
+ IOobject::AUTO_WRITE
),
mesh_,
- dimensionedScalar("rhoPhi", dimVolume/dimTime, 0.0)
+ dimensionedScalar("zero", dimVolume/dimTime, 0.0)
),
rhoPhi_
(
@@ -234,7 +234,7 @@ Foam::phaseSystem::phaseSystem
IOobject::NO_WRITE
),
mesh_,
- dimensionedScalar("rhoPhi", dimMass/dimTime, 0.0)
+ dimensionedScalar("zero", dimMass/dimTime, 0.0)
),
phaseModels_(generatePhaseModels(phaseNames_)),
phasePairs_(),
@@ -667,16 +667,16 @@ Foam::tmp<Foam::volScalarField> Foam::phaseSystem::kappa() const
tmp<volScalarField> tmpkappa
(
- phaseModelIter()()/phaseModelIter()->kappa()
+ phaseModelIter()()*phaseModelIter()->kappa()
);
++phaseModelIter;
for (; phaseModelIter != phaseModels_.end(); ++ phaseModelIter)
{
- tmpkappa.ref() += phaseModelIter()()/phaseModelIter()->kappa();
+ tmpkappa.ref() += phaseModelIter()()*phaseModelIter()->kappa();
}
- return 1.0/tmpkappa;
+ return tmpkappa;
}
@@ -687,7 +687,7 @@ Foam::tmp<Foam::scalarField> Foam::phaseSystem::kappa(const label patchI) const
tmp<scalarField> tmpKappa
(
phaseModelIter()().boundaryField()[patchI]
- /phaseModelIter()->kappa(patchI)
+ *phaseModelIter()->kappa(patchI)
);
++phaseModelIter;
@@ -695,10 +695,10 @@ Foam::tmp<Foam::scalarField> Foam::phaseSystem::kappa(const label patchI) const
{
tmpKappa.ref() +=
phaseModelIter()().boundaryField()[patchI]
- /phaseModelIter()->kappa(patchI);
+ *phaseModelIter()->kappa(patchI);
}
- return 1.0/tmpKappa;
+ return tmpKappa;
}
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/phaseSystem.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/phaseSystem.H
index 97d9bd844b..55f01a638b 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/phaseSystem.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/phasesSystem/phaseSystem/phaseSystem.H
@@ -520,7 +520,7 @@ public:
) = 0;
//- Calculate mass transfer
- virtual void massTransfer(const volScalarField& T) = 0;
+ //virtual void massTransfer(const volScalarField& T) = 0;
//- Calculate mass transfer for species
virtual void massSpeciesTransfer
diff --git a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/setDeltaT.H b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/setDeltaT.H
index 1f1b930639..9a9c4ed929 100644
--- a/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/setDeltaT.H
+++ b/applications/solvers/multiphase/icoReactingMultiPhaseInterFoam/setDeltaT.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2017 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -58,7 +58,6 @@ if (adjustTimeStep)
maxDeltaT
)
);
-
Info<< "deltaT = " << runTime.deltaTValue() << endl;
}
diff --git a/src/thermophysicalModels/radiation/radiationModels/radiationModel/radiationModel.C b/src/thermophysicalModels/radiation/radiationModels/radiationModel/radiationModel.C
index 72e74dce4e..e9021e0355 100644
--- a/src/thermophysicalModels/radiation/radiationModels/radiationModel/radiationModel.C
+++ b/src/thermophysicalModels/radiation/radiationModels/radiationModel/radiationModel.C
@@ -278,6 +278,20 @@ Foam::tmp<Foam::fvScalarMatrix> Foam::radiation::radiationModel::ST
}
+Foam::tmp<Foam::fvScalarMatrix> Foam::radiation::radiationModel::ST
+(
+ const volScalarField& rhoCp,
+ volScalarField& T
+) const
+{
+ return
+ (
+ Ru()/rhoCp
+ - fvm::Sp(Rp()*pow3(T)/rhoCp, T)
+ );
+}
+
+
const Foam::radiation::absorptionEmissionModel&
Foam::radiation::radiationModel::absorptionEmission() const
{
diff --git a/src/thermophysicalModels/radiation/radiationModels/radiationModel/radiationModel.H b/src/thermophysicalModels/radiation/radiationModels/radiationModel/radiationModel.H
index 68ca5dc21d..792022f512 100644
--- a/src/thermophysicalModels/radiation/radiationModels/radiationModel/radiationModel.H
+++ b/src/thermophysicalModels/radiation/radiationModels/radiationModel/radiationModel.H
@@ -252,6 +252,13 @@ public:
volScalarField& T
) const;
+ //- Temperature source term
+ virtual tmp<fvScalarMatrix> ST
+ (
+ const volScalarField& rhoCp,
+ volScalarField& T
+ ) const;
+
//- Access to absorptionEmission model
const absorptionEmissionModel& absorptionEmission() const;
--
GitLab