32#include "phaseCompressibleTurbulenceModel.H"
33#include "virtualMassModel.H"
40namespace diameterModels
42namespace coalescenceModels
59Foam::diameterModels::coalescenceModels::Luo::
82 const sizeGroup& fi = popBal_.sizeGroups()[i];
83 const sizeGroup& fj = popBal_.sizeGroups()[j];
84 const phaseModel& continuousPhase = popBal_.continuousPhase();
91 popBal_.continuousPhase()
99 popBal_.continuousPhase()
107 *
cbrt(popBal_.continuousTurbulence().epsilon()*fi.
d())
112 pi/4.0*
sqr(fi.
d() + fj.
d())*uij
123 continuousPhase.
rho()*fi.
d()*
sqr(uij)
124 /popBal_.sigmaWithContinuousPhase(fi.
phase())
131 <<
"A virtual mass model for " << fi.
phase().
name() <<
" in "
132 << popBal_.continuousPhase().name() <<
" is not specified. This is "
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Base class for coalescence models.
Model of Luo (1993). The coalescence rate is calculated by.
virtual void addToCoalescenceRate(volScalarField &coalescenceRate, const label i, const label j)
Add to coalescenceRate.
Class that solves the univariate population balance equation by means of a class method (also called ...
This class represents a single sizeGroup belonging to a velocityGroup. The main property of a sizeGro...
const dimensionedScalar & d() const
Return representative diameter of the sizeGroup.
const phaseModel & phase() const
Return const-reference to the phase.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
const dimensionedScalar & rho() const
const word & name() const
virtual tmp< volScalarField > Cvm() const =0
Return the virtual mass coefficient.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
constexpr scalar pi(M_PI)
dimensionedScalar exp(const dimensionedScalar &ds)
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar cbrt(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)