Monitoring vibronic coherences and molecular aromaticity in photoexcited cyclooctatetraene with an X-ray probe: a simulation study

Abstract

Understanding conical intersection (CI) dynamics and subsequent conformational changes is key for exploring and controlling photo-reactions in aromatic molecules. Monitoring of their time-resolved dynamics remains a formidable experimental challenge. In this study, we simulate the photoinduced S₃ to S₁ non-adiabatic dynamics of cyclooctatetraene (COT), involving multiple CIs with relaxation times in good agreement with experiment. We further investigate the possibility to directly probe the CI passages in COT by off-resonant X-ray Raman spectroscopy (TRUECARS) and time-resolved X-ray diffraction (TRXD). We find that these signals sensitively monitor key chemical features during the ultrafast dynamics. First, we distinguish two CIs by using TRUECARS signals with their appearances at different Raman shifts. Second, we demonstrate that TRXD, where X-ray photons scatter off electron densities, can resolve ultrafast changes in the aromaticity of COT. It can further distinguish between planar and non-planar geometries explored during the dynamics, as e.g. two different tetraradical-type CIs. The knowledge gained from these measurements can give unique insight into fundamental chemical properties that dynamically change during non-adiabatic passages.

---

1 Introduction

Aromaticity is a property of cyclic (ring-shaped), typically planar (flat) molecular structures with delocalized π electrons that gives increased stability compared to saturated (non-aromatic) compounds having single bonds or other non-cyclic arrangements with the same set of atoms. Aromaticity plays key roles in chemical reactions (electrophilic aromatic substitution¹), molecular physics (organic semiconductors² and aromatic ring currents³), and biochemistry, where amino acids serve as building blocks of proteins. Thus, monitoring conical intersection (CI) dynamics and conformational changes is key for unravelling and controlling photochemical reactions in aromatic molecules.

Cyclooctatetraene (COT) is a conjugated cyclic 4n π-electron system that has a D_8h planar conjugated aromatic π-network in the lowest excited state but is non-aromatic in higher excited states (S_n>1) as well as in the ground state, where it has a non-planar boat-like D_2d structure with localized single and double C–C bonds. Thus, COT may serve as a prototypical photoactive unit where photon absorption is employed as a control knob to switch between non-aromatic and aromatic states. Therefore, photorelaxation through CIs induces strong modification of the planarity and electron density, and thereby aromaticity, but without ring opening like heterocyclic compounds.^4,5

Its thermal and photochemical relaxation pathways have drawn significant experimental and theoretical attention^6–8 (Scheme 1). A photon initially excites COT to the optically allowed (bright) S_2/3 state (1), followed by ultrafast non-radiative decay to the optically forbidden (dark) S₁ state (2). A non-adiabatic transition to S₀ is controlled by two tetraradical-type conical intersections (CI),^8,9 CI_st (3) and CI_b (4). CI_st has a typical out-of-plane triangular –(CH)₃– kink of a triradical nature similar to other unsaturated hydrocarbons,^10,11 leading to three- or four-membered ring formation (5) and cis → trans isomerization (6).⁸ However, the decay via this channel is suppressed by an energy barrier. CI_b holds C_2v symmetry and has two unpaired electrons centered at single carbon atoms and two resonance-stabilized allyl radicals.⁹ Thanks to its lower barrier, it acts as an elective radiation-less channel leading to the formation of semibullvalene (SBV, 7) as a main product,¹² and a double-bond shifted (DBS, 8) or the original COT (1) is formed as byproducts.¹³ We refer to ref. ⁸ for the comprehensive photo/thermal reaction pathways.

Scheme 1 Photochemical and thermal reaction scheme for cyclooctatetraene.

The DBS yields a degenerate product of the parent COT, where the π-electrons migrate within the octagonal perimeter, resulting in a site exchange between singly- and doubly-bonded carbon atoms eventually leading to photoinduced bond order inversion within the ring. Modulation of these pathways in favor of desired channels could be accomplished by chemical modification with π-donor and -acceptor substituents, or by quantum control strategies, making COT a potent and attractive design unit for e.g. molecular photoswitches that allow control of aromaticity and bond order.^14,15

While the S₁ to S₀ relaxation pathway has been extensively studied, the nonradiative S₃ to S₁ decay has remained largely unexplored. It is a key issue prior to unravelling the role of S₁/S₀ CI, since the dark S₁ state is only accessible by relaxation from the higher excited states. The S₃ → S₁ relaxation can further affect the nuclear dynamics on S₁, thus influencing the resulting photoproducts. CIs play important roles in many photophysical and photochemical processes.^16,17 Monitoring CI pathways is important for controlling and achieving the desired photoproducts. Various studies have shown that modifying the initial condition of nuclear wavepackets in the electronic states forming the CIs can have tremendous effects on the photoproducts,^18,19 corroborating the necessity of including the S₃ → S₁ relaxation mechanism in photoexcited COT.

Thanks to the unique temporal, spectral and spatial resolutions provided by free-electron X-ray light sources, many X-ray techniques have been proposed to monitor the CIs in molecules.^20,21 A popular technique is femto/attosecond spectroscopy, which indirectly probes CIs by depletion/appearance/bifurcation of the absorptive lines.^22–24 Polli et al. implemented ultrafast optical spectroscopy to probe light-induced photoisomerization of rhodopsin and mapped out the energy gap between the ground and excited electronic states as a function of time.²⁵ Such evidence is indirect and circumstantial and does not give the direct signature of CI.

To this aim, transient redistribution of ultrafast electronic coherences in attosecond Raman signals (TRUECARS) has been theoretically proposed.²⁶ This technique directly monitors vibronic coherences created during the CI passage with no contributions of populations (background-free), which is more direct evidence or which signifies the presence of CI better than the energy gap between the involved states. A hybrid broadband/narrowband pulse used in TRUECARS can offer a good combination of both spectral and temporal resolutions. This stems from the ultrafast timing of CIs as well as the few tens of eV energy range spanned by the vibronic coherences. It has been theoretically used to monitor the CI passage in photo-relaxation in (4-thio)uracil^27–29 and energy transfer in a heterodimer³⁰ and a triarylamine trimer.³¹ A major difficulty in the implementation of TRUECARS is the precise phase control between two pulses, which is under development. Herein, we focus on how a TRUECARS signal distinguishes two different CI passages during photorelaxation in COT.

On the other hand, time-evolving electronic charge densities ain CI passages can be imaged with subfemtosecond resolution using ultrafast time-resolved X-ray diffraction/scattering (TRXD).³² In ultrafast TRXD experiments, a molecule is prepared in a time-evolving superposition of states by using an optical laser, undergoing non-stationary dynamics, and then a hard X-ray probe pulse is scattered by the excited molecule onto a detector, yielding the three components of the scattering signal: elastic scattering, inelastic scattering, and mixed (in)elastic scattering related to electronic coherence which contribute to the signal.^32,33 The snapshots at different pump-probe delays create a movie presenting temporal evolution of electron densities triggered by the pump pulse. Since the pioneering theoretical work of Wilson et al.,³² there has been an immense development of theory^34–37 to study non-stationary molecular samples in excited electronic states. A novel development of bright XFELs extended the X-ray scattering measurement in a solid crystal into gas or liquid phases,^38–40 which involves real-time monitoring of the coherent vibrational motion of excited N-methylmorpholine⁴¹ and its orientation of the transition dipole moment using gas-phase X-ray scattering.⁴² In the liquid phase, ultrafast hydrogen bond dynamics,⁴³ solvent reorganization coupled to intramolecular charge transfer,⁴⁴ liquid–liquid phase transition,⁴⁵ and structural changes of proteins⁴⁶ have been monitored. Nielsen et al. recorded coherent nuclear dynamics with atomistic resolution on the excited⁴⁷ and ground state⁴⁸ potential energy surfaces for systems in an environment. XFELs are being upgraded to achieve brighter light sources, higher repetition rate, and greater spatiotemporal resolution to expand their applications.

Our previous theoretical work has demonstrated that the TRXD signals can potentially image transient electron transition densities directly associated with CI passages in azobenzene^49,50 and 4-thiouracil,²⁹ exhibiting characteristic positive/negative oscillations due to the formation of electronic coherences. We had further shown that the two-dimensional diffraction pattern, dominated by elastic scattering, can be used to monitor conformational changes e.g., cis to trans photoisomerization in azobenzene.⁴⁹ Such valence electron densities can be used to monitor the electron density (aromaticity) and planarity variation of the molecule.

In the present study, we implement a non-adiabatic excited state molecular dynamics⁵¹ protocol to track the S₃ → S₁ relaxation pathway in optically excited COT. We employ a semiempirical ab initio multiple cloning (AIMC) approach based on Multi-Configurational Ehrenfest (MCE), which provides an accurate description of non-adiabatic molecular dynamics in large conjugated molecules⁵² with affordable cost.

We capture the appearance of vibronic coherences as well as aromaticity changes upon photorelaxation by using an ultrafast X-ray probe. The vibronic coherences generated at the multiple CIs are tracked by the TRUECARS signal. We show highly diverse scenarios for excited state relaxation and record the temporal and the energetic profiles of CIs by using the TRUECARS signal and its spectrogram. We find that explicit use of transition polarizabilities is crucial to assess accurate observation of vibronic coherences evolving during the CI passage. Due to the quantum nature of nuclear motions, the vibronic coherences do not vanish after passing through the CI passage. The TRUECARS signal provides a clear signature of the two CIs (S₃/S₂vs. S₂/S₁) by their different timings and energy splitting distributions between the involved states.

We find that TRXD is a powerful tool for the real-time tracking of the aromaticity and molecular conformation changes in molecules by tracking the evolving valence electron densities. The 2D elastic scattering pattern can differentiate the different CI pathways, CI_st and CI_b, and the photoproduct, SBV from the reactant COT. This helps map the comprehensive relaxation pathways of COT from bright S₃ to S₀via the dark and aromatic S₁ state.

2 Results and discussion

The AIMC approach was employed to describe the S₃ → S₁ electronic transitions in COT. AIMC naturally includes decoherence through cloning events when mean-field theory fails to describe two electronic states evolving on very different surfaces. This approach has been successfully applied to describe photoinduced dynamics in large molecules, such as a dendrimer³¹ and a bichromic molecule.³⁰

The vibronic coherences emerging at the excited state CIs are tracked by the TRUECARS signal and the geometry and aromaticity changes are monitored by the TRXD signal. TRUECARS uses a hybrid field (2 fs), with a central frequency of 200 eV and (500 as) following the pump-probe waiting time T (Fig. 1). The central frequency was chosen to maximize the signal strength by maximizing the polarizability cross section while still staying off-resonant (pre-resonant).⁵³ Otherwise the populations will contribute and dominate the coherences.^54,55 The signal is finite only when there is an overlap of nuclear wavepackets in different electronic states, making this technique free from population background. A single broadband pulse with the same central frequency and bandwidth (Fig. S1†) is used for TRXD.

Fig. 1 (a) Pulse configuration and (b) loop diagram for TRUECARS with cyclooctatetraene (COT) placed on the xy plane (Lewis structure given). The grey area indicates the electronic and nuclear population in the bright S₃ state created by the pump pulse (not considered explicitly in the simulation) and a free evolution period of the molecule. At time delay T, the hybrid (broad) and (narrow) pulse is applied to probe the dynamics. See the ESI† for loop diagram rules.

We display the TRUECARS/TRXD signal as well as evolution of population, coherence, and molecular geometry, averaged over the total 98/57 trajectories with equal weight, in Fig. 2. Note that only the signal averaged over all trajectories is experimentally observable although individual trajectories illustrate different molecular dynamics scenarios. We first note that our AIMC approach well describes the S₃ → S₁ photo-relaxation dynamics of COT. The S₁ population dynamics (Fig. 2e) is fitted with P(t) = Ae^{−(t−t₀)/k}, with t₀ = 26.5 fs, yielding a growth time of 54.5 fs. This is in a good agreement with our previous surface hopping simulations⁵⁶ and time-resolved mass spectroscopy measurement of photoexcited cyclooctatriene and bicyclooctadiene using a near-IR photoionization probe,⁵⁷ where ¹B₂ to ²A₁ (C_2v symmetry, corresponding to S₃ to S₁ in the current study) relaxation was estimated to occur within 67 fs. The average population is distributed with a large fraction (70%) in S₁ and a smaller one in S₂ (20%) and the others for S₃ to S₄ states at 250 fs.

Fig. 2 TRUECARS and TRXD signals and relevant molecular properties in the ensemble of 98/57 trajectories of COT. The averaged frequency-dispersed TRUECARS signal S(ω_r, T), for (a) total, (b) S₂/S₁ coherence, and (c) S₃/S₂ coherence. (d) Combined potential energy surfaces of the electronic states in all trajectories. (e) Averaged population in the electronic states. (f) Averaged coherence magnitude ρ_KL between electronic states, according to eqn (3). (g) Averaged expectation value of the polarizability operator calculated with eqn (5). (h) Averaged nonadiabatic coupling magnitude for each coherence ρ_KL. (i) Averaged FROG spectrogram, according to eqn (7), which is extracted from the TRUECARS signal by integration over the negative Raman shifts (ω_r < 0). (j) Time evolution of the averaged molecular geometry, bond length (top), dihedral angle (C₁–C₄ in Fig. S6a†), and bond alternation (bottom). (k) Averaged two-dimensional TRXD scattering pattern projected on the xy (left), xz (middle), and yz (right) planes at T = 1 fs. The top, middle, and bottom panel shows contributions from elastic scattering, inelastic scattering, and mixed elastic/inelastic scattering (coherence) to the total signal, respectively. (l) Same as (k) but at T = 250 fs.

The TRUECARS signal (Eqn (4) and (5)) is visible over the entire simulation time, with stronger magnitudes at e.g. 0 to 70 fs. The two CI passages are well captured by the TRUECARS signal (S₂/S₁ in Fig. 2b and S₃/S₂ in Fig. 2c). The molecule enters the S₃/S₂ CI region, with strong nonadiabatic coupling (NAC) (Fig. 2h) due to their close spacing in energy (Fig. 2d), creating a vibronic coherence (Fig. 2f), and thereby, the TRUECARS signal shows up from the beginning and maintains its amplitude until 70 fs (Fig. 2c). A delayed S₂/S₁ CI is observed (Fig. 2b) as the major population transfer occurs between 50 and 100 fs, but the relevant TRUECARS signal remains strong until 170 fs. We note that after the second CI, S₂ and S₁ evolve differently, and their energy splitting increases, as the TRUECARS signal is observed at a higher Raman shift ω_r (Fig. 2b). We find that the S₄ state is only slightly affected in all trajectories.

In our previous studies,^30,31 we used a constant polarizability over the nuclear space and all coherences contributed to the TRUECARS signal according to their magnitude with no further selectivity. This approximation holds when only two electronic states are involved, but as soon as more than one electronic transition is involved, the relative transition polarizability strengths determine the magnitude of the individual contributions to the total signal. This means that while the coherences themselves could be equally strong, the coherence associated with the higher transition polarizability will dominate the TRUECARS signal. We display the TRUECARS signal calculated with geometry-independent polarizability over the nuclear space in Fig. S2.† We find that the signal is particularly strong at 0 to 25 fs for S₃/S₂ and 90 to 100 fs for S₂/S₁. This signal looks more sensitive to the timing of the CIs, but the relative transition polarizability strengths determine the observed signal shown in Fig. 2a. This is more evident in trajectory 3 (Fig. S3†). The total coherence magnitude is maintained between 0 and 60 fs (Fig. S3d†), but the TRUECARS signal is strongest between 40 fs and 60 fs (Fig. S3a†) due to the large transition polarizability (Fig. S3c†). If geometry independent polarizabilities were used, the TRUECARS signal would be equally strong between 0 fs and 60 fs as shown in Fig. S3b.† This implies that, for systems undergoing multiple CIs, both the relative strength of the polarizabilities and the topologies of potential energy surfaces influence the signal and need to be properly accounted for in the simulations.

The integrated TRUECARS spectrogram reveals the energy splitting distribution between electronic states involved in the coherence.²⁸ In turn, a transient energy splitting is encoded in the temporal gain/loss oscillations in the TRUECARS signal at a given Raman shift (ω_r). We display the integrated frequency resolved optical-gating (FROG)⁵⁸ spectrogram, given by eqn (7), in Fig. 2i. The spectrogram between 0 and 25 fs is distributed between 0.3 and 0.7 eV, representing the energy splitting between S₂ and S_3/4 during the first CI passage. The S₂/S₁ CI passage is captured from 30 to 50 fs. It is more evident in Fig. S4a,† where the FROG spectrogram reveals the energy splitting of CI at Raman shift ω_r = 0.02 eV. An increasing energy splitting between S₁ and a higher excited state can be observed where the spectrogram evolves from 1 to 2 eV after 50 fs (Fig. S4b and c†). The frequency profile of the spectrogram maps the energy gap between the relevant states in Fig. 2d.

The time-evolving geometric features are displayed in Fig. 2j. These include the bond length (C₁–C₃, top panel), dihedral angle (C₁ to C₄, middle panel), and bond length alternation (0.5(b₁₅ + b₄₈) − b₁₈, bottom panel) over time (the atomic labels are given in Fig. S6†). Earlier, we have used the same parameters to monitor the non-adiabatic passage of COT using semiempirical trajectory surface hopping dynamics.⁵⁶ Based on the Franck-Condon approximation, the molecular geometry starts from the non-planar S₀ minimum conformation with different C₁–C₅versus C₄–C₈ bond lengths, and then approaches those at planar S₁ minimum geometry with equalized bond lengths. Indeed, the dihedral angle approaches 180°, and atomic distance increases to 3.6 Å as the population is transferred to the S₁ state.

Previously, we had demonstrated that the coherence contribution to the TRXD signal can image the evolving electron densities during the CI passage, which is characterized by its phase oscillation between gain and loss along the temporal axis.^29,49 The temporal oscillation showed the strongest intensities during the CI passage and the observed phase change corresponds to real-space phase changes of electron density as the molecule crosses the CI. The coherence contribution involves mixed elastic and inelastic scattering events, where the latter involves only a single active electron transition. Hence it is weak and buried under the stronger state densities, where all electrons contribute to the signal. We had suggested extraction of the information by observation at higher momentum transfer;⁴⁹ however, it only works for systems where the transition density is more localized in real space (high q) compared to the delocalized state densities (low q). The valence excited states in COT, where an electron in an occupied π orbital is promoted to an unoccupied π* orbital, exhibit delocalized electron densities (Fig. S5†), and thus, the hard X-rays are not beneficial in COT. A frequency-resolved diffraction set-up,⁵⁹ can be alternatively used since the coherences oscillate faster than the populations. In those studies, we had assumed very short waves (more than 20 keV), which are currently unavailable, and convolute the temporal resolution, rendering the experiment difficult. Currently, a new superconducting accelerator at the European X-ray Free Electron Laser in Hamburg and the Stanford Linear Coherent Light Source is being developed to provide up to 25 keV of photon energy, which may enable such frequency-resolved diffraction measurement. A development of large free electron laser facilities will enable direct observation of the evolving coherence electron densities in the future.

The state densities are virtually identical for the ground state and excited state, and sensitive to the molecular conformation and electron densities. Hence, we expect that TRXD can be used to monitor the change in the molecular conformation (non-planar to planar) and aromaticity (localized electron densities at the double bond to delocalized densities). We display a two-dimensional (2D) TRXD signal in Fig. 2k and l. Note that the 2D patterns shown in Fig. 2–4 are imaged only with the valence electron densities since the semiempirical AIMC-NEXMD calculations use basis functions composed of only valence electrons. For comparison, we display the 2D pattern of the TRXD signal for S₀ and S₁ optimal geometries in Fig. S7,† calculated with CASSCF(8e/8o), involving all π and π* orbitals in the 6-31G* basis set. The elastic scattering pattern projected on the xy plane shows the localized double bond features in the S₀ minimum. In contrast, the pattern exhibits well delocalized electron density over the entire ring in the S₁ minimum conformation, as the molecular geometry becomes planar and all valence bonds are equalized. On comparison, the 2D pattern projected on the xz (Elas_xz, Fig. S7a/b† top middle panel) and yz pattern (Elas_yz, Fig. S7a/b† top right panel) is less sensitive but we observe an elongated pattern in the S₁ minimum conformation compared to the non-planar conformation. We observe that the signal shows similarity with that of the S₀ minimum in the beginning (at 1 fs, Fig. 2k) but ends with a pattern (at 250 fs, Fig. 2l) similar to that of the S₁ minimum.

The inelastic scattering contribution from electronic populations (middle panels) or mixed elastic/inelastic scattering contribution from electronic coherence (bottom panels) does not exactly match with those in Fig. S7,† because they are calculated without considering the population (wave function coefficients) of the adiabatic states. It is not straightforward to directly compare them. Nevertheless, we note that the coherence contribution exists from 1 fs as the molecule enters the S₃/S₂ CI immediately. Their phase oscillations along the temporal axis could also be used to directly monitor the CI passage.

On a side note, we show that the TRXD signal can also be used to track the following S₁ → S₀ dynamics, owing to the different conformations of two CIs and photoproducts. We adopted their optimized geometry from the previous study⁸ and display their static 2D XRD pattern, projected on the xy plane in Fig. S6.† We find that the CI_b (Fig. S6c†) and its main product SBV (Fig. S6e†) exhibit sufficiently different XRD patterns compared to the CI_st (Fig. S6d†) and original COT (Fig. S6a†). The presented and discussed 2D TRXD signal may only be accessible by simulation, if the alignment of COT is not achievable.

Next, we explore two individual trajectories illustrating very different molecular dynamics scenarios. This is only possible in simulations, since only the ensemble averaged signal is observed in the experiments. However, exploring individual trajectories do help in understanding the entire molecular physics. We first describe trajectory 1, a representative scenario, where most (83%) of the populations ends with S₁ within 100 fs without cloning events. A molecule enters the S₃/S₂ CI region with an immediate population transfer from S₃ to S₂ (Fig. 3e), and thus the TRUECARS signal shows up from the beginning (Fig. 3a and c). S₃ and S₂ then evolve differently (S₂ approaches the second CI), their energy splitting increases, and the TRUECARS signal (Fig. 3c) shows faster oscillation at 20 to 30 fs than 10 to 20 fs. The second CI is reached (by nuclear wavepacket) at 30 fs, when the major population transfer between S₂ and S₁ is facilitated (Fig. 3e) by their large NAC (Fig. 3h). Then, the TRUECARS signal frequencies are shifted to larger Raman shifts (Fig. 3b) as the electronic energy gap between the involved states, S₁vs. S_2/3/4, increases, while S₁ is stabilized. Note that there is no finite population transfer between S₂ and S₁ (after CI), but their coherence (Fig. 3f) is maintained. The signal gets even stronger as the expectation value of the transition polarizability (Fig. 2f) remains strong in this region.

Fig. 3 The TRUECARS and TRXD signals and relevant molecular properties in trajectory 1. The frequency-dispersed TRUECARS signal S(ω_r, T), for (a) total, (b) S₂/S₁ coherence, and (c) S₃/S₂ coherence. (d) Potential energy surfaces of the excited states. (e) Population in the electronic states. (f) Coherence magnitude ρ_KL. (g) The expectation value of the polarizability operator. (h) Nonadiabatic coupling magnitude for each coherence ρ_KL. (i) Integrated FROG spectrogram. (j) Time evolution of the molecular geometry, bond length (top), dihedral angle (middle), and bond alternation (bottom). (k) Two-dimensional TRXD scattering pattern at T = 1 fs. (l) Same as (k) but at T = 83 fs.

The time-evolving geometric features and the 2D TRXD signal displayed in Fig. 2j–l, show that the molecule ends with planar and aromatic geometry at 83 fs, as most of the population ends with S₁. Note that the simulation is terminated at 83 fs due to the fact that the trajectory has reached the region where the energy gap between S₁ and S₀ is smaller than 0.1 eV (assumed to be S₀/S₁ CI, see the Method section below).

Next, we illustrate an extreme opposite scenario, where several cloning events happen and the non-adiabatic simulation is terminated before reaching the S₁ minimum (Fig. 4). In this scenario, the S₃ state, which decays immediately in trajectory 1, survives until 110 fs (Fig. 4e). The major population transfer between S₂ and S₃ states occurs at 10 fs, slower than the immediate occurrence in trajectory 1. The corresponding TRUECARS signal becomes strongest between 20 and 40 fs when the S₂ state is significantly populated and the coherence ρ₂₃ (and thereby the expectation value of transition polarizability) is large (Fig. 4f). After 20 fs, the S₃ and S₂ states evolve differently, ρ₂₃ decreases, and their energy gap increases. Eventually, at 30 fs, mean-field theory breaks down and a cloning event occurs to describe the different evolutions of the S₂ and S₃ states separately.

Fig. 4 Same as Fig. 3 but for trajectory 2. The yellow vertical lines mark the cloning events. The 2D TRXD pattern shown in (l) is plotted at T = 115 fs.

Another cloning event happens at 60 fs, when the major population transfer between the S₂ and S₁ ends (Fig. 4e). Due to the finite energy splitting between them (Fig. 4d) and small NAC (Fig. 4h), the S₁ state is not significantly populated. Hence, the TRUECARS signal augmented by the S₂/S₁ coherence (Fig. 4b) is an order of magnitude weaker than that from the S₃/S₂ coherence (Fig. 4c). We find the latter dominates the total signal as the 〈Ψ₂(t)|α_S2S3|Ψ₃(t)〉 transition polarizability is strongest (Fig. 4g). The integrated FROG spectrogram shown in Fig. 4i well describes the increasing energy splitting between S₂ and S₃ between 20 and 40 fs, and the decreasing pattern thereafter.

We find that the relevant geometric features do not converge to those of the S₁ minimum conformation at the end of the dynamics simulation (115 fs) because the majority of the population stays in the S₂ state, which has a boat-like non-planar conformation, similar to the S₀ minimum. The dihedral angle is larger than 20° and the bond length is maintained shorter than that of the delocalized ring. Thus, the 2D elastic scattering pattern projected on the xy plane at 115 fs shown in Fig. 4l, maintains the localized double bond feature, when compared to Fig. 4k at 1 fs.

Examination of all 98 (for TRUECARS) or 57 (for TRXD) trajectories shows highly diverse scenarios. The two trajectories discussed above are exemplary cases that contribute to the total ensemble. The consequences of the other non-adiabatic dynamics trajectories are placed in between those of the above-mentioned typical cases. Overall, we demonstrated that AIMC dynamics successfully describes the ultrafast S₃ to S₁ relaxation of excited states, despite the diverse evolution of individual trajectories. Note again, that only an ensemble averaged signal can be observed in the experiments.

The implementation of the proposed TRUECARS and TRXD experiments requires precise phase control between narrow and broad pulses, and the alignment of COT molecules perpendicular to the propagation of the X-ray probe pulse, respectively. Extracting coherence information in the TRXD signal requires a very hard X-ray beam, which is under development. Once achieved, the timing of the CI passages and the energetic nature of vibronic coherences can be captured by using the TRUECARS spectrogram at different Raman shifts. The evolving electron densities and molecular conformation could be tracked by using TRXD signals. The signal is off-resonant with any molecular transition; it does not require any specific core property, and directly reveals the coherences between valence electronic states.

3 Conclusions

We have carried out an AIMC molecular dynamics simulation to study non-adiabatic dynamics of photoexcited cyclooctatetraene, tracing its internal conversion from the bright S₃ to the dark S₁ state. The resulting excited state lifetime shows an excellent agreement with experiment, demonstrating an adequate characterization of the molecular photophysics by the semi-empirical multi-configurational Ehrenfest approach. The vibronic coherences created at the two major CI passages persist across the ensemble averaging over 98 trajectories and are well captured by the TRUECARS signal at different Raman shifts. We demonstrated that TRUECARS and TRXD in combination can distinguish between different CIs explored during the photoinduced dynamics. Changes in the molecular aromaticity, as well as non-planar to planar geometrical dynamics, are directly resolved in the signals. The signals combined with a semiempirical nonadiabatic molecular dynamics protocol thereby provide accurate temporal, structural and energetic profiles of the CI pathway that could reveal novel chemical design opportunities and control knobs for photochemical reactions.

4 Methods

The excited state non-adiabatic dynamics of COT is calculated using the ab initio multiple cloning (AIMC)^60,61 approach implemented in the non-adiabatic excited state molecular dynamics (NEXMD) package.⁶² This is an extension of the Multiconfigurational Ehrenfest (MCE)⁶³ method, which follows the spirit of the ab initio multiple spawning (AIMS) approach,⁶⁴ allowing bifurcations of the molecular wave function in the nuclear configuration space thus naturally accounting for decoherences. Details about the connection between these two similar approaches can be found elsewhere.⁶⁵ For AIMC, ensembles of individual Ehrenfest trajectories (clones) are used as basis functions to represent the quantum wave function of electrons and nuclei |Ψ(t)〉:⁶⁶

Each configuration ψ_n(t) is factorized into a nuclear part χ_n(t) and orthonormalized adiabatic multi-configuration electronic eigenfunctions ϕ⁽ⁿ⁾_I:

The nuclear wavepacket χ_n(t) in each configuration is given by a Gaussian function centered at the Ehrenfest trajectory.

Population transfer between electronic states can occur during CI passages, where the Born Oppenheimer approximation breaks down and the motions of electrons and nuclei are strongly coupled. The original nuclear wave packet branches into multiple parts, where the excess energy follows different relaxation pathways, each dominated by a single adiabatic state. In such cases, Ehrenfest trajectories evolving on an average potential energy surface can lead to unphysical nuclear motions. AIMC recognizes these cases and replaces the original Ehrenfest trajectory configuration with two new configurations and coefficients, each evolving along its own distinct mean-field. This splitting is denoted as a cloning event, which allows to naturally account for decoherence of vibronic wavepackets evolving on the sufficiently different potential energy surfaces. More details of the AIMC method and its implementation can be found in ref. ⁶⁷ and ⁶⁸.

The coherences between electronic states are given by:

The phases of both the electronic and the nuclear parts of the molecular wave function are accounted for when calculating the vibronic coherence magnitude ρ_KL.

The TRUECARS signal is finally given by:²⁶

where “Im” denotes the imaginary part, is a hybrid narrow (2 femtosecond)/broadband (500 attosecond) Gaussian pulse envelope (Fig. 1a), ω_s is the central probe frequency, and T is the time delay between the pump and the probe. The expectation value of the transition polarizability 〈Ψ(t)|α_KL|Ψ(t)〉 is given by

The transition polarizability α_KL is calculated from the transition charge density, σ_KL, where

using the state charge density matrices P_rs^ij, and the basis set of atomic orbitals φ_r(r). Populations do not contribute to the signal, since α_KK is zero along the diagonal, and only the transition polarizabilities (off-diagonal elements) between electronic states are finite. The TRUECARS signal is calculated for a randomly oriented ensemble by averaging over the x, y, and z axes. We shall display the frequency resolved optical-gating (FROG) spectrogram of the TRUECARS signal given by ref. ⁵⁸, by convolving a temporal trace S(t) at a constant ω_r, with a Gaussian gating function E_gate(t) with a full width at half-maximum (fwhm) of 0.484 fs,S(T) oscillates with frequencies that correspond to the energy splitting between the relevant vibronic coherences, and the FROG spectrogram reveals the transient energy splitting along the trajectory. The FROG spectrograms are scanned and integrated over a negative Raman shift (ω_r < 0) window to capture the evolution of the signal away from ω_r = 0.

The gas phase (single-molecule) TRXD signal of a sample with N non-interacting molecules reads^33,69

wherewhere 〈ϕ_I^m|σ(−q,t)|ϕ_K^m〉 = (σ_IK^m)^† and 〈ϕ_Kⁿ|σ(q,t)|ϕ_Jⁿ〉 = σ_KJⁿ. We refer the reader to the ESI† for the derivation more detail.

AIMC simulations of COT have been performed at constant energy using a 0.05 fs time step. The initial conformational structures were sampled from a 520 ps ground-state adiabatic molecular dynamics trajectory in a vacuum using a Langevin thermostat at 300 K with a 0.1 fs time step.⁵⁶ Following a 20 ps equilibration period, 100 snapshots of geometries and velocities were harvested every 10 ps and used as the initial conditions for the AIMC non-adiabatic molecular dynamics simulation. Nuclear dynamics were simulated for all 42 nuclear degrees of freedom, from the initial excited S₃ state populated by an impulsive excitation. Note that S₁ is a dark state and S₂ and S₃ are both bright degenerate states, but the latter has 6 times larger oscillater strength than S₂ during the ground state dynamics (see Fig. 1e in ref. ⁵⁶). The S₄ state is located 0.47 eV higher than the S_2/3 state, and thus, one can safely exclude the excitation to the S₄ state by using a spectrally narrow enough optical pump pulse. If we follow what Levine et al. had done in their butadiene work,⁷⁰ our initial condition should shift 14.3% of the population to S₂ and the remaining major population to S₃. We do not expect substantial differences in TRUECARS and TRXD signals, but the relaxation timescale of COT could be a bit faster.

Excited state properties (e.g., energies, gradients, and non-adiabatic couplings) are calculated on-the-fly at the configuration interaction single (CIS) level of theory using the semiempirical Austin model 1 Hamiltonian.⁷¹ The non-adiabatic transitions to the ground state near S₀/S₁ CI have an inherent superposition or multireference character which cannot be properly described with single-reference CIS, time-dependent Hartree Fock or time-dependent density functional theory due to an incorrect description of the topology of the CI near the crossing (phase factors, etc).^72–74 An alternative approach is the “Open-GS” method^75,76 that enforces such a transition to the ground state when the energy gap between the ground and excited states is smaller than a certain threshold. Hence, we set up a 0.1 eV threshold value for S₁/S₀ CI description so that the AIMC simulation is terminated, once one of the clones reaches this point. Finally, 98 independent trajectories were averaged for further TRUECARS signal analysis, whereas 57 trajectories were used for the TRXD signal. We found that TRXD computations require extremely large data storage (e.g., each trajectory occupies around 66 GB, and the largest one takes up 200 GB) but averaging over 57 trajectories is enough to get converged results for TRXD signals. The computations of TRUECARS signals are not subjected to such computational cost, and we averaged over all 98 trajectories.

The use of transition polarizabilities, though computationally demanding, gives more accurate results compared with our previous work, where we set α_KL to be constant (geometry-independent) over the nuclear space, thereby reducing 〈Ψ(t)|α|Ψ(t)〉 to the overlap between the involved electronic states (vibronic coherence magnitude).^30,31

Data availability

All study data are included in the article and/or ESI.†

Author contributions

Y. N. calculated the X-ray signals and wrote the manuscript. H. S. simulated non-adiabatic molecular dynamics and wrote the manuscript. V. M. F. and S. F.-A. derived the TRUECARS and TRXD signals. D. K. helped in computing and analyzing the TRUECARS signal and its spectrogram. J. Y. L. and S. T. supervised the project and wrote the manuscript. M. G. suggested the cyclooctatetraene molecule and wrote part of introduction. S. M. designed and supervised the project and wrote the manuscript.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work was primarily supported by the U.S. Department of Energy, Office of Science, Office of Basic Energy Sciences under Award DE-SC0022225 and partly supported by the Center for Integrated Nanotechnologies (CINT), a U.S. Department of Energy, Office of Science User Facility. The support of the National Science Foundation through Grant No. CHE-2246379 is gratefully acknowledged (Y. N. and D. K.) S. M. is a fellow of Hagler institute for advanced study at Texas A&M University. H.S. and S.T. acknowledge support from the U.S. Department of Energy, Office of Science, Basic Energy Sciences, and Chemical Sciences, Geosciences, and Biosciences Division under Triad National Security, LLC (“Triad”) contract grant # 89233218CNA000001 (FWP: LANLE3T1). D.K. gratefully acknowledges support from the Alexander von Humboldt Foundation through the Feodor Lynen program. V.M.F. and S.F.-A acknowledge support from CONICET, UNQ, and ANPCyT (PICT-2018-2360).

Notes and references

1. G. A. Olah, Acc. Chem. Res., 1971, 4, 240–248.
2. A. Pron and P. Rannou, Prog. Polym. Sci., 2002, 27, 135–190.
3. Y. Nam, J. R. Rouxel, J. Y. Lee and S. Mukamel, Phys. Chem. Chem. Phys., 2020, 22, 26605–26613.
4. M. N. R. Ashfold, M. Bain, C. S. Hansen, R. A. Ingle, T. N. V. Karsili, B. Marchetti and D. Murdock, J. Phys. Chem. Lett., 2017, 8, 3440–3451.
5. A. Bhattacherjee, K. Schnorr, S. Oesterling, Z. Yang, T. Xue, R. de Vivie-Riedle and S. R. Leone, J. Am. Chem. Soc., 2018, 140, 12538–12544.
6. D. A. Hrovat and W. T. Borden, J. Am. Chem. Soc., 1992, 114, 5879–5881.
7. P. G. Wenthold, D. A. Hrovat, W. T. Borden and W. C. Lineberger, Science, 1996, 272, 1456–1459.
8. M. Garavelli, F. Bernardi, A. Cembran, O. Castaño, L. M. Frutos, M. Merchán and M. Olivucci, J. Am. Chem. Soc., 2002, 124, 13770–13789.
9. M. Garavelli, F. Bernardi, V. Moliner and M. Olivucci, Angew. Chem., Int. Ed., 2001, 40, 1466–1468.
10. I. J. Palmer, I. N. Ragazos, F. Bernardi, M. Olivucci and M. A. Robb, J. Am. Chem. Soc., 1993, 115, 673–682.
11. F. Bernardi, M. Olivucci and M. A. Robb, Chem. Soc. Rev., 1996, 25, 321–328.
12. H. E. Zimmerman and H. Iwamura, J. Am. Chem. Soc., 1970, 92, 2015–2022.
13. L. Li, M. Lei, Y. Xie, H. F. Schaefer, B. Chen and R. Hoffmann, Proc. Natl. Acad. Sci. U. S. A., 2017, 114, 9803–9808.
14. I. Aprahamian, ACS Cent. Sci., 2020, 6, 347–358.
15. A. A. S. Briquet, P. Uebelhart and H.-J. Hansen, Helv. Chim. Acta, 1996, 79, 2282–2315.
16. H. J. Wörner, J. B. Bertrand, B. Fabre, J. Higuet, H. Ruf, A. Dubrouil, S. Patchkovskii, M. Spanner, Y. Mairesse, V. Blanchet, E. Mével, E. Constant, P. B. Corkum and D. M. Villeneuve, Science, 2011, 334, 208–212.
17. M. S. Schuurman and A. Stolow, Annu. Rev. Phys. Chem., 2018, 69, 427–450.
18. D. Keefer, S. Thallmair, S. Matsika and R. de Vivie-Riedle, J. Am. Chem. Soc., 2017, 139, 5061–5066.
19. B. Gu, D. Keefer, F. Aleotti, A. Nenov, M. Garavelli and S. Mukamel, Proc. Natl. Acad. Sci. U. S. A., 2021, 118, e2116868118.
20. H. Timmers, X. Zhu, Z. Li, Y. Kobayashi, M. Sabbar, M. Hollstein, M. Reduzzi, T. J. Martínez, D. M. Neumark and S. R. Leone, Nat. Commun., 2019, 10, 3133.
21. J. Yang, X. Zhu, T. J. A. Wolf, Z. Li, J. P. F. Nunes, R. Coffee, J. P. Cryan, M. Gühr, K. Hegazy, T. F. Heinz, K. Jobe, R. Li, X. Shen, T. Veccione, S. Weathersby, K. J. Wilkin, C. Yoneda, Q. Zheng, T. J. Martinez, M. Centurion and X. Wang, Science, 2018, 361, 64–67.
22. C. Bressler and M. Chergui, Chem. Rev., 2004, 104, 1781–1812.
23. E. Goulielmakis, Z.-H. Loh, A. Wirth, R. Santra, N. Rohringer, V. S. Yakovlev, S. Zherebtsov, T. Pfeifer, A. M. Azzeer, M. F. Kling, S. R. Leone and F. Krausz, Nature, 2010, 466, 739–743.
24. A. Stolow, A. E. Bragg and D. M. Neumark, Chem. Rev., 2004, 104, 1719–1758.
25. D. Polli, P. Altoè, O. Weingart, K. M. Spillane, C. Manzoni, D. Brida, G. Tomasello, G. Orlandi, P. Kukura, R. A. Mathies, M. Garavelli and G. Cerullo, Nature, 2010, 467, 440–443.
26. M. Kowalewski, K. Bennett, K. E. Dorfman and S. Mukamel, Phys. Rev. Lett., 2015, 115, 193003.
27. D. Keefer, T. Schnappinger, R. de Vivie-Riedle and S. Mukamel, Proc. Natl. Acad. Sci. U. S. A., 2020, 117, 24069–24075.
28. S. M. Cavaletto, D. Keefer and S. Mukamel, Phys. Rev. X, 2021, 11, 011029.
29. Y. Nam, D. Keefer, A. Nenov, I. Conti, F. Aleotti, F. Segatta, J. Y. Lee, M. Garavelli and S. Mukamel, J. Phys. Chem. Lett., 2021, 12, 12300–12309.
30. D. Keefer, V. M. Freixas, H. Song, S. Tretiak, S. Fernandez-Alberti and S. Mukamel, Chem. Sci., 2021, 12, 5286–5294.
31. V. M. Freixas, D. Keefer, S. Tretiak, S. Fernandez-Alberti and S. Mukamel, Chem. Sci., 2022, 13, 6373–6384.
32. J. Cao and K. R. Wilson, J. Phys. Chem. A, 1998, 102, 9523–9530.
33. K. Bennett, J. D. Biggs, Y. Zhang, K. E. Dorfman and S. Mukamel, J. Chem. Phys., 2014, 140, 204311.
34. S. Bratos, F. Mirloup, R. Vuilleumier and M. Wulff, J. Chem. Phys., 2002, 116, 10615–10625.
35. M. Simmermacher, A. Moreno Carrascosa, N. E. Henriksen, K. B. Møller and A. Kirrander, J. Chem. Phys., 2019, 151, 174302.
36. G. Dixit, O. Vendrell and R. Santra, Proc. Natl. Acad. Sci. U. S. A., 2012, 109, 11636–11640.
37. A. Kirrander, K. Saita and D. V. Shalashilin, J. Chem. Theory Comput., 2016, 12, 957–967.
38. D. Arnlund, L. C. Johansson, C. Wickstrand, A. Barty, G. J. Williams, E. Malmerberg, J. Davidsson, D. Milathianaki, D. P. DePonte, R. L. Shoeman, D. Wang, D. James, G. Katona, S. Westenhoff, T. A. White, A. Aquila, S. Bari, P. Berntsen, M. Bogan, T. B. van Driel, R. B. Doak, K. S. Kjær, M. Frank, R. Fromme, I. Grotjohann, R. Henning, M. S. Hunter, R. A. Kirian, I. Kosheleva, C. Kupitz, M. Liang, A. V. Martin, M. M. Nielsen, M. Messerschmidt, M. M. Seibert, J. Sjöhamn, F. Stellato, U. Weierstall, N. A. Zatsepin, J. C. H. Spence, P. Fromme, I. Schlichting, S. Boutet, G. Groenhof, H. N. Chapman and R. Neutze, Nat. Methods, 2014, 11, 923–926.
39. K. J. Gaffney and H. N. Chapman, Science, 2007, 316, 1444–1448.
40. M. P. Minitti, J. M. Budarz, A. Kirrander, J. S. Robinson, D. Ratner, T. J. Lane, D. Zhu, J. M. Glownia, M. Kozina, H. T. Lemke, M. Sikorski, Y. Feng, S. Nelson, K. Saita, B. Stankus, T. Northey, J. B. Hastings and P. M. Weber, Phys. Rev. Lett., 2015, 114, 255501.
41. B. Stankus, H. Yong, N. Zotev, J. M. Ruddock, D. Bellshaw, T. J. Lane, M. Liang, S. Boutet, S. Carbajo, J. S. Robinson, W. Du, N. Goff, Y. Chang, J. E. Koglin, M. P. Minitti, A. Kirrander and P. M. Weber, Nat. Chem., 2019, 11, 716–721.
42. H. Yong, N. Zotev, B. Stankus, J. M. Ruddock, D. Bellshaw, S. Boutet, T. J. Lane, M. Liang, S. Carbajo, J. S. Robinson, W. Du, N. Goff, Y. Chang, J. E. Koglin, M. D. J. Waters, T. I. Sølling, M. P. Minitti, A. Kirrander and P. M. Weber, J. Phys. Chem. Lett., 2018, 9, 6556–6562.
43. J. Yang, R. Dettori, J. P. F. Nunes, N. H. List, E. Biasin, M. Centurion, Z. Chen, A. A. Cordones, D. P. Deponte, T. F. Heinz, M. E. Kozina, K. Ledbetter, M.-F. Lin, A. M. Lindenberg, M. Mo, A. Nilsson, X. Shen, T. J. A. Wolf, D. Donadio, K. J. Gaffney, T. J. Martinez and X. Wang, Nature, 2021, 596, 531–535.
44. E. Biasin, Z. W. Fox, A. Andersen, K. Ledbetter, K. S. Kjær, R. Alonso-Mori, J. M. Carlstad, M. Chollet, J. D. Gaynor, J. M. Glownia, K. Hong, T. Kroll, J. H. Lee, C. Liekhus-Schmaltz, M. Reinhard, D. Sokaras, Y. Zhang, G. Doumy, A. M. March, S. H. Southworth, S. Mukamel, K. J. Gaffney, R. W. Schoenlein, N. Govind, A. A. Cordones and M. Khalil, Nat. Chem., 2021, 13, 343–349.
45. P. Zalden, F. Quirin, M. Schumacher, J. Siegel, S. Wei, A. Koc, M. Nicoul, M. Trigo, P. Andreasson, H. Enquist, M. J. Shu, T. Pardini, M. Chollet, D. Zhu, H. Lemke, I. Ronneberger, J. Larsson, A. M. Lindenberg, H. E. Fischer, S. Hau-Riege, D. A. Reis, R. Mazzarello, M. Wuttig and K. Sokolowski-Tinten, Science, 2019, 364, 1062–1067.
46. Y. Lee, J. G. Kim, S. J. Lee, S. Muniyappan, T. W. Kim, H. Ki, H. Kim, J. Jo, S. R. Yun, H. Lee, K. W. Lee, S. O. Kim, M. Cammarata and H. Ihee, Nat. Commun., 2021, 12, 3677.
47. E. Biasin, T. B. van Driel, K. S. Kjær, A. O. Dohn, M. Christensen, T. Harlang, P. Vester, P. Chabera, Y. Liu, J. Uhlig, M. Pápai, Z. Németh, R. Hartsock, W. Liang, J. Zhang, R. Alonso-Mori, M. Chollet, J. M. Glownia, S. Nelson, D. Sokaras, T. A. Assefa, A. Britz, A. Galler, W. Gawelda, C. Bressler, K. J. Gaffney, H. T. Lemke, K. B. Møller, M. M. Nielsen, V. Sundström, G. Vankó, K. Wärnmark, S. E. Canton and K. Haldrup, Phys. Rev. Lett., 2016, 117, 013002.
48. K. Haldrup, G. Levi, E. Biasin, P. Vester, M. G. Laursen, F. Beyer, K. S. Kjær, T. Brandt van Driel, T. Harlang, A. O. Dohn, R. J. Hartsock, S. Nelson, J. M. Glownia, H. T. Lemke, M. Christensen, K. J. Gaffney, N. E. Henriksen, K. B. Møller and M. M. Nielsen, Phys. Rev. Lett., 2019, 122, 063001.
49. D. Keefer, F. Aleotti, J. R. Rouxel, F. Segatta, B. Gu, A. Nenov, M. Garavelli and S. Mukamel, Proc. Natl. Acad. Sci. U. S. A., 2021, 118, e2022037118.
50. J. R. Rouxel, D. Keefer, F. Aleotti, A. Nenov, M. Garavelli and S. Mukamel, J. Chem. Theory Comput., 2022, 18, 605–613.
51. T. R. Nelson, A. J. White, J. A. Bjorgaard, A. E. Sifain, Y. Zhang, B. Nebgen, S. Fernandez-Alberti, D. Mozyrsky, A. E. Roitberg and S. Tretiak, Chem. Rev., 2020, 120, 2215–2287.
52. V. M. Freixas, A. J. White, T. Nelson, H. Song, D. V. Makhov, D. Shalashilin, S. Fernandez-Alberti and S. Tretiak, J. Phys. Chem. Lett., 2021, 12, 2970–2982.
53. J. T. O'Neal, E. G. Champenois, S. Oberli, R. Obaid, A. Al-Haddad, J. Barnard, N. Berrah, R. Coffee, J. Duris, G. Galinis, D. Garratt, J. M. Glownia, D. Haxton, P. Ho, S. Li, X. Li, J. MacArthur, J. P. Marangos, A. Natan, N. Shivaram, D. S. Slaughter, P. Walter, S. Wandel, L. Young, C. Bostedt, P. H. Bucksbaum, A. Picón, A. Marinelli and J. P. Cryan, Phys. Rev. Lett., 2020, 125, 073203.
54. D. Cho, J. R. Rouxel and S. Mukamel, J. Phys. Chem. Lett., 2020, 11, 4292–4297.
55. Y. Nam, F. Montorsi, D. Keefer, S. M. Cavaletto, J. Y. Lee, A. Nenov, M. Garavelli and S. Mukamel, J. Chem. Theory Comput., 2022, 18, 3075–3088.
56. H. Song, Y. Nam, D. Keefer, M. Garavelli, S. Mukamel and S. Tretiak, J. Phys. Chem. Lett., 2021, 12, 5716–5722.
57. K. Kosma, S. A. Trushin, W. E. Schmid and W. Fuß, Chem. Phys., 2015, 463, 111–119.
58. R. Trebino, K. W. DeLong, D. N. Fittinghoff, J. N. Sweetser, M. A. Krumbügel, B. A. Richman and D. J. Kane, Rev. Sci. Instrum., 1997, 68, 3277–3295.
59. S. M. Cavaletto, D. Keefer, J. R. Rouxel, F. Aleotti, F. Segatta, M. Garavelli and S. Mukamel, Proc. Natl. Acad. Sci. U. S. A., 2021, 118, e2105046118.
60. V. M. Freixas, S. Fernandez-Alberti, D. V. Makhov, S. Tretiak and D. Shalashilin, Phys. Chem. Chem. Phys., 2018, 20, 17762–17772.
61. H. Song, V. M. Freixas, S. Fernandez-Alberti, A. J. White, Y. Zhang, S. Mukamel, N. Govind and S. Tretiak, J. Chem. Theory Comput., 2021, 17, 3629–3643.
62. W. Malone, B. Nebgen, A. White, Y. Zhang, H. Song, J. Bjorgaard, A. Sifain, B. Rodriguez-Hernandez, V. Freixas, S. Fernandez-Alberti, E. Roitberg, R. Nelson and S. Tretiak, J. Chem. Theory Comput., 2020, 16, 5771–5783.
63. D. V. Shalashilin, J. Chem. Phys., 2009, 130, 244101.
64. M. Ben-Nun, J. Quenneville and T. J. Martínez, J. Phys. Chem. A , 2000, 104, 5161–5175.
65. B. F. E. Curchod and T. J. Martínez, Chem. Rev., 2018, 118, 3305–3336.
66. S. Fernandez-Alberti, D. V. Makhov, S. Tretiak and D. V. Shalashilin, Phys. Chem. Chem. Phys., 2016, 18, 10028–10040.
67. D. V. Makhov, W. J. Glover, T. J. Martinez and D. V. Shalashilin, J. Chem. Phys., 2014, 141, 054110.
68. D. V. Makhov, K. Saita, T. J. Martinez and D. V. Shalashilin, Phys. Chem. Chem. Phys., 2015, 17, 3316–3325.
69. K. Bennett, M. Kowalewski, J. R. Rouxel and S. Mukamel, Proc. Natl. Acad. Sci. U. S. A., 2018, 115, 6538–6547.
70. B. G. Levine and T. J. Martínez, J. Phys. Chem. A, 2009, 113, 12815–12824.
71. J. J. P. Stewart, J. Mol. Model., 2004, 10, 155–164.
72. B. G. Levine, C. Ko, J. Quenneville and T. J. MartÍnez, Mol. Phys., 2006, 104, 1039–1051.
73. L. Stojanović, S. Bai, J. Nagesh, A. F. Izmaylov, R. Crespo-Otero, H. Lischka and M. Barbatti, Molecules, 2016, 21, 1603.
74. S. M. Parker, S. Roy and F. Furche, Phys. Chem. Chem. Phys., 2019, 21, 18999–19010.
75. Y. Zhang, S. Mukamel, M. Khalil and N. Govind, J. Chem. Theory Comput., 2015, 11, 5804–5809.
76. F. Furche and J. P. Perdew, J. Chem. Phys., 2006, 124, 044103.

---

Footnotes

† Electronic supplementary information (ESI) available: Signal derivation, loop diagram rule, and supplementary figures. See DOI: https://doi.org/10.1039/d2sc04335a

‡ These authors contributed equally to this work.

---