chemistrykit.thermo

chemistrykit.thermo hub · History · Examples · API reference

chemistrykit.thermo: chemical thermodynamics.

Equations of state (ideal gas, van der Waals, Redlich-Kwong, Peng-Robinson) and Lewis fugacities; Hess’s-law thermochemistry; Clausius-Clapeyron phase boundaries and the Gibbs phase rule; reaction equilibrium (Kp/Kc, the reaction quotient, the van’t Hoff equation, and a Gibbs-energy-minimization equilibrium-composition solver); and Raoult’s/Henry’s law mixtures with colligative properties (freezing-point depression, boiling-point elevation, osmotic pressure) and the Margules activity-coefficient model of non-ideal solutions.

class chemistrykit.thermo.BinaryIdealSolution(P_A_star, P_B_star)

Bases: object

A binary liquid mixture obeying Raoult’s law for both components.

Combining Raoult’s law for each component with Dalton’s law for the vapor phase gives the total vapor pressure as a function of liquid composition, and the vapor-phase composition in equilibrium with it (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 5.4a):

\[P_{tot}(x_A) = x_A P_A^* + (1-x_A) P_B^*, \qquad y_A = \frac{x_A P_A^*}{P_{tot}(x_A)}\]

Parameters:

• P_A_star (float) – Vapor pressure of pure A.

• P_B_star (float) – Vapor pressure of pure B.

P_A_star: float

P_B_star: float

total_pressure(x_A)

Total vapor pressure above a liquid of composition x_A.

Parameters:

x_A (float or array-like of float) – Liquid-phase mole fraction of A.

Returns:

float or ndarray

vapor_composition(x_A)

Vapor-phase mole fraction of A, \(y_A\), in equilibrium with liquid composition x_A.

Parameters:

x_A (float or array-like of float) – Liquid-phase mole fraction of A.

Returns:

float or ndarray

Examples

For a 1:1 mixture of equally volatile components, the vapor has the same composition as the liquid:

>>> solution = BinaryIdealSolution(P_A_star=100.0, P_B_star=100.0)
>>> round(float(solution.vapor_composition(0.5)), 6)
0.5

The more volatile component (higher pure vapor pressure) is enriched in the vapor relative to the liquid:

>>> solution = BinaryIdealSolution(P_A_star=200.0, P_B_star=50.0)
>>> bool(solution.vapor_composition(0.5) > 0.5)
True

class chemistrykit.thermo.ClausiusClapeyron(delta_h_vap, T_ref, P_ref, R_gas=8.31446261815324)

Bases: object

The (integrated) Clausius-Clapeyron equation for a liquid-vapor (or solid-vapor) boundary.

Assuming the enthalpy of vaporization \(\Delta H_{vap}\) is constant over the temperature range of interest, that the vapor behaves as an ideal gas, and that the molar volume of the condensed phase is negligible compared to the vapor’s (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 4.1, eq. 4.11), the Clausius-Clapeyron differential equation \(dP/dT = \Delta H_{vap}P/(RT^2)\) integrates to

\[\ln\frac{P(T)}{P_{ref}} = -\frac{\Delta H_{vap}}{R} \left(\frac{1}{T} - \frac{1}{T_{ref}}\right)\]

Parameters:

• delta_h_vap (float) – Enthalpy of vaporization, in J/mol (assumed temperature-independent – an approximation whose quality degrades far from T_ref).

• T_ref (float) – A reference temperature at which the vapor pressure P_ref is known, in K.

• P_ref (float) – Vapor pressure at T_ref, in Pa.

• R_gas (float)

Examples

Water: normal boiling point 373.15 K at 101325 Pa, using the (literature) enthalpy of vaporization at the normal boiling point, 40700 J/mol (Atkins & de Paula, Table 4.1), to estimate the vapor pressure at 363.15 K (90 degC). The constant-\(\Delta H_{vap}\) approximation gives 70.6 kPa, close to the accepted steam-table value of about 70.1 kPa – the ~1% discrepancy is exactly the approximation error flagged above (\(\Delta H_{vap}\) actually decreases somewhat as T rises toward the critical point):

>>> water = ClausiusClapeyron(delta_h_vap=40700.0, T_ref=373.15, P_ref=101325.0)
>>> round(float(water.pressure(363.15)) / 1000.0, 1)
70.6

P_ref: float

R_gas: float = 8.31446261815324

T_ref: float

boiling_point(P)

Return the temperature at which the vapor pressure equals P (invert pressure()).

Parameters:

P (float or array-like of float) – Target vapor pressure, in Pa.

Returns:

float or ndarray – Temperature, in K.

Examples

Inverting pressure() recovers the reference point exactly:

>>> water = ClausiusClapeyron(delta_h_vap=40700.0, T_ref=373.15, P_ref=101325.0)
>>> round(float(water.boiling_point(101325.0)), 2)
373.15

delta_h_vap: float

classmethod from_two_points(T1, P1, T2, P2)

Determine \(\Delta H_{vap}\) from two (T, P) points on the phase boundary.

Solving the integrated equation for the slope between two known points:

\[\Delta H_{vap} = \frac{-R \ln(P_2/P_1)}{1/T_2 - 1/T_1}\]

Parameters:

• T1 (float) – First reference point (K, Pa).

• P1 (float) – First reference point (K, Pa).

• T2 (float) – Second point (K, Pa).

• P2 (float) – Second point (K, Pa).

Return type:

ClausiusClapeyron

Returns:

ClausiusClapeyron – With T_ref, P_ref set to (T1, P1).

Examples

Recovering a known enthalpy of vaporization from two synthetic points generated with it:

>>> model = ClausiusClapeyron(delta_h_vap=35000.0, T_ref=300.0, P_ref=1.0e4)
>>> T2 = 320.0
>>> P2 = model.pressure(T2)
>>> fit = ClausiusClapeyron.from_two_points(300.0, 1.0e4, T2, P2)
>>> round(fit.delta_h_vap, 2)
35000.0

pressure(T)

Return the vapor pressure at temperature(s) T.

Parameters:

T (float or array-like of float) – Absolute temperature(s), in K.

Returns:

float or ndarray – Vapor pressure, in Pa (same units as P_ref).

class chemistrykit.thermo.EquationOfState

Bases: ABC

Common interface for a pure-substance PVT equation of state.

Concrete subclasses implement pressure() (an explicit formula) and molar_volume() (generally a root-find, for any EOS cubic or higher in \(V_m\)); compressibility_factor() then follows for every subclass for free.

R: float = 8.31446261815324

compressibility_factor(P, T, **kwargs)

Return the compressibility factor \(Z = PV_m/(RT)\).

\(Z = 1\) exactly for an ideal gas; deviations from 1 quantify real-gas non-ideality (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 1.3).

Parameters:

• P (float) – Pressure, in Pa.

• T (float) – Absolute temperature, in K.

• **kwargs – Forwarded to molar_volume() (e.g. branch for a cubic EOS).

Returns:

float

abstractmethod molar_volume(P, T, **kwargs)

Return the molar volume at pressure P and temperature T.

Parameters:

• P (float) – Pressure, in Pa.

• T (float) – Absolute temperature, in K.

Returns:

float – Molar volume, in m^3/mol.

abstractmethod pressure(Vm, T)

Return the pressure at molar volume Vm and temperature T.

Parameters:

• Vm (float or array-like of float) – Molar volume, in m^3/mol.

• T (float or array-like of float) – Absolute temperature, in K.

Returns:

float or ndarray – Pressure, in Pa.

class chemistrykit.thermo.EquilibriumComposition(species, n, x, extents, converged)

Bases: object

Result of a solve_equilibrium_composition() call.

Parameters:

• species (tuple)

• n (ndarray)

• x (ndarray)

• extents (ndarray)

• converged (bool)

converged: bool

Whether the underlying optimizer reported success.

Type:

bool

extents: ndarray

Equilibrium extent(s) of reaction \(\xi_j\).

Type:

ndarray, shape (n_reactions,)

moles(name)

Return the equilibrium mole amount of a single named species.

Parameters:

name (str)

Return type:

float

Returns:

float

n: ndarray

Equilibrium mole amounts.

Type:

ndarray, shape (n_species,)

species: tuple

Species names, in the order used throughout.

Type:

tuple of str

x: ndarray

Equilibrium mole fractions.

Type:

ndarray, shape (n_species,)

class chemistrykit.thermo.IdealGas

Bases: EquationOfState

The ideal (perfect) gas law: \(PV_m = RT\).

See Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 1.1.

Examples

The molar volume of an ideal gas at STP (1 atm-scale pressure of 101325 Pa, 0 degC) is the textbook 22.4 L/mol:

>>> gas = IdealGas()
>>> round(float(gas.molar_volume(P=101325.0, T=273.15)) * 1000.0, 1)
22.4
>>> round(float(gas.compressibility_factor(P=101325.0, T=273.15)), 6)
1.0

molar_volume(P, T, **kwargs)

Return the molar volume at pressure P and temperature T.

Parameters:

• P (float) – Pressure, in Pa.

• T (float) – Absolute temperature, in K.

Returns:

float – Molar volume, in m^3/mol.

pressure(Vm, T)

Return the pressure at molar volume Vm and temperature T.

Parameters:

• Vm (float or array-like of float) – Molar volume, in m^3/mol.

• T (float or array-like of float) – Absolute temperature, in K.

Returns:

float or ndarray – Pressure, in Pa.

class chemistrykit.thermo.MargulesSolution(A, P1_star, P2_star)

Bases: object

A non-ideal binary liquid described by the one-parameter (two-suffix) Margules model.

Margules expanded the logarithms of the activity coefficients as power series in mole fraction (M. Margules, Sitzungsber. Kais. Akad. Wiss. Wien, Math.-Naturwiss. Kl. 104, 1243-1278 (1895)). Truncated at its first term, the model has a symmetric excess Gibbs energy \(G^E/(RT) = A\,x_1x_2\), giving

\[\ln\gamma_1 = A\,x_2^2, \qquad \ln\gamma_2 = A\,x_1^2\]

and the modified Raoult’s law \(P_i = x_i\gamma_i P_i^*\) (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 5.3). \(A > 0\) gives positive deviations from Raoult’s law, \(A < 0\) negative ones, and \(A = 0\) recovers BinaryIdealSolution.

Parameters:

• A (float) – Dimensionless Margules parameter, \(W/(RT)\).

• P1_star (float) – Vapor pressure of pure component 1.

• P2_star (float) – Vapor pressure of pure component 2.

Examples

At infinite dilution component 1 obeys Henry’s law with \(K_H = P_1^* e^{A}\), while near \(x_1 = 1\) its activity coefficient tends to 1 (Raoult’s law):

>>> import numpy as np
>>> sol = MargulesSolution(A=1.2, P1_star=30.0, P2_star=20.0)
>>> round(sol.henry_constant_1 / 30.0, 6) == round(float(np.exp(1.2)), 6)
True
>>> g1, g2 = sol.activity_coefficients(1.0)
>>> float(g1), round(float(g2), 6) == round(float(np.exp(1.2)), 6)
(1.0, True)

A: float

P1_star: float

P2_star: float

activity_coefficients(x1)

Return \((\gamma_1, \gamma_2)\) at liquid mole fraction x1.

Parameters:

x1 (float or array-like of float) – Liquid-phase mole fraction of component 1.

Returns:

tuple of (float or ndarray)

excess_gibbs(x1, T, R_gas=8.31446261815324)

Molar excess Gibbs energy \(G^E = RT\,A\,x_1x_2\), in J/mol.

Parameters:

• x1 (float or array-like of float) – Liquid-phase mole fraction of component 1.

• T (float) – Absolute temperature, in K.

• R_gas (float) – Gas constant.

Returns:

float or ndarray

property henry_constant_1: float

Henry’s-law constant of component 1 at infinite dilution, \(P_1^* e^{A}\).

partial_pressures(x1)

Return \((P_1, P_2)\) from the modified Raoult’s law \(P_i = x_i\gamma_iP_i^*\).

Parameters:

x1 (float or array-like of float) – Liquid-phase mole fraction of component 1.

Returns:

tuple of (float or ndarray)

total_pressure(x1)

Total vapor pressure above a liquid of composition x1.

Parameters:

x1 (float or array-like of float) – Liquid-phase mole fraction of component 1.

Returns:

float or ndarray

vapor_composition(x1)

Vapor-phase mole fraction of component 1 in equilibrium with liquid x1.

Parameters:

x1 (float or array-like of float) – Liquid-phase mole fraction of component 1.

Returns:

float or ndarray

class chemistrykit.thermo.PengRobinson(Tc, Pc, omega=0.0)

Bases: EquationOfState

The Peng-Robinson equation of state.

\[P = \frac{RT}{V_m - b} - \frac{a\,\alpha(T)}{V_m^2 + 2bV_m - b^2}\]

with \(a = 0.45724\,R^2T_c^2/P_c\), \(b = 0.07780\,RT_c/P_c\), and the temperature-dependent attraction factor

\[\alpha(T) = \left[1 + \kappa\left(1 - \sqrt{T/T_c}\right)\right]^2, \qquad \kappa = 0.37464 + 1.54226\,\omega - 0.26992\,\omega^2\]

where \(\omega\) is Pitzer’s acentric factor. The \(\kappa(\omega)\) correlation was fitted by Peng and Robinson so that the equation reproduces pure-substance vapor pressures, which is why, unlike van der Waals or Redlich-Kwong, it predicts saturated liquid densities and vapor pressures well enough for process design (D.-Y. Peng & D. B. Robinson, Ind. Eng. Chem. Fundam. 15, 59-64 (1976), eqs. 7-9 and 17-18).

Parameters:

• Tc (float) – Critical temperature, in K.

• Pc (float) – Critical pressure, in Pa.

• omega (float) – Pitzer acentric factor (0 for a simple fluid such as argon).

Examples

At the critical temperature \(\alpha = 1\), and the equation’s universal critical compressibility factor is \(Z_c \approx 0.3074\), closer to real fluids’ 0.23-0.31 than van der Waals’s 3/8 or Redlich-Kwong’s 1/3. The critical isotherm passes through \((V_c, P_c)\) with \(V_c = Z_c RT_c/P_c\):

>>> eos = PengRobinson(Tc=304.13, Pc=7.3773e6, omega=0.224)  # CO2
>>> float(eos.alpha(304.13))
1.0
>>> Vc = 0.30740 * eos.R * 304.13 / 7.3773e6
>>> round(float(eos.pressure(Vc, 304.13)) / 7.3773e6, 3)
1.0

alpha(T)

Return the attraction factor \(\alpha(T) = [1 + \kappa(1 - \sqrt{T/T_c})]^2\).

Parameters:

T (float or array-like of float) – Absolute temperature(s), in K.

Returns:

float or ndarray

molar_volume(P, T, branch='vapor')

Solve the Peng-Robinson cubic for Vm at pressure P, temperature T.

In terms of \(A = a\alpha P/(RT)^2\) and \(B = bP/(RT)\), the equation becomes a cubic in \(Z = PV_m/(RT)\) (Peng & Robinson 1976, eq. 5):

\[Z^3 - (1-B)Z^2 + (A - 3B^2 - 2B)Z - (AB - B^2 - B^3) = 0\]

Parameters:

• P (float) – Pressure, in Pa.

• T (float) – Absolute temperature, in K.

• branch (str) – Return the largest (“vapor”) or smallest (“liquid”) admissible root.

Return type:

float

Returns:

float

Examples

>>> eos = PengRobinson(Tc=304.13, Pc=7.3773e6, omega=0.224)
>>> Vm = eos.molar_volume(P=1.0e5, T=300.0)
>>> round(float(eos.pressure(Vm, T=300.0)), 3)
100000.0

pressure(Vm, T)

Return the pressure at molar volume Vm and temperature T.

Parameters:

• Vm (float or array-like of float) – Molar volume, in m^3/mol.

• T (float or array-like of float) – Absolute temperature, in K.

Returns:

float or ndarray – Pressure, in Pa.

class chemistrykit.thermo.RedlichKwong(a, b)

Bases: EquationOfState

The Redlich-Kwong equation of state.

\[P = \frac{RT}{V_m - b} - \frac{a}{\sqrt{T}\,V_m(V_m+b)}\]

An improvement on van der Waals that gives the attractive term an explicit temperature dependence (O. Redlich & J. N. S. Kwong, Chem. Rev. 44, 233 (1949)).

Parameters:

• a (float) – Attraction parameter, in Pa m^6 K^0.5 mol^-2.

• b (float) – Excluded-volume parameter, in m^3/mol.

classmethod from_critical_constants(Tc, Pc)

Build a RedlichKwong EOS from critical temperature and pressure.

\[a = 0.42748\,\frac{R^2 T_c^{2.5}}{P_c}, \qquad b = 0.08664\,\frac{R T_c}{P_c}\]

(Redlich & Kwong 1949; the numeric coefficients are tabulated e.g. in Smith, Van Ness & Abbott, Introduction to Chemical Engineering Thermodynamics, 7th ed., Table 3.1).

Parameters:

• Tc (float) – Critical temperature, in K.

• Pc (float) – Critical pressure, in Pa.

Return type:

RedlichKwong

Returns:

RedlichKwong

Examples

As with van der Waals’s 3/8, Redlich-Kwong’s critical-point construction fixes a universal critical molar volume \(V_c = RT_c/(3P_c)\) and hence a universal critical compressibility factor, here exactly 1/3, again independent of which substance’s (Tc, Pc) is used. (As in VanDerWaals.from_critical_constants(), this checks the closed-form critical relation directly rather than numerically solving the cubic exactly at its ill-conditioned triple root.)

>>> Tc, Pc = 304.13, 7.3773e6  # CO2
>>> eos = RedlichKwong.from_critical_constants(Tc=Tc, Pc=Pc)
>>> Vc = eos.R * Tc / (3.0 * Pc)
>>> round(Pc * Vc / (eos.R * Tc), 6)
0.333333

molar_volume(P, T, branch='vapor')

Solve the Redlich-Kwong cubic for Vm at pressure P, temperature T.

Multiplying \(P = RT/(V_m-b) - a/(\sqrt{T}V_m(V_m+b))\) through by \(\sqrt{T}V_m(V_m+b)(V_m-b)\) and collecting terms in \(V_m\) gives

\[V_m^3 - \frac{RT}{P} V_m^2 + \left(\frac{a}{P\sqrt{T}} - b^2 - \frac{bRT}{P}\right) V_m - \frac{ab}{P\sqrt{T}} = 0\]

Parameters:

• P (float) – Pressure, in Pa.

• T (float) – Absolute temperature, in K.

• branch (str) – Return the largest (“vapor”) or smallest (“liquid”) real positive root.

Return type:

float

Returns:

float

Examples

Solving for Vm and substituting back reproduces the original pressure:

>>> rk = RedlichKwong(a=6.4239, b=2.7143e-5)  # e.g. representative CO2-like values
>>> Vm = rk.molar_volume(P=1.0e5, T=300.0)
>>> round(float(rk.pressure(Vm, T=300.0)), 3)
100000.0

pressure(Vm, T)

Return the pressure at molar volume Vm and temperature T.

Parameters:

• Vm (float or array-like of float) – Molar volume, in m^3/mol.

• T (float or array-like of float) – Absolute temperature, in K.

Returns:

float or ndarray – Pressure, in Pa.

class chemistrykit.thermo.VanDerWaals(a, b)

Bases: EquationOfState

The van der Waals equation of state.

\[P = \frac{RT}{V_m - b} - \frac{a}{V_m^2}\]

where a corrects for intermolecular attraction and b for the finite volume of the molecules themselves (J. D. van der Waals, doctoral thesis, Leiden, 1873; Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 1.3b).

Parameters:

• a (float) – Attraction parameter, in Pa m^6 mol^-2.

• b (float) – Excluded-volume parameter, in m^3/mol.

Examples

With a = b = 0 van der Waals reduces exactly to the ideal gas law:

>>> vdw = VanDerWaals(a=0.0, b=0.0)
>>> ideal = IdealGas()
>>> round(float(vdw.pressure(Vm=0.01, T=300.0)), 6) == round(float(ideal.pressure(Vm=0.01, T=300.0)), 6)
True

classmethod from_critical_constants(Tc, Pc)

Build a VanDerWaals EOS from critical temperature and pressure.

At the critical point the van der Waals cubic has a triple root, which fixes

\[a = \frac{27 R^2 T_c^2}{64 P_c}, \qquad b = \frac{R T_c}{8 P_c}\]

(Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 1.3c).

Parameters:

• Tc (float) – Critical temperature, in K.

• Pc (float) – Critical pressure, in Pa.

Return type:

VanDerWaals

Returns:

VanDerWaals

Examples

A structural consequence of these formulas – independent of which substance’s (Tc, Pc) is used – is that the van der Waals critical molar volume is always \(V_c = 3b\), and hence the critical compressibility factor is always exactly 3/8 (Atkins & de Paula, Table 1.6), unlike the true value for real gases (typically 0.23-0.31). (Solving the cubic numerically exactly at the critical point, where it has a triple root, is a textbook-classic ill-conditioned root-finding problem, so this checks the closed-form critical relation directly rather than round-tripping through molar_volume().)

>>> Tc, Pc = 304.13, 7.3773e6  # CO2
>>> eos = VanDerWaals.from_critical_constants(Tc=Tc, Pc=Pc)
>>> Vc = 3.0 * eos.b
>>> round(Pc * Vc / (eos.R * Tc), 6)
0.375

molar_volume(P, T, branch='vapor')

Solve the van der Waals cubic for Vm at pressure P, temperature T.

Rearranging \(P = RT/(V_m-b) - a/V_m^2\) into a cubic in \(V_m\):

\[V_m^3 - \left(b + \frac{RT}{P}\right) V_m^2 + \frac{a}{P} V_m - \frac{ab}{P} = 0\]

Below the critical temperature this can have three real positive roots (a metastable liquid branch, an unstable middle root, and a vapor branch).

Parameters:

• P (float) – Pressure, in Pa.

• T (float) – Absolute temperature, in K.

• branch (str) – Return the largest (“vapor”) or smallest (“liquid”) real positive root.

Return type:

float

Returns:

float

Examples

Solving for Vm and substituting back reproduces the original pressure:

>>> vdw = VanDerWaals(a=0.1448, b=3.913e-5)  # e.g. representative CO2-like values
>>> Vm = vdw.molar_volume(P=1.0e5, T=300.0)
>>> round(float(vdw.pressure(Vm, T=300.0)), 3)
100000.0

pressure(Vm, T)

Return the pressure at molar volume Vm and temperature T.

Parameters:

• Vm (float or array-like of float) – Molar volume, in m^3/mol.

• T (float or array-like of float) – Absolute temperature, in K.

Returns:

float or ndarray – Pressure, in Pa.

class chemistrykit.thermo.VantHoffFit(delta_h, delta_s, r_squared, R_gas=8.31446261815324)

Bases: object

Result of fitting \(\ln K\) vs. \(1/T\) data to the van’t Hoff equation.

Since \(\ln K = -\Delta H^\circ/(RT) + \Delta S^\circ/R\) (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 6.4, eq. 6.16), the slope of a plot of \(\ln K\) against \(1/T\) gives \(-\Delta H^\circ/R\) and the intercept gives \(\Delta S^\circ/R\).

Parameters:

• delta_h (float)

• delta_s (float)

• r_squared (float)

• R_gas (float)

R_gas: float = 8.31446261815324

Gas constant used in the fit, in J mol^-1 K^-1.

Type:

float

delta_h: float

Fitted standard reaction enthalpy, in J/mol.

Type:

float

delta_s: float

Fitted standard reaction entropy, in J mol^-1 K^-1.

Type:

float

predict(T)

Evaluate the fitted van’t Hoff equation at temperature(s) T.

Parameters:

T (float or array-like of float) – Absolute temperature(s), in K.

Returns:

float or ndarray

r_squared: float

Coefficient of determination of the linear (ln K vs 1/T) fit.

Type:

float

chemistrykit.thermo.boiling_point_elevation(Kb, b, i=1.0)

Boiling-point elevation \(\Delta T_b = i K_b b\).

Parameters:

• Kb (float) – Ebullioscopic constant of the solvent, in K kg/mol (see CRYOSCOPIC_CONSTANTS for common solvents).

• b (float) – Molality of solute, in mol/kg.

• i (float) – Van’t Hoff factor (see freezing_point_depression()).

Return type:

float

Returns:

float – Boiling-point elevation, in K.

Examples

>>> round(boiling_point_elevation(Kb=0.51, b=1.00, i=2.0), 2)
1.02

chemistrykit.thermo.fit_van_t_hoff(T, K, R_gas=8.31446261815324)

Fit equilibrium-constant vs. temperature data to the van’t Hoff equation.

Linearizes \(\ln K = -\Delta H^\circ/R \cdot (1/T) + \Delta S^\circ/R\) and fits by ordinary least squares – the “van’t Hoff plot” (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 6.4), structurally identical to chemistrykit.kinetics.systems.arrhenius.fit_arrhenius()’s Arrhenius plot.

Parameters:

• T (array-like of float) – Absolute temperatures, in K (at least 2 distinct values).

• K (array-like of float) – Equilibrium constants measured at each temperature in T.

• R_gas (float) – Gas constant, in J mol^-1 K^-1.

Return type:

VantHoffFit

Returns:

VantHoffFit

Examples

Generate exact data from a known (ΔH°, ΔS°) and recover ΔH°:

>>> import numpy as np
>>> T = np.array([280.0, 300.0, 320.0, 340.0, 360.0])
>>> K_ref, delta_h = 1.0, 45_000.0
>>> K = van_t_hoff_equilibrium_constant(T, T_ref=300.0, K_ref=K_ref, delta_h=delta_h)
>>> fit = fit_van_t_hoff(T, K)
>>> round(fit.delta_h, 2)
45000.0
>>> round(fit.r_squared, 6)
1.0

chemistrykit.thermo.freezing_point_depression(Kf, b, i=1.0)

Freezing-point depression \(\Delta T_f = i K_f b\).

A colligative property: it depends on the number of dissolved solute particles per kg of solvent, not their identity (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 5.5). See also chemistrykit.solutions, which uses the same van’t Hoff factor i convention for strong-electrolyte dissociation, but does not duplicate this colligative-property machinery – it lives here in chemistrykit.thermo only.

Parameters:

• Kf (float) – Cryoscopic constant of the solvent, in K kg/mol (see CRYOSCOPIC_CONSTANTS for common solvents).

• b (float) – Molality of solute, in mol/kg.

• i (float) – Van’t Hoff factor: the number of particles each formula unit of solute dissociates into (1 for a nonelectrolyte, 2 for e.g. NaCl assuming complete dissociation).

Return type:

float

Returns:

float – Freezing-point depression, in K (the new freezing point is \(T_f^* - \Delta T_f\)).

Examples

A 1.00 molal aqueous NaCl solution (i=2, ideal complete dissociation):

>>> round(freezing_point_depression(Kf=1.86, b=1.00, i=2.0), 2)
3.72

chemistrykit.thermo.fugacity(eos, P, T, branch='vapor')

Fugacity \(f = \phi P\) of a pure fluid, in Pa.

Parameters:

• eos (EquationOfState) – Pure-fluid equation of state.

• P (float) – Pressure, in Pa.

• T (float) – Absolute temperature, in K.

• branch (str) – Root of a cubic equation of state to use.

Return type:

float

Returns:

float

Examples

An ideal gas’s fugacity is its pressure:

>>> from chemistrykit.thermo import IdealGas
>>> round(fugacity(IdealGas(), 5.0e5, 300.0), 6)
500000.0

chemistrykit.thermo.fugacity_coefficient(eos, P, T, branch='vapor')

Fugacity coefficient \(\phi = f/P\) of a pure fluid described by eos.

The residual integral is evaluated in density, \(\rho = 1/V_m\), as \(\int_0^{\rho}(P - RT\rho')/\rho'^2\,d\rho'\) with 96-point Gauss-Legendre quadrature, so any EquationOfState works.

Parameters:

• eos (EquationOfState) – Pure-fluid equation of state.

• P (float) – Pressure, in Pa.

• T (float) – Absolute temperature, in K.

• branch (str) – Root of a cubic equation of state to use (ignored by IdealGas).

Return type:

float

Returns:

float

Examples

For van der Waals the integral has the closed form \(\ln\phi = b/(V_m-b) - 2a/(RTV_m) - \ln[P(V_m-b)/(RT)]\):

>>> import numpy as np
>>> from chemistrykit.thermo import VanDerWaals
>>> vdw = VanDerWaals(a=0.3640, b=4.267e-5)
>>> P, T = 2.0e6, 300.0
>>> V = vdw.molar_volume(P, T)
>>> RT = vdw.R * T
>>> exact = np.exp(vdw.b / (V - vdw.b) - 2 * vdw.a / (RT * V) - np.log(P * (V - vdw.b) / RT))
>>> bool(np.isclose(fugacity_coefficient(vdw, P, T), exact, rtol=1e-10))
True

chemistrykit.thermo.gibbs_energy_of_mixture(n, gibbs_formation, T, P=100000.0, P_standard=100000.0, R_gas=8.31446261815324)

Total Gibbs energy of an ideal-gas mixture, \(G = \sum_i n_i \mu_i\).

Using the ideal-gas chemical potential \(\mu_i = \Delta G_{f,i}^\circ + RT\ln(x_i P/P^\circ)\) (Smith, Van Ness & Abbott, Introduction to Chemical Engineering Thermodynamics, 7th ed., Ch. 13.2), where x_i is species i’s mole fraction. This is the objective function minimized by solve_equilibrium_composition().

Parameters:

• n (array-like of float) – Mole amount of each species (must sum to a positive total).

• gibbs_formation (array-like of float) – Standard Gibbs energy of formation of each species, in J/mol, same order as n.

• T (float) – Absolute temperature, in K.

• P (float) – Total pressure, in Pa.

• P_standard (float) – Standard-state pressure, in Pa.

• R_gas (float) – Gas constant, in J mol^-1 K^-1.

Return type:

float

Returns:

float – Total Gibbs energy, in J.

chemistrykit.thermo.gibbs_phase_rule(n_components, n_phases, reactions=0)

The Gibbs phase rule: \(F = C - P + 2 - r\).

F is the number of intensive degrees of freedom (e.g. temperature, pressure, composition variables) that can be varied independently while the system remains in the same set of coexisting phases at equilibrium; C is the number of independent components, P the number of phases in equilibrium, and r the number of independent reaction equilibria constraining the composition (0 for a non-reactive system) (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 4.5; J. W. Gibbs, Trans. Connecticut Acad. 3, 108 (1876)).

Parameters:

• n_components (int) – Number of independent chemical components, C.

• n_phases (int) – Number of phases in equilibrium, P.

• reactions (int) – Number of independent reaction-equilibrium constraints, r.

Return type:

int

Returns:

int – Degrees of freedom, F.

Examples

A pure substance (C=1) at its triple point (three phases coexisting) has zero degrees of freedom – the triple point is a single fixed (T, P):

>>> gibbs_phase_rule(n_components=1, n_phases=3)
0

A pure substance with two phases in equilibrium (e.g. liquid-vapor) has one degree of freedom – the phase boundary is a curve, P(T):

>>> gibbs_phase_rule(n_components=1, n_phases=2)
1

chemistrykit.thermo.henry_law_pressure(x, K_H)

Henry’s law: \(P_B = x_B K_H\).

For a dilute solute B, the partial vapor pressure is proportional to its mole fraction with the (empirical, solute- and solvent-specific) Henry’s law constant \(K_H\) as the proportionality constant – unlike Raoult’s law, \(K_H \neq P_B^*\) in general, because the solute’s local environment in dilute solution (surrounded by solvent) differs from pure solute (W. Henry, Philos. Trans. R. Soc. 93, 29 (1803); Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 5.4b).

Parameters:

• x (float or array-like of float) – Mole fraction of the (dilute) solute in solution.

• K_H (float) – Henry’s law constant for this solute/solvent pair, same pressure units as the desired result.

Returns:

float or ndarray

Examples

>>> round(float(henry_law_pressure(x=0.001, K_H=1.5e5)), 3)
150.0

chemistrykit.thermo.hess_law_enthalpy(multipliers, step_enthalpies)

Enthalpy of a reaction built as a linear combination of steps, \(\Delta H = \sum_k c_k\,\Delta H_k\).

Parameters:

• multipliers (array-like of float) – Coefficient \(c_k\) applied to each step reaction (negative to reverse a step, 2 to double it, and so on).

• step_enthalpies (array-like of float) – Enthalpy change of each step as written, in J/mol (or any consistent unit).

Return type:

float

Returns:

float

Examples

Burning graphite to CO (-110.5 kJ/mol) and then CO to CO2 (-283.0 kJ/mol) releases the same heat as burning graphite straight to CO2:

>>> hess_law_enthalpy([1.0, 1.0], [-110.5, -283.0])
-393.5

chemistrykit.thermo.kc_from_kp(Kp, delta_n, T, R_gas=8.31446261815324)

Convert \(K_p\) to \(K_c\) via \(K_c = K_p/(RT)^{\Delta n}\) (inverse of kp_from_kc()).

Parameters:

• Kp (float) – Equilibrium constant in terms of pressure.

• delta_n (float) – Change in moles of gas, products minus reactants.

• T (float) – Absolute temperature, in K.

• R_gas (float) – Gas constant.

Return type:

float

Returns:

float

Examples

Round-trips with kp_from_kc():

>>> Kp = kp_from_kc(Kc=1.8, delta_n=1.0, T=350.0)
>>> round(kc_from_kp(Kp, delta_n=1.0, T=350.0), 6)
1.8

chemistrykit.thermo.kp_from_kc(Kc, delta_n, T, R_gas=8.31446261815324)

Convert \(K_c\) to \(K_p\) via \(K_p = K_c(RT)^{\Delta n}\).

Valid for reactions among ideal gases, where \(\Delta n\) is the change in moles of gas (products minus reactants) (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 6.2). Kc must be expressed in concentration units consistent with R_gas (e.g. mol/m^3 with the SI R_gas, or mol/L with R_gas = 0.0831446 L bar / (mol K)).

Parameters:

• Kc (float) – Equilibrium constant in terms of concentration.

• delta_n (float) – Change in moles of gas, products minus reactants.

• T (float) – Absolute temperature, in K.

• R_gas (float) – Gas constant, consistent with the concentration units of Kc.

Return type:

float

Returns:

float

Examples

With delta_n = 0 (no change in moles of gas), Kp = Kc:

>>> round(kp_from_kc(Kc=2.5, delta_n=0.0, T=298.15), 6)
2.5

chemistrykit.thermo.osmotic_pressure(M, T, i=1.0, R_gas=8.31446261815324)

The van’t Hoff equation for osmotic pressure, \(\Pi = i M R T\).

Formally identical to the ideal gas law – van’t Hoff originally noted the (coincidental, but pedagogically useful) analogy (J. H. van’t Hoff, 1887; Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 5.5e).

Parameters:

• M (float) – Molarity of solute, in mol/m^3 (use M_mol_per_L * 1000 to convert from the more commonly tabulated mol/L).

• T (float) – Absolute temperature, in K.

• i (float) – Van’t Hoff factor (see freezing_point_depression()).

• R_gas (float) – Gas constant, in J mol^-1 K^-1; with M in mol/m^3 this gives Pi in Pa.

Return type:

float

Returns:

float – Osmotic pressure, in Pa.

Examples

A 0.100 mol/L (=100 mol/m^3) nonelectrolyte solution at 298.15 K:

>>> round(osmotic_pressure(M=100.0, T=298.15) / 1000.0, 2)
247.9

chemistrykit.thermo.raoult_vapor_pressure(x, P_pure)

Raoult’s law: \(P_A = x_A P_A^*\).

The partial vapor pressure of a component in an ideal mixture is proportional to its mole fraction, with the pure-component vapor pressure as the constant of proportionality (F. M. Raoult, 1887; Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 5.4a).

Parameters:

• x (float or array-like of float) – Mole fraction of the component in the liquid.

• P_pure (float) – Vapor pressure of the pure component, \(P_A^*\).

Returns:

float or ndarray

Examples

>>> round(float(raoult_vapor_pressure(x=0.4, P_pure=100.0)), 6)
40.0

chemistrykit.thermo.reaction_enthalpy_from_formation(stoich_coeffs, formation_enthalpies)

Standard reaction enthalpy from enthalpies of formation, \(\Delta_rH^\circ = \sum_i \nu_i\,\Delta_fH_i^\circ\).

This is Hess’s law applied to the cycle “decompose reactants into their elements, then build the products from them”.

Parameters:

• stoich_coeffs (array-like of float) – Signed stoichiometric coefficients (negative for reactants).

• formation_enthalpies (array-like of float) – Standard enthalpy of formation of each species, same order, in J/mol or kJ/mol (0 for elements in their reference state).

Return type:

float

Returns:

float

Examples

Combustion of methane, \(CH_4 + 2O_2 \to CO_2 + 2H_2O(l)\), with formation enthalpies -74.8, 0, -393.5 and -285.8 kJ/mol:

>>> round(reaction_enthalpy_from_formation([-1, -2, 1, 2], [-74.8, 0.0, -393.5, -285.8]), 1)
-890.3

chemistrykit.thermo.reaction_quotient(activities, stoich_coeffs)

The reaction quotient \(Q = \prod_i a_i^{\nu_i}\).

stoich_coeffs are signed net stoichiometric coefficients (positive for products, negative for reactants), so a reactant with coefficient \(-\nu\) contributes \(a^{-\nu} = 1/a^{\nu}\) to the product, matching the usual “products over reactants” definition (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 6.1). activities may be partial pressures relative to the standard pressure (dimensionless, for \(Q_p\)) or concentrations relative to a standard concentration (for \(Q_c\)); which one determines whether the result should be compared against \(K_p\) or \(K_c\).

Parameters:

• activities (array-like of float) – Activity (or activity-like relative pressure/concentration) of each species.

• stoich_coeffs (array-like of float) – Signed net stoichiometric coefficient of each species, same order as activities.

Returns:

float

Examples

For \(N_2O_4 \rightleftharpoons 2NO_2\), at the equilibrium mole fractions corresponding to \(Q = K = 4\) (i.e. extent \(\xi = \sqrt{K/(4+K)}\), as in solve_equilibrium_composition()’s worked example below):

>>> import numpy as np
>>> xi = np.sqrt(4.0 / 8.0)
>>> x_N2O4, x_NO2 = (1.0 - xi) / (1.0 + xi), 2.0 * xi / (1.0 + xi)
>>> round(float(reaction_quotient([x_N2O4, x_NO2], [-1.0, 2.0])), 3)
4.0

chemistrykit.thermo.saturation_pressure(eos, T)

Vapor pressure of a cubic equation of state from Lewis’s equal-fugacity condition.

Finds the pressure at which the liquid and vapor roots of eos have equal fugacity, \(f^L(P_{sat}) = f^V(P_{sat})\). The search is bracketed between the isotherm’s two spinodal pressures (its local minimum and maximum in \(P(V_m)\)), where three roots exist.

Parameters:

• eos (EquationOfState) – A cubic equation of state with an excluded-volume attribute b (VanDerWaals, RedlichKwong, PengRobinson).

• T (float) – Absolute temperature, below the equation’s critical temperature, in K.

Return type:

float

Returns:

float – Saturation pressure, in Pa.

Examples

The van der Waals coexistence curve at reduced temperature 0.9 lies at reduced pressure 0.647, whatever the substance:

>>> from chemistrykit.thermo import VanDerWaals
>>> Tc, Pc = 304.13, 7.3773e6
>>> vdw = VanDerWaals.from_critical_constants(Tc, Pc)
>>> round(saturation_pressure(vdw, 0.9 * Tc) / Pc, 3)
0.647

chemistrykit.thermo.solve_equilibrium_composition(species, stoich_matrix, n0, gibbs_formation, T, P=100000.0, R_gas=8.31446261815324)

Solve for the equilibrium composition of a reacting ideal-gas mixture.

Rather than solving \(Q(\xi) = K\) algebraically (which gets unwieldy for several simultaneous reactions), this parameterizes the feasible mole amounts by the extent(s) of reaction \(\xi_j\), \(n_i(\xi) = n_{i,0} + \sum_j \nu_{ij}\xi_j\) – which automatically satisfies atomic mass balance for any \(\xi\) – and numerically minimizes the total Gibbs energy gibbs_energy_of_mixture() over the feasible region \(n_i(\xi) \geq 0\), via SLSQP (Smith, Van Ness & Abbott, Introduction to Chemical Engineering Thermodynamics, 7th ed., Ch. 13.7). At the unconstrained minimum, \(dG/d\xi_j = \sum_i \nu_{ij}\mu_i = \Delta G_{rxn,j} = 0\), which is exactly the classical \(Q_j = K_j\) equilibrium condition for every reaction j – Gibbs-energy minimization and root-finding on the equilibrium constant are the same equilibrium condition, just reached by different numerical routes.

Parameters:

• species (sequence of str) – Ordered species names.

• stoich_matrix (array-like, shape (n_species, n_reactions)) – Net stoichiometric coefficient of each species in each reaction.

• n0 (array-like, shape (n_species,)) – Initial (pre-reaction) mole amounts.

• gibbs_formation (array-like, shape (n_species,)) – Standard Gibbs energy of formation of each species, in J/mol.

• T (float) – Absolute temperature, in K.

• P (float) – Total pressure, in Pa.

• R_gas (float) – Gas constant, in J mol^-1 K^-1.

Return type:

EquilibriumComposition

Returns:

EquilibriumComposition

Examples

\(N_2O_4 \rightleftharpoons 2NO_2\), choosing standard Gibbs energies of formation so that \(K = \exp(-\Delta G^\circ_{rxn}/RT) = 4\) exactly at 298.15 K (with \(\Delta G_f^\circ(N_2O_4) = 0\) as the reference), starting from 1 mol of pure \(N_2O_4\). Since \(Q(\xi) = x_{NO_2}^2/x_{N_2O_4} = 4\xi^2/(1-\xi^2)\) at \(P = P^\circ\), setting \(Q = K\) gives the closed-form extent \(\xi = \sqrt{K/(4+K)} = \sqrt{4/8} \approx 0.7071\), which the Gibbs-minimization solver should reproduce:

>>> import numpy as np
>>> from chemistrykit.constants import R
>>> T = 298.15
>>> K_target = 4.0
>>> delta_g_rxn = -R * T * np.log(K_target)
>>> gibbs_formation = [0.0, delta_g_rxn / 2.0]  # [N2O4, NO2]
>>> result = solve_equilibrium_composition(
...     species=("N2O4", "NO2"),
...     stoich_matrix=[[-1.0], [2.0]],
...     n0=[1.0, 0.0],
...     gibbs_formation=gibbs_formation,
...     T=T,
... )
>>> round(float(result.extents[0]), 4)
0.7071
>>> result.converged
True

chemistrykit.thermo.van_t_hoff_equilibrium_constant(T, T_ref, K_ref, delta_h, R_gas=8.31446261815324)

The (integrated) van’t Hoff equation \(\ln(K/K_{ref}) = -\Delta H^\circ/R\,(1/T - 1/T_{ref})\).

Assumes the standard reaction enthalpy \(\Delta H^\circ\) is constant over the temperature range of interest (Atkins & de Paula, Physical Chemistry, 11th ed., Ch. 6.4) – the same approximation, and the same functional form, as chemistrykit.thermo.systems.phase_equilibria.ClausiusClapeyron (a coincidence of both following from \(d\ln(\cdot)/dT \propto 1/T^2\), one for a phase boundary, the other for an equilibrium constant).

Parameters:

• T (float or array-like of float) – Absolute temperature(s) at which to evaluate K, in K.

• T_ref (float) – Reference temperature at which K_ref is known, in K.

• K_ref (float) – Equilibrium constant at T_ref.

• delta_h (float) – Standard reaction enthalpy, in J/mol (positive for endothermic).

• R_gas (float) – Gas constant, in J mol^-1 K^-1.

Returns:

float or ndarray

Examples

For an endothermic reaction, K increases with temperature (Le Chatelier’s principle):

>>> import numpy as np
>>> T = np.array([280.0, 300.0, 320.0])
>>> K = van_t_hoff_equilibrium_constant(T, T_ref=298.15, K_ref=1.0, delta_h=50e3)
>>> bool(np.all(np.diff(K) > 0))
True