Information
• Textbooks
• Media
•
Resources edited by
Molecular Modeling Exercises and Experiments
Alan J. Shusterman Reed College Portland, OR 97202-8199
The Inversion Potential of Ammonia: An Intrinsic Reaction Coordinate Calculation for Student Investigation
W
Arthur M. Halpern,* B. R. Ramachandran,† and Eric D. Glendening Department of Chemistry, Indiana State University, Terre Haute, IN 47809; *
[email protected] The inversion of ammonia represents one of the fundamental examples of intramolecular transformations. In this case pyramidal ammonia, in which the N atom is sp3 hybridized, is raised to the higher-energy, planar structure (in which the N atom is sp2 hybridized) and the molecule then reverts back to a pyramidal form, but with the mirror-image (but identical) structure compared with the original form. An equivalent way to describe the inversion process is to imagine the N atom oscillating along the threefold symmetry axis through the plane defined by the three H atoms. A schematic diagram representing NH3 inversion is shown in Figure 1. This inversion, which transforms one pyramidal (C3v) equilibrium structure of NH3 into the other through a planar (D3h) transition state, is represented by a symmetric double minimum potential. This situation was recognized by Dennison and Uhlenbeck in 1932 (1). Computational as well as experimental studies relating to the double minimum potential of NH3 have been presented in this Journal (2, 3). Such double minimum potentials describe the energetics of a number of different molecular processes such as intramolecular proton transfer in malonaldehyde and tropolone (4, 5). The two key features of this potential are the height of the barrier separating the two equilibrium structures and the distance between the minima (or, alternatively, the width of the barrier). These two properties greatly affect the dynamics of interconversion between the two stable species. In the example of NH3, as we will discuss below, the combination of the small masses of the H atoms, along with the relatively low barrier height, results in a rate of quantum mechanical tunneling through the barrier that is an order of magnitude larger than classical thermal inversion over the barrier. Even in the case of tertiary amines, where the larger masses of the N-bonded substituents essentially eliminate tunneling, the thermal inversion rates are sufficiently large at ordinary temperatures to preclude the isolation of enantiomeric forms of NR1R2R3 compounds (6, 7) despite the larger inversion barriers for these substituted amines, for example, 2900 cm᎑1 for trimethylamine (8) versus 2000 cm᎑1 for NH3 (9), respectively. This article describes how students can be empowered to construct the full, double-minimum inversion potential for ammonia by performing intrinsic reaction coordinate (IRC) calculations exclusively from ab initio methods and to obtain both tunneling and thermal inversion rates from the results. This project, although relatively straightforward to carry out, has all the hallmarks of a rigorous quantum chemical calcu-
lation and its applications. This work can be associated with the third-year physical chemistry lecture–laboratory or an upper-level course in computational chemistry. The project can also be expanded to challenge beginning graduate students. We envision that the entire project can be carried out over a two- or three-week period, depending on the students’ background in quantum chemistry and on whether they acquire the IR spectrum themselves or are given it by the instructor. In the first session the instructor can review needed background information, set up the problem, and then allow students to perform the calculations and assemble the IRC potential. In the second session students, with assistance from the instructor, can present and analyze their IRC results, including the prediction of the IR transitions (and the comparison with experiment). On this occasion, a discussion of the tunneling aspects of the project can also be led by the instructor. A valuable component of this project is the opportunity for students to obtain the IR spectrum of NH3 vapor so they can determine for themselves experimentally the 34.4 cm᎑1 splitting in the ν2 fundamental transition (the “inversion” bending mode), which lies between about 900 and 1000 cm᎑1. (Alternatively, the instructor can provide students with the spectrum.) This splitting arises from the coupling between the two equivalent NH3 C3v structures in the double minimum potential and is a manifestation of tunneling between the zeroth-order degenerate eigenstates. Students will compare this value with that predicted from their IRC potential results. Another useful feature of the project is that it can be readily extended to related systems, such as ND3, phosphine,
† Current address: Department of Science & Mathematics, Saint Mary-of-the-Woods College, Saint Mary-of-the-Woods, IN 47876.
Figure 1. Potential energy diagram of NH3 inversion showing the double minimum potential.
www.JCE.DivCHED.org
•
Vol. 84 No. 6 June 2007
•
Journal of Chemical Education
1067
Information •
Textbooks •
Media
•
Resources
and other group 15 hydrides or the isoelectronic species OH3+ and CH3−. Students will learn and appreciate more about quantum chemistry as a result of applying what they have learned in lecture to a basic and interesting problem. Discussion Although students are exposed to the idea of a “reaction coordinate” in second-year organic chemistry (10), they are not usually called upon to consider its quantitative use. For the present purposes, the intrinsic reaction coordinate corresponds to a set of mass-weighted coordinates that describe the minimum energy potential (MEP) between a transition state (or other high-energy, nonequilibrium state) and an equilibrium species.1 Thus the IRC takes into account the displacements of all atoms (i.e., all internal degrees of freedom) of a system as it follows the MEP from any high-energy state (usually, but not necessarily restricted to the transition state) to an equilibrium structure. The IRC was originally described rigorously by Fukui (11), who, along with Hoffmann, received the Nobel Prize in Chemistry in 1981 for his early work in applying the concept of frontier orbitals to chemical reactions. Fortunately, there are a number of quantum chemistry application packages that perform IRC calculations, including the Gaussian 03 suite of programs (12) (and earlier versions) and GAMESS (13). These applications accept an input file that represents a transition-state structure of a molecule (or other constrained high-energy geometry), and produce a set of structures and energies corresponding to the minimum energy path that ends at the equilibrium state of the molecule. The input files needed to perform these calculations using Gaussian 03 for Windows, as well as the operational procedures used to construct the complete IRC potential, are provided in the Supplemental Material.W Because the geometries along the IRC lie relatively far from the equilibrium structure, it is necessary to use a level of theory that describes accurately the sought-after potential surface. We have determined that a good compromise between rigor and computational efficiency is to perform the calculations at the MP2 level of theory using the cc-pVTZ basis set. MP2 is a perturbative method that approximately accounts for electron correlation effects (14), and cc-pVTZ is the correlation-consistent polarized valence triple zeta basis set (15) used to represent the molecular orbitals of NH3. We outline here the general sequence of steps taken to create the inversion potential for NH3. Detailed information about running the calculations, obtaining the needed information, constructing the IRC potential, and obtaining the split energy levels of the inversion mode is provided in the Supplemental Material.W 1. Students optimize the equilibrium and transition-state structures of NH3 and obtain the respective vibrational frequencies. From these results they will immediately obtain the inversion barrier (excluding zero-point energies) and also ascertain the structural properties of these species. They can also compare their calculated frequencies with experiment. Students will confirm that the optimized D3h structure will have only one imaginary frequency (corresponding to the one degree of freedom that has negative curvature) and that this condition is a necessary criterion for a first-order saddle point, or transition state. 1068
Journal of Chemical Education
•
2. The next objective is to create an NH3 structure that represents a high-energy point on the “wing” of the potential. Students can experiment with several partially optimized C3v geometries with different fixed H⫺N⫺H bond angles until they identify a structure has an energy of 2–3 times the inversion barrier (relative to the equilibrium energy). Such a structure will have an H⫺N⫺H bond angle of ca. 80–85⬚. 3. At this point the student can run the IRC calculations. This is done in two parts: one calculation for the relaxation of the D3h transition-state structure to the C3v equilibrium structure (the barrier), and the second for the relaxation of the high-energy C3v structure (obtained from step 2) to the equilibrium state. The run times of all calculations described in steps 1–3 are less than one hour on a reasonably fast PC. 4. Now the student can construct the full IRC potential from the two IRC calculations. The IRC runs carried out in step 3 are pieced together to form one-half of the full IRC potential; this curve is copied, reflected, and then joined with the first curve to create the potential surface. 5. At this point the student can use the IRC potential to obtain the eigenvalues and eigenfunctions associated with the NH3 inversion mode. To achieve this important objective the student can utilize FINDIF, a LabVIEW module that is discussed below. 6. In this last step, the students will evaluate their results by comparing the predicted fundamental and overtone ∼ transition and compare them with the observed valν 2 ∼ funues. Those who acquire the IR spectrum of the ν 2 damental will be able to make this comparison directly. They will also calculate the tunneling rate and, if desired, the thermally-induced inversion rate.
Structures and Barrier To assess the quality of the results of the calculations and gauge accordingly the learning experience students will have, we continue by presenting the results for the optimized equilibrium and transition-state structures of NH3. This information is summarized in Table 1. These energies lead to an inversion barrier of 2078 cm᎑1, which is somewhat higher than the value obtained from higher-level calculations, that is, 1769 cm᎑1 (16). Students should be made aware of the fact that the value of the barrier height depends considerably on the level of the calculation, as well as the basis set used. To illustrate this point, we present results of barrier height calculations obtained from several methods and basis sets in Table 2. In Table 1 we also report the MP2兾cc-pVTZ frequencies of ammonia in its equilibrium and transition states. Because the calculation obtains the second derivative of the energy with respect to nuclear displacements, that is, the force constants, ki, it provides the harmonic frequencies, ∼ νi, for each vibrational mode i: ∼ νi =
1
ki
2πc
µi
(1)
∼ is in cm᎑1 units, c is the speed of light, and µ is in which ν
Vol. 84 No. 6 June 2007
•
www.JCE.DivCHED.org
Information
• Textbooks
• Media
•
Resources
Table 1. Structural and Vibrational Properties of Optimized Structures of NH3 from MP2/cc-pVTZ Calculations Property
NH3 (C3v)
rNH/Å
NH3 (D3h) a,b
0.9937
105.99 [106.67]a,b
(120.0)
᎑56.4529922
᎑56.4435240
1.0113 [1.0124]
∠HNH/deg c
E/Eh
∼ ᎑1 ν1/cm ᎑1 ∼ ν2/cm ∼ν /cm᎑1
3512 (a1) [3337]b 1086 (a1) [950]
878i (a2´´)
b
3899 (e´)
b
1589 (e´)
3
3656 (e) [3444]
4
1686 (e) [1627]
∼ ᎑1 ν /cm a
3678 (a1´)
b
b
Equilibrium values. Experimental values are shown in brackets (17). The energies reported in the calculations are in atomic units (hartrees); energies of spectroscopic quantities are customarily given in cm᎑1 (1 Eh = 2.194746 x 105 cm᎑1).
c
Figure 2. The IRC potential for NH3 obtained from the MP2/cc-pVTZ calculation. The first three pairs of eigenvalues are shown. The 0⫾ and 1⫾ eigenvalues are nearly coincident on this energy scale.
Table 3. Eigenvalues of the IRC Potential of NH3
Table 2. The NH3 Inversion Barrier Obtained with Several Methods and Basis Sets
State
E/cm−1
0+
0545.738
Method
pVDZ
pVTZ
pVQZ
0−
0546.517
MP2
2881
2078
1879
1+
1508.185
1746
1−
1545.736
2038
+
2
2182.838
2−
2480.646
3+
2993.334
3−
3522.933
B3LYP
2506
MP4
3098
CCSD(T)
1946
1852 2248 1937
1867
NOTE: Values are in cm᎑1. The basis sets refer to double, triple, and quadruple zeta sets, respectively.
the reduced mass corresponding to the particular set of nuclear displacements along i. The harmonic frequencies will differ from the experimental fundamental transitions because of anharmonicity, especially for ∼ ν 2, the bending mode, which corresponds to motion along the inversion coordinate. It should be noted that because the potential energy surface of the D3h transition state has negative curvature with respect to ∼ is imaginary. the inversion coordinate (i.e., k2 < 0), ν 2 In Table 1, we compare the bond lengths and bond angles for ammonia obtained from the calculation with the respective parameters based on the equilibrium structure, namely the bottom of the potential wells—as opposed to the zeropoint-based quantities. These comparisons are very favorable. IRC Potential and Vibrational Eigenstates The Schrödinger equation for NH3 inversion can be expressed as −h 2
d
2
8 π c µeff d x
2
+ V (x ) ψ ( x ) = E ψ( x )
(2)
where µeff is the effective reduced mass, V(x) is the potential energy, x is the displacement of the nuclei, and ψ(x) is the www.JCE.DivCHED.org
•
vibrational eigenfunction. Note that V and E (the eigenvalue) are in cm᎑1 units. Because the IRC coordinate, χ, is mass weighted (with units of Å amu1/2), eq 2 can be recast as −h
d
2
2 2
8 π c dχ
+ V (χ) Ψ( χ) = E Ψ( χ)
(3)
and therefore students do not have to specify a value for µeff. Furthermore, the incorporation of µ into χ effectively allows the coordinate to treat variations in the reduced mass across the potential. Figure 2 depicts the calculated IRC potential for ammonia on which are superimposed the first three pairs of split eigenvalues. These inversion eigenvalues are obtained from FINDIF. This module is a Windows-based application that diagonalizes a one-dimensional numerical potential of equallyspaced points. FINDIF provides the eigenvalues and also displays selected eigenfunctions (or their squares) along with the potential being investigated. The eigenvalues associated with the IRC potential shown in Figure 2 are listed in Table 3. These values represent the ∼ ) and can be energy levels of the inversion mode of NH3 (ν 2 used to obtain the tunneling rate between the invertomers for various n⫾ states, as will be shown below.
Vol. 84 No. 6 June 2007
•
Journal of Chemical Education
1069
Information •
Textbooks •
Media
•
Resources
Figure 3. The IR spectrum of NH3 vapor showing the inversion splitting of the ν∼ 2 fundamental transition.
Figure 4. The symmetric (left) and antisymmetric (right) wavefunctions of the n = 0 and n = 1 inversion states.
∼ Transitions IR Spectrum and ν 2 The inversion eigenvalues can be used to predict the splitting in the ∼ ν 2 fundamental and its overtone transitions. This aspect of the project creates the opportunity to include an experimental component to this project in which students ∼ region (ca. obtain the IR spectrum of NH3 vapor in the ν 2 ᎑1 900–1000 cm ). Figure 3 depicts this spectrum. To interpret the IR data and to compare the eigenvalues in Table 3 with literature values, students will have to understand the selection rule that governs the transitions to the symmetry components of the n states from the 0± ground state. For dipole-allowed transitions, only the fundamental (i.e., n = 0 → n = 1) is allowed in the harmonic oscillator approximation, and students can easily observe this transition directly, but they can obtain experimental values of the first three overtones (0 → 2, 3, 4) from the literature (18). A detailed analy∼ fundamental, not sis of the rotational structure of the ν 2 completely shown in Figure 3, has been presented by David (3). Given the symmetric nature of the IRC potential, the wavefunctions are expected to be either symmetric or antisymmetric. Using FINDIF students can immediately confirm this prediction by noting that the lower of the two eigenfunctions of a pair of split states (e.g., 0+) is symmetric while the upper one (i.e., 0−) is antisymmetric. The pairs of eigenfunctions for the n = 0 and 1 states are shown in Figure 4. ∼ fundamental transition as an example, we Using the ν 2 observe that the probability of a dipole-induced transition between the 0⫾ and 1⫾ states is proportional to the square of the transition moment integral,
µ1,0 =
Ψ1± ˆ µ Ψ0 ± dτ
(4)
where the dipole moment operator, µ, is antisymmetric (since it depends on x, y, or z). Thus for µ1,0 to be nonzero, the wavefunctions Ψ1 and Ψ0 must be of opposite symmetry. The general selection rule based on symmetry can be represented as Ψ0+ → Ψn − and Ψ0− → Ψn+
(5)
∼ transitions can be predicted From this result, values of the ν 2 from the eigenvalues in Table 3, and Table 4 contains the ex-
1070
Journal of Chemical Education
•
Table 4. Values of the Splitting in the 0⫾ → n⫿ Transitions in NH3 This calculation n
0− → n+
0+ → n−
splitting
Experiment 0− → n+
0+ → n− splitting
1 0961.17 1000.00 038.83 0932.43 0968.12 035.69 2 1636.32 1934.91 298.59 1597.47 1882.18 284.71 3 2446.82 2977.20 530.38 2384.15 2895.51 511.36 Note: Values are in cm᎑1. The experimental data are from ref 18.
pected frequencies. These predictions are compared with the experimental assignments. We see in Table 3 that the splitting between the 0± levels obtained from the IRC calculation is in remarkable agreement with results from the microwave spectroscopy (i.e., 0.779 cm᎑1 calculated versus 0.793 cm᎑1 experimental; ref 18). Students should not, however, become complacent. This agreement is somewhat fortuitous since that value depends on the computational method and basis set used (see, e.g., Table 1). While the calculated 0⫾ → n⫿ transition energies in Table 4 are all high by ca. 30–60 cm᎑1 the splitting values are all reasonably close to the experimental data. Inversion Tunneling Rate Students who have had quantum chemistry will be aware that tunneling through a barrier can be a rapid process under conditions when the barrier is relatively low and narrow and when the mass of the system is small (19). In the case of NH3, the conditions are favorable for rapid tunneling between the invertomers. From their IRC results students can readily calculate the tunneling rate, ktunn, for the ∼ ν 2 ground state (0⫾) and compare their results with experiment. The tunneling rate is given by (2, 20) (6)
ktunn = 2 c ∆E0 ᎑1
where ∆E0 is the splitting of the n = 0 level in cm and c is the speed of light in cm s᎑1. From this expression and the IRC results (see Table 3) we obtain a value of ktunn = 4.67 × 1010 s᎑1.
Vol. 84 No. 6 June 2007
•
www.JCE.DivCHED.org
Information
The experimental tunneling frequency, which can be directly obtained from the microwave spectrum of NH3 in the 1.3-cm region, is well known and has a value of 2.37 × 1010 s᎑1.2 This value, when expressed as the tunneling rate,3 that is, 4.74 × 1010 s᎑1, is in excellent agreement with the IRC result. Students, thus encouraged by the success of their calculation, will perhaps be motivated to carry out two other calculations of the ammonia inversion rate. The first, known as the WKB approximation (22), is a semiclassical method that often provides good results. Using this approach the tunneling rate is given by the expression
2ν exp(−θ ) π
k WKB =
(7)
where ν is the frequency with which the system encounters the classical turning point at the barrier of the n = 0 level. In this calculation ν may be assigned as twice the average zeropoint energy,4 which, in s᎑1, is ν =
(E
0+
)
(8)
+ E0− c
where the energy values are in cm᎑1 (see Table 3) and c is the speed of light in cm s᎑1. The value of θ is obtained from the integral θ =
2π h
χ>
∫
2 V ( χ) − E 0
dχ
(9)
χ
