Modeling of microreactors for ethylene epoxidation and total oxidation

Chemical Engineering Science 134 (2015) 563–571

Contents lists available at ScienceDirect

Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces

Modeling of microreactors for ethylene epoxidation and total oxidation V. Russo a,b, T. Kilpiö a, J. Hernandez Carucci a, M. Di Serio b, T.O. Salmi a,n a b

Åbo Akademi University, Department of Chemical Engineering, Laboratory of Industrial Chemistry and Reaction Engineering, FI-20500 Turku/Åbo, Finland University of Naples “Federico II”. Chemical Sciences Department, IT-80126 Naples, Italy

H I G H L I G H T S


Ethylene oxide kinetics was studied.
A kinetic model was developed.
A good description of data was obtained.

art ic l e i nf o

a b s t r a c t

Article history: Received 6 December 2014 Received in revised form 18 March 2015 Accepted 12 May 2015 Available online 9 June 2015

Microreactors are especially well-suited for laboratory-scale studies of rapid exothermic reactions, because they have excellent mass and heat transfer characteristics. Two novel reactor models were used for plate microreactors, one for a washcoated reactor and the other one for a silver plate microreactor. The aim of this modeling study was to precisely explain both the concentration and the temperature dependencies of the reaction rates of ethylene oxide (EO) synthesis over a wide range of operating conditions (p, T, flows) by using all the extensive data generated with our microreactors. Microreactors are especially well suited for producing accurate kinetic data because of the uniformity of reaction conditions in the reactor system. Modeling took into account the reaction and mass transfer effects. The axial and radial concentration profiles and for the washcoated case also the concentration profiles along the coating direction were solved numerically. The models were based on dynamic mass balances for the gas phase, with convection, axial and radial dispersion terms included. From the modeling viewpoint, an interesting phenomenon, the interaction between intrinsic kinetics and diffusion in the porous catalyst layer, was tackled with the aid of mathematical modeling for the washcoated reactor case. The reactor models were solved numerically by using gPROMS software. Ethylene oxidation on silver catalyst was the selected example reaction system, because the main product, ethylene oxide, is a key compound and an important intermediate for chemical industry. The reactor and kinetic models were able to describe all the experimental data with a very satisfactory agreement. The microreactor models developed are generic and applicable to various kinds of heterogeneously catalyzed gas-phase reaction systems. The industrial breakthrough of flow chemistry and catalyzed single phase microreactor technology is expected to take place first in the production of fine and specialty chemicals. For some bulk products (e.g EO), flow chemistry has already arisen interest. The benefits that microreactors can offer are the exact control of process conditions, the compact reactor size and the ease of scalability by simply “numbering up”. For a new product, ease of scalability shortens the time to market, whereas for a bulk product like EO, millireactors (structurally similar to microreactors, just having larger flow channels) may be more potential candidates for on-site production applications. & 2015 Elsevier Ltd. All rights reserved.

Keywords: Microreactor Ethylene oxide Washcoat Kinetics

1. Introduction Microreactors are strong tools in the investigation of the intrinsic kinetics of rapid gas-phase reactions. Two main properties of microreactors make them particularly attractive for this

n

Corresponding author.

http://dx.doi.org/10.1016/j.ces.2015.05.019 0009-2509/& 2015 Elsevier Ltd. All rights reserved.

purpose: the catalyst layers in microreactor channels are very thin, typically 20–30 μm or even less and the heat transfer characteristics of microreactors are excellent. Thin catalyst layers suppress the internal mass transfer resistance in the catalyst pores thus guaranteeing the operation within the regime of intrinsic kinetics, which simplifies the interpretation of experimental data. Mathematical modeling of microreactors is discussed in detail in some monographies (Hessel et al., 2004; Keil, 2007).

564

V. Russo et al. / Chemical Engineering Science 134 (2015) 563–571

The use of simplified description of fluid dynamics (dispersion models and semi-empirical equations for pressure drop) is widely accepted and used way to replace the momentum balances (Navier-Stokes equations, PDEs) in chemical engineering. The accurate solution of momentum balances demands a tight grid especially in places of high velocity gradients. By avoiding the evaluation of velocity profiles with momentum balances PDEs, one reduces the size of the problem essentially and can use also a less tight grid. The merit then is that calculation times shorten an order of magnitude (e.g. from hours to minutes). Parameter estimation task typically demands a multitude of simulations and the benefit of having shorter calculation times cumulates. The disadvantage of the simplification is that local velocities are not evaluated and only the net effect of backmixing is taken into account. On the hand, these systems have been sometimes treated with rough approaches, such as plug flow conditions, in terms of fluid dynamics. In this paper, an attempt in sufficiently describing the fluid dynamics of two types of microreactors for the ethylene oxide synthesis reaction was made. The differences in total calculation times depend on how much more complicated the other modeling options are. Computational Fluid Dynamics (CFD) is a popular modeling tool nowadays. Singlephase CFD models for this kind of simple geometry (when calculating only the local flows) are converging well perhaps within an hour (without reactions). However, it is hard to predict how much more time is consumed when non-linear Langmuir-Hinshelwood-HougenWatson (LHHW) reaction kinetic equations are incorporated. They may retard the calculation progress a lot. Another issue to keep in mind is that although CFD models are strong in calculating vector variables, with scalar variables, numerical diffusion is often faced. When multiple parameter estimations are carried out, one executes the program code tens, even 100 of times, which multiplies the run time difference. When doing model development with rapidly running program codes, a multitude of test runs can be performed each and every day and abundant amount of feedback is received. With CFD, the number of trials becomes limited. From the modeling viewpoint, microchannels and monolith reactor channels resemble each other a lot and both are, thus, shortly addressed. Models for monolith reactors both with and without a washcoat have been developed (Stutz and Poulikakos, 2008; Deutchmann et al., 2001; Pattegar and Kothare, 2004; Carucci et al., 2010; Hayes et al., 1999). All the above-mentioned studies were conducted for various different gas-only reaction systems. Plug flow, axial dispersion and backmixing models for microreactors are discussed in a recent review article of our group (Carucci et al., 2012). An illustrative parameter estimation study for ethylene oxide using a combined radial molecular diffusion and laminar velocity profiles in the model was included. Parameters for ethylene oxidation were estimated for a single temperature (Carucci et al., 2010). The current work is a summary of experimental and modeling activities carried out in our laboratory-scale microreactors for ethylene oxidation. This effort is a logical continuation of the previous work using the entire data sets together with this far the most detailed model for ethylene oxidation.

Fig. 1. Reaction scheme of direct ethylene oxide synthesis from ethylene and oxygen.

The reaction scheme for ethylene oxide synthesis is presented in Fig. 1. The exothermic nature of the reactions is clearly demonstrated there. The generally accepted view on the reaction mechanism is that ethylene epoxidation and total oxidation are essentially parallel processes; thus the role of reaction 3 in the scheme is minor compared to reactions 1 and 2. The high values of the reaction enthalpies imply typically nonisothermal conditions in standard reactors. Typical adiabatic temperature rise curves have been presented for the industrial ethylene oxidation. The point of how much the surface-to-volume ratio improves when scaling down from industrial pipes to microreactor channels has been highlighted (Kursawe, 2009). Selectivity/productivity results obtained by a simplified plug flow model have also been presented. The kinetic investigations pursued for revealing the exact kinetic expression of partial oxidation of ethylene into ethylene oxide using a silver based catalyst are summarized and rate equations for partial and total oxidation of ethylene are presented in a previous paper of our group (Salmi et al., 2013). Through the years, several hypotheses have been proposed for the molecular reaction mechanism of ethylene oxidation on silver surface (Carucci et al., 2010; Kursawe, 2009; Salmi et al., 2013; Borman and Westerterp, 1995; Ghazali et al., 1983; Lafarga et al., 2000; Petrov et al., 1988). In a recent molecular modeling study based on Density Functional Theory (DFT) calculations, a common surface intermediate was proposed for the formation of ethylene oxide and carbon dioxide (Özbek and van Santen, 2013).

2. Modeling approach and methods 2.1. General principles of the microreactor models The level of sophistication in the modeling of microreactors is a challenge. The key issues are the interaction of chemical kinetics, mass and heat transfer effects as well as the flow related phenomena in micro-channels. Rate equations describing the intrinsic kinetics of a particular chemical system can be found in literature or derived from plausible reaction mechanisms. An important issue is to reveal the role of the internal mass transfer resistance in the pores of the washcoat layer. Even though the catalyst layer is thin in microreactor channels, the reaction rate can be retarded by diffusion limitations, if the kinetic process itself is rapid. This has recently been illustrated (Schmidt et al., 2013). The role of flow modeling is somewhat contradictory: on one hand, the most sophisticated approach based on CFD should be used, but on the other hand, particularly if kinetic parameters should be estimated from experimental data, a simplified approach is attractive, based on the concepts of plug flow, laminar flow or dispersion effects. In general it can be stated that the model should not be too heavy from the computational point of view, if it is used for the estimation of kinetic parameters from experimental data. Here we present a parameter estimation approach based on reaction kinetics, internal diffusion and dispersion. The approach allows to describe all the important phenomena in the microreactor channels, still keeping the computational effort within very reasonable limits. The models presented in this work are dynamic, based on radial and axial dispersion and reaction was assumed to take place only where it practically proceeds, within the washcoat or on the silver wall. Two kinds of micro-reactor elements were used, wash-coated elements and silver plates in the experiments reported previously. Two separate models for heterogeneously catalyzed oxidation of ethylene to ethylene oxide were generated in this work. One model was developed for the washcoated reactor and the other one for the silver plate reactor. In the washcoated reactor model,

V. Russo et al. / Chemical Engineering Science 134 (2015) 563–571

both the reactions and internal diffusion of each chemical compound inside the washcoated layer were included. In the model for the silver plate reactor, the reaction was assumed to take place at the radially outermost (annulus shaped) volume element, which in practice represented the solid wall. Both models were based on mass balances with a simplified description of the fluid dynamics (axial and radial dispersion model combined with a semiempirical expression for pressure drop). Energy balances were not included, because microreactors are well known to possess excellent heat transfer characteristics and isothermal conditions prevailed in the experiments.

Dispersion in uF , C i,in rp Dz

565

r

Convective flow in

uF , C i,in

Outermost balance volume element

2.2. Washcoated plate microreactor model Eq. (1) gives the gas-phase mass balances for the washcoated reactor. The model was a dynamic one (including the accumulation terms) and convective and dispersion fluxes (both radial and axial) were included. The reaction terms were not present in the gas-phase mass balances, because the reactions took place within the washcoat layer. Rapid dynamics is characteristic for microreactors, so the use of dynamic models was not absolutely necessary, but they were used for getting a good initialization for the reactor model to enhance the convergence. The gas-phase mass balances were described with parabolic partial differential equations. Axial and radial dispersion terms were included and the form of the radial dispersion terms show how the cylindrical geometry was taken into account. Accumulation

Convection

Axialdispersion

∂C i;F ðt;z;rÞ ¼ ∂t

∂C i;F ðt;z;rÞ þ Dz;F ∂z

U

 uF U

∂2 C i;F ðt;z;rÞ þ Dr;F ∂z2

U




Radialdispersion



∂2 C i;F ðt;z;rÞ 1 ∂C i;F ðt;z;rÞ þ r U ∂r ∂r 2

ð1Þ Dynamic solid-phase mass balances are presented by Eq. (2). The balances included the internal diffusion term for cylindrical geometry

uF U



∂C i;F ðt;zÞ ∂z r ¼ R

¼ Dz;F U



Dispersion out uF , C i,out Dz

Accumulation
Internaldiffusion  Reactionrates  P D ∂C i;s ðt;z;rp Þ ∂2 C i;s ðt;z;r p Þ ∂C ðt;z;r Þ ¼ εefpf ;i U þ rsp U i;s∂rp p þ ε1p U νi;j Ur j ðt; z; rp Þ ∂t ∂rp 2

ð2Þ

The boundary conditions are collected in Eqs. (3)–(6). At the feed entrance, the concentrations of the reactants were used for both gas and solid phases, while the product concentrations were set to zero. At the outlet, the first derivatives of concentrations were set to zero in both gas and solid phases (the Danckwerts' closed boundary condition). In the radial direction, the concentration derivatives of all components in the gas-phase were set to zero at the very center due to the symmetry and on the outer surface of the washcoat concentrations were set to zero. The boundary conditions at the interface between the gas-phase and the washcoat layer were expressed so that the net flow of material by convection and axial dispersion from gas phase into the balance volume element was set equal to the mass flux to the washcoat layer, Fig. 2. This boundary condition arouse from the simplification that in the outermost volume element of the gas side, radial dispersion was not included. This simplification was justified by the small size of the outermost volume element and by the fact that in the calculations of fluid dynamics. Entrance  C i;s ðt; r p Þz ¼ 0 ¼ C i;F IN

ð3Þ

Convective flow out RW R

uF , C i,out

Fig. 2. Fluxes of all components in and out the outermost balance volume. Notation of various radii included.

Outlet  ∂C i;F ðt; rÞ ¼0  ∂z z¼L

 ∂C i;s ðt; r p Þ ¼0  ∂z z¼L

ð4Þ

CenterWall  ∂C i;F ðt; zÞ  ∂r

 ∂C i;s ðt; zÞ ¼0 ∂r p rp ¼ RW

ð5Þ

r¼0

Catalystsurface    C i;s ðt; zÞrp ¼ R ¼ C i;F ðt; zÞr ¼ R

D ;i U s ∂C i;s ðt;zÞ ∂2 C i;F ðt;zÞ þ RefWf   R U ∂r p r ¼ R ∂z2 r¼R p W

and the reaction terms. The active sites in the washcoat layer were assumed to be evenly distributed.

 C i;F ðt; rÞz ¼ 0 ¼ C i;F IN

Flux to the washcoat

¼0

ð6Þ

Generally speaking, whenever using momentum balances or CFD, the boundary condition for the flow on walls is that the axial velocity is set zero. However, it is exactly zero on the wall and normally high axial velocity gradients are present already “near the wall”. Therefore, even for the outermost volume element which represented our boundary condition, the axial convection term was included. For that element also the axial dispersion term was included while the radial dispersion was neglected. The boundary condition stated then: what arrives the outermost balance volume by convection or backmixing (axial dispersion) departs by convection or backmixing or if not, is transferred (to/ from) the washcoat. Fig. 2 above shows for the outermost volume element, the terms that were included. The boundary condition for the boundary at the interface of the reactor and washcoat is illustrated below in Fig. 2. Also the various radii are explained there. r is the radial location within the reactor (excluding the washcoat). rp represent the radial location within the washcoat. R is the outer radius of reactor (excluding the washcoat) and RW is the radius of reactor (including the washcoat), Fig. 2. The used discretization schemes were central difference approximations for axial and radial derivatives (Numerical Method of Lines). In particular, both axial and radial coordinates were approximated with first and second order central finite differences, with 50 discretization points for the reactor length and 10 for the reactor radius.

566

V. Russo et al. / Chemical Engineering Science 134 (2015) 563–571

2.3. Silver plate microreactor model In case of the silver plate reactor, the reactions proceed on the solid wall. Then, the reactions could be modeled to take place at the outermost radial volume element only. Mass balances were required only for the gas-phase, because the silver wall was nonporous. The gas-phase mass balances are actually exactly the same for all the gas-phase volume elements except for the outermost one. The mass balances for the outermost layer are different, because they include only the reaction terms. The gas-phase mass balances are given in Eq. (7). Accumulation ∂C i;F ðt;z;rÞ ¼ ∂t

 uF U

Convection ∂C i;F ðt;z;rÞ þ Dz;F ∂z

U

The experimental data suggest the following selection of the exponents: α1 ¼ 1, α2 ¼ 1/2, β1 ¼ β2 ¼1. In this way a practically applicable kinetic model is obtained, which is able to describe the experimental data. This approach also explains the dependence of the epoxide selectivity on the oxygen partial pressure.For the sake of simplicity, the contribution of dissociatively adsorbed oxygen was ignored in the denominator of Eq. (11). The temperature dependences of the model parameters can be modeled in different ways. The first one is to include it in the reaction rate constant only. The second way is to make temperature adjustment also for all the adsorption constants. Both approaches appear in

Axialdispersion Radialdispersion Reactionrates
2  P
  ∂ C i;F ðt;z;rÞ 1 ∂C i;F ðt;z;rÞ þ U þ r U ∂r νi;j Ur j ðz; rÞr ¼ R ∂r 2

ð7Þ

∂2 C i;F ðt;z;rÞ þDr;F ∂z2

By assuming that the reactions proceed exclusively at the microreactor plate wall, the reaction rate expressions were used only at the cylindrical volume element closest to the wall, i.e. r ¼ R. The boundary conditions are presented in Eqs. (8)–(10) and they are the same as the ones for the gas-phase of the washcoated reactor: for the inlet, concentrations were the feed concentrations, for the outlet, concentration derivatives were zero and for the radial center, concentration derivatives were zero due to the symmetry. Eq. (7) together with the reaction terms was used only for the outmost cylindrical volume element.

literature, but the former one is used more often, because it has far less adjustable parameters. Therefore, the former approach was selected to be tested first to find out, whether it could produce an adequate agreement with the experiments. This was done for both epoxidation and total oxidation. The dependence of the kinetic constants was described by a modified Arrhenius law, especially designed to improve the quality of parameter identification,  

Ea;j 1 1 kj ¼ kj;ref Uexp   ð12Þ U T T ref Rg

Entrance

 C i;F ðt; rÞz ¼ 0 ¼ C i;F IN

ð8Þ

Outlet

 ∂C i;F ðt; rÞ ¼0  ∂z z¼L

ð9Þ

Center

 ∂C i;F ðt; zÞ ¼0  ∂r r¼0

ð10Þ

2.5. Physical properties The molecular diffusivities of each reactive component were calculated in helium (B), which was used for diluting the oxygen/ ethylene mixture, with the correlation (Wilke and Lee, 1955; Reid et al., 1988, Poling et al., 2001) h i 1=2 3:03  ð0:98=M AB Þ UT 3=2 U 10  3 DA;B ¼ ð13Þ 1=2 P U M AB U σ 2AB U ΩD

2.4. Rate equations

The parameters in Wilke and Lee correlations were calculated starting from known physical properties, i.e. normal boiling point and critical volumes, as summarized in Eqs. (14)–(18)

The following rate expression was used for the ethylene epoxidation and total oxidation

 1 M AB ¼ 2 U ð1=M A Þ þ ð1=M B Þ

βj

ð1 þ K E U C E þ K O2 U C O2 Þ2

ð11Þ

where j¼1 and j ¼2 refer to epoxidation and total oxidation, respectively. Our experimental observations suggest that the effective reaction order with respect to oxygen was 0.9 in the epoxidation reaction, while it was less for the total oxidation (Salmi et al., 2013). The importance of each adsorption term in the denominator of rate expressions depends on the concentration and temperature ranges covered by the original data. Theoretically, all the reactants have to adsorb on the active sites before reacting there and all the reaction products should desorb away. The fact that our kinetic expressions were capable to describe accurately the concentration and temperature dependences of the reaction rates over the used operation conditions, gives convincing evidence that the proposed kinetic models were good enough ones. The inclusion of more terms might have improved the model accuracy a little, but it would also have added more parameters in the model, which in general impairs the accuracies of the parameters.

O /He

2

k1 UC E U C Oαj2

DAB [m /s]

rj ¼

ð14Þ

1.4x10

-4

1.2x10

-4

1.0x10

-4

8.0x10

-5

6.0x10

-5

4.0x10

-5

480

500

CO /He

520

H O/He

C H /He

540

C H O/He

560

580

Temperature [K] Fig. 3. Binary diffusion coefficients as a function of temperature. P ¼ 1 bar.

V. Russo et al. / Chemical Engineering Science 134 (2015) 563–571

Table 1 Wilke and Lee parameters for the main components in helium.

σ

Component (i)

a (m2 s  1 K  1)

b (m2 s  1)

Ethylene Oxygen Ethylene oxide Carbon dioxide Water

2.36  10  7 3.87  10  7 1.86  10  7 2.83  10  7 3.27  10  7

 5.51  10  5  8.50  10  5  4.67  10  5  6.74  10  5  8.46  10  5

1=3 AB ¼ 1:18 UV b

¼ 1:18 U ð0:285 UV 1:048 Þ1=3 c

0:193 1:03587 þ ΩD ¼ n 0:1561 þ expð0:47635 U T n Þ expð1:52996 UT n Þ ðT Þ 1:76474 þ expð3:89411 U T n Þ

ð15Þ

1:06036

Tn ¼

εAB k

kUT

εAB

¼ 1:15 U T b

ð16Þ

ð17Þ ð18Þ

By plotting the diffusion coefficients, calculated from the Wilke and Lee correlation, as a function of temperature, a linear trend was always observed as demonstrated by Fig. 3. In order to reduce the number of the equations constituting the numerical problem, the diffusion coefficients were then considered linearly temperature dependent according to the following equation, DAB ¼ ai U T þbi

ð19Þ

The values of the parameters (a and b) for each component are reported in Table 1. The parameters were obtained by applying linear regression to the data points which were obtained by using the Wilke and Lee correlation. Effective diffusivities were calculated from liquid diffusivities by making an order of magnitude correction due to porosity and tortuosity, 1/10, as suggested by (Belfiore, 2003). The thickness of the washcoat layer was 15–20 μm. Due to the preparation method (impregnation on the wall) no precisely regular pore size can be given. The pore channels that can be seen in SEM pictures are roughly 20 μm in diameter and porosity is in the range of 0.3–0.4. As the washcoat model takes into account the internal diffusion of each component in/from the catalytic layer, effective diffusivities must be provided. The mean transport pore model was used (Toppinen, 1997). By considering that the experiments of ethylene oxide synthesis were performed at different temperatures and pressures, gas density depended on the operation conditions. For this purpose, the ideal gas law was applied. The gas density is molecular mass dependent; for the main component, helium, the molar mass is 4 g/mol, whereas for oxygen and ethylene it is nearly 30 g/mol. With 70% of helium, the average molar mass is around 11 g/mol. Finally, the dispersion coefficients were calculated by using an order of magnitude value of 100 for the Peclet number. This represents a typical order of magnitude value for microreactors working in this kind of gas-phase application, Eqs. (20)–(21). The radial dispersion coefficient was taken as one third the axial dispersion coefficient. The order of magnitude is based on earlier work conducted by (Arve (2005)). In their studies, Pe number for hydrocarbons was in the range of 40–60 in the microreactor experiments. For so high Pe numbers, the productivity and selectivity are only very weakly dependent on the Pe number, especially in diluted systems as in the present case. One benefit of the type of microreactors used in our study (multiple narrow channels) is that they can approach the plug flow

567

operation mode pretty well. The Peclet number value of 100 corresponds to this kind of a situation. This fact was verified by repeating some simulations with different Pe numbers within the range, observing no significant changes in the outcome of both models. Dz;F ¼

uF UL Pe

ð20Þ

1 Dr;F ¼ U Dz;F 3

ð21Þ

The value of radial dispersion coefficient was taken as 1/3 of the axial dispersion coefficient. The value is based on a heuristic rule. Radial dispersion is surely weaker than axial but to know it precisely, it would demand local measurements. In larger scale fixed bed reactors, it could be found out by making tracer test and by having a radial and axial dispersion coefficients present to be estimated from the obtained data. For microreactors it is hardly measurable. Simulations were added to illustrate that radial dispersion was not dominant in our experiments and the effect of it in our case was minor. To test the sensitivity of the reaction system on radial dispersion, additional three simulations with 1/5, 1/10 and 1/15 times the axial dispersion coefficient value were conducted. Quite the same conversion was reached (0.18322188, 0.1832213 and 0.18322243) for 1/5, 1/ 10 and 1/15 respectively which means that our heuristic value was within a reasonable range.

2.6. Numerical methods The model equations were solved simultaneously numerically for both the gas-phase and the catalyst washcoat or wall by using gPROMS ModelBuilder software. The sum of squares, i.e. the quadratic sum of the differences between experimental and modeled concentrations was minimized to get the best possible fit of the model to the experimental data. The model equations were solved simultaneously for both the gas-phase and the catalyst washcoat or wall by using gPROMS ModelBuilder software. The PDE systems were solved by using the built-in method of lines. The spatial coordinates were approximated with central differences (50 points for length, 10 for the radius). The sum of squares, i.e. the quadratic sum of the differences between experimental and modeled concentrations was minimized to get the best possible fit of the model to the experimental data. The solvers needed to perform both parameter estimations and simulations were chosen from gPROMS library. The approach that the software uses is the Maximum Likelihood Method. This method allows one to introduce uncertainty and standard deviation for each and every measurement and analysis, thus weighing the data points and the measurements which are most reliable. However in our case, the deviation was set as constant which implied that the measurement data and analysis were regarded as reliable. The target became then to minimize the least squares of the target variables: concentrations. Table 2 Experimental data for washcoated reactor. Ag/Al2O3, (Carucci et al., 2010) Exp.

Cin(E) [mol/m3]

Cin(O2) [mol/m3]

T/[K]

uF/[m/s]

1 2 3 4 5 6 7 8

4.35 4.35 4.35 4.35 4.88 4.69 4.51 4.2

2.17 2.17 2.17 2.17 2.44 2.34 2.26 2.1

553 553 553 553 493 513 533 573

0.1 0.12 0.15 0.17 0.1 0.1 0.1 0.1

568

V. Russo et al. / Chemical Engineering Science 134 (2015) 563–571

0.3

Table 3 Kinetic and adsorption parameters obtained for the washcoat microreactor. Kinetic constants (subscript ‘ref’) are calculated at 553 K. C.I. stands for the confidence interval which is related to a certain percentage.

X Y

X E [%], YEO [%]

Simulation

0.2

0.1

0.0 0.08

0.10

0.12

0.14

0.16

0.18

0.20

Fluid velocity [m/s] Fig. 4. Experiments and modeling of washcoated microreactor. Ethylene conversion and ethylene oxide yield as a function of the fluid velocity. The symbols represent the experimental data, and the continuous lines the simulations. Temperature was 553 K.

0.4 X Y

XE [%], YEO [%]

0.3

Simulation

0.2

0.1

0.0

480 490 500 510 520 530 540 550 560 570 580 590

Temperature [K] Fig. 5. Experiments and modeling of washcoated microreactor. Ethylene conversion and ethylene oxide yield as a function of temperature. The symbols represent the experimental data and the continuous lines the model results. Superficial gas velocity¼ 0.1 m/s.

Parameter

Optimal Estimate C.I. 90%

C.I. 95%

C.I. 99%

Ea,1/[J/mol] Ea,2/[J/mol] k1,ref /[m3/(molns)] k2,ref/[m3/(molns)] KE/[m3/mol] K O2 /[m3/mol]

4.59  104 5.70  104 4.98  10  1 3.73  100 4.10  100 3.16  10  1

1.68  102 1.73  102 3.62  10  1 2.71  100 1.57  100 2.43  10  1

2.26  102 2.33  102 4.89  10  1 3.66  100 2.12  100 3.28  10  1

1.39  102 1.43  102 3.00  10  1 2.25  100 1.30  100 2.01  10  1

oxide yield were obtained by increasing temperature. However, the selectivity suffered from the temperature increase, which indicated that higher temperatures favored the undesirable side reactions. The model results followed precisely the experimental observations. As can be seen, the reaction is strongly dependent on both fluid velocity and temperature. By applying the washcoat model to the experimental data, and the reaction rate expressions reported in Eqs. (11) and (12), a good agreement was obtained as can be appreciated from Figs. 4 and 5. The related kinetic parameters are given in Table 3. Table 3 shows that the adsorption constant for ethylene is roughly one order of magnitude larger than that for oxygen. The activation energy values indicate a rather strong temperature dependence of the rate constants. Numerical simulations revealed that the role of diffusion resistance inside the washcoat was minor in this case. The gPROMS software includes built-in computation of the statistics, where the confidence intervals are calculated automatically. The software calculates the Fischer information matrix which is equal to the inverse of the variance-covariance matrix. It evaluates the 90%, 95% and 99% F-values for this matrix (calculated using internal statistical functions). A confidence region means that if we repeat the experiments (which produces nearly the same measurements, but with slightly different observation values and, therefore, a different distribution of the measurement errors), and estimate the parameters out of the repeated experimental data, the values of the estimated parameters will lie in this confidence region with the announced probability. The confidence ellipsoid is only a linear approximation of the nonlinear confidence region and may not be very accurate for models which are highly non-linear in the parameters. 3.2. Silver plate microreactor

3. Modeling results and discussion 3.1. Washcoated microreactor The experimental data collected by (Carucci et al., 2010), Table 2, were interpreted by using both models (recalling that two kinds of microreactors were used). An effort was made to reveal both the temperature and the flow rate dependences of the ethylene conversion and the ethylene oxide yield. Figs. 4 and 5 present the experimental and parameter estimation results for the wash-coated reactor. In the earlier study (Carucci et al., 2010), results at only a single temperature were used for parameter estimation, whereas now the experimental data used for the parameter estimation covered a range of temperatures and superficial gas velocities. Higher superficial velocities of the fluid decreased both the yield and conversion almost linearly as expected. The average residence time became the shorter the higher the velocity of the fluid was and this was the very reason behind this kind of behavior. When temperature was increased, it enhanced the synthesis reaction but also the other reactions. Higher ethylene conversion and ethylene

Further experimental work was performed previously (Carucci et al., 2010) by using a pure silver plate microreactor. In this case, a higher epoxidation activity was found, and numerous experiments were performed to investigate the effect of the oxygen and ethylene concentrations on the reaction rate. The collected experimental data were been here elaborated with the second model, Eq. (7), considering that the reactions take place at the microplate surface only. The same reaction mechanisms was assumed as for the washcoated microreactor, so Eq. (11) was used without any changes. The trends of the experimental data, with the related simulations, are reported in Figs. 6–9. Table 4 shows the values of each input variables. For the silver plate reactor, the temperature and superficial gas velocity trends were qualitatively similar to the results obtained with the washcoated reactor for the very same reasons (higher temperatures favoring the side reactions, higher velocities reducing the residence time). The effects of feed concentrations (ethylene, and oxygen) were also included, since experimental data for these effects were available. The differing effects of feed

V. Russo et al. / Chemical Engineering Science 134 (2015) 563–571

0.30

0.14

X

X

0.12

Y

Y

0.25

Simulation

Simulation

XE [%], YEO [%]

0.10

XE [%], YEO [%]

569

0.08 0.06 0.04

0.20 0.15 0.10 0.05

0.02 0.00 0.08

0.10

0.12

0.14

0.00

0.16

Fluid velocity [m/s]

0

1

2

3

4

5

6

7

3

Oxygen initial concentration [mol/m ]

Fig. 6. Experiments and modeling of pure silver microreactor. Ethylene conversion and ethylene oxide yield as a function of the fluid velocity. The symbols represent the experimental data, the continuous lines the model simulations.

Fig. 9. Experiments and modeling of pure silver microreactor. Ethylene conversion and ethylene oxide yield as a function of the initial oxygen concentration. The symbols represent the experimental data, the continuous lines the simulations.

0.30 X

XE [%], YEO [%]

Table 4 Experimental data for silver plate reactor.

Y

0.25

Simulation Exp.

0.20 0.15 0.10 0.05 0.00 480

500

520

540

560

580

Temperature [K] Fig. 7. Experiments and modeling of pure silver microreactor. Ethylene conversion and ethylene oxide yield as a function of temperature. The symbols represent the experimental data, the continuous lines the simulations.

0.30

X Y

XE [%], YEO [%]

0.25

Simulation

0.20 0.15 0.10 0.05 0.00

0

1

2

3

4

5

6

7

3

Ethylene initial concentration [mol/m ] Fig. 8. Pure silver microreactor. Ethylene conversion and ethylene oxide yield as a function of the ethylene initial concentration. The symbols represent the experimental data, the continuous lines the simulations.

Cin(E) [mol/m3]

Silver plate 1 4.51 2 4.51 3 4.51 4 4.51 5 4.88 6 4.69 7 4.51 8 4.35 9 4.2 10 2.26 11 2.26 12 2.26 13 2.26 14 2.26 15 1.13 16 2.26 17 3.38 18 4.51 19 5.64 Silver plate, (Salmi et al., 2013) 1 18.4 2 18.4 3 18.4 4 18.4 5 18.4 6 18.4 7 18.4 8 4.6 9 7.36 10 9.2 11 13.8 12 18.4 13 23

Cin(O2) [mol/m3]

T/[K]

uF/[m/s]

2.26 2.26 2.26 2.26 2.24 2.34 2.26 2.17 2.1 1.13 2.26 3.38 4.51 5.64 2.26 2.26 2.26 2.26 2.26

533 533 533 533 493 513 533 533 573 533 533 533 533 533 533 533 533 533 533

0.09 0.11 0.13 0.15 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.09

4.6 7.36 9.2 13.8 18.4 23 27.6 18.4 18.4 18.4 18.4 18.4 18.4

523 523 523 523 523 523 523 523 523 523 523 523 523

0.148 0.148 0.142 0.142 0.142 0.142 0.142 0.142 0.142 0.142 0.142 0.142 0.142

concentrations could be well explained by the differing concentration dependences of LHHW (Langmuir-Hinshelwood-HougenWatson) kinetics for ethylene and oxygen, as the obtained model results illustrate. As for the washcoated microreactor, it is possible to observe a strong dependence of the reaction rates on both temperature and residence time. Moreover, by varying the feed concentration, there is a strong variation of the reaction rates: by fixing the oxygen and

570

V. Russo et al. / Chemical Engineering Science 134 (2015) 563–571

3.0

Table 5 Kinetic and adsorption parameters obtained for the pure silver microreactor. Kinetic constants (subscript ‘ref’) are calculated at 553 K. C.I. stands for the confidence interval and it is related to a certain percentage. Optimal Estimate C.I. 90%

C.I. 95%

C.I. 99%

Ea,1/[J/mol] Ea,2/[J/mol] k1,ref/[m3/(molns)] k2,ref/[m3/(molns)] KE/[m3/mol] KO2/[m3/mol]

3.10  104 5.44  104 1.54  10  4 3.14  10  4 3.59  10  1 2.69  10  1

7.15  101 9.05  101 2.55  10  7 6.55  10  7 7.35  10  4 6.62  10  4

1.10  102 1.20  102 3.39  10  7 8.70  10  7 9.76  10  4 8.79  10  4

5.32  101 7.57  101 2.14  10  7 5.48  10  7 6.15  10  4 5.53  10  4

Calculated data [%]

Parameter

2.5 2.0 1.5

± 20%

1.0 0.5

X Y

0.0

Table 6 Kinetic and adsorption parameters obtained for the pure silver microreactor. Kinetic constants (subscript ‘ref’) are calculated at 553 K. C.I. stands for confidence interval and it is related to a certain percentage. Parameter

Optimal Estimate C.I. 90%

k1,ref/[m3/(molns)] k2,ref/[m3/(molns)] KE/[m3/mol] KO2/[m3/mol]

3.63  10  3 1.08  10  2 9.98  10  1 5.99  10  1

C.I. 95%

1.81  10  5 2.17  10  5 5.34  10  5 6.40  10  5 2.52  10  3 3.02  10  3 1.61  10  3 1.93  10  3

C.I. 99% 2.89  10  5 8.54  10  5 4.03  10  3 2.58  10  3

increasing the ethylene concentration, Fig. 7, a decrease of both the ethylene conversion and ethylene oxide yield can be observed; on the contrary, by increasing the oxygen inlet concentration, there is an increase in the conversion and the ethylene oxide yield. These aspects are properly described by the developed model. The regressed parameters are reported in Table 5. In this case, the two adsorption constants are of the same order of magnitude. A good agreement can be obtained for this second case study too, as it can be appreciated from Figs. 6–9. Finally, the modeling of the experimental work (Salmi et al., 2013) was done using the same silver plate microreactor but an improved pretreatment of the catalyst. In this case, a higher activity was found and the experiments were performed at a high total pressure, higher global concentration, by varying the inlet concentration only. The data have been here interpreted with the same reaction rate law Eq. (11). A new parameter estimation study was required to take into account the pretreatment effect on the catalyst activity and selectivity. The new kinetic parameters are given in Table 6 and the overall agreement is shown in the parity plot, Fig. 10. The higher activity is further reflected in the higher kinetic constants: one order of magnitude for the partial oxidation and two orders of magnitudes for the total oxidation of ethylene. Also in this case the agreement between the model and the experimental data can be considered to be very satisfactory.

0.0

0.5

1.0

1.5

2.0

2.5

3.0

Experimental data [%] Fig. 10. Experiments and modeling of pure silver microreactor. Ethylene conversion and ethylene oxide yield parity plot.

modeling approaches are generic and applicable more broadly to other heterogeneously catalyzed gas-phase reactions. The obtained modeling results were well in-line within the experimental observations. This implied on one hand that the experimental arrangement was capable to generate explainable data and, on the other hand that the proposed physical models (based on mass balances, LHHW kinetics, simple flow description by dispersion models and semi-empiric equations for pressure drop) were capable to punctually describe the concentration and temperature dependences of the reaction rates associated with the EO synthesis. Microreactors proved to be well suited devices for laboratory studies. By ruling out undesirable heat and mass transfer related effects, modeling of their operation became straightforward. The results are also scalable, by “numbering up” for single phase applications. The “numbering up” of microscale reactors will work well if the pressure drop of the microchannel is high compared to the one of feed entrance. The kinetic results are intrinsic. For bulk products, larger channel diameters are to be used to get desired production volumes. Then, the models described here would demand also the energy balances to be included and simultaneously solved. The addition of energy balances would not be a hard challenge, because from the mathematical point of view, they are analogous to mass balances and consequently can be solved with exactly the same method, the Numerical Method of Lines.

Nomenclature

4. Conclusions

ai bi C i;F C i;F IN

Two physical microreactor models, one for a washcoated silver catalyzed reactor and the other one for a silver plate microreactor, were developed and used for ethylene epoxidation and total oxidation. Both models were based on the solution of dynamic mass balances for the gas-phase components. The fluid dynamics in the microreactor elements was expressed in simple terms to develop rapid program codes for parameter estimation. Complete data elaboration using the experimental material previously generated in our laboratory for a precise description of the reaction kinetics was done. The results proved that the models were capable to explain the experimental data very well. The used

C i;S CE C O2 C i;0 C i;F DA;B Def f ;i Dr;F Dz;F Ea,j i

diffusivity parameter, [–] diffusivity parameter, [–] concentration of component i in gas, [mol/m3] concentration of component i in inlet gas stream, [mol/ m3 ] local concentration of component i in solid, [mol/m3] concentration of ethylene, [mol/m3] concentration of oxygen, [mol/m3] initial concentration of component i, [mol/m3] concentration of component i in fluid, [mol/m3] diffusivity of compound A in B , [m2/s] effective diffusivity of component I [m2/s] radial dispersion coefficient, [m2/s] axial dispersion coefficient, [m2/s] activation energy of reaction j, [J/mol] component i

V. Russo et al. / Chemical Engineering Science 134 (2015) 563–571

j k k1 ; k2 kj kj;ref K E , K O2 L MA MB M AB P Pe R RW r rp Re Rg rj s t T Tref Tn Tb uF

νi;j

V_ Vc XE YEO z α, β

εp εAB υi,j

ρ σ AB ΩD

reaction j parameter (Eq. 8) reaction rate constants, in convenient units to produce rate as, [mol/m3s] kinetic constant of reaction j, [(m3/mol)2 s  1)] kinetic constant at reference temperature, [(m3/ mol)2 s  1)] adsorption constant for ethylene and oxygen, [(m3)/mol] reactor length, [m] molecular mass solute [g/mol] molecular mass of solvent, [g/mol] mean molecular mass, [g/mol] pressure, [Pa] peclet number, [–] inner reactor radius, [m] reactor radius, [m] radial location, [m] radial location within the washcoat layer, [m] reynolds number, [–] ideal gas constant, [J/(K mol)] reaction rate expression of reaction j, [mol/(m3∙s)] shape factor [ ] time, [s] temperature, [K] reference temperature, [363 K] boiling point, [K] boiling point, [K] superficial velocity of fluid (gas), [m/s] stoichiometric coefficient, [–] volumetric flow rate, [m3/s] critical molar volume of component i, [cm3/mol] conversion of ethylene [%] yield [%] axial location, [m] parameters in rate equation, [–] porosity, [–] parameter stoichiometric coefficient of component i in reaction j, [–] density, [kg/m3] parameter, [–] parameter, [–]

Acknowledgments Academy of Finland and the University of Naples Federico II are acknowledged for the financial support. Process System Enterprise

571

(PSE) is acknowledged for providing licenses to enable this study. The work is a part of the activities of Johan Gadolin Process Chemistry Center (PCC), a Center of Excellence financed by Åbo Akademi University.

References Arve, K., 2005. Catalytic Diesel Exhaust Aftertreatment: From Reaction Mechanism to Reactor Design Doctoral thesis. Åbo Akademi, Department of Chemical Engineering, Åbo. Belfiore, L., 2003. Transport Phenomena for Chemical Reactor Design. John Wiley & Sons, Hoboken, New Jersey. Borman, P., Westerterp, K., 1995. An experimental study of the kinetics of the selective oxidation of ethene over a silver on α-Al2O3 catalyst. Ind. Eng. Chem. Res. 34, 49–58. Carucci, J., Halonen, V., Eränen, K., Wärnå, J., Ojala, S., Huuhtanen, M., Keiski, R., Salmi, T., 2010. Ethylene oxide formation in a microreactor: from qualitative kinetics to detailed modeling. Ind. Eng. Chem. Res. 49, 10897–10907. Carucci, J., Eränen, K., Salmi, T., Murzin, D., 2012. Gas-phase microreactors as a powerful tool for kinetic investigations. Russ. J. Gen. Chem. 82, 2034–2059. Deutchmann, O., Schwiedernoch, R., Maier, L., Chatterjee, D., 2001. Natural gas conversion in monolithic catalysts: interaction of chemical reactions and transport phenomena. Stud. Surf. Sci. Catal. 136, 251–258. Ghazali, S., Park, D., Gasu, G., 1983. Kinetics of ethylene epoxidation on a silver catalyst. Appl. Catal. 6, 195. Hayes, R., Awdry, S., Kolaczkowski, S., 1999. Catalytic combustion of methane in a monolith washcoat: effect of water inhibition on the effectiveness factor. Can. J. Chem. Eng. 77, 688–697. Hessel, V., Hardt, S., Löwe, H., 2004. Chemical Micro Process Engineering. WileyVCH, Weinheim. Keil, F., 2007. Modeling of Process Intensification. Wiley-VCH, Weinheim. Kursawe, A., 2009. Partial Oxidation of Ethylene to Ethylene Oxide in Microchannel Reactors Doctoral thesis. Technische Universität Chemnitz, Fakultät für Naturwissenschaften, Chemnitz. Lafarga, D., Al-Joaied, M., Bondy, C., Varma, A., 2000. Ethylene epoxidation on Ag– Cs/α-Al2O3 catalyst: experimental results and strategy for kinetic parameter estimation. Ind. Eng. Chem. Res. 39, 2148–2156. Pattegar, A., Kothare, M., 2004. A Microreactor for hydrogen production in micro fuel cell applications. J. Microelectromech. Syst. 13, 7–11. Petrov, L., Eliyas, A., Maximov, C., Shopov, D., 1988. Ethylene oxide oxidation over a supported silver catalyst. Appl. Catal. 41, 23–38. Poling, B., Prausnitz, J.M., O’Connell, B., 2001. The Properties of Gases and Liquids, 5th Ed. McGraw-Hill, New York. Reid, R., Prausnitz, J.M., Poling, B., 1988. The Properties of Gases and Liquids, 4th Ed. McGraw-Hill, New York. Salmi, T., Carucci, J., Roche, M., Eränen, K., Wärnå, J., Murzin, D., 2013. Microreactors as tools in kinetic investigations: ethylene oxide formation on silver catalyst. Chem. Eng. Sci. 87, 306–314. Schmidt, S.A., Kumar, N., Reinsdorf, A., Eränen, K., Wärnå, J., Murzin, D., Salmi, T., 2013. Methyl chloride synthesis on Al2O3 in a microstructured reactor – thermodynamics, kinetics and mass transfer. Chem. Eng. Sci. 95, 232–245. Toppinen, S., 1997. Liquid Phase Hydrogenation of Some Aromatic Compounds Doctoral thesis. Åbo Akademi, Department of Chemical Engineering, Åbo. Stutz, M., Poulikakos, D., 2008. Optimum washcoat thickness of a monolith reactor for syngas production by partial oxidation of methane. Chem. Eng. Sci. 63, 1761–1770. Wilke, C.Y., Lee, C.Y., 1955. Estimation of diffusion coefficients for gases and vapors. Ind. Eng. Chem. Res. 47, 1253–1257. Özbek, M.O., van Santen, R.A., 2013. The mechanism of ethylene epoxidation catalysis. Catal. Lett. 143, 131–141.