Skip to content

polykin.thermo.eos¤

This module implements equations of state (EOS) for gas and liquid mixtures.

Virial ¤

Virial equation of state truncated to the second coefficient.

This EOS is based on the following \(P(v,T)\) relationship:

\[ P = \frac{R T}{v - B_m} \]

where \(P\) is the pressure, \(T\) is the temperature, \(v\) is the molar volume, \(B_m\) is the mixture second virial coefficient. The parameter \(B_m=B_m(T,y)\) is estimated based on the critical properties of the pure components and the mixture composition \(y\).

Important

This equation is an improvement over the ideal gas model, but it should only be used up to moderate pressures such that \(v/v_c > 2\).

References

  • RC Reid, JM Prausniz, and BE Poling. The properties of gases & liquids 4th edition, 1986, p. 37, 40, 80, 82.
PARAMETER DESCRIPTION
Tc

Critical temperatures of all components. Unit = K.

TYPE: float | FloatVectorLike

Pc

Critical pressures of all components. Unit = Pa.

TYPE: float | FloatVectorLike

Zc

Critical compressibility factors of all components.

TYPE: float | FloatVectorLike

w

Acentric factors of all components.

TYPE: float | FloatVectorLike

name

Name.

TYPE: str DEFAULT: ''

Source code in src/polykin/thermo/eos/virial.py 
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
class Virial(GasEoS):
    r"""Virial equation of state truncated to the second coefficient.

    This EOS is based on the following $P(v,T)$ relationship:

    $$ P = \frac{R T}{v - B_m} $$

    where $P$ is the pressure, $T$ is the temperature, $v$ is the molar
    volume, $B_m$ is the mixture second virial coefficient. The parameter
    $B_m=B_m(T,y)$ is estimated based on the critical properties of the pure
    components and the mixture composition $y$.

    !!! important

        This equation is an improvement over the ideal gas model, but it
        should only be used up to moderate pressures such that $v/v_c > 2$.

    **References**

    *   RC Reid, JM Prausniz, and BE Poling. The properties of gases &
        liquids 4th edition, 1986, p. 37, 40, 80, 82.

    Parameters
    ----------
    Tc : float | FloatVectorLike
        Critical temperatures of all components. Unit = K.
    Pc : float | FloatVectorLike
        Critical pressures of all components. Unit = Pa.
    Zc : float | FloatVectorLike
        Critical compressibility factors of all components.
    w : float | FloatVectorLike
        Acentric factors of all components.
    name : str
        Name.
    """

    Tc: FloatVector
    Pc: FloatVector
    Zc: FloatVector
    w: FloatVector

    def __init__(self,
                 Tc: Union[float, FloatVectorLike],
                 Pc: Union[float, FloatVectorLike],
                 Zc: Union[float, FloatVectorLike],
                 w: Union[float, FloatVectorLike],
                 name: str = ''
                 ) -> None:

        Tc, Pc, Zc, w = \
            convert_FloatOrVectorLike_to_FloatVector([Tc, Pc, Zc, w])

        N = len(Tc)
        super().__init__(N, name)
        self.Tc = Tc
        self.Pc = Pc
        self.Zc = Zc
        self.w = w

    def Bm(self,
           T: float,
           y: FloatVector
           ) -> float:
        r"""Calculate the second virial coefficient of the mixture.

        $$ B_m = \sum_i \sum_j y_i y_j B_{ij} $$

        **References**

        *   RC Reid, JM Prausniz, and BE Poling. The properties of gases &
            liquids 4th edition, 1986, p. 79.

        Parameters
        ----------
        T : float
            Temperature. Unit = K.
        y : FloatVector
            Mole fractions of all components. Unit = mol/mol.

        Returns
        -------
        float
            Mixture second virial coefficient, $B_m$. Unit = m³/mol.
        """
        return quadratic_mixing(y, self.Bij(T))

    @functools.cache
    def Bij(self,
            T: float,
            ) -> FloatSquareMatrix:
        r"""Calculate the matrix of interaction virial coefficients.

        The calculation is handled by [`B_mixture`](B_mixture.md).

        Parameters
        ----------
        T : float
            Temperature. Unit = K.

        Returns
        -------
        FloatSquareMatrix
            Matrix of interaction virial coefficients, $B_{ij}$.
            Unit = m³/mol.
        """
        return B_mixture(T, self.Tc, self.Pc, self.Zc, self.w)

    def Z(self,
          T: float,
          P: float,
          y: FloatVector
          ) -> float:
        r"""Calculate the compressibility factor of the fluid.

        $$ Z = 1 + \frac{B_m P}{R T} $$

        where $P$ is the pressure, $T$ is the temperature, and $B_m=B_m(T,y)$
        is the mixture second virial coefficient.

        **References**

        *   RC Reid, JM Prausniz, and BE Poling. The properties of gases &
            liquids 4th edition, 1986, p. 37.

        Parameters
        ----------
        T : float
            Temperature. Unit = K.
        P : float
            Pressure. Unit = Pa.
        y : FloatVector
            Mole fractions of all components. Unit = mol/mol.

        Returns
        -------
        float
            Compressibility factor.
        """
        Bm = self.Bm(T, y)
        return 1.0 + Bm*P/(R*T)

    def P(self,
          T: float,
          v: float,
          y: FloatVector
          ) -> float:
        r"""Calculate the pressure of the fluid.

        Parameters
        ----------
        T : float
            Temperature. Unit = K.
        v : float
            Molar volume. Unit = m³/mol.
        y : FloatVector
            Mole fractions of all components. Unit = mol/mol.

        Returns
        -------
        float
            Pressure. Unit = Pa.
        """
        Bm = self.Bm(T, y)
        return R*T/(v - Bm)

    def phi(self,
            T: float,
            P: float,
            y: FloatVector
            ) -> FloatVector:
        r"""Calculate the fugacity coefficients of all components.

        For each component, the fugacity coefficient is given by:

        $$
        \ln \hat{\phi}_i = \left(2\sum_j {y_jB_{ij}} -B_m \right)\frac{P}{RT}
        $$

        where $P$ is the pressure, $T$ is the temperature, $B_{ij}$ is the 
        matrix of interaction virial coefficients, $B_m$ is the second virial
        coefficient of the mixture, and $y_i$ is the mole fraction.

        **References**

        *   RC Reid, JM Prausniz, and BE Poling. The properties of gases &
            liquids 4th edition, 1986, p. 145.

        Parameters
        ----------
        T : float
            Temperature. Unit = K.
        P : float
            Pressure. Unit = Pa.
        y : FloatVector
            Mole fractions of all components. Unit = mol/mol.

        Returns
        -------
        FloatVector
            Fugacity coefficients of all components.
        """
        B = self.Bij(T)
        Bm = self.Bm(T, y)
        return exp((2*dot(B, y) - Bm)*P/(R*T))

    def DA(self,
           T: float,
           V: float,
           n: FloatVector,
           v0: float
           ) -> float:
        nT = n.sum()
        y = n/nT
        Bm = self.Bm(T, y)
        return -nT*R*T*log((V - nT*Bm)/(nT*v0))

Bij cached ¤

Bij(T: float) -> FloatSquareMatrix

Calculate the matrix of interaction virial coefficients.

The calculation is handled by B_mixture.

PARAMETER DESCRIPTION
T

Temperature. Unit = K.

TYPE: float

RETURNS DESCRIPTION
FloatSquareMatrix

Matrix of interaction virial coefficients, \(B_{ij}\). Unit = m³/mol.
Source code in src/polykin/thermo/eos/virial.py 
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
@functools.cache
def Bij(self,
        T: float,
        ) -> FloatSquareMatrix:
    r"""Calculate the matrix of interaction virial coefficients.

    The calculation is handled by [`B_mixture`](B_mixture.md).

    Parameters
    ----------
    T : float
        Temperature. Unit = K.

    Returns
    -------
    FloatSquareMatrix
        Matrix of interaction virial coefficients, $B_{ij}$.
        Unit = m³/mol.
    """
    return B_mixture(T, self.Tc, self.Pc, self.Zc, self.w)

Bm ¤

Bm(T: float, y: FloatVector) -> float

Calculate the second virial coefficient of the mixture.

\[ B_m = \sum_i \sum_j y_i y_j B_{ij} \]

References

  • RC Reid, JM Prausniz, and BE Poling. The properties of gases & liquids 4th edition, 1986, p. 79.
PARAMETER DESCRIPTION
T

Temperature. Unit = K.

TYPE: float

y

Mole fractions of all components. Unit = mol/mol.

TYPE: FloatVector

RETURNS DESCRIPTION
float

Mixture second virial coefficient, \(B_m\). Unit = m³/mol.
Source code in src/polykin/thermo/eos/virial.py 
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
def Bm(self,
       T: float,
       y: FloatVector
       ) -> float:
    r"""Calculate the second virial coefficient of the mixture.

    $$ B_m = \sum_i \sum_j y_i y_j B_{ij} $$

    **References**

    *   RC Reid, JM Prausniz, and BE Poling. The properties of gases &
        liquids 4th edition, 1986, p. 79.

    Parameters
    ----------
    T : float
        Temperature. Unit = K.
    y : FloatVector
        Mole fractions of all components. Unit = mol/mol.

    Returns
    -------
    float
        Mixture second virial coefficient, $B_m$. Unit = m³/mol.
    """
    return quadratic_mixing(y, self.Bij(T))

DA ¤

DA(T: float, V: float, n: FloatVector, v0: float) -> float

Calculate the departure of Helmholtz energy.

\[ A(T,V,n) - A^{\circ}(T,V,n)\]
PARAMETER DESCRIPTION
T

Temperature. Unit = K.

TYPE: float

V

Volume. Unit = m³.

TYPE: float

n

Mole amounts of all components. Unit = mol.

TYPE: FloatVector

v0

Molar volume in reference state. Unit = m³/mol.

TYPE: float

RETURNS DESCRIPTION
float

Helmholtz energy departure, \(A - A^{\circ}\). Unit = J.
Source code in src/polykin/thermo/eos/virial.py 
231
232
233
234
235
236
237
238
239
240
def DA(self,
       T: float,
       V: float,
       n: FloatVector,
       v0: float
       ) -> float:
    nT = n.sum()
    y = n/nT
    Bm = self.Bm(T, y)
    return -nT*R*T*log((V - nT*Bm)/(nT*v0))

N property ¤

N: int

Number of components.

P ¤

P(T: float, v: float, y: FloatVector) -> float

Calculate the pressure of the fluid.

PARAMETER DESCRIPTION
T

Temperature. Unit = K.

TYPE: float

v

Molar volume. Unit = m³/mol.

TYPE: float

y

Mole fractions of all components. Unit = mol/mol.

TYPE: FloatVector

RETURNS DESCRIPTION
float

Pressure. Unit = Pa.
Source code in src/polykin/thermo/eos/virial.py 
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
def P(self,
      T: float,
      v: float,
      y: FloatVector
      ) -> float:
    r"""Calculate the pressure of the fluid.

    Parameters
    ----------
    T : float
        Temperature. Unit = K.
    v : float
        Molar volume. Unit = m³/mol.
    y : FloatVector
        Mole fractions of all components. Unit = mol/mol.

    Returns
    -------
    float
        Pressure. Unit = Pa.
    """
    Bm = self.Bm(T, y)
    return R*T/(v - Bm)

Z ¤

Z(T: float, P: float, y: FloatVector) -> float

Calculate the compressibility factor of the fluid.

\[ Z = 1 + \frac{B_m P}{R T} \]

where \(P\) is the pressure, \(T\) is the temperature, and \(B_m=B_m(T,y)\) is the mixture second virial coefficient.

References

  • RC Reid, JM Prausniz, and BE Poling. The properties of gases & liquids 4th edition, 1986, p. 37.
PARAMETER DESCRIPTION
T

Temperature. Unit = K.

TYPE: float

P

Pressure. Unit = Pa.

TYPE: float

y

Mole fractions of all components. Unit = mol/mol.

TYPE: FloatVector

RETURNS DESCRIPTION
float

Compressibility factor.
Source code in src/polykin/thermo/eos/virial.py 
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
def Z(self,
      T: float,
      P: float,
      y: FloatVector
      ) -> float:
    r"""Calculate the compressibility factor of the fluid.

    $$ Z = 1 + \frac{B_m P}{R T} $$

    where $P$ is the pressure, $T$ is the temperature, and $B_m=B_m(T,y)$
    is the mixture second virial coefficient.

    **References**

    *   RC Reid, JM Prausniz, and BE Poling. The properties of gases &
        liquids 4th edition, 1986, p. 37.

    Parameters
    ----------
    T : float
        Temperature. Unit = K.
    P : float
        Pressure. Unit = Pa.
    y : FloatVector
        Mole fractions of all components. Unit = mol/mol.

    Returns
    -------
    float
        Compressibility factor.
    """
    Bm = self.Bm(T, y)
    return 1.0 + Bm*P/(R*T)

f ¤

f(T: float, P: float, y: FloatVector) -> FloatVector

Calculate the fugacity of all components.

For each component, the fugacity is given by:

\[ \hat{f}_i = \hat{\phi}_i y_i P \]

where \(\hat{\phi}_i(T,P,y)\) is the fugacity coefficient, \(P\) is the pressure, and \(y_i\) is the mole fraction.

PARAMETER DESCRIPTION
T

Temperature. Unit = K.

TYPE: float

P

Pressure. Unit = Pa.

TYPE: float

y

Mole fractions of all components. Unit = mol/mol.

TYPE: FloatVector

RETURNS DESCRIPTION
FloatVector

Fugacities of all components. Unit = Pa.
Source code in src/polykin/thermo/eos/base.py 
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
def f(self,
      T: float,
      P: float,
      y: FloatVector
      ) -> FloatVector:
    r"""Calculate the fugacity of all components.

    For each component, the fugacity is given by:

    $$ \hat{f}_i = \hat{\phi}_i y_i P $$

    where $\hat{\phi}_i(T,P,y)$ is the fugacity coefficient, $P$ is the 
    pressure, and $y_i$ is the mole fraction.

    Parameters
    ----------
    T : float
        Temperature. Unit = K.
    P : float
        Pressure. Unit = Pa.
    y : FloatVector
        Mole fractions of all components. Unit = mol/mol.

    Returns
    -------
    FloatVector
        Fugacities of all components. Unit = Pa.
    """
    return self.phi(T, P, y)*y*P

phi ¤

phi(T: float, P: float, y: FloatVector) -> FloatVector

Calculate the fugacity coefficients of all components.

For each component, the fugacity coefficient is given by:

\[ \ln \hat{\phi}_i = \left(2\sum_j {y_jB_{ij}} -B_m \right)\frac{P}{RT} \]

where \(P\) is the pressure, \(T\) is the temperature, \(B_{ij}\) is the matrix of interaction virial coefficients, \(B_m\) is the second virial coefficient of the mixture, and \(y_i\) is the mole fraction.

References

  • RC Reid, JM Prausniz, and BE Poling. The properties of gases & liquids 4th edition, 1986, p. 145.
PARAMETER DESCRIPTION
T

Temperature. Unit = K.

TYPE: float

P

Pressure. Unit = Pa.

TYPE: float

y

Mole fractions of all components. Unit = mol/mol.

TYPE: FloatVector

RETURNS DESCRIPTION
FloatVector

Fugacity coefficients of all components.
Source code in src/polykin/thermo/eos/virial.py 
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
def phi(self,
        T: float,
        P: float,
        y: FloatVector
        ) -> FloatVector:
    r"""Calculate the fugacity coefficients of all components.

    For each component, the fugacity coefficient is given by:

    $$
    \ln \hat{\phi}_i = \left(2\sum_j {y_jB_{ij}} -B_m \right)\frac{P}{RT}
    $$

    where $P$ is the pressure, $T$ is the temperature, $B_{ij}$ is the 
    matrix of interaction virial coefficients, $B_m$ is the second virial
    coefficient of the mixture, and $y_i$ is the mole fraction.

    **References**

    *   RC Reid, JM Prausniz, and BE Poling. The properties of gases &
        liquids 4th edition, 1986, p. 145.

    Parameters
    ----------
    T : float
        Temperature. Unit = K.
    P : float
        Pressure. Unit = Pa.
    y : FloatVector
        Mole fractions of all components. Unit = mol/mol.

    Returns
    -------
    FloatVector
        Fugacity coefficients of all components.
    """
    B = self.Bij(T)
    Bm = self.Bm(T, y)
    return exp((2*dot(B, y) - Bm)*P/(R*T))

v ¤

v(T: float, P: float, y: FloatVector) -> float

Calculate the molar volume the fluid.

\[ v = \frac{Z R T}{P} \]

where \(v\) is the molar volume, \(Z(T, P, y)\) is the compressibility factor, \(T\) is the temperature, \(P\) is the pressure, and \(y\) is the mole fraction vector.

PARAMETER DESCRIPTION
T

Temperature. Unit = K.

TYPE: float

P

Pressure. Unit = Pa.

TYPE: float

y

Mole fractions of all components. Unit = mol/mol.

TYPE: FloatVector

RETURNS DESCRIPTION
float

Molar volume of the fluid. Unit = m³/mol.
Source code in src/polykin/thermo/eos/base.py 
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
def v(self,
      T: float,
      P: float,
      y: FloatVector
      ) -> float:
    r"""Calculate the molar volume the fluid.

    $$ v = \frac{Z R T}{P} $$

    where $v$ is the molar volume, $Z(T, P, y)$ is the compressibility
    factor, $T$ is the temperature, $P$ is the pressure, and $y$ is the
    mole fraction vector.

    Parameters
    ----------
    T : float
        Temperature. Unit = K.
    P : float
        Pressure. Unit = Pa.
    y : FloatVector
        Mole fractions of all components. Unit = mol/mol.

    Returns
    -------
    float
        Molar volume of the fluid. Unit = m³/mol.
    """
    return self.Z(T, P, y)*R*T/P

Examples¤

Estimate the molar volume of a 50 mol% ethylene/nitrogen at 350 K and 10 bar.

from polykin.thermo.eos import Virial
import numpy as np

Tc = [282.4, 126.2]    # K
Pc = [50.4e5, 33.9e5]  # Pa
Zc = [0.280, 0.290]
w = [0.089, 0.039]

eos = Virial(Tc, Pc, Zc, w)
v = eos.v(T=350., P=10e5, y=np.array([0.5, 0.5]))

print(f"{v:.2e} m³/mol")
2.87e-03 m³/mol