Go to the documentation of this file.
31 #include "phaseCompressibleTurbulenceModel.H"
38 namespace diameterModels
40 namespace binaryBreakupModels
63 gammaUpperReg2by11_(),
64 gammaUpperReg5by11_(),
65 gammaUpperReg8by11_(),
66 C4_(dimensionedScalar::getOrDefault(
"C4",
dict,
dimless, 0.923)),
67 beta_(dimensionedScalar::getOrDefault(
"beta",
dict,
dimless, 2.05)),
70 dimensionedScalar::getOrDefault(
"minEddyRatio",
dict,
dimless, 11.4)
72 kolmogorovLengthScale_
76 "kolmogorovLengthScale",
77 popBal_.time().timeName(),
83 "kolmogorovLengthScale",
97 for (scalar z = 1
e-2; z <= 10.0; z = z + 1
e-2)
117 gammaUpperReg2by11Table.
append(gamma2by11);
118 gammaUpperReg5by11Table.
append(gamma5by11);
119 gammaUpperReg8by11Table.
append(gamma8by11);
122 gammaUpperReg2by11_.reset
126 gammaUpperReg2by11Table,
127 bounds::repeatableBounding::CLAMP,
132 gammaUpperReg5by11_.reset
136 gammaUpperReg5by11Table,
137 bounds::repeatableBounding::CLAMP,
142 gammaUpperReg8by11_.reset
146 gammaUpperReg8by11Table,
147 bounds::repeatableBounding::CLAMP,
158 kolmogorovLengthScale_ =
178 const phaseModel& continuousPhase = popBal_.continuousPhase();
179 const sizeGroup& fi = popBal_.sizeGroups()[i];
180 const sizeGroup& fj = popBal_.sizeGroups()[j];
184 pow(fi.
x()/fj.
x(), 2.0/3.0) +
pow((1 - fi.
x()/fj.
x()), 2.0/3.0) - 1
189 12.0*cf*popBal_.sigmaWithContinuousPhase(fi.
phase())
191 beta_*continuousPhase.
rho()*
pow(fj.
d(), 5.0/3.0)
192 *
pow(popBal_.continuousTurbulence().epsilon(), 2.0/3.0)
196 const volScalarField xiMin(minEddyRatio_*kolmogorovLengthScale_/fj.
d());
205 2.0*
pow(
b[celli], 3.0/11.0)*tgamma(5.0/11.0)
207 gammaUpperReg5by11_()(
b[celli])
208 - gammaUpperReg5by11_()(tMin[celli])
213 C4_*(1 - popBal_.alphas())/fj.
x()
216 popBal_.continuousTurbulence().epsilon()
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
const phaseModel & phase() const
Return const-reference to the phase.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
static constexpr const zero Zero
Global zero (0)
virtual tmp< volScalarField > rho() const =0
Return the density field.
virtual void addToBinaryBreakupRate(volScalarField &binaryBreakupRate, const label i, const label j)
Add to binary breakupRate.
addToRunTimeSelectionTable(binaryBreakupModel, LehrMilliesMewes, dictionary)
void append(const T &val)
Append an element at the end of the list.
scalar incGammaRatio_Q(const scalar a, const scalar x)
Normalized upper incomplete gamma function.
#define forAll(list, i)
Loop across all elements in list.
dimensionedScalar pow025(const dimensionedScalar &ds)
virtual tmp< volScalarField > nu() const =0
Return the laminar kinematic viscosity.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
This class represents a single sizeGroup belonging to a velocityGroup. The main property of a sizeGro...
LuoSvendsen(const populationBalanceModel &popBal, const dictionary &dict)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Class that solves the univariate population balance equation by means of a class method (also called ...
Macros for easy insertion into run-time selection tables.
const phaseModel & continuousPhase() const
Return continuous phase.
defineTypeNameAndDebug(LehrMilliesMewes, 0)
const dimensionedScalar & d() const
Return representative diameter of the sizeGroup.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const populationBalanceModel & popBal_
Reference to the populationBalanceModel.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Base class for binary breakup models which give the breakup rate between a sizeGroup pair directly,...
const phaseCompressibleTurbulenceModel & continuousTurbulence() const
Return reference to turbulence model of the continuous phase.
const dimensionedScalar e
Elementary charge.
A 2-tuple for storing two objects of dissimilar types. The container is similar in purpose to std::pa...
dimensionedScalar cbrt(const dimensionedScalar &ds)
const dimensionedScalar & x() const
Return representative volume of the sizeGroup.
virtual void correct()
Correct diameter independent expressions.