Ultrafast coherent photoexcited dynamics in a trimeric dendrimer probed by X-ray stimulated-Raman signals

The photoinduced ultrafast coherent inter-chromophore energy redistribution in a triarylamine trimer is explored using nonadiabatic excited state molecular dynamics followed by simulations of X-ray Raman signals. The nitrogencentered system ensures strong interchromophore interactions and, thus, the presence of coherences. Nevertheless, the multitude of non-deterministic photoinduced pathways during the ultrafast inter-branch migration of the excitation results in random confinement on some branches and, therefore, spatial exciton scrambling and loss of phase information at long times. We show that the vibronic coherence dynamics evolving into the incoherent scrambling mechanism on ultrafast 50 fs timescale, is accurately probed by the TRUECARS X-ray stimulated Raman signal. In combination with previous results, where the technique has revealed long-lived coherences in a rigid heterodimer, the signal is most valuable for detecting ultrafast molecular coherences or their absence. We demonstrate that X-ray Raman spectroscopy is a useful tool in the chemical design of functional molecular building blocks.


A Cloning criteria and thresholds
• Criteria #1: Cloning should take pace in situations where more than one excited state is significantly populated. This can be monitored by the distribution width W n [1]: (S1) W n takes values between 1 and the total number of excited states considered N . Values near 1 means that all the population is concentrated in a single state, while values near N means an even distribution. Therefore, the cloning events were restricted to situations where: The time evolution of the expectation value of W n might be misleading in the sense that it has a broad distribution and presents sudden changes after cloning events. Moreover, it only quantifies whether the electronic populations for a given trajectory spread over the different states or not. A more decisive criterion is the breaking angle defined in criteria #2: the populations might be distributed over several states, but if all those states pull the nuclear system in the same direction of motion there is no need for bifurcating the wavefunction. Figure S1 depicts the distribution of values of W n for all time steps and trajectories of the complete ensemble for the current case. The peak at approximately W n = 2 is associated to the threshold δ 1 = 2, while the peak at approximately W n = 1 is associated either to the subsequent reduction of W n after each cloning event or to values acquired after electronic energy relaxation to the S 1 state.  Figure S1: Probability density function (PDF) of W n defined according to (eq. S1) for all trajectories and time steps of the complete ensemble for the current case.
• Criteria #2: Cloning should take place in situations where the wavepacket would split in configuration space. This is monitored by the breaking angle [1]: where F (n) M is the averaged excited state gradient for configuration n: where R n represents the nuclear coordinates for configuration n and V (n) I represents the energy corresponding to the excited state I of configuration n. F (n) max is minus the gradient corresponding to the most populated state. By restricting the cloning events to situations in which: we ensure that the two new configurations produced after a cloning event will split in configuration space providing the desired bifurcation effect and rectifying the physical meaning of the mean field.
• Criteria #3: In order to make the algorithm efficient and avoid a overproduction of cloning events in the strong coupling regions of configuration space, we avoid cloning when the coupling is to high. By this way clones are produced once the system leaves strong coupling regions if the two previous criteria are fulfilled. The coupling is monitored in similar fashion as the quantum transition probability is calculated for the fewest switches surface hopping algorithm [2,3]: where σ I,max is the non-adiabatic coupling between excited states I and the one with the greater population for configuration n: (S7) Therefore, by setting δ 3 = 0.05 we ensure an efficient clone generation away from the strong coupling regions of configuration space.
• Criteria #4: The nuclear amplitude corresponding to each configuration splits after a cloning event. If the total number of clones is not controlled, an exponential growth, with an exponentially decreasing weight for the corresponding configurations, would make the algorithm too inefficient. Therefore, we allowed a total of 4 consecutive clones, leading to a maximum of 16 clones per initial condition. Previous convergence tests have shown that there is no a relevant accuracy improvement after the production of approximately 14 clones per initial condition [2]. In the current case, an average of approximately 6 clones per initial condition were generated. These clones take place through out all the simulation time. Figure S2 shows that cloning events occurs during the simulations in an homogeneous rate.

B Fraction of transition density
Within the Collective Electronic Oscillator (CEO) approach the Configuration Interaction Singlet (CIS) eigenstate I for configuration n, written in Atomic Orbitals (AO) basis, is denoted frequently as transition density matrices [4,5,6]: where |φ (n) g is the ground state wavefunction, andĉ † i andĉ j are the electron creation and annhilation operators with indexes i and j referring to AO basis functions. Diagonal elements ρ are relevant to the changes in the distribution of electronic density in the i th orbital in the case of bound excitonic states caused by excitation [7]. During AIMC simulations, the intramolecular electronic energy redistribution can be followed using the time-dependent spatial localization of ρ (n) I . The fraction of transition density localized on a specific segment X of the molecule can be defined as [1]: In order to calculate the corresponding expectation value for the Multiconfigurational Ehrenfest (MCE) wavefunction we can introduce the operatorρ X such that: which expectation value can be calculated as: where stands for the real part.

S5
C Parameters of the fit to eq. (1) The function in eq. (1) describes the relaxation process to S 1 as a superposition of an exponential decay plus an exponentially damped oscillation, such that both contributions have the same weight at t = 0. Parameters a and b corresponds to the asymptotic value of S 1 population and its initial value, respectively; τ 1 is the S 1 decay time, T is the period of the vibronic oscillations and τ 2 is the vibronic decoherence time. R is the corresponding R-squared.  Table S1: Parameters obtained by fitting the population of S 1 in Figure 2(a) to f (t) given by eq. (1).

D Pure electronic dephasing
Within linear response theory [8], the correlation between excited states I and J is quantified by the autocorrelation function: where δE IJ is the corresponding energy gap and the angular brackets means average over the ensemble. The energy corresponding to a given state K is calculated as the expectation value: where E (n) K is the energy corresponding to the adiabatic state K for the center of configuration n. The pure electronic dephasing D IJ (t) is defined by the second order cumulant approximation [8]: Finally, the pure electronic decoherence time τ IJ can be obtained by fitting the pure dephasing D IJ to a Gaussian function:

E Surface hopping simulations
NEXMD surface hopping (SH) simulations were performed in order to compare the capability of the SH method to reproduce the oscillatory behavior of the electronic populations resulting from vibronic coherences. The same ensemble of initial configurationshas as for AIMC simulations has been considered . The electronic wavefunction was collapsed to S 2 at t = 0 in order to maximize the appearance of any possible oscillations. Figure (S4) shows the electronic populations, calculated as the fraction of trajectories evolving on each state at any given time. NEXMD-SH simulations were performed at a constant energy. Classical time steps of 0.1 fs have been used for the propagation of nuclei in ground state and a quantum time step of 0.025 fs has been used to propagate the electronic degrees of freedom. Corrections for decoherence [9], parameters and the methodology of NEXMD-SH simulations have been extensively discussed elsewhere [10]. Figure S4: Time evolution of the electronic populations for SH dynamics.

F Transition density matrix snapshots
In order to calculate transition density matrix for the MCE wavefunction |Ψ , we first calculate the expectation value of the transition density matrix from the ground state to a given state K: where ρ (n) 0K is the transition density matrix [7] from the ground state to excited state K for the center of configuration n.
The transition density matrix ρ from the ground state to the superposition of excited states is then calculated as: where a K is calculated according to [1]: (S18) S9 Figure S5: Transition density matrices (eq. (S17)) corresponding to the the different snapshots shown in Figure 3(a).

G Non-adiabatic coupling analysis
A further insight into the vibronic dynamics can be achieved by analyzing the nonadiabatic coupling vector (NACR) in regions of strong coupling near conical intersection seams. The strong coupling regime is defined as situations where the absolute value of NACR is higher than its mean value plus 10 times the standard deviation during each AIMC trajectory. Figure S6(a) shows a real-space representation of a typical NACR. We find that NACR involves antisymmetric nuclear motions on two of the three branches. This behavior involving antisymmetric motions on different branches of a dendrimer has been reported previously [11] and was associated to the fact that the molecular wavefunction adopts a standing wave pattern according to the particle (exciton) in a box model and the S 2 → S 1 transition can be associated with a transition where symmetry between states changes. That is, the vibronic excitation has an antisymmetric form where the two branches experience structural deformations with opposite phases (i.e. expansion and compression). Figure  S6(b) depicts the probability density function of the fraction of NACR absolute values localized on the three branches, denoted as I, II, and III according to their corresponding contribution values in a decreasing order. We can observe that the NACRs are not equally distributed among the three branches but localized mainly in two of them. The distribution of the projections of NACR onto the normal modes basis, calculated at the S 1 minimum energy configuration, is shown in Figure  S6(c). The peak at ∼ 19 fs (∼ 1750 cm −1 ) matches the period T of the vibronic oscillations (∼ 20 fs) obtained by fitting the evolution in time of the S 1 population (see Figure 2(a) and Table S1). Vibronic dynamics are thus assisted by a reduced set of high-frequency normal modes.

H Additional Trajectories
Only trajectory 2 is discussed in the main text. This trajectory is representative of a large fraction of trajectories in the ensemble, with only minor differences especially in the signal. Here, we discuss a few interesting trajectories to give a better impression about the range of properties during the dynamics.
• Figure S7 shows a trajectory that is similar to trajectory 2 in the main text, but where the transition density is located in branches A and C instead of B and C. This exemplifies that most trajectories are quite similar and the three branches of the molecule are equal, while the transition density is localized in two of them simultaneously.
• Figure S8 depicts a trajectory where the ρ 13 contribution to the signal is equally strong to the usually dominating ρ 12 contribution. Thus, they are both well visible in the total signal of this trajectory and well distinguishable at different Raman shifts.
• In Figure S9, trajectory 5 exhibits much larger and more frequent population oscillations, especially between S 1 and S 2 . This does not considerably affect the relative magnitudes and persistence of the individual signal contributions, however. Figure S7: Dynamics and TRUECARS signal in trajectory 3. a) Populations in the S 1 (blue), S 2 (red) and S 3 (yellow) states. b) Fraction of transition density in branch A (purple), B (yellow), C (red) and around the central Nitrogen (blue). c) Frequency-dispersed TRUECARS signal (eq. 2) for the coherence ρ 12 between S 1 and S 2 . d-e) same as c) for ρ 13 and ρ 23 . f ) Total TRUECARS signal for trajectory 2 given by the sum of c-e. S14 Figure S8: Dynamics and TRUECARS signal in trajectory 4. a) Populations in the S 1 (blue), S 2 (red) and S 3 (yellow) states. b) Fraction of transition density in branch A (purple), B (yellow), C (red) and around the central Nitrogen (blue). c) Frequency-dispersed TRUECARS signal (eq. 2) for the coherence ρ 12 between S 1 and S 2 . d-e) same as c) for ρ 13 and ρ 23 . f ) Total TRUECARS signal for trajectory 2 given by the sum of c-e. S15 Figure S9: Dynamics and TRUECARS signal in trajectory 5. a) Populations in the S 1 (blue), S 2 (red) and S 3 (yellow) states. b) Fraction of transition density in branch A (purple), B (yellow), C (red) and around the central Nitrogen (blue). c) Frequency-dispersed TRUECARS signal (eq. 2) for the coherence ρ 12 between S 1 and S 2 . d-e) same as c) for ρ 13 and ρ 23 . f ) Total TRUECARS signal for trajectory 2 given by the sum of c-e.