Modeling of a Three-Phase Continuously Operating Isothermal

ARTICLE pubs.acs.org/IECR

Modeling of a Three-Phase Continuously Operating Isothermal Packed-Bed Reactor: Kinetics, Mass-Transfer, and Dispersion Effects in the Hydrogenation of Citral Teuvo Kilpi€o,* P€aivi M€aki-Arvela, Mats R€onnholm, Victor Sifontes, Johan W€arna, and Tapio Salmi

Laboratory of Industrial Chemistry and Reaction Engineering, Process Chemistry Centre, Åbo Akademi, FI-20500 Turku/Åbo, Finland ABSTRACT: A continuously operating isothermal dynamic packed-bed reactor was modeled. The model included chemical reaction, gasliquid mass transfer, convection, axial dispersion, pore diffusion, and catalyst deactivation. The model was solved by using the method of lines. The model was applied on experimental data from citral hydrogenation over a supported nickel catalyst. The experiments had been carried out at 2565 °C and at 6.1 bar in a laboratory-scale trickle-bed reactor (d = 1 cm; L = 5 cm). The parameters in the model were the rate constants, pore diffusivity, coking rate constant, Peclet number, and gasliquid mass-transfer coefficient. A sensitivity study was performed to reveal how much a change in each of these changed the product concentration trend. The simulations revealed that the gasliquid mass-transfer coefficient and effective diffusivity should have been well below expected values to significantly reduce the productivity. The gasliquid mass transfer and pore diffusion were not rate-limiting; because hydrogen was used in excess, particles were small and the system was dilute. The citral concentration-dependent deactivation model based on site competition was able to describe the observed activity decline. Parameter estimation for the reaction rate and coking rate was carried out. A reasonable agreement with the experimental trends was obtained. An estimate of the Peclet number was obtained from step-response measurements with an inert tracer, which revealed that the reactor did not operate completely as a plug-flow unit. The model described here can be extended to be applicable for other hydrogenation and oxygenation reactions of other fine chemicals.

1. INTRODUCTION Chemical reaction engineering specialists are focusing their efforts on searching for various ways to replace petroleum-based feedstocks with renewable raw materials. Another major trend to be seen is to have nontoxic intermediates in the production processes whenever feasible. In this work, the production of citronellal from citral was studied. Citronellal is a nontoxic, renewable main component in citronella oil, giving it the distinctive lemon scent. It can be used in the synthesis of a multitude of organic compounds, the most well-known of them is menthol. Citronellal can be made from citral, another renewable compound in essential oils. Menthol is in the limelight again because the market situation is favorable and the main producers are increasing their production capacity remarkably. The aim of this work was to develop a model for a continuously operating catalytic three-phase reactor by using mathematical modeling and simulation, based on the component mass balances. Most models for three-phase reactors so far are designed for the production of bulk chemicals and fuel components. The lack of models for laboratory-scale trickle beds is evident; therefore, this work was carried out. The experimental study was carried out in a new packed-bed reactor system, designed in our laboratory (Figure 1). Experimental kinetic data were obtained using six reactors in parallel.1 This modeling study concentrates on the first reaction step in the straight synthesis of menthol, the hydrogenation of citral to citronellal. Experimental data from citral hydrogenation was used for model testing and fine-tuning. The phenomena that could be significant for this case were included (chemical reaction, gasliquid mass transfer, pore diffusion, axial dispersion, and r 2011 American Chemical Society

Figure 1. TBR schematically and in practice.

deactivation). Heat-generation and pressure-drop effects were checked also. Models of varying degree of complexity have been reported in numerous articles and textbooks.27 The specific features of the current work were the application for which the model was done, the range of phenomena included, the dynamic nature of the model, and the pore-diffusion model derivation and solution. The calculation method of pore diffusion8 was extended to Special Issue: CAMURE 8 and ISMR 7 Received: September 1, 2011 Accepted: October 28, 2011 Revised: October 26, 2011 Published: October 28, 2011 8858

dx.doi.org/10.1021/ie2019678 | Ind. Eng. Chem. Res. 2012, 51, 8858–8866

Industrial & Engineering Chemistry Research

ARTICLE

Table 1. Parameter Values in Experiments1 variable

values

units

residence time for the liquid

80 and 160

s

velocity of the liquid (interstitial)

0.032 and 0.064

cm/s

superficial velocity of hydrogen

2.1

cm/s

pressure

6.1

bar

temperature

25, 40, 45, 50, 55, 60, and 65

°C

catalyst particle size

6393 and 150250

μm

catalyst mass

0.1

g

silica mass citral concentration in the feed

1.2 0.02

g mol/L

pore volume of silica

0.83

cm3/g

reactor diameter

0.01

m

reactor length

0.05

m

LangmuirHinshelwood (LH) kinetics for a spherical catalyst particle. Analogous derivation can be applied to other kinds of kinetics and catalyst shapes.

2. MATERIALS AND METHODS 2.1. Model for the Citral Hydrogenation Reactor. 2.1.1. Main Hypothesis. The mathematical model has the following

features: (i) The system is described by a dynamic multiphase axialdispersion model. This gives the opportunity to use robust solvers for ordinary differential equations. The steadystate approach model could also have been used as an alternative, but the reactor dynamics and catalyst activity changes have to be described in an appropriate way. Furthermore, a stable model solution algorithm was desired. Therefore, a dynamic model was developed. (ii) The model was based on partial differential equations, which were solved using the method of lines. (iii) The kinetics of the reaction system, citral hydrogenation, was treated as a single main reaction following LH kinetics, with a consecutive side reaction taking place marginally.1 (iv) Catalyst deactivation was modeled as a coke formation reaction with LH kinetics. (v) The Peclet number for axial dispersion was obtained from a separate tracer test in a nonreactive system. An analytical solution for the residence time distribution (RTD) for such a case was available. The particle-diffusion model solved the concentration profile step by step, and then the effectiveness factor was obtained. (vi) Partial wetting was not included in the model because of the use of a fine-sand-bed feed distributor and because of the small reactor volume. Partial wetting is an important issue to be considered especially in process scale-up. Liquid tends to flow as rivulets, and some part of the liquid may stay stagnant in pockets between particles. Therefore, as the results of this study are interpreted, it has to be kept in mind that they were based on complete wetting. 2.2. Experimental Data. A number of experiments had been carried out in the multitubular laboratory-scale reactor system. The variable values in the experimental study are given in Table 1. The solvent was ethanol (Primalco). Its boiling point at ambient

pressure is 78 °C. At 6 bar, it is even higher. Because the vapor pressure depends exponentially on the temperature, the system could be treated as nonvolatile. 2.3. Model Equations. Model equations are collected in Table 2. The liquid-phase material balances are shown in eq 1. The accumulation, convection, axial-dispersion, mass-transfer, and reaction terms were included.7 For citral, all except the mass-transfer term were used. For dissolved hydrogen, also the gasliquid mass transfer was added. The activity and effectiveness corrections were applied to the reaction rates. The boundary conditions were eqs 2 and 3. The dynamic liquid holdup was dependent only on Rel and the geometry because of low flow rates9 (eq 4). The gas mass balances7 in a general form are given in eq 5. In the current case, the balance was simplified with the reasoning found in eqs 69. The gasliquid molar flux is given in eq 10. The gasliquid mass-transfer coefficient was obtained from a power law10 model (eq 11). It had been previously used successfully in hydrogenation of an organic compound.11 A collection of alternative expressions do also exist.12 The saturation concentration of hydrogn was obtained from Henry’s law (eq 12). The heat of reaction was checked using the group contribution method of Benson,13 and the pressure drop was checked with eqs 13 and 14. The kinetics of the main reaction was described with a LHtype expression14 (eq 15). The reaction order with respect to hydrogen corresponds to molecular adsorption. For the side reaction, eq 16 was used.14 Its extent was minor and diminished rapidly. The selectivity was high (98% citronellal and 2% citronellol). Catalyst deactivation was predominant, and it was modeled as coke formation according to eqs 1719. The concentration profiles inside the catalyst particles were obtained from the mass balance (eq 20) and by use of boundary conditions (eq 21). The shape factor of the catalyst was defined in eq 22. The effectiveness factor was calculated from the concentrations inside the particle with eq 23. The Peclet number was defined15 by eq 24. It was calculated from the RTD using the solution of the mass balance (eq 25) with boundary conditions (eq 26). 2.4. Numerical Techniques. Calculation of the reactor mass balances was initiated by having the first step so small that no reaction took place there. At the outlet, the spatial derivatives of the concentrations were set equal in the last two points making the second derivatives vanish. Altogether, 50 constant steps were used to guarantee the accuracy when solving the mass balances with the numerical method of lines. The major steps in the derivation of the step-by-step calculation equation for the concentration profile inside particles were to use dimensionless variables, to collect model constants in a single term, to use the product differentiation rule, to approximate the spatial derivatives with a central difference,8 and to use a constant step length. The concentration in the next step could then be calculated always from the concentration in the previous two steps. The first two concentrations were required. The concentration in the particle center was the initial guess. The concentration after the first step was the same because the derivative at the particle center was set to zero. The calculation could then progress step by step. As the surface was reached, the dimensionless concentration was to be 1. If not, the initial guess had to be adjusted. A gradient-based routine was applied to perform the iteration. A total of 100 constant steps were used. A spherical shape was assumed, and the shape factor was 2. 8859

dx.doi.org/10.1021/ie2019678 |Ind. Eng. Chem. Res. 2012, 51, 8858–8866

wL

∑

N

∂Ci, L ∂ Ci, L þ εL Da, L 2 þ Ni, L av þ εL ηe, i vi, j rj aj ∂l ∂l j¼1

!

ð4Þ

ð1Þ Ci, L ¼ Ci, F , l ¼ L

l ¼ 0

∂Ci, G 1 ¼ εG ∂t

8860

p5 CA ð1 θC Þ 1 þ k2 CA þ k3 CH2 þ k4 CB

Peclet number: wL L Pe ¼ Da, L

Boundary conditions: dCi ¼ 0 at r ¼ 0 and Ci ¼ Ci, S at r ¼ R dr

Activity trend expression: dα αp5 CA ¼ dt 1 þ k2 CA þ k3 CH2 þ k4 CB

rC ¼

Coking reaction rate:

Main reaction rate expression: p1 CA CH2 rA ¼ ð1 þ k2 CA þ k3 CH2 þ k4 CB Þ2

Henry's law: CH , L pH2 ¼ H 2 CTOT, L

Gasliquid mass transfer: NH2 , L av ¼ kgl aðCH2 , L CH2 , L Þ

ð24Þ

ð21Þ

ð19Þ

ð17Þ

ð15Þ

ð12Þ

ð10Þ

Mass balance for the tracer test: ∂CA, L ∂2 CA, L ∂CA, L ¼ Da, L w ∂t ∂l2 ∂l

Shape factor definition: Ap s ¼ R1 Vp

Mass balance for particles: dCi s r d De dr rj ¼ 0 rs dr

a ¼ α2 ¼ ð1 θC Þ2

Activity:

Side reaction rate expression: pffiffiffiffiffiffiffiffi rB ¼ p6 CB CH2

Friction factor: dp 2ff FG wG 2 ¼ dl dK

ð25Þ

ð22Þ

!

rA, x Vx rA, S Vp

x¼1

∑

N1

CA, L ¼ 0 w 1 at l ¼ 0

Boundary conditions:

ηe, A ¼

0

1

ð11Þ

ð6-9Þ

ð5Þ

ð3Þ

ð2Þ

∂CA, L ¼ 0 at l ¼ L ∂l

ð26Þ

ð23Þ

ð20Þ

ð18Þ

ð16Þ

17:3 B C ff ¼ @31:3 þ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiA ð14Þ ½XG ðReL WeL Þ1=4 3=2 XG ðReL WeL Þ1=4 1

Pressure gradient:

Effectiveness factor:

ð13Þ

k ad 2 dp gl p ¼ 2:0js 0:2 ReL 0:73 ReG 0:2 ScL 0:5 εL dpipe DL 1 εP

Gasliquid mass-transfer coefficient: !0:2

∂Ci, G ∂wG ∂εG ∂2 Ci, G wG Ci, G Ci, G þ εG Da, G 2 Ni, L av ∂l ∂l ∂t ∂l

Gas mass balances in general form:

∂ Ci, L ¼ 0, ∂l2 2

Boundary conditions:

Gas mass balance in our simplified case (small Δp, constant T, and single gas component, dilute liquid) _ n p 1 ∂wG ≈ constant, εL, tot ¼ f ðRel Þ ≈ constant w εG ≈ constant w Ci, G ¼ Ni, L av CH2 , G ¼ _H2 ¼ V H2 RT εG ∂L

Dynamic liquid hold-up: !0:33 dp ReL 0:14 εL, dyn ¼ 0:16 dpipe

∂Ci, L 1 ¼ εL ∂t

2

Liquid-side mass balance for component i:

Table 2. Model Equations

Industrial & Engineering Chemistry Research ARTICLE



Figure 2. Effect of the main reaction rate constant p1.

ARTICLE

Figure 5. Effect of the deactivation rate constant p5 [s1].

Figure 3. Effect of p1, deactivation present. Figure 6. Effect of the Pe number, severe deactivation.

Figure 4. Effect of the effective diffusivity. Figure 7. Concentration trends at various locations.

An analytical solution6 was available for the Pe number evaluation from a step-change experiment. The solution was an infinite series, the terms of which decreased rapidly. One had to use gradient-based iteration for every term in the series. Already, the first five terms gave quite accurate results. 2.5. Parameter Evaluation. First, multiple single-parameter change simulations were made to figure out the response behavior of the models. The selected parameters were the reaction rate constant p1, the coke formation constant p5, the pore diffusivity De, the gasliquid mass-transfer coefficient kgla, and the Peclet number Pe. Thereafter, an optimization study was made. In this study, only the main reaction and coking reaction rates were optimized because for the other above-mentioned parameters there were predictions available.

3. RESULTS FROM CITRAL HYDROGENATION 3.1. Sensitivity Study of the Model. The effects of the key parameters on the product concentration were evaluated. The key parameters were the reaction and deactivation parameters p1 and p5, the pore diffusivity De, the gasliquid mass-transfer coefficient kgla, and the Peclet number Pe. Because the selectivity was high, Figures 28 present only the main product. The simulation inputs are given in Table 3, and the results are displayed in Figures 29. Simulations 1 and 2 were performed to reveal how the product concentration was influenced by the changes in the reaction rate constant p1 in the absence and presence of severe deactivation. The effect of the pore diffusivity 8861



ARTICLE

was then investigated in simulation 3 in the presence of catalyst deactivation. Thereafter, the effect of the coking rate constant was checked in simulation 4. The effect of the Peclet number was studied in simulation 5. Simulation 6 was a single simulation showing how the concentration profiles emerge as a function of the location in the catalyst bed. The effect of the mass-transfer parameter (kgla) was examined in simulation 7. The reactors operated almost isothermally because of dilute conditions, i.e., low concentrations of citral in the solvent. The heat of reaction was 28.9 kJ/mol. If the liquid had absorbed all of the heat, the corresponding temperature rise would have been only around 0.2 °C. The pressure drop according to eqs 13 and 14 was with the larger flow only around 1% of the operating pressure. The reactor length was 0.05 m. A longer reactor would have a more profound plug-flow behavior and a larger pressure drop. The gas holdup was around 0.2. For the gasliquid mass-transfer coefficient, the value 0.2 s1 was selected because it represented the special value on the verge when mass transfer just started to restrict the production rate at the highest reaction rate observed in the experiments. The masstransfer coefficients were checked with the equations of Wild and Fukushima,12 and both correlations gave a considerably larger value for the mass-transfer coefficient than the 0.2 s1 used. So, the used value was on the safe side. The mass-transfer coefficient was given as a constant, which was justified by having the physical properties nearly constant because of a narrow temperature range (20 °C) and by having flows almost constant (liquid dilute and gas used in excess). The united gasliquid mass-transfer coefficient kgla was used, and the specific mass-transfer area, consequently, was not required separately. A spherical size was used. The shape factor was then 2.

The extreme effective diffusivity, 1 1012 m2/s, was selected to show how unrealistically small diffusivities are required before particle diffusion started to reduce the productivity. The porosity of a good commercial catalyst is up to 0.7. The tailor-made catalyst in our laboratory had a somewhat lower porosity, 0.65. The adsorption parameters were taken from the literature,14 k2 = 4.54 L/mol, k3 = 12.9 L/mol, and k4 = 4.89 L/mol. Although one study has shown16 them to be temperature-dependent while having platinum as the catalyst, the temperature dependence in this study with nickel as the catalyst was carried out so that the temperature dependence of the main reaction was completely included in p1. The parameters above are for nickel as the catalyst, and the optimization study was also done in a similar way.14 The denominator of the reaction rate expression was D = (1 + k2CA + k3CH2 + k4CB)2. It took into account adsorption of the raw material, hydrogen, and the main product, respectively. The adsorption term of hydrogen was the largest of these. Because of the relatively slow reaction, the liquid was saturated with hydrogen. When using a short reactor with a negligible pressure drop, the concentration of hydrogen in the liquid in each experiment remained almost constant. The solubility of hydrogen depended slightly on the temperature (Henry's law constant: 416 MPa at 40 °C and 387 MPa at 60 °C). This meant that the adsorption term k3CH2 for hydrogen was 0.320.35. The adsorption terms of the raw material and main product put together were around 0.0900.098. The whole denominator varied then between 2.00 and 2.09, therefore being almost constant. Mathematically, the kinetic rate equation thus approached a power law. The rate constant was selected to be adjustable because trickle-bed reactor (TBR) fluid dynamics deviates so much from the operation of a batch reactor. In our model, fluid dynamics was described with an axial-dispersion model, which is a great simplification when, in

Figure 8. Effect of kgla, severe deactivation.

Figure 9. Extent of the side reaction.

Table 3. Sensitivity and Simulation Study Variables and Their Values run

De

studied variable

1 2

reaction rate constant p1 reaction rate constant p1

3a

effective diffusivity De

p1

Pe

kgla [s1]

p5 [s1]

no resistance no resistance

0.252.5 15

3 3

0.2 0.2

0.0005 0.04

10101012 m2/s

0.016

1.22

3

0.2

4

coking rate constant p5 [s ]

no resistance

2.5

3

0.2

0.0010.08

5

Peclet number Pe

no resistance

2.5

110

0.2

0.04

6

location in pipe

no resistance

2.5

3

0.2

0.04

7

kgla [s1]

no resistance

2.5

3

0.0010.1

0.04

1

a The pore diffusivity run was done by using the biggest particle size, 250 μm, and the maximum reaction rate parameter values from the parameter estimation below.

8862



Figure 10. Initial rates of the hydrogenation reaction.1

Figure 11. RTD as a function of the Peclet number: (b) measurements giving the average Pe = 3.

reality, the liquid tends to flow as rivulets and stagnant zones also exist. The rate constant adjustment was used to compensate for concentration nonuniformities caused by complex threephase flow. If the adsorption parameters were not at all known, then the estimation could have been done with a strategy similar to that used now but for a larger number of parameters: by first making a multitude of simulations with a wide range of parameters to get a clear idea on how the system behaves and where to focus the final search. In this special case, because the denominator remained almost constant and the same in all of the runs, a power law could also have given an agreeable result. Concentration Trends of the Main Product. The increase in p1 made the startup curve steeper and increased the steady-state product concentration (Figure 2). The steady state depended nonlinearly on p1. At high p1, the concentration approached total conversion, 0.02 mol/L. Parameter p1 has a similar influence on the product concentration under catalyst deactivation (Figure 3), although deactivation gradually decreases the product concentration. The pore diffusivity depends on the diffusivity in the liquid and pore structure. The structural effects were combined into a single correction factor,5 the magnitude of which was 0.1. Because the pore diffusivity was larger than 1011 m2/s, the effect of pore diffusion was minimal. Below 1011 m2/s, it started to decrease the reaction rate. Pore diffusion changed the slope of the startup curve and the steady-state concentration (Figure 4). Because the effective diffusivity of citral for this case was roughly 1010 m2/s, the diffusion resistance was minimal. The deactivation model had as the key parameter the coking rate constant, p5. The adsorption parameters were the same as those for the main reaction. When p5 was increased, the rate at

ARTICLE

which the activity dropped speeded up. The concentration trend became more nonlinear, as can be seen in Figure 5. The coking effect, being dependent on the raw material concentration, reduced the difference between the results obtained by higher and lower Pe numbers (Figure 6). With a large Pe number, the reactor started to resemble a plug-flow unit, whereas with a small one, it acted more like a continuously stirred tank reactor (CSTR). In a pipe-flow type of operation, the concentration is initially high, causing coking to take place at the inlet of the reactor. In a stirred-tank-like operation, the concentration is constant and coking even. The product concentration obtained with a single CSTR is far less than that of an ideal plug-flow reactor of the same size, when concentrationdependent deactivation is absent. The reactor loses its effectiveness mostly in the beginning, and reaction takes place further and further in the reactor as catalyst deactivation progresses. Figure 7 illustrates how the concentration profile develops within the reactor for a single case. Figure 8 shows how low the value of kgla has to be in order to start to restrict the reaction rate. For the example system, kgla of around 0.2 s1 was the calculated value, which was close to the critical limit. Figure 9 shows the extent of the side reaction, which was really minor in the actual case. The temperature range in our simulations was narrow, 20 °C. Because catalyst deactivation was predominant and present even when the first samples were taken, as Figure 10 clearly indicates, even the initial rates of the main reaction did not follow the exponential behavior of the Arrhenius equation. Deactivation was regarded as the main reason behind the low activation energies. The major thing that explains the difference in the kinetics between a batch and a continuous system is in the fluid dynamics. Flow conditions are very different in a TBR compared to a stirred tank. An axial-dispersion model was used to explain axial mixing. It is still a very simple way of modeling three-phase fluid flow. In a stirred tank, it is easier to approach intrinsic kinetics than it is in a TBR. One can set the stirrer speed high enough to eliminate the effects of nonuniformities. In a TBR, liquid tends to flow in rivulets and some part of it stays more stagnantly in the pockets between particles. Another cause for the difference is that deactivation in the TBR is dependent on the axial location, whereas in an ideally stirred tank, it is uniform. 3.2. RTD. The aim was to get an estimate for the Peclet number from RTD measurements. An example is presented in Figure 11. The catalyst particle size in these experiments was 150250 μm. The result was obtained by taking into account only the points in which the concentration clearly deviated from the initial and final values. Only a few samples from the small reactor could be taken during the signal rise period. Experiments were performed by changing the citral concentration from zero up to 0.02 M. For the current system, the Pe number was in the range of 15, preferably 3, which was used in the parameter estimation task. This implies that axial dispersion was clearly present. 3.3. Parameter Estimation. The LevenbergMarquardt algorithm was used for the final search of the parameters.17 The number of parameters in the estimation task was minimized. The Pe number was obtained from RTD and was set as 3, which demonstrated the axial-dispersion effect. An estimate for the effective diffusivity was obtained by applying a correction factor of 0.1 to the liquid diffusivity;5 this accounts for the effect of the catalyst activity and tortuosity. The adsorption parameters k2, k3, and k4 were taken directly from the literature14 and kept constant 8863



ARTICLE

In a TBR, liquid-flow paths, local flow velocities, and flow directions differ along the way. In this study, the well-established axial-dispersion model was used to take into account flow nonidealities along the main flow direction. In reality, flows presumably vary with some random degree in all directions. The axial-dispersion coefficient was estimated from separate residence time measurements. It is not only the residence time, however, that counts but also the kind of local conditions that are faced along the way. The axial-dispersion coefficient influences the productivity, and in this way, it is linked with kinetics and deactivation in a TBR. Even when the experiments are carried out in the kinetic regime in a TBR, with no mass-transfer resistances, this coupling is still present. Only if the unit worked as a complete ideal plug-flow unit, this coupling would not be there.

Figure 12. Citronellal concentration trends in the reactor: marks = experimental values; lines = model values.

Table 4. Parameter Values (Pe = 3; Effectiveness Factor = 1) T/°C

p1

p5 [s1]

degree of explanation

standard error

40

0.374

0.00814

82.79

0.684 104

45 55

0.704 1.026

0.0114 0.0152

97.97 94.87

0.318 104 0.644 104

60

1.218

0.0160

93.17

0.693 104

during the estimation. Thereafter, in the model for the main reaction, there existed only two parameters, the value of which was to be estimated. These were the reaction rate constant, p1, and the coke formation reaction constant, p5. In this way, the parameters became well identified. As discussed earlier, the whole denominator of the LH rate expression varied between 2.00 and 2.09 in our case. It was almost constant and, mathematically, the kinetic rate equation was very close to a power law, and the kinetic part could be adjusted by changing only one parameter, the rate constant. The parameters were estimated individually for each temperature. The experimental datapoints and produced curves are shown in Figure 12. The agreement between the model and experiments is reasonably good, as can be seen from the degrees of explanation and standard errors in Table 4, although some deviations were observed. The catalyst deactivation was more strongly present in the beginning of the experiments than what was predicted. The parameter values are presented in Table 4. In terms of a more general case, the difference that adsorption parameters can make is that at high concentrations they can change the concentration dependency of the reaction rates on individual reactants. Then, for those kinds of cases, it would be beneficial to have a concentration profile curve from experiments carried out in the kinetic regime over the whole wide concentration range and study the dependency in that way. In continuous TBR, the study would demand a long reactor with multiple sampling points. In a batch operation, it would demand longer reaction times. One alternative for a general case would be to do the estimation the same way as now, but for a larger number of parameters, it becomes increasingly difficult. One should then make a lot of simulations first to get an idea of the sensitive ranges of each parameter and the interactions between the parameters.

4. CONCLUSIONS Modeling of a three-phase continuous catalytic TBR is a highly complicated task. Many of the phenomena produce nearly similar effects. The best would be if all of the phenomena could be studied separately. This is extremely seldom the case for a TBR because many of the phenomena are coupled and linked together by complex fluid dynamics of the three-phase flow. A lot of simplifying assumptions are required to obtain a tractable model. The modeling approach based on the component mass balances was selected in this study. Axial dispersion in the liquid phase was taken into account. Catalyst deactivation was modeled as coking. Pore diffusion was calculated using a numeric step-by-step routine. The main reaction and consecutive side reaction were studied. The modeling work was done to explore what took place in laboratory-scale experiments, which of the phenomena were present, and to what extent they had an impact on the productivity. The second reason for the effort was to use the data for model development, testing, and verification purposes. The model was used for studying citral hydrogenation to citronellal. It was first used for a sensitivity study to reveal how key parameters did change the product concentration. Sensitivity analysis revealed the following things. The effect of the rate constant of the main reaction on the productivity was pronounced far away from the total conversion, whereas it became smaller close to the total conversion. The deactivation rate constant did not remarkably change the achieved maximum product concentration. It dominated all the more later on in the experiments. The effect of the effective diffusivity on the productivity was only minor in the realistic range of values. The Peclet number had a significant effect when working with rather high conversions. The dominating deactivation softened the effect later on in the experiments. The gasliquid mass-transfer coefficient should be unlikely small to start limit the productivity. The parameter estimation of the main reaction and coking rate was finally carried out. A reasonable agreement was achieved. The gasliquid mass transfer and pore diffusion were found not to be rate-limiting in the example case, whereas catalyst deactivation and axial dispersion were very dominant features in the system. TBR is a powerful tool to reveal catalyst deactivation because in it deactivation becomes visible as the conversion decreases during the time on stream. Furthermore, six TBRs in parallel were used to intensify generation of the experimental data.1 When deactivation in a batch reactor is compared to that in a 8864


Industrial & Engineering Chemistry Research TBR, one of the main differences is that, in a TBR, deactivation takes place locally, whereas in an ideally stirred reactor, it is uniformly present. This is because deactivation was raw material concentration-dependent in this case and close to the feed entrance of the trickle bed, the concentration was at its highest and deactivation, consequently, as its strongest. The degree of axial dispersion (value of the axial-dispersion coefficient) also in itself directly influences the productivity of the TBR reactor and in that way also the local concentrations and local deactivation.

’ AUTHOR INFORMATION Corresponding Author

*E-mail: teuvo.kilpio@abo.fi.

’ ACKNOWLEDGMENT This work is a part of the activities at the Åbo Akademi Process Chemistry Centre within the Finnish Centre of Excellence Program (20062011) by the Academy of Finland. ’ NOTATION a aj av Ap CA, CB, CH2

activity at time t, in location l activity of the jth reaction specific mass-transfer area [m2/L] particle surface area [m2] local concentration of raw material, product, and H2 [mol/L] concentration of raw material, main product, CA,L, CB,L, CC,L and side product in the liquid [mol/L] CH2,G, CH2,L, CH2,L* concentration of hydrogen in gas, in liquid, and in saturation [mol/L] total concentration in the liquid [mol/L] CTOT,L concentration of component i [mol/L] Ci concentration of component i in the liquid at Ci,L, Ci,S, Ci,G surface and in gas [mol/L] KrischerKast hydraulic diameter [m], dK dp[16εP3/(9π(1 εP)2)]1/3 dp particle diameter [m] pipe diameter [m] dpipe D denominator of the reaction rate equation diffusivity in the liquid [m2/s] DL axial dispersion coefficient in gas and liquid Da,G, Da,L [m2/s] effective pore-diffusion coefficient [m2/s] De two-phase friction factor ff H Henry's law constant [Pa] i, j indexes reaction rate constants in convenient units kA, kB, kC* to produce rate in units [mol/(L s)] united mass-transfer coefficient for gasliquid kgla mass transfer [s1] adsorption coefficients for the raw material KA, KH2, KB and hydrogen [ mol1] l location [m] L length of the reactor [m] N total number of reactions, eq 1 N total number of points in particle, eq 23 molar flux of component i in liquid Ni,L, NH2,L [mol/(m2/s)] pressure and partial pressure of H2 [Pa] p, pH2

ARTICLE

p1, k2, k3, k4, p5, p6 parameters of reaction rate equations in convenient units to produce reaction rate in units [mol/(L s)], where p1 = kAKAKH2, k2 = KA,k3 = KH2, k4 = kB, p5 = kC*KA dp/dl pressure gradient for the reactor [Pa/m] Pe Peclet number reaction rate of the main, side, and coking rA, rB, rC* reactions [mol/(L s)] reaction rate of the main reaction at the rA,S, rA,x surface and location x [mol/(L s)] reaction rate of the jth reaction [mol/(L s)], rj eq 1 R characteristic mean radius of catalyst particles [m] particle Reynolds number in liquid, ReL FLwLdp/μL ReG particle Reynolds number in gas, FGwGdp/μG S selectivity s shape factor, s = 2 for a sphere, s = 1 for a long cylinder, and s = 0 for a flake Schmidt number, μL/FLDL ScL t time [s] T temperature [K] stoichiometric coefficient for component vi,j i in reaction j particle volume and volume of the xth eleV p, V x ment of the particle [m3] · volume flow of hydrogen [m3/s] VH2 superficial gas velocity [m/s] wG superficial liquid velocity [m/s] wL w wL/εL WeL Weber number in liquid, FLwL2dp/σ x dimensionless distance from center = r/R modified LockhartMartinelli ratio, wGFG1/2/ XG wLFL1/2 α activity coefficient volume fraction of liquid, dynamic volume εL, εL,dyn fraction of liquid porosity in packing εp gas holdup εG surface shape factor for the particle ϕs effectiveness factor for component i and ηe,i, ηe,A component A viscosity of the liquid and gas [kg/(m s)] μL, μG density of the liquid and gas [kg/m3] FL , FG coke coverage θC* σ surface tension [N/m]

’ REFERENCES (1) M€aki-Arvela, P.; Er€anen, K.; Alho, K.; Salmi, T.; Murzin, D. Multi-Tubular Reactor Design as an Advanced Screening Tool for Three-Phase Catalytic Reactions. Top. Catal. 2007, 45, 223. (2) Levenspiel, O. Chemical Reaction Engineering; John Wiley & Sons: New York, 1999. (3) Froment, G.; Bishoff, K. Chemical Reactor Analysis and Design; John Wiley & Sons: New York, 1979. (4) Belfiore, L. A. Transport Phenomena for Chemical Reactor Design; John Wiley& Sons: New York, 2003. (5) Butt, J. B. Reaction Kinetics and Reactor Design; Marcel Dekker Inc.: New York, 2000. (6) Salmi, T.; Mikkola, J. P.; W€arna, J. Chemical Reaction Engineering and Reactor Technology; CRC Press: Boca Raton, FL, 2009. 8865



ARTICLE

(7) R€onnholm, M.; W€arna, J.; Salmi, T. Comparison of Three-Phase Reactor Performances with and without Packing Elements. Catal. Today 2003, 7980, 285. (8) Hoyos, B.; Cadavid, J. G.; Rangel, H. Formulation and Numeric Calculation of Non-Isothermal Effectiveness Factor for Finite Cylindrical Catalysts with Bi-Dimensional Diffusion. Lat. Am. Appl. Res. 2004, 34, 17. (9) Lange, R.; Schubert, M.; Bauer, T. Liquid Holdup in Trickle-Bed Reactors at Very Low Liquid Reynolds Numbers. Ind. Eng. Chem. Res. 2005, 44, 6504. (10) Al-Dahhan, M. H.; Larachi, F.; Dudukovic, M. P.; Laurent, A. High-Pressure Trickle Bed Reactors: A Review. Ind. Eng. Chem. Res. 1997, 36, 3292. (11) Chaudhari, R. V.; Jaganathan, R.; Mathew, S. P.; Julcour, C.; Delmas, H. Hydrogenation of 1,5,9-Cyclododecatriene in Fixed-Bed Reactors: Down- vs. Upflow Modes. AIChE J. 2002, 48, 110. (12) Midoux, N.; Morsi, B.; Purwasasmita, M.; Laurent, A.; Charpentier, J. Interfacial Area and Liquid Side Mass Transfer Coefficient in Trickle Bed Reactors Operating with Organic Liquids. Chem. Eng. Sci. 1984, 39, 781–794. (13) Poling, B.; Prausnitz, J.; O’Connell, J. The Properties of Gases and Liquids; McGraw-Hill Book Co.: New York, 2001. (14) Salmi, T.; M€aki-Arvela, P.; W€arna, J.; Er€anen, K.; Denecheau, A.; Alho, K.; Murzin, D. Y. Modelling of Concecutive Reactions with a Semibatch Liquid Phase: Enhanced Kinetic Information by a New Experimental Concept. Ind. Eng. Chem. Res. 2007, 46, 3912. (15) Cussler, E. L. Diffusion Mass Transfer in Fluid Systems; Cambridge University Press: New York, 1997. (16) Mukherjee, S.; Vannice, M. Solvent Effects in Liquid-Phase Reactions II. Kinetic Modeling for Citral Hydrogenation. J. Catal. 2006, 243, 131. (17) Haario, H. Modest User’s Guide; Profmath Oy: Helsinki, Finland, 1994.

8866


Modeling of a Three-Phase Continuously Operating Isothermal

Recommend Documents