Chemical Modification of Conical Intersections in Photoisomerization

Dec 13, 2016 - Guided by a notion of symmetry-breaking modulation or control of the so-called symmetry-allowed conical intersection by shining laser p...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/JPCA

Chemical Modification of Conical Intersections in Photoisomerization Dynamics of Butadiene Derivatives Published as part of The Journal of Physical Chemistry virtual special issue “Mark S. Gordon Festschrift”. Hiroki Ichikawa Department of Basic Science, The University of Tokyo, Komaba 3-8-1, Tokyo 153-8902, Japan

Kazuo Takatsuka* Fukui Institute for Fundamental Chemistry, Kyoto University, Takano-Nishikiraki-cho 34-4, Sakyo-ku, Kyoto 606-8103, Japan ABSTRACT: Guided by a notion of symmetry-breaking modulation or control of the so-called symmetry-allowed conical intersection by shining laser pulses [Arasaki, Y.; et al. Phys. Chem. Chem. Phys. 2010, 12, 1239], we here explore a possibility of the modulation of the symmetry-allowed conical intersection by chemical substitution with functional groups. As a first case study, we choose photoisomerization dynamics of s-trans-1,3-butadiene H2CCHCHCH2 with one of the terminal hydrogen atoms being replaced by −CF3. The target here is not the control of the rate of nonadiabatic transition but to know which one of the double bonds is more frequently isomerized in the radiationless quenching process on the way back to the ground state. We analyze when and how the symmetry is broken by tracking ab initio molecular dynamics paths, the mean-field paths with use of the nonadiabatic electron wavepacket dynamics, and the associated branching paths.

1. INTRODUCTION Conical intersection (CI) is one of the key notions to determine the rate process in photochemical dynamics.1−7 Among others the so-called symmetry-allowed conical intersection is practically found most common.2 It emerges when two or more eigenfunctions belonging to mutually different irreducible representations of a symmetric group of an electronic Hamiltonian happen to have a common electronic energy at a given nuclear configuration (molecular shape). Let us take an example from the electronic states of NO2 molecule. This molecule, having C2v symmetry, happens to have a one-dimensional manifold of CI, which is made up by the 2A1 ground state and 2B2 state. Because their matrix element always satisfies ⟨Φ(2B2)|Ĥ (el)|Φ(2A1)⟩ = 0, where Ĥ (el) and Φ(1B2) are the electronic Hamiltonian and the electronic wave function of symmetry 2B2, their adiabatic potential energy hypersurfaces can exactly cross each other, and thereby the noncrossing rule leading to the so-called avoided crossing (AX) does not apply. A nuclear wavepacket sitting on the ground-state 2A1 is promoted to the 2B2 state by a short-pulse pump laser and spontaneously comes to the CI region.8−13 The nonadiabatic interaction due to nuclear kinematic interaction then causes nonadiabatic state mixing between the two states, resulting in a wavepacket bifurcation after the passage. A natural idea then is that because the symmetry-allowed conical intersection emerges from symmetry requirement, the energy crossing can be modified by introducing a symmetry-breaking perturbation into the system.9,11,12 In the above NO2 system one can actually © XXXX American Chemical Society

displace and/or resize the CI manifold by shining an external weak laser on it. With an appropriate tuning of the laser parameters, one can achieve a good control of the wavepacket bifurcation. Note that this is a modulation of the nonadiabatic interaction itself with only a little deformation of the original potential surfaces. We have shown that such an effect on the relevant rate process of nonadiabatic transition can be indeed observed in terms of pump−probe photoelectron spectroscopy.12 Nature cannot make use of such a laser control and needs some other efficient alternatives in photochemical processes of organic or bioorganic chemical compounds such as photoisomerization dynamics of retinal. Therefore, we herein explore a possible effect of symmetry breaking in terms of chemical modification such as asymmetric substitution of atoms and/or functional groups on the symmetry-allowed conical intersections. As such a rather simple prototype, we here perform a fulldimensional nonadiabatic electron wavepacket dynamics study on photoisomerization of s-trans-1,3-butadiene H2CCH CHCH2 with one of the terminal hydrogen atom being replaced by −CF3 group as H2CCHCHCHCF3 in the trans position with respect to the group of H2CCH−. It is well-known that pure 1,3-butadiene is photoexcited to its 1Bu state (S2 state) and followed by nonradiative transitions Received: October 23, 2016 Revised: December 12, 2016 Published: December 13, 2016 A

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

C3C4. We hereafter track which one of C1C2 and C3C4 is isomerized.

through conical intersections to the ground state, leaving the isomerization in one or two of the two symmetrical double bonds behind. The probability of which one of the double bonds is isomerized is of course even to them due to the symmetry. The entire mechanism of this series of nonradiative transitions has been extensively studied in the literature. The ́ latest one by Levine and Martinez in terms of the so-called multiple spawning method14 with carefully chosen electronicstate spaces is a landmark in the study of photochemical dynamics,15 their clarified mechanism of which will be summarized later in the context of our study. (This paper15 contains an extensive list of references on photoisomerization of 1,3-butadiene.) The photochemical dynamics of the present system is far more complicated and chemically interesting than that in laser control of NO2. In particular, our interest here is not in the rate of quenching but in identifying which one of the double bonds is to be isomerized with higher probability and the associated mechanism. We hope that the present study may eventually set a theoretical foundation about the selectivity of specific double bonds in photoisomerization of retinal,16 although it is recognized that the photodynamics of retinal is not a simple extrapolation of butadiene.15,17 Chemical modulation over conical intersections is not as simple as pulse-laser control, because substitution or addition of functional groups to break molecular symmetry should more or less modulate the potential energy hypersurfaces, too. In the present system, the effect of the substitution of −CF3 is found to be small in the π−π* (11Ag to 1Bu) photoabsorption process. Therefore, one can pay main attention to the chemical modification of the relevant conical intersections. Because the conical interactions play the critical role in choosing the product channels and yet exist only very locally, modification of CI is expected to be far easier than reshaping the global potential energy surfaces. This paper is organized as follows. We first briefly review the dynamical mechanism of photoinduced isomerization of s-trans1,3-butadiene in section 2, which gives the background about the conical intersections of this molecule. After presenting the theoretical and computational methods used in this paper in section 3, we track in section 4 the dynamics of the CF3substituted butadiene in terms of ab initio molecular dynamics, full-dimensional nonadiabatic electron wavepacket theory with use of the semiclassical Ehrenfest method, and path-branching representation. Shift of the location and conversion of the nature of the conical intersection between the S2 and S1 states are also discussed. The paper concludes in section 5 with remarks with respect to more general aspects about chemical modification of the conical intersections.

Figure 1. Schematic representaion in terms of the Pauling resonance picture of the main electronic configurations dominating the optically dark S1 (21Ag) state and bright S2 (11Bu).

Suppose the molecule is placed on a plane having C2h symmetry. As schematically represented in this figure, the main electronic configurations of the valence structure of 1Bu and 21Ag are biradical, tetraradical, zwitterion (charge-transfer state or polarized state in addition to single and double bonds) in rather classical terminology: The major characteristic of the Bu state is its biradical feature, [C1C2Ċ 3Ċ 4] − [Ċ 1 Ċ 2C3C4]. On the contrary, the major configurations of the 21Ag states are not only of [C1C2Ċ 3Ċ 4] + [Ċ 1Ċ 2 C3C4] but also of ionic character [C1+C2−C3C4] + [C1−C2+C3C4] + [C1C2C3+C4−] + [C1C2 C3−C4+] along with a tetraradical character [Ċ 1Ċ 2Ċ 3 Ċ 4]. Although the wave function of Bu state is antisymmetric with respect to the 180° rotation around the C2 axis perpendicular to the molecular plane, that of Ag is symmetric. These two states can form a manifold of the symmetryallowed conical intersection for the molecular geometry at which the two states happen to have the same energy in this molecular symmetry, because they do not mix each other through the electronic Hamiltonian Ĥ (el) because (el) ⟨Φ(1Bu)|Ĥ |Φ(21A g)⟩ = 0

(1)

On the contrary, the nonadiabatic interaction can cause a mixing among the states. (See ref 18 for an example of spontaneous symmetry breaking induced by the nonadiabatic interactions along with its possible spectral observation.) Symmetry breaking due to a substituent effect would also mix them. It is useful then to recall a very basic fact that the addition of [C1C2Ċ 3Ċ 4] − [Ċ 1Ċ 2C3C4] and [C1C2 Ċ 3Ċ 4] + [Ċ 1Ċ 2C3C4] (notice the difference in the sign) should give a larger component that has biradical character in Ċ 3Ċ 4, whereas the corresponding subtraction should result in Ċ 1Ċ 2 with the large amplitude. In the symmetric situation of the pure butadiene, anything that happens on the C1C2 bond can take place on the C3C4 bond with the equal probability. However, chemical substitution of one of the hydrogen atoms can break this equal opportunity. ́ 2.2. Levine−Martinez (LM) Mechanism of Photoisomerization. Through a very extensive and careful study on the photodynamical process of s-trans-1,3-butadiene starting ́ 15 from 1Bu back to the ground state, Levine and Martinez observed and came to the following conclusions about the mechanism of the radiationless transition (referred to as LM mechanism in this paper): (1) In the initial 5 fs after excitation, the molecule undergoes nonadiabatic transitions that result in

2. BUTADIENE AND ITS SYMMETRY-ALLOWED CONICAL INTERSECTION; BACKGROUND 2.1. Symmetry in the Electronic States in the Valence Bond Representation. The ground state of s-trans-1,3butadiene is photoexcited to the optically bright 1Bu state, which is S2 state having the lower and dark counterpart state S1 of the symmetry 1Ag. The electronic structures of these states are rather complicated in the configuration interaction picture based on the molecular orbitals. We therefore depict them in the valence-bond representation along with their Pauling resonance picture in Figure 1. The carbon atoms consisting of the skeletal structure of butadiene is numbered as C1C2 B

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 2. Schematic representation of the typical histrory of the excited-state paths for s-trans-1,3-butadiene running across the nonadiabatic regions after photoexcitation to the isomerization in, for instance, C3−C4 in the ground state. Twisting motion around a broken double bond and pyramidization like in sp3 conjugation are observed as characteristic deformation in molecular geometry.

55% of the excited population remaining on the bright 11Bu state and the remainder transferring to the dark 21Ag state. In contrast to longer polyenes, which transfer to the dark state without breaking planarity, butadiene twists about a single terminal carbon−carbon bond regardless of whether it is on the bright (1Bu) or dark (21Ag) state. (2) After twisting, the remaining S2 population is quenched to S1 through a seam of conical intersections. (3) The population on S1 is split between two oppositely polarized charge-transfer states. (4) Population in both wells is quenched to the ground state quickly via conical intersections. (5) Half of the population has quenched to the ground state (S0) within 200 fs after the photoexcitation. The successive events in the geometrical deformation and nonadiabatic transitions for the butadiene to undergo right after photoabsorption to the S2 (1Bu) state is schematically displayed in Figure 2, which is drawn on the basis of the mechanism ́ 15 and our own computaestablished by Levine and Martinez tional observation. The substituted butadiene should experience the similar history of electronic and nuclear dynamics, except that the critical symmetry is broken. Note that the lowest energy of S2, if it does not undergo nonadiabatic transition to S1, is mainly of the configuration [C1C2Ċ 3Ċ 4] − [Ċ 1Ċ 2C3C4] with the twisting angle 90° being high enough with a large energy gap from the energy of S0 at the same molecular geometry. Therefore, nonadiabatic transition between the biradical state and the ground state is too improbable to be ascribed to the main path for the direct nonadiabatic quenching from S2 to S0. On the contrary, the major components of the S1 state, which is correlated to 21Ag at the initial geometry, are of charge-transfer (zwitterionic) nature, which guides the excited state to the ground state through conical intersection. As observed in the LM mechanism, the charge-transfer state like, for instance, C3+C4− in the S1 state proceeds with the pyramidization motion at C4. On the contrary, the 90° twisted structure around the biradical pair of carbons like Ċ 3Ċ 4 is obviously not favorable for charge transfer, because the relevant 2p molecular orbitals come to a mutually perpendicular orientation having the zero overlap. Therefore, the pyramidization should take place simultaneously along with the twisting motion. Our less extensive calculations on s-trans-1,3-butadiene with a small number of configurations based on a smaller atomic basis functions qualitatively support the LM mechanism. We therefore continue our study on the substituent effect of one of the terminal hydrogen atoms attached on C4 with CF3 group.

3. THEORETICAL METHODS To make this paper self-contained as much as possible in a concise manner, we very briefly outline the theoretical bases of mixed quantum-classical representation of nonadiabatic electron wavepacket dynamics along with the relevant computational details. 3.1. Theory of Path-Branching Representation and the Semiclassical Ehrenfest Method. The theoretical method for nonadiabatic electron wavepacket dynamics in path-branching representation is first outlined (see refs 6, 19, 20, and 22−27 for details.) The standard nonrelativistic full dimensional electronic and nuclear Hamiltonian is usually written as 1 2

H(r ,Q ) =

∑ Pk̂

2

+ H (el)(r ;Q ) (2)

k

(el)

where H (r;Q) is the electronic Hamiltonian for the electronic r and nuclear Q coordinates. We represent the electronic wavepackets that are to be propagated in time along nuclear paths Q(t) as Ψelec(r ,Q (t )) =

∑ CI(t ) ΦI (r ;Q (t ))

(3)

I

with t being time. {ΦI(r;Q(t))} serves as a basis set, which is composed of the Slater determinants or the so-called configuration-state functions (CSF). The time-dependent coefficients CI(t) are to be determined with respect to the basis functions. The motion of Ψelec is driven by the coupled equations iℏ

∂CI = ∂t



∑ ⎜⎜HIJ(el) − iℏ ∑ Q̇ kXIJk − J



k

ℏ2 4



∑ (Y IJk + Y JIk*)⎟⎟CJ k



(4)

where HIJ(el) = ⟨ΦI |Ĥ Y IJk =

ΦI

(el)

|ΦJ ⟩,

∂2 ΦJ ∂Q k 2

XIJk =

ΦI

∂ ΦJ , ∂Q k

and

(5)

with Qk being the kth component of the vector Q. The secondorder terms YkIJ are neglected in our usual practice, which should be justified by the presence of the factor ℏ2 multiplied to C

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A them. The nuclear path Q(t) is driven by the so-called force matrix19 ⎡ ∂H (el) IJ + FIJk = −⎢ ⎢⎣ ∂R k

mixing coefficients of the electronic configurations thus far given by the SET calculations made up to the exit point. We practically ignore minor sub-branches carrying very small electronic population. 3.2. Computational Details. First the molecular structures of both the nonsubstituent butadiene and the CF3 derivative are optimized in their ground (S0) states. The vibrational analysis is conducted at the minimum energy geometries, assigning the zero-point energy, and giving a phase space point on each (separated) normal coordinate, (Qi, Pi) with Qi and Pi being the coordinate and its associated momentum for the ith normal mode. By writing (Qi, Pi) = (Ri cos θi, Ri sin θi), we randomly choose the phase θi in each mode, the resultant (Qi, Pi) of which are now used as an initial condition of the nuclear dynamics to initiate the dynamics on the S2 (1Bu) state just after photoabsorption for each molecule. Recall that the lowest excited state S1 is optically dark. With these initial conditions, classical paths begin to run on the S2 surfaces according to the ab initio MD method (AIMD) using the energy gradients generated in terms of configuration interaction method, and as soon as the energy difference between two adiabatic potential energy surfaces (S2 and S1) is smaller than 0.02 hartree (0.544 eV), the method of evolving the paths is switched to the SET method. These paths starting from the S2 state are supposed to have the electronic-state coefficients to be used in eq 4 as unity for S2 and zero for all the other adiabatic states. Then after the passage of the strong nonadiabatic regions, where the energy difference exceeds 0.02 hartree, the path is branched into two daughter paths according to the above-described simplified path-branching scheme, each of which resumes to run in the AIMD scheme until it comes to another quasi-degenerating regions between S2 and S1 states or S1 and S0. Therefore, some of them can remain on an adiabatic surface for somewhat a long time even after many passages of the nonadiabatic regions, and the other can pass across quickly. Thus, paths proliferate every time the path-branching procedure is applied. Because in the SET scheme the dynamics of electronic-state mixing by eq 4 is carried out within the quantum chemical scheme (actually we have implemented our own codes into the GAMESS package29,30), the coefficients associated with the individual electronic states, not only the probability but also the phase information, are accurately generated. These coefficients are transformed from the appropriate basis functions like CSF to those for the adiabatic electronic states. Thus, each path running on one of the adiabatic potential energy surfaces in an adiabatic region always retains the information about electronic-state population of this adiabatic state, even after having branched many times. The trajectories are propagated from 30 different initial conditions with the time interval of 0.1 fs for up to 500 fs using velocity-Verlet method. All the nonadiabatic electron wavepacket calculations are performed in the level of CISD/RHF/631G with all occupied orbitals except the HOMO and HOMO−1 are constrained to be doubly occupied in the configurations used in the CI. The number of CSFs is equal to 2278 for the nonsubstituent and 5356 for the CF3 derivative, respectively. The number of thus branched (large component) paths has amounted to 335 for normal butadiene, the electronic-state norm of which covers 74% of the total electronic populations, whereas 180 for CF3-substituted butadiene covering 65% of the total electronic-state norm. Other minor branching paths (carrying a small amount of electronic population) are not to be included into the statistical



∑ (XIKk HKJ(el) − HIK(el)XKJk )⎥ ⎥⎦

K

⎡ ∂X l ∂XIJk ⎤ IJ ⎥ − + iℏ ∑ Q̇ l ⎢ ⎢⎣ ∂Q k ∂Q l ⎥⎦ l

(6)

which is equivalent to (el)

FIJk =

ΦI

∂Ĥ ∂R k

ΦJ

⎡ ∂X l ∂XIJk ⎤ IJ ⎥ − + iℏ ∑ Q̇ l ⎢ ⎢⎣ ∂Q k ∂Q l ⎥⎦ l (7)

when the basis set {ΦI(r;Q(t))} was complete. The second term in eq 6 represents the components analogous to the Lorentz force in electromagnetic theory.19,20 Because FkIJ is a matrix, it provides as many forces to nuclei as the number of involved electronic states at every small time step to induce infinitely many branching paths.20 Technically though, it is virtually impossible to obtain the exact solutions, and instead, we approximately generate a small number of representative paths to substitute a cascade consisting of the infinitely many branching paths.20 A more drastic approximation is to take an average of the force matrix over an electron wavepacket under study in such a way that k FSET =−

∑ CI*(t ) FIJkCJ(t ) I ,J

=−

Ψelec(Q (t ))

∂H (el) Ψelec(Q (t )) ∂R k

r

(8)

which gives an averaged path but not branching paths. This mean-field method is called the semiclassical Ehrenfest theory (SET).28 (Note that the conventional SET does not include the ℏ2

terms of − 4 ∑k(YkIJ + YkJI*) in eq 4.) The equations of motion eqs 4 and 6 or 8 are to be numerically integrated in the on-thefly scheme. It is widely known, however, that semiclassical Ehrenfest theory gives unphysical paths after the passage of nonadiabatic regions, because they keep running on the averaged potential energy surfaces even after nonadiabatic interactions are already switched off. Theory of path-branching representation gives a way for the averaged path to branch into the individual asymptotic regions, thereby representing realistic dynamics of quantum wavepacket bifurcation due to nonadiabatic interactions. The method of path-branching should be made theoretically on the basis of diagonalization of the force matrix, eq 6. However, this procedure is rather time-consuming for the present system in view of the molecular size. Therefore, we apply a simplified version of the path branching scheme,20 in which each SET path once generated in a nonadiabatic region is forced to branch at an exit of it, where the relevant coupling elements become smaller than a predetermined threshold value. Suppose that the end point of this SET path is located at a phase-space point (Q, P). Then a single path is to be generated on each adiabatic potential energy surface (PES) with its own energy gradient as in the ab initio molecular dynamics. The electronic population on each path is assigned according to the D

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 3. Dynamics of the adiabatic electronic-state mixing along the set of branching paths after excitation to the final stage. The population of the three relevant states dominating the paths are plotted against time for (a) s-trans-1,3-butadiene (b) CF3-substituted butadiene.

conditions. However, we did not see a major difference between ours and those of LM in a qualitative and phenomenological level and therefore conceive that this difference does not matter in the arguments on symmetry breaking. 4.2. Tracking the Isomerization Dynamics. 4.2.1. Geometrical Parameters. To illustrate the dynamics of geometrical isomerization, we track the time evolution of selected geometrical parameters. Among the parameters used in the full dimensional calculations, we choose, in harmony with the choice in LM,15 (i) the dihedral angles for the bond of C1C2 (called bond 1) and for C3C4 (bond 2), and (ii) the pyramidal angles associated with bond 1 and bond 2. (See Figure 4 for definitions.) The dihedral angles indicate to what

ensemble, from which the average or statistical distributions are deduced.

4. DYNAMICS MODULATED BY SYMMETRY BREAKING WITH CHEMICAL SUBSTITUTION In this section we study the details about the effects of substitution in the various dynamical aspects by comparing strans-1,3-butadiene and H2CCH−CHCHCF3. The global potential energy surfaces, if available, are certainly helpful to capture the possible chemical events. It is, however, prohibitively time-consuming to generate them, and the information thereof is mostly limited to the static aspects of the electronic states. Instead, we let nonadiabatic electron wavepacket dynamics run along averaged paths or the branching paths with the initial conditions stated above. These methods reveal the electronic-state mixing along the paths. 4.1. Mixing among the Adiabatic Electronic States. Using the computational scheme presented above, we first track the dynamics of electronic-state mixing for about 500 fs right after photoexcitation, which is long enough for these molecules to come to the electronic ground state after geometrical isomerization. Figure 3 shows the time evolution of the populations over the adiabatic states of S2, S1, and S0, each of which has been summed up over the individual branching paths at each time (recall that the initial paths at t = 0 start with AIMD method). The overall successive behaviors of statemixing seem to be similar to each other of the two molecular systems, the normal butadiene in panel a of Figure 3 and the substituted butadiene in panel b. A minor difference can be seen in that the penetration from the S2 state to S1 is a little slower and fewer in the substituted butadiene. On the contrary, CF3 substitution seems to make the final transition from S1 to S0 faster than the nonsubstituted one, once the transition to S1 has been made. We note again that our target in this study is not the rate of the radiationless transition but the way of symmetry breaking in the isomerization. However, this figure alone cannot reveal how the isomerization is hindered by the substitution. Incidentally, the mixing of S1 state to S2 state begins around 50 fs in average after photoexcitation of the normal butadiene. Our estimated time scale is rather longer than the average time of the mixing between 1Bu and 21Ag observed by LM for the pure butadiene.15 This is presumably due to our poorer basis set and smaller electronic configurations limited to single and double excitations and also the difference in the initial

Figure 4. Definitions of the twisting and pyramidal angles.

extent the double bonds are broken by formation of biradicals (like Ċ 3Ċ 4) and/or charge transfers (like C3+C4−). The charge-transfer state forms a zwitterion in the bond(s); for instance, C3+−C4− is necessarily led to reorganization of the s− p hybridization from sp2 in the double bond to sp3 in the single bond in forming a pyramidal structure with C4− being its top. Thus, the relevant electronic change is also indirectly monitored by the dynamics of pyramidal angle. Note, however, the zwitterions should disappear as the S0 state begins to dominate the state, because the double bonds are to be rebuilt. Thus, the pyramidal angles return to zero in average as it was before photoexcitation, although vibrational excitation in the out-of-plane mode may remain. In this paper we take only the absolute values for both the dihedral and pyramidal angles measured, because we are only interested in which one of the double bonds is isomerized in the radiationless transition. (It is of course not possible to recognize the difference in the structures before and after isomerization in nonsubstituted butadiene. Nonetheless, the study of isomerization dynamics is crucial in a wider aspect.) The bond lengths of bond 1 and bond 2 are also crucial E

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

of bond 1, which should be resulted by their twisting vibrational motion. Even after returning to the electronically ground state, the molecule should be in a highly excited state in its vibrational mode, unless it is quenched. Each branching path like this one has a different history of electronic-state transition and geometrical deformation with a different transition probability. We therefore need to scrutinize the entire events in terms of the path ensemble. 4.2.3. Probability of Geometrical Change in the Ensemble of Branching Paths. In Figure 6 we show the time evolution of

geometrical parameters to see. However, we do not show them, because the above two quantities are enough to indicate the evolution of asymmetry. 4.2.2. Snapshots of Isomerization Dynamics along a Branching Path; An Example. To give an idea on how the dynamics proceeds, we take one randomly selected piece of branching path in a set of bifurcating sequence, which passes a nonadiabatic region nearby the conical intersection between S2 and S1 at 46.0 fs after photoexcitation, and another between S1 and S0 at 310.3 fs. Shown in Figure 5 are the geometries on these occasions.

Figure 5. Snapshots of a randomly chosen branching path, which is performing successive isomerization, first in the bond 1 at 46.0 fs after excitation in the transition from S1 to S0, and the second in bond 2 at 310.3 fs in passing from S1 to S0. Black, green, and gray balls represent carbon, florine, and hydrogen atoms, respectively.

To see how radical electrons are distributed over the molecule, which is one of the keys to characterize the excited state having biradicals or more, we display the unpaired electron distribution superposed on the molecular skeleton. Unpaired electron density D(r) is defined as31 D(r) = 2ρ(r,r) −

∫ dr′ ρ(r,r′) ρ(r′,r)

Figure 6. Time evolution of the geometry parameters (dihedral angle (θ), pyramidal angle (ϕ)) in the ensemble of the paths of the nonsubstituted butadiene: (a) dihedral angle for C1−C2 (bond 1); (b) diheral angle for C3−C4 (bond 2); (c) pyramidal angle for bond 1; (d) pyramidal angle for bond 2.

(9)

where ρ(r,r′) is the first-order spinless density matrix in the coordinate representation. The total electron density is just ρ(r) = ρ(r,r). D(r) is rewritten as D(r) =

∑ ni(2 − ni)|λi(r)|2 i

the dihedral angles (θ in the range of 0−180°) for bond 1 (C1− C2, panel a) and bond 2 (C3−C4, panel b), and the pyramidal angle (ϕ in 0−90°) with respect to C1 (panel c) and C4 (panel d for the nonsubstituted butadiene, just for reference. Presented there as a color chart are the total electronic-state populations summed over the S2, S1, and S0 states at each time and given angles. The values of the other geometrical parameters at each instance are ignored. The chart therefore simply represents a kind of probability for the molecule to have a given angle at a given time, the darker place indicating the higher probability. Note also that the present charts alone do not provide the full information about geometrical isomerization. First of all, it is noticed from Figure 6a,b that the twisting motion around the axis of C1−C2 or C3−C4 is very fast, and moreover, the place of large amplitude twisting motion can be frequently relocated between the bonds C1−C2 and C3−C4. This reflects the much faster motion of the electrons, by which the position of the biradicals is shifted (schematically) as

(10)

where λi(r) are natural orbitals with ni their associated occupation numbers. Although natural orbitals in timedependent wavepacket formalism can be complex-valued, their occupation numbers are certainly kept real. Obviously, singly occupied natural orbitals should make the highest contributions to D(r) in an invariant manner, explaining why D(r) can represent unpaired electrons. At time zero in panel a of Figure 5, CF3-substituted butadiene starts with almost a flat plane structure. The πunpaired electrons (radical electrons) are seen to be delocalized over C1, C2, C3, and C4 atoms, of the butadiene skeleton. This delocalization reflects the resonance structure as schematically represented in the upper panel of Figure 1 for the 11Bu state. (See ref 21 about a relationship between the Pauling resonance structure and electron wavepacket dynamics.) Panel b shows the moment (t = 46.0 fs) for bond 1 to be twisted and the electronic-state transition is occurring from S2 and S1. We here observe a distribution of the unpaired electrons localized mainly over the bond 1. About 300 fs later (panel c at 310.3 fs), another nonadiabatic transition from S1 and S0 takes place. The position of the CF3 group demonstrates that isomerization is now taking place in bond 2. The decrease of the unpaired electron population on C4 atom is presumably due to the partial formation of lone-pair on it by charge polarization. Still a large population of the unpaired electrons is seen over the area













C1C2C3C4 ↔ C1C2C3C4 ↔ C1C2C3C4 (11)

Figure 6 reveals that both the distributions of dihedral and pyramidal angles become larger after about 50 fs, indicating that significant geometrical deformation takes place since then. This is consistent with the history of the electronic-state mixing as shown in Figure 3a. The smaller fluctuation in the angles with nonzero values at the outset is due to the initial conditions for the trajectories, which have been chosen to have the zero-point F

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

CF3-substituted butadiene, should be determined on the passage of a path from S2 state to S1. We next confirm this fact. 4.3.1. Location of the Areas of Nonadiabatic Transition. Figure 8 shows the distribution of state passage from S2 to S1 (panel a) and from S1 to S0 (panel b) in the dihedral angle space for bond 1 and bond 2. The paths launched from the S2 state eventually pass across the region near the conical intersection, some of them proceeding to the S1 state. Such nonadiabatic passage of a path can be readily detected by the decrease of the electronic-state population of S2 during a short time-interval along each path. By accumulating the total decrease of the electronic population of S2 at each point of the dihedral angles of bond 1 and bond 2, without respect to the other geometrical parameters, we have attained a distribution function as in Figure 8a. To smooth the distribution functions, we have applied the method of kernel density estimation.32 The distribution function in panel b has been obtained in a similar way for the passage from S1 to S0. These figures suggest the location of the conical intersections projected onto this angle space. A quick inspection over panel b shows a dramatic asymmetry in the distribution for the transition from S1 to S0. More importantly, the asymmetric distribution in the transition from S2 to S1 suggests that bond 2 is twisted far more frequently than bond 1 is. 4.3.2. Statistical Ratio of the Passage through the Electronic States. To make a little more quantitative discussion, we tabulate the ratio of the probability of the passage from S2 to S1 versus that from S1 to S0 between bond 1 and bond 2. Table 1 presents not only the ratios for the

energy of the molecular vibration on the S0 states, as described above. The graphs in this figure for bond 1 and bond 2 look more or less symmetric, although a minor difference between them is noticed. We will see more about the statistical results in the next subsection. Next we show the similar quantities but for CF3-substituted butadiene in Figure 7. The most remarkable difference from the

Figure 7. Same as Figure 6 for the CF3-substituted butadiene.

nonsubstituted butadiene is found in the marked asymmetry between the dihedral angles for bond 1 and bond 2 (compare panels a and b), clearly demonstrating that the isomerization takes place mostly on bond 2, C3−C4, to which CF3 is bonded. Panel b of this figure illustrates that the large amplitude twisting motions eventually leading to internal rotation that begins at around 50 fs: The high probability of finding the twisting angle around 180° after about 50 fs lasting to 500 fs clearly illustrates that the isomerization is dominated by bond 2, whereas such a distribution around 180° for bond 1 is seen with only very low probability. As for the pyramidal angles, a significant difference is not observed except that the extent and probability of the pyramidal angle in bond 2 in CF3-substituted butadiene (panels c and d in Figure 7) are slightly smaller than those of the nonsubstituted butadiene (panels c and d in Figure 6). 4.3. Passage near the Conical Intersections from S2 to S1 and from S1 to S0. It is almost certain that the symmetry breaking in the isomerization dynamics observed as above in

Table 1. Probability (in %) of Nonadiabatic Passage from S2 to S1 and S1 to S0

S2 to S1 S1 to S0 a

bond 1 (nBDEa)

bond 2 (nBDEa)

bond 1 (−CF3b)

bond 2 (−CF3b)

51.4 49.4

48.6 50.6

28.0 10.2

72.0 89.7

Nonsubstituted butadiene. bCF3-substituted butadiene.

substituted butadiene but also those for the nonsubstituted one. We immediately notice that both bond 1 and bond 2 in the nonsubstituted butadiene have almost equal probabilities in the transitions both from S2 to S1 and from S1 to S0, as expected. On the contrary, it is observed that CF3-substituted butadiene by far favors isomerization of bond 2 in both of the transitions from S2 to S1 and from S1 to S0. We are thus able to conclude

Figure 8. Distribution of the probability of transition in CF3-substituted butadiene (a) from S2 to S1 and (b) from S1 to S0, represened in the dihedral angle space (θ-coordinates) for bond 1 and bond 2. See the text for the physical meaning of the distribution. G

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 9. For the nonsubstituted butadiene, (a) the adibatic potential energy surface for the S2 state and (b) difference between the adiabatic electronic energies of the S2 and S1 states, which are projected onto the dihedral (θ) and pyramidal (ϕ) angle space for bond 1. All the other geometrical parametrs are frozen.

Figure 10. For CF3-substituted butadiene, (a) potential energy surface (PES) of bond 1 in S2, (b) PES of bond 2 in S2, (c) energy difference between the S2 and S1 states for bond 1, and (d) energy difference between the S2 and S1 for bond 2, projected in the same manner as in Figure 9.

4.4.1. Deformation of the Conical Intersections. First Figure 9a shows the PES of the S2 state of the nonsubstituted butadiene projected onto the two-dimensional space with respect to bond 1. Due to the symmetry, this PES is exactly the same as that for bond 2. A rather wide range of the potential basin is observed in the area centered at θ = 90°. Notice that this potential basin is also flat and wide in the direction of the pyramidal angle. More interesting is the difference map of the PES between S2 and S1, shown in Figure 9b. The bottom, that is the smallest energy difference, is seen near θ = 105° and ϕ = 20°. Because all the geometrical parameters other than θ and ϕ are frozen, we do not see the exact subset of the conical intersection. Nevertheless, it is expected that the conical intersection can be readily accessed by relaxing the bond length of bond b1 (C3 C4) and the other parameters. Once again it is confirmed that

that bond 2 tends to be isomerized more frequently in CF3substituted butadiene and that this symmetry breaking is already induced by the passage from S2 to S1 state. 4.4. Characteristic Features of the Potential Energy Surfaces Nearby the Conical Intersections. Having presented the dynamical consequence of symmetry breaking, we next survey how the conical intersections between S2 and S1 have been modified by the chemical substitution. To do so, we consider the potential energy surface (PES) of S2 and the energy difference between S 2 and S 1 . Our graphical representation is made only in the two-dimensional space of the dihedral angle (θ) and pyramidal one (ϕ) of the relevant chemical bonds, bond 1 and bond 2, whereas all the other geometrical parameters are frozen to those values attained at the ground-state (S0) optimized geometry to save the computational time. H

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 11. (a) Dynamics of the adiabatic electronic-state mixing along the set of branching paths after excitation to the final stage for a hypothetical molecule, in which only the mass of the terminal hydrogen atom at C4 in the trans position is replaced with the mass of CF3 in the s-trans-1,3butadiene. The plots are made in the same manner as in Figure 3. (b) log plot of the figure in panel (a). The straight line represents a fitted one in the form of exp(−at + b).

and vice versa according to the scheme of eq 11 before the path successfully make a transition to S1 with a large probability. It is rather clear that bond 2 of the CF3-substituted butadiene has a favorable potential energy surface to make a nonadiabatic transition with large probability. The minimum energy in the (θ, ϕ)-space (Figure 10b) is located near (θ, ϕ) = (90°, 0°), which is readily accessed from the initial place at (θ, ϕ) = (0°, 0°) without significant change of the pyramidization motion. Actually, motion in this mode is more stiff and thereby pyramidization to the large angle tends to be suppressed. Besides, this lowest energy is a little lower than the minima of the other potential basins (Figure 9a and Figure 10a). Furthermore, Figure 10d shows that the transition region lies near the axis of ϕ = 0° with a wide range of θ coordinate around 90°. Thus, the paths for the biradical state created on bond 2 need not waste much time to reach the transition region by scanning a larger space. This explains why the isomerization in bond 2 for the CF3-substituted butadiene has much higher possibility than in bond 1. 4.4.3. On the Mass Effect. The above explanation based on the geometry of the potential energy surface alone does not seem to be consistent with a basic fact observed in Figure 3 that transition from S2 to S1 in the CF3-substituted butadiene is slower than the nonsubstituted butadiene. Although the goal of this paper is not to evaluate the absolute rates of the relevant reactions and transitions, we briefly touch on the initial decay rate of the S2 state from photon absorption to the first passage of the nonadiabatic region between S2 and S1. To account for this dynamics, one may need quantitative analysis on the amplitude of the nonadiabatic transitions besides the static electronic-state effects. Also, one of the most likely reasons should be the mass effect in the twisting motion; the CF3 group is far heavier than H. The mass and electrostatic effects of substituents on photoisomerization dynamics have been ́ extensively studied by Martinez and his colleague.14,17,33 We here survey the individual effects of mass and electronic states in this initial stage by approximately decoupling them by means of the method proposed in ref 17. To do so, we consider a hypothetical molecule, in which only the mass of the terminal hydrogen atom at the C4 in the trans position is replaced with the mass of CF3 in the s-trans-1,3-butadiene. No other parameters are changed from those of the pure 1,3-butadiene. Figure 11a shows the population dynamics of the adiabatic states S2, S1, and S0 for this hypothetical molecule as in the study of Figure 3. Several interesting features are noticed upon

the pyramidization should arise as a consequence of the charge polarization changing from Ċ 3Ċ 4 to C3+C4− through reorganization from sp2 (like in the methyl radical) conjugation to sp3 (as in the methyl anion having a lone pair). Therefore, the nonadiabatic transition from S2 to S1 must predominantly proceed when the paths come to this basin area. If so, it is anticipated that CF3 substitution with one of the hydrogen atoms attached on C4− should soften the pyramidization back toward the planar structure, because the strong electron affinity of CF3 group is expected to depopulate the lone pair on C4−. Then we turn to the similar quantities of the CF3-substituted butadiene in Figure 10; panels a and b represent the PES of S2 state for bond 1 and bond 2, respectively. Panels c and d show the difference map of the PES between S2 and S1 again for for bond 1 and bond 2, respectively. As for the PESs of the S2 state, they are both similar to that of the nonsubstituted butadiene. In particular, that of bond 1 (panel a) is particularly akin to (or almost the same as) that of nonsubstituted butadienes. On the contrary, we see in panel b that the potential basin centered at θ = 90° and ϕ = 0° is a little deeper to the extent of 0.5 eV. The PES in panel b is a little steeper toward the direction of the pyramidal angle. This should be a reflection of the effect of charge polarization induced by the CF3. In other words, the necessary charge polarization can be realized with less assistance of the sp3 hybridization. The dramatic breaking of the symmetry is seen in the difference map for bond 2 in panel d, whereas again the chart of the energy difference for bond 1 in panel c seems very similar to that of nonsubstituted butadiene. Although the figure does not specify the exact location of the conical intersections, it is certain that the dramatic modulation of the conical intersection has been induced almost exclusively in bond 2 by the substitution at the C4 atom. 4.4.2. Dynamical Consequences of the Modulation of the Potential Surfaces. Upon excitation to the S2 state in both types of butadiene, the biradical is created on bond 1 and bond 2 in a manner of the Pauling resonance and the large amplitude motion of twisting starts. The pyramidal motion can start simultaneously and the path begins to wander around in the potential basin in each S2 state (Figure 9a, Figure 10a,b) before finding a place of large probability of nonadiabatic transitions. Although the exact manifold on the conical intersection is among the best places to give the highest probability, the biradical can wander about frequently from bond 1 to bond 2 I

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

results in such a study. All the functional groups but CH3 are substituted to the C4 atom, in the trans position, as CF3 has been in the present paper. Only for the CH3 group have we examined the dependence of the position of substitution; on the C3 atom and C4 atom. The table shows that for the NH2substituted butadiene, for instance, bond 1 is predominantly modulated in the passage from its S2 to S1, whereas bond 2 is mainly affected in the transition from S1 to S0. Therefore, the double bonds are likely to undergo sequential isomerization in this order. [CHBH2]+ works in the opposite way to CF3 in that isomerization in bond 1 is mainly driven. On the contrary, normal butyl does not seem to specify the bonds. This table already suggests that there are many other interesting mechanisms behind the substituent effects in the radiationless deactivation for the conjugate molecules. Thus, the study of chemical modulation or the direct control of conical intersection in photoisomerization may need more extensive investigation from a wider scope.

inspection. However, we focus only on the initial stage dynamics and fit the adiabatic population of the S2 state in terms of an exponential function exp(−at + b). The parameters a and b thus obtained are listed in Table 2. Table 2. S2 State Population in the Initial Statge Dynamics As Fitted in exp(−at + b) molecules

a (fs−1)

b

pure 1,3-butadiene mass-only substituted butadiene CF3-substituted butadiene

2.71 × 10−2 2.06 × 10−2 8.21 × 10−3

0.117 0.303 0.376

A large difference in the parameter b between the nonsubstituted butadiene (0.117) and the group of substituted ones with the heavy “elements” (0.376 and 0.303) should reflect that the start of the large scale twisting and pyramidization motions is suppressed by the inertia due to the heavier mass. Likewise, the difference in the a value, a for the mass-only substituted butadiene subtracted from a for the pure butadiene subtracted, which represents essentially the mass effect, is estimated to be −0.65 × 10−2. Further, a for the CF3−butadiene subtracted from a for the mass-only substituted butadiene, which should almost purely reflect the electronicstate effect, is estimated to be −1.24 × 10−2. Therefore, the electronic-state effect is about 2 times larger than that of the mass effect. The role of the CF3 substitution in this stage is summarized as follows. First the heavier mass tends to retard the start of the large amplitude twisting and pyramidization motions. Once the large amplitude motion begins, both the electronic-state modification by CF3 and the mass effect make the population decay slower, but the former effect appears to be predominant in the present system.



Corresponding Author

*K. Takatsuka. E-mail: [email protected]. ORCID

Kazuo Takatsuka: 0000-0002-5624-1810 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was supported in part by a Grant-in-Aid for Basic Science from the Ministry of Education, Culture, Sports, Science, and Technology of Japan with the grant number 15H05752.



5. CONCLUDING REMARKS We have explored the possibility of modulating the manner of electronic-state mixing by introducing asymmetry into the symmetry-allowed conical intersection in terms of chemical substitution. As a first example, we chose photoisomerization dynamics of s-trans-1,3-butadiene H2CCH−CHCH2 with one of the terminal hydrogen atoms being replaced by −CF3. The target here has been to clarify whether we can control which one of the double bonds is more favorable for isomerization in radiationless transitions to the ground state. We have analyzed when and how the symmetry is broken by tracking the sequence of nonadiabatic electron wavepacket dynamics along the AIMD paths, the mean-field paths, and nonadiabatically branching paths. We have chosen −CF3 substitution for a rather obvious reason. However, there can exist other mechanisms to modulate the conical intersections. More systematic studies are worthwhile and necessary. Table 3 shows our preliminary

bond 1 bond 2 both

in transitionb from S1 to S0

CH3, [CHBH2]+, NH2, F, CH3a CF3 n-Bu

[CHBH2]+ CH3, CF3, NH2, CH3a n-Bu, F

REFERENCES

(1) Jasper, A.; Kendrick, B. K.; Mead, C. A.; Truhlar, D. G. In Modern Trends in Chemical Reaction Dynamics Part I Chapter 8; Yang, X., Liu, K., Eds.; World Scientific: Singapore, 2004. (2) Domcke, W., Yarkony, D., Köppel, H., Eds. Conical Intersections, Advanced Series in Physical Chemistry 15; World Scientific: Singapore, 2004. (3) Baer, M. Beyond Born−Oppenheimer: Electronic Nonadiabatic Coupling Terms and Conical Intersections; Wiley: Hoboken, NJ, 2006. (4) Levine, B. G.; Martínez, T. D. Isomerization Through Conical Intersections. Annu. Rev. Phys. Chem. 2007, 58, 613−6334. (5) Domcke, W., D. Yarkony, D., Köppel, H., Eds. Conical Intersections, Advanced Series in Physical Chemistry 17; World Scientific: Singapore, 2011. (6) Yonehara, T.; Hanasaki, K.; Takatsuka, K. Fundamental Approaches to Nonadiabaticity: Towards a Chemical Theory beyond the Born-Oppenheimer Paradigm. Chem. Rev. 2012, 112, 499−542. (7) Takatsuka, K.; Yonehara, T.; Hanasaki, K.; Arasaki, Y. Chemical Theory Beyond the Born-Oppenheimer Paradigm: Nonadiabatic Electronic and Nuclear Dynamics in Chemical Reactions; World Scientific: Singapore, 2015. (8) Arasaki, Y.; Takatsuka, K. Quantum Wavepacket Dynamics for Time-Resolved Photoelectron Spectroscopy of the NO2 Conical Intersection. Chem. Phys. 2007, 338, 175−185. (9) Arasaki, Y.; Takatsuka, K. Communication: Optical Conversion of Conical Intersection to Avoided Crossing. Phys. Chem. Chem. Phys. 2010, 12, 1239−1242. (10) Arasaki, Y.; Takatsuka, K.; Wang, K.; McKoy, V. Time-Resolved Photoelectron Spectroscopy of Wavepackets through a Conical Intersection in NO2. J. Chem. Phys. 2010, 132, 124307. (11) Arasaki, Y.; Wang, K.; McKoy, V.; Takatsuka, K. Monitoring the Effect of a Control Pulse on a Conical Intersection by Time-Resolved

Table 3. Effects of the Substituent Groups in the Isomerization Dynamics of Butadiene in transition from S2 to S1

AUTHOR INFORMATION

a

This methyl group is bound to the C3 atom. bThis transition is made successively after one from S2 to S1. J

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A Photoelectron Spectroscopy. Phys. Chem. Chem. Phys. 2011, 13, 8681− 8689. (12) Arasaki, Y.; Takatsuka, K. Reply to “Comment: On the Persistence of Conical Intersections under Perturbations. Phys. Chem. Chem. Phys. 2011, 13, 4756−4758. (13) Kraus, P. M.; Arasaki, Y.; Bertrand, J. B.; Patchkovskii, S.; Corkum, P. B.; Villeneuve, D. M.; Takatsuka, K.; Wörner, H. J. TimeResolved High-Harmonic Spectroscopy of Nonadiabatic Dynamics in NO2. Phys. Rev. A: At., Mol., Opt. Phys. 2012, 85, 043409. (14) Martínez, T. J. Insights for Light-Driven Molecular Devices from Ab Initio Multiple Spawning Excited-State Dynamics of Organic and Biological Chromophores. Acc. Chem. Res. 2006, 39, 119−126. (15) Levine, B. G.; Martínez, T. J. Ab Initio Multiple Spawning Dynamics of Excited Butadiene: Role of Charge Transfer. J. Phys. Chem. A 2009, 113, 12815−12824. (16) Bassolino, G.; Sovdat, T.; Liebel, M.; Schnedermann, C.; Odell, B.; Claridge, T. D. W.; Kukura, P.; Fletcher, S. P. Synthetic Control of Retinal Photochemistry and Photophysics in Solution. J. Am. Chem. Soc. 2014, 136, 2650−2658. (17) Kuhlman, T. S.; Glover, W. J.; Mori, T.; Møller, K. B.; Martínez, T. J. Between Ethylene and Polyenes - the Nonadiabatic Dynamics of cis-dienes. Faraday Discuss. 2012, 157, 193−212. (18) Li, Z.-W.; Yonehara, T.; Takatsuka, K. Nonadiabatic Electron Wavepacket sStudy on Symmetry Breaking Dynamics of the LowLying Excited States of Cyclic-B4. Chem. Phys. 2016, 464, 14−25. (19) Takatsuka, K. Generalization of Classical Mechanics for Nuclear Motions Nonadiabatically Coupled with Electron Wavepacket Dynamics and in Quantum-Classical Mixed Representation. J. Phys. Chem. A 2007, 111, 10196−10204. (20) Yonehara, T.; Takatsuka, K. Phase-Space Averaging and Natural Branching of Nuclear Paths for Nonadiabatic Electron Wavepacket Dynamics. J. Chem. Phys. 2008, 129, 134109. (21) Nagashima, K.; Takatsuka, K. Electron-Wavepacket Reaction Dynamics in Proton Transfer of Formamide. J. Phys. Chem. A 2009, 113, 15240−15249. (22) Yonehara, T.; Takahashi, S.; Takatsuka, K. Non-BornOppenheimer Dynamics for Electronic and Nuclear Wavepacket Dynamics. J. Chem. Phys. 2009, 130, 214113. (23) Takatsuka, K.; Yonehara, T. Nonadiabatic Chemical Dynamics in Intermediate and Intense Laser Fields. Adv. Chem. Phys. 2010, 144, 93−156. (24) Takatsuka, K.; Yonehara, T. Exploring Dynamical Electron Theory beyond the Born-Oppenheimer Framework: From Chemical Reactivity to Nonadiabatically Coupled Electronic and Nuclear Wavepackets On-the-fly under Laser Field. Phys. Chem. Chem. Phys. 2011, 13, 4987−5016. (25) Nagashima, K.; Takatsuka, K. Early-Stage Dynamics in Coupled Proton-Electron Transfer from π − π* State of Phenol to Solvent Ammonia Clusters: An Electron Dynamics Study. K. J. Phys. Chem. A 2012, 116, 11167−11179. (26) Yamamoto, K.; Takatsuka, K. Electronic Quantum Effects Mapped onto Non-Born-Oppenheimer Nuclear Paths: Nonclassical Surmounting over Potential Barriers and Trapping above the Transition States due to Nonadiabatic Path-Branching. J. Chem. Phys. 2014, 140, 124111. (27) Yamamoto, K.; Takatsuka, K. Dynamical Mechanism of Charge Separation by Photoexcited Generation of Proton-Electron Pairs in Organic Molecular Systems. A Nonadiabatic Electron Wavepacket Dynamics Study. Chem. Phys. 2016, 475, 39−53. (28) Schiff, L. I. Quantum mechanics; McGraw-Hill: New York, 1968. (29) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S. T.; Gordon, M. S.; Jensen, J. H.; Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S.; Windus, T. L.; Dupuis, M.; Montgomery, J. A. General Atomic and Molecular Electronic Structure System. J. Comput. Chem. 1993, 14, 1347−1363. (30) 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: Amsteldam, 2005.

(31) Takatsuka, K.; Fueno, T.; Yamaguchi, K. Distribution of Odd Electrons in Ground-State Molecules. Theor. Chim. Acta 1978, 48, 175−183. (32) Parzen, E. On Estimation of a Probability Density Function and Mode. Ann. Math. Stat. 1962, 33, 1065−1076. (33) Virshup, A. M.; Levine, B. G.; Martínez, T. J. Steric and Electrostatic Effects on Photoisomerization Dynamics Using QM/MM Ab Initio Multiple Spawning. Theor. Chem. Acc. 2014, 133, 1506.

K

DOI: 10.1021/acs.jpca.6b10680 J. Phys. Chem. A XXXX, XXX, XXX−XXX