lopscience

¡opscience.iop.org

Home Search Collections Journals About Contact us My lOPscience

Femtosecond dynamics of correlated many-body states in Cgg fullerenes

This content has been downloaded from lOPscience. Please scroll down to see the full text.

View the table of contents for this issue, or go to the journal homepage for more Download details: IP Address: 91.200.81.30

This content was downloaded on 28/01/2017 at 07:19 Please note that terms and conditions apply.

You may also be interested in:

Multiscale dynamics of C60 from attosecond to statistical physics F Lépine

Theoretical methods for attosecond electron and nuclear dynamics: applications to the H2 molecule Alicia Palacios, José Luis Sanz-Vicario and Fernando Martin

Advances in attosecond science

Francesca Calegari, Giuseppe Sansone, Salvatore Stagira et al.

Charge migration induced by attosecond pulses in bio-relevant molecules Francesca Calegari, Andrea Trabattoni, Alicia Palacios et al.

Ultrafast molecules and clusters I V Hertel and W Radloff

Time-, angle- and kinetic-energy-resolved photoelectron spectroscopy of highly excited states of NO Peter Trabs, Franziska Buchner, Masood Ghotbi et al.

Nonlinear effects in photoionization over a broad photon-energy range within the TDCIS scheme Antonia Karamatskou

Controlling the dissociation dynamics of acetophenone radical cation through excitation of ground and excited state wavepackets

Katharine Moore Tibbetts, Maryam Tarazkar, Timothy Bohinski et al.

Theory of strong-field attosecond transient absorption Mengxi Wu, Shaohao Chen, Seth Camp et al.

New Journal of Physics

The open access journal at the forefront of physics

Dautsdie PhyilbUiicha GwalUdiaft DPG

IOP Institute of Physics

Published in partnership with: Deutsche Physikalische Gesellschaft and the Institute of Physics

Cross Mark

OPEN ACCESS

RECEIVED

5 October 2016

ACCEPTED FOR PUBLICATION

3 November 2016

PUBLISHED

29 November 2016

Original content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence.

Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation andDOI.

Femtosecond dynamics of correlated many-body states in C6q fullerenes

Sergey Usenko1, Michael Schüler2, Armin Azima3,4,5, Markus Jakob1,5, Leslie L Lazzarino3,

Yaroslav Pavlyukh2,6, Andreas Przystawik1, Markus Drescher3,4,5, Tim Laarmann1,5,7 and Jamal Berakdar2,7

1 Deutsches Elektronen-Synchrotron DESY, D-22607 Hamburg, Germany

2 Institut für Physik, Martin-Luther-Universität Halle-Wittenberg, D-06099 Halle, Germany

3 Department of Physics, University of Hamburg, D-22761 Hamburg, Germany

4 Center for Free-Electron Laser Science CFEL, DESY, D-22607 Hamburg, Germany

5 The Hamburg Centre for Ultrafast Imaging CUI, D-22761 Hamburg, Germany

6 Department of Physics and Research Center OPTIMAS, University of Kaiserslautern, D-67653 Kaiserslautern, Germany

7 Authors to whom any correspondence should be addressed.

E-mail: tim.laarmann@desy.de andjamal.berakdar@physik.uni-halle.de Keywords: fullerene, ionization dynamics, superatom molecular orbitals, TDDFT

Abstract

Fullerene complexes may play a key role in the design of future molecular electronics and nanostructured devices with potential applications in light harvesting using organic solar cells. Charge and energy flow in these systems is mediated by many-body effects. We studied the structure and dynamics of laser-induced multi-electron excitations in isolated C60 by two-photon photoionization as a function of excitation wavelength using a tunable fs UV laser and developed a corresponding theoretical framework on the basis of ab initio calculations. The measured resonance line width gives direct information on the excited state lifetime. From the spectral deconvolution we derive a lower limit for purely electronic relaxation on the order of Te\ = 10+ fs. Energy dissipation towards nuclear degrees of freedom is studied with time-resolved techniques. The evaluation of the nonlinear autocorrelation trace gives a characteristic time constant of Tvib = 400 ± 100 fs for the exponential decay. In line with the experiment, the observed transient dynamics is explained theoretically by nonadiabatic (vibronic) couplings involving the correlated electronic, the nuclear degrees of freedom (accounting for the Herzberg-Teller coupling), and their interplay.

1. Introduction

Molecular junctions, molecular transistors and organic solar cells rely on charge transport channels with negligible energy dissipation during the carriers propagation time. In nanostructured materials and molecular complexes the characteristic timescale is determined by the long-range polarization interaction and by the formation and breaking of chemical bonds mediated by the electronic and nuclear motion. Transient structures and dynamics on the femto and sub-femtosecond timescale is the focus of ultrafast spectroscopy. Time-resolved experiments using femtosecond (fs) laser pulses unravel the dynamic response of promising materials that could serve for instance as molecular building blocks for organic photovoltaics. Polymer solar cells are commonly composed of a photoactive film of a conjugated polymer donor and a fullerene derivative acceptor [ 1-3], which makes use of the fullerenes' unique ability to form stable C -0 anions. Electron correlation plays an important role in the formation of four bound states of the fullerene anion [4-6]. In fact, electronic correlations are responsible for the binding of the 2Ag state, whereas the bindings of the states 2T1u, 2T2u and 2T1g are less affected by electronic correlations (cf [4] and further references therein).

With its special structure consisting of 174 nuclear degrees of freedom, 60 essentially equivalent delocalized n electrons, and 180 structure-defining localized a electrons, neutral C60 serves as a model for a large—but still finite—molecular system with many electronic and nuclear degrees of freedom. Because of the large charge conjugation, its finite 'energy gap', and quantum confinement of electronic states, C60 maybe viewed as an

© 2016 IOP Publishing Ltd and Deutsche Physikalische Gesellschaft

interesting intermediate case between a molecule and a condensed matter system. In fact, applying solid-state concepts to the valence 'Bloch electrons' on the C60 sphere results in an 'angular band structure' [7] from which other relevant quantities (such as plasma frequencies and group velocities) can be extracted. Photophysical studies of fullerenes using fs laser fields cover the whole range from atomic through molecular to solid state physics [8-10]. The molecular response is truly a multi-scale phenomenon. It ranges from attosecond dynamics in electronic excitation and ionization to statistical physics describing thermalization processes. So, light-induced processes in fullerenes cover more than 15 orders of magnitude in time [11].

Using low-temperature scanning tunneling microscopy (LT-STM) of C60 molecules deposited on copper surfaces Feng etal observed tunneling through electronic states that possess nearly atom-like character [12]. These 'superatom' molecular orbitals (SAMOs), also discussed below, have a well-defined symmetry and can be characterized by the nodal structure (principle quantum numbers) n and angular momentum quantum numbers L [13]. In addition, these virtual states show a remarkable stability [14], i.e., their initial stage of decay proceeds substantially slower than other states (even LUMO or HOMO) which qualifies them as robust channel for hot electron transport. From a molecular physics point of view, SAMOs can be regarded as low-lying, mixed valence Rydberg states [15] that exhibit substantial electron density inside the hollow sphere. They form chemical bonds affected by hybridization when the system is excited optically or probed with STM in deposited nanostructures. The resulting free-electron bands in self-assembled one-dimensional wires or two-dimensional quantum wells are holding great promise for unique applications in molecular electronics [ 16], but also for new functionalities such as current carrying states and hence nanometer-sized magnetic field generators [17]. We note in passing that recently SAMO states have also been observed in planar, non-fullerene materials [ 18,19],as well as in isolated C60 fullerenes [20-22], where they are energetically located below the known high-lying Rydberg states [23,24].

It is of fundamental importance for designing fullerene derivatives as building blocks for solid state chemistry to go beyond the characterization of static many-body electronic structure [25]. In particular for optimization and control of the charge flow and energy dissipation, rigorous dynamic studies on the fs time scale using ultrafast lasers are indispensable but, still in their infancy [26,27]. Additionally, time-dependent density-functional theory (TDDFT) calculations of the absorption and photoelectron spectra, accounting for full structural analysis, were performed [28,29] and constitute the basis for the further development presented below.

Accessing energy dissipation upon the excitation of correlated many-body states in gas phase C60 became feasible recently which allows to connect the coherent quantum [28,30], classical and statistical mechanisms [31]. It is known from experiments with optical lasers that electron thermalization mediated by inelastic electron-electron collisions takes place on a time scale of50-100 fs (see [11] and references therein). This is where the present experimental and theoretical work comes into play. Here, the objective is to reduce the complexity of laser-induced multi-photon processes by populating highly-excited many-body states in a resonant one-photon transition in the ultraviolet (UV) spectral range at low laser peak intensity on the order of 3.5 x 1010 Wcm-2. This allows for a rather detailed probing ofthe correlated electron dynamics in highly excited states. The study is based on a resonance-enhanced multi-photon ionization (REMPI) scheme, i.e., two-photon photoemission (2PPE) spectroscopy as depicted in figure 1(a). The photoionization yield recorded for resonant excitation is enhanced as compared to an experiment performed in the off-resonance regime. Thereby, we trace the time-dependent electronic structure of intermediate states, which are free from any perturbation caused by metallic substrates affecting the energetics in LT-STM experiments. Furthermore, REMPI on gas phase fullerenes provides information on the neutral molecule whereas 2PPE and LT-STM essentially probe the binding energy and DOS of an anion deposited on a metal surface. Our experiments are compared to ab initio calculations.

This paper is organized as follows. In section 2 many-body states below the C60 ionization threshold are calculated as a guideline for the 2PPE experiments. Section 3 describes some experimental details with a focus on the tunable fs laser system in the UV spectral range used for the time-resolved studies. Excitation energy dependent mass spectra are evaluated in section 4.1 and discussed in terms of resonance-enhanced ionization and excited state lifetimes. Section 4.2 concentrates on the time-resolved experiments. A detailed theoretical analysis of the experimental data is given in section 5 followed by a short summary and outlook. In appendix A we provide details on how the electron-vibron coupling matrix elements were computed and the master equation in Lindblad form is derived in appendix B.

2. Optical excitations

The first step of a REMPI or 2PPE experiment entails the calculation of excited states, which are typically more difficult to describe than ground-state properties. Often exact diagonalization (full configuration interaction) is

(a) 12 10 8

>> 6 2> b o c

0 4 2 0

excitation G (6T1u)^ 3p-SAMO-3S-SAMO-

Figure 1. (a) Excitation spectrum of the C60 molecule starting from the ground state (GS). Optically accessible excited states are denoted as bright states (BS), while excitations with vanishing dipole transition moment are referred to as dark states (DS). The two-photon REMPI experiment is sensitive to BSs between 5 and 6 eV. The first three excitations are labeled according to [39]. (b) Excitation G resolved with respect to the constituting occupied (blue), virtual (orange) andSAMO (purple) electronic states. The thickness of the red arrows is proportional to the corresponding weight \Aa=G \2.

not feasible for large systems in which case one may resort to only a few methods: TDDFT [32], equation-of-motion (EOM) quantum chemistry methods [33], and many-body perturbation theory based on a Green's function formulation [34]. For its reduced computational cost as compared to the other methods, we have employed the linear-response TDDFT approach (Casida's method) [35].

As a first step we calculated the Kohn-Sham (KS) orbitals using the OCTOPUS package [36] with a modified version of the asymptotically corrected functional by Leeuwen and Baerends [37], which was shown to considerably improve excited-state properties [38]. In order to account for a multitude of highly-excited states, we have chosen a relatively large box to which all KS states fi (r) are confined (a sphere with a radius of 12 A with uniform grid spacing of 0.15 A). This ensures that higher virtual orbitals (including the SAMOs) are well represented. After converging the ground-state and computing a sufficient number ofvirtual orbitals, we computed the singly-excited (i.e., single particle-hole excitations) many-body excitations by Casida's method. Formally this procedure amounts to approximating the excited many-body states \ Fa) by

\F«)=EE Apj Q\$o), (1)

i£ occ jGvirt

where \F0) denotes the determinant built by the ground-state KS orbitals and ci (c/) is the annihilation (creation) operator with respect to the KS basis. The particle-hole amplitudes At" are determined by Casida's equation based on linear response. The major approximation hereby is related to the exchange-correlation (xc) kernel fxc (r, rw), defined as the functional derivative of the KS potential with respect to the density. We use the local-density approximation for the xc kernel, as it is local and (within adiabatic TDDFT) frequency-independent. Casida's equation is thus transformed into an eigenvalue problem. We computed the Casida vectors Aiaj with OCTOPUS, taking 75 occupied and 60 virtual orbitals into account, yielding well-converged results for excitation energies up to 10 eV. For testing purposes we also computed the binding energies of the virtual orbitals analogously to [21]. We obtained very similar results for the low-lying states relevant for the present experiments.

The energies of the obtained many-body excitations are shown in figure 1 (a), where we distinguish the states with vanishing dipole transition moment (these we refer to as dark states, DS) from the ground state GS, and optically accessible states (bright states, BS). The onset of the visible to UV (UV-vis) optical absorption is well documented (e.g., cf [39] for a review and a comparison with previous experiments). Three distinct absorption peaks (labeled according to literature as C, E and G bands, respectively) are found in the UV-vis region. An overall agreement between our calculations ofthe optical absorption with the experimental results is found, though the excitation energies are slightly underestimated (which is typical for DFT calculations). Note, the specific experimental conditions may affect the spectral positions of the absorption peaks (as discussed in [39]). In particular, the present experiment probes the optical properties of isolated molecules by ultrafast pulses and thus potentially eliminates energy-loss channels (e.g. collisions and inter- or intra molecular decay) that might shift the absorption peaks to higher energy. For these reasons, we base the subsequent calculations on our DFT

0 2 4 6 8 10 12

angular momentum L

0 2 4 6 8 10 12

angular momentum L

Figure 2. Occupied and virtual KS states close the negative ionization potential (IP) ordered (from bottom to top) according to figure 1 (b). The respective orbital character is illustrated by the relative weight of each angular-momentum component L. Only one representative of degenerate states is shown. (a) Occupied states, (b) virtual states below the SAMOs, and (c) orbitals including the SAMOs and higher states.

Figure 3. Schematic layout of the femtosecond laser pulse generation scheme in the ultraviolet spectral range.

results without any adjustments. As detailed below, very good agreement to the present experiment is achieved, justifying this approach a posteriori.

According to our calculations, excitation G represents one triply degenerate many-body state, which is identified as 6T1u excitation [40]. The projection onto the KS basis is depicted in figure 1 (b), where we illustrate the relative weight of the excitation from occupied (i) to virtual ( j) orbitals by the thickness of the corresponding arrows. The relevant orbitals are also presented in figure 2. As clear from figures 1 (b) and 2, the excitation G is predominantly composed of the transitions (i) from HOMO to virtual states above the first SAMOs with dominant angular momentum of L = 8 and L = 6, and (ii) from HOMO-1 to LUMO+2. The angular-momentum analysis of the individual orbitals (figure 2) clarifies why the transition to the many-body state associated to excitation G is optically allowed. Analogously one can conclude that neither the 3s-SAMO nor the 3p-SAMO can be populated in a direct single-photon dipole transition from HOMO. This is consistent with a full symmetry analysis of the initial orbitals and the ab initio calculations.

3. Experimental setup

In order to study electronic transitions into highly excited many-body states close to the ionization potential, fs pulses at UV frequencies are required [41]. The present time-resolved mass-spectrometric study is based on second-order UV autocorrelation making use of a time-of-flight (TOF) mass spectrometer and state-of-the-art nonlinear optics. The laser setup for generating fs pulses in the range of216-222 nm comprises a Ti:Sa laser system, an optical parametric amplifier (OPA) and several frequency mixing stages. The outline of the system is sketched in figure 3. The Ti:Sa laser (Amplitude Technologies) is the backbone of the overall generation scheme and provides 35 fs (FWHM) pulses with a pulse energy of 12 mJ behind the compressor at a repetition rate of 25 Hz and 800 nm central wavelength. This output is split into several arms. First, the beam is split in a 90:10 ratio. The more intense fraction is sent to a commercial OPA (Light Conversion, TOPAS-C + HE-TOPAS). The OPA

is continuously tunable in the infrared spectral range (1140-3500 nm) and pumped by the Ti:Sa laser. For the subsequent frequency conversion in the UV the output of the OPA is tuned to 1200 nm. The 11 mJ of the 800 nm pump pulse are converted to »0.4 mJ at this wavelength. The low-intensity fraction of the Ti:Sa laser (10%) is equally split into two branches. One half is recombined with the 1200 nm output of the TOPAS in a ^-barium borate (BBO, 0.2 mm thick) crystal to generate 480 nm light by sum frequency generation (SFG). The second half is frequency doubled in another BBO, and then together with the 480 nm beam directed to the third SFG stage (BBO, 40 ,um thick) to finally generate the 5.7 eV (218 nm) photons. The UV pulse energy was measured using a calibrated XUV photodiode and a pyro detector to be »2 ^J.

The output UV beam is directed into a vacuum chamber where it is split into two pulses by a reflective split-and-delay unit (SDU) in order to generate two synchronized pulse replicas. The complete nonlinear optical setup was simulated using the software package LAB2 [42] including dispersion induced by UV pulse propagation in air and through the 2 mm thick entrance window of the vacuum chamber. According to the calculation the 218 nm pulse duration in the interaction region is of the order of 100 fs FWHM with a spectral bandwidth of 2.8 nm. A coarse cross-correlation measurement performed between the 400 nm and 480 nm pulses of 150 fs FWHM supports the derived UV laser beam parameters.

The SDU consists of a Si split-mirror with one half mounted on a delay stage which can displace the mirror along its normal to set the time delay between the pump and the probe pulses. The SDU is followed by a focussing mirror which spatially overlaps the two pulse replicas in the laser-sample interaction area. The laser beams are focused onto the C60 molecular beam with a spherical mirror (f = 300 mm). Its reflectivity is above 80% in the 200-245 nm range. The beam waist in the interaction area is on the order of 150 ^m. The maximum peak intensity reached in the experiments is approximately 3.5 x 1010 W cm-2, which is derived from first principles based on Gaussian beam propagation, pulse energy, pulse duration and the far-field laser profile.

The molecular beam is produced by evaporation of C60 powder (purity >99.5%) in a resistively heated oven at 775 K. The UV laser beams are focused perpendicular to both, the effusive molecular beam and the TOF spectrometer axis. The ions created in the intersection volume are extracted by a static electric field (Wiley-McLaren configuration [43]), directed onto multichannel plates, and finally counted after amplification and discrimination by a digital oscilloscope. The mass resolution of the TOF spectrometer is 0.2% at M/q = 720.

4. Two-photon photoemission

4.1. Excitation energy dependence

The photoionization signal recorded for resonant excitation is enhanced compared to an experiment performed off-resonance. The spectral width of the resonance yields information on the excited state lifetime. A UV wavelength scan was performed by tuning the IR wavelength of the OPA. The populated many-body state in the neutral molecule is subsequently ionized during the pulse duration of 100 fs. A full mass spectrum is accumulated over »250 laser shots for each excitation wavelength. The C+0 yield was normalized to the relative pulse energy monitored by a photodiode. The C +0 ion yield as a function of laser wavelength in the range of 216-222 nm is shown in figure 4. The cut-off at 216 nm corresponds to the OPA's lower wavelength limit of 1140 nm. Relatively large error bars at short wavelengths result from the corresponding low pulse energy and thus poor statistics. The C +0 signal disappears for excitation wavelengths longer than 222 nm, thus representing the low-energy threshold of the resonance. The wavelength scan clearly indicates a resonance-enhanced two-photon photoionization at Aexc = 218 ± 0.5 nm. The observed width oftheREMPI signal is of the order of 3.65 nm (94 meV) pointing towards ultrafast ionization within a lifetime that can be as short as 10+5 fs. The lifetime estimate is derived from the deconvolution ofthe observed resonance with a Gaussian laser pulse spectrum of 2.8 nm (FWHM) and a Lorentzian describing the homogeneous broadening.

4.2. Pump-probe delay dependence

Time-resolved mass spectrometry traces the excited state dynamics in neutral C60 molecules directlyin the time domain. The transient electronic structure is initiated by a UV 100 fs pump pulse and followed (probed) by a delayed pulse replica that photoionizes the molecule (see figure 1 (a)). The single-shot mass spectra are taken at varying delay times between the UV pulses ranging from -60 fs to 900 fs with »3000 laser shots per delay point. The time-dependent on-resonance ion signal shown in figure 5 is derived by taking the average number of C +0 counts for each delay and normalizing it to the relative pulse energy monitored as the UV stray light peak in the TOF spectrum. The pump-probe scan is repeated two times.

In the most general case the ion signal from a three-level system with a transient intermediate electronic state exposed to resonant excitation will result from two excitation pathways: direct two-photon photoionization from the ground state to the continuum and the ionization via the transient state (REMPI). Therefore, the total ion signal Stot can be described as a sum of three components [44]: the coherent term Sac (coherent artifact), the

o '>»

212 214 216 218 220 wavelength [nm]

Figure 4. Normalized C+0 ion yield (black scatter) as a function of excitation wavelength showing the resonance-enhanced two-photon photoionization. The fit to the data points (black line) is a convolution of a Gaussian profile representing the laser pulse spectrum and a Lorentzian profile, representing the natural linewidth of the resonance, respectively. The experimental data is compared with calculations (see subsection 5.3) for different values of the electron-vibron coupling strength g.

pump-probe delay [fs]

Figure 5. Top graph: normalized stray light signal on the TOF detector induced by the excitation pulse as a function of the pumpprobe delay, which monitors the UV laser stability throughout the experiment. The straight dotted line designates 1. Bottom graph: normalized C+0 counts as a function of pump-probe delay for resonant excitation at Aexc = 218 nm. The black scatter are the data points and the curve is a fit obtained from equation (4) (for details see the text). The normalized ion counts derived from theory are compared to the experimental results for different scaling factors g of the electron-vibron interaction strength (see section 5.3) in the limit of long pump-probe delays (^300 fs). The bold value indicates the best agreement between theory and experiment.

incoherent term Sinc and a constant background abg. The coherent artifact reflects the direct nonlinear ionization process and is proportional to the autocorrelation (AC) function of the laser pulse:

I (t - t)I (t) dt, (2)

where I(t) is the laser intensity and t is the delay between the pump and probe pulses. The incoherent term Sinc carries information about the population dynamics of the transient state. It is a convolution of the laser pulse AC with a symmetric decay function. In case of an exponential decay with a characteristic time constant tvib the incoherent term is given by:

Sac (t - t)exp(-|t|/tvib)dt. (3)

The constant background abg consists of contributions from each excitation pulse individually and is independent of the pump-probe delay. The total signal reads as:

— a ac Sac

ainc Sinc (t) + abg (4)

with a; being the relative amplitudes of the different components. In general the amplitudes have a ratio depending on the spatial overlap of the two pulses and the ionization pathway of the system. To extract t^b from

the measurement, the experimental data is fit by the least squares method using expression (4). The background issetto «bg = 1 and the laser pulses are assumed to be Gaussian with tfwhm = 100 fs. Other variables, i.e., aac, ainc, and rvib, are free fit parameters. The best fit curve (black line in the bottom graph of figure 5) yields a time constant rvib = 400 ± 100 fs (95% confidence band) and an amplitude ratio aac : ainc : abg = 0.24 : 1.13 : 1. The laser intensities in the interaction region in the present experiment are as low as 3.5 x 1010 W cm-2 which makes the contribution from the direct (nonlinear) two-photon process small. The ratio ainc : abg is close to 1 as expected for single photon ionization from an occupied transient state.

The observed exponential time constant rvib = 400 ± 100 fs is significantly longer than the characteristic electron-electron interaction time derived from pump-probe spectroscopy [45] and pulse duration dependent studies [46,47] in the optical spectral range. It seems that electron thermalization mediated by inelastic electron-electron collisions on a time scale of50-100 fs does not play a key role when high lying correlated many-body states are excited directly at rather low peak intensity. In the following a detailed theoretical description of the observed resonance and its time-dependent structure is discussed.

5. Simulations

To describe the response of the molecule upon pulsed laser irradiation, the interplay between the electronic excitations and the vibrationally hot molecule has to be taken into account (for an overview on this topic we refer to the book [48] and further references therein). In full generality, an ab initio description for both the electronic and nuclear degrees of freedom and their coupling is not feasible currently. Hence, one has to rely on suitable approximations. In an Ehrenfest approach the electrons are described on the level of TDDFT and the nuclei are subject to classical equations of motion due to the forces exerted by the electron distribution. Despite the success of this molecular dynamics approach [49] for predicting vibrationally-assisted charge-transfer processes [50] in photovoltaics the Born-Oppenheimer (BO) approximation is an inherent limitation. For excited-state properties, where the nuclear and electronic dynamics are strongly mixed, BO-type molecular dynamics is not predictive. Molecular dynamics beyond the BO approximation [51-53] has been employed for the C60 molecule [54]; merging such schemes with a treatment of the electronic excitations beyond the KS level is computationally too demanding for our system. Alternatively, one can treat the electrons in a single-particle atomic basis within a tight-binding model [55], removing the adiabaticity constraint with respect to the many-body states. Besides the inevitable empirical ingredient, such theory is also not directly compatible with the ab initio description of the many-body states in section 2, i.e. electronic correlations can be taken into account only with great efforts.

5.1. Initial laser-induced dynamics

In order to elucidate the laser- and vibration-induced dynamics we take a different angle. Since a considerable amount of energy is stored in thermally activated vibrations, which can only be transferred to the electronic subsystem in small portions, the vibrons can be treated as an effective heat bath for the electrons. A similar model has successfully been employed for incorporating the influence of the vibrations on charge-transfer processes in organic photovoltaic systems based on C60 [56]. To construct an appropriate model for our case, several ingredients are required. For the vibrations we restrict ourselves to the harmonic approximation of the bottom of the BO surfaces. The vibronic eigenmodes along with their eigenfrequencies and reduced masses were computed using the OCTOPUS code, as well. The resulting density of states (DOS) of the vibronic modes is shown in figure 6. Our results compare very well with those tabulated in the literature, for instance table 6.2 in the book [48].

As inferred from figure 6 the high-energy modes (which affect the electrons at most) are only weakly populated for T = 775 K. Therefore, the oscillations of the nuclei around their equilibrium positions can be considered small. Hence, the Herzberg-Teller (HT) expansion [57] of the full Hamiltonian, including electrons and nuclei, yields a reasonable description for both subsystems and their interaction. The first-order HT Hamiltonian amounts to approximating the electron-vibron coupling as linear in the mode amplitudes Qn. On the KS level, the electron-vibron matrix elements are thus given by

K, = (fi\ jQ- if)

Details on the evaluation of equation (5) are provided in appendix A. Since we are opting for a model in the many-body basis of the excitations discussed in section 2, the matrix elements (5) are transformed according to equation (1) (see appendix A). We thus obtain the model Hamiltonian (atomic units are used throughout)

H (t) = Hel (t) + Hel-vib + Hvib, (6)

g 2 ■o

CL * <0 1

Bose distribution

i v J .w

o o o c T3

2 S3-o'

50 100 150

energy [meV]

Figure 6. Vibrational density of states (DOS) along with the occupation according to the Bose distribution for T = 775 K.

Hel (t) = Y,Ea\Fa) (Fa\ + f (t) \Fa) (Fp \, (7)

Hel-vib = gEE^aplFa) (Fp\Qn, (8)

HHvib = -1E

J3 - 2

-7J- + kn Q n

Here, Map are the dipole transition matrix elements from our TDDFT calculations, while f (t) comprises the time-dependent fields. The prefactor gin equation (8) is introduced as an overall scaling factor for the strength of the vibronic coupling. Ideally, g = 1 should be fixed; however, due to the perturbative description resulting from the HT expansion, the electron-vibron interaction might be underestimated. Hence, gis kept as a parameter.

Instead of describing the dynamics of the full density matrix according to Hamiltonian (6) (which is a formidable task), we treat the vibrations as a heat bath. This allows to obtain a master equation for the density matrix in the electronic subspace only. Here we assume Markovian dynamics and thus employ the Lindblad master equation, following the standard derivation and formulation from [58]. More details are presented in appendix B. This procedure requires an additional parameter: the vibronic broadening n, which corresponds to the lifetime ofthe vibrational modes. The laser-induced dynamics is, besides the electron-vibron interaction, treated on ab initio level endorsing the predictive power of the approach. Note that the electron-vibron coupling (8) includes, in principle, Jahn-Teller and Herzberg-Teller couplings, which where identified as the main mechanisms for vibrational coupling in fullerenes [59-61 ].

The time evolution of the occupation of the states depicted in the level scheme figure 1 (a) is presented in figure 7(a). The driving pulse f (t) is chosen as a Gaussian pulse with a FWHM of 100 fs as in the experiment, while the central frequencyis adjusted to the vertical excitation energy DE = 5.66 eVbetween ground state and the BSs corresponding to the G peak. The peak amplitude amounts to the intensity of 3.5 x 1010 Wcm-2.The values for gand n are chosen to match the time-dependent pump-probe signal observed in the experiment (see section 5.3).

As one can infer from figure 7, the population transfer between the ground state and the bright excited states is clearly not in the perturbative regime, as two Rabi cycles are apparent during the 100 fs UV pulse interaction. The relaxation dynamics, transferring part of the excitation to the dark states, takes place on two time scales: for short delays one can observe a rapid energy transfer, while for longer times the distribution thermalizes. The depletion dynamics of the laser-excited BS primarily takes place due to the coupling to two lower-energy states at »5.5 eV (see figure 7(b)). Closer inspection reveals that these DSs involve the excitation of the 3s and the 3p SAMOs; the respective weights are given in figure 7(a). This relaxation mechanism is the dominant consequence of the electron-vibron coupling. This behavior is expected, as dissipation pathways are generally preferred compared to bath-induced excitations. These thermalization processes are hence less pronounced and occur on a longer time scale.

time [fs]

Figure 7. (a) Population dynamics induced by the pump pulse (sketched in the background). The color coding of ground state (GS), bright states (BS) and dark states (DS) is identical to figure 1 (a). The insets show the weight of the 3s-and 3p-SAMOs in the dominantly populated states. The DSs involving the strongest SAMO excitations is highlighted by the purple color. (b) Dominant population mechanisms for the dynamics in (a).

5.2. Pump-probe dynamics

In order to compute the pump-probe signals from the dynamics of the density matrix pab (t), as discussed above, an extension to the scattering states is required. However, a straightforward implementation of the Lindblad equation including both bound and unbound many-body states is not feasible. This is due to the large dimension of the Hamiltonian after the inevitable discretization of the continuum. Hence, we opt for a perturbation description which allows to compute the pump-probe dynamics from the time-dependent density matrix without incorporating the ionization dynamics explicitly. This is achieved by a method known from time-resolved photoelectron spectroscopy [62]. Adopting a straightforward derivation for the many-body case, one obtains

Nk k Re £ p dt ft dt'Fprobe (t)Fprobe (t') e-i(Ej+w-E+-'k)(t-1'^M^J2paa (t') (10)

for the number Nk of released photoelectrons with momentum k and energy ek. The probe laser pulse is assumed as /probe (t) = Fprobe (t) e-iwt with the pulse envelop Fprobe. The matrix element k,a = ($,k |D|4j)

(D denotes the dipole operator) describes the transition from the intermediate states |4j ) to the final states |$b,k) with one photoelectron and the ion state labeled by ft (energy E+ + ek). Note that only the population of excited intermediate states and not the full density matrix enters equation (10). This is an approximation, which relies on the fact that pathway interferences play only a minor role for the considered two-step ionization process.

For the excited state with one electron in the continuum state |k), |$g,k), we write the usual anti-symmetrized product ansatz

\$b,k) = £¿14+), (1)

where |4+) is an eigenstate of the ionized system with energy E+. In this case the photoemission matrix element reduces to

($b,k |D|4j ) = (k|D|fDb). (12)

Here, fJb (r) stands for the corresponding Dyson orbital. As the sum over all excited states ft is implied, a (computationally expensive) precise calculation of the Dyson orbitals can be omitted by approximating them by simple hole states, i.e., by assuming

|4+= (m, J)) » £ml$a) , E+ » Ea - 6m. (13)

Here, em stands for the KS eigenvalue of orbital fm (r). The matrix element (12) can thus be evaluated in terms of a superposition of matrix elements with respect to the KS orbitals. The scattering states |k) were computed with respect to the spherically averaged KS potential [28], which is known to be an adequate approximation for angle-integrated quantities. Asymptotic corrections ensuring the r-1 behavior are incorporated by smoothly interpolating between the short-range and long-range regimes. Further orthogonalization with respect to the bound KS orbitals is performed.

To reflect the experimental situation, the integration over all photoelectron states has to be performed. Furthermore, we note that equation (10) balances spectral resolution versus temporal resolution in terms of the convolution of the phase factor e—lAEt, AE = Ea + w — E+ — ek, with the envelop function. Due to very long pulses as compared to one oscillation period, this convolution practically yields a Dirac ¿-function with respect to the energy balance. Taking this into account, the total ionization pump-probe signal S (t) = Jdk Nk simplifies to

S (t) dt I dt' Fprobe (t — t) Fprobe (t'— t) Paß (Ea + w) Paa (t'), (14)

where Pap (e) = k JdWk\Mbk,a\2 with k = V2e is proportional to the energy-resolved ionization probability with respect to the initial state \Fa) and final state \F+). As in section 4.2, t denotes the pump-probe delay.

5.3. Theory versus experiment

The pump-probe signal based on the laser-driven and vibron-coupled dynamics of the density matrix can now be compared to the experiment. We remark that equation (14) assumes that the two pulses can be separated and does not account for the scenario of overlapping pulses. We thus limit the comparison between the theoretical calculations and the experiment, presented in figure 5, to the region of exponential decay for t ^ 300 fs. Furthermore, the dynamics presented in figure 7 indicates that the laser intensity exceeds the perturbative regime. However, even with stronger pulses, the bright states around 5.66 eV are the only accessible channels, as absorbing another photon leads to immediate ionization. This will, however, only affect the background signal. Therefore, we adjust the background Sbg = S (t ^ <) to the experiment. After solving the Lindblad master equation for various values of the parameters g and n, we found the best fitfor g= 1.128 and h = 0.77 meV.The latter corresponds to a vibrational lifetime of »5.4 ps, which is in accordance with previous experiments [8]. The small deviation of the scaling factor g from unity underpins the predictive power of our treatment. Note that also g = 1 results in a decay dynamics which closely resembles, apart from the background, the experimental data, whereas varyingg to smaller or larger values clearly deviates from the measurements.

We also calculated the ionization signal for varying photon energy u and compared the resulting spectra in figure 4. To match the laser spectral bandwidth to the experiment, the obtained curves were convolved with a Gaussian having a FWHM of 2.8 nm. The theoretical spectra are centered around the vertical excitation energy, which is lower than what is observed in the experiment. This kind of underestimating bandgaps and thus excitation energies is typical for most DFT calculations and for the C60 molecular, in particular [40]. However, the theoretical and the experimental peak differ only by »10 meV. Generally, both the experimental as well as the theoretical spectra are considerably sharper than in optical absorption measurements. This is a major advantage of the current experimental setup: the dipole excitation dominates over possible loss channels, if the pulse strength is increased, which leads to a narrow resonance.

5.4. Steady-state versus time-dependent picture

The interpretation of the present data on the relaxation dynamics requires the full picture including both, correlated electronic and nuclear degrees of freedom and their interactions [63,64]. It is known that highly excited electronic and vibrational C60 states are strongly mixed [55,65]. In turn, relaxation channels open up depopulating the electronically excited states by internal conversion close to conical intersections similar to electron-phonon coupling in solid state materials. We note that characteristic fragmentation patterns observed in TOF mass spectra using optical lasers as a function of pump-probe delay revealed (nonadiabatic) vibronic coupling on a time scale of tvib = 200-300 fs [45], which is in good agreement with the present findings. Furthermore, nonadiabatic coupling of mixed valence and Rydberg states is ubiquitous in polyatomic molecules affecting potential energy surfaces, energy relaxation and dissociation dynamics [21]. Similar processes have been considered to evoke the loss of small neutral fragments from C60 on a picosecond time scale [66].

The identification of vibronic coupling playing a key role in the electronic energy loss of correlated many-body states may open new vistas for optical control of charge-transport phenomena in smart materials containing these nanospheres. For instance, coherently induced radial symmetric 'breathing' motion of the cage atoms strongly impacts the structure and dynamics of the molecule [67]. The carbon cage and the electron system start to exchange many eV in energy periodically on a sub-100 fs period timescale. The coherent oscillation prevails for several cycles [67], which might be interesting for novel ultrafast switching applications in molecular electronics.

6. Summary and outlook

Ultrashort pulses in the UV spectral range excite correlated many-body states of isolated C60 molecules below the ionization continuum. The population is followed by subsequent UV pulses of the same wavelength that ionize the molecules. By recording the C+0 ion yield as a function of time delay between the pump and probe pulses we observe an exponential decaywith a time constant of tvib = 400 ± 100 fs which is explained by the coupling ofelectronic excitation to nuclear motion in the neutral molecule. The initial electronic relaxation can be as fast as tei = 10+f fs according to the evaluation of the resonance line width in single pulse experiments as a function of excitation wavelength. The experimental results are in good agreement with ab initio calculation of structure and dynamics including electronic correlation and vibronic coupling. However, the UV pulse duration of 100 fs did not allow to observe pure electron dynamics in time-resolved experiments. Future experimental work making use of shorter UV pulses shall reveal the predicted laser-driven Rabi oscillation and time-resolved transformation of the electronic orbitals, i.e., the coupling between different electronic states.

Acknowledgments

This research was supported by the Deutsche Forschungsgemeinschaft through the excellence cluster 'The Hamburg Centre for Ultrafast Imaging (CUI)—Structure, Dynamics and Control of Matter at the Atomic Scale' (DFG-EXC1074), the collaborative research centres 'Light induced dynamics and control of correlated quantum systems' (SFB 925), 'Functionality of Oxide Interfaces' (SFB 762), the GrK 1355 and Grant Number PA 1698/ 1-1.

Appendix A. Electron-vibron matrix elements

In this appendix we provide the details on how the first-order (Herzberg-Teller) electron-vibron coupling matrix elements, equation (5), were computed. After the calculation of the vibrational eigenmodes v, one readily obtains the associated deformations of the fullerene cage. Denoting the collection of nuclear coordinates in equilibrium by {R(0)} = (R(0),..., Rgo), the vibrational distortion is characterized by

(R(n)(q)} = (R(n)(q), ...,R60)(q)) with Rm)(q) = R+ qV(A.1)

where V^) is the displacement eigenvector. Based on equation (A.1), we performed a DFT calculation for each vibrational mode, fixing the magnitude of the distortions at Sq = 0.01. Instead of taking derivatives of the KS potential, we employ the equivalent formulation

k" ^ -1 (fh - to, k =-1V2 + vKS (r, (R(n)(Sq)}), (A. 2)

where we have approximated the derivative by finite differences. The KS Hamiltonian describing the molecule in equilibrium is denoted by h0, while hn describes the distorted molecule. For the evaluation of equation (A.2) we insert a (approximate) completeness relation and obtain

kn ^ 1 [(fffk(n)i(fk(nVj>- ,

hnlfr > = if >, h0if!) = f. (A.3)

Note that the overlaps flf^) incorporate the symmetry properties of the vibronically induced transitions.

The transformation into the many-body basis is accomplished by expressing the one-body coupling operator in second quantization, k = Y^ijkijc} Cj, and evaluating (Flk^lFp) according to the algebra of the fermionic creation and annihilation operators.

Appendix B. Lindblad master equation

For the derivation of the Lindblad master equation in the weak-coupling limit, we follow [58]. For the Hamiltonian (7)-(9) one obtains the following EOM for the electronic density matrix

P(t) = T^abPab (t) lFa) l

dp(t) = — i [Hel (t), p(t)] dt

+ EGaßa'ß' ipßß' (t)IFaXFa' I — 1 (|Fß')<Fß IP(t)} } (B.1)

aßdß' V 2 )

The square (curly) brackets denote the commutator (anti-commutator). The vibronic bath enters into

Gaßdß' = Egn (Eß — Ea) dEß—Ea,Eß'—Ea> K0ißKaß', (B2)

where gv (E) = (NB (E) + 1) Av (E). NB (E) denotes the Bose distribution (displayed in figure 6) which accounts for the occupation ofthe vibronic modes for the given temperature. The vibrational frequencies determine the spectral function Av (E) = 2pd (E — Wv), which we replace by the smeared form

An(E) = V2p/h2 exp [ — (E — Wv)2/2V2] (B.3)

to account for the finite lifetime ofthe vibrations.

References

[1] Sariciftci N S, Smilowitz L, Heeger A J and Wudl F 1992 Science 258 1474-6

[2] Yu G, Gao J, Hummelen J C, Wudl F and Heeger A J1995 Science 270 1789-91

[3] Scholes G D and Rumbles G 2006 Nat. Mater. 5 683-96

[4] Klaiman S, Gromov EVandCederbaumLS2013 J. Phys. Chem.Lett. 4 3319-24

[5] KlaimanS,GromovEVandCederbaumLS2014Phys. Chem. Chem.Phys. 16 13287-93

[6] VooraVK and Jordan KD 2014 NanoLett. 14 4602-6

[7] Pavlyukh Y and Berakdar J 2009 Chem. Phys. Lett. 468 313-8

[8] Hertel IV, Laarmann T and Schulz C P 2005 Adv. At. Mol. Opt. Phys. 50 219-86

[9] Jacquemin R, Kraus S and Eberhardt W 1998 SolidState Commun. 105 449-53

[10] Link S, Scholl A, Jacquemin R and Eberhardt W 2000 Solid State Commun. 113 689-93

[11] LepineF2015 J. Phys. B 48 122002

[12] Feng M, Zhao J and PetekH 2008 Science 320 359-62

[13] Feng M, Zhao J, Huang T, Zhu X and Petek H 2011 Acc. Chem. Res. 44 360-8

[14] Pavlyukh Y and Berakdar J 2011 J. Chem.Phys. 135 201103

[15] Reisler H and Krylov A12009 Int. Rev. Phys. Chem. 28 267-308

[16] Zhao J and Petek H 2014 Phys. Rev. B 90 075412

[17] Watzel J, Pavlyukh Y, Schäffer A and Berakdar J 2016 Carbon 99 439-43

[18] ZoppiL,Martin-SamosLandBaldridgeKK2015Phys. Chem. Chem.Phys. 176114-21

[19] Dougherty D B, Feng M, Petek H, Yates J T and Zhao J 2012 Phys. Rev. Lett. 109 266802

[20] Johansson J O, Henderson G G, Remacle F and Campbell E E B 2012 Phys. Rev. Lett. 108 173401

[21] Mignolet B, Johansson JO, Campbell EE B and Remacle F 2013 Chem. Phys. Chem 14 3332-40

[22] Bohl E, SokölKP, Mignolet B, Thompson J O F, Johansson J O, Remacle F and Campbell EE B 2015 J. Phys. Chem. A119 11504-8

[23] Boyle M, Hoffmann K, Schulz C P, Hertel IV, Levine R D and Campbell E E B 2001 Phys. Rev. Lett. 87 273401

[24] Boyle M, Laarmann T, Hoffmann K, Heden M, Campbell E E, Schulz C P and Hertel IV 2005 Eur. Phys. J. D 36 339-51

[25] RoyX etal2013 Science 341 157-60

[26] Shchatsinin I, Laarmann T, Zhavoronkov N, Schulz C P and Hertel IV 2008 J. Chem. Phys. 129 204308

[27] LiHetal2015Phys. Rev. Lett. 114 123004

[28] Schüler M, Berakdar J and Pavlyukh Y 2015 Phys. Rev. A 92 021403

[29] WoppererP etal2015 Phys. Rev. A 91 042514

[30] Schüler M, Pavlyukh Y, Bolognesi P, Avaldi L and Berakdar J 2016 Sci. Rep. 6 24396

[31] Garcia de AbajoF J 2004 Phys. Rev. B 70 115422

[32] Marques MAL, Maitra N T, Nogueira F M S, Gross E K U and Rubio A 2012 Fundamentals of Time-Dependent Density Functional Theory (Berlin: Springer)

[33] Bartlett Rand MusiaA M 2007 Rev. Mod. Phys. 79 291-352

[34] Aryasetiawan F and Gunnarsson O 1998 Rep. Prog. Phys. 61237

[35] Chong D P 1995 Recent Advances in Density Functional Methods (Singapore: World Scientific)

[36] AndradeXetal2015Phys. Chem. Chem.Phys. 1731371-96

[37] van Leeuwen R and Baerends E J1994 Phys. Rev. A 49 2421-31

[38] Schipper P R T, Gritsenko O V, Gisbergen S J A V and Baerends E J 2000 J. Chem. Phys. 112 1344-52

[39] Smith A L 1996 J. Phys. B 29 4975

[40] Bauernschmitt R, Ahlrichs R, Hennrich F H and Kappes M M 1998 J. Am. Chem. Soc. 120 5052-9

[41] Johansson J O, Henderson G G and Campbell E E B 2013 EPJ Web Conf. 4102015

[42] Schmidt B, Hacker M, Stobrawa G and Feurer T LAB2-A Virtual Femtosecond Laser Lab (http://lab2.de/)

[43] Wiley W C and McLaren IH 1955 Rev. Sci. Instrum. 26 1150-7

[44] Knoesel E, Hotzel A and Wolf M 1998 Phys. Rev. B 5712812-24

[45] Shchatsinin I, Laarmann T, StibenzG, Steinmeyer G, StalmashonakA, ZhavoronkovN, Schulz C P and Hertel IV 2006 J. Chem. Phys. 125194320

[46] Campbell E E B, Hansen K, Hoffmann K, Korn G, Tchaplyguine M, Wittmann M and Hertel IV 2000 Phys. Rev. Lett. 84 2128-31

[47] Hansen K, Hoffmann K and Campbell E E B 2003 J.Chem.Phys. 119 2513-22

[48] Chancey C C and O'Brien M C M 1997 The Jahn-Teller Effect in C60 and Other Icosahedral Complexes (Princeton, NJ: Princeton University Press)

[49] Leszczynski J 2012 Handbook of Computational Chemistry (Berlin: Springer)

[50] Pittalis S, Delgado A, Robin J, Freimuth L, Christoffers J, Lienau C and Rozzi C A 2015 Adv. Funct. Mater. 25 2047-53

[51] Ben-Nun M, Quenneville J and Martinez T J 2000 J. Phys. Chem. A104 5161-75

[52] Nyman G and Yu H G 2000 Rep. Prog. Phys. 63 1001-59

[53] Schlegel H B 2003 J. Comput. Chem. 24 1514-27

[54] Fischer M, Handt J, Seifert G and Schmidt R 2013 Phys. Rev. A 88 061403

[55] Zhang G P, Sun X and George T F 2003 Phys.Rev. B 68 165410

[56] Tamura H, Martinazzo R, Ruckenbauer M and Burghardt 12012 J. Chem. Phys. 137 22A540

[57] Köppel H, Yarkony D Rand Barentzen H 2009 The Jahn-Teller Effect: Fundamentals and Implications for Physics and Chemistry (Berlin: Springer)

[58] Breuer H P and Petruccione F 2002 The Theory of Open Quantum Systems (Oxford: Oxford University Press)

[59] Sassara A, Zerza G, Chergui M, Negri F and Orlandi G 1997 J. Chem. Phys. 107 8731-41

[60] Menéndez J and Page J B 2000 Vibrational spectroscopy of C60 Light Scatteringin Solids VIII (Berlin: Springer) 27-95

[61] Orlandi G and Negri F 2002 Photochem. Photobiol. Sci. 1289-308

[62] Kemper A F, SentefM A, Moritz B, Freericks J K and Devereaux T P 2014 Phys. Rev. B 90 075126

[63] Lezius M, Blanchet V, Rayner D M, Villeneuve D M, Stolow A and Ivanov M Y 2001 Phys. Rev. Lett. 86 51-4

[64] Lezius M, Blanchet V, Ivanov M Y and Stolow A 2002 J. Chem. Phys. 117 1575-88

[65] Torralva B, Niehaus T A, Elstner M, Suhai S, Frauenheim T and Allen R E 2001 Phys. Rev. B 64153105

[66] Boyle M, Laarmann T, Shchatsinin I, Schulz C P and Hertel IV 2005 J. Chem. Phys. 122 181103

[67] Laarmann T, Shchatsinin I, Stalmashonak A, Boyle M, Zhavoronkov N, Handt J, Schmidt R, Schulz C P and Hertel IV 2007 Phys. Rev. Lett. 98 058302