Go to the documentation of this file.
39 namespace diameterModels
57 tmp<volScalarField> tInvDsm
71 const sizeGroup& fi = sizeGroups_[i];
81 Foam::diameterModels::velocityGroup::fSum()
const
83 tmp<volScalarField> tsumSizeGroups
97 sumSizeGroups += sizeGroups_[i];
100 return tsumSizeGroups;
104 void Foam::diameterModels::velocityGroup::renormalize()
106 Info<<
"Renormalizing sizeGroups for velocityGroup "
113 sizeGroups_[i] *=
pos(sizeGroups_[i]);
118 sizeGroups_[i] /= fSum_;
124 Foam::diameterModels::velocityGroup::mvconvection()
const
132 phase_.alphaRhoPhi(),
133 phase_.mesh().divScheme
135 "div(" + phase_.alphaRhoPhi()().
name() +
",f)"
153 popBalName_(diameterProperties.
lookup(
"populationBalance")),
167 phase.time().timeName(),
174 formFactor_(
"formFactor",
dimless, diameterProperties),
177 diameterProperties.
lookup(
"sizeGroups"),
193 phase.time().timeName(),
203 phase.time().timeName(),
216 phase.time().timeName(),
225 phase_.mesh().solverDict(popBalName_).getOrDefault<
Switch>
227 "renormalizeAtRestart",
231 phase_.mesh().solverDict(popBalName_).getOrDefault<
Switch>
245 mag(1 - fSum_.weightedAverage(fSum_.mesh().V()).value()) >= 1
e-5
246 ||
mag(1 -
max(fSum_).value()) >= 1
e-5
247 ||
mag(1 -
min(fSum_).value()) >= 1
e-5
251 <<
" Initial values of the sizeGroups belonging to velocityGroup "
252 << this->phase().
name()
253 <<
" must add to" <<
nl <<
" unity. This condition might be"
254 <<
" violated due to wrong entries in the" <<
nl
255 <<
" velocityGroupCoeffs subdictionary or bad initial conditions in"
256 <<
" the startTime" <<
nl
257 <<
" directory. The sizeGroups can be renormalized at every"
258 <<
" timestep or at restart" <<
nl
259 <<
" only by setting the corresponding switch renormalize or"
260 <<
" renormalizeAtRestart" <<
nl
261 <<
" in the fvSolution subdictionary " << popBalName_ <<
"."
262 <<
" Note that boundary conditions are not" <<
nl <<
"renormalized."
268 fields_.add(sizeGroups_[i]);
286 mvConvection_ = mvconvection();
294 Info<< this->phase().
name() <<
" Sauter mean diameter, min, max = "
295 << d_.weightedAverage(d_.mesh().V()).value()
296 <<
' ' <<
min(d_).value()
297 <<
' ' <<
max(d_).value()
302 Info<< phase_.name() <<
" sizeGroups-sum volume fraction, min, max = "
303 << fSum_.weightedAverage(phase_.mesh().V()).value()
304 <<
' ' <<
min(fSum_).value()
305 <<
' ' <<
max(fSum_).value()
310 phase_.mesh().solverDict(popBalName_).getOrDefault<
Switch>
defineTypeNameAndDebug(constant, 0)
velocityGroup(const dictionary &diameterProperties, const phaseModel &phase)
Construct from components.
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,...
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
A class for managing temporary objects.
const dimensionSet dimDensity
virtual bool read(const dictionary &diameterProperties)
Read diameterProperties dictionary.
static tmp< convectionScheme< Type > > New(const fvMesh &mesh, const surfaceScalarField &faceFlux, Istream &schemeData)
Return a pointer to a new convectionScheme created on freestore.
Helper class to manage multi-specie phase properties.
addToRunTimeSelectionTable(diameterModel, constant, dictionary)
A2stract base-class for dispersed-phase particle diameter models.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual ~velocityGroup()
Destructor.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &ds, const word &patchFieldType=fvPatchField< scalar >::calculatedType())
Return tmp field from name, mesh, dimensions and patch type.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
void postSolve()
Corrections after populationBalanceModel::solve()
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fv::convectionScheme< scalar > > mvConvection(fv::convectionScheme< scalar >::New(mesh, fields, phi, mesh.divScheme("div(phi,Yi_h)")))
GeometricField< scalar, fvPatchField, volMesh > volScalarField
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
ITstream & lookup(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const word & name() const
const phaseModel & phase_
void preSolve()
Corrections before populationBalanceModel::solve()
const word & name() const
Return the name of this phase.
virtual tmp< volScalarField > d() const
Return the Sauter-mean diameter.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Return a pointer to a new sizeGroup created on freestore.
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar e
Elementary charge.
virtual bool read(const dictionary &phaseProperties)=0
Read phaseProperties dictionary.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
dimensionedScalar pos(const dimensionedScalar &ds)