Unraveling the Reactions that Unravel Cellulose - American Chemical

Jun 11, 2012 - production.1 Oxygen-rich liquid bio-oil can be more easily stored and ..... based on bond dissociation energies at 0 K.61 ..... Pyrolys...
0 downloads 0 Views 2MB Size
Article pubs.acs.org/JPCA

Unraveling the Reactions that Unravel Cellulose Heather B. Mayes and Linda J. Broadbelt* Department of Chemical and Biological Engineering, Northwestern University, Evanston, Illinois, United States S Supporting Information *

ABSTRACT: For over 90 years, researchers have postulated mechanisms for the cleavage of cellulose’s glycosidic bonds and resulting formation of levoglucosan without reaching consensus. These reactions are key primary reactions in thermal processes for the production of carbon-neutral, renewable transportation fuels. Previous literature reports have proposed a variety of mainly heterolytic and homolytic mechanisms, but there has been insufficient evidence to settle the debate. Using density functional theory (DFT) methods and implicit solvent, we compared the likelihood of forming either radical or ionic intermediates. We discovered a concerted reaction mechanism that is more favorable than previously proposed mechanisms and is in better alignment with experimental findings. This new understanding of the mechanism of cellulose thermal decomposition opens the door to accurate process modeling and educated catalyst design, which are vital steps toward producing more cost-efficient renewable energy.



INTRODUCTION Biomass fast pyrolysis is a promising technology for the production of carbon-neutral, drop-in transportation fuels.1 Using less severe temperatures than in gasification, biomass fast pyrolysis is typically run near 500 °C to maximize liquid production.1 Oxygen-rich liquid bio-oil can be more easily stored and transported than a gas product2 and can be upgraded to produce carbon-neutral substitutes for conventional transportation fuels.3 Biomass fast pyrolysis can also be used to produce specialty chemicals and electricity and can use feeds from a variety of carbon-neutral, nonfood sources including agricultural waste and energy crops.1 Despite growing interest in biomass fast pyrolysis, kinetic models of this process are still rudimentary, relying on lumped components and lumped reaction steps.4,5 Even models for the pyrolysis of pure cellulose, a polymer of β-D-glucose and generally the dominant component of biomass, rely on lumped components or a fixed-percent conversion and predict very few individual products.5−14 A sampling of these kinetic models is shown in Figure 1. Such models cannot predict product speciation, limiting their usefulness as tools for process design. A mechanistic model capturing individual elementary reactions would predict product speciation and be sensitive to differences in feed and process conditions, allowing for process optimization to produce components of maximum value. Building a mechanistic kinetic model requires understanding of each elementary step and the associated kinetic parameters. However, there is uncertainty even in dominant reaction mechanisms. Levoglucosan is known to be a main product of cellulose pyrolysis and can account for up to 60 wt % yield.15−19 Despite investigation of this reaction over the past 90 years, there is controversy over whether the reaction proceeds via heterolytic or homolytic glycosidic bond © 2012 American Chemical Society

Figure 1. Several proposed kinetic models for cellulose pyrolysis, by (a) Broido and Nelson, 1975,7 (b) Bradbury, Sakai, and Shafizadeh (the “Broido-Shafizadeh Model”), 1979,8 and (c) Ranzi et al., 2008.13

cleavage.4,20−26 A simplified version of these competing hypotheses is shown in Scheme 1. Figure 2 illustrates proposed reaction mechanisms for the formation of levoglucosan from cellulose. An early hypothesis by Irvine and Oldham in 1921 put forth a two-step process by which acids catalyze starch hydrolysis to form glucose which then undergoes dehydration to levoglucosan (furanose form), shown in Figure 2a.20 To test this hypothesis, Ivanov et al. impregnated cellulose with glucose and found that the yield of levoglucosan decreased by half.15 Other studies found a more than 5-fold decrease in levoglucosan yield from glucose pyrolysis compared to cellulose pyrolysis.19,27 The production of levoglucosan from glucose pyrolysis confirms that it can be formed from glucose. However, if the main mechanism for levoglucosan formation from cellulose were through a glucose Received: January 12, 2012 Revised: June 7, 2012 Published: June 11, 2012 7098

dx.doi.org/10.1021/jp300405x | J. Phys. Chem. A 2012, 116, 7098−7106

The Journal of Physical Chemistry A

Article

Scheme 1. Simplified Representation of Transformation of a Cellulose Monomer (Left, in Brackets) to Levoglucosan (Right) via Heterolytic Cleavage (Top) or Homolytic Cleavage (Bottom)

The multitude of reactions in cellulose pyrolysis and the fast reaction times have stymied definitive experimental discrimination of the reaction scheme. This work aimed to elucidate this key step in cellulose pyrolysis using quantum mechanics (QM) and density functional theory (DFT). We discovered a concerted mechanism for glycosidic bond cleavage and levoglucosan formation that offers a lower-energy path than either the radical or ionic mechanisms. We used methyl-cellobiose, shown in Figure 3, as a model molecule for cellulose. Cellobiose, a glucose dimer, is the smallest subunit of cellulose that contains the glycosidic bond. A methyl group was added to the O4 of the nonreducing end, where the rest of a cellulose chain would be, to prevent hydrogen bonding at that location. Such hydrogen bonding could alter the conformation of the model molecule in a manner inconsistent with cellulose. The nonreducing end of cellobiose is of special interest since its O6 forms a bond with its C1 to create the second cyclic bond found in levoglucosan. Cellulose remains a solid at pyrolysis conditions and is a polar molecule. To capture the potential of the polar environment to stabilize reactions, we used implicit solvent to better model the electrostatic environment.

intermediate, glucose addition to cellulose pyrolysis or pure glucose pyrolysis would be expected to maintain or increase levoglucosan yield. The decreased levoglucosan yield indicates that the dominant mechanism for the formation of levoglucosan does not include a glucose intermediate. Several groups have reported that levoglucosan forms through an SN2 mechanism in which an intermediate epoxide conformer is formed.23,28,29 This mechanism can only proceed if the position of the C-1 hydroxyl group is in a conformation which does not sterically prohibit the SN2 attack by the C-2 hydroxyl. Later studies with sugars with hydroxyl groups in various positions showed that the formation of levoglucosan was insensitive to the orientation of the hydroxyl groups and thus disproved this mechanism.25 Other groups have proposed that cellulose decomposes to levoglucosan through a diradical mechanism as shown in Figure 2c.4,22−24 This theory arose by analogy with other macromolecules to explain observed bond rearrangements. The mechanism for simultaneous creation of two radical ends and subsequent rearrangement to stable products is unclear. Kislitsyn et al. presented a more detailed mechanism, showing each atom movement and conformation change (Figure 2d).21 After initial homolytic cleavage, three additional steps are required to form a molecule of levoglucosan. More recently, Shen and Gu proposed a mechanism in which an acetal reaction between C-1 and C-6 creates a hydroxyl radical from C-6, followed by the free hydroxyl radical bonding with a radical at C-4 to form levoglucosan, as shown in Figure 2e.26 For these multistep radical mechanisms, the same path must be consistently followed to create levoglucosan in the amounts experimentally observed. A less selective path leading to a variety of product species might be expected due to the reactive nature of radicals. Ponder et al. theorized a two-step mechanism with an ionic intermediate to yield levoglucosan, depicted in Figure 2f.25 The large influence of organic salts on cellulose pyrolysis makes it plausible that ions are also important in the pyrolysis of clean cellulose. However, the path of migration of the hydrogen from the O-6 on the reducing end to the glycosidic oxygen to form levoglucosan in the final step is not evident. Assary and Curtiss recently proposed a two-step mechanism for formation of levoglucosan from cellobiose.30 In the first step, the glycosidic bond is broken and a carbon−carbon double bond is formed. It is not clear how the hydrogen is transferred from O6 to C2 in the second step. Hosoya et al. reported a concerted reaction for thermal degradation of methyl β-glycoside.31 However, this smaller model and in vacuo study does not capture the full steric and electrostatic character of the glycosidic bond.



COMPUTATIONAL METHODS All QM calculations used Gaussian 09.32 Calculations were carried out in the gas phase as well as using the polarized continuum model (PCM) for implicit tetrahydrofuran (THF) or water solvent.33 It is important to note, however, that explicit water or THF molecules were not used in the calculations as they are not relevant to cellulose pyrolysis conditions. Consideration of the electrostatic environment engendered by cellulose is important, especially when comparing alternatives to the ionic mechanism, as ionic intermediates are significantly stabilized by a polar environment. Cellulose has a dielectric constant of 6.0 at 25 °C and shows a slight rise at higher temperature.34 While the dielectric constant of cellulose has not specifically been reported at fast pyrolysis conditions of 500 °C and near atmospheric pressure, THF, with its moderate dielectric constant of 7.43, was chosen to represent the cellulose environment at pyrolysis conditions. Calculations in gas and implicit water solvent were used to estimate how significantly the kinetic parameters would change due to any uncertainty in the estimate of the dielectric constant of cellulose at 500 °C. Only the electrostatic contribution from an implicit solvent was included, to better model the environment’s electrostatic potential. Cavitation and dispersion interaction energies were not included. Based on the comparison of a series of different basis sets, we selected the 6-31+G(2df,p) basis set for geometry optimiza7099

dx.doi.org/10.1021/jp300405x | J. Phys. Chem. A 2012, 116, 7098−7106

The Journal of Physical Chemistry A

Article

Figure 2. Several proposed mechanisms for the production of levoglucosan from cellulose.

functions on the non-hydrogen atoms were important especially for modeling intramolecular hydrogen bonding. Similar conclusions have been reported by other research

tions. Our calculations showed that diffuse functions on hydrogen atoms in this system did not appreciably affect results, while diffuse functions and additional polarized 7100

dx.doi.org/10.1021/jp300405x | J. Phys. Chem. A 2012, 116, 7098−7106

The Journal of Physical Chemistry A

Article

expensive to explore, encumbering identification of the lowestenergy reactant conformation. While little has been published about the structure of methyl-cellobiose, both experimental and theoretical studies have been completed to elucidate the structure of cellobiose. Therefore, we first found a low-energy structure for cellobiose based on reported structures and used it to build methyl-cellobiose. We began with the cellobiose structure reported by Chu and Jeffrey from their X-ray diffraction study.56 One-dimensional dihedral scans were performed in gas and implicit H2O environments to search for lower-energy conformations. Performing the dihedral scans with the PCM implicit water solvent did not uncover different low-energy conformations from those found by performing the dihedral scans in vacuo. These studies led to the low-energy conformer shown in Figure 4a, which was compared with other low-energy conformers from QM studies previously reported in the literature.

Figure 3. Methyl-cellobiose.

groups.35,36 Single point energies were computed with the 6311+G(3df,2p) basis set. Some calculations were completed with different basis sets in order to compare results with other studies. The M05-2X37 functional has been recommended for carbohydrate QM studies.35,38 We implemented the newer M06-2X39 functional. For selected molecules and reactions, we compared the M06-2X results with results from B3LYP,40 BMK,41 and MP2,42 as they have been used by other groups for other biomass conversion reactions. Frequency calculations used the same functional as was employed for geometry optimizations. Frequencies were scaled by the most appropriate factors reported by Merrick et al.43 The factors applied to frequencies from M06-2X/6-31+(2df,p) are 0.9657 for ZPE, 0.9450 for fundamental frequencies, 0.9194 for low frequencies (less than 260 cm−1), 0.9440 for enthalpy calculations, and 0.9355 for entropy calculations. Factors used for other functionals are enumerated in the Supporting Information. All minima have zero imaginary frequencies. One-dimensional dihedral scans were performed to search for lower-energy local minima. Transition states were verified to have exactly one imaginary frequency. The frequencies were visualized and intrinsic reaction coordinates (IRCs) were followed in both directions to ensure that the transition state was matched with the correct reaction. The Hessian-based Predictor-Corrector integration method was used for the IRC calculation.44−47 The Hessian was recomputed analytically every three predictor steps. Additionally, the transition states were confirmed to have no internal instabilities in the single-determinant wave function.48,49 Radical and ionic bond dissociation energies and formation of the products from concerted routes include counterpoise corrections calculated in the gas phase.50,51 Enthalpies and entropies at multiple temperatures were calculated from the Gaussian output using CalcTherm52 and the specified frequency scaling factors. Transition state theory53 was used to calculate rate coefficients at 100 K intervals between 300 and 1500 K. The results were fitted to an Arrhenius equation to yield frequency factors and activation energies. The minimum value for the coefficient of determination was 0.9999. Corrections for tunneling and internal rotation were not included. At 500 °C, the transmission coefficient is generally near unity.54,55 Accurate treatment of internal modes of rotation is impeded by coupled motion and by large contributions to the low frequencies from ringpuckering. Thus, the frequencies were calculated based on the harmonic oscillator assumption. These simplifications affect the estimation of the pre-exponential factor more than the activation energy.

Figure 4. Three low-energy conformations of cellobiose.

Strati et al. performed a study of 27 cellobiose conformers at the B3LYP/6-311++G** level of theory and found low-energy syn and anti conformers in vacuo.57 The syn structures correspond to glycosidic dihedral angles similar to those experimentally obtained for native cellulose,58 in which the C-6 hydroxymethyl “arms” of the two monomers are on opposite sides of the glycosidic bond. The syn structure with the lowest electronic energy is shown in Figure 4b. This structure had the second-lowest (0.24 kcal/mol higher than the lowest) Gibbs free energy at 298 K of the syn conformers and corresponds to their conformer X. Strati’s global minimum energy structure was an anti conformer (conformer XIX) and is shown in Figure 4c. The relative energies of these three conformers are compared in Table 1, using the conformer in Figure 4c as the reference point. Comparison of the electronic energies in vacuo is included for the B3LYP/6-311++G** level of theory to aid comparison with Strati et al.’s work.57 The rest of the calculations include the electrostatic contribution from PCM implicit THF solvent. The total energies are included in the Supporting Information. While an anti conformer was found to be the globalminimum energy structure in both gas and implicit solvent, it was not used for this study as this orientation is not representative of cellulose. While it is possible that an end group could rotate to the anti conformation, the bulk of the glycosidic linkages are expected to be found in the



RESULTS AND DISCUSSION It is vital to use low-energy conformers of reactant molecules in QM studies to ensure that the full energy barrier between the reactant and a transition state is captured. Methyl-cellobiose’s 24 non-hydrogen atoms and 138 degrees of freedom result in a multitude of low-energy conformers, which are computationally 7101

dx.doi.org/10.1021/jp300405x | J. Phys. Chem. A 2012, 116, 7098−7106

The Journal of Physical Chemistry A

Article

Table 2. Comparison of Estimated Energy Requirementsa for Glycosidic Cleavage of Methyl-Cellobiose via Homolytic, Heterolytic, or Concerted Mechanisms In Vacuo or Including Electrostatic Effects of Implicit THF or Water Solvent

Table 1. Comparisons of Relative Energies for the Three Conformers in Figure 4, in kcal/mol level of theory

solvent

ΔE 0 K B3LYP/6-311++G** none B3LYP/6-311++G** THF MP2/6-311+G(3df,2p)//B3LYP/ THF 6-31+G(2df,p) M06-2X/6-311+G(3df,2p)//M06− THF 2X/6-31+G(2df,p) ΔG 500 °C B3LYP/6-311++G** none B3LYP/6-311++G** THF MP2/6-311+G(3df,2p)//B3LYP/ THF 6-31+G(2df,p) M06-2X/6-311+G(3df,2p)//M06THF 2X/6-31+G(2df,p)

Figure 4a

Figure 4b

Figure 4c

5.56 2.76 4.25

3.32 2.87 5.45

0.00 0.00 0.00

4.51

5.20

0.00

6.05 4.48 6.06

4.20 5.00 7.68

0.00 0.00 0.00

4.00

5.86

0.00

phase

radicalb

ionicb

concertedc

gas implicit THF implicit water

97.8 94.5 94.6

167.2 90.3 78.6

55.2 53.9 53.5

a

Values are in kcal/mol and were calculated at the M06-2X/6311+G(3df,2p)//M06-2X/6-31+G(2df,p) level of theory. bDissociation enthalpies at 500 °C. cActivation energy.

expected, the implicit solvent has a significant stabilizing effect on the ionic mechanism and little effect on the radical mechanism. The results with implicit water at 500 °C provide a lower bound of 78.6 kcal/mol for the cleavage of this bond into two ions and a lower bound of 94.6 kcal/mol for the cleavage into two radicals. The gas phase results are consistent with a recent DFT study using the B3LYP/6-31+G* level of theory in gas phase estimating the energy barrier at 178 kcal/mol for heterolytic cleavage and 79 kcal/mol for homolytic cleavage based on bond dissociation energies at 0 K.61 As an alternative to the ionic or radical mechanisms, we propose a concerted reaction in which the glycosidic bond is cleaved at the same time that the 1,6-ring of levoglucosan is formed. The isodesmic reaction proceeds via a four-centered transition state, with both a carbon−oxygen bond and an oxygen−hydrogen bond being both broken and formed. Kinetic parameters were determined for several low-energy conformations for this transition state, using the lowest-energy conformation of methyl-cellobiose as the reactant. The Cartesian coordinates for the atoms in three low-energy conformations are included in the Supporting Information. The transition state which resulted in the highest reaction rate coefficient at 500 °C is shown in Figure 6. From the frequency

experimentally observed syn conformation. The vast majority of the bonds will be interior bonds in the solid, without freedom to rotate. Furthermore, follow-up studies to Strati et al.’s work done with DFTMD show that when solvent and entropy effects are better captured, the resulting simulations find the syn conformation to be lower in energy.59 Of the syn conformers found, the conformer in Figure 4a has the lowest Gibbs free energy at 500 °C and implicit THF. The hydroxymethyl groups are in the gauche−gauche (GG) orientation. The GG orientation has been shown to be an accessible cellulose structure at higher temperature.60 A methyl group was added to the nonreducing end of this conformer to build the methyl-cellobiose conformation used as a reactant in this study, and one-dimensional dihedral scans were performed to ensure the resulting structure was also a low-energy conformation. The resulting structure is shown in Figure 5. The Supporting Information lists the energies and Cartesian coordinates of several low-energy conformations of methylcellobiose.

Figure 6. Transition state for the concerted mechanism initiation step from geometry optimization at the M06-2X/6-31+G(2df,p) level of theory with implicit THF solvent. Distances shown are in Å.

Figure 5. Geometry for methyl-cellobiose used as the reactant to model the initiation step of the conversion of cellulose to levoglucosan. The geometry was optimized using the M06-2X/6-31+G(2df,p) level of theory and implicit THF.

calculation at the M06-2X/6-31+G(2df,p) level of theory with implicit THF, the unscaled wavenumber of the imaginary frequency is 90i cm−1, and the smallest real unscaled wavenumber is 20.2 cm−1. In the gas phase, the unscaled imaginary frequency is 137i cm−1, and the distances between the carbon and oxygens in the transition state center are slightly decreased from those shown in Figure 7 for the structure optimized in implicit THF. In implicit water, the unscaled imaginary frequency is 82i cm−1, and the distances between the carbon and oxygens in the transition state center are slightly increased. The gas phase imaginary frequency is similar to those reported for other transition states for converted bond cleavage

The energy requirements for heterolytic and homolytic glycosidic cleavage were investigated by comparing the enthalpy for dissociation into either a cation and anion or two radicals, as shown in Scheme 1. Table 2 compares the calculated bond dissociation energies for the initial glycosidic cleavage for methyl-cellobiose in gas, implicit THF, and implicit water. Enthalpy differences between the reactant and radical or ionic fragments provide lower-bound estimates for the activation energy of homolytic or heterolytic cleavage. As 7102

dx.doi.org/10.1021/jp300405x | J. Phys. Chem. A 2012, 116, 7098−7106

The Journal of Physical Chemistry A

Article

Figure 7. Geometry for methyl-cellobiosan used as the reactant to model depropagation of cellulose to levoglucosan, with the methyl group in place of the rest of the chain. The geometry was optimized using the M06−2X/6-31+G(2df,p) level of theory and implicit THF.

Figure 9. Free energy diagram for the glycosidic cleavage of methylcellobiose at 500 °C, modeling the initiation step of the concerted mechanism which unravels cellulose and leads to levoglucosan. The Gibbs free energies were calculated at the M06-2X/6-311+G(3df,2p)//M06-2X/6-31+G(2df,p) level of theory with implicit THF.

and formation processes in the gas phase.62−64 The IRC connecting the transition state to the reactant and product is included in the Supporting Information. The flat transition state involves large geometric changes, separating the dimer into parts. The kinetic parameters for this transition state are shown in Table 3. There is little difference in the energy barriers between

Table 4 includes the calculated energy requirements using B3LYP and MP2 to aid comparison of the results reported here Table 4. Gibbs Free Energy Differences between the Transition State and Reactant and Activation Energies for the Glycosidic Cleavage of Methyl-Cellobiose Calculated with Implicit THF Solvent, in kcal/mol

Table 3. Arrhenius Parameters and Reaction Rate Coefficient at 500 °C for Concerted Glycosidic Cleavage of Methyl-Cellobiose (Initiation) at the M06-2X/6311+(3df,2p)//M06-2X/6-31+(2df,p) Level of Theory environment

A s−1

gas implicit THF implicit water

4.7 × 10 1.2 × 1013 9.4 × 1012 13

Ea kcal/mol

k at 500 °C s−1

55.2 53.9 53.5

1.2 × 10−2 6.9 × 10−3 6.9 × 10−3

the gas phase, implicit THF, and implicit water for the concerted mechanism. Figure 9 displays the free energy diagram for glycosidic cleavage of methyl-cellobiose calculated with M06-2X/6-311+G(3df,2p)//M06-2X/6-31+G(2df,p) and implicit THF solvent.

level of theory

ΔG‡, 500 °C

Ea

B3LYP/6-31+G(2df,p) BMK/6-31+G(2df,p) M06-2X/6-31+G(2df,p) MP2/6-311+G(3df,2p)//B3LYP/6-31+G(2df,p) M06-2X/6-311+G(3df,2p)//M06-2X/6-31+G(2df,p)

46.3 53.5 54.3 54.5 54.4

48.6 53.3 53.8 56.8 53.9

with previously reported energies for alternate mechanisms and similar systems.30,57,61 BMK was also included to provide an additional point of reference as it was optimized specifically for kinetic parameters.41 The M06-2X functional has been recently used in many biomass decomposition studies with different basis sets than the ones used in this work.65−67 There are only small differences in our results for M06-2X/6-31+G(2df,p) and M06-2X/6-311+G(3df,2p)//M06-2X/6-31+G(2df,p), as shown in Table 4. As previously noted, M06-2X has been shown to provide reliable estimates of energy barriers. The significantly lower energy barrier from B3LYP is consistent with previous studies that noted its underestimation of energy barriers.35,38,66,68−70 The MP2 barrier is the highest, which is consistent with other studies that noted a tendency for MP2 to overestimate energy barriers.66,69,71 For cellulose, a midchain initial concerted glycosidic cleavage yields two polymeric fragments: a cellulose-like polymer with a terminal levoglucosan-like end and a shorter cellulose chain. Thus, the initiation step can be followed by depropagation steps in which a molecule of levoglucosan is released from the chain end for each subsequent scission of a glycosidic bond. The depropagation step can be repeated until the cellulose chain has fully decomposed into levoglucosan molecules, as shown in Scheme 2. We investigated the depropagation step in Scheme 2 using cellobiosan with an additional methyl group to represent the rest of the chain, as shown in Figure 7. The smallest unscaled wavenumber for this structure is 36.4 cm−1. The methyl group is included in place of the rest of the chain to prevent hydrogen bonding at that position that would not be available in the chain. While the middle of a cellulose chain is not expected to have freedom to rotate to the anti conformer,

Figure 8. Transition state geometry for methyl-cellobiosan representing the depropagation of cellulose to levoglucosan. The geometry was optimized using the M06-2X/6-31+G(2df,p) level of theory and implicit THF. Distances shown are in Å.

The reaction rate coefficients for ionic or radical scission were estimated to allow a comparison of reaction rates for the initial step of all three mechanisms. Using the lower bounds for activation energy reported above and a typical upper bound value of 1 × 1016 s−1 for the pre-exponential factor, the upper bound of the reaction rate coefficient is 1.9 × 10−11 s−1 for initiation via the homolytic mechanism and 5.9 × 10−7 s−1 for initiation via the heterolytic mechanism. Thus, initiation via the concerted mechanism is at least 4 orders of magnitude faster than either of these previously proposed mechanisms. 7103

dx.doi.org/10.1021/jp300405x | J. Phys. Chem. A 2012, 116, 7098−7106

The Journal of Physical Chemistry A

Article

calculated kinetic parameters for initiation and depropagation in gas, implicit THF, and implicit water. The values calculated with implicit THF solvent are believed to best represent the reaction under the conditions of cellulose pyrolysis. This concerted mechanism is clearly more favorable than either homolytic or heterolytic cleavage as indicated by the calculated reaction rate coefficients. It is also lower in energy than the two-step mechanism for levoglucosan formation from cellobiose without radical or ionic intermediates proposed by Assary and Curtiss.30 They reported a ΔH‡ at 298 K in vacuo of 59.1 kcal/mol at the MP2/6-311++G(3df,3pd)//B3LYP/631G(2df,p) level of theory for the first step of the mechanism. They assumed that the second step is fast, despite the simultaneous requirements of dissociation of an oxygen− hydrogen bond, dissociation of a carbon−carbon double bond, formation of an oxygen−carbon bond, and formation of a carbon−hydrogen bond. No transition state or rate coefficient was provided to support this assumption. We calculated ΔH‡ at 298 K in vacuo for the one-step mechanism for formation of levoglucosan from cellobiose using the same level of theory and without frequency scaling factors to allow direct comparison of calculated energies. The mechanism presented here has a lower ΔH‡ at 298 K in vacuo of 52.7 kcal/mol. While it is experimentally difficult to validate this finding because efforts to observe the kinetics of this process rely on mass-loss data, which provides convoluted values for the production of condensable vapor and is often further confounded by heat-lag effects, some comparisons are possible. Antal et al. set out to control heat-lag effects and estimated the activation energy for pyrolytic weight loss of cellulose to be 54.5 kcal/mol.72 Because the primary condensable product is levoglucosan, it is reasonable to compare this value to the activation energy for concerted glycosidic cleavage. It is in excellent agreement to the activation energy of 53.9 kcal/mol presented here for the concerted initiation step and reasonable agreement with the activation energy of 58.5 kcal/mol for depropagation. Further confirmation has been provided by kinetic models our group is developing. When the concerted mechanism and its calculated activation energy were used, our model predicts a product distribution that is in good agreement with the detailed speciation results experimentally determined by Patwardhan et al.19

Scheme 2. Concerted Mechanism for Levoglucosan Formation via Initiation and Depropagation Steps

as discussed above, we expect ends to have freedom to rotate. In the absence of experimental studies to help determine the lowest-energy structure of methyl-cellobiosan or its simpler surrogate, cellobiosan, we did not restrict our search for a lowenergy reactant to syn conformations of the levoglucosan-like end. The anti-conformation of the reactant and transition state were lower in energy than syn conformations. The structure of the transition state for the depropagation reaction is shown in Figure 8. Cartesian coordinates and relative energies of several low-energy conformations of the reactant and transition state are included in the Supporting Information. The transition state shown in Figure 8 corresponds to the conformation which yielded the highest rate coefficient as determined using transition state theory. The unscaled wavenumber of the imaginary frequency for this conformation is 156.5i cm−1 and the smallest real unscaled wavenumber is 25.6 cm−1. A free energy diagram for the depropagation step is shown in Figure 10. Table 5 shows the



CONCLUSIONS This work unraveled the reaction mechanism for the primary reaction in fast pyrolysis of clean cellulose. The one-step mechanism shows a highly plausible route for the high yield of levoglucosan from cellulose. The concerted glycosidic cleavage and levoglucosan formation were modeled in a PCM THF solvent to capture the electrostatic effects of the cellulose environment. The energy barriers for this concerted reaction are approximately 25 kcal/mol lower and the rate coefficient is at least 4 orders of magnitude higher than for ionic cleavage, the previous paradigm for cellulose glycosidic cleavage. Incorporation of the concerted mechanism and its associated kinetic parameters into a kinetic model our group is developing results in product speciation and yields that are in very good agreement with experimental data. In addition to providing data for mechanistic modeling of cellulose pyrolysis, this unraveled mechanism for pure cellulose pyrolysis provides a basis for studying catalytic interventions in cellulose pyrolysis, including the action of the naturally occurring inorganic salts. Efforts to quantify these effects at the mechanistic level are underway.

Figure 10. Free energy diagram for the glycosidic cleavage of methylcellobiosan at 500 °C, modeling the depropagation step of the concerted mechanism which unravels cellulose chains that have been initiated and forms levoglucosan. The Gibbs free energies were calculated at the M06-2X/6-311+G(3df,2p)//M06-2X/6-31+G(2df,p) level of theory with implicit THF.

Table 5. Arrhenius Parameters and Reaction Rate Coefficient at 500 °C for Concerted Glycosidic Cleavage of Methyl-Cellobiosan (Depropagation) at the M06-2X/6311+(3df,2p)//M06-2X/6-31+(2df,p) Level of Theory environment

A s−1

gas implicit THF implicit water

9.7 × 10 1.2 × 1014 2.5 × 1013 13

Ea kcal/mol

k at 500 °C s−1

60.7 58.5 56.3

6.4 × 10−4 3.4 × 10−3 3.1 × 10−3

7104

dx.doi.org/10.1021/jp300405x | J. Phys. Chem. A 2012, 116, 7098−7106

The Journal of Physical Chemistry A



Article

(23) Gardiner, D. J. Chem. Soc. C 1966, 1473−1476. (24) Kato, K. Agric. Biol. Chem. 1967, 31, 657−663. (25) Ponder, G. R.; Richards, G. N.; Stevenson, T. T. J. Anal. Appl. Pyrolysis 1992, 22, 217−229. (26) Shen, D. K.; Gu, S. Bioresour. Technol. 2009, 100, 6496−6504. (27) Golova, O. P.; Andrievskaya, E. A.; Pakhomov, A. M.; Merlis, N. M. Russ. Chem. Bull. 1957, 6, 399−401. (28) Kilzer, F. J.; Broido, A. Pyrodynamics 1965, 2, 151−163. (29) Molton, P. M.; Demmitt, T. F. Reaction Mechanisms in Cellulose Pyrolysis: A Literature Review; Battelle, Pacific Northwest Laboratories: Richland, WA, 1977 (30) Assary, R. S.; Curtiss, L. A. ChemCatChem 2012, 4, 200−205. (31) Hosoya, T.; Nakao, Y.; Sato, H.; Kawamoto, H.; Sakaki, S. J. Org. Chem. 2009, 74, 6891−6894. (32) Frisch, M. J. et al. Gaussian 09, Revision B.01; Gaussian, Inc.: Wallingford, CT, 2010. (33) Tomasi, J.; Mennucci, B.; Cammi, R. Chem. Rev. 2005, 105, 2999−3093. (34) Calkins, C. R. TAPPI 1950, 33, 278−285. (35) Csonka, G. I.; French, A. D.; Johnson, G. P.; Stortz, C. A. J. Chem. Theory Comput. 2009, 5, 679−692. (36) Curtiss, L. A.; Redfern, P. C.; Raghavachari, K. Chem. Phys. Lett. 2010, 499, 168−172. (37) Zhao, Y.; Schultz, N. E.; Truhlar, D. G. J. Chem. Theory Comput. 2006, 2, 364−382. (38) Wodrich, M. D.; Corminboeuf, C.; Schreiner, P. R.; Fokin, A. A.; von Ragué Schleyer, P. Org. Lett. 2007, 9, 1851−1854. (39) Zhao, Y.; Truhlar, D. G. Theor. Chem. Acc. 2008, 120, 215−241. (40) Becke, A. D. J. Chem. Phys. 1993, 98, 5648−5652. (41) Boese, A. D.; Martin, J. M. L. J. Chem. Phys. 2004, 121, 3405− 3416. (42) Pople, J. A.; Binkley, J. S.; Seeger, R. Int. J. Quantum Chem. 1976, 10, 1−19. (43) Merrick, J. P.; Moran, D.; Radom, L. J. Phys. Chem. A 2007, 111, 11683−11700. (44) Hratchian, H. P.; Schlegel, H. B. J. Chem. Phys. 2004, 120, 9918−9924. (45) Hratchian, H. P.; Schlegel, H. B. In Theory and Applications of Computational Chemistry: The First Forty Years; Dykstra, C. E., Frenking, G., Kim, K. S., Scuseria, G., Eds.; Elsevier Science: New York, 2005; pp 195−249. (46) Hratchian, H. P.; Schlegel, H. B. J. Chem. Theory Comput. 2005, 1, 61−69. (47) Collins, M. A. Theor. Chem. Acc. 2002, 108, 313−324. (48) Seeger, R.; Pople, J. A. J. Chem. Phys. 1977, 66, 3045−3050. (49) Bauernschmitt, R.; Ahlrichs, R. J. Chem. Phys. 1996, 104, 9047− 9052. (50) Boys, S. F.; Bernardi, F. Mol. Phys. 1970, 19, 553−566. (51) Simon, S.; Duran, M.; Dannenberg, J. J. J. Chem. Phys. 1996, 105, 11024−11031. (52) Pfaendtner, J.; Yu, X.; Broadbelt, L. J. Theor. Chem. Acc. 2007, 118, 881−898. (53) Moore, J. W.; Pearson, R. G. Kinetics and Mechanism, 3rd ed.; John Wiley and Sons: New York, 1981. (54) Truong, T. N.; Duncan, W. T.; Tirtowidjojo, M. Phys. Chem. Chem. Phys. 1999, 1, 1061−1065. (55) Hong, X.; Sun, H.; Law, C. K. Comput. Theor. Chem. 2011, 963, 357−364. (56) Chu, S. S. C.; Jeffrey, G. A. Acta Crystallogr., Sect. B 1968, 24, 830−838. (57) Strati, G. L.; Willett, J. L.; Momany, F. A. Carbohydr. Res. 2002, 337, 1833−1849. (58) Nishiyama, Y.; Langan, P.; Chanzy, H. J. Am. Chem. Soc. 2002, 124, 9074−9082. (59) Momany, F. A.; Schnupf, U. Carbohydr. Res. 2011, 346, 619− 630. (60) Matthews, J. F.; Bergenstråhle, M.; Beckham, G. T.; Himmel, M. E.; Nimlos, M. R.; Brady, J. W.; Crowley, M. F. J. Phys. Chem. B 2011, 115, 2155−2166.

ASSOCIATED CONTENT

S Supporting Information *

Scaling factors for all levels of theory used; comparison of computational methods for predicting geometries; cellobiose conformation total electronic and free energies; and Cartesian coordinates for atoms in key molecular structures. This material is available free of charge via the Internet at http://pubs.acs.org.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was supported by the Department of Energy (DOE), Office of Energy Efficiency and Renewable Energy (EERE), through the Office of Biomass Program, Grant Number DEEE0003044, and by the DOE Computational Science Graduate Fellowship (CSGF), which is provided under Grant Number DE-FG02-97ER25308. This research used resources of the National Energy Research Scientific Computing Center, which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231. The work was also supported by an award from the ARCS Foundation Inc., Chicago Chapter.



REFERENCES

(1) Huber, G. W.; Iborra, S.; Corma, A. Chem. Rev. 2006, 106, 4044− 4098. (2) Bridgwater, A. V.; Meier, D.; Radlein, D. Org. Geochem. 1999, 30, 1479−1493. (3) Elliott, D. C. Energy Fuels 2007, 21, 1792−1815. (4) Mohan, D.; Pittman, C. U., Jr.; Steele, P. H. Energy Fuels 2006, 20, 848−889. (5) Di Blasi, C. Prog. Energy Combust. Sci. 2008, 34, 47−90. (6) Byrne, G. A.; Gardiner, D.; Holmes, F. H. J. Appl. Chem. 1966, 16, 81−88. (7) Broido, A.; Nelson, M. A. Combust. Flame 1975, 24, 263−268. (8) Bradbury, A. G. W.; Sakai, Y.; Shafizadeh, F. J. Appl. Polym. Sci. 1979, 23, 3271−3280. (9) Agrawal, R. K. Can. J. Chem. Eng. 1988, 66, 403−412. (10) Diebold, J. P. Biomass Bioenergy 1994, 7, 75−85. (11) Antal, M. J., Jr.; Varhegyi, G. Ind. Eng. Chem. Res. 1995, 34, 703−717. (12) Várhegyi, G.; Antal, M. J., Jr.; Jakab, E.; Szabó, P. J. Anal. Appl. Pyrolysis 1997, 42, 73−87. (13) Ranzi, E.; Cuoci, A.; Faravelli, T.; Frassoldati, A.; Migliavacca, G.; Pierucci, S.; Sommariva, S. Energy Fuels 2008, 22, 4292−4300. (14) Authier, O.; Ferrer, M.; Khalfi, A.-E.; Lédé, J. Int. J. Chem. React. Eng. 2010, 8, A78. (15) Ivanov, V. I.; Golova, O. P.; Pakhomov, A. M. Russ. Chem. Bull. 1956, 5, 1295−1296. (16) Shafizadeh, F. J. Anal. Appl. Pyrolysis 1982, 3, 283−305. (17) Lakshmanan, C. M.; Gal-or, B.; Hoelscher, H. E. Ind. Eng. Chem. Prod. Res. Dev. 1969, 8, 261−267. (18) Lomax, J. A.; Commandeur, J. M.; Arisz, P. W.; Boon, J. J. J. Anal. Appl. Pyrolysis 1991, 19, 65−79. (19) Patwardhan, P. R.; Satrio, J. A.; Brown, R. C.; Shanks, B. H. J. Anal. Appl. Pyrolysis 2009, 86, 323−330. (20) Irvine, J. C.; Oldham, J. W. H. J. Chem. Soc., Faraday Trans. 1921, 119, 1744−1759. (21) Kislitsyn, A. N.; Rodionova, Z. M.; Savinykh, V. I.; Gusev, A. V. Zh. Prikl. Khim. (Leningrad) 1971, 44, 2518−2524. (22) Pakhomov, A. M. Russ. Chem. Bull. 1958, 6, 1525−1527. 7105

dx.doi.org/10.1021/jp300405x | J. Phys. Chem. A 2012, 116, 7098−7106

The Journal of Physical Chemistry A

Article

(61) Zhang, X.; Li, J.; Yang, W.; Blasiak, W. Energy Fuels 2011, 25, 3739−3746. (62) Kraka, E.; Wu, A.; Cremer, D. J. Phys. Chem. A 2003, 107, 9008−9021. (63) Tranter, R. S.; Klippenstein, S. J.; Harding, L. B.; Giri, B. R.; Yang, X.; Kiefer, J. H. J. Phys. Chem. A 2010, 114, 8240−61. (64) Bianciotto, M.; Barthelat, J.-C.; Vigroux, A. J. Am. Chem. Soc. 2002, 124, 7573−7587. (65) Beste, A.; Buchanan, A. C., III J. Org. Chem. 2009, 74, 2837− 2841. (66) Kim, S.; Chmely, S. C.; Nimlos, M. R.; Bomble, Y. J.; Foust, T. D.; Paton, R. S.; Beckham, G. T. J. Phys. Chem. Lett. 2011, 2, 2846− 2852. (67) Parthasarathi, R.; Bellesia, G.; Chundawat, S. P. S.; Dale, B. E.; Langan, P.; Gnanakaran, S. J. Phys. Chem. A 2011, 115, 14191−14202. (68) Nimlos, M. R.; Blanksby, S. J.; Ellison, G. B.; Evans, R. J. J. Anal. Appl. Pyrolysis 2003, 66, 3−27. (69) Temelso, B.; Sherrill, C. D.; Merkle, R. C.; Freitas, R. A., Jr. J. Phys. Chem. A 2006, 110, 11160−11173. (70) Zhao, Y.; Truhlar, D. G. J. Chem. Theory Comput. 2008, 4, 1849−1868. (71) Beste, A.; Buchanan, A. C.; Britt, P. F.; Hathorn, B. C.; Harrison, R. J. J. Mol. Struct.: THEOCHEM 2008, 851, 232−241. (72) Antal, M. J., Jr.; Várhegyi, G.; Jakab, E. Ind. Eng. Chem. Res. 1998, 37, 1267−1275.

7106

dx.doi.org/10.1021/jp300405x | J. Phys. Chem. A 2012, 116, 7098−7106