Subscriber access provided by NEW MEXICO STATE UNIV
Article
Nonadiabatic Transition State Theory and Trajectory Surface Hopping Dynamics: Intersystem Crossing Between B and A States of SiH 3
1
1
1
2
Ryan R Zaari, and Sergey A Varganov J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/jp509515e • Publication Date (Web): 30 Jan 2015 Downloaded from http://pubs.acs.org on February 1, 2015
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
The Journal of Physical Chemistry A is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 23
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
Nonadiabatic Transition State Theory and Trajectory Surface Hopping Dynamics: Intersystem Crossing Between 3B1 and 1A1 States of SiH2 Ryan R. Zaari and Sergey A. Varganov∗ Department of Chemistry, University of Nevada, Reno, Nevada, U.S.A., 89557 E-mail:
[email protected] Phone: +1 (775) 784-1406
∗
To whom correspondence should be addressed
1
ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Abstract The ability of time-independent nonadiabatic transition state theory (NA-TST) to reproduce intersystem crossing dynamics obtained from the more computationally demanding Tully’s fewest switches trajectory surface hopping method is investigated. The two approaches are applied to the intersystem crossing between the ground 1 A1 state and lowest energy 3 B1 state of SiH2 , coupled through the spin-orbit interaction. For NA-TST, the transition probabilities are calculated using the Landau-Zener formula and the Delos formula which accounts for tunneling. The fewest switches method produces rate constants of 7.6×1011 s−1 for the triplet to singlet transition and 5.2×1011 s−1 for the singlet to triplet transition, using a first-order kinetics model. This corresponds to a triplet electronic state lifetime of 781 fs. The NA-TST predicted rate constants are 1-2 orders of magnitude smaller, leading to a larger triplet state lifetime, as compared with the fewest switches method. This discrepancy cannot be explained by the difference in transition probabilities obtained from NA-TST and fewest switches molecular dynamics, and is believed to be a result of the NA-TST semilocal description of nonadiabatic transitions in the vicinity of the intersystem crossing. Also, the larger triplet state lifetime obtained from NA-TST could be a result of the quantum sampling of rovibrational states, which is missing in classical trajectories traversing the crossing seam.
Keywords spin-orbit coupling, fewest switches, Landau-Zener formula, transition probability
Introduction Providing theoretical methods for physical descriptions of underlying phenomenon occurring during excited state nonadiabatic dynamics is challenging. 1–5 The use of full quantum descriptions is limited to relatively small systems or small time durations. Thus in order to 2
ACS Paragon Plus Environment
Page 2 of 23
Page 3 of 23
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
study larger systems or longer dynamics, semi-classical descriptions to the quantum phenomenon may be necessary. We focus on the specific area of study within nonadiabatic dynamics between states of different spin multiplicity and the resulting dynamics near an intersystem crossing between states. This is in contrast to interactions between states of the same spin occurring near a conical intersection. For the majority of systems, the spin-orbit coupling for states of differing spin is weak and the excitations are thus typically known as spin-forbidden. For this reason the effects of spin-orbit coupled states are frequently neglected. The spin-orbit coupling can become large enough for the spin-forbidden pathway to be substantially important in systems containing heavy nuclei, transition metals within biological systems, 6 or even those without heavy nuclei, such as benzene and cyclobutadiene, 7 and carbon nanotubes. 8 The influence of spin-forbidden processes on reaction and dissociation dynamics has also been studied. 9–13 The nonadiabatic dynamics at an intersystem crossing between two potential energy surfaces (PES) can be described by introducing the concept of a minimum energy crossing point (MECP). The MECP is the lowest energy geometry on the seam of intersection between the two potential energy surfaces. The Landau-Zener formula provides a simple time-independent approach to calculate nonadiabatic transition probabilities between states of different spin for energies at and above the MECP. 14,15 To describe the transition probability for energies below the MECP, a method that includes features of tunneling was developed by Delos. 16 Furthermore, nonadiabatic transition state theory (NA-TST) can be used in conjunction with the transition probability equations to produce rate constants for reactions proceeding through an intersystem crossing. 17,18 In NA-TST, the MECP is treated analogously to the transition state in transition state theory, and the transition probability acts as a transmission coefficient for crossing the transition state. Nonadiabatic transitions can also be described by time-dependent molecular dynamics approaches. One such approach is the direct dynamics Ab Initio Multiple Spawning (AIMS) method 19,20 which is currently implemented 21,22 within the Molpro 23,24 and GAMESS 25,26
3
ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 4 of 23
electronic structure packages. The nuclei are treated as frozen Gaussian wavepackets propagating on an electronic PES as calculated on-the-fly using an electronic structure method. As the wavepacket propagates it can separate (spawn) depending on the nonadiabatic coupling between the electronic states. Another approach is the Tully’s fewest switches (FS) trajectory surface hopping method. 27 In this method, the nuclear dynamics is described by the classical motion of the nuclei on an electronic potential energy surface. The probability of hopping is determined by the coupling between electronic states (spin-orbit coupling or nonadiabatic coupling between nuclear-electronic degrees of freedom) and the relative electronic energy between the two states. A number of additions to the FS approach have been suggested to introduce quantum effects in classical dynamics, such as the time-uncertainty (TU) 28,29 and the gradV (∇V) 30 methods. In this work we investigate how the transition probabilities and rate constants obtained from the more computationally demanding FSTU-∇V nonadiabatic molecular dynamics compare to rate constants and transition probabilities obtained from NA-TST utilizing the analytical Landau-Zener and Delos formulas. As a model, the intersystem crossing between the lowest energy triplet 3 B1 state and the ground singlet 1 A1 state of SiH2 , coupled through the spin-orbit interaction, is used.
Theoretical Methods The simplest method for describing nonadiabatic dynamics is based on the Landau-Zener formula for the probability of transitioning (P LZ ) between two states, 14,15 with the probability of surface hopping (pSH ) expressed here as a double passage, 31 LZ pLZ )(1 + P LZ ) SH = (1 − P
P
LZ
−2π |HSO |2 (E − EM ECP ) = exp ~ |∆G|
4
s
!
µH , 2(E − EM ECP )
ACS Paragon Plus Environment
(1)
Page 5 of 23
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
which describes the probability of hopping along the direction orthogonal to the crossing seam at the MECP. The Landau-Zener transition probability is a function of the spinorbit coupling between the spin states HSO , the reduced mass µH for the degree of freedom orthogonal to the crossing seam, the relative energy of the MECP, EM ECP , and the absolute difference of the gradients |∆G| between the two PES, |∂V3 /∂q − ∂V1 /∂q|. To account for tunneling or transition probabilities for energies less than EM ECP , the transition probability for a weak-coupling double passage is used as previously formulated by Delos, 16 2 2 pDelos SH (E − EM ECP ) = 4π HSO
2µH ~2 G∆G
2/3
Ai2 −(E − EM ECP )
2µH ∆G2 ~2 G
4
!1/3 ,
(2)
where Ai is the Airy function and G is the geometric mean of the norms of the two gradients at the MECP,
q
|∂V3 /∂q| × |∂V1 /∂q|.
The rate constant can be determined by adopting the unimolecular microcanonical form of NA-TST for spin-forbidden reactions with internal energy E, 17,18 Z E 1 ρM ECP (E − EH )pSH (EH − EM ECP )dEH . k(E) = hρR (E) 0
(3)
The density of rovibrational states for the degrees of freedom parallel to the crossing seam ρM ECP (E − EH ) is convoluted with the Landau-Zener or Delos transition probability expressions for the degree of freedom orthogonal to the crossing seam pSH (EH − EM ECP ), for energies (velocities) along the reaction coordinate orthogonal to the crossing seam EH . This resulting ‘integrated density of states’ is divided by the density of rovibrational states for the triplet or singlet state minima, ρR (E). Harmonic frequencies at the MECP for degrees of freedom parallel to the crossing seam were obtained by first constructing an effective Hessian, applying a projection operator that removes the translational, rotational and degree of freedom orthogonal to the crossing seam, and then diagonalizing the resulting Hessian. 18,32 Rovibrational density of states were determined by convoluting the density of vibrational 5
ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 6 of 23
states, calculated through direct counting, 33 with the density of rotational states expressed for an asymmetric top 34 using rotational moments of inertia. The reduced mass µH , in equations (1) and (2), corresponds to the degree of freedom orthogonal to the crossing seam that was projected out of the effective Hessian. In the FSTU-∇V method the probability of hopping from state l to state k is implemented for the case of spin-orbit coupling, 35 FS Pl→k (t)
"
#
2∆t Im(ρkl ) = max 0, HSO . ~ ρll
(4)
The FS method describes both hopping from the triplet to singlet state (l=3, k=1) and FS rehopping from the singlet to triplet state (l=1, k=3). The FSTU-∇V probability, Pl→k (t),
is dependent on the electronic density matrix elements, ρkl and ρll . The off-diagonal elements are defined by ρkl = ck (t)c∗l (t) and the diagonal probabilities as ρll = |cl (t)|2 , with electronic time-dependent coefficients c(t). The coefficients, c(t), are determined by solving the timedependent Schr¨odinger equation with off-diagonal spin-orbit coupling values (HSO ) in the electronic Hamiltonian and wavefunction Ψ = c1 (t)ψ1 (r) + c3 (t)ψ3 (r). A hop occurs if a FS randomly generated number is less than the calculated FSTU-∇V probability (Pl→k ), and
as long as the total energy after the hop does not increase, otherwise the hop is defined as frustrated. In contrast to NA-TST, where nonadiabatic properties are defined at the MECP, the transitions in FSTU-∇V dynamics can occur at any molecular geometry along the trajectory. The fraction of trajectories in the triplet state is given by, F3 (t) = N3 (t)/NT , where N3 (t) is the number of trajectories in the triplet state at any time, t, and NT is the total number of trajectories. The system can be described by first-order kinetics, k3
F3 ⇄ F1
(5)
k1
The solution of the rate equation describing the above process for the fraction of trajectories 6
ACS Paragon Plus Environment
Page 7 of 23
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
in the triplet state, 36 F3 (t), is given by
F3 (t) = (F3 (t0 ) − k1 τ ) exp −
t − t0 + k1 τ, τ
(6)
with initial time t0 and F3 (t0 )=1. The triplet state lifetime is defined in terms of the rate constants as τ =
1 . k1 +k3
Computational Details The MECP geometry, spin-orbit coupling, potential energies and energy gradients required for NA-TST and FSTU-∇V molecular dynamics were calculated using the GAMESS electronic structure package. GAMESS was also used to calculate reactant harmonic frequencies, and MECP harmonic frequencies through application of the effective Hessian projection operator and subsequent diagonalization of the effective Hessian. In order to achieve a balance between accuracy and long trajectory times, the complete-active space self-consistent field (CASSCF) method 37,38 consisting of a full valence active space of 6 electrons in 6 orbitals was used, with the cc-pVTZ basis set. Individual state specific energies and gradients were calculated for the singlet and triplet states. The spin-orbit couplings were determined using the CAS-CI method with singlet-triplet state averaged CASSCF orbitals, and a partial 2-electron and full 1-electron spin-orbit Hamiltonian (HSO2P). 39,40 To perform the molecular dynamics calculations, GAMESS was linked to the direct nonadiabatic trajectory (DiNT) software package 41 which provides the FSTU-∇V method. The molecular dynamics calculations were carried out without any symmetry restrictions on the SiH2 molecule. All integration was carried out using the Runge-Kutta fourth-order method with a time step of 0.5 fs. Initially 96 trajectories were run, all starting on the triplet state PES, but two failed to complete due to SCF convergence issues, and were neglected from the data set. As a consequence, all reported data was taken from the 94 converged trajectories, each with a total propagation time of 2 ps (4000 steps in total). The maximum difference 7
ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
between the initial (t=0) and final (t=2 ps) total energies of the 94 trajectories, indicating the extent of energy conservation, is 50 cm−1 , with an average of 3 cm−1 . The transition probability (eq. (4)) was calculated and a hop occurs for probabilities greater than a randomly generated number between 0 and 1. Hops between the triplet and singlet states can thus occur multiple times (hopping and rehopping) if the conditions allow. The initial conditions (coordinates and momenta) for each of the trajectories were randomly sampled from 3 B1 state-selected normal modes using the harmonic approximation. 41,42 The harmonic frequencies and zero-point energy (ZPE) of each normal mode (bending, symmetric and asymmetric stretching) were determined from a Hessian diagonalization. To generate initial conditions, the ZPE of each vibrational mode was used to determine the corresponding maximum harmonic displacements. An initial distribution of geometries was generated around the triplet state minimum by assuming a cosine distribution of random displacements. The remaining energy, in the form of kinetic energy, is used to generate initial momentum. The true potential energy of each initial geometry is evaluated through an electronic structure calculation at the beginning of each trajectory. The initial triplet state A, whereas the initial bond A to 1.68 ˚ bond lengths, from all 94 trajectories, range from 1.33 ˚ angle ranges from 105◦ to 128◦ .
Results and Discussion Relevant data pertaining to the 1 A1 ground state and 3 B1 excited state, and to the MECP, are contained in Table 1. The total vibrational zero-point energy (ZPE) for the triplet state is 2570 cm−1 . The calculated 1 A1 state and 3 B1 state vibrational frequencies, equilibrium and MECP geometries, are in very close agreement with theoretical 43–49 and experimental 50,51 values reported for SiH2 . Our calculated MECP is located only 816 cm−1 (0.10 eV) above the 3
B1 state minimum. Of the 94 trajectories originating in the 3 B1 state of SiH2 , 75 trajectories
hopped to the 1 A1 state and of those, 48 trajectories eventually rehopped back into the 3 B1
8
ACS Paragon Plus Environment
Page 8 of 23
Page 9 of 23
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
state. Over the 2 ps duration, there were a total of 141 hops from the triplet state and a total of 87 rehops from the singlet state, thus for some trajectories there were multiple hops between the singlet and triplet states. There were 19 trajectories that did not hop at all during the 2 ps simulation. Analysis showed that these trajectories were unable to reach regions of large transition probability as determined by the FSTU-∇V method. Table 1: Geometric parameters (r, Θ), relative energies (∆E), spin-orbit coupling (HSO ) and vibrational harmonic frequencies (ω) for the 1 A1 minimum, 3 B1 minimum and MECP geometry. Energies are relative to the 3 B1 minimum. 1
r (˚ A) Θ (deg) ∆E (eV) ∆E (cm−1 ) HSO (cm−1 ) ω (cm−1 )
3
A1 min. 1.54 94 -0.81 -6533 57.0 997 1972 1978
B1 min. MECP 1.50 1.48 118 133 0.00 0.10 0 816 54.6 50.7 877, 2106 1835, 2254 2157
Figure 1 illustrates the geometries on the 3 B1 potential energy surface where hops to the 1
A1 state occur. In contrast to NA-TST, where the transition probability and rate constant
are calculated using data at the MECP geometry (Figure 1; filled square box), the geometries that result in hops according to FSTU-∇V dynamics are not localized at the MECP. Also shown is the average geometry of all hops (open square box) for the data plotted in Figure 1. The average hopping bond length (1.5 ˚ A) is near the MECP bond length (1.48 ˚ A), but the average hopping bond angle is 120◦ , compared to the MECP bond angle of 133◦ . Figure 2 shows the combined distribution of both bond lengths (r1 , r2 ), the distribution of bond angles and distribution of spin-orbit couplings between the 141 geometries that resulted in hops from the triplet state (dark color) and the 87 geometries that rehopped from the singlet state (light color). The bond lengths are grouped in 0.04 ˚ A increments, the bond angles in 5◦ increments and the spin-orbit couplings in 1 cm−1 increments. According to Figure 2a, the majority of hopping and rehopping occurs for bond lengths in the range A. The majority of hopping and rehopping occurs for bond angles in A to 1.64 ˚ from 1.40 ˚ the range of 105◦ to 125◦ , as shown in Figure 2b. Though the bond length of the MECP 9
ACS Paragon Plus Environment
The Journal of Physical Chemistry
135 1.7
130
1.65
120
1.55 115 1.5 110 1.45
O (degrees)
125
1.6
r2 (Å)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 10 of 23
105 1.4 100
1.35 1.35 1.4 1.45 1.5 1.55 1.6 1.65 1.7
95
r1 (Å)
Figure 1: SiH2 geometries for each occurrence of hopping from the 3 B1 PES. The x- and y-axis denote the bond lengths in ˚ A. The colored z-axis is used to denote the SiH2 bond angle from 95◦ (magenta) to 135◦ (dark blue). Filled circles characterize the 85 hops from the triplet to singlet state. The filled square box marks the MECP geometry and the open square box is the average geometry of all the hops. geometry, r=1.48 ˚ A, is located within the central region of the distribution in Figure 2a, the MECP bond angle of Θ=133◦ is not. The harmonic approximation used in generating initial conditions produces an average bond angle of 116◦ , with the maximum reaching only 128◦ . Thus the MECP region of the PES corresponding to smaller angles is better sampled than the MECP region at larger angles. According to Figure 2c, the majority of hopping and rehopping occurs for spin-orbit couplings between 52 cm−1 and 56 cm−1 . The spin-orbit coupling at the MECP (50.7 cm−1 ) is slightly outside of this range. The maximum difference in spin-orbit couplings associated with hops is less than 10 cm−1 . The minimum and maximum spin-orbit coupling sampled among all trajectories are 46.8 cm−1 and 59.5 cm−1 , respectively. The fraction of all trajectories in the 3 B1 state at any time, F3 (t), during the 2 ps simulation is described in Figure 3. Also plotted is the first-order decay fit to the data according to eq. (6), where the first trajectory hops at t0 =37 fs and F3 (t0 )=1. From the figure nearly 43% of the trajectories are located in the 3 B1 state after 2 ps. From the firstorder decay fit of FSTU-∇V transitions, the predicted rate constants are k3 =7.6×1011 s−1
10
ACS Paragon Plus Environment
Page 11 of 23
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment
Page 12 of 23
Page 13 of 23
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment
The Journal of Physical Chemistry
zero for energies below the first vibrational state, ω = 0.23 eV (1835 cm−1 ). In combination with the convolution integral of eq. (3), this presents a threshold of internal energy for which the rate constant is nonzero. When the Landau-Zener formula is used, the rate constant is nonzero for internal energies corresponding to E > EM ECP +ω, thus for triplet internal energies E > 0.33 eV or singlet internal energies E > 1.14 eV. When the Delos tunneling formula is used, the rate constant is nonzero for internal energies corresponding to E ≥ ω= 0.23 eV. !"#$")
,-#./01#.!231#!456713618!7ï"
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 14 of 23
!"#$"(
!"#$""
!"#$"%
!"#$%'
!"#$%& !%
!%*+
!"
!"*+
!(
!(*+
961#-63:!;6#-
Figure 5: Predicted rate constants, k3 (blue) and k1 (black), for the spin-forbidden transition from the triplet 3 B1 to singlet 1 A1 state of SiH2 , as calculated using the Landau-Zener formula (solid line) and the Delos formula (dashed line). The internal energy corresponding to k3 is relative to the triplet minimum and the internal energy corresponding to k1 is relative to the singlet minimum. The range of internal energies sampled from the 94 trajectories that resulted in hops is used to compare the rate constants obtained using FSTU-∇V surface hopping and NA-TST. The range of internal energies sampled relative to the triplet state minimum is 0.26 eV – 0.38 eV. Relative to the singlet state minimum this corresponds to a range of internal energies of 1.07 eV – 1.19 eV. From rate constant curves in Figure 5, these internal energies corresponds to a k3 range of 1.3×1010 s−1 – 9.3×1010 s−1 and a k1 range of 6.3×109 s−1 – 2.1×1010 s−1 , when using the Delos formula. The corresponding rate constants when using Landau-Zener formula are zero for the minimum energy, though the maximum internal energy has associated rate constants of k3 =1.2×1011 s−1 and k1 =1.2×1010 s−1 . As previously indicated, the rate constants obtained from the first-order decay fit of the FSTU-∇V tra14
ACS Paragon Plus Environment
Page 15 of 23
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
jectory data are k3 =7.6×1011 s−1 and k1 =5.2×1011 s−1 . Thus the lifetime predicted using NA-TST is more than an order of magnitude greater than the calculated FSTU-∇V dynamics lifetime of 781 fs. One may hypothesize that the larger transition probabilities obtained using Landau-Zener or Delos formulas (Figure 4) would result in NA-TST predicting larger rate constants and smaller triplet state lifetime compared to the FSTU-∇V method, but this is contrary to what was found. The order of magnitude smaller triplet state lifetime obtained from FSTU-∇V dynamics could be a result of the classical nature of the trajectories traversing the crossing seam in comparison to the quantum sampling of rovibrational states implemented in NATST. The trajectories encounter a crossing seam that is lower in energy than predicted by NA-TST, within the low energies that are sampled around the ZPE. This could outweigh the effects of transition probability, thus reducing the obtained FSTU-∇V state lifetime substantially.
Conclusion Spin-orbit coupling mediated nonadiabatic transitions between the triplet 3 B1 state and the ground singlet 1 A1 state of SiH2 were studied using nonadiabatic transition state theory and trajectory surface hopping molecular dynamics. The MECP between the states is located 816 cm−1 above the 3 B1 minimum. The transition probabilities were calculated using the Landau-Zener and Delos formulas using properties calculated at the MECP including spinorbit coupling, PES gradients, and the reduced mass associated with the velocity orthogonal to the crossing seam. The rate constants of triplet-singlet and singlet-triplet transitions were obtained using Landau-Zener and Delos transition probabilities in combination with the density of rovibrational states at the MECP, at the triplet and at the singlet minima. Trajectory surface hopping calculations were carried out using Tully’s fewest switches with time uncertainty and gradV (FSTU-∇V) method using 94 trajectories, each 2 ps long,
15
ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 16 of 23
originating in the 3 B1 state. There is a distribution of geometries for which FSTU-∇V hopping occurs. The MECP bond length (1.48 ˚ A) is located near the maximum of the bond length distribution corresponding to FSTU-∇V hops. Whereas, the MECP bond angle (133◦ ) and spin-orbit coupling (50.7 cm−1 ) do not match the corresponding maxima of the FSTU-∇V distributions. This can be attributed to the algorithm used to generate initial conditions and to the position of the intersystem crossing seam between the two PES. The rate constants of triplet-singlet and singlet-triple transitions, and triplet state lifetime, were estimated by fitting the fraction of trajectories in the triplet state to first-order decay kinetics. The FSTU-∇V nonadiabatic dynamics method produces a triplet-singlet rate constant of 7.6×1011 s−1 and a singlet-triplet rate constant of 5.2×1011 s−1 , resulting in a triplet state lifetime of 781 fs. NA-TST predicts rate constants that are an order of magnitude smaller than those obtained using FSTU-∇V dynamics. Consequently, the triplet state lifetime obtained using NA-TST is an order of magnitude greater than obtained using FSTU-∇V. This discrepancy cannot be explained by the difference in transition probabilities obtained from NA-TST and FSTU-∇V molecular dynamics, and could be a result of the NA-TST semilocal description of nonadiabatic transitions in the vicinity of the intersystem crossing. Also, the smaller triplet state lifetime obtained from FSTU-∇V dynamics could be a result of the classical nature of the trajectories traversing the crossing seam in comparison to the quantum sampling of rovibrational states implemented in NA-TST. At the low energies sampled here this effect could result in a large discrepancy in the reported lifetimes between methods.
Acknowledgement The authors are grateful to the University of Nevada, Reno for financial support in the form of S.A.V. startup funding.
16
ACS Paragon Plus Environment
Page 17 of 23
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
References (1) Tully, J. C. Perspective: Nonadiabatic Dynamics Theory. J. Chem. Phys. 2012, 137, 22A301. (2) Barbatti, M. Nonadiabatic Dynamics with Trajectory Surface Hopping Method. WIREs Comp. Molec. Sci. 2011, 1, 620–633. (3) Tully, J. C. Mixed Quantum-classical Dynamics. Faraday Discuss. 1998, 110, 407–419. (4) Yarkony, D. R. Current Issues in Nonadiabatic Chemistry. J. Phys. Chem. 1996, 100, 18612–18628. (5) Persico, M.; Granucci, G. An overview of nonadiabatic dynamics simulations methods, with focus on the direct approach versus the fitting of potential energy surfaces. Theor. Chem. Acc. 2014, 133, 1526. (6) Yson, R. L.; Gilgor, J. L.; Guberman, B. A.; Varganov, S. A. Protein Induced SingletTriplet Quasidegeneracy in the Active Site of [NiFe]-hydrogenase. Chem. Phys. Lett. 2013, 577, 138–141. (7) Penfold, T. J.; Worth, G. A. The Effect of Molecular Distortions on Spin-Orbit Coupling in Simple Hydrocarbons. Chem. Phys. 2010, 375, 58–66. (8) Habenicht, B. F.; Prezhdo, O. V. Ab Initio Time-Domain Study of the Triplet State in a Semiconducting Carbon Nanotube: Intersystem Crossing, Phosphorescence Time, and Line Width. J. Am. Chem. Soc. 2012, 134, 15648–15651. (9) Remmert, S. M.; Banks, S. T.; Harvey, J. N.; Orr-Ewing, A. J.; Clary, D. C. Reduced Dimensionality Spin-Orbit Dynamics of CH3 + HCl ⇋ CH4 + Cl on Ab Initio Surfaces. J. Chem. Phys. 2011, 134, 204311. (10) Favero, L.; Granucci, G.; Persico, M. Dynamics of Acetone Photodissociation: A Surface Hopping Study. Phys. Chem. Chem. Phys. 2013, 15, 20651–20661. 17
ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 18 of 23
´ Ray´on, V. M.; Largo, A. Computational Study of the Reaction of P+ with (11) Cimas, A.; Acetylene: Does Spin-Crossing Play a Significant Role? J. Chem. Phys. A 2012, 116, 3014–3022. (12) Hack, M. D.; Jasper, A. W.; Volobuev, Y. L.; Schwenke, D. W.; Truhlar, D. G. Quantum Mechanical and Quasiclassical Trajectory Surface Hopping Studies of the Electronically ˜ State of NaH2 . J. Chem. Phys. A 1999, 103, Nonadiabatic Predissociation of the A 6309–6326. (13) Besora, M.; Harvey, J. N. Understanding the Rate of Spin-forbidden Thermolysis of HN3 and CH3 N3 . J. Chem. Phys. 2008, 129, 044303. (14) Zener, C. Non-Adiabatic Crossing of Energy Levels. Proc. R. Soc. London A 1932, 137, 696. (15) Wittig, C. The Landau-Zener Formula. J. Phys. Chem. B 2005, 109, 8428–8430. (16) Delos, J. B. On the Reactions of N2 with O. J. Chem. Phys. 1973, 59, 2365. (17) Lorquet, J. C.; Leyh-Nihant, B. Nonadiabatic Unimolecular Reactions. 1. A Statistical Formulation for the Rate Constants. J. Phys. Chem. 1988, 92, 4778–4783. (18) Harvey, J. N.; Aschi, M. Spin-forbidden Dehydrogenation of Methoxy Cation: A Statistical View. Phys. Chem. Chem. Phys. 1999, 1, 5555–5563. (19) Ben-Nun, M.; Mart´ınez, T. J. Nonadiabatic Molecular Dynamics: Validation of the Multiple Spawning Method for a Multidimensional Problem. J. Chem. Phys. 1998, 108, 7244. (20) Ben-Nun, M.; Quenneville, J.; Mart´ınez, T. J. Ab Initio Multiple Spawning: Photochemistry from First Principles Quantum Molecular Dynamics. J. Chem. Phys. A 2000, 104, 5161–5175.
18
ACS Paragon Plus Environment
Page 19 of 23
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
(21) Levine, B. G.; Coe, J. D.; Virshup, A. M.; Mart´ınez, T. J. Implementation of Ab Initio Multiple Spawning in the Molpro Quantum Chemistry Package. Chem. Phys. 2008, 347, 3–16. (22) Gaenko, A.; DeFusco, A.; Varganov, S. A.; Mart´ınez, T. J.; Gordon, M. S. Interfacing the Ab Initio Multiple Spawning Method with Electronic Structure Methods in GAMESS: Photodecay of trans-Azomethane. J. Phys. Chem. A 2014, 118, 10902– 10908. (23) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Sch¨ utz, M.; Celani, P.; Korona, T.; Lindh, R.; Mitrushenkov, A.; Rauhut, G. et al. MOLPRO, version 2012.1, a package of ab initio programs. 2012; see http://www.molpro.net. (24) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Sch¨ utz, M. Molpro: A General Purpose Quantum Chemistry Program Package. WIREs Comput. Mol. Sci. 2012, 2, 242–253. (25) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S. T.; S.Gordon, M.; Jensen, J. H.; Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S. J. et al. General Atomic and Molecular Electronic Structure System. J. Comput. Chem. 1993, 14, 1347–1363. (26) Gordon, M. S.; Schmidt, M. W. In Theory and Applications of Computational Chemistry: The First Forty Years; Dykstra, C. E., Frenking, G., Kim, K. S., Scuseria, G. E., Eds.; Elsevier: Amsterdam, 2005. (27) Tully, J. C. Molecular Dynamics with Electronic Transitions. J. Chem. Phys. 1990, 93, 1061. (28) Jasper, A. W.; Stechmann, S. N.; Truhlar, D. G. Fewest-switches with Time Uncertainty: A Modified Trajectory Surface-hopping Algorithm with Better Accuracy for Classically Forbidden Electronic Transitions. J. Chem. Phys. 2002, 116, 5424.
19
ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
(29) Jasper, A. W.; Stechmann, S. N.; Truhlar, D. G. Erratum: “Fewest-switches with Time Uncertainty: A Modified Trajectory Surface-hopping Algorithm with Better Accuracy for Classically Forbidden Electronic Transitions”. J. Chem. Phys. 2002, 117, 10427. (30) Jasper, A. W.; Truhlar, D. G. Improved Treatment of Momentum at Classically Forbidden Electronic Transitions in Trajectory Surface Hopping Calculations. Chem. Phys. Lett. 2003, 369, 60–67. (31) Harvey, J. N. Understanding the Kinetics of Spin-forbidden Chemical Reactions. Phys. Chem. Chem. Phys. 2007, 9, 331. (32) Miller, W. H.; Handy, N. C.; Adams, J. E. Reaction Path Hamiltonian for Polyatomic Molecules. J. Chem. Phys. 1980, 72, 99. (33) Baer, T.; Hase, W. L. Unimolecular Reaction Dynamics: Theory and Experiments.; Oxford university press: New York (N.Y.), 1995. (34) Green, N. J. B. In Unimolecular Kinetics, Part 1. The Reaction Step.; Green, N. J. B., Ed.; Elsevier Science: Amsterdam, 2003. (35) Granucci, G.; Persico, M.; Spighi, G. Surface Hopping Trajectory Simulations with Spin-Orbit and Dynamical Couplings. J. Chem. Phys. 2012, 137, 22A501. (36) Granucci, G.; Persico, M. Critical Appraisal of the Fewest Switches Algorithm for Surface Hopping. J. Chem. Phys. 2007, 126, 134114. (37) Roos, B. O.; Szulkin, M.; Jaszunski, M. Dynamic Correlation for MCSCF Wave Functions: An Effective Potential Method. Theor. Chim. Acta 1987, 71, 375–384. (38) Ivanic, J.; Ruedenberg, K. Identification of Deadwood in Configuration Spaces Through General Direct Configuration Interaction. Theor. Chem. Acc. 2001, 106, 339–351.
20
ACS Paragon Plus Environment
Page 20 of 23
Page 21 of 23
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
(39) Fedorov, D. G.; Koseki, S.; Schmidt, M. W.; Gordon, M. S. Spin-Orbit Coupling in Molecules: Chemistry Beyond the Adiabatic Approximation. Int. Rev. Phys. Chem. 2003, 22, 551–592. (40) Fedorov, D. G.; Gordon, M. S. A Study of the Relative Importance of One and Twoelectron Contributions to Spin-Orbit Coupling. J. Chem. Phys. 2000, 112, 5611–5623. (41) Jasper, A. W.; Oana, C. M.; Truhlar, D. G. Direct Nonadiabatic Trajectories (DiNT), version 1.1. 2013. (42) Bonhommeau, D.; Valero, R.; Truhlar, D. G.; Jasper, A. W. Coupled-surface Inves˜ Effect of Exciting the Symmetric and tigation of the Photodissociation of NH3 (A): Antisymmetric Stretching Modes. J. Chem. Phys. 2009, 130, 234303. (43) Matsunaga, N.; Koseki, S.; Gordon, M. S. Relativistic Potential Energy Surfaces of XH2 (X=C, Si, Ge, Sn, and Pb) Molecules: Coupling of 1 A1 and 3 B1 States. J. Chem. Phys. 1996, 104, 7988. ˜ 1 A1 , a (44) Yamaguchi, Y.; Huis, T. J. V.; Sherrill, C. D.; Schaefer III, H. F. The X ˜ 1 B1 , ˜ 1 A1 Electronic States of SiH2 . Theor. Chem. Acc. 1997, 97, 341–349. A˜ 1 B1 and B (45) Bauschlicher, C. W.; Langhoff, S. R.; Taylor, P. R. On the 1 A1 -3 B1 Separation in CH2 and SiH2 . J. Chem. Phys. 1987, 87, 387. (46) Balasubramanian, K.; McLean, A. D. The Singlet-Triplet Energy Separation in Silylene. J. Chem. Phys. 1986, 85, 5117. (47) Pettersson, L. G. M.; Siegbahn, P. E. M. Accurate Effective Core Potential for Germanium. Application to the Singlet-Triplet Splitting in GeH2 . Chem. Phys. 1986, 105, 355. (48) Balasubramanian, K. J. Relativistic Configuration Interaction Calculations for Polyatomics: Applications to PbH2 , SnH2 , and GeH2 . J. Chem. Phys. 1988, 89, 5731. 21
ACS Paragon Plus Environment
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
(49) Dyall, K. G. All-electron Molecular Dirac-Hartree-Fock Calculations: Properties of the XH4 and XH2 Molecules and the Reaction Energy XH4 → XH2 +H2 , X=Si, Ge, Sn, Pb. J. Chem. Phys. 1992, 96, 1210. (50) Dubois, I. The Absorption Spectrum of the Free SiH2 Radical. Can. J. Phys. 1968, 46, 2485. (51) McKellar, A. R. W.; Bunker, P. R.; Sears, T. J.; Evenson, K. M.; Saykally, R. J.; Langhoff, S. R. Far Infrared Laser Magnetic Resonance of Singlet Methylene: SingletTriplet Perturbations, Singlet-Triplet Transitions, and the Singlet-Triplet Splitting. J. Chem. Phys. 1983, 79, 5251. (52) Borsella, E.; Caneve, L. Frequency and Time Resolved Luminescence of Intermediate Reaction Products in IR Laser Decomposition of Silane. Appl. Phys. B. 1988, 46, 347– 355. (53) Santos, M.; Diaz, L.; Torresano, J. A.; Pola, J. Visible Luminiscence Study in the Infrared Multiphoton Dissociation of Beta-chloroethenylsilane. J. Photochem. Photobiol., A 1997, 104, 19–23. (54) Thoman Jr., J. W.; Steinfeld, J. I.; McKay, R. I.; Knight, A. E. W. Wide Fluctuations in Fluorescence Lifetimes of Individual Rovibronic Levels in SiH2 (A˜ 1 B1 ). J. Chem. Phys. 1987, 86, 5909. (55) Koseki, S.; Gordon, M. S. Potential Energy Surfaces and Dynamical Properties of Three Low-Lying States of Silylene. J. Molec. Spec. 1987, 123, 392–404.
22
ACS Paragon Plus Environment
Page 22 of 23
Page 23 of 23
The Journal of Physical Chemistry
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment