33 template<
class CloudType>
42 const scalar Yliq = massliq/(massliq + masssol);
43 const scalar Ysol = 1 - Yliq;
44 Xliq = Yliq/liquids_.properties()[liqToLiqMap_].W();
45 Xsol = Ysol/this->owner().thermo().solids().
properties()[solToSolMap_].W();
46 Xliq /= (Xliq + Xsol);
51 template<
class CloudType>
62 this->owner().thermo().carrier().Y()[i][celli]
63 /this->owner().thermo().carrier().W(i);
70 template<
class CloudType>
81 template<
class CloudType>
93 <<
"Activity coefficient UNIFAC is not implemented " <<
nl 99 const scalar ic = this->coeffDict().getScalar(
"ic");
100 return inv((1 + ic*Xsol/(Xliq + ROOTVSMALL)));
110 template<
class CloudType>
119 gamma_(this->coeffDict().getScalar(
"gamma")),
120 alpha_(this->coeffDict().getScalar(
"alpham")),
121 liquids_(owner.
thermo().liquids()),
122 solution_(this->coeffDict().
lookup(
"solution")),
123 liqToCarrierMap_(-1),
130 <<
"Solution is not well defined. It should be (liquid solid)" 135 Info<<
"Participating liquid-solid species:" <<
endl;
142 const label idLiquid =
owner.composition().idLiquid();
147 const label idSolid =
owner.composition().idSolid();
152 const word activityCoefficienType
154 this->
coeffDict().getWord(
"activityCoefficient")
157 if (activityCoefficienType ==
"Hoff")
161 else if (activityCoefficienType ==
"UNIFAC")
168 <<
"activityCoefficient must be either 'Hoff' or 'UNIFAC'" 175 template<
class CloudType>
178 const LiquidEvapFuchsKnudsen<CloudType>& pcm
182 method_(pcm.method_),
185 liquids_(pcm.owner().
thermo().liquids()),
186 solution_(pcm.solution_),
187 liqToCarrierMap_(pcm.liqToCarrierMap_),
188 liqToLiqMap_(pcm.liqToLiqMap_),
189 solToSolMap_(pcm.solToSolMap_)
195 template<
class CloudType>
215 const scalar
rhog = this->owner().thermo().thermo().rho()()[celli];
217 const label gid = liqToCarrierMap_;
218 const label lid = liqToLiqMap_;
219 const label sid = solToSolMap_;
221 const scalar W = liquids_.properties()[lid].W();
223 const scalar YeInf = this->owner().thermo().carrier().Y()[gid][celli];
225 const scalar
sigma = liquids_.properties()[lid].sigma(pc, Ts);
231 const scalar Dab = liquids_.properties()[lid].D(pc, Ts);
234 const scalar pSat = liquids_.properties()[lid].pv(pc,
T);
236 scalar Xliq(0), Xsol(0);
237 calcXcSolution(liqMass[lid], solMass[sid], Xliq, Xsol);
240 const scalar
gamma = activityCoeff(Xliq, Xsol);
243 const scalar Rliq =
RR/W;
246 const scalar Kn = 2*gamma_/d;
248 (1+Kn)/(1+ (4/(3*alpha_) + 0.377)*Kn +
sqr(Kn)*4/(3*alpha_));
251 const scalar Sc =
nu/(Dab + ROOTVSMALL);
254 const scalar Sherwood = Sh(
Re, Sc);
257 const scalar Ni = (
rhog*Sherwood*Dab*Cm/d)*
log((1 - YeInf)/(1 - YeSurf));
261 dMassPC[lid] += Ni*As*dt;
265 template<
class CloudType>
277 switch (parent::enthalpyTransfer_)
279 case (parent::etLatentHeat):
284 case (parent::etEnthalpyDifference):
286 scalar hc = this->owner().composition().carrier().Ha(idc,
p,
T);
287 scalar hp = liquids_.properties()[idl].h(
p,
T);
303 template<
class CloudType>
313 template<
class CloudType>
323 return liquids_.pvInvert(
p, X);
dimensionedScalar rhog("rhog", dimDensity, transportProperties)
const dictionary & properties() const
Return const access to the properties dictionary.
dimensionedScalar Pr("Pr", dimless, laminarTransport)
void size(const label n)
Older name for setAddressableSize.
DSMCCloud< dsmcParcel > CloudType
Templated phase change model class.
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
activityCoeffMethodType method_
Method used.
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual void calculate(const scalar dt, const label celli, const scalar Re, const scalar Pr, const scalar d, const scalar nu, const scalar rho, const scalar T, const scalar Ts, const scalar pc, const scalar Tc, const scalarField &X, const scalarField &Xsol, const scalarField &liqMass, scalarField &dMassPC) const
Update model.
label liqToLiqMap_
Mapping between local and global liquid species.
Lookup type of boundary radiation properties.
virtual scalar Tvap(const scalarField &X) const
Return vapourisation temperature.
virtual scalar dh(const label idc, const label idl, const scalar p, const scalar T) const
Return the enthalpy per unit mass.
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
label liqToCarrierMap_
Mapping between liquid and carrier species.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const CloudType & owner() const
Return const access to the owner cloud.
#define forAll(list, i)
Loop across all elements in list.
scalar activityCoeff(const scalar Xliq, const scalar Ysol) const
Return activity coefficient.
psiReactionThermo & thermo
LiquidEvapFuchsKnudsen(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedScalar cbrt(const dimensionedScalar &ds)
virtual scalar TMax(const scalar p, const scalarField &X) const
Return maximum/limiting temperature.
tmp< scalarField > calcXc(const label celli) const
Calculate the carrier phase component volume fractions at celli.
void calcXcSolution(const scalar massliq, const scalar masssol, scalar &Xliq, scalar &Xsol) const
Calculate volumetric fractions of components in the solution.
constexpr scalar pi(M_PI)
scalarField Re(const UList< complex > &cf)
Extract real component.
const dictionary & coeffDict() const
Return const access to the coefficients dictionary.
errorManip< error > abort(error &err)
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
label solToSolMap_
Mapping between local and global solid species.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
PtrList< volScalarField > & Y
messageStream Info
Information stream (stdout output on master, null elsewhere)
const scalar RR
Universal gas constant: default in [J/(kmol K)].
A class for managing temporary objects.
Templated base class for dsmc cloud.
List< word > solution_
List of active liquid names i.e (liquidName solidName)
scalar Sh() const
Sherwood number.
static constexpr const zero Zero
Global zero (0)