## Jacketed Vessel Heat Transfer (Half Pipe Coil)

Agitator equipped vessels with half pipe coil jackets are widely used in variety of process applications. This article shows how to calculate heat transfer in an agitated vessel provided with an external half pipe coil jacket.

Overall heat transfer coefficient, U is defined as

``1/U = 1/hi + ffi + x/k + ffo + 1/ho``

where,

• hi : film coefficient process side
• ho : film coefficient coil side
• ffi : fouling factor process side
• ffo : fouling factor coil side
• x : vessel wall thickness
• k : vessel wall thermal conductivity

Featured Resources :

Hydrocarbon Processing the premier magazine providing job-help information to technical and management personnel in petroleum refining, gas processing, petrochemical/chemical and engineer/constructor companies throughout the world – since 1922.

### Process Side, hi

Process side film coefficient, hi depends upon type of impellor, Reynold’s number (Re) and Prandtl number (Pr).

````Re = D².N.ρ / μ`
`Pr = Cp.μ/k````

where, D is impellor diameter, N is impellor rpm, ρ is fluid density, μ is fluid viscosity, Cp is fluid specific heat and k is fluid thermal conductivity at bulk fluid temperature. hi is defined as following :

````Nui = C.Rea. Prb. (μ/ μw)c. Gc `
`hi.DT/k = Nui````

where constants C, a, b & c are available in literature for different type of impellors. DT is vessel diameter. Gc is geometric correction factor for non-standard geometries. μ/ μw is viscosity correction factor due to difference in viscosities at bulk fluid and wall temperatures. These constants are available in references mentioned below.

### Coil Side, ho

Pipe coils are made with a 180° central angle or a 120° central angle. Equivalent diameter (De) and Flow area (Ax) is
defined as following.

#### 180° Coil

````De = (π/2).dci`
`Ax = (π/8).dci²````

#### 120° Coil

````De = 0.708 dci`
`Ax = 0.154 dci²````

where, dci is inner diameter of pipe.

Reynold’s and Prandtl number are calculated based on jacket fluid properties and velocities.

````Re = De.v.ρ / μ`
`Pr = Cp.μ / k````

where, ρ is coil fluid density, μ is coil fluid viscosity and k is coil fluid thermal conductivity. v is fluid velocity in coil.

Featured Resources :

Serves chemical engineering professionals in the chemical process industry including manufacturing, engineering, government, academia, financial institutions and others allied to the field serving the global chemical process industry.

#### Turbulent Flow

For Re > 10000

``Nuc = 0.027 Re0.8 Pr0.33 (μ/ μw)0.14 (1 +3.5 De/Dc)``

where Dc is the mean or centerline diameter of the coil. Coil outer diameter (Do) is determined as following.

````180° Coil : Do = DT + 2(dci/2) + 2.x`
`120° Coil : Do = DT + 2(dci/4) + 2.x````
``Dc = ( Do + DT) / 2``

#### Laminar Flow

For Re < 2100

``Nuc = 1.86 [ Re.Pr.De/L ]0.33 (μ/ μw)0.14``

where, L is coil length along the vessel.

#### Transient Flow

For 2100 < Re < 10000
Calculate NuLaminar based on Re = 2100 and NuTurbulent based on Re = 10,000.

``NuTransient = NuLaminar + (NuTurbulent - NuLaminar).(Re - 2100)/(10000 - 2100)``

ho is determined as following

``ho.De/k = Nuc``

ho and hi thus calculated are used to get value of overall heat transfer coefficient U.

Spreadsheet for Half Pipe Coil Agitated Vessel Heat Transfer

### References

• Heat Transfer Design Methods – John J. McKetta (1992)
• Heat Transfer in Agitated Jacketed Vessels – Robert F. Dream, Chemical Engineering, January 1999

## Estimating Binary Interaction Parameter for Peng Robinson EOS

This article shows how to estimate binary interaction parameter (BIP) used in Peng Robinson (PR) Equation of State (EOS) from experimental data by regression in excel spreadsheet.

Example
Determine binary interaction parameter used in Peng Robinson EOS for a mixture of acetone and chloroform from Tx-y experimental data available at 101.33 kPa.

Obtain pure component properties like critical temperature (Tc), critical pressure (Pc) and accentric factor ω from literature. Experimental K1 is determined from x, y data.

``K1Experimental = y1 / x1``

Peng Robinson EOS parameters are calculated.

````κi = 0.37464 + 1.54226ω - 0.26992ω²`
`αi = [ 1 + κi (1 - (T/Tc)0.5)]²`
`ai = 0.45724 (RTc)²α / Pc`
`bi = 0.07780 RTc / Pc````

Mixture parameters are calculated next

````a12= [(a1.a2)0.5(1 - k12)] = a21`
`a = a1.x1² + 2.a12.x1.x2 + a2.x2²`
`b = b1.x1 + b2.x2`
`A = aP/(RT)²`
`B = bP/RT````

where k12 is Binary Interaction Parameter (BIP) to be estimated. For first run, use a guess value to start calculation.

Featured Resources :

Hydrocarbon Processing the premier magazine providing job-help information to technical and management personnel in petroleum refining, gas processing, petrochemical/chemical and engineer/constructor companies throughout the world – since 1922.

Following cubic equation is solved to get ZL.

``Z³ + (B-1)Z² + (A-3B² -2B)Z + (B³+B²-AB) = 0``

Above equation can be written as

``Z³ + C2.Z² + C1.Z + C0 = 0``

### Solving Cubic Equation

Cubic equation is solved using following procedure. Calculate Q1, P1 & D.

````Q1 = C2.C1/6 - C0/2 - C2³/27`
`P1 = C2²/9 - C1/3`
`D = Q1² - P1³````

If D >= 0, then equation has only one real root provided by

``Z1 = (Q1 + D0.5)1/3 + (Q1 - D0.5)1/3 - C2/3``

If D < 0, then equation has 3 real roots, following parameters are calculated

````t1 = Q1² / P1³`
`t2 = (1 - t1)0.5 / t10.5. Q1/abs(Q1)`
`θ = atan(t2)````

Roots are calculated as following –

````Z0 = 2.P10.5.cos(θ/3) - C2/3`
`Z1 = 2.P10.5.cos((θ + 2*Π)/3) - C2/3`
`Z2 = 2.P10.5.cos((θ + 4*Π)/3) - C2/3````

Roots thus calculated are arranged in descending order, highest root gives ZV and lowest root gives ZL.

### Fugacity

Based on ZL, liquid fugacity φiL is calculated for each data set.

In next step, Vapor phase fugacity is calculated. Mixture properties are estimated as following –

````a = a1.y1² + 2.a12.y1.y2 + a2.y2²`
`b = b1.y1 + b2.y2`
`A = aP/(RT)²`
`B = bP/RT````

Cubic equation is solved as per method shown above to get ZV.

``Z³ + (B-1)Z² + (A-3B² -2B)Z + (B³+B²-AB) = 0``

Based on ZV, vapor fugacity φiV is calculated for each data set.

Featured Resources :

Serves chemical engineering professionals in the chemical process industry including manufacturing, engineering, government, academia, financial institutions and others allied to the field serving the global chemical process industry.

K1 is calculated as:

``K1Calculated = φ1L/φ1V``

Difference of K1Experimental and K1Calculated is obtained for all data points. An objective function is defined as summation of all these differences. Click on Solver in Data Ribbon (Excel 2010) to open dialog box for Solver parameters and input data as shown below.

Minimize the objective function by changing values of kij. Uncheck Make Unconstrained Variables Non-Negative, as k12 can take negative value. Click solve to start regression and new value k12 is calculated.

`` k12 = -0.0605``

Above value is then used to calculate y1 values and results are plotted to check the deviation.

Example
Determine temperature dependent binary interaction parameter used in Peng Robinson EOS ( k12 = k1 + k2.T + k3/T ) for a mixture of acetone and n-Hexane from P x-y experimental data available at 318.15 °K.

Follow above steps to calculate pure component and mixture parameters. Change formula for kij to make it temperature dependent. As an initial guess put values for 3 parameters (k1,k2,k3) which gives a numeric difference value. In solver select these 3 variables to be varied for minimizing the difference value and then solve it. After doing regression binary interaction parameters are obtained and result is plotted as following.

````k1 = -0.0218`
`k2 = 3.142x10-4`
`k3 = 10.155````

Note : In case of non-convergence, try different initial guess values for kij and then let solver obtain optimized values.

Spreadsheet for Regressing BIP for Peng Robinson EOS

## 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````

Featured Resources :

Hydrocarbon Processing the premier magazine providing job-help information to technical and management personnel in petroleum refining, gas processing, petrochemical/chemical and engineer/constructor companies throughout the world – since 1922.

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)co.bd(i)ca``

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.

### References

• Fundamentals of Multicomponent Distillation – C.D. Holland

## Shortcut Sizing for Air Cooled Heat Exchanger

An air-cooled exchanger is used to cool fluids with ambient air. 1 inch OD tube is the most popular diameter and the most common fins are 1/2 inch or 5/8 inch high. Tube configuration used in this guideline : 1 inch OD tube, 5/8 inch fin height, 10 fins per inch and 2.5 inch triangular pitch.

APSF – External area in ft²/ft² of bundle face area.

Rows APSF Face Velocity, ft/min
3 80.4 700
4 107.2 660
5 134.0 625
6 160.8 600

APF – Total external area/ ft of fintube in ft²/ft ~ 5.58.

Face Velocity – Typical air face velocities (VFace ) used in design are tabulated above, these value result in optimum cost of exchanger.

Obtain Process Parameters
Obtain process duty (Q), hot process side inlet (T1) and outlet (T2) temperature. Select an overall heat transfer coefficient (U) from literature based on type of fluids. Select an air inlet temperature (t1) that is not exceeded during a certain percentage of time over the year (e.g. 95% of the time).

Featured Resources :

Hydrocarbon Processing the premier magazine providing job-help information to technical and management personnel in petroleum refining, gas processing, petrochemical/chemical and engineer/constructor companies throughout the world – since 1922.

Calculate Air Density
Density at air inlet temperature and site elevation

````ρo = 14.696 x 29 /(10.7316 x (t1 + 459.67))`
`ρAir / ρo = exp(-29 x z/ (1545 x (t1 + 459.67)))````

where,
t1: Air Inlet Temperature, °F
z : Elevation above sea level, feet

Calculate MTD
Assume an Air outlet temperature (t2) and calculate LMTD.

````LMTD = ((T1-t2) -(T2-t1))/ ln((T1-t2)/(T2-t1))`
`R = (T1 - T2)/(t2 - t1)`
`S = (t2 - t1)/(T1 - t1)````

LMTD Correction factor (F) is estimated based on following graphs

Cross Flow Single Pass

Cross Flow Two Pass

Calculate Air Flowrate
Finned area is estimated using following equation :

``AFinned = Q / (U * F * LMTD)``

Bundle Face area is estimated using following equation :

``AFace = AFinned / APSF``

Air Flow is estimated using following relation :

``VAir = AFace * VFace``

Air Mass flowrate is estimated :

``MAir = VAir * ρAir``

Air temperature rise is calculated :

``ΔT = Q / (MAir * CpAir )``

Revised air outlet temperature is calculated :

``t2 = t1 + ΔT``

This temperature is again used in above steps to re-estimate air outlet temperature. Above steps are iterated till there is no change in air outlet temperature.

Air Cooler Dimension
Air cooler width is calculated :

``W = AFace / LTube``

Number of Tubes are calculated :

``NTube = AFinned / (APF * LTube)``

Number of Tubes per Row are calculated :

``Nr Tube = NTube / No of Rows``

where,
LTube : Length of Tube

Air Side Pressure Drop
Air Side pressure drop is calculated as following :

``ΔP Total = ΔPStatic + ΔPVelocity``

Static Pressure Drop
Pressure drop across tube bundle.

````ΔP Static = FP * No of Rows / DR`
`FP = 6*10-8 * ( GFace )1.825````

where,
DR : ρAir / ρAir at seal level and 70°F
GFace : Air face mass velocity in lb/h.ft² face area
ΔPStatic : Static pressure drop in inch of H2O

Velocity Pressure Drop
Typically 2 fans are used in air cooler. Fan area per fan (FAPF) is calculated as following :

``FAPF = 0.4 * AFace / No. of Fans``

Fan Diameter is calculated :

````D = (4 * FAPF / π )0.5`
`ΔPVelocity = (ACFM / (4005 * (π* D2/4)) )2* DR````

where,
ACFM : Air flowrate in Actual Cubic Feet per Minute
ΔPVelocity : Velocity pressure drop in inch of H2O

Power Calculation
Break power for fan is calculated as following :

``BHP = ΔP Total * ACFM / 6356 / ηFan``

Motor power is calculated as following :

``Power = BHP / ηMotor``

Spreadsheet for Shortcut Air Cooled Heat Exchanger Design

### References

1. Process Heat Transfer: Principles and Applications,2007, Robert W Serth
2. Handbook of Chemical Engineering Calculations, Nicholas P Chopey
3. Rules of Thumb for Chemical Engineers, Carl R Branan
4. GPSA, Engineering Databook, 12th Edition FPS

## Double Pipe Heat Exchanger Design

This article shows how to do design for Double Pipe Heat Exchanger and estimate length of double pipe required.

Obtain flowrate (W ), inlet, outlet temperatures and fouling factor for both hot and cold stream. Calculate physical properties like density (ρ), viscosity (μ), specific heat (Cp) and thermal conductivity (k) at mean temperature. Determine heat load by energy balances on two streams.

````Q = mH.CpH(THot In - THot Out)`
`  = mC.CpC(tCold Out - tCold In)````

where,
mH , mC: Mass flow rate of Hot and Cold Stream
CpH , CpC: Specific Heat of Hot and Cold Stream
THot In , THot Out: Inlet and outlet temperature of Hot Stream
tCold In , tCold Out: Inlet and outlet temperature of Cold Stream

Calculate Logarithmic Mean Temperature Difference (LMTD)

``LMTD = (ΔT1 - ΔT2)/ln( ΔT1 / ΔT2)``

For Counter-current flow

````ΔT1 = THot In - tCold Out`
`ΔT2 = THot Out - tCold In````

For Co-current flow

````ΔT1 = THot In - tCold In`
`ΔT2 = THot Out - tCold Out````

Featured Resources :

Serves chemical engineering professionals in the chemical process industry including manufacturing, engineering, government, academia, financial institutions and others allied to the field serving the global chemical process industry.

Calculate Film Coefficient
Allocate hot and cold streams either in inner tube or annular space. General criteria for fluid placement in inner tube is corrosive fluid, cooling water, fouling fluid, hotter fluid and higher pressure stream. Calculate equivalent diameter (De) and flow area (Af) for both streams.

Inner Tube

````De = Di`
`Af = π Di²/4````

Annular Space

````De = D1 - Do`
`Af = π (D1² - Do²)/4````

where,
Di : Inside Pipe Inner Diameter
Do : Inside Pipe Outer Diameter
D1 : Outside Pipe Inner Diameter

Calculate velocity (V), Reynolds No. (Re) and Prandtl No. (Pr) number for each stream.

````V = W / ( ρ Af )`
`Re = De V ρ / μ`
`Pr = Cp μ / k````

For first iteration a Length of double pipe exchanger is assumed and heat transfer coefficient is calculated. Viscosity correction factor (μ / μw)0.14 due to wall temperature is considered 1.

For Laminar Flow (Re <= 2300), Seider Tate equation is used.

``Nu = 1.86 (Re.Pr.De/L )1/3(μ/ μw)0.14``

For Transient & Turbulent Flow (Re > 2300), Petukhov and Kirillov equation modified by Gnielinski can be used.

````Nu = (f/8)(Re - 1000)Pr(1 + De/L)2/3/[1 + 12.7(f/8)0.5(Pr2/3 - 1)]*(μ/μw)0.14`
`f  = (0.782* ln(Re) - 1.51)-2````

where,
L : Length of Double Pipe Exchanger
μw : Viscosity of fluid at wall temperature
Nu : Nusselts Number (h.De / k)

Estimate Wall Temperature

Wall temperature is calculated as following.

``TW = (hitAve + hoTAveDo/Di)/(hi + hoDo/Di)``

where,
hi : Film coefficient Inner pipe
ho : Film coefficient for Annular pipe
tAve : Mean temperature for Inner pipe fluid stream
TAve : Mean temperature for Annular fluid stream

Viscosity is calculated for both streams at wall temperature and heat transfer coefficient is multiplied by viscosity correction factor.

Overall Heat Transfer Coefficient
Overall heat transfer coefficient (U) is calculated as following.

``1/U = Do/hi.Di + Do.ln(Do/Di)/2kt + 1/ho+ Ri.Do/Di + Ro``

where,
Ri : Fouling factor Inner pipe
Ro : Fouling factor for Annular pipe
kt : Thermal conductivity of tube material

Calculate Area and length of double pipe exchanger as following.

````Area = Q / (U * LMTD )`
`L = Area / π * Do````

Compare this length with the assumed length, if considerable difference is there use this length and repeat above steps, till there is no change in length calculated.

Number of hair pin required is estimated as following.

``N Hairpin = L / ( 2 * Length Hairpin )``

Calculate Pressure Drop
Pressure drop in straight section of pipe is calculated as following.

``ΔPS = = f.L.G²/(7.5x1012.De.SG.(μ/ μw)0.14)``

where,
ΔP : Pressure Drop in PSI
SG : Specific Gravity of fluid
G : Mass Flux ( W / Af ) in lb/h.ft²

For Laminar flow in inner pipe, friction factor can be computed as following.

``f = 64/Re``

For Laminar flow in annular pipe.

````f = (64 / Re) * [ (1 - κ²) / ( 1 + κ² + (1 - κ²) / ln κ) ]`
`κ = Do / D1````

For turbulent flow in both pipe and annular pipe

``f = 0.3673 * Re -0.2314``

Pressure Drop due to Direction Changes

For Laminar Flow

``ΔPR = 2.0x10-13. (2NHairpin - 1 ).G²/SG``

For Turbulent Flow

``ΔPR = 1.6x10-13. (2NHairpin - 1 ).G²/SG``

Total Pressure Drop

``ΔPTotal = ΔPS + ΔPR``

Spreadsheet for Double Pipe Exchanger Design

## Flash Calculation (Raoult’s Law)

This article shows step by step procedure to do Bubble Point, Dew Point and Flash Calculation based on Raoult’s Law.

### Bubble Point Calculation

Bubble point of a system is the temperature at which liquid mixture begins to vaporize.

Obtain process parameters
Get liquid mixture molar composition ( Xi ) and Pressure (P) of the system. Obtain equilibrium ratios ( Ki ) for the components. Ki can be calculated from Antoine equation.

````Yi = Ki.Xi`
`Ki(T) = (e A - B / (T + C) )/ P````

where Yi is vapor phase molar composition in equilibrium with liquid and A,B,C are Antoine equation constants.

Calculate Bubble Point
At Bubble Point temperature summation of vapor phase molar fraction should be 1.

``Σ Ki(T).Xi = 1``

Above equation can be solved iteratively using Newton Raphson method. An initial temperature T is assumed. Function F(T) is calculated as following.

``F(T) = Σ Ki(T).Xi - 1``

Derivative of F(T) is calculated as following.

``F'(T) = Σ (B.Ki(T)/(T+C)² ).Xi``

New estimate of temperature is calculated as following.

``TNew = T - F(T)/F'(T)``

Function F(T) and F'(T) are calculated based on new temperature and this process is repeated till there is negligible difference in between T and TNew. Bubble point temperature thus obtained is then used to calculate vapor phase composition based on above equations.

Featured Resources :

Serves chemical engineering professionals in the chemical process industry including manufacturing, engineering, government, academia, financial institutions and others allied to the field serving the global chemical process industry.

### Dew Point Calculation

Dew point is the temperature at which liquid begins to condense out of the vapor.

Obtain process parameters
Get vapor mixture molar composition ( Yi ) and Pressure (P) of the system. Obtain equilibrium ratios ( Ki ) for the components.

````Yi = Ki.Xi`
`Ki(T) = (e A - B / (T + C) )/ P````

Calculate Dew Point
At Dew Point temperature summation of liquid phase molar fraction should be 1.

``Σ Yi / Ki(T) = 1``

Above equation can be solved iteratively using Newton Raphson method. An initial temperature T is assumed. Function F(T) is calculated as following.

``F(T) = Σ Yi / Ki(T) - 1``

Derivative of F(T) is calculated as following.

``F'(T) = Σ -Yi.( B/(Ki.(T+C)² ))``

New estimate of temperature is calculated as following.

``TNew = T - F(T)/F'(T)``

Function F(T) and F'(T) are calculated based on new temperature and this process is repeated till there is negligible difference in between T and TNew. Dew point temperature thus obtained is then used to calculate liquid phase composition based on above equations.

### Flash Calculation

A mixture when flashed to conditions between bubble and dew point separates in vapor and liquid phases. Flash calculation is done to determine vapor fraction and composition of liquid, vapor formed when a mixture is flashed at a given pressure and temperature.

Obtain process parameters
Get molar composition ( Zi )of the mixture and flash conditions mainly pressure (P) and temperature (T) of the system. Obtain equilibrium ratios ( Ki ) for the components.

``Yi = Ki.Xi``

Solve Flash Equations
Based on material balance on the system

````1 = V + L`
`Zi = V.Yi + L.Xi````

where V & L are vapor and liquid fractions. Solving above equations for Xi gives

``Xi = Zi / ( V.( Ki - 1) + 1 )``

At Flash conditions

``0 = Σ Yi - Σ Xi``

Above equation can be solved iteratively using Newton Raphson method. An initial vapor fraction V is assumed. Function F(V) is calculated as following.

````F(V) = Σ Yi  - Σ Xi`
`     = Σ [Zi.(Ki - 1)/(V.(Ki - 1) + 1)]````

Derivative of F(V) is calculated as following.

``F'(V) = Σ -[Zi.(Ki - 1)² /( V.(Ki - 1) + 1)²]``

New estimate of vapor fraction is calculated as following.

``VNew = V - F(V)/F'(V)``

Function F(V) and F'(V) are calculated based on new vapor fraction and this process is repeated till there is negligible difference in between V and VNew. Vapor fraction thus obtained is then used to estimate vapor and liquid molar composition based on above equations.

Spreadsheet for Flash Calculation based on Raoult’s Law

## Heat Exchanger Analysis

Heat Exchanger Analysis based on Effectiveness (ε) – NTU method is done when inlet temperatures are known and outlet temperatures are to be determined.

#### Obtain Process Parameters

Get process stream mass flowrate (M), specific heat (Cp) and inlet temperature (T). Obtain the heat transfer area (A) and overall heat transfer coefficient (U) for the given dimensions of heat exchanger.

Calculate heat capacities and obtain the minimum heat capacity.

````CH = MH * CpH`
`CC = MC * CpC`
`CMin = Minimum (CH , CC)`
`CR = CMin / CMax````

where MH, MC are hot and cold fluid mass flowrate; CpH, CpC are hot and cold fluid specific heat.

#### Calculate NTU and QMax

Number of transfer units ( NTU ) is calculated using following equation :

``NTU = U.A/ CMin``

Maximum heat transfer rate ( QMax ) is calculated using following equation :

``QMax = CMin.(THot In - TCold In)``

Featured Resources :

Serves chemical engineering professionals in the chemical process industry including manufacturing, engineering, government, academia, financial institutions and others allied to the field serving the global chemical process industry.

#### Determine Effectiveness

Based on NTU and CR (Ratio of heat capacities) determine heat exchanger effectiveness (ε) from Effectiveness – NTU curves available in literature.

ε – NTU Curve for Cross Flow exchanger Both stream unmixed

#### Calculate Outlet Temperature

Heat exchanger duty is calculated as:

``Q = ε * Q Max``

Outlet temperature are estimated as following :

````THot Out = THot In - Q /( MH.CpH)`
TCold Out = TCold In + Q /( MC.CpC)```

Example
Hot exhaust gases 1.5 kg/s, enter a finned-tube, cross-flow heat exchanger at 250°C and is used to heat pressurized water at a flow rate of 1 kg/s available at 35°C. Exhaust gas specific heat is 1000 J/kg.K and water specific heat is 4197 J/kg.K. Overall heat transfer coefficient is 100 W/m².K and area is 40 m². What is rate of heat transfer by exchanger and gas and water outlet temperatures ?

Calculate Heat Capacities

````CH = 1.5 kg/s * 1 kJ/kg.°K `
`   = 1.5 kW/ °K`
`CC = 1.0 kg/s * 4.197 kJ/kg.°K `
`   = 4.197 kW/ °K`
`CR = 0.36````

Determine NTU and Q Max

````NTU = 2.67`
`QMax = 322.5 kW````

Determine Effectiveness
Based on ε – NTU curves for cross flow heat exchanger, determine ε

``ε = 0.845``

Calculate Outlet Temperature
Heat duty is calculated as following –

````Q = 272.36 kW`
`THot Out = 68.43 ° C`
`TCold Out = 99.89 ° C````

All above calculation for some common exchanger geometries have been provided in below spreadsheet.

## PT Flash Calculation using PR EOS

PT Flash calculation determines split of feed mixture F with a molar composition Zi, into Vapor V and Liquid L at pressure P and temperature T. These calculations can be done in a excel spreadsheet using Peng Robinson Equation of State (PR EOS). To start with bubble point pressure (PBubble) and dew point pressure (PDew) are determined for feed mixture.

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

Initial guess of vapor fraction V and Ki is made as following.

````V = (PBubble - P)/(PBubble - PDew)`
`Ki = exp[ ln(Pc/P) + ln(10)(7/3)(1 + ω )(1-T/Tc)]````

Based on initial Ki values, iteration is done to get value of V which satisfies material balance on system.

````Yi = Ki.Xi`
`1 = V + L`
`Zi = V.Yi + L.Xi````

where V & L are vapor and liquid fractions. Solving above equations for Xi gives :

``Xi = Zi / ( V.( Ki - 1) + 1 )``

Featured Resources :

Hydrocarbon Processing the premier magazine providing job-help information to technical and management personnel in petroleum refining, gas processing, petrochemical/chemical and engineer/constructor companies throughout the world – since 1922.

At Flash conditions

``Σ Yi - Σ Xi = 0``

Above equation can be solved by iteration using Newton Raphson method. Function F(V) is defined as:

````F(V) = Σ Yi  - Σ Xi`
`F(V) = Σ [Zi (Ki - 1)/( V.(Ki - 1) + 1)]````

Derivative of F(V) is calculated as:

``F'(V) = Σ -[Zi(Ki - 1)² /( V.(Ki - 1) + 1)²]``

New estimate of vapor fraction is calculated as:

``V New = V - F(V)/F'(V)``

Function F(V) and F'(V) are calculated based on new vapor fraction and this process is repeated till there is negligible difference in between V and VNew. Vapor fraction thus obtained is then used to estimate vapor and liquid molar composition (Yi & Xi).

### Iteration for Ki

Vapor (Yi) and Liquid (Xi) mol fractions estimated above are used to generate values for Ki. Parameters for Peng Robinson EOS are calculated for each component i.

````κi = 0.37464 + 1.54226ω - 0.26992ω²`
`αi = [ 1 + κi (1 - (T/Tc)0.5)]²`
`ai = 0.45724 (RTc)²α / Pc`
`bi = 0.07780 RTc / Pc````

#### φiL Calculation

Mixture parameters are calculated.

````aij = [(ai.aj)0.5(1 - kij)] = aji`
`a = ΣiΣj aij.Xi.Xj`
`b = Σi bi.Xi`
`A = aP/(RT)²`
`B = bP/RT````

where, kij’s are Binary Interaction Parameter available from literature. Following cubic equation is solved to get ZL.

``Z³ + (B-1)Z² + (A-3B² -2B)Z + (B³+B²-AB) = 0``

Roots calculated are arranged in descending order, highest root gives ZV and lowest root gives ZL.

Based on ZL, liquid fugacity φiL is calculated for each component.

Featured Resources :

Serves chemical engineering professionals in the chemical process industry including manufacturing, engineering, government, academia, financial institutions and others allied to the field serving the global chemical process industry.

#### φiV Calculation

Mixture parameters are calculated.

````a = ΣiΣj aij.Yi.Yj`
`b = Σi bi.Yi`
`A = aP/(RT)²`
`B = bP/RT````

Cubic equation is solved to get ZV.

``Z³ + (B-1)Z² + (A-3B² -2B)Z + (B³+B²-AB) = 0``

Based on ZV, vapor fugacity φiV is calculated for each component.

Ki is calculated as:

``Ki = φiL/φiV``

New values of Ki thus calculated are again used to estimate V and thereafter Xi & Yi. Iteration is repeated till there is no further change in Ki values. Typically, in 10 iterations change in Ki values become negligible.

All above calculations along with iterative procedure for flash calculation have been provided in below spreadsheet.

Spreadsheet for PT Flash calculation using PR EOS

## Estimating Binary Interaction Parameter by Regression

Wilson equation is commonly used to predict non ideality in binary mixture vapor liquid equilibrium. This article shows how to estimate binary interaction parameters used in wilson equation from experimental data by regression in excel spreadsheet.

Example
Determine binary interaction parameters used in wilson equation for a mixture of acetone and chloroform from T x-y experimental data available at 101.33 kPa.

Obtain pure component properties of acetone and chloroform from literature mainly vapor pressure data and liquid molar volume.

Featured Resources :

Hydrocarbon Processing the premier magazine providing job-help information to technical and management personnel in petroleum refining, gas processing, petrochemical/chemical and engineer/constructor companies throughout the world – since 1922.

Based on modified Raoult’s Law following relation is obtained.

```` yi.P = xi.γi.Pisat `
` yi   = xi.γi.Ki````

From above equation experimental γ1 is obtained for acetone.

```` γ1(Exp) = y1/ ( x1.K1 )`
` K1 = [ e(A - B/ (T + C))] / P````

An initial value of binary interaction parameter is assumed.

```` A12 = 200 cal/gmol`
` A21 = 200 cal/gmol````

Liquid phase γ1 is obtained from above interaction parameter and using wilson equation.

```` Y12 = V2/V1.e-A12/RT`
` Y21 = V1/V2.e-A21/RT`
` lnγ1 = -ln(x1 + (1-x1)*Y12) + (1-x1)[ Y12/(x1 + Y12.(1-x1)) - Y21/(1 - x1 + Y21.x1) ]````

Square of difference of γ1_experimental and γ1_calculated is obtained for all data points. An objective function is defined as summation of all these differences.

Click on Solver in Data Ribbon (Excel 2010) to open dialog box for Solver parameters and input data as shown below.

Minimize the objective function by changing values of A12, A21. Uncheck Make Unconstrained Variables Non-Negative, as these variables can take negative values. Click solve to start regression and new values of A12 and A21 are calculated.

```` A12 = 157.9 cal/gmol`
` A21 = -570.3 cal/gmol````

Above values are then used to calculate y1 values and results are plotted to check the deviation.

Example
Determine binary interaction parameters used in wilson equation for a mixture of benzene and acetonitrile from P x-y experimental data available at 318.15 °K.

Use above steps and change formula for Y12, Y21 as temperature is fixed. After doing regression binary interaction parameters are obtained and result is plotted as following.

## Heat Loss from Insulated Pipe

Heat loss/gain takes place from a pipe carrying hotter/ colder fluid than ambient temperature. Insulation reduces the heat loss to surroundings. Heat loss depends upon number of factors like insulation thickness, ambient temperature, wind speed etc. This article shows how to calculate heat loss from an insulated pipe and a bare pipe to surroundings.

Example
A 3″ Carbon steel pipe is carrying hot oil at 180°C and insulated with 50 mm thick calcium silicate insulation. Insulation is cladded with a sheet with surface emissivity of 0.9. Ambient temperature is 28°C and wind velocity is 3.5 m/s. Calculate surface temperature and heat loss from insulated and bare pipe.

Overall heat transfer coefficient of an insulated pipe is defined as following.

where, kPIPE, kINSULATION are thermal conductivity of pipe and insulation. hin is heat transfer coefficient for fluid flowing in pipe and hair is heat transfer coefficient due to air flowing outside the pipe. The first two terms of denominator in above equation are generally smaller compared to remaining terms and can be neglected. For this example first term due to pipe fluid is ignored.

Featured Resources :

Hydrocarbon Processing the premier magazine providing job-help information to technical and management personnel in petroleum refining, gas processing, petrochemical/chemical and engineer/constructor companies throughout the world – since 1922.

### Air Side Heat Transfer Coefficient, hAIR

Air side heat transfer is due to combined effect of convection and radiation. Assume a temperature at cladding surface t_surface and steel pipe surface t_interface. Calculate average air film temperature as following.

`` t_average = ( t_surface + t_ambient )/ 2``

Estimate thermodynamic properties of air like thermal conductivity (k), viscosity (μ), expansion coefficient (β = 1/t_average), air density (ρ), kinematic viscosity (ν), specific heat (Cp) and thermal diffusivity (α) at average air film temperature. These properties are available in literature in form of tables, these can be fitted into a polynomial form using excel’s LINEST function. Reynolds’s number (Re), Prandtl number (Pr) and Rayleigh number (Ra) are calculated based on above properties.

Heat transfer coefficient due to radiation is calculated using following relation.

`` h_radiation = σ ε (t_surface4 - t_ambient4)/ (t_surface - t_ambient)``

where σ is Stefan Boltzmann coefficient and ε is emissivity for cladded surface.

### h_convection

Convective heat transfer coefficient comprises of forced and free convection. Forced convection can be modeled based on correlation by Churchill and Bernstein.

`` h_forced = Nu.k_air / D3``

Free convection is calculated based on correlation by Churchill and Chu.

`` h_free = Nu.k_air / D3``

Combined heat transfer coefficient due to forced and free convection is calculated using following relation.

```` Nu_combined = ( Nu_forced 4 + Nu_free 4) 0.25`
` h_convection = Nu_combined.k_air / D3````

Air side heat transfer coefficient is calculated as following.

`` h_air = h_radiation + h_convection``

### Overall Heat Transfer Coefficient, U

Thermal conductivity for insulation material and pipe is available in literature and depends upon temperature. It can be fitted into a polynomial equation using LINEST function in excel. Heat transfer resistance due to pipe and insulation is calculated using following relation.

```` r_pipe = D3.ln(D2/D1) / 2.k_pipe`
` r_insulation = D3.ln(D3/D2) / 2.k_insulation````

Overall heat transfer coefficient is calculated as.

```` r_overall =  r_pipe + r_insulation + 1/h_air`
` U = 1/r_overall````

Heat flowing through insulation is estimated.

`` Q = (t_operating - t_ambient)/r_overall``

A revised estimate for interface and surface temperature is made.

```` t_interface = t_operating - Q.r_pipe`
` t_surface   = t_interface - Q.r_insulation````

Above steps are repeated with these new estimates till there is negligible difference in temperature.

Heat loss per unit length of pipe is estimated as following.

`` HeatLoss = πD3 Q``

### Bare Pipe

For heat loss from bare pipe all above steps are repeated with resistance due to insulation not considered.

```` r_pipe = D2.ln(D2/D1) / 2.k_pipe`
` r_overall =  r_pipe + 1/h_air````

For this example surface temperature and heat loss are as following.

Spreadsheet for Heat Loss from Insulated Pipe