Bubble Point and Dew Point Calculations for Non-Ideal Mixtures Using the Gamma-Phi Method.

The purpose of this article is to provide a practical, engineer-ready workflow for computing bubble point and dew point conditions of non-ideal mixtures using the gamma-phi approach, including robust iteration strategies, model choices, and implementation details for real VLE work.

1. Why bubble point and dew point matter in real process design.

Bubble point and dew point calculations are core steps in vapor–liquid equilibrium (VLE) work for distillation, flash drums, condensers, reboilers, and phase stability checks in process simulation.

For ideal mixtures at low pressure, Raoult’s law with ideal gas behavior can be adequate, but many real mixtures deviate due to non-ideal liquid behavior, non-ideal vapor behavior, or both.

The gamma-phi method directly accounts for non-ideality in the liquid phase via activity coefficients and in the vapor phase via fugacity coefficients.

Note : If your mixture is strongly associating, highly non-ideal, near an azeotrope, or at elevated pressure, the gamma-phi approach and a careful selection of activity and EOS models are necessary to avoid large equilibrium errors.

2. Definitions and problem statements.

2.1 Bubble point definition.

At the bubble point, a liquid mixture of known overall liquid composition first begins to form an infinitesimal amount of vapor at equilibrium.

Common specifications are bubble pressure at fixed temperature and liquid composition, or bubble temperature at fixed pressure and liquid composition.

2.2 Dew point definition.

At the dew point, a vapor mixture of known overall vapor composition first begins to form an infinitesimal amount of liquid at equilibrium.

Common specifications are dew pressure at fixed temperature and vapor composition, or dew temperature at fixed pressure and vapor composition.

2.3 Notation used in this post.

xi is the liquid mole fraction of component i.

yi is the vapor mole fraction of component i.

P is system pressure, and T is system temperature.

Pisat(T) is pure-component saturation pressure at temperature T.

γi(T, x) is the activity coefficient of component i in the liquid mixture.

φiV(T, P, y) is the vapor-phase fugacity coefficient of component i at system conditions.

φisat(T) is a fugacity-coefficient-like correction at saturation, often needed when vapor non-ideality is treated consistently.

Fugacity equality at equilibrium is enforced component-by-component.

3. The gamma-phi equilibrium condition.

The equilibrium condition is equality of component fugacities between vapor and liquid phases.

For each component i, fiV = fiL is imposed.

In the gamma-phi framework, liquid non-ideality is represented by activity coefficients, while vapor non-ideality is represented by fugacity coefficients from an equation of state.

3.1 Practical K-value form.

A common engineering form is written as yi = Ki xi.

In gamma-phi, Ki depends on T, P, x, and y through γi and φiV.

A widely used expression is.

K_i(T,P,x,y) = (γ_i(T,x) * P_i^sat(T) * Φ_i^corr(T,P)) / (P * φ_i^V(T,P,y)).

Φicorr is a compact way to represent optional saturation and pressure corrections that improve consistency, such as a Poynting correction and a saturated-vapor fugacity correction.

At low pressure, Φicorr and φiV often approach 1, and the equation reduces to the classic γ-Raoult relation.

3.2 Optional corrections that improve accuracy.

When pressure is not low, or when high accuracy is required, include pressure and vapor non-ideality corrections explicitly.

A common corrected form is.

K_i = (γ_i * P_i^sat / P) * (φ_i^sat / φ_i^V) * Π_i.

Here φisat is the fugacity coefficient of pure i at (T, Pisat) in the vapor, computed with the same EOS used for the mixture when consistency is desired.

Πi is a Poynting factor, typically written as exp(∫(v_i^L dP)/(RT)), often approximated using liquid molar volume over a modest pressure range.

Note : If you omit corrections, state the assumptions clearly, such as low-to-moderate pressure and negligible liquid compressibility, because these assumptions determine whether K-values remain physically consistent.

4. Choosing models for γ and φ in real VLE calculations.

4.1 Activity coefficient models for γ.

Activity coefficient models represent excess Gibbs energy behavior in the liquid phase.

Common choices include Wilson, NRTL, and UNIQUAC.

Wilson is often effective for fully miscible systems but can struggle with liquid–liquid equilibria.

NRTL is flexible and widely used for strongly non-ideal systems, including partial miscibility.

UNIQUAC is structurally grounded and often robust across diverse chemistries when parameters are available.

4.2 Equation of state choices for φiV.

Vapor fugacity coefficients are computed from an EOS such as Peng–Robinson or Soave–Redlich–Kwong.

At low pressure, ideal gas behavior may be acceptable and φiV can be set to 1 as a simplifying assumption.

At elevated pressure, EOS-based φiV becomes important for correct dew point and bubble point predictions.

4.3 Model matching guidance.

Situation. Liquid model for γ. Vapor model for φ. Practical note.
Low pressure, mild non-ideality. Wilson or NRTL or UNIQUAC. Ideal gas with φiV = 1. Often adequate for screening calculations.
Low pressure, strong non-ideality or azeotrope risk. NRTL or UNIQUAC. Ideal gas may still be adequate. Prioritize good binary parameters for γ.
Moderate to high pressure. NRTL or UNIQUAC. Peng–Robinson or SRK EOS. Include φiV and consider correction factors.
Near-critical or highly compressible vapor. As appropriate for liquid region. EOS mandatory. Convergence can be sensitive, so use damping.

5. Bubble pressure at fixed temperature and liquid composition.

Given T and x, the unknown is P, and the equilibrium vapor composition y is also unknown.

In the gamma-phi method, K-values depend on y through φiV, so a coupled iteration is used.

5.1 Governing condition for bubble point.

At equilibrium, yi = Ki xi, and vapor-phase composition must sum to 1.

The bubble point condition is therefore.

Σ (x_i * K_i(T,P,x,y)) = 1.

5.2 Robust iteration strategy for bubble pressure.

A practical approach is an outer loop on pressure and an inner loop on y for the EOS fugacity coefficients.

One common workflow is.

1) Inputs: T, x_i, property models for Psat(T), γ(T,x), EOS for φ^V(T,P,y). 2) Initialize P using an idealized estimate such as P ≈ Σ (x_i * γ_i * P_i^sat). 3) Initialize K_i ignoring φ effects, then set y_i ∝ x_i * K_i and normalize y. 4) Repeat until convergence. a) Compute γ_i(T,x) from the activity model. b) Compute P_i^sat(T) from vapor pressure correlation. c) With current P and y, compute φ_i^V(T,P,y) from EOS. d) Compute K_i from gamma-phi expression. e) Update y_i = x_i * K_i, then normalize y so Σ y_i = 1. f) Update pressure using bubble condition, for example: P_new = Σ (x_i * γ_i * P_i^sat * Φ_i^corr / φ_i^V) / Σ y_i_target, or more directly enforce Σ (x_i*K_i) = 1 by a scalar root solve on P. g) Apply damping: P := λ*P_new + (1-λ)*P, with 0<λ≤1. 5) Stop when |Σ(x_i*K_i) - 1| and |P_new - P|/P are below tolerance.
Note : If the EOS computation of φiV is sensitive, update pressure slowly using damping, because aggressive pressure jumps can destabilize the y-dependent fugacity coefficients.

6. Bubble temperature at fixed pressure and liquid composition.

Given P and x, the unknown is T, and the solution is found by enforcing Σ(xiKi) = 1 at the correct temperature.

Because Pisat(T) changes strongly with temperature, a temperature root solve is typically efficient.

6.1 Root function for bubble temperature.

F(T) = Σ (x_i * K_i(T,P,x,y(T))) - 1.

Here y depends on T through K-values and EOS fugacity coefficients.

6.2 Practical temperature solving approach.

Use a bracketing method when possible, because it is robust for highly non-ideal mixtures.

A typical approach is.

1) Bracket T_low and T_high such that F(T_low) < 0 and F(T_high) > 0. 2) For each trial T in the bracket: a) Compute Psat_i(T) and γ_i(T,x). b) Iterate y and compute φ_i^V if vapor non-ideality is included. c) Evaluate F(T). 3) Apply bisection or safeguarded secant until |F(T)| is within tolerance.
Note : For mixtures that form azeotropes, F(T) can be flat or have multiple near-roots in some regions, so bracketing and safeguards are preferable to pure Newton updates.

7. Dew pressure at fixed temperature and vapor composition.

Given T and y, the unknown is P, and the incipient liquid composition x is unknown.

The dew point condition is derived from xi = yi/Ki and Σ xi = 1.

7.1 Governing condition for dew point.

Σ (y_i / K_i(T,P,x,y)) = 1.

This condition is coupled because γi depends on x, and K-values can depend on y through φiV.

7.2 Robust iteration strategy for dew pressure.

A practical approach alternates between updating x from current K-values and updating K-values from updated x.

1) Inputs: T, y_i, models for Psat(T), γ(T,x), EOS for φ^V(T,P,y). 2) Initialize P using an idealized estimate such as: P ≈ 1 / Σ (y_i / (γ_guess * P_i^sat / P_ref)). 3) Initialize x_i ∝ y_i / K_i and normalize x. 4) Repeat until convergence. a) Compute γ_i(T,x) from activity model. b) Compute Psat_i(T). c) Compute φ_i^V(T,P,y) from EOS using known y. d) Compute K_i. e) Update x_i = y_i / K_i, normalize x so Σ x_i = 1. f) Enforce dew condition by scalar root solve on P so that Σ(y_i/K_i) = 1. g) Apply damping to P and optionally to x updates. 5) Stop when |Σ(y_i/K_i) - 1| and composition changes are below tolerance.
Note : Dew point calculations are often more sensitive than bubble point calculations because γ depends on x, which is itself computed from K-values, so apply damping to x updates when the mixture is strongly non-ideal.

8. Dew temperature at fixed pressure and vapor composition.

Given P and y, solve for T such that Σ(yi/Ki) = 1 with x updated self-consistently at each trial temperature.

8.1 Root function for dew temperature.

G(T) = Σ (y_i / K_i(T,P,x(T),y)) - 1.

8.2 Practical solving approach.

Use a bracketing method if possible, and at each trial temperature iterate x until Σ x = 1 is satisfied by normalization and K-value consistency.

Safeguarded secant is often efficient when vapor pressure correlations are smooth and monotonic across the bracket.

9. Implementation details that prevent common errors.

9.1 Unit consistency checklist.

Ensure P, Pisat, and EOS calculations use consistent pressure units.

Ensure the vapor pressure correlation returns absolute pressure, not gauge pressure.

Ensure the gas constant matches the unit system used for any exponential correction factors.

9.2 Numerical stability and convergence criteria.

Use both a scalar equilibrium residual and composition change criteria.

Typical scalar residuals are |Σ(xK)−1| for bubble problems and |Σ(y/K)−1| for dew problems.

Track maximum changes in y or x, because small scalar residual does not guarantee stable compositions when K-values are strongly composition-dependent.

9.3 Handling near-azeotrope behavior.

Near an azeotrope, x and y can become similar and K-values approach 1 for some components, which can slow convergence.

Use damping, and prefer bracketing methods for temperature solves.

Consider generating a good initial guess from a nearby known solution point to improve robustness.

9.4 When the gamma-phi approach may be insufficient.

If the mixture is at very high pressure and the liquid phase is significantly compressible, or if you need consistent treatment across all phases and wide conditions, an EOS-based phi-phi approach may be more appropriate.

Even then, gamma-phi remains a practical workhorse when liquid non-ideality dominates and reliable activity parameters exist.

Note : Do not mix parameter sets across models, because activity parameters are model-specific, and using mismatched parameters can produce physically inconsistent phase behavior.

10. Engineer-ready summary table of bubble and dew equations.

Case. Given. Unknown. Core condition. Key update relations.
Bubble pressure. T, x. P, y. Σ(xK)=1. y_i = x_i K_i, normalize y.
Bubble temperature. P, x. T, y. Σ(xK)=1. y_i = x_i K_i, normalize y, root solve in T.
Dew pressure. T, y. P, x. Σ(y/K)=1. x_i = y_i / K_i, normalize x.
Dew temperature. P, y. T, x. Σ(y/K)=1. x_i = y_i / K_i, normalize x, root solve in T.

11. Practical pseudo-code template for gamma-phi bubble and dew point.

The following template shows a consistent structure you can adapt into spreadsheets or code.

Common functions: - Psat_i(T): vapor pressure correlation for each component. - gamma_i(T, x): activity coefficient model (NRTL, Wilson, UNIQUAC). - phiV_i(T, P, y): EOS fugacity coefficients for vapor mixture (PR, SRK). - K_i(T,P,x,y) = (gamma_i * Psat_i(T) * Phi_corr_i) / (P * phiV_i).
BubbleP(T, x):

initialize P

initialize y from y_i ∝ x_i * (gamma_i*Psat_i/P)

loop:
compute gamma_i(T,x)
compute Psat_i(T)
compute phiV_i(T,P,y)
compute K_i
S = Σ(x_i * K_i)
update y_i = x_i * K_i; normalize y
update P using scalar root solve so S → 1
check convergence

DewP(T, y):

initialize P

initialize x from x_i ∝ y_i / (gamma_guess*Psat_i/P)

loop:
compute gamma_i(T,x)
compute Psat_i(T)
compute phiV_i(T,P,y)
compute K_i
D = Σ(y_i / K_i)
update x_i = y_i / K_i; normalize x
update P using scalar root solve so D → 1
check convergence

FAQ

What is the main advantage of the gamma-phi method for bubble point and dew point calculations.

The gamma-phi method separates liquid and vapor non-ideality into two physically meaningful factors, where γ captures excess Gibbs energy behavior in the liquid and φ captures real-gas behavior in the vapor.

This separation makes it practical for many chemical systems where liquid non-ideality is strong and high-quality activity parameter sets are available.

How should I choose between NRTL, Wilson, and UNIQUAC for activity coefficients.

Choose the model based on the availability of validated binary interaction parameters for your mixture and the expected behavior of the liquid phase.

NRTL is commonly selected for strongly non-ideal and partially miscible systems, Wilson is often used for fully miscible systems, and UNIQUAC is broadly applicable when parameters exist and structural consistency is desired.

When do I need to include vapor fugacity coefficients from an EOS.

Include EOS-based φiV when pressure is moderate to high, when the vapor phase is non-ideal, or when dew point predictions are sensitive to vapor non-ideality.

At low pressure, setting φiV to 1 can be an acceptable approximation for many systems, but you should verify that the approximation does not distort your design margins.

Why do dew point calculations sometimes converge more slowly than bubble point calculations.

Dew point calculations solve for an unknown liquid composition x, and γ depends on x, which is computed from K-values that depend on γ.

This feedback can create stronger nonlinearity than bubble calculations, so damping and safeguarded root solvers often improve stability.

What is the most common reason bubble point and dew point results look inconsistent in spreadsheets.

The most common reasons are unit inconsistencies for pressure and vapor pressure, incorrect normalization of x or y, and mixing parameter sets across incompatible thermodynamic models.

Another common issue is forgetting that K-values can depend on composition through φiV and γ, which requires iteration rather than a single-pass calculation.