ARTICLE pubs.acs.org/IECR
Reaction and Deactivation Rates of Methane Catalytic Cracking over Nickel Ashraf Amin, William Epling,* and Eric Croiset Department of Chemical Engineering, University of Waterloo, 200 University Avenue West, Waterloo, Ontario N2L 3G1 Canada
bS Supporting Information ABSTRACT: The kinetics of methane catalytic cracking over nickel supported on porous and nonporous aluminas was modeled using a separable kinetics approach in order to develop initial rate and activity decay equations. The model parameters were estimated using a set of experiments conducted in an electrobalance. The experimental work covered the 500650 °C temperature range, using pure methane, as well as different partial pressures of CH4/N2 and CH4/H2 mixtures at atmospheric pressure. The model results showed a good match with the experimental data, and the estimated kinetic parameters agreed well with those reported in the literature. The morphology of the support affected the initial reaction rate and catalyst deactivation. The methane cracking activation energy was estimated to be 88 and 75 kJ/mol for the porous and nonporous catalysts, respectively. The activation energy for the encapsulating carbon formation was estimated to be 147 and 149 kJ/mol for the porous and nonporous catalysts, respectively. The deactivation reaction was found to be half-order in surface carbon. The model was expanded to include cracking/ regeneration cycles. The model showed good agreement with the experimental data at different experimental conditions and up to 39 cycles. Cracking/regeneration cycles suggest that the porous catalyst can be used for conducting continuous methane cracking.
1. INTRODUCTION Currently, CO-free hydrogen production is of significant commercial and therefore research interest. CO-free hydrogen has wide application in industry and fuel cell applications. Methane catalytic cracking (MCC) is a potential process for the production of CO-free hydrogen. MCC has the advantage over the conventional steam reforming process by sequestering the carbon, which eliminates any cross contamination of the products with COx and reduces the evolution of greenhouse gases.15 Methane decomposition without a catalyst requires very high temperatures, >1300 °C. An active catalyst is required to obtain high methane conversion at reasonable temperatures: 500700 °C for nickel-based catalysts, 700950 °C for iron-based catalysts, 850 950 °C for carbon-based catalysts, and 7001000 °C for Co, Pd, Pt, Cr, Ru, Mo, and W catalysts.3,68 Besides its being active in a relatively low temperature range compared to other catalysts, nickel is studied due to its high activity in methane cracking and has the highest carbon capacity among the iron group components.9 Also, the support material and its textural properties affect methane conversion. Minimal interaction between the catalyst and support is important to ensure a highly active catalyst for the reaction.9,10 The textural properties of the support affect the reactant and product species diffusion to the active sites, which affects the methane catalytic cracking rate.11 Takenaka et al.11 reported that nickel supported on nonporous silica shows the highest activity among different surface area and textural structure silicas. A major problem encountered with methane cracking is rapid catalyst deactivation. Catalyst deactivation occurs as carbon accumulates on the catalyst, with carbon of course being one of the reaction products.1,7,12 The deposited carbon can diffuse through and around nickel to form carbon filaments. However, some of the carbon accumulates on the external surface of the r 2011 American Chemical Society
nickel sites, gradually blocking the active nickel surface until it is completely encapsulated, causing complete deactivation.5,6 According to this deactivation scenario, the reaction rate starts with a maximum value and then decreases gradually until it reaches zero. The filamentous carbon formation does not affect the catalyst activity, since the carbon filament is formed by diffusion of carbon through nickel keeping the nickel on the tip of the filament, and therefore does not affect either the nickel surface area subjected to the reaction or the nickel activity. On the other hand, the encapsulating carbon is formed as a layer on the catalyst surface preventing any contact between the catalyst and the reactants, and it is therefore expected that catalyst deactivation is a function of the encapsulating carbon concentration. The deposited carbon is therefore an indication of reaction extent, via carbon filament formation, but it is also a measure of deactivation as encapsulating carbon forms. Any change in reaction conditions that causes an increase in the reaction rate will lead to an increase in the deactivation rate and vice versa. To achieve maximum hydrogen production and to extend the catalyst lifetime, a balance should be achieved between reaction rate and deactivation rate. The optimum operating conditions can be explored using a well-developed reaction rate model that takes into account the deactivation mechanism and the effect of the operating parameters. Since the first methane catalytic cracking mechanism was proposed by Grabke,13,14 different models have been proposed.1319 The deactivation behavior has been described through equations derived through mathematical fits in different studies.1519 Chen Received: June 3, 2011 Accepted: September 19, 2011 Revised: September 12, 2011 Published: September 19, 2011 12460
dx.doi.org/10.1021/ie201194z | Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research Table 1. Comparison of Literature and Thermodynamic Values for KP temperature, K
thermodynamics
Zavarukhin and Kuvshinov16
773
0.50
0.42
823 873
1.08 2.33
0.83 2.05
923
4.47
3.52
et al.20 proposed two reaction steps to explain the dependency of the catalyst activity on encapsulating carbon formation. Hazra et al.4 used Chen’s steps to develop a kinetic model for methane catalytic cracking using nonseparable kinetics. In this model, deactivation was evaluated based on the reaction mechanism, quantifying the rate of deactivation through the separable kinetic form, which helped provide an explanation for the deactivation phenomenon. The reaction equilibrium constant, KP, values used in the model were calculated from thermodynamics and were within 30% of the values used by Zavarukhin and Kuvshinov16 as shown in Table 1. To conduct methane cracking in a continuous process, the catalyst is expected to go through cracking/regeneration cycles. Modeling catalyst behavior during cyclic operation is critical in understanding/predicting the catalytic performance over long periods of operation. To develop such a model, two nickel catalysts, supported on porous γ-alumina (porous catalyst) and nonporous α-alumina (nonporous catalyst), were prepared in this study. Experiments were conducted using an electrobalance, with effects of particle size, flow rate, temperature, and methane partial pressure with nitrogen as the inert gas or with different partial pressures of hydrogen on the main reaction evaluated. The experimental study was used to investigate the model parameters for each catalyst individually. The model parameters considered here are the pre-exponential factor and activation energy of the rate and adsorption constants when expressed in an Arrhenius form.
2. EXPERIMENTAL METHODS The catalyst was prepared by wet impregnation of alumina using an aqueous solution of Ni(NO3)2 3 6H2O (99.99%, Alfa Aesar). Two different supports were used: (1) a porous γalumina support [99.97% metal basis, 3 μm APS powder, Alfa Aesar] and (2) a nonporous α-alumina support [99.99% metal basis, 0.92.2 μm APS powder, Alfa Aesar]. The catalyst preparation method consisted of drying the support overnight at 150 °C, preparing a solution of nickel nitrate hexahydrate in deionized water, and then adding the alumina to the solution (the nickel and alumina amounts used were in a 10/100 weight ratio). The slurry was then stirred for 3 h at 80 °C for the porous catalyst and in an oil bath at 100 °C for the nonporous catalyst. The activity experiments were performed in a Cahn TG 151 electrobalance, manufactured by Thermo Cahn. The electrobalance is composed of a quartz tube that can withstand high temperature and pressure applications (up to 1000 °C at 60 bar). The sample holder is placed in the middle of a quartz tube hanging down from a wire attached to a balance that can record any change in the sample weight, with 1 μg resolution. The temperature is controlled by a K-type thermocouple placed just below the sample holder. The catalysts were calcined in air at 600 °C inside the thermal gravimetric analyzer (TGA) for 30 min. The reduction step was then also carried out in the TGA, at the reaction temperature, for 30 min using a 10 vol % H2 and 90
ARTICLE
vol % N2 mixture.5 The weight of the fresh catalyst was kept around 10 mg for each experiment to avoid spatial limitations. Also, a flat bottom sample holder was used to ensure that each catalyst particle behaved independently, as shown in Figure 1. All of the gases used in this study were 99.99% pure, supplied by Praxair. Except when the flow rate effect and particle diameter effect were studied, 120 mL/min was used as a default flow rate and 725 μm particles were used as the default particle size. The activity experiments included experiments at different partial pressures of CH4/ N2 and CH4/H2 mixtures. The BET surface areas were measured with a Micromeritics Gemini III 2375 surface area analyzer. The nickel dispersion was measured using a Hiden Catlab Microreactor with a QIC-20 mass spectrometer. A 5% H2 in He mixture was used, with chemisorption measurements obtained at 25 °C, after the sample had been pretreated at 550 °C using 10% H2 in N2 mixture, and then cooled in pure He to 25 °C.
3. RESULTS AND DISCUSSION 3.1. Typical Experimental Results. Figure 2 shows typical experimental results from the electrobalance, i.e., catalyst weight change, in grams of carbon per gram of nickel (gC/gNi), as a function of time. By differentiating the catalyst weight change with time, the rate of methane decomposition was calculated in terms of mmolCH4/ gNi/min. The maximum rate was used as the initial catalyst rate and the time was set at zero when the electrobalance control valve was switched from nitrogen (inert gas) to the reactant mixture. As seen in Figure 2, the maximum rate occurs some period of time after switching the valve to the reactant mixture. The time to reach the maximum rate was further investigated by simulating the residence time distribution (RTD) of the gases inside the electrobalance and around the catalyst particles. 3.2. RTD. The RTD of the electrobalance reactor vessel was studied using COMSOL. Two built-in models, the NavierStokes and the Convection and Diffusion models, were used to fully simulate momentum and mass transfer inside the reactor tube and around the catalyst particles. Figure 3 shows the simulation results for pure methane at 550 °C and 1 atm at different flow rates. Figure 3 shows the change in bulk gas methane concentration surrounding the catalyst particle over time. The results indicate that the maximum rate occurs around the same time methane reaches its maximum concentration in the gas bulk around the catalyst particle. These RTD results show that attributing the maximum rate as the initial rate is acceptable. 3.3. Controlling Regime. Before kinetic data were collected, the reaction was studied at different flow rates and with different particle diameters to examine the regime that controls the reaction, i.e., if the reaction is controlled by internal diffusion, external diffusion, or kinetics.21 The particle diameter was varied between 300 and 1000 μm, and the gas flow rate was varied between 72 and 240 mL/min. The BET surface areas for each catalyst with different particle sizes are listed in Table 2. Figure 4 shows the reaction rate versus time for different porous catalyst particle diameters. These data indicate that the initial rate (represented as the maximum rate) is not affected by catalyst particle size. Figure 4 also shows that there is no clear effect of particle size on the time required for full deactivation. Experiments with the nonporous catalyst resulted in trends similar to those shown in Figure 4. The effect of flow rate on the initial rate was studied at different temperatures, and the results are illustrated for the porous and nonporous catalysts in Figures 5 and 6, respectively. The flow 12461
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research
ARTICLE
Figure 1. Catalyst sample before reaction and after complete deactivation.
Table 2. BET Surface Areas for Porous and Nonporous, Reduced Catalysts surface area, m2/g
Figure 2. Typical experimental results at 550 °C and 1 atm for 100% methane, 120 mL/min, using Ni/γ-Al2O3.
particle size, μm
porous
nonporous
300
70.4
18.4
725
67
17.9
1000
59
15.9
H2, the initial rate can be expressed as initial reaction rate ¼ R0 ¼ f1 ðPCH4 , PH2 , TÞ
ð1Þ
Once catalyst deactivation occurs, the instantaneous reaction rate is a function of species partial pressures, temperature, and time (time accounts for the activity decay over time as a result of active site deactivation): instantaneous reaction rate ¼ R1 ¼ f2 ðPCH4 , PH2 , T, tÞ
ð2Þ
It is conventional to define another variable called “activity”, which relates the rate at any time (t) to the initial rate, R0. The activity can be expressed by the following equation: a¼ Figure 3. Simulation of the RTD in the case of pure methane inside the electrobalance at different inlet flow rates at 550 °C and 1 atm for 100% methane. The concentration corresponds to the methane concentration surrounding the catalyst particle.
rate has no effect on the initial reaction rate at flow rates at or above 120 mL/min for both the porous and nonporous catalysts. As the change in initial rate was insignificant above 120 mL/min, this was chosen for initial rate tests. Since, under the operating conditions used in the present study, the reaction rate is independent of particle diameter and gas flow rate (at least for flow rates above 120 mL/min), it is concluded that the reaction is kinetically controlled. 3.4. Reaction Rate. The main problem in modeling the kinetics of methane catalytic cracking is taking into account rapid catalyst deactivation. Usually, the rate of a chemical reaction is a function of chemical species concentrations (or partial pressures) and temperature. Since the gaseous species involved in methane cracking are CH4 and
R1 ðtÞ R0
ð3Þ
A practical approach toward studying deactivating reactions is to use the separable kinetics technique, and it was used in the present study. Using this method, the reaction rate is separated into two terms: the first term accounts for reaction kinetics effect and is time independent; the second term accounts for the deactivation and is time dependent. A complete description for the separable kinetic approach is available in Butt and Petersen.22 3.4.1. Reaction Mechanism. Methane cracking occurs on the catalyst surface in a multistep process. The reaction starts with the adsorption of methane, followed by a series of dehydrogenation reactions until it ends as adsorbed carbon on the top of the nickel surface.17,23 The adsorbed carbon segregates/dissolves on/through the nickel surface based on the reacting mixture concentration and the active sites available for the reaction.18,23 The adsorbed carbon forms a solution with the nickel and diffuses through the nickel from the top to the bottom (initially in 12462
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research
ARTICLE
Figure 4. Particle diameter effect on the reaction rate for the porous catalyst at 550 °C, 100% methane, 120 mL/min, and 1 atm.
Figure 6. Flow rate effect on the initial reaction rate for the nonporous catalyst using 100% methane, 725 μm particles, and 1 atm.
step 6. 2H 3 I S 2H2 þ 2I
1=KH
ð9Þ
step 7. C 3 I S CNi, f þ I
1=KC
ð10Þ
carbon diffusion process ðcarbon filament formationÞ : step 8: Figure 5. Flow rate effect on the initial reaction rate for the porous catalyst using 100% methane, 725 μm particles, and 1 atm.
contact with the support), where carbon forms carbon filaments. The carbon diffusion process occurs simultaneously with encapsulating carbon formation. The encapsulating carbon is adsorbed carbon on the top of the nickel surface but remains and blocks the nickel site from the reaction and causes deactivation.4,20,23 Based on this scheme, methane cracking can be divided into the following three processes: methane cracking process: step 1. CH4 þ I S CH4 3 I
ð4Þ
K1
step 2. CH4 3 I þ I S CH3 3 I þ H 3 I
rds=K2
ð5Þ
K3
ð6Þ
step 3. CH3 3 I þ I S CH2 3 I þ H 3 I
CNi;r f Cfil kp ðcarbon filament formationÞ deactivation (encapsulating carbon formation): step 10. nC 3 I f nCP 3 I
kd
ð13Þ
In this scheme, rds is the rate determining step, n is the deactivation order, I represents a vacant site, and X 3 I represents adsorbed X species. CNi,f, CNi,r, and CP 3 I are carbon adsorbed on the nickel top, carbon adsorbed on the nickel bottom, and encapsulating carbon adsorbed on the nickel top, respectively. 3.4.2. Initial Rate Equation. The reaction step represented by eq 5, the second step in the cracking sequence, is assumed to be the rate-limiting step.2325 Therefore, it is assumed that all methane cracking reaction steps, eqs 410, are in equilibrium except step 2 (eq 5), and that the catalyst is uniformly distributed over the support. Based on the previous assumptions, the following equation was developed to represent the methane cracking rate: " # dθCH3 3 I θCH3 3 I θH 3 I ¼ C0 k2 θCH4 3 I θI ð14Þ dt K2 where θj is the fraction of sites occupied by species j, and is defined as
step 4. CH2 3 I þ I S CH 3 I þ H 3 I
K4
ð7Þ
θj ¼
step 5. CH 3 I þ I S C 3 I þ H 3 I
ð12Þ
K5
ð8Þ
Cj C0
ð15Þ
C0 is the total site concentration (mol/cm2) and Cj is the concentration of sites occupied by species j (mol/cm2). 12463
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research
ARTICLE
Replacing θCH4 3 I, θI, θCH3 3 I, and θH 3 I by their expressions from the reaction steps at equilibrium, eq 14 becomes " # dθCH3 3 I C0 k2 K1 θC 3 I 2 KH 2 KC CCNi, f 2 ¼ rCH4 ¼ PCH4 PH dt K1 K2 K3 K4 K5 2 kC 2 CCNi, f 2
Table 3. Parameter Estimation for the Porous Catalyst constant
KH2
where rCH4 is the overall methane cracking rate, in mmolCH4/gNi/ min, and is equal to the rate of the rate-limiting step (dθCH3 3 I/dt). Using the overall site balance to replace θC 3 I in eq 16
KCH4
ð17Þ
1
θC 3 I ¼
1 þ þ
1 KCH4 PCH4 KH2 3=2 PH2 3=2 KH2 P H2 þ þ þ KC CCNi, f KC CCNi, f K3 K4 K5 K4 K5
KH2 1=2 PH2 1=2 KH2 1=2 PH2 1=2 þ K5 KC CCNi, f
θC 3 I ¼
1 1 þ KC CCNi, f
88 ( 6
AH2 (atm3/2)
2 108 ( 5 109
7.8
ΔHH2 (kJ/mol)
144 ( 31
9.18
28
ACH4 (atm1)
3.75 105 ( 4.68 106 16
ΔHCH4 (kJ/mol)
56 ( 8
13.5
a
2834
number of samples
30
F tabulated (0.05,2,28)
4.2
R2
0.995
t critical value (0.05/2,2).
where k is the specific reaction rate for methane cracking, in mmolCH4/gNi/min/atm. C0 k2 KCH4
k¼ 2
kC CCNi, sat
ð18Þ
If we consider that the surface concentrations of CH2, CH, and H, θCH2 3 S, θCH 3 S, and θH 3 S, are negligible, then eq 18 becomes
39
E (kJ/mol)
F observed
By replacing the value of each θX from the individual reaction mechanism rate equations and by assuming that all methane cracking reaction steps, eqs 410, are in equilibrium except step 2 (eq 5), the following can be derived:
t value
4.3a
θI þ θCH4 3 I þ θCH3 3 I þ θCH2 3 I þ θCH 3 I þ θC 3 I þ θH 3 I ¼ 1
estimate
A (mmolCH4/gNi/min/atm) 4.64 107 ( 2.36 106
k
ð16Þ
parameter
2
1 1 þ KC CCNi, sat
!2
KCH4 and KH2 are the overall adsorption constants for methane and hydrogen, in atm1 and atm3/2, respectively. KCH4
KCH4 ¼
1 KCH4 PCH4 KH2 3=2 PH2 3=2 þ þ KC CCNi, f K3 K4 K5
KC CCNi, sat
!
1 1 þ KC CCNi, sat
ð19Þ Replacing θC 3 I in eq 16 by its value from eq 19 dθCH3 3 I dt
C0 k2 KCH4
¼ 2
kC CCNi, f
" PCH4
2
1 KCH4 PCH4 KH2 3=2 PH2 3=2 1 þ þ þ KC CCNi, f KC CCNi, f K3 K4 K5 2
KH2 KC CCNi, f PH 2 KCH4 K2 K3 K4 K5 2
K3 K4 K5
!2
!
1 1 þ KC CCNi, sat
KP is the equilibrium constant of the overall reaction with units of atm and its values, from thermodynamic calculations, are listed in Table 1.
#
KP ¼
ð20Þ Snoeck et al.23 observed that the carbon concentration is uniform over the nickel particle. The authors23 assumed that CCNi, f ≈ CCNi, r ≈ CCNi, sat With this assumption, CCNi,f can be replaced with CCNi,sat. By dividing the denominator and numerator by !2 1 1 þ KC CCNi, sat and incorporating it into the other constants: ! P H2 2 k PCH4 KP rCH4 ¼ ð1 þ KH2 PH2 3=2 þ KCH4 PCH4 Þ2
KH2 3=2
K H2 ¼
KH2 2 KC CCNi, sat KCH4 K2 K3 K4 K5
3.4.3. Model Discrimination. The model parameters were evaluated statistically using a three-step method. In the first step, the kinetic and adsorption constants were estimated at four temperatures (500, 550, 600, and 650 °C). In the second step, from the estimated values of the constants in the first step at different temperatures, the pre-exponential factor and activation energy were determined (by plotting ln k vs 1/T). In the third step, the estimated parameters in the second step were used as initial values for an overall fitting exercise that took into account all experimental results for a particular catalyst. The parameter estimation was achieved by minimizing the sum of squares function:26,27 f ¼
ð21Þ 12464
n
ðrexp rcalc Þ2 ∑ p¼1
ð22Þ
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research
ARTICLE
Table 4. Parameter Estimation for the Nonporous Catalysta constant k KH2
parameter
A (mmolCH4/gNi/min/atm) 1 107 ( 13 105
t value 15
E (kJ/mol)
75 ( 15
AH2 (atm3/2)
7.6 108 ( 2.4 108 6.3
ΔHH2 (kJ/mol) KCH4
estimate
133 ( 36
10 7.3
ACH4 (atm1)
7.8 106 ( 4.4 105 0.36
ΔHCH4 (kJ/mol)
69 ( 34
4 4.3b
F observed
914
number of samples
31
F tabulated (0.05,2,29)
4.2
R2
0.984
A, E, and ΔH are the pre-exponential factor for the Arrhenius equation, the activation energy in J/mol, and the heat of adsorption in J/mol, respectively. b t critical value (0.05/2,2).
a
rexp, rcalc, and n are the reaction rate obtained experimentally, the rate calculated from the model, and the number of experiments, respectively. The results of the parameter estimation and analysis of variance are given in Tables 3 and 4 for the porous and nonporous catalysts, respectively. In the statistical analysis, three different tests were used: the Fischer (F) test was used to check the adequacy of the model to represent the data, the correlation coefficient was used to check if the model output fits the data well, and finally the t-distribution test was used to examine if each term in the model contributes significantly to the model. The critical F tabulated value is 4.2 according to Montgomery and Runger.27 If the calculated F test value is greater than 4.2, it implies that the model is adequate to represent the data, and the higher the F test value the better the model. From Table 3, the calculated F values are 2834 and 914 for the porous and nonporous catalysts, indicating that the model is adequate to represent the kinetic data for both catalysts. It also indicates that the model is more adequate for the porous catalyst data than for the nonporous catalyst data. Also, the 0.995 and 0.984 correlation coefficient values for the porous and nonporous catalysts indicate a good fit. Parity plots are provided in the Supporting Information (Figures S1 and S2) and show that most of the data points are distributed symmetrically around the line with the majority of the points falling on the diagonal. The parity plots represent the accuracy of the prediction of the model for the kinetic data which agree with the results of the F test and R2 test. Before conducting the t test (Student t test), the model coefficients were examined to check that all are positive, that the rate constant increased with temperature, and that the adsorption constants decreased with increasing temperature. For a 95% confidence interval, the t value must be greater than 4.3 in order to accept that the term contributes significantly to the model. The true mean of the parameters is estimated within an interval of 95% confidence. For the porous catalyst (see Table 3), all model parameters contribute significantly to the model. For the nonporous catalyst (see Table 4), except for the methane adsorption parameters, all other parameters show significant contributions to the model. This can be explained by the poor fit of the methane adsorption constant to the Arrhenius equation. However, if the methane
adsorption term is removed from the model for the nonporous catalyst, the F value becomes 220 and R2 is reduced to 0.93, and the sum of squares error is increased by 5-fold, which suggests the necessity of keeping the methane adsorption term. 3.4.4. Comparison with Literature Data. The estimated parameters were compared to literature values. The measured activation energy is 88 ( 6 kJ/mol for the porous catalyst and 75 ( 15 kJ/mol for the nonporous catalyst. The lower activation energy is associated with the higher activity observed with the nonporous catalyst; however, the calculated error suggests that the activation energies are too similar to discriminate or conclude that they differ. In the literature, reported activation energies over nickel catalysts are 88,17 90,18 97,19 96,28 and 75 kJ/mol.29 The heats of hydrogen adsorption are estimated as 144 ( 31 and 133 ( 36 kJ/mol for the porous and nonporous catalysts, respectively. For nickel supported on alumina, Bartholomew30 reported a value of 125 kJ/mol (wet-impregnation preparation method), and Germer and MacRae31 reported 117 kJ/mol (using nickel rods). Bartholomew30 mentioned that the heat of adsorption of hydrogen on nickel is affected by the catalyst’s preparation method and the catalyst’s structure. The values in the literature are in the same range but are slightly lower than the estimated value in this work. The heats of methane adsorption are estimated at 56 ( 8 and 69 ( 34 kJ/mol for the porous and nonporous catalysts, respectively. Beebe et al.32 found that the heat of adsorption of methane over the nickel (110) surface, which is the most active surface for methane cracking, was 55.6 kJ/mol, very close to the estimated value for the porous catalyst. The higher value observed with the nonporous catalyst is attributed to the lower specific surface area available for the reaction, which increases competition over active sites. 3.5. Deactivation. Once the initial rate, as calculated from eq 21, is achieved, deactivation starts and the rate decreases with time until the catalyst is completely deactivated. The actual rate at any time can be expressed using the activity term defined in eq 3. The initial rate is a function of methane and hydrogen partial pressures, and temperature, and the rate decreases due to encapsulation of the active sites by carbon. Therefore, the activity term should be a function of the fraction of encapsulated sites. The overall site balance before deactivation is expressed by the following equation: θI þ θCH4 3 I þ θCH3 3 I þ θCH2 3 I þ θCH 3 I þ θC 3 I þ θH 3 I ¼ 1
ð23Þ When taking into consideration encapsulated sites, whose fraction is expressed as θCP 3 I, the overall site balance is expressed as follows: θI þ θCH4 3 I þ θCH3 3 I þ θCH2 3 I þ θCH 3 I þ θC 3 I þ θ H 3 I þ θ CP 3 I ¼ 1
ð24Þ
Replacing θI, θCH4 3 I, θCH3 3 I, θCH2 3 I, θCH 3 I, and θH 3 I in eq 24 by their expressions from the reaction steps at equilibrium, the following relation between θC 3 I and θCP 3 I is obtained: θC 3 I ¼ 12465
1 θCP 3 I k0
ð25Þ
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research
ARTICLE
From step 10 in the reaction mechanism (eq 13), we have
where k0 ¼ 1 þ þ
dθCP 3 I
1 KCH4 PCH4 KH2 3=2 PH2 3=2 þ þ KC CCNi, f KC CCNi, f K3 K4 K5
K H2 P H2 KH2 1=2 PH2 1=2 KH2 1=2 PH2 1=2 þ þ K4 K5 K5 KC CCNi, f
dt ð26Þ
rCH4
!2 " PCH4
dα kd ¼ 0n αn dt k
KH2 2 KC CCNi, f PH 2 KCH4 K2 K3 K4 K5 2
#
ð27Þ A characteristic of the deactivation is the reduction in the total number of active sites. Therefore, it is important to estimate the total number of active sites at any time during the deactivation process. To do so, we introduce the fraction, α, of active sites not deactivated relative to the initial total number of active sites:22 α¼
Nt N0
ð29Þ
Ultimately, the goal is to derive an expression of the activity as a function of time. The activity is related to the fraction α, and thus to 1 θCP 3 I. From eqs 26 and 27, it can be seen that the cracking rate is a complex function of temperature, PCH4, PH2, CCNi,f, and θCP 3 I. The variables CCNi,f and θCP 3 I are related to each other, but they are difficult to estimate. Because α is equal to 1 θCP 3 I, it will be assumed that the cracking rate is proportional to (1 θCP 3 I)P, where P is a fitting parameter: ðrCH4 Þt ¼ ðrCH4 Þt ¼ 0 ð1 θCP 3 I ÞP
ð30Þ
Note that P is very likely different from 2, despite the dependence of the cracking rate on (1 θCP 3 I)2 as seen in eq 27. Indeed, CCNi,f is a function of θCP 3 I. Equation 30 satisfies the fact that, at no deactivation, θCP 3 I is zero and the reaction rate is equal to the initial rate. After deactivation occurs, the reaction rate decreases as θCP 3 I increases until the rate reaches zero, when all the active sites are encapsulated. It should be clear that eq 30 is equivalent to defining the activity, a, as a ¼ ð1 θCP 3 I ÞP
ð31Þ
Therefore, the relation between the ratio α and the activity a is simply a ¼ αP
ð32Þ
Differentiating eq 29, we get dθCP 3 I dα ¼ dt dt
ð33Þ
ð35Þ
By integration 1=ðn 1Þ 1 α¼ 1 þ ðn 1ÞCd t
ð36Þ
Because Cd depends on k0 , the constant Cd is therefore a function of T, PCH4, PH2, and CCNi,f. Cd is fit to the following equation, which is a simpler form of k0 , where the four terms that depend on PH2 are combined in one term for which the power of hydrogen partial pressure must range between 0.5 and 1.5: Cd ¼ kd ðkd, C þ kd, CH4 PCH4 þ kd, H2 PH2 b Þ
ð37Þ
From eqs 32 and 36, we finally have P=ðn 1Þ 1 a¼ 1 þ ðn 1ÞCd t
ð28Þ
where N0 is the initial total number of active sites per mass of catalyst and Nt is the total number of active sites per mass of catalyst at time t. Because N0 = Nt + NCP 3 I, α becomes α ¼ 1 θ CP 3 I
ð34Þ
Combining eqs 25, 29, 33, and 34, we obtain a differential equation in α:
By incorporating eq 25 in eq 16, it becomes C0 k2 KCH4 1 θCP 3 I ¼ 2 k0 kC CCNi, f 2
¼ kd θC 3 I n
ð38Þ
To calculate the activity as a function of time, three parameters need to be determined: Cd, n, and P. These three parameters are determined by fitting each experiment independently. Note that the parameters n and P must be the same for all experiments, but Cd is dependent. Figures 7 and 8 show the fitting results for experiments containing 100% methane at different temperatures for the porous and nonporous catalysts, respectively (comparisons between the experimental and model results for 90% methane and 10% nitrogen at 550 °C for the porous and nonporous catalysts are shown in the Supporting Information, Figures S3 and S4). Figures 7 and 8 show the effect of the value of encapsulating order, n, on the modeling results. The best fits were obtained for the following values for P = 0.4 and n = 0.5. The order, n, is dependent on the carbon surface coverage. Using n = 0.5 and P = 0.4, the equation for deactivation becomes a¼
1
!0:8
1 0:5kd ðkd, C þ kd, CH4 PCH4 þ kd, H2 PH2 0:83 Þt
ð39Þ The parameter estimation results for eq 39 for both catalysts are listed in Table 5. The calculated parameters listed in Table 5 demonstrate the adequacy of the model to represent the deactivation behavior data for both catalysts (F test) at reasonable correlation coefficient values. The true mean of the parameters is estimated within a 95% confidence interval. The estimated values for the activation energies and adsorption heats are similar for porous and nonporous catalysts, which indicates that the deactivation mechanism is the same for both catalysts. The activation energy for encapsulating carbon formation is 147149 kJ/mol, and values obtained from the literature range from 9919 to 135 kJ/mol,16 and to 178 kJ/mol.15 12466
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research
ARTICLE
Figure 7. Deactivation order for the porous catalyst at different temperatures for 100% methane, 120 mL/min, 725 μm particles, and 1 atm; P = 0.4.
Figure 8. Deactivation order for the nonporous catalyst at different temperatures for 100% methane, 120 mL/min, 725 μm particles, and 1 atm; P = 0.4.
The value of Cd decreased significantly when CH4/H2 mixtures were used, compared to pure methane, at temperatures between 500 and 550 °C. The value of Cd is an indication of the deactivation rate, so a higher value of Cd leads to lower activity, as shown by eq 38. Hydrogen adsorbs more on the active nickel sites at lower temperatures,30,33 which decreases the availability or reactivity of the active sites, decreasing the reaction rate and the deactivation rate. However, at higher temperatures, the rapid kinetics shift the reaction to produce more carbon, which increases the fraction of active sites covered by carbon and increases the deactivation rate. The effect of adsorption on activity may explain conclusions drawn by Villacampa et al.,34
where any factor that increases the rate of methane cracking also increases the rate of deactivation and vice versa. Increasing temperature leads to increasing the fraction of adsorbed carbon sites, which in turn increases encapsulating carbon formation.33 For hydrogen, increasing PH2 decreases the reaction rate, but also occupies a fraction of the active sites, which slows the formation of encapsulating carbon. Since hydrogen adsorption is reversible, the longer deactivation time observed with CH4/H2 mixtures is explained. 3.6. Cracking Cycles. After studying deactivation of the fresh catalyst, regenerated catalyst deactivation was also examined. A series of cracking/regeneration cycles was 12467
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research
ARTICLE
Table 5. Parameter Estimation for the Activity Equation constant kd kd,C kd,CH4
parameter
porous
Ad
4904 ( 172
2214 ( 118
Ed (kJ/mol)
147 ( 0.26
149 ( 0.41
Ad
313 ( 6
531 ( 22
ΔHd (kJ/mol)
26 ( 0.15
26 ( 0.33
Ad,CH4 (atm1)
4082 ( 234 1845 ( 830
ΔHd,CH4 (J/mol) 3.56 ( 0.43 kd,H2
Ad,H2 (atmb)
n
ΔHd,H2 (kJ/mol) 81 ( 1 0.5
P
nonporous
3.33 ( 3.02
0.34 ( 0.03 0.39 ( 0.11
0.4
81 ( 2 0.5 0.4
b
0.83 ( 0.05
0.83 ( 0.11
F observed
301
131
number of experiments
30
31
F tabulated (0.05,2,29)
4.2
R2
0.98
Figure 9. Experimental activity data for different cracking cycles using the porous catalyst.
0.96
performed for porous and nonporous catalysts at 550 °C, 120 mL/min, and 100% methane, using 725 μm particles, and at 1 atm. The regeneration was conducted at the same conditions but using 100% air. The porous catalyst sustained its activity for 24 cycles, at which point the experiment was ended, while the nonporous catalysts showed deterioration in performance in the second cycle and activity was significantly decreased by the fifth cycle (Figures 9 and 10). The initial reaction rate varied with the number of cracking cycles, and an equation was developed for each catalyst to predict the initial rate for each catalyst. The initial reaction rate for the porous catalyst as a function of cycle number starting from the second cycle can be written as
Figure 10. Experimental activity data for different cracking cycles using the nonporous catalyst.
rCH4 ðcycleÞ ¼ rCH4 ðcycle ¼ 1Þ 3 1:275 expð 0:004 3 cycleÞ
ð40Þ The initial reaction rate for the nonporous catalyst as a function of cycle number can be calculated from the following equation: rCH4 ðcycleÞ ¼ rCH4 ðcycle ¼ 1Þ 3 1:5181 expð 0:426 3 cycleÞ
ð41Þ The activity changes with time for different cycles using the porous and nonporous catalysts are shown in Figures 9 and 10, respectively. The porous catalyst shows slower deactivation as the number of cycles increases, while the nonporous catalyst shows a remarkable increase in the deactivation rate as the number of cycles increases. The rapid decrease in activity for the nonporous catalyst implies that a change in the active site concentration or distribution leads to a change in the initial rate and deactivation. Sintering was observed as indicated by the reduction of nickel dispersion from 6.3% for the fresh, nonporous catalyst to 3% for the regenerated catalyst, as determined using hydrogen chemisorption. Accordingly, catalyst activity deteriorated quickly. The porous catalyst showed excellent stability over 24 cycles, which indicates that the porous catalyst could be applicable in cyclic hydrogen production via methane cracking. Equation 38 was applied for the different cracking cycles for the porous catalyst. The value of Cd was estimated based on the
Figure 11. Carbon deposited on the porous catalyst at 550 °C during different cycles.
best fit of the data. Due to the slower deactivation observed with the second cycle for the porous catalyst, Cd decreased as the number of cycles increased. An equation similar to eq 39 is developed to predict the activity at any cycle, but the value of Cd is corrected by multiplying it by a factor, f, to account for the change in deactivation rate over different cycles, as shown in the following equation: a¼
1 þ ðn 1Þkd f ðkd, C
1 þ kd, CH4 PCH4 þ kd, H2 PH2 0:83 Þt
!0:8
ð42Þ 12468
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research
ARTICLE
maintained good activity after 24 cycles. The model predicts for the porous catalyst a steady decrease in the maximum initial rate: this rate reaches 25% of its initial rate after 350 cycles.
’ ASSOCIATED CONTENT
bS Figure 12. Model prediction of initial reaction rate of the porous catalyst at 550 °C for different cycles.
For the porous catalyst, the factor f is calculated from the following polynomial: f ¼ 7E 5ðcycleÞ3 þ 0:0036ðcycleÞ2 0:0546ðcycleÞ þ 0:7829
Supporting Information. Parity plots for the porous and nonporous catalysts; comparsion beween experimental results and model predicion for nickel supported on porous and nonporous aluminas at 550 °C for 90% methane (balance nitrogen). This material is available free of charge via the Internet at http://pubs.acs.org.
’ AUTHOR INFORMATION Corresponding Author
*E-mail:
[email protected].
ð43Þ Figure 11 compares the experimental data and model results for different cycles for the porous catalyst. The model results agree well with the experimental data. Figure 12 shows the model prediction for the initial rate of the porous catalyst as a function of the number of cycles. The model prediction shows a steady decrease in the initial rate as the number of cycles increases until the initial rate reaches ∼50 mmolCH4/gNi/min after 170 cycles, which is equal to half the initial cracking rate calculated for the second cycle. After 380 cycles, the initial rate dropped to 25 mmolCH4/gNi/min. Equation 42 can be applied to predict the activity for cycles 139.
4. CONCLUSIONS Methane cracking kinetics were investigated using nickel supported on porous and nonporous aluminas. Models for predicting the initial reaction rate and activity decay were developed for the 500650 °C temperature range, at atmospheric pressure. The reaction mechanism begins with molecular adsorption of methane, with the release of the first hydrogen atom as the rate-limiting step. At constant temperature, the initial rate is a function of methane partial pressure and can be described by the LangmuirHinshelwood type equation that accounts for the competition between hydrogen and methane adsorption over the active sites. The activation energy for methane cracking is estimated at 88 and 75 kJ/mol for the porous and nonporous catalysts, respectively. The activation energy for encapsulating carbon formation, which is responsible for deactivation, was similar for both catalysts, 147 and 149 kJ/mol for the porous and nonporous catalysts, respectively. Catalyst deactivation depends on the temperature, methane and hydrogen partial pressures, and carbon surface coverage. The deactivation rate is half-order in surface carbon. Deactivation during successive cracking cycles was studied. The change in the active site distribution and concentration from one cycle to another altered the deactivation behavior of the catalyst. The support type affected the activation energy for methane cracking. The presence of all the nickel on the outer surface of a low surface area support increased the sintering possibility in successive cracking/regeneration cycles. While the nonporous catalyst significantly deactivated after five cycles, the porous catalyst
’ ACKNOWLEDGMENT The authors wish to thank the Natural Sciences and Engineering Research Council of Canada and the Ministry of Higher Education and Scientific Research of Egypt for financial support. ’ NOMENCLATURE a = activity coefficient A = Arrhenius constant b = order of hydrogen adsorption effect on a, atmb CNi,f = concentration of carbon dissolved in nickel on the nickel/ gas interface, molC/m3Ni CNi,r = concentration of carbon dissolved in nickel on the nickel/ support interface, molC/m3Ni Cfil = carbon deposited as carbon filaments E = activation energy, J/mol h = order of methane adsorption effect on activity, atm I = vacant site k = specific rate constant for methane cracking, mmolCH4/gNi/ min/atm K1, K2, K3, K4, K5 = equilibrium constants for steps 15 (eqs 48 KC = reciprocal of the equilibrium constant for step 7 (eq 10) KCH4 = adsorption constant for methane, atm1 KH = reciprocal of the equilibrium constant for step 6 (eq 9) KH2 = adsorption constant of hydrogen, atm3/2 KP = equilibrium constant for methane cracking, atm kd = specific rate constant for encapsulating carbon formation (deactivation) kd,H2 = specific rate constant for hydrogen adsorption effect on a kd,CH4 = specific rate constant for methane adsorption effect on a n = order of deactivation PCH4, PH2 = methane and hydrogen partial pressures, atm rexp = methane cracking rate from experimental work rcalc = methane cracking rate calculated from the developed model RTD = residence time distribution t = time, min X 3 I = active site occupied by intermediate X θX = fractional coverage of active sites with intermediate X α = fraction of active sites not deactivated ΔH = heat of adsorption, J/mol 12469
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470
Industrial & Engineering Chemistry Research
’ REFERENCES (1) Abbas, H. F.; Daud, W. M. A. W. Hydrogen production by thermocatalytic decomposition of methane using a fixed bed activated carbon in a pilot scale unit: Apparent kinetic, deactivation and diffusional limitation studies. Int. J. Hydrogen Energy 2010, 35 (22), 12268–12276. (2) Muradov, N. In Thermo-catalytic CO2-free Production of Hydrogen from Hydrocarbon Fuels; DOE hydrogen program review; The National Renewable Energy Laboratory for the U.S. Department of Energy, Baltimore, MD, USA, 2001; pp 271296. (3) Pinilla, J. L.; Suelves, I.; Lazaro, M. J.; Moliner, R.; Palacios, J. M. Parametric study of the decomposition of methane using a NiCu/Al2O3 catalyst in a fluidized bed reactor. Int. J. Hydrogen Energy 2010, 35 (18), 9801–9809. (4) Hazra, M.; Croiset, E.; Hudgins, R. R.; Silveston, P. L.; Elkamel, A. Experimental investigation of the catalytic cracking of methane over a supported Ni catalyst. Can. J. Chem. Eng. 2009, 87 (1), 99–105. (5) Rahman, M.; Croiset, E.; Hudgins, R. Catalytic Decomposition of Methane for Hydrogen Production. Top. Catal. 2006, 37 (2), 137–145. (6) Pinilla, J. L.; Suelves, I.; Lazaro, M. J.; Moliner, R.; Palacios, J. M. Activity of NiCuAl catalyst in methane decomposition studied using a thermobalance and the structural changes in the Ni and the deposited carbon. Int. J. Hydrogen Energy 2008, 33 (10), 2515–2524. (7) Muradov, N. Z.; Veziroglu, T. N. From hydrocarbon to hydrogen-carbon to hydrogen economy. Int. J. Hydrogen Energy 2005, 30 (3), 225–237. (8) Pinilla, J. L.; Moliner, R.; Suelves, I.; Lazaro, M. J.; Echegoyen, Y.; Palacios, J. M. Production of hydrogen and carbon nanofibers by thermal decomposition of methane using metal catalysts in a fluidized bed reactor. Int. J. Hydrogen Energy 2007, 32 (18), 4821–4829. (9) Ermakova, M. A.; Ermakov, D. Y.; Kuvshinov, G. G. Effective catalysts for direct cracking of methane to produce hydrogen and filamentous carbon: Part I. Nickel catalysts. Appl. Catal., A: Gen. 2000, 201 (1), 61–70. (10) Ermakova, M. A.; Ermakov, D. Y. Ni/SiO2 and Fe/SiO2 catalysts for production of hydrogen and filamentous carbon via methane decomposition. Catal. Today 2002, 77 (3), 225–235. (11) Takenaka, S.; Ogihara, H.; Yamanaka, I.; Otsuka, K. Decomposition of methane over supported-Ni catalysts: effects of the supports on the catalytic lifetime. Appl. Catal., A: Gen. 2001, 217 (12), 101–110. (12) Guo, J.; Lou, H.; Zheng, X. The deposition of coke from methane on a Ni/MgAl2O4 catalyst. Carbon 2007, 45 (6), 1314–1321. (13) Grabke, H. Evidence on the surface concentration of carbon on gamma iron from the kinetics of the carburization in CH4H2. Metall. Mater. Trans. B 1970, 1 (10), 2972–2975. (14) Grabke, H. The kinetics of decarburization and carburization of gamma-iron in methane-hydrogen mixtures. Phys. Chim. 1965, 69 (5), 14. (15) Borghei, M.; Karimzadeh, R.; Rashidi, A.; Izadi, N. Kinetics of methane decomposition to COx-free hydrogen and carbon nanofiber over Ni-Cu/MgO catalyst. Int. J. Hydrogen Energy 2010, 35 (17), 9479–9488. (16) Zavarukhin, S. G.; Kuvshinov, G. G. The kinetic model of formation of nanofibrous carbon from CH4-H2 mixture over a highloaded nickel catalyst with consideration for the catalyst deactivation. Appl. Catal., A: Gen. 2004, 272 (12), 219–227. (17) Demicheli, M. C.; Ponzi, E. N.; Ferretti, O. A.; Yeramian, A. A. Kinetics of carbon formation from CH4-H2 mixtures on nickel-alumina catalyst. Chem. Eng. J. 1991, 46 (3), 129–136. (18) Alstrup, I.; Tavares, M. T. Kinetics of Carbon Formation from CH4 + H2 on Silica-Supported Nickel and Ni-Cu Catalysts. J. Catal. 1993, 139 (2), 513–524. (19) Kuvshinov, G. G.; Mogilnykh, Y. I.; Kuvshinov, D. G. Kinetics of carbon formation from CH4-H2 mixtures over a nickel containing catalyst. Catal. Today 1998, 42 (3), 357–360. (20) Chen, D.; Lødeng, R.; Anundskas, A.; Olsvik, O.; Holmen, A. Deactivation during carbon dioxide reforming of methane over Ni catalyst: microkinetic analysis. Chem. Eng. Sci. 2001, 56 (4), 1371–1379. (21) Fogler, H. S. Elements of Chemical Reaction Engineering, 3rd ed.; Prentice Hall International Series: New York, 2002.
ARTICLE
(22) Butt, J. B.; Petersen, E. E. Activation, Deactivation, and Poisoning of Catalysts; Academic Press: San Diego, 1988. (23) Snoeck, J. W.; Froment, G. F.; Fowles, M. Kinetic Study of the Carbon Filament Formation by Methane Cracking on a Nickel Catalyst. J. Catal. 1997, 169 (1), 250–262. (24) Lee, M. B.; Yang, Q. Y.; Tang, S. L.; Ceyer, S. T. Activated dissociative chemisorption of CH4 on Ni(111): Observation of a methyl radical and implication for the pressure gap in catalysis. J. Chem. Phys. 1986, 85 (3), 1693–1694. (25) Hamza, A. V.; Madix, R. J. The activation of alkanes on Ni(100). Surf. Sci. 1987, 179 (1), 25–46. (26) Englezos, P.; Kalogerakis, N. Applied Parameter Estimation For Chemical Engineers; Marcel Dekker: New York, 2001. (27) Montgomery, D. C.; Runger, G. C. Applied Statistics and Probability for Engineers, 3rd ed.; John Wiley & Sons: New York, 2003. (28) Gilliland, E. R.; Harriott, P. Reactivity of Deposited Carbon. Ind. Eng. Chem. 1954, 46 (10), 2195–2202. (29) Chesnokov, V. V.; Zaikovskii, V. I.; Buyanov, R. A.; Molchanov, V. V.; Plyasova, L. M. Formation of Carbon Morphological Structures from Hydrocarbons over Nickel Containing Catalysts. Kinet. Catal. 1994, 35 (1), 6. (30) Bartholomew, C. H. Hydrogen adsorption on supported cobalt, iron, and nickel. Catal. Lett. 1990, 7 (1), 27–51. (31) Germer, L. H.; MacRae, A. U. Adsorption of Hydrogen on a (110) Nickel Surface. J. Chem. Phys. 1962, 37 (7), 1382–1386. (32) Beebe, J. T. P.; Goodman, D. W.; Kay, B. D.; Yates, J. J. T. Kinetics of the activated dissociative adsorption of methane on the low index planes of nickel single crystal surfaces. J. Chem. Phys. 1987, 87 (4), 2305–2315. (33) Schouten, F. C.; Kaleveld, E. W.; Bootsma, G. A. AES-LEEDellipsometry study of the kinetics of the interaction of methane with Ni(110). Surf. Sci. 1977, 63, 460–474. (34) Villacampa, J. I.; Royo, C.; Romeo, E.; Montoya, J. A.; Del Angel, P.; Monzon, A. Catalytic decomposition of methane over NiAl2O3 coprecipitated catalysts: Reaction and regeneration studies. Appl. Catal., A: Gen. 2003, 252 (2), 363–383.
12470
dx.doi.org/10.1021/ie201194z |Ind. Eng. Chem. Res. 2011, 50, 12460–12470