Microsoft Word - PAPER-V4-N2-2 NEW.doc IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 13 MODELING OF A REACTIVE DISTILLATION COLUMN: METHYL TERTIARY BUTYL ETHER (MTBE) SIMULATION STUDIES MUHAMAD NAZRI MURAT, ABDUL RAHMAN MOHAMED AND SUBHASH BHATIA∗ School of Chemical Engineering, Engineering Campus, Universiti Sains Malaysia, Seri Ampangan, 14300 Nibong Tebal, Penang, Malaysia Email:chbhatia@eng.usm.my Abstract: A process simulation stage-wise reactive distillation column model formulated from equilibrium stage theory was developed. The algorithm for solving mathematical model represented by sets of differential-algebraic equations was based on relaxation method. Numerical integration scheme based on backward differentiation formula was selected for solving the stiffness of differential-algebraic equations. Simulations were performed on a personal computer (PC Pentium processor) through a developed computer program using FORTRAN90 programming language. The proposed model was validated by comparing the simulated results with the published simulation results and with the pilot plant data from the literature. The model was capable of predicting high isobutene conversion for heterogeneous system, as desirable in industrial MTBE production process. The comparisons on temperature profiles, liquid composition profile and operating conditions of reactive distillation column also showed promising results. Therefore the proposed model can be used as a tool for the development and simulation of reactive distillation column. Key words: Modeling, simulation, reactive distillation, relaxation method, equilibrium stage, heterogeneous, MTBE 1. INTRODUCTION Reactive distillation process has been given special attention in the past two decades because of its potential for process intensification for certain types of chemical reactions [1]. The most important benefit of reactive distillation technology is a reduction in capital investment, because two process steps can be carried out in the same device. Such integration leads to lower costs in pumps, piping and instrumentation. For exothermic reaction, the reaction heat can be used for vaporisation of liquid. This leads to savings of energy costs by the reduction of reboiler duties. Reactive distillation process is also advantageous when the reactor product is a mixture of species that can form several ∗ To whom the correspondence should be addressed IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 14 azeotropes with each other. Reactive distillation conditions can allow the azeotropes to be “reacted away” through reaction [2, 3]. Although the advantageous of reactive distillation process was known since 1920, even until 1980, the technology was only utilised for homogeneous esterification reaction and was under utilised in other areas [4]. Developing reactive distillation column is a challenging task because of the complexities in column design, process synthesis and operability of reactive distillation processes resulting from the interaction of reaction and distillation [5]. With the aid of computer simulation, it was only recently been shown that difficult or unsolvable reaction-separation modeling problems are indeed solvable [6]. In more recent development, many reactive distillation processes such as esterification, hydration, alkylation have been utilised for reactions that rely on solid heterogeneous catalyst, often termed as catalytic distillation process. In the last two decades, research on heterogeneously catalysed reactive distillation has been extended to various industrial processes, including a well-known process for the production of octane booster methyl tert-butyl ether (MTBE), methyl acetate process and Nylon 6,6 process [5, 7]. Despite the recent advances by many researchers, it can be claimed that there is still no generally accepted design methods for reactive distillation. The main reason is the difficulty associated with reactive distillation processes and the parametric uncertainty of the model parameters [6]. With proper understanding of the thermodynamic and kinetic theories, several approaches and techniques have been proposed and claimed to be able to handle these difficult problems. Most of reactive distillation mathematical models are originally derived from conventional distillation calculations, which are based on equilibrium stage model. The liquid and vapour phase are assumed to be in phase equilibrium, and the stage is associated with the MESH equations (Material balance, Equilibrium relationship, Summation (constraint) equation, Heat (energy) balance). Relaxation method can be used as a solving tool for reactive distillation model with the MESH equations written in unsteady state form and thus utilising numerical integration to find the steady-state solution. This method is closely related to dynamic model; however, it is not often used because it was previously considered to be too demanding of computer time [2]. Today, with the high performance personal computer (PC) available in the market at a very reasonable price, the issue of too much of processing time related to PC requirement for relaxation method is not a critical issue. One of the numerical integration schemes that have been implemented in reactive distillation is the linear multistep method based on Backward Differentiation Formula [8]. As being known to author, so far, only the solution algorithm of Chen et al. [5] used this integration method coupled with relaxation method for solving homogeneous reactive distillation problem. The following are the advantages of this method [5]: • It is robust for arbitrary initialisation, e.g. convergence is almost always guaranteed for steady state solution. • Equilibrium stage models often have multiple steady state solutions, and the method greatly help the iteration to converge towards stable state, not to unstable state. IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 15 • Model is very stiff especially when the system approaches equilibrium reactive limit and the method is chosen to maximise accuracy and stability in this limit. This paper demonstrates the applicability and reliability of relaxation method that can be extended to heterogeneous system of reactive distillation process. The numerical method scheme is based on Backward Differentiation Formula linear multistep method for stiff system of equations [7]. The MTBE system is chosen for the verification of the proposed mathematical model. 2. MODELING The model assumes adiabatic operation of the column and the reaction takes place in the liquid phase. Each reactive stage is treated as a perfectly mixed continuous stirred tank reactor (CSTR). The vapour and liquid leaving any stage are in phase equilibrium with negligible heat of mixing of liquid and vapour mixtures. The vapour hold-up is assumed to be negligible. The hydrodynamic effects are neglected in the present modeling work in order to simplify the modelling complexities. With reference to the column configuration given in Figure 1, the model equations are presented in this section. They include mass and energy balances, vapour-liquid equilibrium (VLE) model and summation (constraint) equation. The kinetic model is shown separately although it can be considered as part of the material balance. The overall material balance for equilibrium stage j: Fj + Lj–1 + Vj+1 + δjRj = Lj + Vj (mol/s) (1) Unsteady state component i material balance: Hj t x ij d d , = (zj,iFj + xj–1,iLj–1 + yj+1,iVj+1) – (xj,iLj + yj,iVj) + δj ∑ = R r rjir rv 1 ,, )( (mol/s) (2) where j = 1, 2… N-1; i = 1, 2... c; r = 1, 2,..., R j is the stage number, i is the component number, r is the specific reaction number. The following definitions are useful in understanding and clarifying the derived material balance: • Liquid hold-up on stage j, Hj is defined as the molar quantity (or volumetric quantity) of liquid mixture that remains or held at certain level on stage j • Ri is the total numbers of moles generated or disappear through reaction on stage j. • ri,r is the rate of reaction r on stage j (mol/s). Parameter δj (0 or 1) refers to reaction occurrence on stage j. When reaction occurs on stage j, δj is set to unity but when there is no reaction; δj is set to zero value. IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 16 • Symbols zj,i, xj,i, yj,i are component i mole fractions of feed flow Fj (mol), liquid flow Lj (mol), and vapour flow Vj (mol) on stage j respectively. • The vr,i term is the stoichiometric coefficient of component i for reaction r. Stage N-1 Stage j Hj, Rj Vj Lj-1 Vj+1 Lj Stage j +1 Stage j -1 Stage 1 Lj+1 VN L1 Vj-1 V1 Stage 0 D = LD + VD L0 LD VD (V0) V1 VN LN-1 VN LN-1 Stage N B = LN F1 Fj-1 Fj Fj+1 FN-1 Stage R1 Reactive/Catalytic Distillation Zone Stage R2 Fig. 1: Schematic Representation of a Reactive Distillation Column IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 17 Energy Balance: ( ) ( ) ( )∑∑∑ = +++ = −−− = −+−+− c i L ij V ijijj c i L ij L ijijj c i L ij f ijijj hhyVhhxLhhzF 1 ,,1,11 1 ,,1,11 1 ,,, ( ) ( )∑∑ == −=∆− c i L ij V ijijj R r rj R rjjj hhxVrHW 1 ,,, 1 ,,δ (J/s), (3) where hL is the partial molar enthalpy of liquid (J/mol), hV is the partial molar enthalpy of vapour (J/mol), ∆HR is the heat of reaction (J/mol), and W is the weight of catalyst (kg). Phase Equilibrium: x P P y φ γ 0 = (4) The vapour phase is assumed to be ideal so that the entire fugacity coefficients φ for the system are equated to unity. The activity coefficients γ that characterise liquid phase non- ideality are calculated from the UNIQUAC method. The saturated vapour pressure Po is calculated from the Antoine equation and P is the total pressure. Summation (constraint) equations: 0.1 1 , =∑ = c i ijx (Liquid phase) (5) 0.1 1 , =∑ = c i ijy (Vapour phase) (6) Kinetic Model: The principal reaction between methanol(MeOH) and isobutene (IB) for MTBE synthesis is represented as: CH3OH + (CH3)2C=CH2 ⇔ (CH3)3COCH3 (7) MeOH IB MTBE The activity based rate model is:[9] IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 18         −= 2 MeOH MTBE MeOH IB aK a a a Wqkrate eq f (mol/s) (8) where W is the weight of catalyst (kg) and q is the amount of acid groups on the resin per unit mass (4.9 equiv/kg). The forward rate constant, kf is given by [9, 11]       −×= T k f 11110 exp1067.3 12 mol/(s equiv) (9) The rate constant units are based on equiv of hydrogen ions present in the catalyst rather than weight of the catalyst. Since hydrogen ions present in the catalyst are responsible for the catalytic reaction. The equilibrium constant Keq is given by [9, 11] Keq = 284 exp[f(T)] (10) and ( ) ( )+−+−+      +      −= 20 2 403 0 2 0 1 ln 11 )( TTATTA T T A TT ATf ( ) ( )40463035 TTATTA −+− (11) where T0 = 298.15 K, A1 = -1.49277×103 K, A2 = -77.4002, A3 = 0.507563 K-1, A4 = - 9.12739×10-4 K-2, A5 = 1.10649×10-6 K-3, A6 = -6.27996×10-10 K-4 [9, 11]. Figure 2 shows the sequence of calculations and algorithm for solving reactive distillation simulation model. The following are the explanation of the steps involved [5, 10]: 1. Specify all the specifications for a reactive distillation column to operate, and read all the constants and parameters related to the model equations 2. Initialise liquid compositions xj,i for the column, or use the last converged simulation. 3. Calculate vapour compositions yj,i and temperature Tj on stage j through bubble point flash calculation. Calculate liquid enthalpy and vapour enthalpy on stage j.Calculate the liquid and vapour rates throughout the column. Calculate all the mass balance derivatives, then calculate the residual from the average of the absolute value of all derivatives. Check if the residual is less than the specified error tolerance (e.g. 1.0×10 -4 ). If YES, converged steady state solution is obtained, then go to next step. If NO, predict next xj,i by calling the integration routine, then go to step 3. Output the composition and temperature profiles of the column. IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 19 Specify all column specifications and read all the required constants and parameters Calculate all mass balance derivatives, and then calculate the residual Initialise liquid composition throughout column Output the temperature and composition profiles of the column Calculate liquid and vapour flow rates throughout column Calculate liquid and vapour enthalpy throughout column Calculate vapour composition and temperature throughout column by bubble point calculation Check the specified tolerance for residual, and are new temperature and composition profiles near the previous values? YES NO Call integration routine Fig. 2: Algorithms for Solving Reactive Distillation Simulation Model. IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 20 3. SIMULATION AND RESULTS 3.1 Simulation Basis As a basis for simulations, it will be assumed that 197 mol/s (500 000 metric tons/year) MTBE is produced. Total mass of catalyst is 6400 kg for 8 reactive stages (800 kg catalyst per reactive tray), which is located in the middle of the column from stage 3 to stage 10. The total numbers of stages are 17 (numbered top-down) with total condenser as stage 0 and partial reboiler as stage 16. The column pressure is 11 bar (10.87 atm.) and the reflux ratio is set to 7. The column has two feed streams: a methanol feed and a mixed butenes (isobutene and 1-Butene) feed; a small stoichiometric excess of methanol is used. Table 1 and Figure 3 summarise simulation basis for heterogeneously catalysed MTBE synthesis [11]. Catalytic reaction zone Methanol Feed: Liquid at 320 K, 11 bar Stage 10 Flow rate: 206 mol/s Butenes Feed: Vapor at 350 K, 11 bars Stage 10 Flow rate: 549 mol/s Mole fraction Isobutene = 0.356 1-Butene = 0.644 Total condenser Partial reboiler Bottom Product Distillate Stage 0 Stage 16 Fig. 3: Column Configuration and Feed Specifications for MTBE Synthesis 3.2 Simulation Assumptions The MTBE system will be assumed to consist of 4 components: (1) methanol, (2) isobutene, (3) MTBE and (4) 1-Butene. The industrial MTBE process involves a minimum of 12 components with reactive components of methanol, isobutene and MTBE. The inert (C4 components) in MTBE system are quite similar to each other in physical and chemical properties, and thus assumed to be represented by 1-Butene. The by products are very small in quantity and can be regarded as negligible [1]. IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 21 Table 1 Feed Specifications, Column Configuration and Operating Conditions for MTBE Synthesis Quantity Units Methanol Butenes Feed mol/s 206 549 Phase Liquid Vapour Temperature K 320 350 Pressure bar 11 11 Feed stage 10 10 xMeOH 1.0 0.000 xIB 0.0 0.356 xMTBE 0.0 0.000 x1-Butene 0.0 0.644 Number of stages, N 17 Column Pressure, P bar 11 Catalyst Weight, W kg 6400 (total) Reflux ratio, r 7 Bottom flow, B mol/s 197 3.3 Simulation Results and Comparison Table 2 shows the results comparison between the proposed model and with published simulation results [11]. As it can be seen from Table 2, it is clear that the simulation results produced are quite close to the published simulation results [11]. This similarity verifies that the proposed model is capable of successfully simulating heterogeneously catalysed MTBE synthesis. It has been reported that a process utilising reactive distillation technology can boost the isobutene conversion up to 99% and higher [12], and possibly up to 99.99% isobutene conversion [13]. With agreement to the reported isobutene conversion, simulation results show an isobutene conversion of 99.5%, which is as desirable in real pilot plant reactive distillation column. Figure 4 and 5 show the temperature and liquid composition profiles obtained from simulation. Similar trend of temperature and composition profiles was reported[11] which is based on the UNIQUAC model and Redlich-Kwong equations for the non-ideality in liquid and vapour phases respectively. In the present research, the simulation results presented are based on the same non-ideal liquid phase using UNIQUAC model but with IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 22 the assumption of ideal vapour phase. The similarity in the simulation results indicate that the non-ideality in liquid phase is dominant, whereas the non-ideality of vapour phase does not significantly contribute to the differences in simulation results [10]. Table 2: Comparison of Simulation Results for Operating Conditions, Feed and Product Flow compositions Figure 6 shows the liquid phase activity coefficient profiles corresponding to the temperature and liquid composition profiles. It can be seen that quite large activity coefficients of methanol were obtained throughout the column. These large activity coefficients (as high as value of 11.6) express the highly non-ideal liquid behaviour of methanol in MTBE system. Other components, however, behave almost ideally as individual activity coefficient shows values of nearly unity throughout the column. The highly non-ideal behaviour of methanol in liquid phase for MTBE system has been reported by the experimental work [9]. The calculation based on UNIQUAC model has been reported to fit well the experimental data, by which the reaction kinetics has been formulated in activity basis [9]. Nijhuis et al. (1993) Proposed Model Quantity Units Top Bottom Top Bottom xMeOH 0.04 0.00 0.032 0.000 xIB 0.01 0.00 0.003 0.000 xMTBE 0.00 0.98 0.000 0.986 x1-Butene 0.95 0.02 0.965 0.014 Temperature, T K 348 424 348.8 423.8 Product Flow mol/s 366 197 363.7 197.0 Heat Duty, Q MW 49.7 35.1 49.7 39.3 IB conversion Mole % 98.5 99.5 MTBE Purity Mole % 98 98.6 Reflux Flow, L0 kmol/s 2.6 2.55 IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 23 TEMPERATURE PROFILE 340 350 360 370 380 390 400 410 420 430 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 Stage Number T em pe ra tu re (K ) Fig. 4: Temperature Profile throughout Column for MTBE Synthesis COMPOSITION PROFILE 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 Stage Number L iq ui d M ol e Fr ac tio n MeOH IB Fig. 5: Composition Profile throughout Column for MTBE Synthesis IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 24 ACTIVITY COEFFICIENT PROFILE 0 1 2 3 4 5 6 7 8 9 10 11 12 13 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 Stage Number A ct iv ity C oe ff ic ie nt γ MeOH IB MTBE 1-Butene Figure 6: Activity Coefficient Profile throughout Column for MTBE synthesis Figure 7 shows the reaction rate profiles calculated from activity basis rate. In continuous reactive distillation operation, there are many factors influencing the reaction conversion. The main parameters influencing reaction rate namely the temperature, liquid composition and the activity coefficients are also subjected to changes due to distillation effect in the column. For a non-ideal system, the reaction rate kinetic is always expressed in activity basis for best fitting of the experimental data. Therefore, it is very important to use reliable activity coefficient model since it is highly desirable to have correct activity coefficient values for the use in vapour-liquid equilibrium relationship as well as in activity based reaction rate kinetics. REACTION RATE PROFILE 0 10 20 30 40 50 60 70 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 Stage Number R at e (m ol s -1 ) Figure 7: Reaction Rate Profile throughout Column for MTBE Synthesis IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 25 A proposed model for reactive distillation system cannot be considered a practical model until it is verified by the experimental or industrial operation data. The comparison with other researchers’ simulation works is only a base platform for model verification. In literature, there are very limited industrial data for reactive distillation system especially for MTBE system. Recent research and simulation work [14] did include some of pilot plant data for MTBE system. All the data given are in weight (mass) basis, thus conversion to molar basis has to be made using relative molecular weight Mr. Table 3 shows the flows and liquid compositions of distillate and bottom product given in weight basis and corresponding calculated molar basis. Table 4 shows the comparison of MTBE plant data [14] and simulation results from proposed model. From Table 4, the overall comparison of reported MTBE plant data and simulation data show a reasonable agreement, but the simulation values for the C4 mole fraction in the product stream are not as close as to the reported plant data. This is because the typical industrial process for the MTBE production involves at least 12 components.[14] Besides the reactive components, which are methanol, isobutene (isobutylene) and MTBE, the inert in C4 mixture and the by-products involve a minimum of 5 other components namely trans-2-butene, cis-2-butene, diisobutene, isobutane, 1-butane, tert butyl alcohol (TBA) and propylene. Table 3: Flowrate and Liquid Composition of Distillate and Bottom Product in MTBE Plant [14] Given Data Calculated Component xi (Wt. %) Flowrate (kg/h) Rel. Mol. Wt., Mr Flow rate (mol/s) xi (mole frac.) MeOH 4.61 253.27 32.04 2.1957 0.0780 IB 0.24 13.19 56.11 0.0653 0.0023 MTBE 0 0 88.15 0.0000 0.0000 C4 95.15 5227.54 56.11 25.8794 0.9197 Distillate 5494 28.1404 MeOH 0.02 1.0857 32.04 0.0094 0.0005 IB 0 0 56.11 0.0000 0.0000 MTBE 98.64 5354.60 88.15 16.8734 0.9786 C4 1.34 72.74 56.11 0.3601 0.0209 Bottom Pdt. 5428.43 17.2429 IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 26 Table 4: MTBE Plant Data [14] and Simulation Results Comparison for MTBE Synthesis Distillate Bottom Product Composition Plant Model Plant Model MeOH 0.0780 0.1128 0.0005 0.0112 IB 0.0023 0.0000 0.0000 0.0000 MTBE 0.0000 0.0309 0.9786 0.9888 C4 0.9197 0.8563 0.0209 0.0000 Operating Parameters Plant Model % error Distillate (mol/s) 28.14 26.06 7.39 Bottom flow (mol/s) 17.24 17.24 (specified) Top temp. (K) 336 343 2.08 Bottom temp. (K) 417 417 0.00 Reflux flow, L0 (mol/s) 14.13 13.03 7.78 IB conversion (%) 99.62 99.99 0.37 Small amounts of water is also present in the hydrocarbon (C4 mixture) and in methanol feedstock but the reported data from product stream in the pilot plant [14] did not show any trace of water. This is because water is consumed in pre-reactor before it is fed into the reactive distillation column. Although the inert in C4 mixture have very similar physical and chemical properties, the by-products formed have somewhat different properties. The main side reactions involved are the dimerisation of isobutene to diisobutene (2,4,4-trimethyl-1-pentene) and the hydration of isobutene to tert-butyl-alcohol (TBA). The existence of the by-products is the main reason why the simulation results show some differences in the amount of C4 mixture. However, one advantage of MTBE reaction is that the by-products formed are very small in quantity [11] so that many researchers have assumed negligible side- reactions in MTBE synthesis. The hydration of isobutene to TBA can be neglected since kinetic study shows that the reaction can only take place when excess of methanol is insufficient [14]. Another reason that may contribute to the differences in simulation results is that the MTBE plant [14] uses structured catalyst packing, which has specific reaction efficiency and hydrodynamic conditions. Modeling the catalyst reaction efficiency and hydrodynamic conditions of specific structure packing is a challenging task due to insufficient information and research on the structured packing itself [2, 14]. Therefore the catalyst reaction efficiency and the hydrodynamic conditions are usually ignored especially in equilibrium stage modeling of reactive distillation column. IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 27 In the present simulation of MTBE synthesis, 1-Butene (an inert component) has been assumed to be the representative of the entire inert components in C4 mixture including that of small quantity of by-products. This assumption is justified by the fact that the C4 mixture in the distillate stream in MTBE plant [14] contains the largest fraction of 1- Butene. All the inert in C4 mixture are light components, thus only small amount of inert will be recovered at the bottom of the column. 4. CONCLUSION A steady state process simulation reactive distillation model was developed from unsteady state material balance equations based on equilibrium stage model. The model was verified and utilised for the simulation of heterogeneously catalysed MTBE synthesis. The proposed mathematical model describing reactive distillation process consists of mixed sets of ordinary differential equations (ODE) and algebraic equations. The model equations were solved numerically using backward differentiation formula linear multistep method for stiff system of equations [8] based on relaxation method. The main advantage of the proposed solution strategies is that it is robust for arbitrary initialisation. Convergence towards steady state is almost always guaranteed though it was considered slow when approaching the solution. The results obtained from the simulation of the proposed model have been shown to be capable of giving results as good as Aspen Plus simulation results [11] and MTBE plant data [14]. The simulation of MTBE synthesis shows the capability of the proposed model to give good prediction and thus capable of giving the insights of the performance of actual reactive distillation column. The temperature profiles, liquid composition profile and operating conditions of reactive distillation column showed promising results when compared to the published simulation results and the real MTBE plant data. Therefore the proposed model can be used as a tool for the development and simulation of reactive distillation column. In addition, the optimisation of reactive distillation can also be carried out much easier and safer, with time efficient and reduced costs. ACKNOWLEDGEMENT The financial support by the Ministry of Science, Technology & Environment, Malaysian Govt. under long term IRPA grant (Project no. 03-02-05-0007) is gratefully acknowledged. REFERENCES [1] T. Pöpken, S. Steinigeweng and J. Gmehling, “Synthesis and hydrolysis of methyl acetate by reactive distillation using structured catalytic packings: Experimental and simulation,” Ind. and Eng. Chemistry Res., 40, 1566-1574, 2001. [2] R. Taylor and R. Krishna, “Review: Modelling reactive distillation,” Chem. Eng. Science, 55, 5183-5229, 2000. IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 28 [3] A. Tuchlenski, A. Beckmann, D. Reusch, R. Düssel, U. Wiedlich and R. Janowsky, “Reactive distillation – industrial applications, process design & scale-up,” Chem. Eng. Science, 56, 387-394, 2001. [4] M.G. Sneesby, M.O. Tade, and T.N. Smith, “Steady-states transitions in the reactive distillation of MTBE,” Computers and Chem. Eng., 22(7-8), 879-892, 1998. [5] F. Chen, R.S. Huss, M.F. Malone and M.F. Doherty, “Simulation of kinetics effect in reactive distillation,” Computers and Chem. Eng., 24, 2457-2472, 2000 [6] P.A. Pilavachi, M. Schenk, E. Perez-Cisneros and R. Gani, “Modeling and simulation of reactive distillation operations,” Ind. and Eng. Chemistry Res., 36, 3188-3197, 1997. [7] J. Wang, X. Ge, Z. Wang and Y. Jin, “Development of a novel catalytic distillation column,” Chem. Eng. and Processing, 41, 115-121, 2002. [8] C.W. Gear, “Numerical Initial Value Problems in Ordinary Differential Equations,” Prentice Hall Inc., 1971. [9] A. Rehfinger and U. Hoffmann, “Kinetics of methyl tertiary butyl ether liquid phase synthesis catalyzed by ion exchange resin-I. Intrinsic rate expression in liquid phase activities,” Chem. Eng. Science, 45(6), 1605-1617, 1990. [10] M.N. Murat, “Modeling and simulation of reactive distillation column for the production of methyl tertiary butyl ether (MTBE),” M.Sc. Thesis, Universiti Sains Malaysia, 2002. [11] S.A. Nijhuis, F.P.J.M. Kerkhof and A.N.S. Mak, “Multiple steady states during reactive distillation of methyl tert-butyl ether,” Ind. and Eng. Chemistry Res., 32, 2767-2774, 1993. [12] J.L. DeGarmo, V.N. Parulekar and V. Pinjala, “Consider reactive distillation,” Chem. Eng. Progress, 88(3), 43-50, 1992. [13] G.G. Podrebarac, F.T.T. Ng and G.L. Rempel, “More uses for catalytic distillation,” Chemtech, 27(5), 37-45, 1997. [14] J. Bao, B. Gao, X. Wu, M. Yoshimoto, and K. Nakao, “Simulation of industrial catalytic- distillation process for production of methyl tert-butyl ether by developing user’s model on Aspen Plus platform,” Chem. Eng. J., 90(3), 253-266, 2002. NOTATION a Activity - c Total number of components - Fj Feed flowrate of stage j mol/s hL Partial molar enthalpy of liquid J/mol hV Partial molar enthalpy of vapour J/mol ∆HR Heat of reaction J/mol Keq Thermodynamic reaction equilibrium constant - kf Forward rate constant mol/(s equiv) Lj Liquid flowrate from stage j mol/s IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 29 Mr Relative molecular weight g/mol N Total number of stages - P Pressure Pa Po Saturated vapour pressure Pa q Ion exchange capacity equiv/kg Rj Total rate of generation of moles on stage j mol rj Rate of reaction on stage j mol/s T Temperature K t Time s v Stoichiometric coefficient - Vj Vapour flow rate from stage j mol/s W Weight of catalyst kg x Liquid composition - y Vapour composition - z Feed composition - 2.1.1. Greek letters δj Parameter for reaction occurrence on stage j (0 or 1) - φ Fugacity coefficient - γ Activity coefficient - 2.1.2. Subscript i Component number j Stage number r Reaction number 2.1.3. Superscript f Feed L Liquid R Reaction(s) V Vapour IIUM Engineering Journal, Vol. 4, No. 2, 2003 M.N. Murat et al. 30 BIOGRAPHIES Muhamad Nazri Murat received his B.Sc. (Hons) degree in Chemical Engineering from the University of Manchester Institute of Science & Technology (UMIST) in 1999. He obtained his M.Sc. degree in 2002 from the School of Chemical Engineering, Universiti Sains Malaysia. His research interests are in modeling and simulation, design and control aspect of reaction distillation column. Dr. Abdul Rahman Mohamed obtained his Ph.D in Chemical Engineering from the University of New Hampshire. His research interest is mainly in the environmental issues with emphasis in air pollution control. He is currently the Dean and Associate Professor of the School of Chemical Engineering, Unversiti Sains Malaysia. Dr. Abdul Rahman has conducted and completed numerous funded research projects and has published more than 80 papers at various international and national conferences and journals. Prof. Subhash Bhatia joined the School of Chemical Engineering, Universiti Sains Malaysia in 1995. He was a full Professor at the Department of Chemical Engineering, Indian Institute of Technology, Kanpur (India). Professor Bhatia was a visiting faculty at the University of Queensland, Australia from 1988 – 1989 and 1994 – 1995. His research interests are zeolite catalysis, chemical reaction engineering and environmental catalysis. He has written a book on Zeolite Catalysis that was published by CRC Press, USA and has published more than 100 papers at national and international journals.