35const Foam::scalar Foam::liquidMixtureProperties::TrMax = 0.999;
49 properties_.setSize(components_.
size());
78 components_(lm.components_),
79 properties_(lm.properties_.size())
83 properties_.set(i, lm.properties_(i)->clone());
109 scalar x1 = X[i]*properties_[i].Vc();
111 vTc += x1*properties_[i].Tc();
114 return vTc/(vc + ROOTVSMALL);
124 Tpt += X[i]*properties_[i].Tt();
142 if (
p >= pv(
p, Thi, X))
146 else if (
p < pv(
p, Tlo, X))
149 <<
"Pressure below triple point pressure: "
150 <<
"p = " <<
p <<
" < Pt = " << pv(
p, Tlo, X) <<
nl <<
endl;
155 scalar
T = (Thi + Tlo)*0.5;
157 while ((Thi - Tlo) > 1.0e-4)
159 if ((pv(
p,
T, X) -
p) <= 0.0)
181 Tpc += X[i]*properties_[i].Tc();
195 Vc += X[i]*properties_[i].Vc();
196 Zc += X[i]*properties_[i].Zc();
199 return RR*Zc*Tpc(X)/Vc;
209 omega += X[i]*properties_[i].omega();
230 scalar Ti =
min(TrMax*properties_[i].Tc(), Tl);
231 Xs[i] = properties_[i].pv(
p, Ti)*Xl[i]/
p;
244 W += X[i]*properties_[i].W();
258 Y[i] = X[i]*properties_[i].W();
262 Y /= (sumY + ROOTVSMALL);
275 X[i] =
Y[i]/properties_[i].W();
279 X /= (sumX + ROOTVSMALL);
299 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
300 scalar
rho = properties_[i].rho(
p, Ti);
304 scalar Yi = X[i]*properties_[i].W();
311 return sumY/(v + ROOTVSMALL);
329 scalar Yi = X[i]*properties_[i].W();
332 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
333 pv += Yi*properties_[i].pv(
p, Ti);
337 return pv/(sumY + ROOTVSMALL);
355 scalar Yi = X[i]*properties_[i].W();
358 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
359 hl += Yi*properties_[i].hl(
p, Ti);
363 return hl/(sumY + ROOTVSMALL);
381 scalar Yi = X[i]*properties_[i].W();
384 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
385 Cp += Yi*properties_[i].Cp(
p, Ti);
389 return Cp/(sumY + ROOTVSMALL);
408 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
409 scalar Pvs = properties_[i].pv(
p, Ti);
415 Xs /= (XsSum + ROOTVSMALL);
421 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
422 sigma += Xs[i]*properties_[i].sigma(
p, Ti);
443 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
444 mu += X[i]*
log(properties_[i].
mu(
p, Ti));
465 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
467 scalar Vi = properties_[i].W()/properties_[i].rho(
p, Ti);
472 phii /= (pSum + ROOTVSMALL);
478 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
482 scalar Tj =
min(TrMax*properties_[j].Tc(),
T);
487 1.0/properties_[i].kappa(
p, Ti)
488 + 1.0/properties_[j].kappa(
p, Tj)
490 K += phii[i]*phii[j]*Kij;
512 scalar Ti =
min(TrMax*properties_[i].Tc(),
T);
513 Dinv += X[i]/properties_[i].D(
p, Ti);
517 return 1/(Dinv + ROOTVSMALL);
CGAL::Exact_predicates_exact_constructions_kernel K
const vector & W() const
Return the linear acceleration of the reference frame.
scalar Tc() const
Return the continuous phase temperature.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
void size(const label n)
Older name for setAddressableSize.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
bool isDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Check if entry is found and is a sub-dictionary.
wordList toc() const
Return the table of contents.
scalar Tpc(const scalarField &X) const
Return pseudocritical temperature according to Kay's rule.
scalarField Xs(const scalar p, const scalar Tg, const scalar Tl, const scalarField &Xg, const scalarField &Xl) const
Return the surface molar fractions.
scalar sigma(const scalar p, const scalar T, const scalarField &X) const
Estimate mixture surface tension [N/m].
scalar Tpt(const scalarField &X) const
Return pseudo triple point temperature (mole averaged formulation)
scalar pv(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture vapour pressure [Pa].
scalar pvInvert(const scalar p, const scalarField &X) const
Invert the vapour pressure relationship to retrieve the boiling.
scalar Ppc(const scalarField &X) const
Return pseudocritical pressure (modified Prausnitz and Gunn)
scalar hl(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture latent heat [J/kg].
const volScalarField & omega() const
Access functions.
Base-class for thermophysical properties of solids, liquids and gases providing an interface compatib...
PtrList< volScalarField > & Y
const volScalarField & mu
const volScalarField & Cp
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
const dimensionedScalar mu
Atomic mass unit.
const scalar RR
Universal gas constant: default in [J/(kmol K)].
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedScalar log(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
constexpr char nl
The newline '\n' character (0x0a)
const dimensionedScalar & D
#define forAll(list, i)
Loop across all elements in list.