The ground state dynamics of s-trans-1,3-butadiene cation: An ab initio quantum dynamical study

Journal of Electron Spectroscopy and Related Phenomena 237 (2019) 146899

Contents lists available at ScienceDirect

Journal of Electron Spectroscopy and Related Phenomena journal homepage: www.elsevier.com/locate/elspec

The ground state dynamics of s-trans-1,3-butadiene cation: An ab initio quantum dynamical study

T

Behnam Nikoobakht School of Chemistry, The University of Sydney, Sydney, Australia

ARTICLE INFO

ABSTRACT

Keywords: Nuclear quantum dynamics Wavepacket propagation method Potential energy surfaces Electronic structure calculation Hamiltonian model and ionization potentials

We used high levels of ab initio electronic structure methods and nuclear quantum dynamical approach for studying the ground state dynamics of the s-trans-1,3-butadiene cation. To benchmark the quality of the electronic structure calculations, potential energy surfaces for the energetically lowest lying doublet X2Bg, A2Au and B2Ag states along seven vibrational modes with considering mode S8 (asymmetric stretching C]C bond) were calculated at the MRCI/CAS(10,7) and CAS(8,8) levels of theories. Considering contribution of S8 mode in the nuclear quantum dynamical investigation leads to a more realistic description of the ground state dynamics of strans-1,3-butadiene, which was not addressed before. Our calculation shows that the population transfer takes place in a time scale 50-60 fs. We also calculated the photoelectron spectrum of the molecular system under study. The excellent agreement of spectrum with the experimental one leads us the conclusion that the potential energy surfaces of three lowest cationic X2Bg, A2Au and B2Ag states and time-dependent population analysis of X2Bg were accurately determined.

1. Introduction The importance of linear polyenes has been related to their pivotal roles in understanding of the wide range of biological systems, material science, and also other fields [1–5]. Even though linear polyenes have been subject of experimental and theoretical studies, important aspects remain unanswered even for the smallest systems (such as butadienes and hexatrienes). For instance, we refer to a few selective review articles [6–9] and a few critical problems relevant to the current work. The latter comprises the potential energy surfaces (PESs) relevant to the system dynamics after photoionization, the vibrational structure of the photoelectron spectral bands and the population transfer between the relevant electronic states important in the system dynamics. In the previous work [10], we started an ab initio based quantum dynamical investigation of the vibronic structure of the X2Bg photoelectron spectral band of s-trans-1,3-butadiene (TBD). PESs along the six internal displacement coordinates corresponding to the electronic states X2Bg, A2Au and B2Ag were computed and it was shown that the coupling modes S17, S19 and S20 play a dominant role in the generation of the vibronic structures of the X2Bg photoelectron spectrum. The resulting complex shape of the X2Bg photoelectron spectrum in 9 eV energy range has been investigated employing various number of vibrational modes and compared to experiment. Most likely, due to the limited number of modes (five) included in

our previous study, the population transfer from ground to the excited states was rather incomplete and was found within 30-40 fs. Another relevant mode, which might be interesting to consider its contribution in the dynamics, is S8 (asymmetric stretching C]C bond). This mode is actually a coupling mode between the X2Bg and A2Au states and thus its contribution on both the population transfer from the ground state to the excited states and the first band of the photoelectron spectrum should be investigated. The contribution of the S8 mode was not considered in our previous work [10] and it thus leaves room for finding the better understanding of the ground state dynamics of TBD cation. Thus, our aim is to reach a reliable description of the ground state dynamics of TBD with considering the additional mode S8 (asymmetric stretching C]C bond) in the quantum dynamical calculation. This requires to modify the Hamiltonian model described in our previous work [10]. In this investigation, we therefore continue and extend our previous ab initio quantum dynamical study of the X2Bg photoelectron spectrum of this prototypical short polyene. This comprises the electronic structure calculations as well as dynamical aspects concerning the strongly coupled X2Bg, A2Au and B2Ag doublet states with considering more vibrational modes in our investigation in comparison with our previous work. In more details, the calculation includes (i) accurate investigations of the vertical ionization potentials (IPs) and PESs using various relevant ab initio electronic structure method, e.g. MRCI using different active spaces (ii) the vibronic coupling strength

E-mail address: [email protected]. https://doi.org/10.1016/j.elspec.2019.146899 Received 17 April 2019; Received in revised form 9 October 2019; Accepted 12 October 2019 Available online 28 October 2019 0368-2048/ © 2019 Elsevier B.V. All rights reserved.

Journal of Electron Spectroscopy and Related Phenomena 237 (2019) 146899

B. Nikoobakht

computed between these cationic states along the relevant vibrational modes at the different levels of theories. (iii) performing the quantum dynamical calculation and realizing that it is important to establish the impact of different vibrational modes and vibronic coupling strengths on the photoelectron spectrum of TBD and especially on the ultrafast population transfer between the X2Bg, A2Au and B2Ag states. In this work, we also present the derivation of the kinetic operator using the Gmatrix method, which has not been addressed in more details in our previous works. This paper is organized as follows. The modified Hamiltonian model, which was originally presented in our previous work [10] for describing the photoionization process in TBD, will be discussed in Section 2.1. The computational details including the electronic structure methods and nuclear quantum dynamics treatment are discussed in Section 2. Section 3 is devoted to the discussion on the ab initio electronic structure and the quantum dynamical results. Section 4 concludes with the summary of the main results.

W11 (S ) W12 (S ) W13 (S ) W (S ) = W12 (S ) W22 (S ) W23 (S ) . W13 (S ) W23 (S ) W33 (S )

In Eq. (3), the vector S refers to the relevant internal coordinates defined in Table 1. W11(S), W22(S) and W33(S), which are diagonal elements of potential energy operator, are corresponding to the PESs of the X2Bg, A2Au and B2Ag electronic states, along relevant internal coordinates of Table 1. They read [10],

Wii (S ) = Vii (S3) + Vii (S8) + Vii (S13) + Vii (S17 ) + Vii (S19) + Vii (S20),

i (4)

= 1, 2, 3.

The index i in Vii(Sn) (n = 3, 8, 13, 17, 19 and 20) refers to the electronic state and Vii(Sn) are written as a function of the internal coordinates as follows [10],

Vii (S3) = a1i + b1i S3 + c1i S32 + d1i S33 + e1i S34, Vii (S8) = a4i + b4i S82 + c4i S84, 2 3 4 Vii (S13) = a3i + b3i S13 + c3i S13 + d3i S13 + e3i S13 ,

2. Theoretical methodology

2 4 Vii (S17) = a5i + b5i S17 + c5i S17 ,

Vii (S19) = a6i + b6i cos(S19) + c6i cos(2S19) + d6i cos(3S19) + e6i cos(4S19)

2.1. Construction of the Hamiltonian model

+ f6i cos(5S19),

In this work, we used the Hamiltonian model, which was constructed in Ref. [10] and it could describe the ionization process corresponding to the first band of the photoelectron spectrum of TBD assigned to the X2Bg ⟵ 11Ag transition. Here, we describe briefly the Hamiltonian model introduced in Ref. [10]. This model has been constructed within a diabatic basis including the X2Bg, A2Au and B2Ag electronic states, where the totally symmetric vibrational modes S3 and S13 and the coupling modes S17, S19 and S20 were treated. These modes are important for the simulation of the X2Bg ⟵ 11Ag photoelectron spectrum in TBD [10,11]. Furthermore, in this work, we add the coupling mode S8 (asymmetric stretching C]C bond) to the Hamiltonian model. Originally we also considered the mode S12 but our investigation showed this mode does not play an important role in the evaluation of the first band of the photoelectron spectrum of TBD. Thus, in this work, we did not treat it in our simulation. These symmetry-adapted internal coordinates are defined in Table 1. The Hamiltonian model with dimension 3 × 3 reads [10]

Vii (S20) = a7i + b7i cos(S20) + c7i cos(2S20) + d7i cos(3S20) + e7i cos(4S20) + f7i cos(5S20). (5) All off-diagonal matrix elements in Eq. (3) representing the vibronic coupling among all electronic states involved in the dynamics are [10]

W12 (S8) = W12 (S17) = W13 (S20) = W23 (S19) =

S8 S8, S17 S17, S20 S20,

(6)

S19 S19.

2.2. Quantum dynamical treatment The multi-configuration time-dependent Hartree (MCTDH) approach has been employed to study the nuclear quantum dynamics of TBD [14,15]. The nuclear wavefunction Ψ in the MCTDH, which solves the time-dependent Schrödinger equation, can describe the dynamics of a molecule possessing several electronic states and f nuclear degrees of freedom. The wavefunction Ψ reads [16,17],

(1)

H = TI3 × 3 + W (S ).

In Eq. (1), I is the unit matrix and T refers to the kinetic operator, which is written as a diagonal 3 × 3 matrix, whose elements can be obtained using

TN =

(3)

( )|

(Q1, Q2, …, Qf ) =

,

=1

1 T p Gp, 2

(2)

( ) (Q , 1

where p refers to the vector of momenta conjugate to the symmetry coordinates of Table 1. For the evaluation of the kinetic operator, the Gmatrix technique is used (for more information see Refs. [12,13]). The detailed information on the derivation of G-matrix required for the calculation of the kinetic operator is mentioned in Appendix A. In Eq. (1), W(S) refers to the potential energy operator and reads [10],

Q2, …, Qf ) =

n1

nf

j1 = 1

jf = 1

A j(1 …) jf (t )

f, =1

( , ) (Q j

, t ), (7)

where Q1, Q2, …, Qf are the nuclear coordinates and α refers to the electronic states. In Eq. (7), A j1 … jf (t ) are called the MCTDH expansion

coefficients. The time-dependent basis functions j( , ) (Q , t ) are the socalled single-particle functions (SPFs). In the MCTDH, one can study non-adiabatic systems using the multi-set formulation [16,17]. In this

Table 1 Definition of the symmetry-adapted internal coordinates of TBD. Modes

Characteristic

Definition

Symmetry

S3 S8 S13 S17 S19 S20

C]C stretch Asymmetric stretching C]C bond CeCeC bending Antisymmetric CeCeC bending Disrotatory twisting of CH2 Conrotatory twisting of CH2

Δr4,3 + Δr1,2 Δr4,3 − Δr1,2 Δϕ4,3,2 + Δϕ3,2,1 Δϕ4,3,2 − Δϕ3,2,1 Δτ5,1,2,3 + Δτ6,1,2,7 + Δτ5,1,2,7 + Δτ6,1,2,3 + Δτ2,3,4,9 + Δτ8,3,4,10 + Δτ2,3,4,10 + Δτ8,3,4,9 Δτ5,1,2,3 + Δτ6,1,2,7 + Δτ5,1,2,7 + Δτ6,1,2,3 − Δτ2,3,4,9 − Δτ8,3,4,10 − Δτ2,3,4,10 − Δτ8,3,4,9

ag bu ag bu au bg

2

Journal of Electron Spectroscopy and Related Phenomena 237 (2019) 146899

B. Nikoobakht

formulation, the wavefunction Ψ is expanded in the set {|α〉} of electronic states, and hence different sets of SPFs can be used for each electronic state. The SPFs, j( , ) (Q , t ) , are written by using a set of primitive basis functions of the discrete variable representation (DVR) type, which are time-independent basis functions. In this work, we chose an appropriate HO-DVR (a harmonic oscillator), resulted in an efficient and accurate form of the wavefunction Ψ. From Eq. (7), it is possible to compute an electronic spectrum and electronic populations corresponding to each electronic states. Here we explain briefly how to calculate the electronic spectrum within the MCTDH context. The autocorrelation function C(t) = 〈Ψ(0)|Ψ(t)〉, the overlap between the initial and the propagated wavefunction, is employed to calculate the photoelectron spectrum. Assuming that the initial state Ψ(0) is real and the Hamiltonian is symmetric (H = HT = H†*) leads the following formula for the autocorrelation function [16]

C (t ) =

(t /2)*| (t /2) .

Fig. 1. Definition of atom numbering, bond lengths and angles and torsional angles for TBD.

(8)

This equation indicates that the autocorrelation function can be evaluated from the wavefunction with only half the propagation time. The Fourier transform of the autocorrelation function C(t) within the framework of Fermi's golden rule gives the electronic spectrum,

P (E )

e iEtC (t )dt.

electronic state. For CAS(8,8) it includes two valence bonding π orbitals, two valence antibonding π* orbitals, two highest occupied σ orbitals at the equilibrium geometry of the ground electronic state, and two valence antibonding σ* orbitals. The π orbitals corresponding to the HOMO and HOMO-1 together with the π* orbitals corresponding to the LUMO and LUMO+1 are depicted in Fig. 2. These CAS spaces were also employed in our earlier works [19,20,24]. For evaluation of PESs using distortions along the coordinates S3 and S13, the molecular system under study remains in the ground state equilibrium geometry point group symmetry, which is C2h. For computations of PESs using distortions along the coordinates S19, S20, S17 and S8, we carried out in the point groups C2, Ci, Cs, and Cs, respectively. All MRCI calculations were carried out using the MOLPRO 2008.1 suite of programs [25]. The EOMIP and ADC(3) calculations were carried out by employing the CFOUR [26] and Dirac [27] program packages.

(9)

Based on Eq. (9), we assume a molecule, which is prepared initially in the ground state, is excited vertically by an ideally short laser pulse into the manifold of vibronically interacting, final electronic states. We have employed the relaxation scheme in an imaginary time to obtain the ground state function Ψ(0) [16,18], where the ground state PES of the neutral molecular system is utilized as described in more details in our previous works [19,20]. Having obtained the ground state wave function, it is lifted to the X2Bg ground state of the cation. Since we perform the propagation within the finite time T, the Fourier transformation of the autocorrelation function generates artifacts [16], which is called the Gibbs phenomenon. In order to fix this issue, the autocorrelation function is multiplied by some appropriate damping function cos 2(πt/ 2T) [16]. Additionally, the autocorrelation function is multiplied by an 2 extra Gaussian function e (t / d) allowing to simulate the experimental line broadening, where τd is the damping parameter (dephasing time). This multiplication is equivalent to the convolution of the spectrum with a Gaussian with a full width at half maximum (fwhm) of 4(ln 2)1/ 2 /τd. Thus, the experimental resolution and line broadening effects can be evaluated.

3. Results and discussion 3.1. Potential energy surfaces The PESs along the coordinates introduced in Table 1 were computed at the MRCI/CAS(10,7) and CAS(8,8) levels of theories. In addition to these methods, we employed the EOMIP [28] and ADC(3) [29] methods to calculate the ionization potentials (IPs). In Table 2, we presented the IPs for the three lowest electronic states X2Bg, A2Au and B2Ag in TBD, where the IPs of the previous theoretical and experimental studies were systematically compared with this work. The first IP corresponds to the removing of an electron from the HOMO (see Fig. 2) with considerable contribution from the p orbital of the first carbon in TBD (see Fig. 1 for numbering of atoms in TBD) [10]. The closest IP value for the X2Bg to the experimental value obtained from the RS2C/ CAS(10,7), which is 9.30 eV. For the electronic state A2Au, the best value for the IP obtained from the MRCI/CAS(10,7) and CAS(8,8), which are 11.49 and 11.48 eV, respectively. For the electronic states B2Ag, we observe that the MRCI/CAS(10,7) and CAS(8,8) produce results in a good agreement with the corresponding experimental values. Overall, there is a fair agreement between our results and the previous theoretical and experimental studies as inferred from Table 2. In Figs. 3 and 4, we calculated the PESs along the S3, S8, S13, S17, S19 and S20 coordinates using different levels of theory as indicated in Section 2.3. In these figures, the points and solid lines refer to the ab initio data and the fitted model to data, respectively. Here, we present for instance the corresponding coefficients of the Hamiltonian model evaluated at the MRCI/CAS(10,7) level of theory in Tables 3 and 4. The

2.3. Electronic structure calculations The ground state equilibrium geometry is obtained by employing the MP2/aug-cc-pVTZ level of theory and imposing the C2h point group, which was reported in our earlier work [10,19,20] (see Table 3 of Ref. [19]). The optimized structure of TBD is depicted in Fig. 1, in which we have specified some internal coordinates. Suitable linear combinations of these internal coordinates generate the 24 symmetry-adapted internal coordinates. Here we only use the most relevant coordinates S3, S8, S13, S17, S19 and S20 (see Table 1 and Ref. [21] for a complete list of all coordinates). We employed the multireference configuration interaction (MRCI) method to compute the PESs of the X2Bg, A2Au and B2Ag electronic states along the totally symmetric modes S3 and S13 and as well as the coupling modes S19, S20, S17 and S8. In this calculation, the cc-pCVTZ and cc-pVTZ basis sets were utilized for atoms C and H, respectively [22,23]. In these calculations, two complete active spaces CAS(10,7) and CAS(8,8) have been used. For CAS(10,7), it contains two valence bonding π-orbitals, two valence anti-bonding π*-orbitals and the three highest σ-orbitals at the equilibrium geometry of the ground

3

Journal of Electron Spectroscopy and Related Phenomena 237 (2019) 146899

B. Nikoobakht

Fig. 2. The four π and π* orbitals at the ground state equilibrium geometry of TBD.

coupling between the two states (X2Bg and A2Au), (A2Au and B2Ag) and (X2Bg and A2Ag), respectively and enter the off-diagonal elements in the Hamiltonian model of Eq. (1). The coupling parameters were determined according to the following equation for the adiabatic potentials V ± , [30]

Table 2 Three lowest vertical ionization potentials (in eV) of TBD at various levels of theory. Statea

X2Bg A2Au B2Ag

Exp. [11]

9.26 11.48 12.2

EOMIPCCSD

CAS(10,7)

CAS(8,8)

MRCI

RS2C

MRCI

RS2C

9.13 11.82 12.55

9.02 11.49 12.49

9.30 11.76 12.63

8.99 11.48 12.70

9.18 11.63 13.05

ADC(3)

OVGFb

9.00 11.42 12.51

8.91 – 12.2

V± =

Va (Si ) + Vb (Si ) ± 2

(

Va (Si )

Vb (Si ) 2

)2 +

2S 2 i

i = 8, 17, 19 and 20.

(10) In Table 5, the coupling mode values computed at the MRCI level of theory are almost close to each other especially for S17 and S19 .

a For the carbon and hydrogen atoms, we used the cc-pCVTZ and cc-pVTZ basis sets, respectively. b Taken from Ref. [11].

3.2. Time-dependent electronic population of X2Bg state

root mean square deviation between the model surfaces and ab initio points is 0.004 eV. In Fig. 3, we observe the computed PESs along the totally symmetric modes S3 and S13 show nearly the same behavior for the two different CAS(10,7) and CAS(8,8). It can be seen that the PES corresponding to the electronic ground cationic states X2Bg is almost isolated in the neighborhood of the ground state equilibrium structure but in far away from the origin, the PESs corresponding to the excited states A2Au and B2Ag cross the ground state PES at the specific molecular geometries. The last statement is especially true for PESs along S13. The bu modes S17 and S8, au mode S19 and bg mode S20 lead to a

For investigation of the time-dependent electronic populations of X2Bg state in the diabatic representation, we employed the Hamiltonian model in Section 2.1 and used the MCTDH algorithm (see Section 2.2). Due to the similarity of PESs and vibronic coupling constants computed using the different CAS(10,7) and CAS(8,8), we thus limit ourselves to only the PESs evaluated employing CAS(10,7) in the nuclear quantum dynamics calculation. We performed the wavepacket propagation up to 500 fs, where time step of 0.5 fs was employed, and thus the electronic population of the X2Bg state was computed. The primitive and the timedependent single particle function (SPF) provided in Table 6. The timedependent population of the X2Bg states in the diabatic representation

Fig. 3. The calculated PESs along S3 and S13 at different levels of theories. The red, green and black lines refer to the PESs corresponding to the electronic states X2Bg, A2Au and B2Ag, respectively. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) 4

Journal of Electron Spectroscopy and Related Phenomena 237 (2019) 146899

B. Nikoobakht

Fig. 4. The calculated PESs along S8, S17, S19 and S20 at different level of theories. The red, green and black lines refer to the PESs corresponding to the electronic states X2Bg, A2Au and B2Ag, respectively. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

at the MRCI/CAS(10,7) level of theory is shown in Fig. 5, where we used the 6D and 5D Hamiltonian models with vibrational modes (S3, S8, S13, S17, S19 and S20) and (S3, S13, S17, S19 and S20), respectively. In Fig. 5, the diabatic populations of X2Bg state using the 5D and 6D Hamiltonian models are presented at the MRCI/CAS(10,7) level of theory. The computed population using the 5D Hamiltonian model shows a

decrease in the population on time scale 30-40 fs and the population transfer from the ground states to the excited state is about 8% at the end of propagation time but this population shows a lot of fluctuations in the transfer between the ground state and the two excited states within the propagation time (the red curve in Fig. 5). After inclusion the coupling mode S8 in the 6D Hamiltonian model, 5

Journal of Electron Spectroscopy and Related Phenomena 237 (2019) 146899

B. Nikoobakht

Table 3 Coefficients of polynomial functions obtained by fitting the analytical expressions of Section 2.1 to the ab initio data along S3 and S13. Coefficients are given in eV. S3

X2Bg

A2Au

B2Ag

a b c d e

9.01634 −5.11814 64.2884 −100.440 60.8636

11.4973 −6.91843 60.0202 −95.7983 61.1717

12.4453 3.11358 50.1343 −127.976 102.509

S13

X2Bg

A2Au

B2Ag

a b c d e

9.01634 −5.11814 64.2884 −100.440 60.8636

11.4973 −6.91843 60.0202 −95.7983 61.1717

12.4453 3.11358 50.1343 −127.976 102.509

Table 6 Number of basis functions for the primitive harmonic and SPF basis in the MCTDH calculations. Modes

Primitive basis

SPF basis

S3 S8 S13 S17 S19 S20

70 30 70 70 68 75

9, 6, 9 9, 6, 9 9, 6, 9 9, 6, 9 9,10, 9 9,10, 9

Table 4 The same as Table 3 but for S17, S19 and S20 coordinates. S17

X2Bg

A2Au

B2Ag

a b

9.016 0.002

11.49 0.002

12.44 0.0007

S8

X2Bg

A2Au

B2Ag

a b c

9.016 70.928 6.142

11.49 101.061 −161.938

12.44 45.019 185.410

S19

X2Bg

A2Au

B2Ag

a b c d e f

16.716 −10.830 4.871 −2.451 0.840 −0.129

−90.632 178.749 −117.794 55.685 −17.535 3.013

233.910 −387.454 258.917 −127.555 42.135 −7.524

S20

X2Bg

A2Au

B2Ag

a b c d e f

15.802 −9.137 3.679 −2.013 0.791 −0.110

−26.751 68.089 −46.500 23.165 −7.885 1.380

−130.612 253.493 −170.225 83.390 −28.558 4.9452

Fig. 5. The computed diabatic electronic populations for the X2Bg state. The diabatic populations are shown by red and blue solid lines using the 5D and 6D Hamiltonian models, respectively. For 5D model, modes S3, S13, S19, S20, S17 are treated while for 6D model, the mode S8 is added. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

blue with red curves). In our previous work [10], the contribution of S8 in the dynamics was not considered. 3.3. The X2Bg ⟵ 11Ag photoelectron spectrum For the evaluation of the photoelectron spectrum of TBD, we compute the Fourier transform of the autocorrelation function C(t) which is obtained from overlapping between the initial and the propagated wavefunction as indicated in Section 2.2. In Ref. [11], authors measured a high resolution photoelectron spectrum of TBD within the energy of range of 8–35 eV. They showed that the vibrational structure was observed in several of the photoelectron bands excited with He I radiation especially in the first band of the photoelectron spectrum, which is the focus of the current investigation. In the first band of the experimental spectrum, four progressions with a common spacing of ≈0.2 eV are visible, which fits very well to the S3 mode, combined with the other two modes S13 and S12. In our simulation, these modes are considered together with the coupling modes S17, S19 and S20. As indicated in Section 2.1, we also add the coupling mode S8 to the Hamiltonian model. We checked that the role of the mode S12 was negligible in our nuclear quantum dynamics calculation and we thus did not consider it [10]. For computation of the first band of the photoelectron spectrum, we considered the 5D and 6D Hamiltonian models as indicated in Section 3.2. The outcomes of these two calculations are shown in Fig. 6. In Fig. 6a and b, the low and high resolution photoelectron spectra of TBD using 5D and 6D models, respectively are only shown. In all calculations, the dephasing time 100 fs for the low-resolution spectrum is used, while the dephasing time 1000 fs for the

Table 5 Computed coupling constants relevant to the modes S8, S17, S19 and S20 at different active spaces. For the modes S17, S19 and S20, the unit is eV/deg while for the mode S8, the unit is eV/Å. Coupling constant S8

S17

S19 S20

MRCI/CAS(10,7)

MRCI/CAS(8,8)

6.873

6.901

0.051

0.051

0.071 0.048

0.071 0.051

the behavior of population of ground state changes considerably, reflecting the important role of vibrational mode S8 in the dynamics. One can observe that the population transfer is slightly increased to 10% within the propagation time and the largest initial decay in the population on time scale 50-60 fs can be seen. This initial decay is followed by the oscillatory motions for the remaining of the propagation time. However, the damping behavior in this case is faster than the one in the populations computed using the 5D Hamiltonian model (comparing the

6

Journal of Electron Spectroscopy and Related Phenomena 237 (2019) 146899

B. Nikoobakht

It can be seen that both computed spectra in Fig. 6a and b are in good agreement with the corresponding experimental spectrum in Fig. 6c. It seems that the energetic position of each peak in Fig. 6b in comparison with Fig. 6a is closer to the corresponding peak in the experimental spectrum in Fig. 6c. Also, in Fig. 6b, the onset of the experimental photoelectron spectrum is nicely reproduced by employing the 6D model, while this is not well reproduced using the 5D model calculation. Furthermore, the splitting in the second peak (at ≈9.12 eV) appearing pronouncedly in Fig. 6a is reduced considerably in Fig. 6b and it is almost similar to the corresponding peak in the experimental spectrum of Fig. 6c. The most striking difference between Fig. 6a and b is related to the high-resolution spectra. It can be seen that the highresolution spectrum shows a lot of vibrational structure in the inset of Fig. 6b, while considerable reduction in the vibrational structures in Fig. 6a can be observed. This demonstrates that with a larger number of modes a better microscopic description of the complex processes underlying the experimental line broadening is obtained. Thus, the role mode S8 in regeneration of the experimental photoelectron spectrum is important. 4. Summary and conclusions In this work, we provided comprehensive theoretical investigations for better understanding of the ground state dynamics of TBD with focusing on the first band of the photoelectron electron spectrum using the high levels of the electronic structure methods and the nuclear quantum dynamical approach. The potential energy curves computed using the MRCI method with CAS(8,8) and CAS(10,7) together with the MCTDH algorithm yield a fully quantal description of the nuclear motion of X2Bg ⟵ X1Ag TBD. PESs calculated using different levels of theories are quite similar along the modes S3, S8, S13, S17, S19 and S20. Based on ab initio electronic results, we are certain that the shape and behavior of PESs and coupling constants calculated at the MRCI/CAS (10,7) and CAS(8,8) are accurate. This finding has been confirmed by the nuclear dynamics simulations using PESs computed at the MRCI/ CAS(10,7) level of theory, where the energetic position of the peaks and vibronic structures of the computed spectra were close to the experimental spectrum. Especially, the photoelectron spectrum generated using the 6D model including the coupling mode S8 is in good agreement with the corresponding experimental spectrum of TBD cation. This finding leads to this conclusion that the time-dependent population analysis considering the seven vibrational modes considering the S8 provides an accurate description of the ground state dynamics of TBD. In fact, after taking into account the mode S8 in the dynamics, it is found that the population transfer occurs within 50-60 fs, while this population transfer using the 5D Hamiltonian model lacking of S8 takes place within 30-40 fs (see Fig. 5). We demonstrated that the contribution of vibrational mode S8 is vital to reach a reliable description of the ground state dynamics of TBD.

Fig. 6. The simulated photoelectron spectrum of TBD corresponding to the X2Bg cationic state obtained using the S3, S13, S17, S19 and S20 modes in (a) while in (b) the coupling mode S8 included. The spectra in a and b were calculated at the MRCI/CAS(10,7) level of theory, respectively. The dephasing time is 100 and 1000 fs for the low and high resolution spectra, where the latter ones are shown in insets. In (c), the first band of the experimental spectrum. The experimental spectrum has been reprinted with permission from [11].

Acknowledgements We are grateful to the fruitful discussion with Profs. H. Köppel and Andreas Dreuw and their continuous supports.

high-resolution spectrum is employed (all high-resolution spectra are shown in insets of Fig. 6). Appendix A. G-matrix elements

In this appendix, we demonstrate the detail derivation of the G-matrix elements which are used for the evaluation of the kinetic operator in the Hamiltonian model of Eq.(1). The parameters p are designated for the inverse bong lengths (the atom numbering in Fig. 1):

p1 =

1 , rC1 C2

p2 =

1 , rC2 C3

p4 =

1 , rC1 H5

p5 =

1 , rC1 H6

p6 =

1 . rC1 H7

(A.1)

The required sine and cosine terms of the bond angles are defined as follows:

7

Journal of Electron Spectroscopy and Related Phenomena 237 (2019) 146899

B. Nikoobakht

sa1 = sin

C1 C2 C3,

sa2 = sin

ca1 = cos

C1 C2 C3,

ca2 = cos

H5 C1 C2 ,

sa3 = sin

H5 C1 C2 ,

H7 C2 C1,

ca3 = cos

H7 C2 C1.

(A.2)

The G-matrix elements are computed by symmetry adaptation from Tables III–IV of Ref. [31], where they are written for localized internal coordinates which are shown in Table 1. Thus, all possible G-matrix elements for all relevant vibrational modes computed using the method described in Refs. [13,12] are written as follows (mC and mH refer to the mass of carbon and hydrogen, respectively and cta is the cosine of the torsion angle for the CeCeCeC skeleton in Eqs. (A.3) and (A.5)):

G S3 S3 = 0.5mC , G S3 S13 = 0.5. p2 . mC . sa1 0.5. mC . p12

G S13 S13 =

0.5. p2 . mC . cta. sa1, 0.5. mC . p22 + 0.5. p2 . mC . p1 . ca1

+ 0.5. mC . cta. p22 G S17 S17 =

0.5.

mC . p12

0.5. mC . p2 . cta. p1 . ca1,

0.5. mC . p22 + 0.5. p2 . mC . p1 . ca1 p22

+ 0.5. mC . p2 . cta. p1 . ca1, 0.5. mC . cta. G S8 S8 = 0.5mC , G S17 S8 = 0.5p2 mC sa1 + 0.5p2 mC sa1cta, G S19 S19 = 0.046875mH p4 sa 2 2 + 0.03125mC p2 sa1 2 0.0625ca 22mC p12 sa 2 2

(A.3)

0.046875mC p42 sa 2 2 + 0.03125ca2mC p1 p4 sa 2 2

0.0625ca12 mC p22 sa 1 2

0.078125ca1mC p1 p4 sa1 1sa2 + 0.03125ca2mC p1 p5 sa 2

2

1

1

0.03125ca1mC p1 p2 sa1 2

+ 0.25ca1ca2mC p12 sa1 1sa 2 1

0.03125mC p2 p6 sa1 1sa 3 1 + 0.03125mC ca3p1 p2 sa1 1sa 3 1 1

0.0625ca1mC p6 p1 sa1 sa 3 + 0.109375mC ca1ca3p12 sa1 1sa 3 1 + 0.0625ca2ca3mC p12 sa 2 1sa 3 1

0.046875ca1mC p1 p5 sa 2 1sa1 1

0.03125ca12mC p12 sa1 2 + 0.015625ca3mC p4 p1 sa 2 1sa 3 1,

G S20 S20 = 0.046875mH p52 sa 2 2

(A.4)

0.046875mH p62 sa 3 2

0.046875mC p52 sa 2 2 2 0.046875p6 mC sa 3 2 + 0.09375ca3mC p6 p1 sa 3 2 0.046875ca3mC p1 p5 sa 2 1sa 3 1 0.46875mH p42 sa 2 2 0.21875mC p22 sa1 2 0.046875mC p42 sa 2 2 + 0.03125ca2mC p1 p4 sa 2 2 0.0625ca22mC p1 2sa 2 2 0.09375ca 23mC p12 sa 3 2

0.0625ca12mC p22 sa1 2 + 0.21875ca1mC p1 p2 sa1 2 + 0.25ca1ca2mC p12 sa1 1sa 2 1

0.09375ca1mC p1 p4 sa1 1sa 2 1

0.1171875ca2mC p2 p1 + 0.0625mC p4 p5 sa 2 2

2

+ 0.03125ca2 mC p1 p5 sa 2 + 0.140625mC p2 p6 sa1 1sa 3 1 1

0.15625mC ca3p1 p2 sa1 1sa 3 1

0.0625ca1mC p6 p1 sa1 sa 3 + 0.109375mC ca1ca3p12 sa1 1sa 3 1 0.03125mC ca2p1 p6 sa 2 1sa 3 1 + 0.0625ca2ca3mC p12 sa 2 1sa 3 1 0.03125ca1mC p1 p5 sa 2 1sa1 1 + 0.015625ca3mC p4 p1 sa 2 1sa 3 1 0.046875mH p62 sa 3 2 0.046875p62 mC sa 3 2

1

0.046875mC p52 sa 2 2

0.03125ca12 mC p12 sa 1 2

0.09375ca 32mC p12 sa 3 2

+ 0.09375ca3mC p6 p1 sa 3 2

0.046875ca3mC p1 p5 sa 2 1sa 3 1.

(A.5)

Inserting Eqs.(A.3)–(A.5) into Eq. (2) results in the kinetic operator.

[17] H.-D. Meyer, G.A. Worth, F. Gatti (Eds.), Multidimensional Quantum Dynamics: MCTDH Theory and Applications, Wiley-VCH, Weinheim, 2009. [18] R. Kosloff, H. Tal-Ezer, Chem. Phys. Lett. 127 (1986) 223. [19] A. Komainda, B. Ostojić, H. Köppel, J. Phys. Chem. A 117 (2013) 8782. [20] B. Nikoobakht, H. Köppel, Chem. Phys. Lett. 651 (2016) 221. [21] P. Pulay, G. Fogarasi, G. Pongar, J. Boggs, A. Vargha, J. Am. Chem. Soc. 105 (1983) 7037. [22] T.H.J. Dunning, Chem. Phys. 90 (1989) 1007. [23] D.E. Woon, T.H.J. Dunning, Chem. Phys. 103 (1995) 4572. [24] A. Komainda, D. Lefrancois, A. Dreuw, H. Köppel, Chem. Phys. 485 (2017) 27. [25] H.-J. Werner, et al., Molpro, Version 2012.1, A Package of Ab Initio Programs, (2012). [26] J.F. Stanton et al., CFOUR, Coupled-Cluster Techniques for Computational Chemistry, A Quantum-Chemical Program Package. For the current version, see http://www.cfour.de. [27] L. Visscher, et al., Dirac, Version 2014, A Package of Ab Initio Programs, (2014). [28] J.F. Stanton, R.J. Bartlett, J. Chem. Phys. 98 (1993) 7029. [29] M. Pernpointner, Chem. Phys. 338 (2007) 44. [30] H. Köppel, W. Domcke, L.S. Cederbaum, Adv. Chem. Phys. 57 (1984) 59. [31] J.H. Frederick, C.J. Woywood, Chem. Phys. 111 (1999) 7255.

References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13]

T. Polivka, V. Sundstrom, Chem. Rev. 104 (2004) 2021. C. Dugave, L. Demange, Chem. Rev. 103 (2003) 2475. M.A. van der Horst, K.J. Hellingwerf, Acc. Chem. Res. 37 (2004) 13. D. Schuetzle, et al., Environ. Health Perspect. 102 (1994) 3. N.L. Morrow, Chem.-Biol. Interact. 135 (2001) 137. B. Hudson, B. Kohler, K. Schulten, Excited States 6 (1982) 1. G. Orlandi, F. Zerbetto, M.Z. Zgierski, Chem. Rev. 91 (1991) 867. W. Fu, Y. Haas, S. Zilberg, Chem. Phys. 259 (2000) 273. V. Vaida, Acc. Chem. Res. 19 (1986) 114. B. Nikoobakht, A. Dreuw, H. Köppel, Chem. Phys. 515 (2018) 654. D.M.P. Holland, et al., J. Phys. B: At. Mol. Opt. Phys. 29 (1996) 3091. E.B.J. Wilson, Chem. Phys. 7 (1939) 1047. E.B.J. Wilson, J. Decius, P.C. Cross, Molecular Vibrations – The Theory of Infrared and Raman Vibrational Spectra, Dover Publications, New York, 1980. [14] H.D. Meyer, U. Manthe, L.S. Cederbaum, Chem. Phys. Lett. 165 (1990) 73. [15] M.H. Beck, et al., The Heidelberg MCTDH Package, Version 8.4.10, (2015) See http://www.pci.uni-heidelberg.de/tc/usr/mctdh/doc/index.html. [16] M.H. Beck, A. Jäckle, G.A. Worth, H.-D. Meyer, Phys. Rep. 324 (2000) 1.

8