Ind. Eng. Chem. Res. 1993,32, 2147-2158
2147
Characterization of Multicomponent Diffusion Effects in MTBE Synthesis David A. Berg and Thomas J. Harris' Department of Chemical Engineering, Queen's University, Kingston, Ontario, Canada, K7L 3N6
The application of multicomponent diffusion theory t o characterize intracatalyst mass-transfer resistances in the liquid-phase synthesis of methyl tert-butyl ether using acidic ion exchange resins is discussed. By formulation of invariant solutions to the multicomponent conservation equations, the system of coupled differential equations is reduced to a single ordinary differential equation and solved using orthogonal collocation. The parameters characterizing the invariant solutions are shown t o yield significant insight into the nature of the diffusion processes and illustrate the importance of accounting for multicomponent diffusion effects. A comparison between model predictions and published experimental data indicates good agreement between model and experiment in the limit where diffusion of methanol is rate limiting.
Introduction With the phaseout of leaded gasoline additives and the increased focus on development of cleaner burning automotive fuels, there has been an increasing demand for methyl tert-butyl ether (MTBE) as a gasoline additive. Some estimates are that by 1995demand for MTBE could increase at a rate of over 25%/year (Ainsworth, 1991). In the United States, recent amendments to the Clean Air Act require specific areas of the country that have not met carbon monoxide emission standards to utilize reformulated gasoline having a minimum oxygen content of 2.7 % Relative to other oxygenates such as alcohols, MTBE has superior blending properties and is more widely accepted by refiners. Commercially, MTBE is typically produced via direct addition of methanol to isobutylene using sulfonated ionexchange resins as catalysts:
.
CH,OH
+ (CH,),C=CH,
F?
(CH,),COCH,
(1)
In practice, selectivities approaching 100% can be achieved with dimerization of isobutylene the only significant side reaction. Hutchings et al. (1992) have reviewed many of the attributes of MTBE synthesis that must be taken into consideration in the design of commercial processes including a summary of existing technology. An increasingly popular route is reactive distillation using a solid catalyst contained within a distillation column. In this manner both chemical reaction and fractionation of products can proceed simultaneously. While commercial reactive distillation processes utilize structured packing in which the catalyst is containedinside either fabric pouches (Smith, 1982) or a reticulated wire envelope (DeGarmo et al., 19921, developmental work utilizing ion-exchange resins manufactured in the form of raschig rings has produced promising results (Flato and Hoffmann, 1992). A comprehensive review of the literature pertaining to MTBE kinetics is given by Hutchings et al. (1992). The most complete results appear to be those of Rehfinger and Hoffmann (199Oa),who completed an exhaustive study of the effect of both catalyst design parameters and reactant concentrations. The observed data were found to correlate well with a Langmuir-Hinshelwood model in which the reaction of adsorbed species is rate limiting: * To whom correspondence should be addressed.
where ai is the liquid phase activity of component i. In characterizing the above kinetics, Rehfhger and Hoffmann observed that diffusional limitations in the macropores of the catalyst could have significant effects upon the observed rate of reaction. This effect was quantified by developing a pseudobinary diffusion model in which the intracatalyst concentration of all reactants except methanol were assumed constant (Rehfinger and Hoffmann, 1990b). This simplification is valid for processes in which isobutylene is present in large excess and the model was consequently able to predict many characteristics of the experimental data. For commercial processes in which the ratio of methanol to isobutylene is approximately stoichiometric, the above assumptions are invalid and it is necessary to include multicomponent diffusion effects. The treatment of molecular diffusion processes in multicomponent mixtures is typically based upon the generalized Maxwell-Stefan (GMS) equations. For mixtures of ideal gases, the GMS equations can be reformulated and used to compute an effective multicomponent diffusivity for each species in the mixture (Froment and Bischoff, 1979). In liquid mixtures, the treatment of molecular diffusion is considerably more complex due both to the effect of thermodynamic nonidealities and the variation of the liquid phase GMS diffusion coefficients with mixture composition (Krishna and Taylor, 1986). While many researchers utilize an effective mixture diffusivity in modeling liquid-phase diffusion and reaction problems, little consideration is generally given as to whether this approach is consistent with the GMS equations. Sundmacher and Hoffmann (1992) extended the previous work of Rehfinger and Hoffmann (1990b)to include multicomponent diffusion effects by utilizing the GMS equations to develop species conservation equations in terms of the chemical potential of each species. The resulting system of differential equations was reduced to an equivalent binary form using an eigenvector transformation to diagonalize the matrix of diffusion coefficients (see Stewart and Prober, 1964). By utilizing an approximate rate expression, the problem was then cast into the form of the Frank-Kamenetzki equation which can be solved analytically for a cylindrical catalyst or numerically for a spherical catalyst.
oaaa-~a~~19312~32-2i47$04.00/00 1993 American Chemical Society
2148 Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993
In the present research, we consider the general problem of accounting for multicomponent diffusion effects in liquid-filled heterogeneous catalysts and apply the results to MTBE synthesis. It will be shown that by representing molecular diffusion processes within the catalyst using the GMS equations that the resulting system of species conservation equations can be reduced to a single ordinary differential equation through the development of invariant solutions. This approach is preferable to an eigenvector transformation in that the solution is obtained directly in terms of the independent variables rather than the transformed quantities. Additionally the parameters describing the invariant solutions are also shown to yield significant insight into the nature of intracatalyst diffusion processes. Last, the application of the above methodology to characterize the experimental data of Rehfinger and Hoffmann(1990b) is discussed. Multicomponent Diffusion i n Liquid Systems Diffusion processes in multicomponent liquid systems can be described using either the generalized MaxwellStefan (GMS) equations or a generalization of Fick's law. For a n-component system, Fick's law may be represented in matrix form as J = -c,[D]VX
(3)
where J is an n - 1 column vector of component diffusive fluxes (relative to a molar average reference velocity) and [D] is an ( n- 1) X ( n- 1) matrix of multicomponent Fick diffusion coefficients. In practical applications, a significant limitation of the Fickian formulation is that the Fick diffusion coefficients are strongly composition dependent even for simple binary liquid mixtures (Krishna and Taylor, 1986). Because of this, most attempts to correlate or predict multicomponent Fick diffusion coefficients have been relatively unsuccessful. Instead, most methods for correlating multicomponent diffusivities are based upon the Maxwell-Stefan equations. In the absence of thermal and pressure gradients, the GMS equations may be expressed as Xi -vTPi
t(xgJj-xJJi)
= c(xiNj-xPi)
and [I'] is an ( n - 1) X ( n - 1) matrix of thermodynamic factors with elements
In comparing (6) to (31, it can be readily seen that the GMS equations can be expressed in a form similar to Fick's law with [D] = [Bl-l[I']. For mixtures of ideal gases, the treatment is simplified considerably as [I'l reduces to the identity matrix and [D] = [Bl-l. In contrast, for highly nonideal liquid mixtures, the cross-term elements of [I'l can be relatively large in magnitude and thermodynamic effects can be very important. To explicitly relate the n independent Ni to the n - 1 independent Ji , it is necessary to consider the nature of convective transport within the mixture. This problem is discussed in detail by Krishna and Taylor (1986) and solutions for several types of diffusion problems are presented. For mixtures in which one component (component n) can be considered nontransferring or stagnant, (5) can be further simplified to N = [@I J = -c,[@][B]-'[I']Vx = -c,[D]VX where [@]is an ( n - 1)
X
(10) ( n - 1) matrix with elements
pij = 6 ,
+ Xj/Xn
(11)
[Dl = [@l[BI-'[I'l (12) Alternatively, a somewhat simpler representation can also be obtained by reformulating the GMS equations (4) directly and incorporating the constraint N,, = 0 for nontransferring n: N = -c,[A]-'[I']Vx (13) where [AI is an ( n - 1) X ( n - 1) matrix with elements
(4)
iz i
=
RT ]=I c&j c&j where Bij is the multicomponent Maxwell-Stefan diffusion coefficient and Ni is the total molar flux of component i relative to a fixed coordinate frame. For each component, the molar flux Ni can be computed as the sum of the diffusive and convective fluxes: Ni = Ji
+
t
Nj
(5)
J=1
While Sundmacher and Hoffmann (1992) have characterized diffusion processes in MTBE catalysis using the GMS equations expressed in terms of gradients in chemical potential, it is generally more convenient to develop solutions in terms of mole fraction gradients. The GMS equations can be reduced to (Krishna and Taylor, 1986) J = -c,[B]-'[I']Vx (6) where [Bl is an ( n - 1) X ( n - 1) matrix with elements given by (7)
Aij = - x i J q j (15) It can be shown by inverting 181 that [AI-' = [@I [Bl-l so that both (10) and (13) are in fact numerically equivalent (Krishna and Taylor, 1986). While (13) is somewhat easier to evaluate, (10) is utilized in the following analysis as it is a more useful representation for modeling of external transport resistances (Berg, 1993). Further reference will be made to (10) in developing conservation equations to describe multicomponent reaction and diffusion processes in MTBE synthesis. Development of Conservation Equations In developing conservation equations to describe multicomponent diffusion and reaction in MTBE synthesis, it is desirable to simplify the analysis by considering the minimum number of components required to adequately characterize the mixture. To account for multicomponent diffusion effects, a minimum of four components is required to adequately characterize the mixture composition: (1) methanol, (2) isobutylene, (3) MTBE, and (4) a pseudo compound to represent nonreacting components. In commercial processes in which isobutylene is obtained from a mixed Cq raffinate stream from a catalytic cracking
Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993 2149 unit, the nonreacting components will consist primarily of mixed butanes, l-butene, and 2-butene (DeGarmo et al., 1992). In the present research l-butene has been utilized to represent nonreacting components as this is also the solvent utilized in the experimental studies performed by Rehfinger and Hoffmann (1990b). The application of multicomponent diffusion theory to describe transport processes in porous catalysts requires consideration of pore geometry effects. The acidic ionexchange resins used for catalysis of MTBE typically consist of a macroreticular structure in which large macropores are interdispersed throughout a catalytically active gel phase. Reactants diffuse into the interior of the catalyst primarily through the macropores and subsequently react in the gel phase. Assuming Knudson diffusion effects can be neglected due to the large size of the macropores, the equations describing species transport within the catalyst may be formulated as N = -c~(c/T)[D]VX (16) where e is the voidage and T is the catalyst tortuosity. In accounting for the convective contribution to the total molar flux, l-butene can be consiiered nonreacting and hence stagnant. In this manner [Dl is represented by a 3 X 3 matrix. Restricting the analysis to steady-state conditions, the species conservation equations describing diffusion and reaction within the catalyst pellet can be formulated as
V-(C~(C/~)[D]VX) = -ucLr, (17) where u is a column vector of the stoichiometric coefficients, CL is the mean concentration of acid groups within the catalyst, and rv is the rate of reaction. To formulate the boundary conditions for (171, it is necessary to consider how the geometry of the problem is affected by the reversibility of the reaction kinetics. For diffusionally controlled systems the geometry will be described by a shell-core model in which the diffusion of reactants occurs in a shell region surrounding a central core in which reaction equilibrium is achieved. For a nonzero core volume the appropriate boundary conditions are
x i = xi,
at r = R,
xi = xis, dxJdr = 0
at r = R,
(18)
where R, is the radius of the catalyst particle and R, is the radial position of the shell-core boundary. The mixture composition at the shell-core boundary, xieq, will be dictated by the constraint of reaction equilibrium. For systems in which the diffusional resistances are small enough such that reaction equilibrium is not obtained in the center of the catalyst, R, = 0 and the constraint on mixture composition does not apply. It is readily apparent from the structure of (17) that if ct[Dl is considered variable over the diffusion path, the complexity of the species conservation equations is increased considerably. This contrasts with the problem of accounting for variable [Dl in nonreacting transport problems for which well-established methods are available (see Krishna & Taylor, 1986). In the following analysis it will be assumed that ct[Dl is constant over the diffusion path and can be evaluated using the mixture composition at the catalyst surface. The implications of this assumption will be discussed further in evaluating the results of the model predictions. For systems such as reactive distillation processes where the concentration of reactants is typically very small, it will be shown that this is a reasonable approximation.
For constant ct[Dl it is possible to employ reduction strategies to reduce the dimensionality of the species conservation equations prior to considering numerical solution techniques. Sundmacher and Hoffmann (1992) utilized an eigenvector transformation (Stewart and Prober, 1964) to diagonalize [Dl and reduce the coupled system of species conservation equations (17) to an equivalent series of equations which are coupled only through the reaction rate expression. Such systems of equations can be easily reduced to a single differential equation using invariant solutions which can be readily derived for systems with a diagonal matrix of diffusion coefficients (e.g., see Aris, 1975). A consequence of this transformation is that the solution is obtained not in terms of the original variables but in terms of the transformed quantities. As discussed by Taylor (1982), these transformed variables have no physical basis and solutions must be transformed back into the original variables to facilitate meaningful interpretation of the results. An alternative reduction strategy which avoids the above problems involves the development of invariant solutions directly from the multicomponent diffusion relationships. For systems in which a single chemical reaction occurs, invariant solutions for multicomponent diffusion problems can be formulated from relationships of the form
- tipi= 0
(19) where ( i j is a constant which can be computed from the stoichiometric coefficients of the reaction. The formulation of reaction invariants for multicomponent diffusion problems including the necessary and sufficient conditions for their existence is discussed in detail by Stewart (1978). Invariant solutions to the present problem can be shown to satisfy these conditions and are developed in the Appendix. The results can be expressed as N j
(xm - xm,) = Qm(xM&H - xM&H,) ('MTBE
- xMTBE,) = QMTBE(XMeOH - 'MeOH,)
(20)
(21)
where Qm and CMTBE are constants defined in the Appendix. The above invariant solutions indicate that for constant ct[Dl, the relative changes in the mole fractions of the reacting species will be linearly related. In addition to reducing the dimensionality of the species conservation equations, the above relationships provide a means of computing the equilibrium composition of the mixture and specifyingthe appropriate boundary condition at the shell-core boundary. Substitution of the invariant solutions into the rate expression (2) yields a simple nonlinear equation which can be solved to yield the equilibrium composition of the mixture. A simple bisection method proved suitable for this purpose in the present research. Substitution of (20) and (21) into the system of ordinary differential equations (17) yields (Appendix)
where D e f can be interpreted as the effective multicomponent diffusivity for methanol in the system and is defined in the Appendix. For a diagonal matrix [Dl in which there are no multicomponent diffusion effects, D,f reduces to the expected value &. The invariant solutions can also be utilized to express the reaction rate expression (2) explicitly in terms Of x M ~ H and , in this manner the problem is reduced to a single ordinary differential equation. The numerical solution of this problem is considerably simpler than the solution of the original system of differential
2150 Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993 IB-1B
Table I. Infinite Dilution Binary Diffusion Coefficients (%to) at T = 333 K MeOH MeOH MeOH IB IB MTBE
IB MTBE 1B MTBE 1B 1B
17.19 9.604 14.74 7.033 10.80 9.370
3.157 2.739 3.157 10.92 12.59 7.033
equations (17) and is preferred to the eigenvector transformation employed by Sundmacher and Hoffmann (1992). Prediction of Multicomponent Diffusion Coefficients Estimation_ of the modified matrix of Fick diffusion coefficients [Dl requires computation of the matrices [Bl-l and [I"]. While the elements of [I"] can be computed directly using models to predict the liquid phase activity coefficients, prediction of the elements of [ B Frequires a method to compute the GMS diffusion coefficients. Most methods for predicting the composition dependency of the GMS diffusion coefficients rely on correlations of the form
qj= f ( X , B $ ,
~~~
MeOH-IB
MTBE-IB
Figure 1. Contour plot of constant 6m for the quaternary mixture MeOH-IB-MTBElB. XIB = 0.5, T 333 K. IB-IB
qio)
where Byorepresents the diffusion coefficient for i in an infinitely dilute mixture of j . Both Kosanovich and Cullinan (1976) and Bandrowski and Kubaczka (1982) developed complexrelationships for r>, requiring estimates of several mixture properties and parameters obtained from experimental data. Kooijman and Taylor (1991) developed several simpler models for &j derived by extending the Vignes relationship (Vignes, 1966) to multicomponent mixtures. In the absence of experimental data to test the validity of more rigorous models such as that due to Bandrowski and Kubaczka (1982), the use of
is recommended which is the approach used in the present research. Values for the infinite dilution diffusion coefficients for use in (23) were predicted using the method developed by Hayduk and Minhas (1982) with the required pure component viscosity data obtained from experimental correlations in Reid et al. (1987). Liquid molar volumes for both the pure components and the mixture were estimated using the Hankinson-Brobst-Thomson correlation (Thomson et al. (1982). Tabulated values for the infinite dilution diffusion coefficients at 333 K are summarized in Table I. The elements of the matrix [I'] were computed using liquid-phase activity coefficients obtained from the UNIQUAC model (Prausnitz et al., 1986) using numerical estimates of the composition derivatives. The UNIQUAC binary interaction coefficients for the system MeOH-IBMTBE were obtained from Rehfinger and Hoffmann (199Oa)who noted that these coefficients predict complete miscibility of the liquid phase over all composition ranges which is in agreement with experimental observations. This is an important consideration as [I'l can be related to the criteria for thermodynamic stability of multicomponent liquid mixtures (Krishna and Taylor, 1986). For 1-butene, it is assumed that 1-butene and isobutylene behave similarly so that the binary interaction coefficients
MTBE-IB
MeOH-IB
Figum2. Contourplotofconstant tZmEforthequaternarymixture MeOH-€B-MTBE-1B. X ~ B= 0.5, T 333 K. IB-IB
MTBE-IB
MeOH-IB
Figure 3. Contour plot of constant& (XloB m2/s)for the quaternary mixture MeOH-IB-MTBE-1B. XIB 0.5, T = 333 K.
reported for isobutylene can be utilized for 1-butene. This is in agreement with the conclusions of Colombo et al. (1983) who observed that binary mixtures of Cr hydrocarbons with methanol or MTBE behave similarly. The effects of mixture composition on the computed , D,f at a temperature of 333 K values for 6m, ~ M T B E and are summarized in Figures 1-3, respectively. To illustrate the resulta graphically it has been necessary to fix the composition of 1-butene in the mixture ( X ~ B= 0.5) so that the effect of the remaining three cornponenta can be illustrated using ternary plots. In this manner the vertices of the plots represent binary mixtures of 1B-MeOH, 1BIB, and 1B-MTBE. The orientation of the contour lines with respect to the
Ind. Eng. Chem. Res., Vol. 32, NO. 9,1993 2151 coordinate axis provide an indication of the relative effect of various composition changes on the system parameters. At low X M ~ H the , lines of constant Def are approximately parallel to the lines of constant XMeOH indicating that the effect of XMeOH is considerably more important than that of either X ~ or Xm. E The plots of &m and & m exhibit ~ ~ similar behaviour at low XWBE though the effect Of X appears to become more significant with increasing X M ~ E and XMeOH. From Figure 1 it is apparent that &m can assume either positive or negative values. This has some very important implications regarding the nature of the diffusion process within the catalyst as the convective transport of both isobutylene and methanol necessarily occurs in the same direction. In the absence of multicomponent diffusion effects it would be expected that &IB> 0 and that the convective transport of both species would occur in the same direction as their respective composition gradients. With strong multicomponent diffusion effects though, &IB I0 is possibledue to the coupling effect of multicomponent diffusion in which the flux of IB is related to the composition gradients of all species in the mixture. For &m = 0, it follows from (20) that X B is constant so that mass transport occurs even in the absence of a gradient in the composition of IB. For &m < 0, the composition gradients of IB and MeOH are opposite in sign so that for Def> 0, IB will diffuse in a direction against its composition gradient. This type of diffusion process has been termed reverse diffusion in the literature and is characteristic of the diffusional coupling that can exist in multicomponent systems. Similar multicomponent diffusion effects have been observed in other chemical systems (Krishna and Taylor, 1986). The magnitude of the predicted values of &m, &;MTBE, and Defwere observed to be highly sensitive to the method used to estimate the GMS diffusion coefficients (Berg, 1993). The application of other models given by Bandrowski and Kubaczka (1982) or Kooijman and Taylor (1991) yielded results which differed by as much as 40%. The replication of Figures 1-3 using these methods indicated that the general orientation of the contour lines remained similar though, and all methods indicated a sign change in &IB in the same general region of the composition space. Experimental determination of the true multicomponent GMS diffusion coefficients is clearly required before a more rigorous comparison of the above models can be conducted. Numerical Solution of the Conservation Equations The reduced form of the species conservation equations (22) is similar to that of a binary diffusion and reaction problem described by a single ordinary differential equation with split boundary conditions. The application of orthogonal collocation to solve such problems is well documented in the literature (e.g., Villadsen and Michelsen, 1978) and is consequently utilized here. Discretization of the conservation equation (22) using global collocation yields a set of nonlinear algebraic equations which can be solved in combination with the boundary conditions (18) to approximate the solution to (22). Due to the highly nonlinear nature of the reaction rate expression (21, solution of the discretized equations proved to be extremely difficult. Global convergence using Marquardt’s method (Dennis & Schnabel, 1983) to solve the nonlinear equations could only be obtained for problems with R, = 0 and insignificant diffusion resistances. For problems with increasing diffusion limitations
and R, # 0, convergence using Marquardt’s method could not always be obtained and it was necessary to employ a homotopy-continuation algorithm. Homotopy methods have proven to be very successful in solving highly nonlinear boundary value problems (Ascher et al., 1988) providing the solution to the differential equations can be E readily parametrized in terms of a key variable describing the nonlinearity of the problem. In the present research, solution strategies were successfully developed employing both the surface methanol concentration and the pellet radius as continuation parameters. The choice of surface methanol concentration as a continuation parameter proved useful in replicating the results of Rehfhger and Hoffmann (199Ob) while continuation using the pellet radius was required to develop plots of effectivenessfactor versus Thiele modulus. In both cases, arc length continuation was utilized to avoid the problem of continuation past turning points. As an example, in using the pellet radius as a continuation parameter, avery small radius was initially selected so that R, = 0 and convergence using Marquardt’s method is globally convergent. Using arc-length continuation in the direction of increasing R,, the solution to the equations as a function of R, is readily mapped out until the point at which the reaction equilibrium a t the center of the catalyst isobtained. At that point, the fundamental nature of the problem changes and it is necessary to include R, as an additional variable (together with the constraint on reaction equilibrium a t the shell-core boundary) to further continue the solution. In this manner the resulting core volume fraction can also be mapped out along the solution path. While the above continuation strategy was generally quite robust, some numerical difficulties were encountered in generating solutions to problems in which the diffusion of methanol was rate limiting and a large number of collocation points were used. These difficulties appeared to be due to the highly nonlinear nature of the rate expression (2) near the reaction equilibrium point and the increasing importance of roundoff error. Rehfhger and Hoffmann (1990b) encountered similar numerical problems using a generalized shooting routine to solve the species conservation equation and were able to effectively resolve the problem by modifying the constraint imposed by reaction equilibrium at the shell-core boundary to nMeOH
-- X X M d H ,
at r = R,
(24)
where x is a constant which was assigned a value of 2.02. Physically, this change is consistent with modifying the reaction kinetics to include a discontinuity at the shellcore boundary:
r, = 0
for r < R,
rv = r, for r > R, (25) Using a similar approach in the present research with x = 1.001-1.20 proved sufficient to eliminate most numerical difficulties and the effect of this change on the predicted results was generally insignificant (Berg, 1993). The effect of the number of collocation points on the accuracy of the predicted results was investigated for several problems in which diffusion was strongly rate limiting. In general, the number of collocation points required to achieve a given level of accuracy varied significantly depending upon the nonlinearity of the reaction rate r, over the diffusion path. While 6-8 collocation points were sufficient for many problems, as many as 16-20 points were required for highly nonlinear
2152 Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993 problems. This relatively large number of points is quite surprising considering the solution profiles in all cases were observed to be extremely well behaved. The above observations, combined with the results of Rehfinger and Hoffman (1990b), appear to suggest that conventional methods such as orthogonal collocation or nonlinear shooting methods are not well suited to problems with highly nonlinear reversible reaction kinetics. This is not an area that has been adequately explored in the literature and requires further research. One possible solution may be the use of alternative collocation schemes employing B-spline basis functions with adaptive grid selection strategies (Ascher et al., 1988).
10
Comparison with Experimental Data The experimental data reported by Rehfinger and Hoffmann (1990b) were used to validate the multicomponent diffusion model for MTBE synthesis. The catalyst tortuosity, 7,was treated as a variable parameter and used to adjust the model predictions to obtain agreement with the experimental data. All model predictions were obtained using 16 collocation points with x = 1.20. The was matrix of modified Fick diffusion coefficients, [DI, assumed constant and evaluated using conditions at the catalyst surface. The effect of this assumption on the model predictions will be considered later. In characterizing the model predictions, it is convenient to present the results in terms of the overall catalyst effectiveness factor:
By using orthogonal collocation to solve the species conservation equations and evaluating the local reaction rate at the internal collocation points, it was possible to compute the above integrand directly using Gaussian quadrature. A comparison between the optimized model predictions and the experimental data of Rehfinger and Hoffmann (1990b) for R, = 0.247 X 103m and R, = 0.419 X le3 m is given in Figures 4 and 5, respectively. Both plots were constructed using a value of 7 = 0.76. Increasing or decreasing T shifted the curve to the right or left respectively and had little effect on the maximum value of q observed. The optimized value 7 = 0.76 was determined then by forcing agreement between the model and the data in the region of low XM~OH,. In this region, the core volume fraction (VJ V) approaches unity and the system is primarily limited by diffusional resistances. For R, = 0.247 X 10-3 m, both the model predictions and the experimental data show a negligible effect of XIB, upon q at low XM~OH, though a more significant effect is observed for R, = 0.419 X m. The present results can be compared with those of Sundmacher and Hoffmann (1992), who determined an optimum value of 7 = 1.25. The differences between these results can probably be attributed to differences between the methods used to predict the GMS diffusion coefficients. In both cases though the magnitude of T is substantially smaller than values typically observed in heterogeneous catalysts. Such low values of 7 are usually associated with surface diffusion effects, which can reduce the observed value of 7 when not directly included in the transport equations (Smith, 1981). While surface diffusion in heterogeneous liquid systems remains poorly understood, it can be a contributing factor in some systems (Yoshida et al., 1991). For MTBE synthesis, it is possible that a similar phenomena is occurring within the gel phase of
Figure 4. Comparison of model predictions with the experimental data of Rehfhger and Hoffmann (199Ob). Effectiveneaa factor and core volume fraction versus X M . O ~ T = 333 K,R, = 0.247 X 10-9 m, (O-,XIB, = O . ~ , X W I B ~ =0.058;O-- - , X I B , = ~ . ~ ~ , X W=I0.025). B&
0.1
LL'
,
I I
0.75
-
0.50
-
0.25
-
0.001
0 01
0 1
1
=hfe0HS
Figure 5. Comparison of model predictions with the experimental data of Rehfinger and Hoffmann (199Ob). Effectivenesa factor and core volume fraction versus x w w T = 333 K,R, = 0.419 X 1O-a m, (0-,ZIB, = 0.3, x m q = 0.015;0- - -,XIB, = O.Gl,zmp, = 0.015).
the ion-exchange resin. Further experimentation including measurement of the true intraparticle multicomponent diffusion coefficientswould be required to verify this effect. With increasing XM~OK, the core volume fraction decreases and the system changes from being primarily diffusion controlled to one in which reaction kinetics become increasingly important. At intermediate X M : ~ ) H , , the agreement between the model predictions and the experimental data deteriorates with the model significantly underpredicting the value of q. Effectiveness factors greater than 1.0 are possible in this region because the rate of reaction increases with decreasing XMeOH,. At high XMeOH,, both the model predictions and the experimental data approach q = 1.0. The values of Def, 6m, and & M ? ~ Eutilized to construct the results in Figures 4 and 5 are themselves functions of the mixture composition a t the catalyst surface and hence vary with XM~OK. While it was shown previously that 6~
Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993 2153 I
2.01
;
I
.........
12
MTBE \,
1
/
1I
LS O o l
I
1
0.001 0.00001 1
0.0001
0.001
0
I
I
0.01
0.1
1
xMeOH
Figure 7. Effectof a onreaction rate using modified reaction kinetics. T 333 K,xm = 0.30,X ~ = 0E.01. Figure 6. Normalized radial profiles of species composition and reaction rate. T = 333 K, R, = 0.247 X 10-9 m, X M ~ K= 0.02, xm, 0.30, X ~ F , 0.058.
can assume either positive or negative values (Figure 11, within the region 0.001 < XM~OH, < 0.2 it was observed that 6~ > 0.0. At larger values of XM~OH,, negative values for 6 m were computed indicating very significant multicomponent diffusion effects. As the reaction is kinetically controlled at large X M ~ H , though, diffusion effects have no influence on 7 in this region. Typical radial profiles of the reaction rate and species composition normalized with respect to the values at the catalyst surface are given in Figure 6. To illustrate the highly nonlinear nature of the reaction rate expression, a value of x = 1.01 was used to construct these profiles. The gradients in species composjtion are approximately linear near the surface of the catalyst, decreasing to zero at the s h e l l - ~ o & b r yBoth ~ MeOH and MTBE show* significant radial variations in composition while the composition of IB, which is present in large excess, is approximately constant. The radial profile of the normalized reaction rate is typical of the extreme nonlinearity of this system. At its maximum, the rate of reaction within the catalyst is approximately9 times the rate at the catalyst surface. Despite the increased rigor of this model, the model predictions are not substantially improved over the pseudobinary model of Rehfinger and Hoffmann (1990b). This can be attributed largely to the nature of the experimental conditions. Because IB is present in large excess,there will be little radial variation in IB composition (as predicted by the present model) and the simplified model of Rehfiiger and Hoffmann will yield similar results. For processes in which IB and MeOY are present in stoichiometric proportions, multicomponent diffusion effects will be more significant, and the two models could be expected to yield substantially different results. A sensitivity analysis of the model predictions to changes in the GMS diffusion coefficients is given by Berg (1993). Changes which had the effect of reducing the magnitude yielded a slight improvement in the model of 6~ or 6-E predictions however the differences between the model predictions and the data were still substantial. The nature of the model mismatch though suggests a deficiency in the kinetic model as the disagreement occurs in a region in which the process is both kinetically and diffusionally controlled. The experimental data utilized by Rehfinger and Hoffmann (1990a) to validate the rate expression (2) is limited to X M ~ O H> 0.02. From the solution of the species conservation equations though the concentration of meth-
c
1 -
I
tI
0.25 0.00 0 001
\ )
I
0.1
0.01
1
=M~OH,
Figure 8. Comparisonof model predictions obtained usingmodified reaction kinetics with the experimental data of Rehfhger and Hoffmann (1990b). Effectiveness factor and core volume fraction versus XWH,, T = 333 K, R, = 0.247 X 10-9 m, ( 0 -, xnj, = 0.3, x ~ F , =, 0.058; 0 - - -,xnj, = 0.63, X ~ E =, 0.025).
an01within the center of the catalyst is predicted to be as low as X M ~ O H = 0.0001. As the reaction rate increases significantly with decreasing XM~OH, the reaction kinetics in this region can have a significant effect upon 7. Ancillotti et al. (1978) have hypothesized that the actual kinetic constants in the rate expression may increase with decreasing XM~OH due to the resulting increased acidity of the catalyst. To evaluate this possibility, an empirical modification of the rate expression (2) has been tested
where a is a constant. Without changing the conditions for reaction equilibrium, the above expression is the simplest means of characterizing the possible effect of variable reaction kinetics. The effect of nonzero a on the rate of reaction, rv, is characterized in Figure 7. The effect of a is only significant for X M ~ O H< 0.05 which is outside the range of experimental data investigated by Rehfinger and Hoffmann. For a = 0.15, a comparison between the model predictions and the experimental data of Rehfinger and Hoffmann (1990b) for R, = 0.247 X m is illustrated in Figure 8. The experimental data points in Figure 8 reflect the increased
2154 Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993 10
0.1
+ Figure 9. Influence of x a on effectiveness factor and core volume , 0.10,X ~ E =, 0.01. fraction for variable @. T = 333 K,X M ~ H =
reaction rate at the catalyst surface (from the modified kinetic model (27)) and hence are different from those in Figure 4. In this case forcing agreement between the model predictions and the experimental data a t low X M ~ H ,yielded an optimum value of T = 1.69. In comparison to the results obtained using (21, the model predictions obtained using the modified reaction kinetics (27) are significantly improved a t intermediate XM~OH,. While the use of larger values of a further improved the agreement between model and experiment in this region, this was at the expense of reducing the agreement between model and experiment at large X M ~ H . The above discussion highlights the importance of obtaining reliable kinetic data at low XM&H in order to properly resolve the observed differences between model predictions and experimental data. While experimental data obtained for XMeOH > 0.02 is adequately characterized by the rate expression (21, the modeling of diffusion and reaction processes in MTBE synthesis requires reaction kinetics in the range 0.0001 < XMeOH < 0.01. The results of the present research lend support to the hypothesis of Ancillotti et al. (1978)that there is a change in the reaction kinetics a t low XM~OH.
Generalized Model Predictions In characterizing the effect of diffusional resistances upon the overall catalyst effectiveness factor, it is useful to parametrize the results in terms of a generalized Thiele modulus, 4. Here we define 4 as
where s is a geometrical shape factor with s = 2 for a spherical catalyst. As discussed by Aris (19751, there is considerable flexibility in defining the diffusion term in (28) for multicomponent diffusion problems. In the present example, it follows from the reduced species conservation equations (22) that the effective diffusivity of MeOH is the most appropriate term. The effect of surface composition on the relationship between 7 and 4 is illustrated in Figure 9. For fixed surface composition, 4 is a function only of R, so that increasing q5 corresponds to increasing R,. The general shape of the plots in Figure 9 are typical of reacting systems exhibiting
Langmuir-Hinshelwood kinetics. For low 4, all systems approach a limiting value of q = 1.0 corresponding to a kinetically controlled state. At intermediate 4, steadystate multiplicity is observed with q > 1.0. This type of behavior is characteristic of negative order reaction kinetics where the reaction rate can increase with decreasing reactant concentration. For large 4, all systems are limited by diffusional resistances and the core volume fraction approaches unity. The significant effect of XIB upon the nature of the plots in Figure 9 can be attributed directly to the form of the rate expression (2). Since the reaction is first order in XD and negative order in XM~OH, the variation of the local rate of reaction over the catalyst diameter will be dependent upon which species is rate limiting. For systems in which IB is in large excess, the diffusion of MeOH will be rate limiting and XIB will be approximately constant over the catalyst cross-section. In this case the localrate of reaction, rv,can be greater than the rate at the catalyst surface, rv, (as in Figure 61, and q > 1.0 is possible. For systems in which MeOH is in large excess, the diffusion of IB will be can be considered approximately rate limiting. As constant in this case, rv will be less than rv, and 7 > 1.0 will not be possible. While Figure 9 illustrates the possibility of steady-state multiplicity for intermediate values of 4, typically not all of these states will be stable with respect to small disturbances and therefore will not be observable experimentally. With the exception of the work by Luss and Lee (1971)who restricted their analysis to a diagonal matrix of diffusion coefficients [DI,the general problem of characterizing the relative stability of solutions to multicomponent diffusion and reaction problems does not appear to have been considered. As discussed by Luss and Lee, the analysis of such problems is necessarily more complex than analogous binary diffusion problems as invariant solutions can only be developed for steady state solutions. This problem will be investigated further in future research. For large 4, all plots approach a limiting asymptotic solution which can be described by (Froment and Bischoff, 1979) 7 = l/$Ja (29) where
The above expression represents an exact solution to the conservation equation (22) for slab catalysts though geometrical effects limit its accuracy for cylindrical or spherical catalyst geometries. For large 4, X M ~ H ,will be constrained by the boundary condition at the shell-core boundary (24)and can be estimated using the methodology described previously. The integrand in (30) cannot be evaluated analytically due to the complex form of the rate expression (2). In the present research the orthogonal collocation techniques developed previously have been utilized and the integrand computed using Gaussian quadrature. A comparison of the predictions obtained using (30) with those obtained from numerical solution of the conservation equations indicated that (30) should be sufficiently accurate for processes with strong masstransfer limitations. For q < 0.1 the differences between the two methods were generally insignificant.
Effect of Variable [D] In the previous analysis all modeling results-have been obtained by assuming that the variation of cJD1 over the
Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993 2156 I
125i t
5
O
i
A
1.00
0.01 0.00
1
0.00
1
1
I
I
0 25
0.50
0.75
1.00
XYeo/XYeoH,
Figure 10. Influence of variable [D]on the relative variation of xm and %M&H within the catalyst. T = 333 K, X M ~ = K 0.05, XIB, = 0.05, XMTBP, = 0.001 (-, constant [blevaluated at catalyst surface; variable [b]).
- .,
150
125
'( H
v,
I
1
I
- :ool
075
0'50 0.25
I
0 00 0 00
0 25
0 50
0 75
100
xMeOH/xMeOH,
Figure 11. Influence of variable [blon the relative variation of xm and XM&H within the catalyst. T = 333 K, X M ~ H , = 0.20, XIB, = 0.20, x-b = 0.001 (-, constant [b]evaluated at catalyst surface; ., variable [b]).
.
diffusion path is negligible and that values for D,f, GIB, and ~ M T B Ecan be computed using the mixture composition at the catalyst surface. The error introduced in this assumption will be directly related to how the above parameters vary over the diffusion path. To characterke this effect, consider first the influence of variable ct [Dl upon the relative variation of X M ~ Hxm, , and X ~ within E the catalyst. Rearranging (A6) yields dX&dxM,oH = 6,
dx-ddx,,
=G
~
(31) E
where 6~ and €-E are now local mixture properties and will vary along the diffusion path. For a specified mixture composition at the catalyst surface, (31) can be characterized as an initial value problem and integrated in terms of X M ~ O H . In the present research an explicit VSVO (variable step, variable order) Adams-Bashforth method with local extrapolation (Lambert, 1991)was used to obtain the solution to (31). To illustrate the effect of reactant concentration upon the variability of ct[Dl over the diffusion path, results are presented in Figures 10 and 11 for two stoichiometric mixtures with X M ~ H ,= 0.05 and XM~OH, = 0.20 respectively. The figures illustrate the variation of XIB r_elativeto XM,OH obtained both by assuming constant cJD1 evaluated at the surface composition and also variable ct[Dl. In both
0.1
1
lo
$ Figure 12. Effectiveness factor versus Thiele modulus for [bl evaluated using both the averagemixture composition and the mixture composition at the catalyst surface T = 333 K, XIB, = 0.20, x m b = 0.01, (- 2j = xi,; * * 2i = (xi, + ~i)/2).
cases, the termination point of the curves at low XM,OH represents the equilibrium composition of the mixture (Le., the resulting mixture composition at the shell-core boundary). The variation of XMTBE relative to XMeOH is characteristically similar and hence is not shown. The results computed for variable ct[Dl show an approximate linear relationship between XIB and XM,OH for low reactant concentrations though the results for high reactant concentrations exhibit significant nonlinearity. In both cases the profiles calculated assuming constant ct[Dl represent straight lines tangent to the true solutions (variable ct[bl) at the catalyst surface. The assumption of constant ct[Dl appears to be reasonable then for relatively dilute systems in which the nonlinearity introduced through composition effects on cJD1 is negligible. To characterize the effect of variable cJd1 upon the computed catalyst effectiveness factor it would be necessary to solve the unreduced species conservation equations (17). While this is beyond the scope of this work, an estimate of this effect can be obtained by characterizing the effect of the mixture composition ( R i ) used to evaluate ct[Dl. In the following, results have been obtained using In the second case, Ri is the arithmetic average of the compositions at the shell-core boundary and the catalyst surface. As the composition at the shell-core boundary is itself a function of ct[bl, an iterative method is required to estimate ct[Dl. In the present research a simple successive substitution routine was employed to solve this problem and convergence was generally achieved within a few iterations. Because the property estimation routines involve substantial computations though, even this relatively simple change increased the required computing time by a factor of 3. A comparison between the model predictions obtained using the above methods is illustrated in Figure 12 for the two systems characterized previously in Figures 10 and 11. For low 4, the effect of ct[D] is insignificant as both systems are kinetically limited. With increasing 4, the effect of 32i is clearly dependent upon the concentration of reactants. For relatively dilute solutions, the results using either the surface composition or the arithmetic average composition in the catalyst are very similar. For even lower reactant concentrations (XM~OH, < 0.02), the differences between the results become essentially insignificant. With increasing reactant concentration though, the parameters D,f, GIB, and &-E exhibit greater variability over the diffusion path and the effect of using an average
2156 Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993 mixture composition to evaluate cJD1 is much more significant. It is interesting to note that the effect of this change is to lower the maximum value of q , an effect opposite to that required to improve the agreement between model predictions and experimental data. The implication of the above results to modeling of industrial systems will be highly dependent upon the concentration of reactants in the system. For systems C 0.02, the parameters Def, &m,and & ‘ ; M ~ E where X M ~ H will not vary significantly over the diffusion path and can be evaluated using conditions at the surface of the catalyst. With increasing concentration of reactants though, the assumption of constant parameters likely introduces significant error into the computations. With regard to the model predictions summarized in Figures 4 and 5, the results obtained at h ’ X M e O H , which were used to optimize the value of T should not be affected by variation of the parameters over the diffusion path. At intermediate X M ~ H , though there is the possibility that errors introduced through the assumption of constant parameters are more significant. Conclusions Most treatment of liquid-phase diffusion and reaction processes in heterogeneous catalysis assume that the diffusion process can be characterized using an effective diffusion coefficient without consideration of multicomponent diffusion effects. The preceding analysis has demonstrated though that the treatment of multicomponent diffusion in heterogeneous reactions does not necessarily lead to increased complexity. Through the development of invariant solutions the multicomponent conservation equations can be reduced to a single ordinary differential equation which can be solved using conventional techniques. This reduction process is preferable to an eigenvector transformation in that the solution can be obtained directly in terms of the original variables rather than the transformed quantities. Additionally, the parameters describing the invariant solutions also provide significant insight into the nature of the diffusion process. The application of multicomponent diffusion theory to the liquid-phase synthesis of MTBE illustrates the importance of accounting for multicomponent diffusion effects. Within some regions of the composition space &‘IB can assume either positive or negative values, a direct result of the coupling effects present in multicomponent diffusion. Application of the model to characterize the experimental data obtained by Rehfinger and Hoffmann (1990b) indicated good agreement between model and experiment in the limiting region in which intracatalyst diffusion is rate limiting. In the region of intermediate 9 in which both diffusion and kinetic effects are important, the disagreement between model and experiment increases, and it is hypothesized that this is due to an deficiency in the kinetic model (2) in describing reaction kinetics at low XM~OH. Experimental verification of reaction kinetics a t H measurement of liquid-phase GMS diffusion low X M ~ and coefficients would be required to further investigate this result. Acknowledgment The authors would like to thank F. Keil (Technical University of Hamburg-Harburg, Germany) for helpful discussion. The first author would like to gratefully acknowledge the financial support provided by DuPont Canada Inc. and the Natural Sciences and Engineering Research Council of Canada.
Nomenclature ai = liquid phase activity of component i [Bl = inverted matrix of diffusion coefficients, s/m2 CL = concentration of acid groups within ion exchange resin, equiv/m3 ct = liquid molar density, km0l/m3 [D] = matrix of Fick diffusion coefficients, m2/s [Dl = matrix of modified Fick diffusion coefficients (= [81CBI-Wl), m2/s Der = effective multicomponent diffusivity for MeOH in mixture, m2/s &, = Maxwell-Stefan diffusion coefficient for i-j pair, m2/s &jo = diffusion coefficientfori in an infinitelydilute mixture of j, m2/s 6 = coefficient in (A6) Ji = molar diffusion flux of component i relative to molar average reference velocity, kmol/m% k = reaction rate constant, kmol/equiv s K, = reaction equilibrium constant Ni = molar flux of component i relative to a fixed coordinate frame, kmol/m% r = radial distance, m rv = reaction rate, kmol/equiv s R = ideal gas constant, 8414.4 J/kmol K R, = radius of shell-core boundary, m R, = radius of catalyst particle, m s = geometrical factor T = temperature, K xi = mole fraction of component i Ri = mixture composition used to evaluate [b] Greek Letters a = coefficient in (27) 181 = matrix of convective transport parameters Ti = liquid phase activity coefficient of component i [TI = matrix of thermodynamic factors 6ij = Kronecker delta (6ij = 1 for i = j, = 0 for i # j ) t = catalyst void fraction q = catalyst effectiveness factor ~.ri= chemical potential of component i vi = stoichiometric coefficient of component i [ i j = coefficient in (19) 7 = catalyst tortuosity 9 = Thiele modulus 9, = effective Thiele modulus defining asymptotic solution x = coefficient in (24) Subscripts 1B = l-butene c = shell-core boundary eq = equilibrium value z, j , k = component number IB = isobutylene MeOH = methanol MTBE = methyl tert-butyl ether n = number of components in mixture s = catalyst surface Appendix. Development of Invariant Solutions Here we consider the specific problem of modeling molecular diffusion processesoccurringin a porous catalyst with a single chemical reaction of the form s1s,
+ SZS,
F?
sss,
with component Sq assumed nonreacting and hence stagnant. This is representative of a general class of problems for which the specific case of MTBE catalysis is included. The species conservation equations for the above reaction scheme can be formulated as
Ind. Eng. Chem. Res., Vol. 32, No. 9,1993 2157
V-N, = vicLrv (A2) where vi is the stoichiometric coefficient for component i and r, is the rate of reaction. The boundary conditions for the solution of (A2) are identical to those illustrated previously (18). Applying the results of Stewart (19781, it follows directly from (A2) and the zero flux boundary condition at the center of the catalyst that the following invariant solutions exist: N 1 (vI/v2)N2 N 1 = (~1/vs)N3 (A3) The above invariant solutions are valid throughout the solution domain and apply irrespective of the mechanism of species transport within the catalyst. If it is assumed that surface and pore diffusion effects can be neglected, transport will occur by molecular diffusion and the total molar flux of each component, N i , can be computed using the results obtained previously (16):
where the modified matrix of Fick diffusion coefficients, [bl, includes the effect of convective transport and is described by (12). Utilization of the invariant solutions (A3) then yields a system of linear equations describing the relationship between the composition gradients of the reacting species: 3
- c t ( € / r ) p i j V x j = (Vi/V1)Nl 11 '
Providing [blis nonsingular, application of Cramer's rule (Anderson, 1986) to express the solution to this system of equations in terms of the composition gradients, V x i , then yields the result V X =~ 6 2 V x 1
V X =~ G 3 V x 1
(A6)
where
6,= det
[
v1
v2
v3
v1
v2
v3
D13
D23
D33
D13
D23
D33
v2
v3
D22
D32
D23
D33
v1
v2
v3
DII
D21
D31
4 2 D22
]bet
[
D32
v1 4
2
D13
1
The above results apply to each point along the diffusion path and indicate that a linear relationship exists between the species composition gradients a t all points. If the modified matrix of Fick diffusion coefficients, [DI, can be evaluated at some average composition and assumed constant over the diffusion path, 6 2 and 6 3 will be constant and the above equations can be integrated relative to the surface conditions to obtain ( ~ -2 ~ 2 , )= 6 2 ( ~ 1 -xi,)
( ~ -3~ 3 , )
6 3 ( ~ 1 -XI.)
647)
The above invariant solutions then provide a direct relationship between the relative changes in composition of the reacting species and effectively eliminate the need to solve two out of the three species conservation equations. Using component 1 as a basis, the conservation equations can be further simplified to -c,(t/r)DefV2xl = vlcLrv (A81 where Def can be interpreted as the effective multicom-
ponent diffusivity for component 1 in the mixture and is defined as
1
De, = v,det[DI det
[
D12 ;13
1, I] D22
D32
(A91
To characterize the range of possible values for Def,it is useful to expand the numerator in (A9). From basic matrix algebra, it follows from (12) that det[Dl = detU31 det1Bl-' det[I'l
(A101
Since component Sq is assumed nonreacting, the elements of 181can be computed from (11). Evaluating the determinant of the matrix then yields det1bI = l/x,
(All)
It can also be shown (Krishna and Taylor, 1986) that det[Bl and det[I'l are both greater than zero for mixtures which are thermodynamicallystable with respect to p h q e splitting. It follows from matrix algebra then that det[Dl > 0 for thermodynamicallystable liquid mixtures and that invariant solutions can always be developed for (A5). The locus of mixture compositions for which d e t [ r l = 0 defiies the spinodal curve which separates the metastable and unstable regions. In this case, (A91 reduces to Der= 0 and no mass transport occurs within the catalyst. Effectively this implies a core volume fraction approachingunity with an overall catalyst effectiveness factor of zero. For situations when the denominator of (A9) becomes singular, the flux of component 1 and ita respective composition gradient become uncoupled and mass transport occurs even though Vxl = 0. Such problems cannot be solved using component 1 as a basis and must be reformulated using either component 2 or 3 as a basis. The above result can be compared to the more simplified case in which the reacting system can be characterized as a mixture of ideal gases. In this case, the effective diffusivity of componentj in the mixture can be obtained directly from the Maxwell-Stefan equations as (Froment and Bischoff, 1979)
The analogous result for liquid systems (A91 is considerably more complexthan the above and can be attributed to the effect of thermodynamicnonidealitieson the driving force for mass transport.
Literature Cited Ainsworth,S.J. BoomingMTBE DemandDraws IncreasingNumber of Producers. Chem. Eng. News. 1991,69, 13-16. Ancillotti, F.; Massi Mauri, M.; Pescarollo, E.; Romagnoni, L. Mechanisms in the Reaction Between Olefins and Alcohols Catalyzedby Ion Exchange Resins. J. Mol. Catal. 1978,4,37-48. Anderson, R. F. V. Introduction to Linear Algebra; Holt, Rinehart and Winston: Toronto, 1986; pp 349-352. Aris, R. The Mathematical Theory of Diffusion and Reaction in Permeable Catalysts. The Theory of the Steady State; Oxford University Press: London, 1975; Vol. I, pp 64-71. Ascher,U. M.; Mattheij,R.M. M.; Russell,R. D. Numerical Solution ofBoundary ValueProblems for OrdinaryDifferentialEquations; Prentice Hall: Englewood Cliffs, NJ, 1988; pp 244-266,490-500. Bandrowski, J.; Kubaczka, A. On the Prediction of Diffusivities in MulticomponentLiquid Systems. Chem. Eng. Sci. 1982,37,13091313.
2158 Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993
Berg, D. A. Transport Resistances in MTBE Reactive Distillation Processes. Ph.D. Dissertation, in progress. Queen’s University, Kingston, ON, 1993. Colombo, F.; Corl, L.; Dalloro, L.; Delogu, P. Equilibrium Constant for the Methyl tert-Butyl Ether Liquid-Phase Synthesis by Use of UNIFAC. Znd. Eng. Chem. Fundam. 1983,22,219-223. DeGarmo, J. L.;Parulekar, V. N.; Pinjala, V. Consider Reactive Distillation. Chem. Eng. Progr. 1992,88, 43-50. Dennis, J. E.; Schnabe1,R.B. Numerical Methods for Unconstrained Optimization and Nonlinear Equations; Prentice-Hall: Englewood Cliffs, NJ, 1983; pp 147-152. Flato, J.; Hoffmann, U. Development and Start-up of a Fixed Bed Reaction Column for Manufacturing Antiknock Enhancer MTBE. Chem. Eng. Technol. 1992,15,193-201. Froment, G. F.; Bischoff, K. B. Chemical Reactor Analysis and Design; Wiley: New York, 1979; pp 141-200. Hayduk, W.; Minhas, B. S. Correlations for Prediction of Molecular Diffusivities in Liquids. Can. J. Chem. Eng. 1982,60, 295-299. Hutchings, G. J.; Nicolaides, C. P.; Scurrell, M. S. Developments in the Production of Methyl tert-Butyl Ether. Catal. Today 1992, 15,23-49. Kooijman, H. A.; Taylor, R. Estimation of Diffusion Coefficients in Multicomponent Liquid Systems. Znd. Eng. Chem. Res. 1991,30, 1217-1222. Kosanovich,G. M.; Cullinan, H. T. Jr. A Study of Molecular Transport in Liquid Mixtures Based on the Concept of Ultimate Volume. Znd. Eng. Chem. Fundam. 1976,15,41-45. Krishna, R.; Taylor, R. Multicomponent Mass Transfer-Theory and Applications. In Handbook of Heat and Mass Transfer; Cheremisinoff, N. C., Ed.; Gulf Publishing Co.: Houston, TX, 1986; Vol. 2,Chapter 7. Lambert, J. D. Numerical Methods for Ordinary Differential Systems: The Initial Value Problem; John Wiley & Sons: Chichester, England, 1991;pp 80-87. Luss, D.; Lee, J. C. M. Stability of an Isothermal Catalytic Reaction with Complex Rate Expression. Chem. Eng. Sci. 1971,26,14331443. Prausnitz, J. M.; Lichtenthaler, R. N.; Gomes de Azevedo, E. Molecular Thermodynamics of Fluid-Phase Equilibria; Prentice-Hall: Englewood Cliffs, NJ, 1986;pp 262-265.
Rehfhger, A.;Hoffmann,U. KineticsofMethylTertiaryButylEther Liquid Phase Synthesis Catalyzed by Ion Exchange Resin. I. Intrinsic Rate Expression in Liquid Phase Activities. Chem. Eng. Sci. 1990a,45,1605-1617. Rehfhger, A.; Hoffmann, U. Kinetics of MethylTertiary Butyl Ether Liquid Phase Synthesis Catalyzed by Ion Exchange Resin. 11. MacroporeDiffusion of Methanol as Rate-ControllingStep. Chem. Eng. Sci. 1990b,45,1619-1626. Reid, R. C.; Prausnitz, J. M.; Poling, B. E. The Properties of Gases and Liquids; McGraw-Hi& New York, 1987;pp 441455. Smith, J. M. Chemical Engineering Kinetics; McGraw H i k New York, 1981;pp 462-473. Smith, L. A. Jr. Catalytic Distillation Process. U.S. Patent 4,336,407, June 1982. Stewart, W. E. Invariant Solutions for Steady Diffusionand Reaction in Permeable Catalysts. Chem. Eng. Sci. 1978,33,547-553. Stewart, W. E.; Prober, R. Matrix Calculation of Multicomponent Mass Transfer. Znd. Eng. Chem. Fundam. 1964,3,224-235. Sundmacher, K.; Hoffmann, U. Importance of Irreversible Thermodynamics for Liquid Phase Ion Exchange Catalysis: Experimental Verification for MTBE-Synthesis. Chem. Eng. Sci. 1992, 47,2733-2738. Taylor, R. Solution of the Linearized Equations of Multicomponent Mass Transfer. Znd. Eng. Chem. Fundam. 1982,21,407-413. Thomson, G. H.; Brobst, K. R.; Hankinson, R. W. An Improved Correlation for Densities of Compressed Liquids and Liquid Mixtures. AZChE J. 1982,28,671-676. Vignes, A. Diffusion in Binary Solutions. Znd. Eng. Chem. Fundam. 1966,5,189-199. Villadsen, J.; Micheleen, M. L. Solution of Differential Equation Models by Polynomial Approximation; Prentice-Hall: Englewood Cliffs, NJ, 1978. Yoshida, H.; Maekawa, M.; Nango, M. Parallel Transport by Surface and Pore Diffusion in a Porous Membrane. Chem. Eng. Sci. 1991, 46,429-438. Received for review February 1, 1993 Accepted May 14, 1993
