Using Electronic Structure Calculations To Investigate the

Chapter 3

Using Electronic Structure Calculations To Investigate the Kinetics of Gas-Phase Ammonia Synthesis Downloaded via AUBURN UNIV on August 20, 2019 at 19:37:10 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

Kelsey M. Stocker* Department of Chemistry and Biochemistry, Suffolk University, 8 Ashburton Place, Boston, Massachusetts 02108, United States *E-mail: [email protected].

Computational chemistry techniques are valuable tools for teaching concepts in thermodynamics and chemical kinetics. In this experiment, undergraduate physical chemistry students gain valuable, authentic experience with the tools of computational chemistry and a more “hands-on” interaction with the energy landscape of a chemical reaction. Students use electronic structure calculations to determine the geometry, vibrational frequencies, and energy of the reactants, products, intermediates, and transition states in a four-step gas-phase ammonia synthesis reaction. Using transition state theory, they construct the reaction coordinate diagram and calculate activation energies, rate constants, and equilibrium constants.

Introduction Computational methods can play a valuable role in undergraduate chemistry courses and have been successfully integrated into the curriculum in a variety of courses (1–4). Electronic structure calculations in particular can contribute to student learning by allowing them to gain “hands on” experiences at the atomic level in ways that traditional undergraduate experiments cannot. Electronic structure methods have been used in excellent computational experiments on reaction mechanisms for use in undergraduate organic or physical chemistry courses (5–11). I developed this lab to fit into the Advanced Theories of Reaction Rates module of my secondsemester 400-level physical chemistry course. At this point in the course, my students have learned about the laws of thermodynamics, physical and chemical equilibria, reaction mechanisms, rate laws, and potential energy surfaces. I wrote this lab to give my students the experience of studying the kinetics of a realistic, multi-step reaction mechanism using reaction coordinate diagrams and transition state theory. There are five primary learning goals associated with this experiment: © 2019 American Chemical Society

Grushow and Reeves; Using Computational Methods To Teach Chemical Principles ACS Symposium Series; American Chemical Society: Washington, DC, 2019.

• • • • •

Perform authentic electronic structure calculations Construct a reaction coordinate diagram for a multi-step reaction Apply the thermodynamic formulation of transition state theory Calculate activation energies, rate constants, and equilibrium constants Use kinetics data to evaluate plausibility of reaction mechanism

Throughout this chapter I have noted possible adaptations and extensions that can be used to make the experiment appropriate for a variety of courses. I have also attempted to make the description of the procedure general enough that it can be used with various electronic structure programs and interface packages. I teach a “quantum first” approach, so by the time we conduct this experiment my students have been exposed to partition functions, statistical mechanics, and computational chemistry techniques. However, the experiment utilizes the thermodynamic formulation of transition state theory so it could easily be used in a “thermo first” course design (12).

Why Ammonia Synthesis? The formation of ammonia is a vitally important process made up of structurally simple species, which makes it an ideal reaction for a computational kinetics experiment (13). After investigating the kinetics of the multi-step gas-phase mechanism, students have a deeper appreciation for the importance of heterogeneous catalysis in industrial processes. The reaction may look harmless enough on paper, but it turns out to be so devilishly difficult to carry out that the first industrially viable process yielded two separate Nobel Prizes in Chemistry: Fritz Haber in 1918 and Carl Bosch in 1931 (14, 15). I typically can’t resist including a brief synopsis of Fritz Haber’s story in the rest of my pre-lab lecture (16). The synthesis of ammonia has history. But it’s not old news by a longshot considering that new, more efficient catalysts remain an area of active research (17). Even the gas-phase reaction, while not industrially viable, was still being investigated nearly a century after the Haber-Bosch process was developed. The energetics of the gasphase mechanism used in this lab weren’t fully mapped until 2003 (18). Reaction Mechanism In their computational study of several potential gas-phase mechanisms, Hwang and Mebel identified the following four-step pathway as the most favorable (18):

The overall reaction is

Conveniently, this mechanism does not include any ions or radical species. The reaction coordinate diagram includes nine distinct stages: five minima (reactants, intermediates, or products) and four maxima (transition states). 22 Grushow and Reeves; Using Computational Methods To Teach Chemical Principles ACS Symposium Series; American Chemical Society: Washington, DC, 2019.

Experimental Procedure The calculations for this lab follow the general procedure shown in Figure 1 and are completed within a single four-hour laboratory period.

Figure 1. Simplified workflow for calculations. To cover all species included in the reaction mechanism, students must perform six geometry optimizations (N2, H2, NNH2, H2NNH2, HNNH3, and NH3), four transition state optimizations (TS1, TS2, TS3, and TS4), and ten vibrational frequency analyses (all species). Computational Details Our lab computers are equipped with Gaussian 09, which is linked to an installation of WebMO version 16.1 (19, 20). All calculations were done using the PBE0 level of theory and the “Routine: 6-31G(d)” basis set. I selected the combination of theory and basis set that gave results for ΔH°rxn (-90.35 kJ mol-1), ΔS°rxn (-197.77 J mol-1 K-1), and ΔG°rxn (-31.40 kJ mol-1) that were closest to the values computed from standard back-of-the-textbook thermodynamic tables (-91.88 kJ mol-1, 198.11 J mol-1 K-1, and -32.80 kJ mol-1, respectively) (21). A thorough discussion of the shortcomings of density functional theory is beyond the scope of my course so I prioritized matching experimental data when choosing a level of theory. I find it also has the added benefit of increasing students’ confidence in computational methods, which can be helpful if you plan to do additional computational labs later in the course. However, the vast majority of other levels of theory and basis sets will still produce reaction coordinate diagrams that are qualitatively similar to the results shown here, so you can use whatever works best for your class. Optimization of Reactants, Products, and Intermediate Species Students build all stable species (non-transition states) in the WebMO molecule editor. I assign pre-lab work that includes drawing Lewis structures of all reactants, products, and intermediates. Students are then equipped to determine the total charge and spin multiplicity of all species when setting up their WebMO jobs. They all turn out to be neutral molecules in the singlet state, which are 23 Grushow and Reeves; Using Computational Methods To Teach Chemical Principles ACS Symposium Series; American Chemical Society: Washington, DC, 2019.

the default settings for WebMO jobs, so you could omit any discussion of charge or spin multiplicity and have students ignore those options when submitting jobs. WebMO has the option to perform a geometry optimization and frequency analysis in a single job. However, there is not an option for a combined transition state optimization and frequency analysis. For ease of instruction and to avoid confusion, I have students do all structure optimizations separately from frequency analyses. Once the structure optimization completes successfully, students can view the results from the Job Manager and select “New Job Using This Geometry” to continue on to the frequency analysis. Optimization of Transition States There are several options for obtaining the transition state structures. You can adapt the transition state optimizations to suit available computational resources, length of lab period, and course level. Import Coordinates Importing coordinates is the fastest and most frustration-free option for students. You can provide the transition state coordinates in XYZ format or as a Gaussian input file. The XYZ coordinates of all four transition states are provided in Tables 1-4. To avoid a “plug and chug” effect, consider asking students to predict the general structure of the transition state based on the reactants and products; this could be assigned as pre-lab work or as part of an in-class discussion. A transition state optimization should still be performed on each structure (which should finish very quickly) before moving on. Table 1. Transition state 1 XYZ coordinates Atom

X (Å)

Y (Å)

Z (Å)

N

0.000000

0.000000

0.000000

N

0.000000

-1.141238

0.000000

H

0.019343

1.011376

0.000000

H

-1.489645

0.955365

0.000000

Table 2. Transition state 2 XYZ coordinates Atom

X (Å)

Y (Å)

Z (Å)

N

0.000000

0.000000

0.000000

N

1.242512

0.347714

-0.117104

H

1.767858

0.109792

-0.953584

H

1.689979

1.010050

0.512958

H

-0.914832

1.024244

-0.816374

H

-0.730331

1.000708

0.092843

24 Grushow and Reeves; Using Computational Methods To Teach Chemical Principles ACS Symposium Series; American Chemical Society: Washington, DC, 2019.


X (Å)

Y (Å)

Z (Å)

N

0.000000

0.000000

0.000000

N

-1.603532

-0.136164

-0.083991

H

-1.844272

0.863625

-0.139143

H

0.415124

-0.758790

-0.527432

H

0.404776

0.906456

-0.205203

H

-0.556551

-0.186787

0.911039


X (Å)

Y (Å)

Z (Å)

N

0.000000

0.000000

0.000000

N

1.911251

-0.071020

0.060975

H

2.072585

0.925068

-0.018935

H

2.326204

-0.540188

-0.737908

H

2.331892

-0.411310

0.920272

H

-0.053464

-1.022976

0.140798

H

-1.686696

-0.084442

-0.203073

H

-1.316616

0.156347

0.502845

Draw Structures A more advanced option is to have students draw the transition states in the WebMO molecular editor. This requires you to provide images of the transition states with bond lengths, angles, and/or dihedrals labeled. An example of this is shown in Figure 2. This approach is much more successful if students have a moderate amount of experience with the WebMO interface. The optimization of TS4 in particular is very sensitive to dihedral angles, so ideally students are comfortable manipulating molecular structures in three dimensions. Transition State Search The most authentic method for locating the transition states is performing a transition state search, or saddle calculation. This calculation requires two structures on either side of the saddle point to be specified and interpolates between them to produce a structure that is close to the transition state. The resulting structure should then be used as input for a transition state optimization before continuing. I have not yet implemented this procedure in my class, primarily due to concerns about time constraints. Additional considerations to keep in mind are that the two input structures must have the same atomic numbering scheme, so students must be especially careful when creating the inputs. 25 Grushow and Reeves; Using Computational Methods To Teach Chemical Principles ACS Symposium Series; American Chemical Society: Washington, DC, 2019.

Figure 2. Transition state structures with interatomic distances in Angstrom. For some transition states, you can approximate a saddle calculation by performing a coordinate scan. To reach TS1, start with the NNH2 structure and perform a coordinate scan over one of the HN-N bond angles, with a maximum angle of 179.9° (extending to 180° will cause the Gaussian job to fail). Optimize the highest energy structure along the coordinate scan as a transition state. Verifying Transition States Regardless of the method chosen to locate and optimize a transition state, the structures should be verified as such. This can be accomplished by a vibrational frequency analysis or an intrinsic reaction coordinate calculation. Aside from verifying a transition state, the results of the vibrational frequency analyses are the source of all thermodynamic data that students record and use later in their calculations. Vibrational Frequencies For each transition state, I have students perform a vibrational analysis to confirm that the structure represents a saddle point on the potential energy surface. The simplest analysis is to record any imaginary (negative) frequencies in the results. There should not be any such frequencies for reactants, products, or intermediates, but there should be exactly one for a transition state. In particular, TS3 and TS4 can be tricky to locate and my students appreciate the straightforward validation that they “got it right”. WebMO also has the ability to animate each vibrational mode. Animation of the transition state imaginary frequency is a helpful tool for demonstrating the simultaneous bond formation and dissociation that occurs along the reaction coordinate. Students also find it fairly entertaining, especially if the structure resembles a stick figure human in any way. Intrinsic Reaction Coordinate An intrinsic reaction coordinate (IRC calculation in WebMO) takes the input transition state structure and moves along the reaction coordinate. The IRC is calculated in both the forward and reverse directions by default. If the transition state structure is correct, the IRC should arrive at the structure of the expected products when moving in one direction, and should arrive at the structure of the expected reactants when moving in the opposite direction. 26 Grushow and Reeves; Using Computational Methods To Teach Chemical Principles ACS Symposium Series; American Chemical Society: Washington, DC, 2019.

The reaction coordinate path can be animated in WebMO and very clearly demonstrates reactants progressing through the transition state to become products. Students can optimize the structures that result from the forward and reverse IRC calculation and compare their geometries to the expected reactant and product species that they have already optimized individually.

Data Analysis and Results The internal energy, enthalpy, and entropy of each species should be recorded from the vibrational frequencies calculation. I structure the data analysis prompts in the lab handout to guide them through a sensible spreadsheet design; it also happens to make it much easier to track down errors when the intermediate values are tabulated and turned in. While I provide the four-step mechanism shown in Reactions 1-4, I ask students to re-write each reaction stage to be stoichiometrically balanced with the end product (2 NH3). The stoichiometrically balanced mechanism steps are

I ask students to get at least this far in their analysis and check in with me before leaving lab. The most common error is adding too many H2 molecules to the transition state stages, essentially forgetting that one has been “used up” to form the transition state. Once all reaction stages have been balanced correctly, they use their data for the individual species to create a table of calculated thermodynamic quantities, as shown in Table 5. The convention of reporting values relative to the first stage is not mandatory, but can make it easier to get an initial sense of the energy landscape. Table 5. Calculated thermodynamic quantities relative to first reaction stage Reaction Stage

Enthalpy (kJ mol-1)

Entropy (J mol-1 K-1)

Internal Energy (kJ mol-1)

0.000

0.00

0.000

TS1 + 2 H2

500.491

-94.73

502.967

NNH2 + 2 H2

261.581

-98.13

264.057

TS2 + H2

369.571

-213.10

374.528

65.094

-222.34

70.051

TS3 + H2

347.593

-211.07

352.549

HNNH3 + H2

265.740

-216.58

270.694

TS4

354.232

-311.81

361.665

2 NH3

-90.354

-197.77

-85.397

N2 + 3 H2

H2NNH2 + H2


Reaction Coordinate Diagram The reaction coordinate diagram is constructed using the relative internal energy of each reaction stage, as shown in Figure 3.

Figure 3. Reaction coordinate diagram based on internal energy of each reaction stage relative to reactants. I ask students to label each stage in the diagram and discuss the major features of the energy landscape. This part of the analysis helps students to grasp the concept of a potential energy surface and its connection to reaction rate. Students’ discussion of their reaction coordinate diagrams are often phrased in terms of hikers, hills, and valleys. This type of framework helps students to rationalize how something that would provide a less hilly path from reactants to products (a catalyst) would be advantageous. Transition State Theory My students make tables showing the calculated enthalpy of activation (ΔH°‡), entropy of activation (ΔS°‡), and change in molecularity (Δn‡) for each step in the forward and reverse directions. These quantities are shown in Tables 6 and 7, respectively. Table 6. Transition state theory values for forward reactions ΔH°‡ (kJ mol-1)

(J mol-1 K-1)

ΔS°‡

Δn‡

1

500.491

-94.73

-1

2

107.989

-114.96

-1

3

282.499

11.28

0

4

88.492

-95.23

-1

Reaction Step

Table 7. Transition state theory values for reverse reactions ΔH°‡ (kJ mol-1)

(J mol-1 K-1)

ΔS°‡

Δn‡

1

238.910

3.40

0

2

304.477

9.25

0

3

81.853

5.51

0

4

444.586

-114.03

-1

Reaction Step


Interestingly, the change in molecularity has been the most error-prone quantity in students’ lab reports. Even when they’re clear on the concept, students often easily default to “products minus reactants” when computing a delta quantity. Perhaps because the number of reactant and product molecules is easier to visualize, they struggled with this concept more than with the enthalpy or entropy of activation. Activation Energy The activation energy, Ea, of each reaction step can be calculated from the enthalpy of activation and change in molecularity:

where R is the gas constant and T is temperature. The activation energies in the forward direction are 505.449, 112.947, 284.977, and 93.450 kJ/mol. I ask for tabulated values, but the activation energy of each step can also be labeled on the reaction coordinate diagram. Students should be able to easily identify the step with the largest activation energy and predict that this step should also be rate-limiting (i.e., have the smallest forward rate constant). You can also provide students with the mechanism and ask them to predict the rate-limiting step as part of the pre-lab assignment. This provides an opportunity for them to apply their previous knowledge of the stability and bond strength of molecular nitrogen in order to predict that the activation of nitrogen is the rate-limiting step. Rate Constants Forward and reverse rate constants are calculated using the thermodynamic formulation of transition state theory (12):

where kB is the Boltzmann constant, R is the gas constant, T is temperature, h is Planck’s constant, and c° is standard concentration. Calculated values for the rate constant of each step in the forward and reverse direction are shown in Table 8. Table 8. Calculated rate constants in the forward (kf) and reverse (kr) directions Step

kf

kr

1

a1.438 × 10-80

b1.298 × 10-29

2

a7.378 × 10-13

b8.537 × 10-41

3

b7.725 × 10-37

b5.500 × 10-2

4

a2.063 × 10-8

a8.790 × 10-72

a Units of M-1 s-1.

b Units of s-1.

Note that I have conserved step numbers in the forward and reverse directions: Step 1 is always step 1, the only difference is whether it proceeds as N2 + 3 H2 → [TS1]‡ + 2 H2 → NNH2 + 2 H2 (forward) or NNH2 + 2 H2 → [TS1]‡ + 2 H2 → N2 + 3 H2 (reverse). Several students inverted the 29 Grushow and Reeves; Using Computational Methods To Teach Chemical Principles ACS Symposium Series; American Chemical Society: Washington, DC, 2019.

entire mechanism when looking at the reverse direction, so what was originally step 4 was relabeled step 1. This doesn’t affect their raw rate constant values, but the equilibrium constants (ratios of forward and reverse rate constants) are affected. In the next iteration of this lab, I plan to be more explicit about the step numbering scheme in my pre-lab lecture, or perhaps adopt a different labeling scheme (e.g., step A, B, C, and D). Equilibrium Constants Once the system has established a dynamic equilibrium, the equilibrium constant is the ratio of the forward and reverse rate constants:

The forward equilibrium constants are 1.108 × 10-51, 8.643 × 1027, 1.405 × 10-35, and 2.348 × 1063. The relatively large equilibrium constant of step 2 makes for an interesting discussion question. The product of step two is the rocket propellant hydrazine (H2NNH2), which is a comparatively stable intermediate compared to its precursor 1,1-diazine (NNH2). Perceptive students will draw connections between the equilibrium constant values and the relative depths of the valleys in the reaction coordinate diagram by drawing upon their prior knowledge of chemical equilibria.

Conclusions In my experience, the major challenges when creating an experiment using electronic structure calculations are maintaining an authentic experience without being too tedious, and communicating the relevance of the experiment to the students. If students have not been previously exposed to computational methods, they may not immediately view their results as “real data” or fully understand how computational tools contribute to our understanding of chemical systems. This experiment uses realistic electronic structure calculations to study the mechanism of a reaction that is connected to one of the most important industrial processes in modern life. Students use their knowledge of thermodynamics, equilibria, and kinetics to discuss their findings, effectively retracing a path throughout the physical chemistry thermodynamics curriculum.

References 1. 2. 3. 4.

5.

Martini, S. R.; Hartzell, C. J. Integrating Computational Chemistry into a Course in Classical Thermodynamics. J. Chem. Educ. 2015, 92, 1201–1203. Esselman, B. J.; Hill, N. J. Integration of Computational Chemistry into the Undergraduate Organic Chemistry Laboratory Curriculum. J. Chem. Educ. 2016, 93, 932–936. Johnson, L. E.; Engel, T. Integrating Computational Chemistry into the Physical Chemistry Curriculum. J. Chem. Educ. 2011, 88, 569–573. Karpen, M. E.; Henderleiter, J.; Schaertel, S. A. Integrating Computational Chemistry into the Physical Chemistry Laboratory Curriculum: A Wet Lab/Dry Lab Approach. J. Chem. Educ. 2004, 81, 475–477. Halpern, A. M. Computational Studies of Chemical Reactions: The HNC-HCN and CH3NHCH3CN Isomerizations. J. Chem. Educ. 2006, 83, 69–76.


6. 7. 8.

9. 10.

11. 12. 13. 14. 15. 16.

17.

18. 19.

Hessley, R. K. Computational Investigations for Undergraduate Organic Chemistry: Predicting the Mechanism of the Ritter Reaction. J. Chem. Educ. 2000, 77, 202–203. Montgomery, C. D. Factors Affecting Energy Barriers for Pyramidal Inversion in Amines and Phosphines: A Computational Chemistry Lab Exercise. J. Chem. Educ. 2013, 90, 661–664. Csizmar, C. M.; Daniels, J. P.; Davis, L. E.; Hoovis, T. P.; Hammond, K. A.; McDougal, O. M.; Warner, D. L. Modeling SN2 and E2 Reaction Pathways and Other Computational Exercises in the Undergraduate Organic Chemistry Laboratory. J. Chem. Educ. 2013, 90, 1235–1238. Marzzacco, C. J.; Baum, J. C. Computational Chemistry Studies on the Carbene Hydroxymethylene. J. Chem. Educ. 2011, 88, 1667–1671. Galano, A.; Alvarez-Idaboy, J. R.; Vivier-Bunge, A. Computational Quantum Chemistry: A Reliable Tool in the Understanding of Gas-Phase Reactions. J. Chem. Educ. 2006, 83, 481–487. Albrecht, B. Computational Chemistry in the Undergraduate Laboratory: A Mechanistic Study of the Wittig Reaction. J. Chem. Educ. 2014, 91, 2182–2185. Chang, R. Physical Chemistry for the Chemical and Biological Sciences, 3rd ed.; University Science Books: Mill Valley, CA, 2000; pp 476−479. United States Geological Survey. Minerals Information: Nitrogen Statistics and Information. https://minerals.usgs.gov/minerals/pubs/commodity/nitrogen/ (accessed July 22, 2018). Nobelprize.org. Fritz Haber – Facts. https://www.nobelprize.org/nobel_prizes/chemistry/ laureates/1918/haber-facts.html (accessed July 22, 2018). Nobelprize.org. Carl Bosch – Facts. https://www.nobelprize.org/nobel_prizes/chemistry/ laureates/1931/bosch-facts.html (accessed July 22, 2018). Everts, S. Who Was Fritz Haber? Chem. Eng. News 2015, 93 (8), 18–23. http://chemicalweapons.cenmag.org/who-was-the-father-of-chemical-weapons/ (accessed July 22, 2018). Kitano, M.; Inoue, Y.; Sasase, M.; Kishida, K.; Kobayashi, Y.; Nishiyama, K.; Tada, T.; Kawamura, S.; Yokoyama, T.; Hara, M.; Hosono, H. Self-Organized Ruthenium-Barium Core-Shell Nanoparticles on a Mesoporous Calcium Amide Matrix for Efficient LowTemperature Ammonia Synthesis. Angew. Chem. Int. Ed. 2018, 57, 2648–2652. Hwang, D.- Y.; Mebel, A. M. Reaction Mechanism of N2/H2 Conversion to NH3: A Theoretical Study. J. Phys. Chem. A. 2003, 107, 2865–2874. Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, Jr., J. A.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg,


J. J.; Dapprich, S.; Daniels, A. D.; Farkas, Ö.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09, revision C.01; Gaussian, Inc.: Wallingford, CT, 2009. 20. Schmidt, J. R.; Polik, W. F. WebMO, version 16.1; WebMO LLC: Holland, MI, 2013. 21. Ball, D. W. Physical Chemistry, 2nd ed.; Cengage Learning: Stamford, CT, 2015; p 808.


Using Electronic Structure Calculations To Investigate the

Recommend Documents