CoolProp  4.2.5
An open-source fluid property and humid air property database
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Public Attributes | List of all members
Mixture Class Reference

#include <Mixtures.h>

Public Member Functions

 Mixture (std::vector< Fluid * > pFluids)
 
 Mixture (std::string FluidsString)
 
 ~Mixture ()
 
void initialize ()
 
double Rbar (const std::vector< double > &x)
 
double Wilson_lnK_factor (double T, double p, int i)
 
double phir (double tau, double delta, const std::vector< double > &x)
 
double d2phir_dDelta_dTau (double tau, double delta, const std::vector< double > &x)
 
double d2phir_dTau2 (double tau, double delta, const std::vector< double > &x)
 
double dphir_dDelta (double tau, double delta, const std::vector< double > &x)
 
double d2phir_dDelta2 (double tau, double delta, const std::vector< double > &x)
 
double dphir_dTau (double tau, double delta, const std::vector< double > &x)
 
double dphir_dxi (double tau, double delta, const std::vector< double > &x, int i)
 
double d2phir_dxi_dTau (double tau, double delta, const std::vector< double > &x, int i)
 
double d2phir_dxi_dDelta (double tau, double delta, const std::vector< double > &x, int i)
 
double d2phirdxidxj (double tau, double delta, const std::vector< double > &x, int i, int j)
 
double fugacity (double tau, double delta, const std::vector< double > &x, int i)
 Returns the fugacity for the given component for the given total reduced density and reciprocal reduced temperature. More...
 
double rhobar_Tpz (double T, double p, const std::vector< double > &z, double rhobar0)
 
void TpzFlash (double T, double p, const std::vector< double > &z, double &rhobar, std::vector< double > &x, std::vector< double > &y)
 
double g_RachfordRice (const std::vector< double > &z, const std::vector< double > &lnK, double beta)
 
double dgdbeta_RachfordRice (const std::vector< double > &z, const std::vector< double > &lnK, double beta)
 
double ndphir_dni__constT_V_nj (double tau, double delta, const std::vector< double > &x, int i)
 
double partial_molar_volume (double tau, double delta, const std::vector< double > &x, int i)
 
double dpdT__constV_n (double tau, double delta, const std::vector< double > &x, int i)
 
double ndpdV__constT_n (double tau, double delta, const std::vector< double > &x, int i)
 
double ndpdni__constT_V_nj (double tau, double delta, const std::vector< double > &x, int i)
 
double ln_fugacity_coefficient (double tau, double delta, const std::vector< double > &x, int i)
 
double dln_fugacity_coefficient_dT__constrho (double tau, double delta, const std::vector< double > &x, int i)
 
double dln_fugacity_coefficient_dT__constp_n (double tau, double delta, const std::vector< double > &x, int i)
 
double dnphir_dni__constT_V_nj (double tau, double delta, const std::vector< double > &x, int i)
 
double dln_fugacity_coefficient_dp__constT_n (double tau, double delta, const std::vector< double > &x, int i)
 
double ndln_fugacity_coefficient_dnj__constT_p (double tau, double delta, const std::vector< double > &x, int i, int j)
 
double dln_fugacity_coefficient_dxj__constT_p_xi (double tau, double delta, const std::vector< double > &x, int i, int j)
 Gernert Equation 3.115. More...
 
double dpdxj__constT_V_xi (double tau, double delta, const std::vector< double > &x, int j)
 
double d2nphir_dni_dxj__constT_V (double tau, double delta, const std::vector< double > &x, int i, int j)
 Gernert Equation 3.117. More...
 
double dphir_dxj__constT_V_xi (double tau, double delta, const std::vector< double > &x, int j)
 
double d_ndphirdni_dxj__constT_V_xi (double tau, double delta, const std::vector< double > &x, int i, int j)
 
double d_dphirddelta_dxj__constT_V_xi (double tau, double delta, const std::vector< double > &x, int j)
 
double ddelta_dxj__constT_V_xi (double tau, double delta, const std::vector< double > &x, int j)
 
double dtau_dxj__constT_V_xi (double tau, double delta, const std::vector< double > &x, int j)
 
double d2nphir_dni_dT (double tau, double delta, const std::vector< double > &x, int i)
 
double d_ndphirdni_dTau (double tau, double delta, const std::vector< double > &x, int i)
 
double d_ndphirdni_dDelta (double tau, double delta, const std::vector< double > &x, int i)
 
double nd2nphirdnidnj__constT_V (double tau, double delta, const std::vector< double > &x, int i, int j)
 
double nddeltadni__constT_V_nj (double tau, double delta, const std::vector< double > &x, int i)
 
double ndtaudni__constT_V_nj (double tau, double delta, const std::vector< double > &x, int i)
 
double d_ndphirdni_dxj__constdelta_tau_xi (double tau, double delta, const std::vector< double > &x, int i, int j)
 
double saturation_p_preconditioner (double p, const std::vector< double > &z)
 
double saturation_p_Wilson (double beta, double p, const std::vector< double > &z, double T_guess, std::vector< double > &K)
 
double saturation_p (double beta, double p, const std::vector< double > &z, std::vector< double > &x, std::vector< double > &y)
 
double rhobar_pengrobinson (double T, double p, const std::vector< double > &x, int solution)
 
void x_and_y_from_K (double beta, const std::vector< double > &K, const std::vector< double > &z, std::vector< double > &x, std::vector< double > &y)
 
void load_excess_values (int i, int j)
 
std::map< std::string, double > load_reducing_values (int i, int j)
 
void check_MethaneEthane ()
 
void check_WaterEthanol ()
 
void test ()
 

Public Attributes

unsigned int N
 
std::vector< Fluid * > pFluids
 
ReducingFunctionpReducing
 
ExcessTermpExcess
 
ResidualIdealMixturepResidualIdealMix
 
PhaseEnvelope Envelope
 
SuccessiveSubstitutionVLE SS
 
NewtonRaphsonVLE NRVLE
 

Detailed Description

This is the class that actually implements the mixture properties

Definition at line 441 of file Mixtures.h.

Constructor & Destructor Documentation

Mixture::Mixture ( std::vector< Fluid * >  pFluids)

Definition at line 319 of file Mixtures.cpp.

Mixture::Mixture ( std::string  FluidsString)

Definition at line 329 of file Mixtures.cpp.

Mixture::~Mixture ( )

Definition at line 919 of file Mixtures.cpp.

Member Function Documentation

void Mixture::check_MethaneEthane ( )

Definition at line 691 of file Mixtures.cpp.

void Mixture::check_WaterEthanol ( )

Definition at line 462 of file Mixtures.cpp.

double Mixture::d2nphir_dni_dT ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

\[ \left(\frac{\partial^2n\alpha^r}{\partial T\partial n_i} \right)_{V,n_j} = \left( \frac{\partial}{\partial T}\left(\frac{\partial n \alpha^r}{\partial n_i}\right)_{T,V,n_j} \right)_{V,\bar n} \]

\[ \left(\frac{\partial^2n\alpha^r}{\partial T\partial n_i} \right)_{V,n_j} = -\frac{\tau}{T}\left[\alpha_{\tau}^r +\left( \frac{\partial}{\partial \tau}\left(n\left(\frac{\partial \alpha^r}{\partial n_i}\right)_{T,V,n_j}\right)\right)_{\delta,\bar x}\right] \]

GERG 2004 Monograph, equations 7.44 and 7.51

Definition at line 965 of file Mixtures.cpp.

double Mixture::d2nphir_dni_dxj__constT_V ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i,
int  j 
)

Gernert Equation 3.117.

Definition at line 1027 of file Mixtures.cpp.

double Mixture::d2phir_dDelta2 ( double  tau,
double  delta,
const std::vector< double > &  x 
)

Definition at line 1218 of file Mixtures.cpp.

double Mixture::d2phir_dDelta_dTau ( double  tau,
double  delta,
const std::vector< double > &  x 
)

Definition at line 1222 of file Mixtures.cpp.

double Mixture::d2phir_dTau2 ( double  tau,
double  delta,
const std::vector< double > &  x 
)

Definition at line 1226 of file Mixtures.cpp.

double Mixture::d2phir_dxi_dDelta ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

Definition at line 1210 of file Mixtures.cpp.

double Mixture::d2phir_dxi_dTau ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

Definition at line 1206 of file Mixtures.cpp.

double Mixture::d2phirdxidxj ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i,
int  j 
)

Definition at line 1202 of file Mixtures.cpp.

double Mixture::d_dphirddelta_dxj__constT_V_xi ( double  tau,
double  delta,
const std::vector< double > &  x,
int  j 
)

Gernert Equation 3.134 Catch test provided

Definition at line 1020 of file Mixtures.cpp.

double Mixture::d_ndphirdni_dDelta ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

The derivative term

\begin{eqnarray*} \left(\frac{\partial }{\partial \delta} \left( n\left(\frac{\partial \phi^r}{\partial n_i} \right)_{T,V,n_j} \right)\right)_{\tau,\bar x} &=& (\alpha_{\delta}^r+\delta\alpha_{\delta\delta}^r)\left[1-\frac{1}{\rho_r}\cdot n\left(\frac{\partial \rho_r}{\partial n_i}\right)_{n_j} \right] \\ &+&\tau\alpha^r_{\delta\tau}\frac{1}{T_r}\cdot n\left(\frac{\partial T_r}{\partial n_i}\right)_{n_j}\\ &+&\phi^r_{\delta x_i}-\sum_{k=1}^{N}x_k\phi^r_{\delta x_k} \end{eqnarray*}

GERG 2004 Monograph Equation 7.50 and Table B4, Kunz, JCED, 2012

Definition at line 1154 of file Mixtures.cpp.

double Mixture::d_ndphirdni_dTau ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

The derivative term

\begin{eqnarray*} \frac{\partial }{\partial \tau} \left( n\left(\frac{\partial \phi^r}{\partial n_i} \right)_{T,V,n_j} \right) &=& \delta \phi^r_{\delta\tau}\left[ 1-\frac{1}{\rho_r}\left[\left(\frac{\partial \rho_r}{\partial x_i}\right)_{x_j} - \sum_{k=1}^N x_k\left(\frac{\partial \rho_r}{\partial x_k}\right)_{x_j} \right]\right]\\ && +(\tau \phi^r_{\tau\tau}+\phi^r_{\tau})\frac{1}{T_r}\left[\left(\frac{\partial T_r}{\partial x_i}\right)_{x_j} - \sum_{k=1}^N x_k\left(\frac{\partial T_r}{\partial x_k}\right)_{x_j} \right]\\ && +\phi^r_{x_i\tau}-\sum_{k=1}^{N}x_k\phi^r_{x_k\tau} \end{eqnarray*}

GERG 2004 Monograph Equation 7.51 and Table B4, Kunz, JCED, 2012

Definition at line 1174 of file Mixtures.cpp.

double Mixture::d_ndphirdni_dxj__constdelta_tau_xi ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i,
int  j 
)

GERG 2004 Monograph equation 7.52:

\begin{eqnarray*} \left( \frac{\partial}{\partial x_j}\left(n\left(\frac{\partial\alpha^r}{\partial n_i}\right)_{T,V,n_j}\right)\right)_{\delta,\tau,x_i} &=& \delta\alpha_{\delta x_j}^{r}\left[ 1-\frac{1}{\rho_r}\cdot n\left(\frac{\partial \rho_r}{\partial n_i}\right)_{n_j}\right] \\ &-& \delta\alpha_{\delta}^{r}\frac{1}{\rho_r}\left[ \left(\frac{\partial}{\partial x_j}\left(n\left(\frac{\partial \rho_r}{\partial n_i}\right)_{n_j}\right)\right)_{x_i}-\frac{1}{\rho_r}\left(\frac{\partial \rho_r}{\partial x_j}\right)_{x_i}\cdot n\left(\frac{\partial \rho_r}{\partial n_i}\right)_{n_j}\right] \\ &+& \tau\alpha_{\tau x_j}^r\frac{1}{T_r}\cdot n\left(\frac{\partial T_r}{\partial n_i}\right)_{n_j}\\ &+& \tau\alpha_{\tau}^{r}\frac{1}{T_r}\left[ \left(\frac{\partial}{\partial x_j}\left(n\left(\frac{\partial T_r}{\partial n_i}\right)_{n_j}\right)\right)_{x_i}-\frac{1}{T_r}\left(\frac{\partial T_r}{\partial x_j}\right)_{x_i}\cdot n\left(\frac{\partial T_r}{\partial n_i}\right)_{n_j}\right] \\ &+& \alpha_{x_ix_j}^r-\alpha_{x_j}^r-\sum_{m=1}^Nx_m\alpha_{x_jx_m}^r \end{eqnarray*}

Definition at line 1124 of file Mixtures.cpp.

double Mixture::d_ndphirdni_dxj__constT_V_xi ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i,
int  j 
)

Gernert Equation 3.118 Catch test provided

Definition at line 1037 of file Mixtures.cpp.

double Mixture::ddelta_dxj__constT_V_xi ( double  tau,
double  delta,
const std::vector< double > &  x,
int  j 
)

Gernert Equation 3.121 Catch test provided

Definition at line 1044 of file Mixtures.cpp.

double Mixture::dgdbeta_RachfordRice ( const std::vector< double > &  z,
const std::vector< double > &  lnK,
double  beta 
)

Objective function from Rachford-Rice (Not to be confused with the Gibbs function)

\[ \frac{dg}{d\beta} = \sum_i z_i\frac{(K_i-1)^2}{(1-\beta+\beta K_i)^2}; \]

Parameters
zBulk mole fractions [-]
lnKLogarithm of the K factor [-]
betaMolar fraction in the gas phase [-]

Definition at line 1558 of file Mixtures.cpp.

double Mixture::dln_fugacity_coefficient_dp__constT_n ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

The derivative term

\[ \left(\frac{\partial \ln \phi_i}{\partial p} \right)_{T,\bar n} = \frac{\hat v_i}{RT}-\frac{1}{p} \]

GERG 2004 Monograph Eqn. 7.30

Definition at line 982 of file Mixtures.cpp.

double Mixture::dln_fugacity_coefficient_dT__constp_n ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

The derivative term

\[ \left(\frac{\partial \ln \phi_i}{\partial T} \right)_{p,\bar n} = \left(\frac{\partial^2n\alpha^r}{\partial T\partial n_i} \right)_{V,n_j} + \frac{1}{T}-\frac{\hat v}{RT}\left(\frac{\partial p}{\partial T}\right)_{V,\bar n} \]

GERG 2004 Monograph Eqn. 7.29

Definition at line 971 of file Mixtures.cpp.

double Mixture::dln_fugacity_coefficient_dT__constrho ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

Derivative of the natural logarithm of the fugacity coefficient with respect to T

Definition at line 953 of file Mixtures.cpp.

double Mixture::dln_fugacity_coefficient_dxj__constT_p_xi ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i,
int  j 
)

Gernert Equation 3.115.

Definition at line 999 of file Mixtures.cpp.

double Mixture::dnphir_dni__constT_V_nj ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

Definition at line 960 of file Mixtures.cpp.

double Mixture::dpdT__constV_n ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

The derivative term

\[ \left(\frac{\partial p}{\partial T} \right)_{V,\bar n} = \rho R(1+\delta \alpha_{\delta}^r-\delta \tau \alpha^r_{\delta\tau}) \]

GERG 2004 Monograph equation 7.61

Definition at line 1060 of file Mixtures.cpp.

double Mixture::dpdxj__constT_V_xi ( double  tau,
double  delta,
const std::vector< double > &  x,
int  j 
)

Gernert Equation 3.130 Catch test provided

Definition at line 1009 of file Mixtures.cpp.

double Mixture::dphir_dDelta ( double  tau,
double  delta,
const std::vector< double > &  x 
)

Definition at line 1214 of file Mixtures.cpp.

double Mixture::dphir_dTau ( double  tau,
double  delta,
const std::vector< double > &  x 
)

Definition at line 1230 of file Mixtures.cpp.

double Mixture::dphir_dxi ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

Definition at line 1198 of file Mixtures.cpp.

double Mixture::dphir_dxj__constT_V_xi ( double  tau,
double  delta,
const std::vector< double > &  x,
int  j 
)

Gernert Equation 3.119 Catch test provided

Definition at line 1032 of file Mixtures.cpp.

double Mixture::dtau_dxj__constT_V_xi ( double  tau,
double  delta,
const std::vector< double > &  x,
int  j 
)

Gernert Equation 3.122 Catch test provided

Definition at line 1050 of file Mixtures.cpp.

double Mixture::fugacity ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

Returns the fugacity for the given component for the given total reduced density and reciprocal reduced temperature.

Definition at line 938 of file Mixtures.cpp.

double Mixture::g_RachfordRice ( const std::vector< double > &  z,
const std::vector< double > &  lnK,
double  beta 
)

Objective function from Rachford-Rice (Not to be confused with the Gibbs function)

\[ g(\beta) = \sum_i(y_i-x_i) = \sum_i z_i\frac{K_i-1}{1-\beta+\beta K_i} \]

Parameters
zBulk mole fractions [-]
lnKLogarithm of the K factor [-]
betaMolar fraction in the gas phase [-]

Definition at line 1547 of file Mixtures.cpp.

void Mixture::initialize ( )

END OF INITIALIZATION END OF INITIALIZATION END OF INITIALIZATION END OF INITIALIZATION

Definition at line 350 of file Mixtures.cpp.

double Mixture::ln_fugacity_coefficient ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

Natural logarithm of the fugacity coefficient

Definition at line 949 of file Mixtures.cpp.

void Mixture::load_excess_values ( int  i,
int  j 
)

Load the excess parameters (departure function parameters)

Parameters
i0-based index of first component
j0-based index of second component

Definition at line 184 of file Mixtures.cpp.

std::map< std::string, double > Mixture::load_reducing_values ( int  i,
int  j 
)

Load the reducing parameters

Parameters
i0-based index of first component
j0-based index of second component

Definition at line 78 of file Mixtures.cpp.

double Mixture::nd2nphirdnidnj__constT_V ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i,
int  j 
)

GERG 2004 Monograph equation 7.41:

\[ n\left(\frac{\partial^2n\alpha^r}{\partial n_i \partial n_j} \right)_{T,V} = n\left( \frac{\partial}{\partial n_j}\left(\frac{\partial n\alpha^r}{\partial n_i}\right)_{T,V,n_j}\right)_{T,V,n_i} \]

and GERG 2004 Monograph equation 7.46:

\[ n\left( \frac{\partial}{\partial n_j}\left(\frac{\partial n\alpha^r}{\partial n_i}\right)_{T,V,n_j}\right)_{T,V,n_i} = n\left( \frac{\partial \alpha^r}{\partial n_j}\right)_{T,V,n_i} + n\left( \frac{\partial}{\partial n_j}\left(n\left(\frac{\partial \alpha^r}{\partial n_i}\right)_{T,V,n_j} \right) \right)_{T,V,n_i} \]

GERG 2004 Monograph equation 7.47:

\begin{eqnarray*} n\left( \frac{\partial}{\partial n_j}\left(n\left(\frac{\partial \alpha^r}{\partial n_i}\right)_{T,V,n_j} \right) \right)_{T,V,n_i} &=& \left( \frac{\partial}{\partial \delta}\left(n\left(\frac{\partial\alpha^r}{\partial n_i}\right)_{T,V,n_j}\right)\right)_{\tau,\bar x}\cdot n\left(\frac{\partial\delta}{\partial n_j}\right)_{T,V,n_i}\\ &+& \left( \frac{\partial}{\partial \tau}\left(n\left(\frac{\partial\alpha^r}{\partial n_i}\right)_{T,V,n_j}\right)\right)_{\tau,\bar x}\cdot n\left(\frac{\partial\tau}{\partial n_j}\right)_{T,V,n_i}\\ &+& \left( \frac{\partial}{\partial x_j}\left(n\left(\frac{\partial\alpha^r}{\partial n_i}\right)_{T,V,n_j}\right)\right)_{\delta,\tau,x_i}-\sum_{k=1}^{N}x_k \left( \frac{\partial}{\partial x_k}\left(n\left(\frac{\partial\alpha^r}{\partial n_i}\right)_{T,V,n_j}\right)\right)_{\delta,\tau,x_i}\\ \end{eqnarray*}

Definition at line 1141 of file Mixtures.cpp.

double Mixture::nddeltadni__constT_V_nj ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

\[ n\left(\frac{\partial \delta}{\partial n_i} \right)_{T,V,n_j} = \delta - \frac{\delta}{\rho_r}\cdot n\left(\frac{\partial \rho_r}{\partial n_i} \right)_{n_j} \]

GERG 2004 Monograph equation 7.48

Definition at line 1114 of file Mixtures.cpp.

double Mixture::ndln_fugacity_coefficient_dnj__constT_p ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i,
int  j 
)

\[ n\left(\frac{\partial \ln \phi_i}{\partial n_j}\right)_{T,p} = n\left(\frac{\partial^2n\alpha^r}{\partial n_j \partial n_i} \right)_{T,V}+1+\frac{n}{RT}\frac{\left(\frac{\partial p}{\partial n_j}\right)_{T,V,n_i}\left(\frac{\partial p}{\partial n_i}\right)_{T,V,n_j}}{\left(\frac{\partial p}{\partial V}\right)_{T,\bar n}} \]

which is also equal to

\[ n\left(\frac{\partial \ln \phi_i}{\partial n_j}\right)_{T,p} = n\left(\frac{\partial^2n\alpha^r}{\partial n_j \partial n_i} \right)_{T,V}+1-\frac{\hat v_i}{RT}\left[n\left(\frac{\partial p}{\partial n_j}\right)_{T,V,n_i}\right] \]

GERG 2004 Monograph Equation 7.31

Definition at line 1108 of file Mixtures.cpp.

double Mixture::ndpdni__constT_V_nj ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

The derivative term

\[ n\left(\frac{\partial p}{\partial n_i} \right)_{T,V,n_j} = \rho RT\left[1+\delta\alpha_{\delta}^r\left[2- \frac{1}{\rho_r}\cdot n\left( \frac{\partial \rho_r}{\partial n_i}\right)_{n_j}\right] +\delta\cdot n\left(\frac{\partial\alpha_{\delta}^r}{\partial n_i}\right)_{T,V,n_j}\right] \]

GERG 2004 Monograph equation 7.63

Definition at line 1074 of file Mixtures.cpp.

double Mixture::ndpdV__constT_n ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

The derivative term

\[ n\left(\frac{\partial p}{\partial V} \right)_{T,\bar n} = -\rho^2 RT(1+2\delta \alpha_{\delta}^r+\delta^2\alpha^r_{\delta\delta}) \]

GERG 2004 Monograph equation 7.62

Definition at line 1066 of file Mixtures.cpp.

double Mixture::ndphir_dni__constT_V_nj ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

The derivative term

\[ n\left(\frac{\partial \phi^r}{\partial n_i} \right)_{T,V,n_j} \]

which is equal to

\begin{eqnarray*} n\left(\frac{\partial \phi^r}{\partial n_i} \right)_{T,V,n_j} &=& \delta \phi^r_{\delta}\left[ 1-\frac{1}{\rho_r}\left[\left(\frac{\partial \rho_r}{\partial x_i}\right)_{x_j} - \sum_{k=1}^N x_k\left(\frac{\partial \rho_r}{\partial x_k}\right)_{x_j} \right]\right]\\ && +\tau \phi^r_{\tau}\frac{1}{T_r}\left[\left(\frac{\partial T_r}{\partial x_i}\right)_{x_j} - \sum_{k=1}^N x_k\left(\frac{\partial T_r}{\partial x_k}\right)_{x_j} \right]\\ && +\phi^r_{x_i}-\sum_{k=1}^{N}x_k\phi^r_{x_k} \end{eqnarray*}

See Table B4, Kunz, JCED, 2012 for the original term and the subsequent substitutions

Definition at line 1092 of file Mixtures.cpp.

double Mixture::ndtaudni__constT_V_nj ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

\[ n\left(\frac{\partial \tau}{\partial n_i} \right)_{T,V,n_j} = \frac{\tau}{T_r}\cdot n\left(\frac{\partial T_r}{\partial n_i} \right)_{n_j} \]

GERG 2004 Monograph equation 7.49

Definition at line 1119 of file Mixtures.cpp.

double Mixture::partial_molar_volume ( double  tau,
double  delta,
const std::vector< double > &  x,
int  i 
)

The partial molar volume

\[ \hat v_i = \left( \frac{\partial V}{\partial n_i}\right)_{T,p,n_j} = \frac{-\left(\dfrac{\partial p}{\partial n_i}\right)_{T,V,n_j}}{\left(\dfrac{\partial p}{\partial V}\right)_{T,\bar n}} \]

from GERG monograph eqn 7.32

Definition at line 977 of file Mixtures.cpp.

double Mixture::phir ( double  tau,
double  delta,
const std::vector< double > &  x 
)

Definition at line 1194 of file Mixtures.cpp.

double Mixture::Rbar ( const std::vector< double > &  x)

Definition at line 310 of file Mixtures.cpp.

double Mixture::rhobar_pengrobinson ( double  T,
double  p,
const std::vector< double > &  x,
int  solution 
)

Calculate the mixture molar density based on the use of the Peng-Robinson equation of state

Definition at line 1493 of file Mixtures.cpp.

double Mixture::rhobar_Tpz ( double  T,
double  p,
const std::vector< double > &  z,
double  rhobar0 
)

Density as a function of T,p,z

Parameters
TTemperature [K]
pPressure [Pa]
zBulk mole fractions [-]
rhobar0Guess value for molar density [mol/m^3]

Definition at line 1277 of file Mixtures.cpp.

double Mixture::saturation_p ( double  beta,
double  p,
const std::vector< double > &  z,
std::vector< double > &  x,
std::vector< double > &  y 
)

Definition at line 1363 of file Mixtures.cpp.

double Mixture::saturation_p_preconditioner ( double  p,
const std::vector< double > &  z 
)

Definition at line 1329 of file Mixtures.cpp.

double Mixture::saturation_p_Wilson ( double  beta,
double  p,
const std::vector< double > &  z,
double  T_guess,
std::vector< double > &  K 
)

Definition at line 1348 of file Mixtures.cpp.

void Mixture::test ( )

Definition at line 391 of file Mixtures.cpp.

void Mixture::TpzFlash ( double  T,
double  p,
const std::vector< double > &  z,
double &  rhobar,
std::vector< double > &  x,
std::vector< double > &  y 
)

Temperature-pressure-bulk mole fraction flash calculation

Parameters
TTemperature [K]
pPressure [Pa]
zBulk mole fractions [-]
rhobarMolar density [mol/m^3]
xLiquid mole fractions [-] (if saturated)
yVapor mole fractions [-] (if saturated)

Definition at line 1400 of file Mixtures.cpp.

double Mixture::Wilson_lnK_factor ( double  T,
double  p,
int  i 
)

Returns the natural logarithm of K for component i using the method from Wilson as in

\[ \ln K_i = \ln\left(\frac{p_{c,i}}{p}\right)+5.373(1+\omega_i)\left(1-\frac{T_{c,i}}{T}\right) \]

Parameters
TTemperature [K]
pPressure [Pa]
iIndex of component [-]

Definition at line 931 of file Mixtures.cpp.

void Mixture::x_and_y_from_K ( double  beta,
const std::vector< double > &  K,
const std::vector< double > &  z,
std::vector< double > &  x,
std::vector< double > &  y 
)

Definition at line 1482 of file Mixtures.cpp.

Member Data Documentation

PhaseEnvelope Mixture::Envelope

Definition at line 459 of file Mixtures.h.

unsigned int Mixture::N

Definition at line 446 of file Mixtures.h.

NewtonRaphsonVLE Mixture::NRVLE

Definition at line 462 of file Mixtures.h.

ExcessTerm* Mixture::pExcess

Definition at line 456 of file Mixtures.h.

std::vector<Fluid *> Mixture::pFluids

Definition at line 454 of file Mixtures.h.

ReducingFunction* Mixture::pReducing

Definition at line 455 of file Mixtures.h.

ResidualIdealMixture* Mixture::pResidualIdealMix

Definition at line 457 of file Mixtures.h.

Definition at line 461 of file Mixtures.h.


The documentation for this class was generated from the following files: