Cubic equations of state¶

Redlich Kwong EOS¶

The base-class is specialised into two sub-classes to address in particular liquid and gas phases. As such, the base-class is disabled from direct instantiation.

Base class¶

class simu.app.thermo.contributions.cubic.rk.RedlichKwongEOS¶

Bases: ThermoContribution, ABC

This contribution implements a general Redlich-Kwong equation of state with Peneloux volume translation:

\[p = \frac{N\,R\,T}{V - B + C} - \frac{A}{(V + C)\,(V + B + C)}\]

The following properties need to be provided by upstream contributions:

Property	Symbol	UOM
ceos_a	\(A\)	J * m³
ceos_b	\(B\)	m³
ceos_c	\(C\)	m³	(optional)

Care is to be taken when utilising a temperature-dependent \(C\) contribution, as doing so can have significant effects on the calorimetric properties. If ceos_c is not provided, the contribution is assumed zero.

As such, there are no further model parameters to be provided at this point. The residual Helmholtz function is

\[A^\mathrm{res} = \int\limits_V^\infty p - \frac{N\,R\,T}{V} \mathrm{d}V = N\,R\,T\,\ln \frac{V}{V + C - B} + \frac{A}{B}\,\ln \frac{V + C}{V + C + B}\]

The explicitly implemented temperature derivative is

\[\begin{split}-S^\mathrm{res} &= N\,R\,\left [ \ln\frac{V}{V - B + C} + T\,\frac{B_T - C_T}{V - B + C} \right ]\\& + \frac1B\left (A_T - \frac{A}{B}\,B_T\right )\, \ln \frac{V + C}{V + B + C} + \frac{A}{B}\left [ \frac{C_T}{V + C} - \frac{B_T + C_T}{V + B + C} \right ]\end{split}\]

The volume derivative is the negative residual pressure:

\[-p^\mathrm{res} = N\,R\,T\, \left [ \frac1V - \frac1{V - B + C}\right ] + \frac{A}{(V + C)\,(V + B + C)}\]

The derivative with respect to molar quantities is

\[\begin{split}\mu_i^\mathrm{res} &= R\,T\,\left [ \ln\frac{V}{V - B + C} + N\,\frac{B_{n,i} - C_{n,i}}{V - B + C} \right ]\\& + \frac1B\left (A_{n,i} - \frac{A}{B}\,B_{n,i}\right )\, \ln \frac{V + C}{V + B + C} + \frac{A}{B}\left [ \frac{C_{n,i}}{V + C} - \frac{B_{n,i} + C_{n,i}}{V + B + C} \right ]\end{split}\]

The contribution updates are

\[\begin{split}S &\leftarrow S + S^\mathrm{res}\\ p &\leftarrow p + p^\mathrm{res}\\ \mu_i &\leftarrow \mu_i + \mu_i^\mathrm{res}\\\end{split}\]

Important

The original publication on the volume correction of EOS [PRF82] suggests to shift both V and B, not only V. The derivation of this proposal is not given, and a numerical experiment shows that doing so drastically impacts the equilibrium and the density predictions in a negative way. We therefore stick to the approach taken by [Aspen Technology, Inc.01], namely a pure volume translation. For this translation, it can easily be shown, as by [PRF82] that equilibrium conditions are not affected:

Consider the chemical potential derived from the Helmholtz function

\[A = A(T, V, \vec n)\quad\Rightarrow\quad \vec \mu = \left . \frac{\partial A}{\partial \vec n} \right |_{T, V}\]

If we translate \(\tilde V = V - \vec c\cdot\vec n\), we can derive

\[\tilde A = \tilde A(T, \tilde V, \vec n)\quad\Rightarrow\quad \vec{\tilde \mu} = \left . \frac{\partial \tilde A}{\partial \vec n} \right |_{T, \tilde V} = \left . \frac{\partial \tilde A}{\partial \vec n} \right |_{T, V} + \left . \frac{\partial \tilde A}{\partial \tilde V} \right |_{T, n}\cdot \left . \frac{\partial \tilde V}{\partial \vec n} \right |_{T} = \vec \mu - p\,\vec c\]

Hence at given temperature and pressure, the chemical potential only shifts by a constant, not affecting any equilibrium constraints.

Liquid root¶

class simu.app.thermo.contributions.cubic.rk.RedlichKwongEOSLiquid¶

Bases: RedlichKwongEOS

As a subclass of RedlichKwongEOS, this entity specialises on describing liquid (and super-critical) phases. The distinct element is the initialisation.

Gas root¶

class simu.app.thermo.contributions.cubic.rk.RedlichKwongEOSGas¶

Bases: RedlichKwongEOS

As a subclass of RedlichKwongEOS, this entity specialises on describing gas (and super-critical) phases. The distinct elements is the initialisation.

A-function¶

class simu.app.thermo.contributions.cubic.rk.RedlichKwongAFunction¶

Bases: ThermoContribution

Given critical temperature T_c (\(T_{c,i}\)) and pressure p_c (\(p_{c,i}\)), this contribution scales the \(\alpha\)-function alpha (\(\alpha_i\)) to define the \(a\)-contribution rk_a_i (\(a_i\)) for the individual species. It is

\[a_i = \alpha_i\,\Omega_a\,\frac{R^2\,T_{c,i}^2}{p_{c,i}} \quad\text{with}\quad \Omega_a = \frac19\,(2^{1/3} - 1)^{-1}\]

B-function¶

class simu.app.thermo.contributions.cubic.rk.RedlichKwongBFunction¶

Bases: ThermoContribution

Given critical temperature T_c (\(T_{c,i}\)) and pressure p_c (\(p_{c,i}\)), this contribution calculates the \(b\)-contribution ceos_b_i for the individual species. It is

\[b_i = \Omega_b\,\frac{R\,T_{c,i}}{p_{c,i}} \quad\text{with}\quad \Omega_b = \frac13\,(2^{1/3} - 1)\]

RedlichKwongMFactor¶

class simu.app.thermo.contributions.cubic.rk.RedlichKwongMFactor¶

Bases: ThermoContribution

This contribution calculates the Redlich Kwong m-factor that is used in various alpha-functions. Based on provided acentric factors omega (\(\omega_i\)), it calculates m_factor (\(m_i\)) as

\[m_i = 0.48 + (1.57 - 0.176\,\omega_i)\,\omega_i\]

BostonMathiasMFactor¶

class simu.app.thermo.contributions.cubic.rk.BostonMatthiasMFactor¶

Bases: ThermoContribution

This contribution calculates the Boston-Matthias m-factor that is used in various alpha-functions. Based on provided acentric factors omega (\(\omega_i\)), it calculates m_factor (\(m_i\)) as

\[m_i = 0.48508 + (1.55171 - 0.15613\,\omega_i)\,\omega_i\]

Critical parameters¶

class simu.app.thermo.contributions.cubic.core.CriticalParameters¶

Bases: ThermoContribution

This class does not perform any calculations, but provides the basic critical parameters as a basis for the typical equation of state contributions.

The following parameters need to be provided (all as species vector):

Property	Symbol	Description
T_c	\(T_c\)	Critical temperatures [K]
p_c	\(p_c\)	Critical pressure[Pa]
omega	\(\omega\)	Acentric factor [-]

The same symbols will just be published as intermediate results for the actual model contributions to be consumed.

Non-symmetric mixing rule¶

class simu.app.thermo.contributions.cubic.core.NonSymmetricMixingRule¶

Bases: ThermoContribution

The mixing rule combines the pure species a-contributions a_i (\(a_i\)) into a lumped a property.

The contribution requires an option dictionary with the following entries:

target: The name of the property \(a\).

source: The name of the property \(a_i\). If omitted, it will be generated as the target name with _i appended.

Implemented with sparse binary interaction coefficient structure, the contribution includes a symmetric interaction (\(k\)) and an antisymmetric contribution (\(l\)).

\[a = \sum_{ij} \sqrt{a_i(T)\,a_j(T)} \left [ n_i\,n_j - 2\,n_i\,n_j\,k_{ij}(T) - \frac2N\,(n_i^2\,n_j - n_i\,n_j^2)\,l_{ij}(T) \right ]\]

These dimensionless interaction coefficients can be a function of temperature according to

\[\begin{split}k_{ij}(T) &:= k_{0,ij} + k_{1,ij}\,\left (1 - \frac{T}{T_\mathrm{ref}} \right ) + k_{2,ij}\,\left (1 - \frac{T_\mathrm{ref}}{T} \right )\\ l_{ij}(T) &:= l_{0,ij} + l_{1,ij}\,\left (1 - \frac{T}{T_\mathrm{ref}} \right ) + l_{2,ij}\,\left (1 - \frac{T_\mathrm{ref}}{T} \right )\end{split}\]

Each pair of species can only be defined once. The pre-factor 2 is used to yield the same parameter values as if \(k\) was a symmetric matrix and \(l\) was antisymmetric.

Despite what above equation suggests based on the double sum, the complexity of this contribution in terms of both memory and runtime is linear in the number of species and linear in the number of non-zero interaction parameters. This is achieved using the relationship

\[\sum_{ij} \sqrt{a_i(T)\,a_j(T)} n_i\,n_j = \left (\sum_i \sqrt{a_i(T)}\,n_i \right )^2\]

Linear mixing rule¶

class simu.app.thermo.contributions.cubic.core.LinearMixingRule¶

Bases: ThermoContribution

This linear mixing rule represents any contribution that computes a lumped quantity as weighted sum over the molar quantities:

\[x = \sum_i x_i\,n_i\]

The contribution requires an option dictionary with the following entries:

target: The name of the property \(x\).

source: The name of the property \(x_i\). If omitted, it will be generated as the target name with _i appended.

SRKAlphaFunction¶

class simu.app.thermo.contributions.cubic.core.SRKAlphaFunction¶

Bases: ThermoContribution

This contribution represents the SRK alpha function without any extrapolation for super-critical temperatures.

The following properties need to be provided upstream:

Property

Symbol

Description

T

\(T\)

Actual temperatures [K]

T_c

\(T_c\)

Critical temperatures [K]

m_factor

\(m_i\)

m-factor as function of acentric factor [-]

\[\alpha_i^{\frac12} = 1 + m_i\,(1 - \sqrt{T/T_{c,i}})\]

The calculated vector is provided as a property called alpha

Boston-Mathias alpha-function¶

class simu.app.thermo.contributions.cubic.core.BostonMathiasAlphaFunction¶

Bases: ThermoContribution

This contribution represents the Mathias alpha function with the Boston-Mathias extrapolation for supercritical temperatures.

The following properties need to be provided upstream:

Property

Symbol

Description

T

\(T\)

Actual temperatures [K]

T_c

\(T_c\)

Critical temperatures [K]

m_factor

\(m_i\)

m-factor as function of acentric factor [-]

Additionally, the contribution requires a polar parameter \(\eta_i\), named eta. We define the root of the reduced temperature as \(\tau_i := \sqrt{T/T_{c,i}}\). Then, we define for \(\tau_i \le 1\):

\[\alpha_i^{\frac12} = 1 + m_i\,(1 - \tau_i) - \eta_i\,(1 - \tau_i)(0.7 - \tau_i^2)\]

Important

While the paper [Mat83] is currently unavailable to me, I must recognize that all but [Aspen Technology, Inc.01] quote the expression differently, namely:

\[\alpha_i^{\frac12} = 1 + m_i\,(1 - \tau_i) - \eta_i\, (1 - \boxed{\tau_i^2})(0.7 - \tau_i^2)\]

Likely, there was a misprint in the original paper [Mat83], as the parameter \(\eta_i\) shall not alter the impact of the acentric factor \(\omega_i\) at \(T = 0.7\, T_{c,i}\), as this point is defined to match the saturation pressure via

\[\omega_i = -\log_{10} \frac{p^\mathrm{sat}_i(0.7\cdot T_{c,i})}{p_{c, i}}\]

The reason why AspenTech could and did use the correct form is that the author P. Mathias was working for the company at the time of these developments, and the co-workers were not dependent on the journal publication.

As described in Appendix (Extensions of alpha-functions to super-critical regions), the extrapolation into the super-critical region is implemented as

\[\alpha_i^{\frac12} = \left [\frac{c}{d}(1-\tau^{d})\right ] \quad\text{with}\quad c = m + 0.3\eta \quad\text{and}\quad d = 1 + \frac{4\,\eta}{c} + c\]

The calculated vector is provided as a property called alpha

Appendix¶

Extensions of alpha-functions to super-critical regions

Property	Symbol	Description
T	\(T\)	Actual temperatures [K]
T_c	\(T_c\)	Critical temperatures [K]
m_factor	\(m_i\)	m-factor as function of acentric factor [-]