An Ab Initio Potential Surface for the Reaction N++H,-NH++H BY MARTIN A. GITTINS AND DAVID M. HIRST* Department of Molecular Sciences, University of Warwick, Coventry CV4 7AL MARTYN F. GUEST Science Research Council Atlas Computer Laboratory, Chilton, Didcot, Oxon OX11 OQY AND Received 3rd May, 1976 Ab initio coiifiguration interaction calculations using Dunning’s gaussian basis set are reported for the potential energy surface for the reaction N+ + H2 --t NH+ + H. For the collinear approach of N+ to H2 the 3X- surface has a shallow minimum and the 311 surface is repulsive. For CZu geometries the 3Bz surface is strongly repulsive and the 3A2 surface has a shallow minimum. The 3B1 surface has a deep well which is not adiabatically accessible at low relative energies. An avoided crossing in C, symmetry between the 3A” surfaces correlating with the 3Az and 3B1 surfaces gives an adiabatic route to the deep potential well, but trajectory hopping is thought to be significant except at the lowest relative energies.Preliminary work on the fitting of an analytic function to attractive surfaces is discussed. 1. INTRODUCTION Ion-molecule reactions are of particular interest for several reasons. For the primary ion beam one has a greater degree of control over the translational energy than is the case for a neutral beam. This means that it is possible to study ion- molecule reactions over a much larger range of relative translational energies than is possible for neutral reactions. A second feature of interest is that many ion-molecule reactions show a change in mechanism as the relative translational energy is changed. Potential energy surfaces for ground and excited states are often comparatively close so that, even at the lowest energies, non-adiabatic processes may compete effectively with adiabatic reactions.Finally, because most ion-molecule reactions do not have an energy barrier, reactivity and dynamics can be studied in the absence of the domin- ating influence of the activation barrier. The reaction N+ +H2+NH++ H (1) has been studied recently by Fair and Mahanl in an ion-beam collision apparatus in the range of initial relative energies from 0.79 to 2.8 eV. Above 2 eV asymmetric product distributions were obtained, indicating that the reaction proceeds by a direct mechanism. As the relative energy is lowered, the intensity peak moves towards the centre-of-mass velocity, and the intensity contours become more symmetric.At 0.79 eV the distribution is quite symmetric, suggesting that at this relative energy the reaction proceeds through a long lived complex. It is surprising that this behaviour is only observed below 2 eV because the 3B1 ground state of symmetric NH2+ lies68 AN A B ZNZTZO POTENTIAL SURFACE approximately 6 eV below the energies of reactants and products. This is in contrast to the reaction C++H2+CH++H (2) which shows considerable symmetry up to 3.5 eV,2 despite the fact that the potential well of CH2+ is only 4.4 eV deep. Fair and Mahan interpreted their observations in terms of qualitative correlation diagrams constructed with the aid of ab initio calculations for NH,3 NH+ and NH2+ and the ab initio potential surfaces of Liskow et aL6 for reaction (2).For the col- linear approach of N+(3P) to H2(lZ;) they concluded that the potential surface of lowest energy ("-) would have neither a deep well nor a high potential barrier. For configurations of C2" symmetry, surfaces of 3A2, 3B1, and 3B2 symmetry correlate with the reactants. It was considered unlikely that complex formation would result from the 3B2 surface, which gives highly excited configurations. The 3B1 surface does have a deep well, but this is thought to be inaccessible at low relative energies because of tlie existence of a substantial potential barrier arising from an avoided crossing between two 3B1 surfaces, one correlating with the reactants and the other with N(2D) + H2+(2Z3.The 3A2 surface, on the other hand, was thought to be rather flat. However, in C, symmetry the 3A2 and 3B1 states relax to 3A" symmetry, and a conical intersection gives a route to the potential well which does not involve any substantial activation barrier. The main features of the collinear and C2, cases have been confirmed by prelimin- ary calculations ' using the unrestricted Hartree-Fock method. We report here a more thorough investigation of the potential surfaces thought to be important in this reaction. 2. A B INITIO CONFIGURATION INTERACTION CALCULATIONS The calculations reported here were made using the ATMOL 3 integral and SCF packages and the MUNICH configuration interaction programg*lo as implemented at the Atlas Computer Laboratory with Dunning's l1 (9s5p/4s2p) and (4s/3s) Gaussian basis sets for nitrogen and hydrogen respectively.The MUNICH system employs the bonded function formalism for configuration interaction. In the configuration interaction calculations we included all single and double excitations, with the exception that the lowest orbital was always doubly occupied and that the highest virtual orbital was excluded. For the collinear approach of N+(3P) to H2('Z:) to give 3Z-(NHH)+ we consider the 1293 configurations derived from the configurations and la2 2a2 3a2 ln,l lnyl la2 202 3 d 401 lzx1 lnyl, (3) (4) It was necessary to include root function (4) to describe correctly the dissociation of 3Z- (NHH)+to NH+(4C-) + H because the restricted Hartree-Fock function dis- sociates to NH(3C-) + H+.A limited number of points on the 311 surface were cal- culated using the SCF orbitals from the 3Xc- calculations with the root functions la2 202 3a2 In1 401 (5) and la2 202 30' In1 402 (6) to allow for the virtual orbitals for the 3C- state not being optimum for the 311 state.MARTIN A . GITTINS, ET AL. 69 For geometries of C,, symmetry we considered the 3B2, 3B1 and 3Az surfaces. The 3B2 surface was calculated using the 546 configurations derived from the function la12 2a12 3a12 lb,' 4a11. la12 2a12 3a12 lbll 1b21. la12 2a12 1 bZ2 3a11 lbll (7) (8) (9) For the surface of 3A2 symmetry 505 configurations were derived from the function The 3B1 ground state of (HNH)+ has the configuration and correlates, at large N-H2 distances, with N(2D) + H2+(2Z:g+) whereas the reactants N+ (3P) +.H2('Eg+) correlate with the excited configuration la12 2a12 3aI2 1 bll 4a11. (10) Thus the transition from reactants to 3B1 (HNH)+ proceeds through an avoided intersection between surfaces corresponding to configurations (9) and (10). In order to describe properly the lowest 3B, surface we included the 859 terms derived from configurations (9) and (10). The case of C, symmetry is more complicated. We considered the 2556 configura- tions derived from the functions and Configurations (11), (13) and (14) are the A" equivalents in C, symmetry to functions (8), (9) and (10) respectively and configuration (12) is included to allow for correct dissociation to NH+(4C-) + H. Fig. 1, 2 and 3 illustrate the coordinate systems used for linear, C,, and C, geo- metries respectively.FIG. 1 i" R ( N-MI + H*-*H r (HH) FIG. 2 f" s--:'") H*-r(HHl-) ' FIG. 370 A N AR I N I T I O POTENTIAL SURFACE -54.9 We also considered some geometries of C, symmetry with the configuration illustrated in fig. 4. * N \ - c H H* f-) t---+ 0 '1 '2 FIG. 4 u) 0) ? f -55.0- !? - L Q ). C a, 3. POTENTIAL ENERGY SURFACES Fig. 5 shows a cut through the potential surfaces for the collinear case as N+ approaches H2 with the H-H distance fixed at 1.5 a.u. (1 a.u. = 5.291 77 x m). The 3Xc- surface has no activation barrier, and has a shallow well of depth approxi- - 55.2 I I I I t -55.1 \ mately 1.5 eV; it therefore is unlikely to contribute to complex formation. The 311 surface is repulsive and need not be considered further.Cuts through the three surfaces for C,, geometry with H-H distances of 1.5 and 2.0 a.u. are shown in fig. 6 and 7. Motion on the 3B2 surface clearly does not contri- bute to the dynamics as this surface is initially strongly repulsive. There is evidence for the existence of an avoided crossing with a bound surface of high energy. The 3A2 surface is similar to the collinear 3 X - surface in having no activation barrier and a shallow well of depth of approximately 1.7 eV. Thus the entrance channel for collisions of C,, symmetry is slightly more favourable energetically than is the collinear approach. The 3B1 surface clearly exhibits the avoided crossing between surfaces correlating with N+(3P) + H2(%:) and with N(2D) + H,+(2C:).There is a substantial energy barrier of height approximately equal to 3.6 eV in the entrance channel making theMARTIN A . GITTINS, ET A L . 71 -54.60 -5L.70 -5L.80 Ln 4, 2 b ‘-54.90 B c -c x c -55.00 -55 10 - 55.20 10 2 0 3.0 4.0 5.0 6 0 R (N-MI FIG. 6.-Energy (in hartrees) of C2,(HNH)+ as a function of distance R(N - M) between N and mid-point M of HI - H2 for r(HH) = 1.5 a.u. 3B1 ground state of (NHH)+ inaccessible, in CZU symmetry, to collisions with low relative energy. Fig. 8 shows the variation of energy for the 3B1 and 3A2 states as a function of angle for the case where the N-H distance is 2.0 a.u. The energy of the 3B1 state falls rapidly as the bond angle is increased from 60” to give a well depth of approximately 5.66 eV for a bond angle of 150-1 55”.The potential well is flat with res- pect to variation in bond angle as observed by other worker^.^*^^-^^ The intersections of the 3A, and 3B1 curves in fig. 6 and 7 occur on the repulsive part of the 3A2 potential. It is not until the H-H distance has been stretched to about 3.0 a.u. that the 3B1 and 3A2 curves are attractive at the point of intersection. An avoided crossing in C, symmetry giving a low energy path to the potential well would be expected to occur in this region of the surface when the geometry is distorted slightly. Fig. 9 shows sections of the 13A” and 23At’ &Is surfaces for the configuration of fig. 3 with 8 = 90” and r(HH) = 1.5 and 3.0 a.u. The minima on the 13A” surface are shallow, and the avoided intersections between the two surfaces occur in the repulsive region.Thus, reactive collisions with these sorts of geometries will proceed by a direct mechanism. Calculations with 8 = 120” and r(HH) = 3.0 a.u. gave results similar to those with 8 = 90”. It is for geometries having the configuration of fig. 4 that the avoided crossing between the 3A” surfaces correlating with the 3B1 and 3A2 states in Czu symmetry gives a low energy path to the deep potential well of 3B1 (NH,)’. Fig. 10 shows cuts through the 1 3Att and 2 3A” surfaces for r(HH) = 4.0 a.u. with rl = 2.5 and r2 = 1.5. The avoided crossing now occurs on the attractive portions of the surfaces. Similar potential curves were obtained for r(HH) = 3.0 a.u. (rl = 1.0, r2 = 2.0) and r(HH) =72 -54.60 -54.70 -54.80 IA aJ b- c ;i P c .;-54.9 a c -55.00 -55.1C -55.2( AN A B INITIO POTENTIAL SURFACE I I ~~ ~~ ~~ O 2 .o 3.0 L.0 5.0 6'0 R (N-M)/a.u.FIG. 7.-Energy (in hartrees) of CZD (HNH)+ as a function of distance R(N - M) between N and mid-point M of HI - H2 for r(HH) = 2.0 a.u. -5L.S v) -55.c E aJ - L CI Jz -. 2. a C u k-55.1 -55. z -55.: I I I I I 0 63 90 120 150 180 8 /degrees Frc,. &-Variation of energy (in hartrees) of C2,(HNH+) as a function of angle 0 for r(NH) = 2.0 a.u.MARTIN A . GITTINS, ET AL, 73 -54.90 -55 00 In W L c L 0 L % -. ? c -55.10 I_-__- - .A- I 1 5 2 0 2 5 3 0 R (NHl/a u FIG. 9.-Energy (in hartrees) of C,(HNH)+ for 6 = 90" as a function of R(NH) - r(HH) = 1.5 a.u.; - - - - - - - - r(HH) = 3.0 a.u. -54.0 - 54.9 VI -55.0 23 L \ x 9 5 -55.1 - 55.2 -55.3 .I 1 2 3 L 5 R3/a.u FIG.lO.-Energy (in hartrees) of C,(HNH)+ for perpendicular approach for r1 = 2.5, r2 = 1.5 as a function of R3.74 AN A B ZNZTIO POTENTIAL SURFACE 3.5 a.u. (rl = 1.5, r2 = 2.0). Table 1 summarises the separation between the potential surfaces. It can be seen that the splitting is very small for r = 4.0 a.u. (8.4 kJ mol-l) and increases as r(HH) decreases from 4.0 a.u. Thus, for slightly off-centre perpen- dicular collisions, motion on the 1 3A’r surface gives a route, which does not involve an activation barrier, to the deep potential well. TABLE 1 r(HH)/a.u. R3/a.u. 2.0 2.25 2.5 2.75 3.0 E/eV 0.5 1 .o 2.0 (4 (6) .-DETAILS OF AVOIDED INTERSECTION OF 3A” SURFACES 3 .OO 3.50 4.00 E(13A”) - E(2 3A”)/kJ mol-1 11 8.9 127.2 28.5 54.5 127.8 48.3 8.4 107.7 55.8 241 153.6 92.0 surface hopping probability 0.53 0.90 0.64 0.93 0.73 0.95 Our ab initio potential surfaces thus confirm the qualitative conjectures of Fair and Mahan.’ Neither CZc nor linear geometries give a sufficiently deep potential well to account for the symmetric product distributions they obtained at low relative energies.However, our calculations for geometries of C, symmetry clearly demon- strate the existence of a conical intersection which gives an adiabatic route to a deep well. Thus, trajectories which pass through this intersection and into the potential well will contribute to a symmetric product distribution. In order to get a rough idea of the importance of transitions between the 3A‘‘ surfaces at the avoided crossing, we made estimates of the probability of potential surface hopping as a function of relative energy, using the approximations of Bausch- licher et aZ.15 to the Landau-Zener-Stuckelberg forniulation.16-18 Using their approximations, the hopping probability P for a relative energy E is given by P = e~p[-(E,/E)~-l (15) where with AE, being the separation between the surfaces and AE,” the curvature of the separation AE at the avoided intersection.The results are included in table 1. They are very approximate because we did not explore the region of the avoided inter- sections in sufficient detail to obtain accurate values for the curvature AE,”. Even at low energies, there is a very high probability of hopping from the lower to the upper surface for r(HH) = 4.0 a.u. Thus the majority of reactive trajectories passing through this part of the potential surface will not visit the potential well.For the case of r(HH) = 3.5 a.u., the hopping probabilities are lower and fall off more rapidly as the energy decreases. Surface hopping is clearly still important here, but at lower relative energies a significant proportion of trajectories will remain on the lower surface. These rough estimates of the hopping probabilities give us some understanding of why the direct mechanism giving an asymmetric product distribution persists to very low relative energies, despite the fact that the potential well is so deep.MARTIN A . GITTINS, ET A L . 75 4. ANALYTIC FITTING OF ATTRACTIVE POTENTIAL SURFACES Before an ab initio potential energy surface obtained from a limited number of calculations for a series of geometries can be used for dynamical calculations, it has to be expressed in a form such that either the energy or the gradient of the energy can be obtained over a continuum of nuclear geometries.Many workers have considered this problem and some of their efforts have been reviewed briefly.19 The two basic approaches are either interpolation or fitting the set of ab initio points to some analytic function. Sathymurthy and Raff l9 investigated the use of three dimensional spline inter- polation and concluded, from a comparison with an analytic function, that although the method was not sufficiently accurate to produce a point-by-point match in classical trajectory calculations, the total cross sections, energy partitioning and spatial dis- tribution were in good agreement with those obtained with the analytical potential.Thus spline interpolation methods may be useful for dynamical calculations but in our experience their application is not straightforward. Spline interpolation in three dimensions requires a complete rectilinear grid of points. Therefore, a coordinate system is required in which the entire three dimensional surface is represented in a rectilinear form. Most coordinates do not do this. Those that do turn out to involve the calculation of a large number of points well away from any reaction co- ordinate and careful selection of grid intervals. Finding a suitable analytic function can be a very difficult procedure. Yarkony et aL2* were unable to obtain a satisfactory fit to their surface for the (HF)2 system.The extended London-Eyring-Polyani-Sato21 surface has been successfully used by Polyani and Schreiber 22 to fit the ab initio collinear surface of Bender et aZ.23 for FHH and this function may well be suitable for fitting potential surfaces having an activation barrier. However, the LEPS function is somewhat inflexible, and Bowman and K~pperrnann~~ have proposed a method involving the rotation of a Morse function and cubic spline interpolation which they believe overcomes the deficiencies of the LEPS function. Little seems to have been done for surfaces which do not have an activation barrier. However, Sorbie and Murrell 25 have recently proposed, for stable triatomic species, an analytic potential which is capable of reproducing both the equilibrium and asymptotic properties of the molecule.The parameters in their potential were obtained by fitting it to the experimental force constants. Their function may be equally useful for fitting a set of ab initio points for an attractive surface. The poten- tial for the molecule ABC is of the form (17) V(R1R2R3) = ‘C/AB(R1) + VBC(R2) + vAC(R3) + VABdRlRZR3) where VAB etc. are the asymptotic diatomic potentials and VABC(RiR2R3) is a suitable three-body function for which they proposed the form (1 - tanh 3y2S2)(1 - tanh 3y3S3). (18) A is a constant, P is a polynomial containing up to quartic terms and y i are variable parameters. The potential is expressed in terms of Si which are displacements from the triatomic equilibrium configuration.In the case of NH2+, if in addition to the yi and the coefficients in the polynomial P, A is varied, there are 24 parameters. In our experience these can best be optimized by using Fletcher’s modification26 to the algorithm of Marquardt.”76 AN A B INITIO POTENTIAL SURFACE Using this function we were able to fit 302 points obtained in previous unrestricted Hartree-Fock calculations for this system with a least squares deviation of 21.4 kJ mol-l. While this is not particularly good, it does indicate that the potential of Sorbie and M~rre11~~ is a useful starting point for the fitting of a set of ab initio points for attractive potential energy surfaces. It should be possible to refine the function to obtain analytic surfaces sufficiently good for dynamical calculations.We are grateful to the Science Research Council for a Research Studentship to one of us (M. A. G.) and for the provision of computer time at the Atlas Computer Laboratory. J. A. Fair and B. H. Mahan, J. Chem. Phys., 1975,62, 515. B. H. Mahan and T. M. Sloane, J. Chem. Phys., 1973,59, 5661. H. P. D. Liu and G. Verhaegen, Int. J. Quant. Chem., 1971,5, 103. H. P. D. Liu and G. Verhaegen, J. Chem. Phys., 1970,53, 735. S . Y. Chu, A. K. Q. Siu and E. F. Hayes, J. Amer. Chem. Soc., 1972, 94,2969. ti D. H. Liskow, C. F. Bender and H. F. Schaefer 111, J. Chem. Phys., 1974,61, 2507. ’ M. A. Gittins and D. M. Hirst, Chem. Phys. Letters, 1975, 35, 534. M. F. Guest and V. R. Saunders, ATMOL 3, Science Research Council, Atlas Computer Laboratory, Chilton, Didcot, Oxon. G. H. F. Diercksen and B. T. Sutcliffe, Theor. Chim. Acta, 1974,34, 105 lo M. F. Guest and W. R. Rodwell, ATLAS C.I., Science Research Council, Atlas Computer Laboratory, Chilton, Didcot, Oxon. l1 T. M. Dunning, J. Chem. Phys., 1970,53,2823. l2 S . T. Lee and K. Morokuma, J. Amev. Chem. Sac., 1971,93, 6863. l3 J. F. Harrison and C. W. Eakers, J. Amer. Chem. Soc., 1973,95, 3467 l4 M. J. S. Dewar, R. C. Haddon, W. K. Li, W. Thiel and P. K. Weiner, J. Amer. Chem. Soc., l5 C. W. Bauschlicher, S. V. O’Neil, R. K. Preston, H. F. Schaefer I11 and C. F. Bender, J. Chem. l6 L. D. Landau, Phys. 2. Sowjetunion, 1932, 2,46. l7 C. Zener, Proc. Roy. Sac. A, 1932, 137, 696. l8 E. C. G. Stuckelberg, HeIv. Phys. Actu, 1932, 5, 369. l9 N. Sathymurthy and L. M. Raff, J. Chem. Phys., 1975, 63,464. *O D. Yarkony, S. V. O’Neil, H. F. Schaefer 111, C. P. Baskin and C. F. Bender, J. Chem. Phys., 21 P. J. Kuntz, E. M. Nemeth, J. C. Polyani, S. D. Rosner and C. E. Young, J. Chem. Phys., 1966, 22 J. C. Polyani and J. L. Schreiber, Chem. Phys. Letters, 1974, 29, 319. 23 C. F. Bender, S. V. O’Neil, P. K. Pearson and H. F. Schaefer 111, Science, 1972,176, 1412. 24 J. M. Bowman and A. Kuppermann, Chem. Phys. Letters, 1975,34,523. 25 K. S. Sorbie and J. N. Murrell, Mol. Phys., 1975,29, 1387. 26 R. Fletcher, U.K.A.E.A. Research Group Report, No. AERE-R6799 (1971). 27 D. W. Marquardt, J. SOC. Zndust. AppI. Math., 1963, 11,431. 1975,97,4540. Phys., 1973, 59, 1286. 1974, 60, 855. 44, 1168.