The article information
 Haimei Shi, Guanghai Guo, Zhigang Sun
 史海梅, 郭广海, 孙志刚
 Numerical Convergence of the Sinc Discrete Variable Representation for Solving Molecular Vibrational States with a Conical Intersection in Adiabatic Representation
 用SinDVR方法求解在绝热表象下锥形交叉分子振动本征态
 Chinese Journal of Chemical Physics, 2019, 32(3): 333342
 化学物理学报, 2019, 32(3): 333342
 http://dx.doi.org/10.1063/16740068/cjcp1812275

Article history
 Received on: December 6, 2018
 Accepted on: January 14, 2019
b. State Key Laboratory of Molecular Reaction Dynamics and Center for Theoretical Computational Chemistry, Dalian Institute of Chemical Physics, Chinese Academy of Sciences, Dalian 116023, China
Nonadiabatic couplings are ubiquitous in polyatomic molecular excited state dynamics, such as in the photodissociation, photochemsitry, and isomerization processes [15]. With developments of the computational methods and computer resources, more and more attention has been payed to nonadiabatic molecular dynamics involving more than one electronic states [1, 613]. Conical intersection (CI) is one of the important nonadiabatic interactions, where two or more potential energy surfaces (PESs) are degenerate and the nonadiabatic couplings between these states cannot be removed simply in adiabatic representation. In the vicinity of the CIs, the BornOppenheimer (BO) approximation breaks down and the couplings between electronic and nuclear motion become important.
The CIs of multielectronic states play a role as "funnels" which convert rapidly the extra electronic energy into nuclear motions [14, 15]. In the adiabatic picture, the presence of a CI results in two additional operators into the Hamiltonian: one is the geometric phase (GP) and the other is the diagonal BO correction (DBOC) [1618]. The introduction of the GP explains a sign change of adiabatic electronic wave functions when it is conveyed along a closed path of nuclear configurations encircling the CI seam, which makes the electronic wave function double valued [19, 20]. Since the total wave function must be single valued, the double valued electronic wave function renders double valuedness of the nuclear wave function as well. In the adiabatic representation, the GP is associated with nuclear motion, which has important consequences in spectroscopy and dynamics for molecules influenced by the CIs, even in situations when the nuclear wave function is localized far from the region of the CI [3, 18, 2049]. The introduced DBOC term is singular and repulsive around the CI, which accounts for the little populations of the adiabatic states in the region around the CIs.
Recently, it was found that during the photodissociation process the GP leads to an additional phase accumulation for fragments of the nuclear wave packet that transports around the CI to opposite sides [18, 28, 31, 50, 51]. This results in destructive interference, resulting in either spontaneously localizing the wave function or slowing down the nuclear dynamics compared with that in the case where the GP is not considered, which alters the lifetime of the dissociated state much [28, 31, 50, 51]. And, the GP shifts the spectrum of a bound state system by altering the pattern of nodes in the nuclear wave function [16, 22, 40, 4345]. An understanding of the GP effects in a chemical reaction has developed only recently, owing mainly to a series of theoretical investigations and experiments on the hydrogenexchange H+H
For understanding the role of CIs in a molecular dynamics, both in diabatic and adiabatic representation, numerical quantum dynamics investigations are essential. A popular numerical method for solving Schrödinger equation is the discrete variable representation (DVR) method. The DVR method based on classical polynomial or triangular functions using Gauss quadrature exhibits exponential convergence when the PESs are smooth, since essentially it is a spectral or pesudospectral method. However, such a spectral or pesudospectral method converges slowly when there are singularities or cusps, as those arise in a molecule with CIs in adiabatic representation. Thus it is a question that if the DVR method is capable of giving accurate results when it is simply applied with usual coordinates. In the work by Guo and his coworkers, cylindrical coordinates are applied to deal with the singularities and the cusps arising in the molecular Hamiltonian with a CI in the adiabatic presentation, and the DVR method was not applied [18].
In this work, the accuracy of the Sinc DVR method for dealing with the CI in normal coordinates both in adiabatic and diabatic representation is investigated, with and without inclusion of the GP and DBOC terms. It is found that the Sinc DVR method for solving the Schrödinger equation with the CI in diabatic representation gives exponential convergence, as expected. Although the Sinc DVR method cannot give numerical results with exponential convergence for solving the Schrödinger equation with the CI in adiabatic representation, especially for those states of energy around the CI, it does give the results with convergence of high order finite difference method. With the usual grid density of the Sinc DVR method, it is also found that the errors of the calculated vibrational energies in the adiabatic representation, which includes the GP and DBOC operators, are less than that introduced by the function forms of the vector potential (or the "mixing" angle for switching the diabatic and adiabatic representation) accounting for the GP. Since accurate function form of the mixing angle or the vector potential for a polyatomic molecule usually is unknown, in a practical calculation in adiabatic representation for accounting for the GP, an arbitrary function, which is simply capable of making the nuclear wavefunction doubleness when transported along the CI one circle, is applied. We conclude that the Sinc DVR method is accurate enough to investigate the molecular quantum dynamics with CIs in the adiabatic representation.
Ⅱ. THEORYThe wellknown twodimensional (2D) JahnTeller model [18, 19, 31, 55] was adopted in our present work, where the CI between two interacting electronic states is simply a point. The CI model of 2D represents a special case, where adiabatic nonadiabatic couplings can be completely removed by the adiabatictodiabatic transformation. The diabatic Hamiltonian assumes the following form (
$ \begin{equation} \hat{H}^\textrm{d}=\hat{H}_N+\hat{H}_e=\hat{T}_N \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix}+ \begin{pmatrix} V_{11} & V_{12}\\ V_{21} & V_{22} \end{pmatrix} \label{dia} \end{equation} $  (1) 
where
$ \hat{T}_N=\frac{1}{2}\frac{\partial^2}{\partial x^2}\frac{1}{2}\frac{\partial^2}{\partial y^2} $ 
which is the kinetic energy operator while the potential energy operator is a 2
$ \begin{equation} \hat{H}^{\textrm{a}}=\mathbf{U}^{\dagger}\hat{H}^\textrm{d}\mathbf{U}=\hat{T}_N \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix}+\mathbf{U}^{\dagger}\cdot[\hat{T}_N, \mathbf{U}]+ \begin{pmatrix} V_{} & 0\\ 0 & V_{+} \end{pmatrix} \label{ad} \end{equation} $  (2) 
which is obtained by a unitary adiabatictodiabatic transformation
$ \begin{matrix} \mathbf{U}= \begin{pmatrix} \cos\alpha & \sin\alpha\\ \sin\alpha & \cos\alpha \end{pmatrix} \end{matrix} $ 
where
$ \alpha=\frac{1}{2}\arctan\frac{2V_{12}}{V_{22}V_{11}} $ 
and
In the adiabatic representation,
$ V_{\pm}=\frac{1}{2}(V_{11}+V_{22})\pm \frac{1}{2}\sqrt{(V_{11}V_{22})^2+4V{_{12}}^2} $ 
which is diagonal. The nonadiabatic couplings
$ \begin{eqnarray} \mathbf{U}^{\dagger}\cdot[\hat{T}_N, \mathbf{U}]=\begin{matrix} \begin{pmatrix} \hat{\tau}_{11} & i\hat{\tau}_{12}\\ i\hat{\tau}_{21} & \hat{\tau}_{22} \end{pmatrix} \end{matrix}\nonumber \end{eqnarray} $ 
with the DBOC as
$ \hat{\tau}_{11}=\hat{\tau}_{22}=\frac{1}{2}\nabla\alpha\cdot\nabla\alpha $ 
and the derivative couplings as
$ \hat{\tau}_{12}=\hat{\tau}_{21}=\frac{1}{2}\left[(i\nabla)^{\dagger}\cdot\nabla\alpha+\nabla\alpha\cdot(i\nabla)\right] $ 
In the presence of CI, although the adiabatic and diabatic representation are related by a unitary transformation, the diabatic Hamiltonian in Eq.(1) and the adiabatic Hamiltonian in Eq.(3) required different boundary conditions, i.e., singlevalued vs. doublevalued boundary conditions. To impose the doublevalued boundary condition of the adiabatic wave function, the vector potential method proposed by Mead and Truhlar may be adopted [20]. In the method, a positiondependent phase factor
$ \begin{array}{l} \hat H_{{\rm{GP}}}^{\rm{a}} = {{\rm{e}}^{  il\eta }}{{\hat H}^{\rm{a}}}{{\rm{e}}^{il\eta }} = {{\hat T}_N}\left( {\begin{array}{*{20}{c}} 1&0\\ 0&1 \end{array}} \right) + \left( {\begin{array}{*{20}{c}} {{{\hat \zeta }_{11}}}&{i{{\hat \zeta }_{12}}}\\ {  i{{\hat \zeta }_{21}}}&{{{\hat \zeta }_{22}}} \end{array}} \right) + \\ \left( {\begin{array}{*{20}{c}} {{V_  }}&0\\ 0&{{V_ + }} \end{array}} \right) \end{array} $  (3) 
where
$ \begin{eqnarray} \hat{\zeta}_{11} = \hat{\zeta}_{22}\hspace{0.05cm}=\hspace{0.05cm}\hat{\tau}_{11}\hspace{0.05cm}+\hspace{0.05cm}\frac{l^2}{2}\nabla\eta\cdot\nabla\eta\hspace{0.05cm}+\hspace{0.05cm}\frac{l}{2}[(i\nabla)^{\dagger}\cdot\nabla\eta\hspace{0.05cm}+\hspace{0.05cm}\nabla\eta\cdot (i\nabla)] \nonumber \\ \hat{\zeta}_{12} = \hat{\zeta}_{21}\hspace{0.05cm}=\hspace{0.05cm}\hat{\tau}_{12}+\frac{l}{2}(\nabla\alpha\cdot\nabla\eta+\nabla\eta\cdot\nabla\alpha)\nonumber \end{eqnarray} $ 
The form of function
Sometimes the adiabatic excited state can be well ignored, then often only the adiabatic ground state need to be considered. The Hamiltonian operators for the adiabatic ground state can be written as
$ \begin{eqnarray} \hat{H}^{(\textrm{a})}_1 = \hat{T}+V_+\zeta_{11}\nonumber\\ = \hat{T}+V_+\hat{\tau}_{11}+\frac{l^2}{2}\nabla\eta\cdot\nabla\eta+\nonumber\\ \frac{l}{2}[(i\nabla)^{\dagger}\cdot\nabla\eta+\nabla\eta\cdot (i\nabla)] \label{adia1} \end{eqnarray} $  (4) 
where the GP effects are included for
$ \begin{equation} \hat{H}^{(\textrm{a})}_1=\hat{T}+V_+\hat{\tau}_{11} \label{adia1a}\end{equation} $  (5) 
In the following numerical calculations, the results calculated with the above Hamiltonian, i.e., twostate diabatic, onestate and twostate adiabatic with/without the GP and the DBOC operators, are compared with each other to illustrate their roles.
In practical calculations, Hamiltonian operators in the adiabatic presentation involves divergent derivatives of
The Sinc DVR is based on a uniform grid with spacing
$ \begin{eqnarray} S_{j, h}(x) = \frac{{\rm Sinc}(xjh)}{h}\nonumber\\ = \begin{cases} \displaystyle\frac{\sin[\pi(xjh)/h]}{\pi(xjh)/h} & \quad \text{for } x\ne x_j\\ 1 & \quad \text{for } x = x_j \end{cases}\nonumber \end{eqnarray} $ 
The grid interval between
The matrix form of first and second derivative operators on this uniform grid are given by
$ \frac{\partial}{\partial x} = \begin{cases} \displaystyle\frac{(1)^{ij}}{(ij)h} & \quad \text{for } i\ne j\\ 0 & \quad \text{for } i = j \end{cases} $ 
and
$ \frac{\partial^2}{\partial x^2} = \begin{cases} \displaystyle\frac{2(1)^{ij}}{(ij)^2h^2} & \quad \text{for } i\ne j\\ \displaystyle\frac{\pi^2}{3h^2} & \quad \text{for } i = j \end{cases} $ 
The Sinc functions are a set of basis functions orthogonal with unit weight function. With the first and second derivative operators given above, the Schrödinger equation can be written in matrix form directly and solved.
Ⅲ. RESULTS AND DISCUSSIONWe first examine the wellstudied 2D JahnTeller model [18, 19, 31, 55]:
$ V_{11} = \frac{\omega{_1}^2}{2}\left(x+\frac{a}{2}\right)^2+\frac{\omega{_2}^2}{2}y^2 $  (6) 
$ V_{22} = \frac{\omega{_1}^2}{2}\left(x\frac{a}{2}\right)^2+\frac{\omega{_2}^2}{2}y^2\Delta $  (7) 
$ V_{12} = cy $  (8) 
with
The lowest two eigenfunctions and their corresponding eigenenergies are presented in FIG. 2. The plots in the left column are calculated with Eq.(1), the plots in the middle column are calculated with Eq.(5), whereas the plots in the right column are calculated with Eq.(4). In the calculations,
$ \begin{eqnarray} \alpha=\eta = \frac{1}{2}\arctan\frac{2V_{12}}{V_{22}V_{11}}\nonumber\\ = \frac{1}{2}\arctan\left(\frac{2cy}{\Delta+a\omega_1x}\right) \end{eqnarray} $  (9) 
Comparing the plots in the middle and right columns of FIG. 2, we see that the patterns of the eigenfunctions have been altered much with inclusion of the GP effects. There is an extra node for the ground state, but one node is missed for the first excited state. Only with inclusion of the GP effects, the calculated wave functions agree with those calculated in the "exact" diabatic representation.
The convergence of typical vibrational states as a function of number of the grid points (
$ {\rm{Error}} = {\log _{10}}\left( {\frac{{{E_n}  E_n^0}}{{E_n^0}}} \right) $  (10) 
where
For understanding their different convergence behaviours, the eigenfunctions of
For obtaining clearer impression of the numerical convergence of all these kinds of calculations, numerical results with many digits calculated using
Since the CI of model Ⅰ is of low energy, most of the low lying states are influenced by the CI strongly. For a model with CI of much higher energy, one may expect that the numerical convergence for the Sinc DVR would be much better. In the following calculations, we will examine the numerical performance of the Sinc DVR for the model with
Similar to FIG. 2, two typical eigenfunctions and their corresponding eigenenergies are presented in FIG. 6. Similarly, the plots in the left column are calculated with Eq.(1), the plots in the middle column are calculated with Eq.(5), whereas the plots in the right column are calculated with Eq.(4). In the calculations,
Similar to FIG. 3, the convergence of typical vibrational states as a function of number of the grid points (
Quite often, in a practical calculation, only the ground PES in adiabatic representation is known. Thus,
Typical eigenenergies of model Ⅰ and Ⅱ in twostate and onestate adiabatic representation for the GP and the DBOC, with
Finally, we want to emphasise that the GP influences the states of a molecular quantum system much, regardless of the energy of the CI. However, different energies of the CI, do result in different effects on the numerical convergence of the Sinc DVR method. The Sinc DVR method converges fast for the states far from the CI. For the states of energies around the CI, the best choice is to set up the "exact" diabatic representation to describe. Anyway, in the case where the diabatic representation is unavailable, the numerical uncertainty introduced by the Sinc DVR is not larger than that introduced by the arbitrary function forms of the vector potential accounting for the GP effects.
Ⅳ. CONCLUSIONIn this work, we examine the numerical convergence of the Sinc DVR method for calculating vibrational eigenstates of 2D JahnTeller model both in diabatic and adiabatic representation. In the adiabatic representation, models both with inclusion of the GP and DBOC Hamiltonian operators and without inclusion of them were considered. At the same time, the effects of the different function forms of the mixing angle for transforming the diabatic and adiabatic representation, and the effect of the resulted vector potential for including the GP, on the eigenstates in the adiabatic representation are examined. It was found that the Sinc DVR method works well for 2D JahnTeller model not only in diabatic representation, but also in the adiabatic representation with inclusion of the GP and the DBOC, even there is cusp and divergent terms in the Hamiltonian operators. At the same time, different function forms of the mixing angle for transforming the diabatic and adiabatic representation, which also gives the vector potential accounting the GP and the DBOC operator, do have small influence on the eigenstate of the system, which is even larger than the numerical error introduced by the Sinc DVR method with the modest number of grid points.
Ⅴ. ACKNOWLEDGMENTSThis work was supported by the National Natural Science Foundation of China (No.21733006 and No.21825303) and NSFC Center for Chemical Dynamics (No.21688102), the Strategic Priority Research Program of Chinese Academy of Sciences (No.XDB17000000), and the Chinese Academy of Sciences, and the Key Research Program of the Chinese Academy of Sciences.
[1] 
W. Domcke, D. R. Yarkony, and H. Köppel, Conical Intersections: Theory, Computation and Experiment, Singapore: World Scientific, (2011).

[2] 
M. Barbatti, and H. Lischka, J. Am. Chem. Soc.
130, 6831(2008).
DOI:10.1021/ja800589p 
[3] 
M. N. R. Ashfold, B. Cronin, A. L. Devine, R. N. Dixon, and M. G. D. Nix, Science
312, 1637(2006).
DOI:10.1126/science.1125436 
[4] 
W. K. Chen, X. Y. Liu, W. H. Fang, P. O. Dral, and G. L. Cui, J. Phys. Chem. Lett.
9, 6702(2018).
DOI:10.1021/acs.jpclett.8b03026 
[5] 
X. Y. Liu, G. L. Cui, and W. H. Fang, Theo. Chem. Acc.
136, 8(2017).
DOI:10.1007/s002140162029z 
[6] 
J. W. Zwanziger, and E. R. Grant, J. Chem. Phys.
87, 2954(1987).
DOI:10.1063/1.453083 
[7] 
H. Koizumi, and S. Sugano, J. Chem. Phys.
102, 4472(1995).
DOI:10.1063/1.469495 
[8] 
W. Domcke, and G. Stock, Adv. Chem. Phys.
100, 1(1997).

[9] 
G. A. Worth, and L. S. Cederbaum, Ann. Rev. Phys. Chem.
55, 127(2004).
DOI:10.1146/annurev.physchem.55.091602.094335 
[10] 
K. Kleinermanns, N. Nachtigallová, and M.S. de Vries, Int. Rev. Phys. Chem.
32, 1(2013).
DOI:10.1080/0144235X.2012.729331 
[11] 
B. G. Levine, and T. J. Martínez, Ann. Rev. Phys. Chem.
57, 613(2007).

[12] 
D. R. Yarkony, Rev. Mod. Phys.
68, 985(1996).
DOI:10.1103/RevModPhys.68.985 
[13] 
A. V. Akimov, and O. V. Prezhdo, J. Chem. Theo. Comput.
9, 4959(2013).
DOI:10.1021/ct400641n 
[14] 
B. Balzer, S. Hahn, and G. Stock, Chem. Phys. Lett.
379, 351(2003).
DOI:10.1016/j.cplett.2003.08.052 
[15] 
S. Hahn, and G. Stock, J. Phys. Chem. B
104, 1146(2000).
DOI:10.1021/jp992939g 
[16] 
A. F. Izmaylov, J. R. Li, and L. JoubertDoriol, J. Chem. Theo. Comput.
12, 5278(2016).
DOI:10.1021/acs.jctc.6b00760 
[17] 
S. Mahapatra, H. Koppel, and L. S. Cederbaum, J. Phys. Chem. A
105, 2321(2001).
DOI:10.1021/jp003784b 
[18] 
C. J. Xie, D. R. Yarkony, and H. Guo, Phys. Rev. A
95, 022104(2017).
DOI:10.1103/PhysRevA.95.022104 
[19] 
H.C. LonguetHiggins, U. Opik, M. H. L. Pryce, and R. A. Sack, Proc. R. Soc. London Ser. A
244, 1(1958).
DOI:10.1098/rspa.1958.0022 
[20] 
C. A. Mead, and D. G. Truhlar, J. Chem. Phys.
70, 2284(1979).
DOI:10.1063/1.437734 
[21] 
V. J. Barclay, C. E. Dateo, I. P. Hamilton, B. Kendrick, R. T. Pack, and D. W. Schwenke, J. Chem. Phys.
103, 3864(1995).
DOI:10.1063/1.470039 
[22] 
C. A. Mead, Rev. Mod. Phys.
64, 51(1992).
DOI:10.1103/RevModPhys.64.51 
[23] 
S. C. Althorpe, J. Chem. Phys.
124, 084105(2006).
DOI:10.1063/1.2161220 
[24] 
Z. G. Lan, W. Domcke, V. Vallet, A. L. Sobolewski, and S. Mahapatra, J. Chem. Phys.
122, 224315(2005).
DOI:10.1063/1.1906218 
[25] 
C. A. Mead, and D. G. Truhlar, J. Chem. Phys.
77, 6090(1982).
DOI:10.1063/1.443853 
[26] 
B. K. Kendrick, J. Hazra, and N. Balakrishnan, Phys. Rev. Lett.
115, 153201(2015).
DOI:10.1103/PhysRevLett.115.153201 
[27] 
B. K. Kendrick, J. Hazra, and N. Balakrishnan, Nat. Comm.
6, 7918(2015).
DOI:10.1038/ncomms8918 
[28] 
F. Bouakline, Chem. Phys.
442, 31(2014).
DOI:10.1016/j.chemphys.2014.02.010 
[29] 
I. G. Ryabinkin, L. JoubertDoriol, and A. F. Izmaylov, J. Chem. Phys.
140, 214116(2014).
DOI:10.1063/1.4881147 
[30] 
L. JoubertDoriol, I. G. Ryabinkin, and A. F. Izmaylov, J. Chem. Phys.
139, 234103(2014).

[31] 
I. G. Ryabinkin, and A. F. Izmaylov, Phys. Rev. Lett.
111, 220406(2013).
DOI:10.1103/PhysRevLett.111.220406 
[32] 
R. N. Dixon, T. A. A. Oliver, and M. N. R. Ashfold, J. Chem. Phys.
134, 194303(2011).
DOI:10.1063/1.3585609 
[33] 
M. G. D. Nix, A. L. Devine, R. N. Dixon, and M. N. R. Ashfold, Chem. Phys. Lett.
463, 305(2008).
DOI:10.1016/j.cplett.2008.08.085 
[34] 
F. Bouakline, S. C. Althorpe, and D. Pelaéz Ruiz, J. Chem. Phys.
128, 124322(2008).
DOI:10.1063/1.2897920 
[35] 
J.C. JuanesMarcos, and S. C. Althorpe, J. Chem. Phys.
122, 204324(2005).
DOI:10.1063/1.1924411 
[36] 
D. Babikov, B. K. Kendrick, P. Zhang, and K. Morokuma, J. Chem. Phys.
122, 044315(2005).
DOI:10.1063/1.1824905 
[37] 
Z. Xu, M. Baer, and A. J. C. Varandas, J. Chem. Phys.
112, 2746(2000).
DOI:10.1063/1.480848 
[38] 
B. K. Kendrick, J. Chem. Phys.
112, 5679(2000).
DOI:10.1063/1.481143 
[39] 
A. J. C. Varandas, and H. G. Yu, J. Chem. Soc. Faraday Trans.
93, 819(1997).
DOI:10.1039/a605777b 
[40] 
B. K. Kendrick, Phys. Rev. Lett.
79, 2431(1997).
DOI:10.1103/PhysRevLett.79.2431 
[41] 
B. K. Kendrick, and R. T. Pack, J. Chem. Phys.
104, 7475(1996).
DOI:10.1063/1.471460 
[42] 
R. Baer, D. M. Charutz, R. Kosloff, and M. Baer, J. Chem. Phys.
105, 9141(1996).
DOI:10.1063/1.472748 
[43] 
B. E. Applegate, T. A. Barckholtz, and T. A. Miller, Chem. Soc. Rev.
32, 38(2003).
DOI:10.1039/A910269H 
[44] 
M. S. Child, Adv. Chem. Phys.
124, 1(2002).

[45] 
H. von Busch, Dev Vas, H. A. Eckel, S. Kasahara, J. Wang, W. Demtröder, P. Sebald, and W. Meyer, Phys. Rev. Lett.
81, 4584(1998).
DOI:10.1103/PhysRevLett.81.4584 
[46] 
F. Bouakline, S. C. Althorpe, P. Larregaray, and L. Bonnet, Mol. Phys.
108, 969(2010).
DOI:10.1080/00268971003610218 
[47] 
B. J. Rao, R. Padmanaban, and S. Mahapatra, Chem. Phys.
333, 135(2007).
DOI:10.1016/j.chemphys.2007.01.012 
[48] 
J.C. JuanesMarcos, S. C. Althorpe, and E. Wrede, J. Chem. Phys.
126, 044317(2007).
DOI:10.1063/1.2430708 
[49] 
J.C. JuanesMarcos, S. C. Althorpe, and E. Wrede, Science
309, 1227(2005).
DOI:10.1126/science.1114890 
[50] 
J. Schön, and H. Köppel, J. Chem. Phys.
103, 9292(1995).
DOI:10.1063/1.469988 
[51] 
C. J. Xie, J. Y. Ma, X. L. Zhu, D. R. Yarkony, D. Q. Xie, and H. Guo, J. Am. Chem. Soc.
138, 7828(2016).
DOI:10.1021/jacs.6b03288 
[52] 
B. K. Kendrick, J. Phys. Chem. A
107, 6739(2003).

[53] 
B. K. Kendrick, J. Hazra, and N. Balakrishnan, New J. Phys.
18, 123020(2016).
DOI:10.1088/13672630/aa4fd2 
[54] 
D. F. Yuan, Y. F. Guan, W. T. Chen, H. L. Zhao, S. R. Yu, C. Luo, Y. X. Tan, T. Xie, X. A. Wang, Z. G. Sun, D. H. Zhang, and X. M. Yang, Science
362, 1289(2018).
DOI:10.1126/science.aav1356 
[55] 
H. Köppel, W. Domcke, and L. S. Cederbaum, Adv. Chem. Phys.
57, 59(1984).

[56] 
D. T. Colbert, and W. H. Miller, J. Chem. Phys.
96, 1982(1992).
DOI:10.1063/1.462100 
[57] 
P. Amore, J. Phys. A: Math. Gen.
39, L349(2006).
DOI:10.1088/03054470/39/22/L01 
[58] 
J. Lo, and B. D. Shizgal, J. Chem. Phys.
125, 194108(2006).
DOI:10.1063/1.2378622 
[59] 
B. Kendrick, and R. T. Pack, J. Chem. Phys.
104, 7502(1996).
DOI:10.1063/1.471461 
[60] 
B. Kendrick, and R. T. Pack, J. Chem. Phys.
106, 3519(1997).
DOI:10.1063/1.473449 
b. 中国科学院大连化学物理研究所，分子反应动力学国家重点实验室和理论计算中心，大连 116023