Process Process Modelli Modelling, ng, Simula Simulatio tion n and Contro Controll for Chemic Chemical al Engi Engine neer erin ing. g. Solv Solved ed prob proble lems ms.. Chap Chapte ter r 4: Nume Numeri rica call Methods. This document contains my own solutions to the problems proposed at the end of each chapter of the book ”Process Modelling, Simulation and Control for Chemica Chemicall Engine Engineers ers”” Second Second Edition Edition,, by William William L. Luyben. Luyben. At such, such, I can’t guarantee guarantee that the proposed proposed solutions solutions are free from errors. errors. Think about them as a starting point for developing or as a means of checking your own soluti solutions ons.. Any Any com commen ments ts or correc correctio tions ns will be apprec appreciat iated. ed. Contac Contactt me at
[email protected] The computer programs developed for this chapter (Matlab) are available at https://www.dropbox.com/sh/38yoxzta1h0wwda/AAC0B0sx15wOj9QkvsM2WAVxa?dl=0
Problem 1 Write a subroutine for bubble point calculation that uses false position convergence. Solution
The example used will be the same presented earlier in the textbook, that is,find the bubble point temperature for a mixture of 50% benzene and 50% toluen toluenee under under a pressu pressure re of 760 mmHg. mmHg. Also the followin followingg values alues of vapor pressures are known: Compo Compone nent nt Benzene Toluene
◦
100 100 C 1360 550
◦
125 C 2600 1140
Assume Assume that the vapor pressure pressure dependence on temperatur temperaturee for every component can be modeled using the following equation: P js =
exp
Aj + Bj T
The false false positio position n alg algori orithm thm is carrie carried d out in two two stages stages.. A first first stage, stage, which is similar to the interval halving algorithm, during which the temperature estimate is updated by a fixed amount until it crosses the correct value, and a second stage, during which the temperature estimate is updated using the recursive algorithm: T n+1 = T = T n − f (T n )
(T n − T n 1 ) f ( f (T n ) − f ( f (T n 1 ) −
−
With T 0 and T 1 beign the first pair of values that are at opposite sides of the correct value. A computer computer program was developed to determine determine the bubble point temperature, ature, the number of iterations iterations required for different different starting points and initial initial step sizes are shown in Figure 1. 1
Figure 1: Results for different starting points and initial step sizes. For the examples analyzed, it is observed that for a given initial step size, the convergence time is smaller for the initial temperature guess that is closer to the correct value. Also for a given temperature guess, the convergence time is smaller for the initial temperature step that is smaller. So, although a greater initial step size crosses the correct solution in fewer steps (in other words the second part of the algorithm can be started earlier), the greater deviation from the correct value that is produced, results in an overall greater convergence time. Code(s) used: P1 false position.m
Problem 2 Compare convergence times, using internal halving, Newton-Raphson, and false position, for an ideal, four-component, vapor liquid equilibrium system. The pure component vapor pressures are (P js (psia)):
2
Component 1 2 3 4
◦
150 F 25 14.7 4 0.5
◦
200 F 200 60 14.7 5
Calculate the correct temperature and vapor compositions for a liquid at 75 psia with a composition x 1 = 0.1, x 2 = 0.54, x 3 = 0.30, and x 4 = 0.06. Solution
The results for interval halving algorithm are shown in Figure 2, for NewtonRaphson in Figure 3, and for false position in Figure 4. For the examples analyzed, it is observed that for every algorithm, the convergence time is smaller for the initial temperature guess that is closer to the correct value. The Newthon-Raphson algorithm converges faster than both interval halving and false position algorithms. For false position algorithm, as discussed in Problem 1, the convergence is faster for smaller initial step sizes, the opposite is observed for the interval halving algorithm, this characteristics results in smaller convergence times for the false position algorithm (compared to interval halving) at smaller initial step sizes, whereas the opposite is true at greater initial step sizes.
Figure 2: Results for interval halving algorithm. Code(s) used: P2 false position.m P2 interval halving.m P2 newton raphson.m 3
Figure 3: Results for Newton-Raphson algorithm.
Problem 3 The design of ejectors requires trial and error to find the ”motive” pressure P m that, with a fixed motive flow rate of gas W m , will suck a design flow rate W s of suction gas from a suction pressure P s and discharge against a higher pressure P d . The motive gas is at 300 F and has a molecular weight of 60. Its flow rate is 5000 lbm /h. The suction gas is at 400 F and 150 psia and has a molecular weight of 50. A suction flow of 7000 lbm /h must be ejected into a discharge pressure of 160 psia. Assume perfect gases and frictionless, reversible, adiabatic operation of the jet; i.e., the expansions and contractions into and out of the throat are isentropic. The ratio of C p to C v for all gasses is 1.2 and C p heat capacities are constant and equal to 0.6 BT U/lbm R. Find the motive pressure P m required and the areas in the throat on the motive suction sides, A m and A s . ◦
◦
◦
Solution
A sketch of the ejector is shown in Figure 5. The system of equations that describe the ejector are: three mass fluxes (considering A s = Ad ): 7000 = A s ρs vs 5000 = A m ρm vm 12000 = A s ρd vd
(1) (2) (3)
Ideal gas expressions for calculating the density: P m M m RT m P d M d ρd = RT d
ρm =
4
(4) (5)
Figure 4: Results for false position algorithm. Isoentropic compression expressions for calculating ”intermediate temperatures” for use at the energy balance:
T s = T s ∗
P d P s
P d P m
T m = T m ∗
γ −1 γ
(6) γ −1 γ
(7)
An energy balance: 12000T d = 5000T m + 7000T s ∗
∗
(8)
and finally a momentum balance: As P s + 7000vs + Am P m + 5000vm = As P d + 12000vd
(9)
This leaves us with 11 unknowns (As , vs , Am , ρm , vm , ρd , vd , P m , T d , T s , T m ) and 9 equations. A design procedure as the one suggested in Perry’s Chemical Engineer’s Handbook relies on determining an optimum discharge over ∗
∗
5
Figure 5: Ejector. motive area ratio (Ad /Am ) and suction flow over motive flow ratio ( ws /wm ), given a discharge over suction pressure ratio (P d /P s ) and suction over motive pressure ratio (P s /P m ) utilizing a design curve (Figure 6). Here two additional relations will be used to complete the system of equations: P s = 0.25 P m As = 100 Am
(10) (11)
Figure 6: Design curves for optimum single-stage ejectors. The system of equations is solved utilizing the succesive substitution method, first the variables that can be solved directly: •
P m from (10).
•
T s from (6). ∗
6
•
T m from (7).
•
ρm from (4).
•
T d from (8).
•
ρd from (5).
∗
second, an iterative procedure is carried out: •
Guess A s .
•
Calculate v s from (1).
•
Calculate A m from (11).
•
Calculate v m from (2).
•
Calculate v d from (3).
•
Calculate A s from (9) (Ascalc ).
•
Reguess A s : As = A s + β (Ascalc
−
As ).
The results of the iteration results are shown in Figure 7 (The starting value for As and β where found by trial and error to assure convergence). The iteration procedure stops when A scalc and A s differ by less than 1%. The final results are shown in Figure 8. Code(s) used: P3 ejector design.m
Problem 4 Find the optimum liquid concentration of the propane-isobutane mixture in an autorefrigerated alkylation reactor. The exothermic heat QR (106 BTU/h) of the alkylation reaction is removed by vaporization of the liquid in the reactor. The vapor is compressed, condensed and flashed back in the reactor through a pressure letdown valve. The reactor must operate at 50 F and the compressed vapors must be condensed at 110 F . Find the liquid mole fraction x of propane that minimizes the compressor horsepower requirements for a given QR . Assume the compressor adiabatic efficiency is 100 percent. ◦
◦
Solution
The compressor power requirements are given by: γRT in wcomp = F γ − 1
P out P in
γ −1 γ
1
−
Where F is the molar flow through the compressor (calculated dividing the reaction heat by the heat of vaporization of the mixture), and γ will be assumed
7
Figure 7: Iteration results (Problem 3). equal to 1.2. The pressure at the inlet is given by the liquid composition and temperature inside the reactor. s s P in = P C xC + P iC 3
3
4
∗
xiC
4
The outlet pressure corresponds to the bubble point pressure at the temperature of the condenser (because we expect a that temperature a total condensation to occur). s s P out = P C yC + P iC 3
3
4
∗
yiC
4
Where yj correspond to the gas composition at the reactor (equal to the liquid composition at the condenser). Data for vapor pressure dependency on temperature as according to (T[K], P[Pa]) 1 : ln(P ) = A +
B + Cln(T ) + DT E T
is given: 1
Perry’s Chemical Engineer’s Handbook. 8th Edition.
8
Figure 8: Final results (Problem 3). Component Propane Butane
A 59.078 66.343
B -3492.6 - 4363.2
C -6.0669 -7.046
D 1.0919 ∗ 10 9.4509 ∗ 10
5
−
6
−
E 2 2
Data for entalphy of vaporization dependency on temperature according to (T[K],∆H v [J/Kmol]): 2
∆H v = C 1 (1 − T r )C +C T r +C T r 2
3
4
with T r = T /T c , is given: Component Propane Butane
C 1 2.9209 ∗ 107 3.6238 ∗ 107
C 2 0.78237 0.8337
C 3 -0.77319 -0.82274
C 4 0.39246 0.39613
T c 369.83 425.12
The different quantities involved in the calculation, and the power requirement are plotted against the liquid composition of propane in the reactor ( a base of Q=1055 J/h was used): It can be observed that the minimum work requirement corresponds to x 1 = 0, meaning that the reactor is filled only with butane. Although the pressure ratio is greater using only butane, as compared to using only propane, the greater molar flow due to the smaller heat of vaporization of propane increases the work done by the compressor. Code(s) used: P4 optimization.m
9
Figure 9: Inlet pressure at compressor.
Figure 10: Vapor mole fraction of propane at reactor.
10
Figure 11: Outlet pressure at compressor.
Figure 12: Flow through the compressor.
11
Figure 13: Power consumption at compressor.
12
