Control of nuclear dynamics in the benzene cation by electronic wavepacket composition

The study of coupled electron-nuclear dynamics driven by coherent superpositions of electronic states is now possible in attosecond science experiments. The objective is to understand the electronic control of chemical reactivity. In this work we report coherent 8-state non-adiabatic electron-nuclear dynamics simulations of the benzene radical cation. The computations were inspired by the extreme ultraviolet (XUV) experimental results in which all 8 electronic states were prepared with significant population. Our objective was to study the nuclear dynamics using various bespoke coherent electronic state superpositions as initial conditions in the Quantum-Ehrenfest method. The original XUV measurements were supported by Multi-configuration time-dependent Hartree (MCTDH) simulations, which suggested a model of successive passage through conical intersections. The present computations support a complementary model where non-adiabatic events are seen far from a conical intersection and are controlled by electron dynamics involving non-adjacent adiabatic states. It proves to be possible to identify two superpositions that can be linked with two possible fragmentation paths.


Supplementary Note 1
The Qu-Eh method solves the time-dependent Schrödringer equation for both CI coefficients and nuclear wavefunction expansion coefficients. The equations for the time-evolution are variational 1 but the orbitals are only propagated to second order. The electronic wavefunction (see equation 1 in main text) is combined with the nuclear wavefunction which is expanded in gaussian wavepackets (see equation S1 below) where is the nuclear coordinate and , ! , ! and ! are respectively, the width, position of centre, momentum of centre and phase of gaussian wavepacket . In the derivation in ref. 1, the time-dependent parameters of the gaussian wavepackets are collected together in the vectors Λ ! .
By inserting the ansatz into the Dirac-Frenkel variational principle, we obtain equation of motions for the expansion coefficients , the gaussian wavepacket parameters Λ, and the electronic state coefficient The parameters are defined as follows: and are the indices of the nth parameter of gaussian wavepacket and , and are the indicies of the gaussian wavepacket. is the overlap matrix between gaussian wavepackets, ! the nuclear Hamiltonian, the time-derivative overlap matrix, the DD-vMCG C matrix (see ref. 2 main paper), !" the electronic hamiltonian and the Ehrenfest potential. The exponent 0 means derivative of the gaussian function on the left side with respect to the th parameter.
In Qu-Eh, the nuclear Hamiltonian is evaluated using a local harmonic expansion around the centre of the individual gaussian wavepacket. The 2 nd order expansion is done by evaluating the energy, gradient and Hessian at each gaussian centre for each nuclear time step. The main steps in the algorithm for Qu-Eh for each time step are: 1. The geometry for the centre of each gwp in cartesian coordinates corresponds to the usual electronic structure input 2. Integrate equations of motion for electronic state coefficients with electronic structure method (electron dynamics): 2.1. Propagate the configuration interaction coefficients with CAS-CI 2.2. Evaluate the Ehrenfest energy 2.3. Evaluate the Ehrenfest gradient and Hessian (including all off diagonal terms and solving the coupled perturbed equations for orbial and CI parameters 3. Build nuclear Hamiltonian in the gaussian wavepacket basis using a local harmonic approximation 4. Integrate the equation of motions for gaussian parameters and expansion coefficients 5. Repeat from 1

Supplementary Note 2
The initial "stationary" nuclear wavepacket is built by using a Gaussian function centred at the intial geometry using a local harmonic approximation (from the 1 st and 2 nd derivatives obtained from a frequency calculation at that geometry) of the neutral ground state S 0 potential energy surface. The 25 gwp are distributed in phase space with a momentum distribution where the gwp form a shell around the centre position. The gwp 1 has no initial momentum and the gwp 2 to 12 are placed in phase space by initially exciting one degree of freedom in a sequential manner, i.e. they have an initial momentum put in one specific normal mode (+ and -to give 25 gwp). The width of the gwp is 0.1 (in dimensionless unit).
In our computations, the expectation value of the momentum for each normal mode, <p> is zero. However, each individual gwp has a momentum and each gwp is initially associated with a normal mode. There are 2 gwp associated with related each vibrational coordinate, one with +p and one -p where p =11.77. Note that p is related to the dimensionless mass-frequency scaled normal coordinates and has units of inverse time. The value is chosen so that the overlap with the central gwp is 0.8. The value of p for each gwp is identical and the initial momentum for gwp 1 is 0.
The initial weight of each gwp is defined by fitting the initial "stationary" nuclear wavepacket as accurately as possible. By using a smaller width (.1 rather than .707) for the gwp, the initial nuclear wavepacket is also narrowed.

Symmetries of the MO of the 8 cationic states
Figure S1 Symmetries of the MO of the 8 cationic states of benzene cation : Shown are the MO from which the electron is removed and the corresponding symmetries in D 6h and D 2h (we choose the principal axis in both groups to be the z-axis). The numbers 1-8 correspond to energy ordering. (Except that the pairs of states D 6 and D 7, B 3 B 4 and X 2 X 1 , have the same energy and belong to the degenerate E 1u or E 1g representations) Supplementary Note 4 Normal mode symmetries and direct products Figure S2 Normal mode symmetries and direct products of electronic states that yield these symmetries. The nm are ordered by energy (from ground state B3LYP frequency computation). The numbers in brackets correspond the modes used in supplementart reference 2.

Electron dynamics E B 3
In figure S3 we illustrate the spin density oscillations E + B 3 on the 6 C atoms together with a schematic illustration of the electron density (parts a b and c). Note the node through atoms 1 and 4 . Once can observe that the spin density on atoms 3 and 5 are out of phase Notice the spin density on atoms 2 and 3 is in phase and out of phase between atoms 2 and 6 or 3 and 5 in agreement with part c.