Trajectory Calculations by the Rolling Ball Model G. M. Fernandez, J. A. Sordo, and T. L. Sordo Universidad de Oviedo, Oviedo, Spain
One of the most important problems in chemistry is understanding chemical reaction rates. T o this end, one of the main objectives in chemical kinetics is to make calculations of the rate of a reaction from first principles. Current experimental and theoretical methods of molecular dynamics allow us to understand a t a microscopic level the mechanism of an elementary chemical process and make possible the interpretation of macroscopic kinetics in molecular terms. A chemical reaction is a process in which the interactions among the reactants determine the formation of the products. The nature of these interactions is reflected by the potential energy surface for the system. Thus, all the geometrical rearrangements and energy variations that a supermolecule undergoes in a reactive collision can he represented hy a trajectory on this potential surface. The theoretical methods of molecular dynamics allow us to study these movements of systems on potential surfaces (14).One of the procedures for dealing with this problem is to choose a set of initial conditions and solve the classical equations of motion for the nuclei of the colliding system. For each initial condition one trajectory on the potential surface is obtained. The analysis of these classical trajectories allows us to interoret and discuss enerev ".disoosal . in exoereic reactions, the relative importance of translational and internal barrier. and collision details. e n e-. m in overcomine -a ~otential . In addition, by averaging trajectories over the initial conditions, quantities such as cross sections and rates can be obtained. The objective of this paper is to propose as an exercise in chemical kinetics the construction of a set of programs that perform the calculation and graphical representation of classical trajectories for a triatomic colinear system on a LEPS surface (6,7)by using the model of the rolling ball. According to this model, for a colinear collision one can simulate exactly the actual trajectory by the process of rolling a hall on the corresponding potential energy surface ( I ). This model represents, for the student, a relatively simple access to the subject of molecular trajectories. The Classical Equations of Motion and the Analogy of the Rolling Ball For a colinear triatomic svstem ABC. in which onlv the relative motion is considered, the Hamiltonian is a function of onlv four variahles-two coordinates and two momenta (8). Under these rondirions, the n~orionoirhe nuclei is governed 11sa set of four differential equations
where i = 1,2, the q; are the generalized coordinates, the p; are conjugate momenta, and q i and p i are their time derivatives. Trajectories for the system are determined by integrating numerically these four simultaneous differential equations
\",.
Hamiltonian is formally identical to a monoparticle Hamiltonian
In this way, solving the classical equations of motion for the actual system is identical to solving the equation for the motion of a article on the ootential surface. The Hamiltonian for the triatomic system may he written (10):
H = T +V =H(rm,r~c,i~~,igc)
(6)
where mA, mB, and mc are the masses of atoms A, B, and C; r m and rBc the distances between atoms A and B, and B and C; and ~ A Band i g the ~ corresponding time derivatives. In order to get T in the form of eqn. (3) the new coordinates x and y are defined:
where Cmc
sin 8 = -mn + m~
and
c=
+ mc) (mn(mg mc(mn + me)1
112
(10)
From V = V ( x , y )and substitution of eqns. (7) and (8) into eqn. (4) we have
where
+
m = m ~ ( m mc) ~
(12)
mn + ma + mc Now we have x = p,lm and y = p,lm, thus the only differential equations to solve are
In this paper, given the initial positionand momentum, the numerical integration of eqns. (13) and (14) has been carried out over small time increments, a t , using the equations
19)
An alternative approach to this problem is to simulate the actual trajectory by a ball rolling on the potential energy surface (I).In order to carry out this analogy i t is necessary to transform the coordinates of the actual system so that the Volume 62
Number 6 June 1985
491
The trajectory is given by the equations x(t
+ at) = r ( t ) + @m! a t
(17)
WBC
1+ 3 K
C +-3 + K DBC@BC -eXP cod
Trajectories on a LEPS Surface For a triatomic system ABC the LEPS potential energy expression is (6)
where K is an adjustable constant, the Sato constant. The result (7) expressions for QMNand JMN
DQAC -=--
ax
(-BBC
3+K 2 DACBAC exp (-PAC
(= CY
- i~co))
(29)
(* - y[tgo - L) C O S ~- rAco))
+-1 +23 K DACBAC e x (~- @ A ~ ( ~- y ( t g ~- L) coso - ~ WAC --1+ 3K -= 2 DAC@AC erp (-WAC ax
+-3 +2 K DACBAC e x (-PAC ~
~
~
(30) (.x
- y (QR
(r - y (t& -
z) -
- coso
rAc~))
5) -
r&)) (31) (32) (33)
we can evaluate eqns. (14) and (15) and, finally, obtain through eqns. (16) and (17) the trajectory on the LEPS surface. where DMNis the classical dissociation energy, r i d is the equilihrium bond length and is a constant for the molecule MN. According to eqns. (7) and (a), QMNand JMNcan be exp~essedas functions of x and y, and so, V = V(x 2).Then, the derivative aHlax in eqn. (15) and b H b y in eqn. (16) take the
492
Journal of Chemical Education
Trajectories for the H, System
As an illustration of the capahilities of the set of programs mentioned in this paper we now present some trajectories for the Hg system on a LEPS surface, calculated by the procedure described in the previous sections. The LEPS surface has been calculated by the program LEPS for a value K = 0.15 (7), DHH= 109.4 kcallmol, rHH0 = 0.74 A, and PHH= 1.94 A-l(l1). Trajectories are calculated hy the program TRAJE and represented on the energy surface V = V(rAB, r ~ cby ) the program REPRE in the way shown in Figures 1,2, and 3. Table 1shows the energy intervals associated with each symbol displayed in this graphical representation of the trajectories. A time interval Bt = 0.5 X 10-16s has been taken determining a maximum axillation in the total energy along the trajectories which is always less than 5 X eV. The TRAJE program performs graphical representations of the values of TAB, TBC,and rAc versus t as well, as shown in Fieures 4 and 5. Trajectories have heen chosen for different experimentally ncccssible ranges of energy (9,121. The range of vihri~tional energy runs fn1m0.25to 1.3eV.rorrespondinl: ru the first few vibrational levels of H,. The initial rclnrive velocitv has heen selected from the range 0.5 X lo4 to 2.0 X io4 m-s-I (-0.09-1.38 eV). For all the trajectories calculated, the starting position of the hall is the bottom of the entrance or reactive valley. In Figures 1-5we present a selection of the most interesting trajectories. For the trajectory in Figure 1the relative velocities are i A B = -10,254 m-s-1 and i ~ =c0 m-s-1 (total energy = 0.363 eV). As can be observed, under these conditions the collision between HAand HB-Hc is nonreactive and inelastic. The hall has not translational energy enough to roll over the saddle-point and returns to the starting position with a transfer between translational and vibrational energy. A reactive trajectory results from ~ A B= -3775 m-s-l and f ~ =c -8815 m-s-I, this last value corresponding to the zero-point vibrational energy of Hz (total energy = 0.433 eV) (Fig. 2). This shows plainly the importance of vibrational energy in overcoming a potential harrier. Figure 3 shows the studied reactive trajectory of minimum pure translational energy without vibration of HB-He. For this trajectory i A B = -10,255m-s-I and bc = Om-s-I (total energy = 0.363 eV).
a
)
In Figures 4 and 5 the values of ~ A B rBc, , and rAc versus time are represented. The trajectory corresponding to Figure 4 has been achieved with ~ A B= -7500 m-s-I and rBc = -8815 m-s-' (total energy = 0.692 eV). I t is clearly a reactive collision by a direct mechanism, the lifetime of the collision complex being L E P S SC=.15 43 TO 97
1.21
7.51
1.78
2.32
205
2.59
about 0.2 X 10-l4 s. Note that the collision is inelastic-the product HA-HB molecule has less vibrational energy than the reactant HB-Hc molecule. In this same figure we can see L E P S SC=.15 -43 10 .9i
1.24
1.51
1.78
2.05
2.32
2.69
2.85
RBC
RBC
2.86
WE
Figure 3. Trajectwy for a reactive collision obtained wilh iM = -10255 and igC= 0 rm-'.
ms?
Figure 1. Trajectory for a nonreactive inelastic collision obtained with im = -10254 m-s-' and iBC= 0 m-s-'.
Table 1. L E P S SC=.15 .43 70 97
13-
.
1.24
1;51
1.78
,
2.32
2;05
I
~
,
. ..
2.59
2;86 RBC
8
.
Symbols Used In the Trajectory Representation and Energy Intervals Associated with Them.
SYMBOL
~~~8
2.602907 2.703015 2.803123 GREATER THAN Figure 2. Trajectory for a reactive Mllision obtained with iM = -3775 and iw = -8815 ms-'.
ms-'
The zero of energy is taken to be me limiling energy fa HZ + H (-109 4 Kcal-i-' = -4.744 ev)
Volume 62
Number 6
June 1985
493
Table 2. Initial Velocltles for Some Trajectories
t
zoo-:
.
,8
,
.
3001
. .'
400:
-
:,,
500-
idm-s-')
Total Energy (eV)
Reactive ?
-10255 -10254 -12024 -19012 -20827 -19012 -17005 -12024 -19673 -19674 -20827
0 0 -8503 -8503 -8503 -12024 -12024 -12024 0 0 0
0.363 0.363 1.104 2.059 2.363 2.541 2.207 1.500 1.337 1.337 1.500
yes
no yes yes no no yes Yes YES
no no
Tha starting polnt f a all of mem is me bonom of me reactive valley.
600700 +.2E16
iA8(m-5-')
(S)
Figure 4. Reactive inelastic trajectory obtained wtlh i m = -7500
ma-' and
iE= -8815 ms-'.
-='
.64 $.,
WC=
o
W\C=
,
+
,
,
".89 l l d 1.39 1.64 1.89 2.14 .2.39 2.84 2.69 3.14 3.39
',
jectory in Tahle 2 shows the minimum translational energy required to obtain the products. This energy is about 3.1 X evgreater than the potential harrier (0.360 eV). The last trajectory in Table 2 shows clearly that, for high-energy collisions. nonreactive traiectories are nossihle. All these calculatio~sand graphical representations have been carried out on an H P 9835B desk computer and a 2631A printer with the ahove-mentioned programs written in Basic laneuaee. A Fortran 77 version of these oroerams . .. for a HP lOdOc&puter is also availahle from theauthurs.'l'he I.KI'S 9 X 97 mid and took surface used hi^^ heen constructed with n : about 35 min on the H P 9835B and about 3 m k on the H P 1000. The calculation and printout of each trajectory takes about 15 min on the H P 9835B and 2 min on the H P 1000. Conclusions This paper reports a very simple alternative for trajectory calculations. All types of mathematical complexities have been avoided trying to concentrate attention on chemical features. The cmistruc&n of these programs is, then. a good exercise in physical chemistry yielding a considerable physical insight into thr dvnnmirs of a chemical urocess. In addition. these programs could he used in freshman chemistry courses to get the student acouainted with verv imoortant conceots in chemistry: the s;bjacent collision hec
