Permalink
Browse files

reactingEulerFoam: improvements to population balance modeling

Removed possibility for the user to specify a driftRate in the constantDrift
model which is independent of a fvOptions mass source. The driftRate must be
calculated from/be consistent with the mass source in order to yield a particle
number conserving result.

Made calculation of the over-all Sauter mean diameter of an entire population
balance conditional on more than one velocityGroup being present. This diameter
field is for post-processing purposes only and would be redundant in case of one
velocityGroup being used.

Solution control is extended to allow for solution of the population balance
equation at the last PIMPLE loop only, using an optional switch. This can be
beneficial in terms of simulation time as well as coupling between the
population balance based diameter calculation and the rest of the equation
system.

Patch contributed by Institute of Fluid Dynamics, Helmholtz-Zentrum Dresden - Rossendorf
(HZDR) and VTT Technical Research Centre of Finland Ltd.
  • Loading branch information...
Henry Weller
Henry Weller committed Jan 4, 2018
1 parent 79207b6 commit 9f54506fbf7432dd7fe58824dbee8f9a88299417
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -51,8 +51,6 @@ Foam::diameterModels::driftModels::constantDrift::constantDrift
)
:
driftModel(popBal, dict),
calculateFromFvOptions_(dict.lookup("calculateFromFvOptions")),
rate_("rate", dimVolume/dimTime, dict.lookupOrDefault<scalar>("rate", 1.0)),
N_
(
IOobject
@@ -64,12 +62,7 @@ Foam::diameterModels::driftModels::constantDrift::constantDrift
popBal.mesh(),
dimensionedScalar("Sui", inv(dimVolume), Zero)
)
{
if (calculateFromFvOptions_ == false)
{
rate_.read(dict);
}
}
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@@ -97,15 +90,7 @@ void Foam::diameterModels::driftModels::constantDrift::driftRate
phaseModel& phase = const_cast<phaseModel&>(fi.phase());
volScalarField& rho = phase.thermo().rho();
if (calculateFromFvOptions_ == true)
{
driftRate +=(popBal_.fluid().fvOptions()(phase, rho)&rho)
/(N_*rho);
}
else
{
driftRate += rate_;
}
driftRate += (popBal_.fluid().fvOptions()(phase, rho)&rho)/(N_*rho);
}
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -25,10 +25,10 @@ Class
Foam::diameterModels::driftModels::constant
Description
Constant growth rate. Used for verification and validation of the drift
formulation implemented in the populationBalanceModel class. Rate is either
given in the corresponding dictionary or calculated from fvOptions mass
source, depending on Switch calculateFromFvOptions.
Constant drift rate within all classes. Used for verification and
validation of the drift formulation implemented in the
populationBalanceModel class. Rate is calculated from fvOptions mass source.
SourceFiles
constant.C
@@ -59,12 +59,6 @@ class constantDrift
{
// Private data
//- If true, calculate drift rate from fvOptions mass source
Switch calculateFromFvOptions_;
//- Optional driftRate
dimensionedScalar rate_;
//- Total number concentration
volScalarField N_;
@@ -969,6 +969,7 @@ Foam::diameterModels::populationBalanceModel::populationBalanceModel
(
fluid.subDict("populationBalanceCoeffs").subDict(name_)
),
pimple_(mesh_.lookupObject<pimpleControl>("solutionControl")),
continuousPhase_
(
mesh_.lookupObject<phaseModel>
@@ -1044,24 +1045,7 @@ Foam::diameterModels::populationBalanceModel::populationBalanceModel
Zero
)
),
d_
(
IOobject
(
IOobject::groupName("d", this->name()),
fluid.time().timeName(),
fluid.mesh(),
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
fluid.mesh(),
dimensionedScalar
(
IOobject::groupName("d", this->name()),
dimLength,
Zero
)
)
d_()
{
this->registerVelocityGroups();
@@ -1083,6 +1067,30 @@ Foam::diameterModels::populationBalanceModel::populationBalanceModel
dimensionedScalar("Sui", inv(dimTime), Zero)
)
);
d_.reset
(
new volScalarField
(
IOobject
(
IOobject::groupName("d", this->name()),
fluid.time().timeName(),
fluid.mesh(),
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
fluid.mesh(),
dimensionedScalar
(
IOobject::groupName("d", this->name()),
dimLength,
Zero
)
)
);
}
if (coalescence_.size() != 0)
@@ -1201,92 +1209,110 @@ bool Foam::diameterModels::populationBalanceModel::writeData(Ostream& os) const
void Foam::diameterModels::populationBalanceModel::solve()
{
const dictionary& populationBalanceControls = mesh_.solverDict(name_);
label nCorr(readLabel(populationBalanceControls.lookup("nCorr")));
scalar tolerance(readScalar(populationBalanceControls.lookup("tolerance")));
calcAlphas();
calcVelocity();
const dictionary& solutionControls = mesh_.solverDict(name_);
bool solveOnFinalIterOnly
(
solutionControls.lookupOrDefault<bool>
(
"solveOnFinalIterOnly",
false
)
);
if (nCorr > 0)
if (!solveOnFinalIterOnly || pimple_.finalIter())
{
preSolve();
}
int iCorr = 0;
scalar initialResidual = 0;
scalar maxInitialResidual = 1;
calcAlphas();
calcVelocity();
while
(
maxInitialResidual > tolerance
&&
++iCorr <= nCorr
)
{
Info<< "populationBalance "
<< this->name()
<< ": Iteration "
<< iCorr
<< endl;
label nCorr(readLabel(solutionControls.lookup("nCorr")));
scalar tolerance
(
readScalar(solutionControls.lookup("tolerance"))
);
sources();
if (nCorr > 0)
{
preSolve();
}
dmdt();
int iCorr = 0;
scalar initialResidual = 0;
scalar maxInitialResidual = 1;
forAll(sizeGroups_, i)
while
(
maxInitialResidual > tolerance
&&
++iCorr <= nCorr
)
{
sizeGroup& fi = *sizeGroups_[i];
const phaseModel& phase = fi.phase();
const volScalarField& alpha = phase;
const dimensionedScalar& residualAlpha = phase.residualAlpha();
const volScalarField& rho = phase.thermo().rho();
Info<< "populationBalance "
<< this->name()
<< ": Iteration "
<< iCorr
<< endl;
fvScalarMatrix sizeGroupEqn
(
fvm::ddt(alpha, rho, fi)
+ fi.VelocityGroup().mvConvection()->fvmDiv
sources();
dmdt();
forAll(sizeGroups_, i)
{
sizeGroup& fi = *sizeGroups_[i];
const phaseModel& phase = fi.phase();
const volScalarField& alpha = phase;
const dimensionedScalar& residualAlpha = phase.residualAlpha();
const volScalarField& rho = phase.thermo().rho();
fvScalarMatrix sizeGroupEqn
(
phase.alphaRhoPhi(),
fi
)
- fvm::Sp
fvm::ddt(alpha, rho, fi)
+ fi.VelocityGroup().mvConvection()->fvmDiv
(
phase.alphaRhoPhi(),
fi
)
- fvm::Sp
(
fvc::ddt(alpha, rho) + fvc::div(phase.alphaRhoPhi())
- fi.VelocityGroup().dmdt(),
fi
)
==
Su_[i]*rho
- fvm::SuSp(SuSp_[i]*rho, fi)
+ fvc::ddt(residualAlpha*rho, fi)
- fvm::ddt(residualAlpha*rho, fi)
);
sizeGroupEqn.relax
(
fvc::ddt(alpha, rho) + fvc::div(phase.alphaRhoPhi())
- fi.VelocityGroup().dmdt(),
fi
)
==
Su_[i]*rho
- fvm::SuSp(SuSp_[i]*rho, fi)
+ fvc::ddt(residualAlpha*rho, fi)
- fvm::ddt(residualAlpha*rho, fi)
);
fi.mesh().equationRelaxationFactor("f")
);
sizeGroupEqn.relax
(
fi.mesh().equationRelaxationFactor("f")
);
initialResidual = sizeGroupEqn.solve().initialResidual();
initialResidual = sizeGroupEqn.solve().initialResidual();
maxInitialResidual = max
(
initialResidual,
maxInitialResidual
);
}
}
maxInitialResidual = max
(
initialResidual,
maxInitialResidual
);
if (nCorr > 0)
{
forAll(velocityGroups_, i)
{
velocityGroups_[i]->postSolve();
}
}
}
if (nCorr > 0)
{
forAll(velocityGroups_, i)
if (velocityGroups_.size() > 1)
{
velocityGroups_[i]->postSolve();
d_() = dsm();
}
}
d_ = dsm();
}
// ************************************************************************* //
@@ -136,6 +136,7 @@ SourceFiles
#include "sizeGroup.H"
#include "phasePair.H"
#include "pimpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -184,6 +185,9 @@ class populationBalanceModel
//- Dictionary
dictionary dict_;
//- Reference to pimpleControl
const pimpleControl& pimple_;
//- Continuous phase
const phaseModel& continuousPhase_;
@@ -245,7 +249,7 @@ class populationBalanceModel
volVectorField U_;
//- Mean Sauter diameter
volScalarField d_;
autoPtr<volScalarField> d_;
// Private member functions
@@ -122,9 +122,6 @@ populationBalanceCoeffs
(
densityChange{}
);
zeroethOrderModels
();
}
}
@@ -104,9 +104,6 @@ populationBalanceCoeffs
(
densityChange{}
);
zeroethOrderModels
();
}
}
@@ -25,9 +25,10 @@ solvers
bubbles
{
nCorr 1;
tolerance 1e-4;
renormalize false;
nCorr 1;
tolerance 1e-4;
renormalize false;
solveOnFinalIterOnly true;
}
p_rgh

0 comments on commit 9f54506

Please sign in to comment.