Category: Distillation

Shortcut Distillation

Shortcut Distillation

Distillation is the most common unit operation for separating liquid mixtures into valuable and high purity products. The normal procedure for solving a typical multi-component distillation problem is to solve the MESH (Material balance, Equilibrium, Summation and Heat) balance equations stage-by-stage.

Though computer programs are normally available for the rigorous solution of the MESH equations, short-cut methods are still useful in the preliminary design work, and as an aid in defining problems for computer simulation.

This article describes a widely used shortcut distillation method commonly referred to as the Fenske-Underwood-Gilliland (FUG) method.

  1. The components that have their distillate and bottoms fractional recoveries specified are called key elements. The most volatile of the keys is called the Light Key (LK) and the least volatile is called the heavy key (HK). The other components are called non-keys(NK). First step is to choose LK and HK and fix their distribution in top and bottom product.
  2. Next step is to estimate overall top and bottom flow rate compositions assuming only light non-key components in top and heavy non-key component in bottoms as first attempt.
  3. Estimate Dew point TDew for top composition and Bubble Point TBubble for bottom composition.
  4. Estimate relative volatility with respect to the HK of all components at the TDew (top) and TBubble (bottom).
    αi = Ki / KHK

    Calculate average volatilities for all components.

    αi,Avg = (αi, Top * αi, Bottom)0.5
  5. Recalculate overall top and bottom flow rate compositions on the basis of the Hengsteback and Geddes equation.
    log(di/bi) = A + B*log(αi)

    Recovery of heavy key in bottoms is defined as bHK/fHK, recovery of light key in distillate is defined as dLK/fLK. Constants A & B are defined as:

    A = log((1-bHK/fHK)/(bHK/fHK))
    B = log[(dLK/fLK)/(1-dLK/fLK)*(bHK/fHK)/(1-bHK/fHK)]/log(αLK)

    Recovery of ith component in distillate is defined as:

    (di/fi) = (10AiB)/(1 + 10AiB)

    Recovery of ith component in bottom is defined as:

    (bi/fi) = 1 - (di/fi)

    After calculating top and bottom compositions, go back to Step (3) and update TDew and TBubble and recalculate relative volatilites.

  6. Estimate minimum number of stages, Nm using Fenske equation.
    Nm = log[(xLK/xHK)d.(xHK/xLK)b]/log(αLK)

    where, (xLK/xHK)d is the ratio of mol fraction of light key and heavy key in distillate and
    (xHK/xLK)b is the ratio of mol fraction of heavy key and light key in bottoms.

  7. Estimate minimum reflux ratio, Rm using Underwood’s equation.
    Rm + 1 = Σ(αi.xi,d)/(αi - θ)

    where Underwood constant, θ is calculated as the root of following equation:

    1 - q = Σ(αi.xi,f)/(αi - θ)

    where, 1 < θ < αLK and q is the quality of feed.

  8. Select an operating reflux ratio, R as a factor of minimum reflux ratio, typically a factor of 1.2 is considered. Number of stages, N is calculated based on Gilliland’s correlation.
    X = (R - Rm)/(R + 1)
    Y = (N - Nm)/(N + 1)
    Y = 1 - exp[(1 + 54.4*X)/(11 + 117.2*X)*(X - 1)/X0.5]
  9. Feed stage location is estimated based on Kirkbride’s equation.
    m/p = [(B/D)*(xHK/xLK)F*(xLK,B/xHK,D)²]0.206
    N = m + p

    where, m is number of stages above feed plate and n is number of stages below feed plate.

Web based calculation available at


Spreadsheet for Shortcut Multi-Component Distillation based on FUG method


Simple Multicomponent Distillation

Simple Multicomponent Distillation

This article shows application of Thiele and Geddes method to solve simple multi-component distillation problem in excel spreadsheet based on θ method of convergence.

Antoine equation is used to calculate equilibrium constant Ki.

Ki = exp(Ai - Bi / (T + Ci))/P

where, Ai, Bi & Ci are Antoine equation constants for a component i, T & P are temperature and pressure at a stage.

Feed Analysis

To start with bubble point pressure (PBubble) and dew point pressure (PDew) are determined for feed mixture.

  • P < PDew, Mixture exists as superheated vapor.
  • P > PBubble, Mixture exists as sub-cooled liquid.
  • PDew < P < PBubble, mixture exist in vapor and liquid phase.

For last case, PT Flash calculation is performed to determine liquid (Lf) & vapor (Vf) fraction and molar composition.

Stage-wise Calculation

Liquid and Vapor flow across column are determined as:

L (Above Feed Stage) = R.D
V (Above Feed Stage) = (1+R).D
L (Below Feed Stage) = R.D + Lf
V (Below Feed Stage) = (1+R).D - Vf

Start of Iteration
A Linear temperature gradient is assumed across the column. Ki values are calculated for all component on each stage based on Antoine equation.

θ Method of Convergence
L/D ratios are calculated for all components from top to stage above feed stage.

Top Stage
Total Condenser

LD(1,i) = L(1)/D

Partial Condenser

LD(1,i) = L(1)/(K(1,i).D)

where 1 indicates top stage, i denotes the component and LD denotes L/D ratio.

Above Feed Stage

LD(j,i) = (LD(j-1,i) + 1).L(j)/(K(j,i).V(j))

where j represents the stage and i represents a single component.

Similarly L/B ratios are calculated from bottom upto feed stage

Bottom Reboiler

LB(bottom) = 1

L/B ratios from bottom upto feed tray is calculated as:

LB(j,i) = LB(j+1,i).K(j+1,i).V(j+1)/L(j+1)  + 1 

For Feed Stage

Vfb(i) = V(f).K(f,i)/( L(f).LB(f,i) )

B/D ratio for feed plate is calculated as:

bd(i) = ( LD(nf-1,i) + Lf.Xfi/(F.Zfi))/(Vfb(i) + Vf.Yfi/(F.Zfi))

where bd represent B/D ratio at feed plate.

Based on overall component material balance.

F.Xi = D.Xdi + B.Xbi
d(i) = D.Xdi
b(i) = B.Xbi
F.Xi = d(i) + b(i)
Σd(i) = D

A multiplier θ is defined as:

[b(i)/d(i)]corrected = θ[b(i)/d(i)]calculated

Above equation is rewritten as:

bd(i)co = θ.bd(i)ca

Function g(θ) is defined as:

d(i)co = FXi/(1+θ.bd(i)ca)
g(θ)  = Σd(i)co - D
      = Σi FXi/(1+θ.bd(i)ca) - D

Above equation is solved using Newton Raphson method.

g'(θ) = -Σi bd(i)ca.FXi/[1+θ.bd(i)ca

New estimate of θ is made as following:

θnew = θ - g(θ)/g'(θ)

Iteration is done till there is negligible change in value of θ. Converged value of θ is used to estimate d(i)co. Value of b(i)co is calculated as following:

b(i)co = θ.d(i)

Liquid mol fraction from top to stage above feed stage are calculated:

X(j,i) = LD(j,i)ca.d(i)co / Σi LD(j,i)ca.d(i)co

Liquid mol fraction from feed to bottom tray are calculated:

X(j,i) = LB(j,i)ca.b(i)co / Σi LB(j,i)ca.b(i)co

Bubble point temperature is calculated for each stage. With this first iteration completes. In next iteration new set of Ki values are calculated for each tray based on updated temperature profile and all of the above steps are repeated till there is no change in temperature profile and θ value becomes 1.

Spreadsheet for Simple Multicomponent Distillation


  • Fundamentals of Multicomponent Distillation – C.D. Holland
Binary Batch Distillation

Binary Batch Distillation

Batch distillation is widely used in chemical processing industries when high value added, low volume chemicals must be separated. This article shows how to model and solve a multi-stage binary batch distillation problem in excel.

The 50 mol feed comprising 70% mols of A and 30% mols of B is to be distilled in the multi-stage batch distillation with 5 equilibrium stages on top of the reboiler (still pot). Reflux is returned to the column as a saturated liquid with the constant reflux ratio of 1.5. Relative volatility of A is 1.8. It is desired to reduce the mols of A in still to 10%. Determine average distillate composition, final amount of liquid in still pot and total amount of distillate collected.

A material balance is done across the column and following relationship commonly known as Rayleigh Equation is derived.

Rayleigh Equation

Where F (Initial mols fed in still), WF (Mols left in still at any time), xF (Initial mol fraction in feed), xW (Mol fraction in still at any time) and xD (Distillate composition at any particular time).

Equilibrium curve is obtained from relative volatility.

 y = αx / (1+ x (α - 1))

Operating line is obtained as following.

 y = Rx /(R+1) + xD /(R+1)

Relationship between xD and xW is obtained by performing stage-by-stage calculations. A value of xD is selected and stages are stepped off (for a given number of equilibrium stages) to find the value of xW. It is done for all possible values of xD ranging from 0 to 1. A plot of 1/(xD – xW) vs xW is made.

Plot for 1/(y-x) vs x Rayleigh Equation

Area under curve is calculated from xW to xF by calculating area of trapezoid formed in small intervals.

Area under curve in Excel

 Area (Cell P59) = (N60-N59)*(O60+O59)/2)

Area of all trapezoids in the range from xW to xF are added to get area under curve.

 Area = 2.153

Amount left in still, WF is calculated as following.

 WF = F.exp(- Area under the curve)
 WF = 5.80 mols

Distillate mols collected is obtained as following.

 D = F - WF
 D = 44.20 mols

Average composition of distillate collected xD is calculated as following.

 xD = (xF F - xW WF)/D
 xD = 0.779 

Spreadsheet for Multi-stage Binary Batch Distillation

McCabe Thiele Diagram

McCabe Thiele Diagram

This article shows how to draw McCabe Thiele diagram for doing binary distillation analysis in a spreadsheet.

Determine total number of stages and feed stage for separation of a binary mixture with relative volatility (α) 2.5, Feed composition (zF) 0.36, Distillate Composition (xD) 0.915, Bottom composition (xB) 0.05 and feed quality (q) 1.5. Reflux ratio to be 1.5 times the minimum reflux ratio (Rm).

Minimum Reflux Ratio (Rm)

Calculate the intersection of q-Line and equilibrium curve.

 y = qx /(q-1) - zF /(q-1)
 y = αx / (1+ x (α - 1))

Solving above equations gives a quadratic expression, which is solved for positive roots.

 x² (q(α-1)) + x (q - zF (α -1) - α(q-1)) - zF = 0
 x1 = 0.470
 y1 = 0.689

Draw a line from ( xD, yD) passing through ( x1, y1) with an intercept xD / (Rm + 1) on Y axis.

Minimum Reflux Ratio intercept on Y axis

 xD / (Rm + 1) = 0.450
 Rm = 1.032

Reflux ratio is calculated as 1.5 times the minimum reflux ratio.

 R = 1.5 x Rm
 R = 1.55

Operating Lines

Intersection of q-Line and rectification section operating line is calculated.

 y = qx /(q-1) - zF /(q-1)
 y = Rx /(R+1) + xD /(R+1)
 x2 = 0.451
 y2 = 0.633

Draw rectification section operating line from (xD, yD) to (x2, y2) and stripping section operating line from (x2, y2) to (xB, yB).

Operating lines on McCabe Thiele Diagram

Draw Stages

Start drawing steps from (xD, yD). Draw a horizontal line by calculating x value from equilibrium curve based on yD. Check whether this x value is greater than x2. As long as it is greater than x2 rectification section operating line is used for calculation of y value. For next step this new value of y is used and rest of the steps are repeated. If x becomes less than x2 stripping section operating line is used for calculation of y value. It is continued till it is greater than xB.

Stage Calculation on McCabe Thiele Diagram

Feed stage is determined by checking at each stage when x value becomes less than x2, similarly total number of stages are determined by checking at each stage when x value become less than xB. At last stage value of x is calculated and difference with xB is added to give total number of stages.

 Feed Stage   = 5
 Total Stages = 11.26

McCabe Thiele Diagram constructed in excel spreadsheet

Spreadsheet for Binary Distillation analysis based on McCabe Thiele method