ceos.rst

Cubic Equations of State

This property package implements a general form of a cubic equation of state which can be used for most cubic-type equations of state. This package supports phase equilibrium calculations with a smooth phase transition formulation that makes it amenable for equation oriented optimization. The following equations of state are currently supported:

• Peng-Robinson
• Soave-Redlich-Kwong

Flow basis: Molar

Units: SI units

State Variables:

The state block uses the following state variables:

Inputs

When instantiating the parameter block that uses this particular state block, 1 optional argument can be passed:

The valid_phase argument denotes the valid phases for a given set of inlet conditions. For example, if the user knows a priori that the it will only be a single phase (for example liquid only), then it is best not to include the complex flash equilibrium constraints in the model. If the user does not specify any option, then the package defaults to a 2 phase assumption meaning that the constraints to compute the phase equilibrium will be computed.

Degrees of Freedom

In general, the general cubic equation of state has a number of degrees of freedom equal to 2 + the number of components in the system (total flow rate, temperature, pressure and N-1 mole fractions). In some cases (primarily inlets to units), this is increased by 1 due to the removal of a constraint on the sum of mole fractions.

General Cubic Equation of State

All equations come from "The Properties of Gases and Liquids, 4th Edition" by Reid, Prausnitz and Poling. The general cubic equation of state is represented by the following equations:

0 = Z³ − (1 + B − uB)Z² + (A − uB − (u − w)B²)Z − AB − wB² − wB³

$$A = \frac{a_mP}{R^2T^2}$$

$$B = \frac{b_mP}{RT}$$

where Z is the compressibility factor of the mixture, a_m and b_m are properties of the mixture and u and w are parameters which depend on the specific equation of state being used as show in the table below.

Equation	u	w	ΩA	ΩB	κj
Peng-Robinson	2	-1	0.45724	0.07780	(1 + (1 − Tr2)(0.37464 + 1.54226ωj − 0.26992ωj2))2
Soave-Redlich-Kwong	1	0	0.42748	0.08664	(1 + (1 − Tr2)(0.48 + 1.574ωj − 0.176ωj2))2

The properties a_m and b_m are calculated from component specific properties a_j and b_j as shown below:

$$a_j = \frac{\Omega_AR^2T_{c,j}^2}{P_{c, j}}\kappa_j$$

$$b_j = \frac{\Omega_BRT_{c,j}}{P_{c,j}}$$

a_m = ∑_i∑_jy_iy_j(a_ia_j)^1/2(1 − k_ij)

b_m = ∑_iy_ib_i

where P_c, j and T_c, j are the component critical pressures and temperatures, y_j is the mole fraction of component j, k_ij are a set of binary interaction parameters which are specific to the equation of state and Ω_A, Ω_B and κ_j are taken from the table above. ω_j is the Pitzer acentric factor of each component.

The cubic equation of state is solved for each phase via a call to an external function which automatically identifies the correct root of the cubic and returns the value of Z as a function of A and B along with the first and second partial derivatives.

VLE Model with Smooth Phase Transition

The flash equations consists of the following equations:

Fⁱⁿ = F^liq + F^vap

z_iⁱⁿFⁱⁿ = x_i^liqF^liq + y_i^vapF^vap

At the equilibrium condition, the fugacity of the vapor and liquid phase are defined as follows:

ln f_i^vap = ln f_i^liq

f_i^phase = y_i^phaseϕ_i^phaseP

$$\ln{\phi_{i}} = \frac{b_i}{b_m}(Z-1) - \ln{(Z-B)} + \frac{A}{B\sqrt{u^2-4w}}\left(\frac{b_i}{b_m}-\delta_i\right)\ln{\left(\frac{2Z+B(u+\sqrt{u^2-4w})}{2Z+B(u-\sqrt{u^2-4w})}\right)}$$

$$\delta_i = \frac{2a_i^{1/2}}{a_m}\sum_j{x_ja_j^{1/2}(1-k_{ij})}$$

The cubic equation of state is solved to find Z for each phase subject to the composition of that phase. Typically, the flash calculations are computed at a given temperature, T. However, the flash calculations become trivial if the given conditions do not fall in the two phase region. For simulation only studies, the user may know a priori the condition of the stream but when the same set of equations are used for optimization, there is a high probability that the specifications can transcend the phase envelope and hence the flash equations included may be trivial in the single phase region (i.e. liquid or vapor only). To circumvent this problem, property packages in IDAES that support VLE will compute the flash calculations at an "equilibrium" temperature T_eq. The equilibrium temperature is computed as follows:

T₁ = max(T_bubble, T)

T_eq = min(T₁, T_dew)

where T_eq is the equilibrium temperature at which flash calculations are computed, T is the stream temperature, T₁ is the intermediate temperature variable, T_bubble is the bubble point temperature of mixture, and T_dew is the dew point temperature of the mixture. Note that, in the above equations, approximations are used for the max and min functions as follows:

$$T_{1} = 0.5{[T + T_{bubble} + \sqrt{(T-T_{bubble})^2 + \epsilon_{1}^2}]}$$

$$T_{eq} = 0.5{[T_{1} + T_{dew} - \sqrt{(T-T_{dew})^2 + \epsilon_{2}^2}]}$$

where ϵ₁ and ϵ₂ are smoothing parameters (mutable). The default values are 0.01 and 0.0005 respectively. It is recommended that ϵ₁ > ϵ₂. Please refer to reference 4 for more details. Therefore, it can be seen that if the stream temperature is less than that of the bubble point temperature, the VLE calculations will be computed at the bubble point. Similarly, if the stream temperature is greater than the dew point temperature, then the VLE calculations are computed at the dew point temperature. For all other conditions, the equilibrium calculations will be computed at the actual temperature.

Other Constraints

Additional constraints are included in the model to compute the thermodynamic properties based on the cubic equation of state, such as enthalpies and entropies. Please note that, these constraints are added only if the variable is called for when building the model. This eliminates adding unnecessary constraints to compute properties that are not needed in the model.

All thermophysical properties are calculated using an ideal and residual term, such that:

p = p⁰ + p^r

The residual term is derived from the partial derivatives of the cubic equation of state, whilst the ideal term is determined using empirical correlations.

Enthalpy

The ideal enthalpy term is given by:

h_i⁰ = ∫_298.15^T(A + BT + CT² + DT³)dT + Δh_form^298.15K

The residual enthalpy term is given by:

$$h_{i}^{r}b_m\sqrt{u^2-4w} = \left(T\frac{da}{dT}-a_m\right)\ln{\left(\frac{2Z+B(u+\sqrt{u^2-4w})}{2Z+B(u-\sqrt{u^2-4w})}\right)} +RT(Z-1)b_m\sqrt{u^2-4w}$$

$$\frac{da}{dT}\sqrt{T} = -\frac{R}{2}\sqrt{\Omega_A}\sum_i{\sum_j{y_iy_j(1-k_{ij})\left(f_{w,j}\sqrt{a_i\frac{T_{c,j}}{P_{c,j}}}+f_{w,i}\sqrt{a_j\frac{T_{c,i}}{P_{c,i}}}\right)}}$$

Entropy

The ideal entropy term is given by:

$$s_{i}^{0} = \int_{298.15}^{T}\frac{(A+BT+CT^2+DT^3)}{T}dT + \Delta s_{form}^{298.15K}$$

The residual entropy term is given by:

$$s_{i}^{r}b_m\sqrt{u^2-4w} = R\ln{\frac{Z-B}{Z}}b_m\sqrt{u^2-4w} + R\ln{\frac{ZP^{ref}}{P}}b_m\sqrt{u^2-4w} + \frac{da}{dT}\ln{\left(\frac{2Z+B(u+\sqrt{u^2-4w})}{2Z+B(u-\sqrt{u^2-4w})}\right)}$$

Fugacity

Fugacity is calculated from the system pressure, mole fractions and fugacity coefficients as follows:

Fugacity Coefficient

The fugacity coefficient is calculated from the departure function of the cubic equation of state as shown below:

$$\ln{\phi_{i}} = \frac{b_i}{b_m}(Z-1) - \ln{(Z-B)} + \frac{A}{B\sqrt{u^2-4w}}\left(\frac{b_i}{b_m}-\delta_i\right)\ln{\left(\frac{2Z+B(u+\sqrt{u^2-4w})}{2Z+B(u-\sqrt{u^2-4w})}\right)}$$

$$\delta_i = \frac{2a_i^{1/2}}{a_m} \sum_j{x_j a_j^{1/2}(1-k_{ij})}$$

Gibbs Energy

The Gibbs energy of the system is calculated using the definition of Gibbs energy:

g_i = h_i − TΔs_i

List of Variables

Variable Name	Description	Units
flow_mol	Total molar flow rate	mol/s
mole_frac_comp	Mixture mole fraction indexed by component	None
temperature	Temperature	K
pressure	Pressure	Pa
flow_mol_phase	Molar flow rate indexed by phase	mol/s
mole_frac_phase_comp	Mole fraction indexed by phase and component	None
pressure_sat	Saturation or vapor pressure indexed by component	Pa
cp_mol_phase	Isobaric molar heat capacity by phase	J/mol/K
cv_mol_phase	Isochoric molar heat capacity by phase	J/mol/K
dens_mol_phase	Molar density indexed by phase	mol/m3
dens_mass_phase	Mass density indexed by phase	kg/m3
enth_mol_phase	Molar enthalpy indexed by phase	J/mol
enth_mol	Molar enthalpy of mixture	J/mol
entr_mol_phase	Molar entropy indexed by phase	J/mol.K
entr_mol	Molar entropy of mixture	J/mol.K
fug_phase_comp	Fugacity indexed by phase and component	Pa
fug_coeff_phase_comp	Fugacity coefficient indexed by phase and component	None
gibbs_mol_phase	Molar Gibbs energy indexed by phase	J/mol
heat_capacity_ratio_phase	Heat capcity ratio by phase	-
isothermal_speed_sound_phase	Isothermal speed of sound by phase	m/s
isentropic_speed_sound_phase	Isentropic speed of sound by phase	m/s
mw	Molecular weight of mixture	kg/mol
mw_phase	Molecular weight by phase	kg/mol
temperature_bubble	Bubble point temperature	K
temperature_dew	Dew point temperature	K
pressure_bubble	Bubble point pressure	Pa
pressure_dew	Dew point pressure	Pa
_teq	Temperature at which the VLE is calculated	K

List of Parameters

Parameter Name	Description	Units
cubic_type	Type of cubic equation of state to use, from CubicEoS Enum	None
pressure_ref	Reference pressure	Pa
temperature_ref	Reference temperature	K
omega	Pitzer acentricity factor	None
kappa	Binary interaction parameters for EoS (note that parameters are specific for a given EoS	None
mw_comp	Component molecular weights	kg/mol
cp_ig	Parameters for calculating component heat capacities	varies
dh_form	Component standard heats of formation (used for enthalpy at reference state)	J/mol
ds_form	Component standard entropies of formation (used for entropy at reference state)	J/mol.K
antoine	Component Antoine coefficients (used to initialize bubble and dew point calculations)	bar, K

Config Block Documentation

idaes.models.properties.cubic_eos.cubic_prop_pack

CubicParameterData

CubicStateBlock

CubicStateBlockData