Few-femtosecond passage of conical intersections in the benzene cation

Observing the crucial first few femtoseconds of photochemical reactions requires tools typically not available in the femtochemistry toolkit. Such dynamics are now within reach with the instruments provided by attosecond science. Here, we apply experimental and theoretical methods to assess the ultrafast nonadiabatic vibronic processes in a prototypical complex system—the excited benzene cation. We use few-femtosecond duration extreme ultraviolet and visible/near-infrared laser pulses to prepare and probe excited cationic states and observe two relaxation timescales of 11 ± 3 fs and 110 ± 20 fs. These are interpreted in terms of population transfer via two sequential conical intersections. The experimental results are quantitatively compared with state-of-the-art multi-configuration time-dependent Hartree calculations showing convincing agreement in the timescales. By characterising one of the fastest internal conversion processes studied to date, we enter an extreme regime of ultrafast molecular dynamics, paving the way to tracking and controlling purely electronic dynamics in complex molecules. Attosecond science is beginning to provide the tools to study the previously unattainable crucial first few femtoseconds of photochemical reactions. Here, the authors investigate extremely rapid population transfer via conical intersections in the excited benzene cation, both by experiment and computation.

T he advent of attosecond technology 1 has opened new frontiers in contemporary physics and chemistry aimed at observing the motion of electrons in atoms, molecules and solids on their natural timescale 2 . In molecular physics, a lot of attention is currently centered around hole migration in excited molecular cations -an ultrafast propagation of charge through a molecular structure induced by the attosecond creation of a coherent superposition of suitable electronic states [3][4][5][6][7][8] . These dynamics are often considered as purely electronic, however in a real molecular system the electronic processes are not independent from accompanying nuclear dynamics. Recent theoretical studies have focused on how the motion of a nuclear wavepacket 9, 10 and nonadiabatic effects 11 lead to decoherence of hole wavepackets in the benzene molecule and benzene derivatives. Motivated by this work, the present study investigates ultrafast coupled electronic-nuclear dynamics in the excited benzene cation with unprecedented time resolution both experimentally and theoretically.
Near conical intersections the natural timescales for electronic and nuclear motion become comparable, leading to a breakdown of the Born-Oppenheimer approximation. In this regime electronic and nuclear motion are strongly coupled 12 . The experimental realisation of extremely short extreme ultraviolet (XUV) pulses by means of high-order harmonic generation (HHG) using few-cycle visible/near-infrared (VIS/NIR) pulses sparked recent experimental and theoretical investigations of such dynamics in cationic molecular systems, for instance in CO þ 2 13 and small hydrocarbon ions [14][15][16] . An accurate description of coupled electronic and nuclear dynamics remains however a formidable theoretical challenge, as it requires all nonadiabatic couplings between potential energy surfaces to be taken into account. The benzene cation is a prototypical complex system in which nonadiabatic dynamics are expected on the few-femtosecond timescale, as a number of conical intersections are located close to the Franck-Condon region. Its dynamics thus became a focal point for modelling, in particular for the multi-configuration timedependent Hartree (MCTDH) method 17 . In a series of publications, one of the co-authors and his coworkers have theoretically investigated the vibronic coupling in the benzene cation [18][19][20] and its substitutes [21][22][23][24][25] . This benchmark multistate, multimode treatment represented the first quantum dynamical investigation of a molecule exceeding three coupled potential energy surfaces.
Electronic structure calculations within the linear vibronic coupling model were used to evaluate the potential energy surfaces and relevant curve crossings for the lowest five cationic states e X-e E in the benzene cation, as schematically depicted in Fig. 1 19 . The e X 2 E 1g , e B 2 E 2g , and e D 2 E 1u states are doubly degenerate, leading to eight component states in total. All of the electronic states investigated are interconnected by low-energy conical intersections. MCTDH simulations established that population in the e E 2 B 2u state is transferred to the e D state and subsequently to the e B state through two conical intersections that can be reached without significant energy barriers. Note, the e Dand e C 2 A 2u -state potentials intersect at high energies and far away from the equilibrium geometry of the cation (cf. Fig. 1) and thus e D-e C internal conversion is not expected to be a relevant decay channel. The lower lying e C → e B → e X states show a similar internal conversion pathway upon initial population of theC state 18,20 .
These in-depth theoretical studies are contrasted by a lack of validating time-resolved experiments. Until recently such studies were experimentally out of reach. To prepare the lowest cationic states from the neutral molecule, photon energies in the XUV frequency range are required. In addition, for experiments to be meaningful, a time resolution of a few femtoseconds is necessary, due to the ultrashort excited-state lifetimes involved.
Here we present a time-resolved experimental study combined with an extended theoretical treatment, which validates the application of the MCTDH theory for the complex dynamics of ultrafast relaxation in the benzene cation. We prepare a superposition of the lowest five electronic states with a spectrally filtered few-femtosecond XUV pulse. The internal conversion dynamics are probed by further exciting the prepared cations with a delayed few-cycle VIS/NIR probe pulse and by measuring specific fragment yields as a function of the XUV-VIS/NIR delay. MCTDH calculations of the internal conversion are performed that are tailored to the experimental conditions, with wavepackets initiated in the relevant e E and e D cationic states, and are quantitatively compared to the experimental results.

Results
Preparation of benzene cations in excited states. The XUV pulses in the experiment were produced by HHG in atomic xenon and filtered with a thin indium metal foil (see Methods section). They consist mainly of the 9th harmonic. Its spectrum is centered around a photon energy of ħω = 16.0 eV as displayed in Fig. 2a Also shown in this figure are partial photoionisation cross sections of benzene, as previously determined by angle-resolved photoelectron spectroscopy 27 . The lowest five cationic electronic states are prepared with significant population in our experiment. For preparation of the lower-lying e X, e B, and e C electronic states, the photoelectron carries away a significant part of the photon energy 28 . The highest-lying states, e E and e D, are prepared with moderate vibrational excitation 28 . Note that compared to the   26 by a two-photon process (orange arrows). The dashed green and light blue curves are a cartoon drawing of the time-evolution of a cation originally transferred to the e E 2 B 2u and to the e D 2 E 1u state, respectively. Internal conversion processes via the conical intersections indicated in the figure lead to population of the e B 2 E 2g state. The potential energy surfaces are reproduced from ref. 19 ionisation cross sections, the cross sections for producing superexcited neutral states are negligible 29 .
Ionisation of benzene with the present XUV spectrum leads to fragmentation as the photon energy exceeds the appearance energy for unimolecular dissociation of ħω > 13.8 eV 26 . Five fragments are observed in the mass spectrum (see Supplementary Note 2), of which the most relevant section is displayed in Fig. 2b (black line), in agreement with tabulated ion appearance energies 26 . The fragmentation processes in benzene cations have previously been characterised as statistically driven, occurring after relaxation to the electronic ground state (see refs 26,[30][31][32] and references therein). Therefore the observed fragmentation patterns are determined by the total internal energy of the cations upon photoionisation.
Time-resolved probing of dynamics in benzene. Probing of the XUV-prepared benzene cations with a VIS/NIR pulse at timeoverlap leads to the appearance of two additional fragments, C 4 H þ 2 and C 4 H þ 3 (see mass spectrum in Fig. 2b (red line)). These two fragments also originate from unimolecular dissociation of benzene cations and have appearance energies of ħω = 17.5 eV 26 , i.e. above the upper edge of the XUV spectrum in Fig. 2a. Note that double ionisation does not play a role in these experiments as the second ionisation potential is 24.9 eV 33 . We thus conclude that the VIS/NIR pulse induces cation-cation transitions from electronic states initially prepared by the XUV. This allows us to investigate the relaxation dynamics of cations produced by the XUV pulse by recording two-colour mass spectra as a function of the pump-probe time delay Δt. The yield of C 4 H þ 3 fragments as a function of time delay is presented in Fig. 3 A positive delay corresponds to the XUV pulse arriving first at the sample. An increase of the fragment yield around time-overlap is observed, followed by an ultrafast decay. At long delays (≥500 fs) the yield is found to decrease to the level observed before Δt = 0 (cf. inset in Fig. 3). The small contribution observed at negative time delay in Fig. 3 probably originates from imperfect dispersion compensation of the VIS/NIR pulse, resulting in a small postpulse.
The transient in Fig. 3 cannot be described by a single exponential decay. Instead, it is found that a fit with a biexponential decay function, convoluted with a Gaussian function which represents the finite duration of the laser pulses, resembles the data excellently (black line). The two contributions to the fit have similar amplitudes and are shown as dashed lines in the figure. Time constants of τ 1 = 11 ± 3 fs and τ 2 = 110 ± 20 fs for the two decay processes are derived from the fit. The full width at half maximum of the Gaussian function is also found from the fit and amounts to τ IRF = 10 ± 1 fs. As discussed below this instrument response function reflects a process involving the absorption of one XUV and two VIS/NIR photons dominant in our experimental signals. The fragment C 4 H þ 2 exhibits similar biexponential dynamics, but with an inferior signal-to-noise level. The decay times are consistent with those for the C 4 H þ 3 fragment. Furthermore, we find the decay times not to vary with the VIS/ NIR probe intensity.
It is instructive to determine the number of photons that are absorbed in the probe step. For this purpose two-colour mass spectra were recorded as a function of the probe intensity. As the pump-probe signals of C 4 H þ 2 and C 4 H þ 3 consist of contributions featuring different decay rates and comparable amplitudes, measurements were carried out both at time-overlap (Δt = 0 fs) as well as at Δt = 80 fs. Linear fits of the probe power-dependent fragment yields in a double logarithmic plot result in slopes of m Δt=0 fs = 2.0 ± 0.3 and m Δt=80 fs = 2.1 ± 0.1 for the fragment C 4 H þ 3 (see also Supplementary Note 3). These two values show that both contributions of the biexponential decay originate from a two-photon probing mechanism.

Discussion
A benzene cation needs to possess a total energy of at least 17.5 eV above the neutral ground state in order for fragmentation producing C 4 H þ 3 and C 4 H þ 2 to be possible (see Fig. 1) 26 . As indicated by the experiment, the probe laser contributes the energy of two VIS/NIR photons (≈3.2 eV) to this cationic energy. Thus, it follows that the initial energy of the cation produced by the XUV ionisation must be at least 14 40 45 50 population, as it depends on different VIS/NIR-induced cationcation transitions (see Supplementary Note 4). Note that in contrast to time-resolved photoelectron spectroscopy 34 , the VIS/ NIR probe step has to be resonant in the present case, which is facilitated by the broad frequency spectrum of the probe pulse. When the XUV-excited cations relax by internal conversion, energy is transferred from electronic to vibrational degrees of freedom. In particular internal conversion from the e D to the low-lying e B state is expected to sharply reduce the probe-induced signal as a lot of energy is transferred from electronic to nuclear degrees of freedom (Fig. 1), leading to a delocalization of the vibrational wavepacket. Following the Condon principle, cation-cation probe transitions would then lead to highly vibrationally excited levels. The limited availability of suitable electronic states in the relevant energy range (i.e. with vibronic ground states at ≈11.6 eV + 2 × 1.6 eV = 14.8 eV) renders such transitions unlikely (cf. Fig. 1 main text and Supplementary  Fig. 1a). Calculated transition probabilities from the e B to the e D, e E, e F, and e G states are found to be all negligible (see Supplementary Note 5).
Note that in our experiment the lower-lying e D state is both directly populated by the XUV and indirectly through decay of the XUV-populated e E state. The rate equations for relaxation are solved by a double-exponential decay with the same time constants that would result from a pure preparation of the higherlying e E state. We thus can attribute the two exponential time constants observed in the experiment as being associated with relaxation dynamics via the e E ! e D and e D ! e B conical intersections.
The assignments made above can be tested against MCTDH calculations capable of tackling the complex problem of benzene cation relaxation. Following the method described in ref. 18 (see Methods section for a brief description) MCTDH calculations were performed taking into account the e E 2 B 2u state, the degenerate e D 2 E 1u states and the degenerate e B 2 E 2g states. Other lowlying states are expected to be irrelevant in the present conditions due to the two-photon probing mechanism established above. The previous results of ref. 18 were extended in two different ways. First, the initial electronic state was chosen both as the e E state and as the e D state of the cation. This is necessary because the pump laser creates both states in the experiment. Second and more important, in the previous work only diabatic electronic populations were considered because the use of diabatic electronic wavefunctions greatly facilitates quantum-dynamical computations; they were thus underlying the Hamiltonian of ref. 18 . Conceptually, however, adiabatic populations are more instructive because they underlie common thinking of how molecular motion proceeds and how molecular structure is defined. Adiabatic electronic wavefunctions are uniquely characterised through eigenvalues of the molecular potential energy operator (fixednuclei electronic hamiltonian) and are naturally obtained in the first place by quantum-chemical computations. We have hence also extracted adiabatic populations from our MCTDH calculations, which would have been impossible at the time of the earlier work reported in 18,19 due to the great numerical cost. We will publish a more detailed description elsewhere and confine ourselves here to describing the comparison with the experimental data. A wavepacket initiated in the e E state is found to undergo internal conversion to the e D state within a fraction of the vibrational period (≤20 fs), whereas cations prepared in the e D state decay towards the e B state within a few hundreds of fs (see Supplementary Note 6). The calculated population dynamics thus further substantiate our interpretation that the ultrafast timescale of τ 1 = 11 ± 3 fs, observed in the experiment, originates from the decay of e E state population via theẼ !D conical intersection while the observed slower decay process with τ 2 = 110 ± 20 fs originates from the decay ofD state population via theD !B conical intersection. In the calculation, a considerable fraction of population remains trapped in theD state, in contrast to the known absence of fluorescence in the benzene cation 35 . This trapped population is an artefact both of the employed diabatic state basis and of unaccounted vibrational modes. The trapped population gets smaller in the adiabatic representation and has been removed before further analysis.
It should be noted that the link between the electronic state populations and the experimental observables is indirect. Nonetheless we argue that the transient behaviour of the experimental signals is mainly influenced by population dynamics and thus the experimental and theoretical timescales can be compared. The XUV pulse prepares both theẼ andD states of the cation. Therefore we construct a superposition of the two aforementioned wavepacket calculations and compare it with the experiment in the following way: the fractional initial populations of thẽ D andẼ cation states upon XUV ionisation are obtained as p 0 ðẼÞ = 36% and p 0 ðDÞ = 64% from Fig. 2b by convoluting the respective partial photoionisation cross sections of these two states with the measured XUV pump spectrum. The total time-dependentẼ andD state populations are hence calculated by adding the respective populations from theẼ-andD-state-initiated wavepacket calculations, using weights of 36% and 64%, respectively. Finally, the totalẼ andD state populations are added, weighted by a relative probing efficiency ofẼ andD, EffẼ/ EffD, to obtain a theoretical 'fragment yield' curve. EffẼ/EffD, which incorporates the overall relative efficiency to form C 4 H þ 3 fragments in the probe step, is the only adjustable parameter. In our analysis, we ignore any potential influence of vibrational dynamics within theẼ andD states on the probing efficiencies. We justify this approximation by the extensive averaging over the coordinate-dependence of Franck-Condon factors due to the many internal coordinates involved. The fragment formation efficiency is determined by the two-photon cation-cation transition probabilities, weighted with the spectral intensity of the probe pulse, and by the final state-dependent relative probability of C 4 H þ 3 formation. A detailed analysis of the probe step, including computed cation-cation transition probabilities, is presented in Supplementary Note 5 (see Supplementary Fig. 1). We find that while more final states are being accessed with significant probability from theD state, the transitions originating from theẼ state lead to higher-lying final states, from which the probability for dissociation towards C 4 H þ 3 is higher. The amplitudes of the two contributions in the experimentally measured pump-probe trace indicate a value of EffẼ/EffD = 2. This ratio is consistent with our computation of the relevant cationcation transition moment amplitudes (see Supplementary Note 5).
The theoretical 'fragment yield' curve resulting from this summation is depicted in Fig. 4. The solid light blue line represents the result from the adiabatic calculation, while the dashed dark blue line represents the result from the diabatic calculation. Adiabatic and diabatic results are similar and well described by a biexponential fit, which yields time constants of τ 1 = 8/12 fs and τ 2 = 170/190 fs (adiabatic/diabatic). These results are weakly dependent on the assumed relative probing efficiency of theD andẼ states. Variation of EffẼ/EffD within limits that reasonably resemble the experimental behaviour leads to decay times in the range of τ 1 = 6-11/11-15 fs and τ 2 = 160-180/180-200 fs (adiabatic/diabatic). An oscillation of the yield derived from the diabatic calculation with a period of about 35 fs is visible, which is absent in the adiabatic result, as well as a pronounced maximum at about 50 fs. These features likely emerge from variations of the adiabatic-to-diabatic mixing angle and are considered to be an artefact of the diabatic theory, observed also in previous two-state ARTICLE NATURE COMMUNICATIONS | DOI: 10.1038/s41467-017-01133-y calculations 36 . The theoretical 'yield' curve for the adiabatic calculation is further convoluted with the Gaussian instrument response function of τ IRF = 10 fs and plotted on top of the experimental result in the inset of Fig. 4. The agreement between experiment and theory is striking, strongly supporting the quantitative validity of the MCTDH calculations. We note that a larger variation at longer time scales is not unexpected since only a selection of vibrational states is included in the simulation. This observation might also hint at a potential coordinatedependence of the transition operator associated with the experimental probe step. We emphasize that this is the first experimental verification of these large-scale computational results, originally published as benchmark calculations and as a purely theoretical prediction.
The comparative study presented here demonstrates that complex nonadiabatic relaxation processes inevitably accompany -already at very short times -the purely electronic dynamics sought after in charge-migration experiments. These relaxation processes may span several electronic states, coupled by many vibrational degrees of freedom. In the present case of the benzene cation, theẼ !D internal conversion limits the purely electronic part of the dynamics to a timescale of about 10 fs, after which time nuclear dynamics and the associated dephasing of the electronic wavepackets cannot be neglected (see also refs 13,37 ). To separate these two effects, experimental and theoretical access to the time-evolution of the involved electronic states is crucial. We demonstrated that monochromatized short XUV laser pulses combined with few-cycle VIS/NIR laser pulses allow for a stateby-state analysis of the population evolution. Moreover, we established that MCTDH calculations are suitable to predict these population evolutions. Our experimental and theoretical approach combined with the ongoing development of frequencytunable few-femtosecond XUV sources 38,39 should serve as a benchmark for future experiments towards the observation of attosecond dynamics.
In conclusion, we have experimentally determined decay rates of cationic states of benzene that are produced by XUV ionisation of the neutral molecule. Two timescales of τ 1 = 11 ± 3 fs and τ 2 = 110 ± 20 fs are observed and interpreted asẼ !D andD !B internal conversion timescales, respectively. TheẼ state relaxation timescale marks one of the fastest internal conversion processes observed to date in a time-resolved experiment.

Methods
Experimental methods. The experiments were carried out with a two-colour XUV-VIS/NIR pump-probe setup, which has already largely been described elsewhere [40][41][42] . The experimental apparatus was upgraded for the present study by compressing the output pulses of a 1 kHz Ti:Sa amplifier by means of a hollow-core fibre setup and a set of chirped mirrors. Few-cycle pulses with a bandwidth of 640-920 nm (FWHM) and a pulse duration of τ = 6 fs, determined with a SEA-F-SPIDER measurement 43 , were obtained. The beam is split into two parts. The first part is used to generate high harmonic radiation in a gas cell containing xenon. The resulting broadband XUV spectrum is then filtered with a 200 nm thick indium foil. This limits the spectral bandwidth 44 and results in the spectrum depicted in Fig. 2a. The XUV photon energy was calibrated by measuring the absorption spectrum of molecular nitrogen (see Supplementary Note 7). The XUV beam was recombined with the second part of the VIS/NIR beam, which acts as a probe, by use of a cored mirror. Using glass wedges for dispersion control, the two arms of the setup were then fine-tuned: (i) the VIS/NIR arm was optimised using argon strong field ionisation yield, and (ii) the XUV arm was optimised using HHG efficiency and spectrum. This allowed us to reach a situation were essentially a single harmonic was transmitted through the Indium metal filter, which was not significantly spectrally cut. The two recombined, collinear beams are then focused by a toroidal mirror into the center of a velocity-map-imaging spectrometer, operated in time-of-flight mode 45 . The benzene sample was provided by a thin supply pipe, with its orifice in close proximity to the interaction region. The VIS/ NIR probe beam intensity was set to approximately 10 12 W cm −2 . Due to the multiphoton nature of the molecular response and potentially not fully compensated pulse propagation effects, we determined a Gaussian instrument response function from a fit of the delay-time dependent ion yield. This resulted in a value of τ IRF = 10 ± 1 fs. Note that our data cannot be fitted convincingly by the response function and a single exponential decay.
Theory methods. The time-dependent electronic populations used in this work are obtained in a fully quantal way with the aid of the Multiconfigurational Time-Dependent Hartree (MCTDH) method. This is a highly efficient tool for solving the time-dependent Schrödinger equation for several (or many) degrees of freedom. The standard method of solution would be the numerically exact propagation of a wavepacket represented in a time-independent product basis set. In the MCTDH scheme the wave function is instead written as a linear combination of Hartree products as follows: Here f is the number of degrees of freedom, q 1 ,…,q f are a set of nuclear coordinates, A j1 jf denote the time-dependent expansion coefficients, and ϕ ðkÞ jk ðq k ; tÞ are a set of time-dependent single-particle functions (SPFs) combined to the Hartree product Φ J . In turn, the SPFs are expressed in a time-independent basis set: where χ ðkÞ ik ðq k Þ are primitive basis functions of the k-th degree of freedom that depend on the particle coordinate q k . The accuracy of MCTDH depends on the number of primitive functions N 1 ,…,N f and the number of SPFs n 1 ,…,n f used. Since both the coefficients and the basis functions in Equation (1) are timedependent, they both are optimised using a variational principle. The equations of motion are then derived from the Dirac-Frenkel variational principle 46,47 . This choice keeps the product representation of the wavepacket optimally short and reduces the length of the state vector by about six orders of magnitude in typical applications (MCTDH contraction effect).
The implementation proceeds along the lines of ref. 18 where a model Hamiltonian for theB-D-Ẽ states of C 6 H þ 6 was set up within a linear vibronic coupling scheme and the ionisation potentials, vibrational frequencies, and coupling constants were determined from ab initio calculations. Up to nine nuclear degrees of freedom have been considered in the calculation, including the strongly Condon-active C-C stretching mode ν 2 (a 1g symmetry) and the Jahn-Teller active modes ν 16 -ν 18 (e 2g symmetry). In view of the degeneracy of theB andD states five component electronic states (five potential energy surfaces) have been treated. All basis set details have been collected in Tab. IV.b of ref. 18 , which shows that 10 11 primitive basis functions are reduced to 10 6 time-dependent basis functions by the MCTDH contraction. The resulting electronic populations of the interacting states are depicted in Fig. 8b of ref. 18   The wave packet has been propagated for 500 fs (i.e. much longer that the 200 fs of ref. 18 ) and also for theD component states as initial state (while in ref. 18 only thẽ E state was considered as initial state). In addition to the diabatic populations also the adiabatic ones have been computed (see Supplementary Note 6 for further details): the adiabatic representation is the more natural one to be used for the interpretation of the nuclear dynamics (for reasons given above). Our simulations allow to obtain a suitably weighted mean of the detectedD andẼ state populations (see above) and hence to implement the proper initial conditions in line with the present experiment.
Data availability. The data that support the findings of this study are available from the corresponding author upon request.