35 const Foam::scalar Foam::liquidMixtureProperties::TrMax = 0.999;
48 components_ =
dict.toc();
49 properties_.setSize(components_.
size());
79 const liquidMixtureProperties& lm
82 components_(lm.components_),
83 properties_(lm.properties_.clone())
108 scalar x1 = X[i]*properties_[i].Vc();
110 vTc += x1*properties_[i].Tc();
113 return vTc/(vc + ROOTVSMALL);
123 Tpt += X[i]*properties_[i].Tt();
141 if (
p >= pv(
p, Thi, X))
145 else if (
p < pv(
p, Tlo, X))
148 <<
"Pressure below triple point pressure: " 149 <<
"p = " <<
p <<
" < Pt = " << pv(
p, Tlo, X) <<
nl <<
endl;
154 scalar
T = (Thi + Tlo)*0.5;
156 while ((Thi - Tlo) > 1.0e-4)
158 if ((pv(
p,
T, X) -
p) <= 0.0)
180 Tpc += X[i]*properties_[i].Tc();
194 Vc += X[i]*properties_[i].Vc();
195 Zc += X[i]*properties_[i].Zc();
198 return RR*Zc*Tpc(X)/Vc;
208 omega += X[i]*properties_[i].omega();
229 scalar Ti =
min(TrMax*properties_[i].Tc(), Tl);
230 Xs[i] = properties_[i].pv(
p, Ti)*Xl[i]/
p;
243 W += X[i]*properties_[i].W();
257 Y[i] = X[i]*properties_[i].W();
261 Y /= (sumY + ROOTVSMALL);
274 X[i] =
Y[i]/properties_[i].W();
278 X /= (sumX + ROOTVSMALL);
298 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
299 scalar
rho = properties_[i].rho(
p, Ti);
303 scalar Yi = X[i]*properties_[i].W();
310 return sumY/(v + ROOTVSMALL);
328 scalar Yi = X[i]*properties_[i].W();
331 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
332 pv += Yi*properties_[i].pv(
p, Ti);
336 return pv/(sumY + ROOTVSMALL);
354 scalar Yi = X[i]*properties_[i].W();
357 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
358 hl += Yi*properties_[i].hl(
p, Ti);
362 return hl/(sumY + ROOTVSMALL);
380 scalar Yi = X[i]*properties_[i].W();
383 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
384 Cp += Yi*properties_[i].Cp(
p, Ti);
388 return Cp/(sumY + ROOTVSMALL);
407 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
408 scalar Pvs = properties_[i].pv(
p, Ti);
414 Xs /= (XsSum + ROOTVSMALL);
420 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
421 sigma += Xs[i]*properties_[i].sigma(
p, Ti);
442 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
443 mu += X[i]*
log(properties_[i].
mu(
p, Ti));
464 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
466 scalar Vi = properties_[i].W()/properties_[i].rho(
p, Ti);
471 phii /= (pSum + ROOTVSMALL);
477 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
481 scalar Tj =
min(TrMax*properties_[j].Tc(),
T);
486 1.0/properties_[i].kappa(
p, Ti)
487 + 1.0/properties_[j].kappa(
p, Tj)
489 K += phii[i]*phii[j]*Kij;
511 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
512 Dinv += X[i]/properties_[i].D(
p, Ti);
516 return 1/(Dinv + ROOTVSMALL);
scalar Cp(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture heat capacity [J/(kg K)].
Base-class for thermophysical properties of solids, liquids and gases providing an interface compatib...
void size(const label n)
Older name for setAddressableSize.
dimensionedScalar log(const dimensionedScalar &ds)
A list of keyword definitions, which are a keyword followed by a number of values (eg...
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
liquidMixtureProperties(const dictionary &dict)
Construct from dictionary.
scalar D(const scalar p, const scalar T, const scalarField &X) const
Vapour diffusivity [m2/s].
CGAL::Exact_predicates_exact_constructions_kernel K
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
scalar mu(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture viscosity [Pa s].
scalar kappa(const scalar p, const scalar T, const scalarField &X) const
Estimate thermal conductivity [W/(m K)].
#define forAll(list, i)
Loop across all elements in list.
scalarField Xs(const scalar p, const scalar Tg, const scalar Tl, const scalarField &Xg, const scalarField &Xl) const
Return the surface molar fractions.
scalar rho(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture density [kg/m^3].
scalar Tpc(const scalarField &X) const
Return pseudocritical temperature according to Kay's rule.
static autoPtr< liquidMixtureProperties > New(const dictionary &)
Select construct from dictionary.
scalar Ppc(const scalarField &X) const
Return pseudocritical pressure (modified Prausnitz and Gunn)
dimensionedScalar exp(const dimensionedScalar &ds)
scalarField X(const scalarField &Y) const
Returns the mole fractions corresponding to the given mass fractions.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
scalar sigma(const scalar p, const scalar T, const scalarField &X) const
Estimate mixture surface tension [N/m].
scalar Tc(const scalarField &X) const
Calculate the critical temperature of mixture.
scalarField Y(const scalarField &X) const
Returns the mass fractions corresponding to the given mole fractions.
scalar Tpt(const scalarField &X) const
Return pseudo triple point temperature (mole averaged formulation)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
const volScalarField & Cp
static autoPtr< liquidProperties > New(const word &name)
Return a pointer to a new liquidProperties created from name.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar hl(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture latent heat [J/kg].
const dimensionedScalar mu
Atomic mass unit.
#define WarningInFunction
Report a warning using Foam::Warning.
PtrList< volScalarField > & Y
scalar pv(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture vapour pressure [Pa].
const scalar RR
Universal gas constant: default in [J/(kmol K)].
scalar W(const scalarField &X) const
Calculate the mean molecular weight [kg/kmol].
static autoPtr< T > New(Args &&... args)
Construct autoPtr with forwarding arguments.
scalar omega(const scalarField &X) const
Return mixture accentric factor.
scalar pvInvert(const scalar p, const scalarField &X) const
Invert the vapour pressure relationship to retrieve the boiling.