34 template<
class CloudType>
45 this->owner().thermo().carrier().Y()[i][cellI]
46 /this->owner().thermo().carrier().W(i);
53 template<
class CloudType>
66 template<
class CloudType>
74 liquids_(owner.thermo().liquids()),
75 activeLiquids_(this->coeffDict().lookup(
"activeLiquids")),
76 liqToCarrierMap_(activeLiquids_.size(), -1),
77 liqToLiqMap_(activeLiquids_.size(), -1)
79 if (activeLiquids_.size() == 0)
83 "Foam::LiquidEvaporationBoil<CloudType>::LiquidEvaporationBoil" 85 "const dictionary& dict, " 88 ) <<
"Evaporation model selected, but no active liquids defined" 93 Info<<
"Participating liquid species:" <<
endl;
98 Info<<
" " << activeLiquids_[i] <<
endl;
100 owner.composition().carrierId(activeLiquids_[i]);
104 const label idLiquid = owner.composition().idLiquid();
108 owner.composition().localId(idLiquid, activeLiquids_[i]);
114 template<
class CloudType>
121 liquids_(pcm.
owner().thermo().liquids()),
130 template<
class CloudType>
137 template<
class CloudType>
155 if ((liquids_.Tc(X) -
T) < SMALL)
161 "void Foam::LiquidEvaporationBoil<CloudType>::calculate" 173 "const scalarField&, " 176 ) <<
"Parcel reached critical conditions: " 177 <<
"evaporating all avaliable mass" <<
endl;
182 const label lid = liqToLiqMap_[i];
183 dMassPC[lid] = GREAT;
190 scalar ps = liquids_.pv(pc, Ts, X);
193 scalar rhos = ps*liquids_.W(X)/(
RR*Ts);
205 scalar Yc = this->owner().thermo().carrier().Y()[i][cellI];
206 Hc += Yc*this->owner().thermo().carrier().Ha(i, pc, Tc);
207 Hsc += Yc*this->owner().thermo().carrier().Ha(i, ps, Ts);
208 Cpc += Yc*this->owner().thermo().carrier().Cp(i, ps, Ts);
209 kappac += Yc*this->owner().thermo().carrier().kappa(i, ps, Ts);
215 const label gid = liqToCarrierMap_[i];
216 const label lid = liqToLiqMap_[i];
219 const scalar TBoil = liquids_.properties()[lid].pvInvert(pc);
222 const scalar Td =
min(T, 0.999*TBoil);
225 const scalar
pSat = liquids_.properties()[lid].pv(pc, Td);
228 const scalar Xc = XcMix[gid];
238 const scalar Dab = liquids_.properties()[lid].D(ps, Ts);
241 const scalar Sc = nu/(Dab + ROOTVSMALL);
244 const scalar
Sh = this->
Sh(Re, Sc);
251 const scalar deltaT =
max(T - TBoil, 0.5);
254 const scalar hv = liquids_.properties()[lid].hl(pc, Td);
260 alphaS = 760.0*
pow(deltaT, 0.26);
262 else if (deltaT < 25.0)
264 alphaS = 27.0*
pow(deltaT, 2.33);
268 alphaS = 13800.0*
pow(deltaT, 0.39);
272 const scalar Gf = alphaS*deltaT*
pi*
sqr(d)/hv;
276 const scalar A = (Hc - Hsc)/hv;
277 const scalar B =
pi*kappac/Cpc*d*
Sh;
286 for (
label i=0; i<50; i++)
290 G = B/(1.0 + Gr)*
log(1.0 + A*(1.0 + Gr));
293 if (
mag(Gr - GrDash)/GrDash < 1
e-3)
300 dMassPC[lid] += (G + Gf)*dt;
307 const scalar Xs = X[lid]*pSat/pc;
310 const scalar Xr = (Xs - Xc)/
max(SMALL, 1.0 - Xs);
315 dMassPC[lid] +=
pi*d*Sh*Dab*rhos*
log(1.0 + Xr)*dt;
323 template<
class CloudType>
335 if (liquids_.properties()[idl].pv(p, T) >= 0.999*
p)
337 TDash = liquids_.properties()[idl].pvInvert(p);
341 switch (parent::enthalpyTransfer_)
343 case (parent::etLatentHeat):
348 case (parent::etEnthalpyDifference):
350 scalar hc = this->owner().composition().carrier().Ha(idc, p, TDash);
351 scalar hp = liquids_.
properties()[idl].h(p, TDash);
360 "Foam::scalar Foam::LiquidEvaporationBoil<CloudType>::dh" 375 template<
class CloudType>
381 return liquids_.Tpt(X);
385 template<
class CloudType>
392 return liquids_.pvInvert(p, X);
dimensionedScalar sqrt(const dimensionedScalar &ds)
const CloudType & owner() const
Return const access to the owner cloud.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
virtual scalar Tvap(const scalarField &X) const
Return vapourisation temperature.
List< word > activeLiquids_
List of active liquid names.
dimensioned< scalar > mag(const dimensioned< Type > &)
const dimensionedScalar & pSat
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
LiquidEvaporationBoil(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
virtual ~LiquidEvaporationBoil()
Destructor.
virtual scalar dh(const label idc, const label idl, const scalar p, const scalar T) const
Return the enthalpy per unit mass.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedScalar cbrt(const dimensionedScalar &ds)
Liquid evaporation model.
tmp< scalarField > calcXc(const label cellI) const
Calculate the carrier phase component volume fractions at cellI.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Templated phase change model class.
dimensionedScalar log(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
#define WarningIn(functionName)
Report a warning using Foam::Warning.
PtrList< volScalarField > & Y
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
List< label > liqToLiqMap_
Mapping between local and global liquid species.
const dictionary & properties() const
Return const access to the properties dictionary.
List< label > liqToCarrierMap_
Mapping between liquid and carrier species.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
errorManip< error > abort(error &err)
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
const scalar RR
Universal gas constant (default in [J/(kmol K)])
virtual scalar TMax(const scalar p, const scalarField &X) const
Return maximum/limiting temperature.
Templated base class for dsmc cloud.
psiReactionThermo & thermo
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
scalar Sh() const
Sherwood number.
virtual void calculate(const scalar dt, const label cellI, const scalar Re, const scalar Pr, const scalar d, const scalar nu, const scalar T, const scalar Ts, const scalar pc, const scalar Tc, const scalarField &X, scalarField &dMassPC) const
Update model.
A class for managing temporary objects.
const dimensionedScalar G
Newtonian constant of gravitation.