Photoinduced non-adiabatic energy transfer pathways in dendrimer building blocks.

The efﬁciency of the intramolecular energy transfer in light harvesting dendrimers is determined by their well-deﬁned architecture with high degree of order. After photoexcitation, through-space and through-bond energy transfer mechanisms can take place, involving vectorial exciton migration among different chromophores within dendrimer highly branched structures. Their inherent intramolecular energy gradient depends on how the multiple chromophoric units have been assembled, subject to their inter-connects, spatial distances, and orientations. Herein, we compare the photoinduced nonadiabatic molecular dynamics simulations performed on a set of different combinations of a chain of linked dendrimer building blocks composed of two-, three-, and four-ring linear polyphenylene chromophoric units. The calculations are performed with the recently developed ab initio multiple cloning-time dependent diabatic basis implementation of the Multiconﬁgurational Ehrenfest (MCE) approach. Despite differences in short time relaxation pathways and different initial exciton localization, at longer time scales, electronic relaxation rates and exciton ﬁnal redistributions are very similar for all combinations. Unlike the systems composed of two building blocks, considered previously, for the larger 3 block systems here we observe that bifurcation of the wave function accounted by cloning is important. In all the systems considered in this work, at the time scale of few hundreds of femtoseconds, cloning enhances the electronic energy relaxation by ∼ 13% compared to that of the MCE method without cloning. Thus, accurate description of quantum effects is essential for understanding of the energy exchange in dendrimers both at short and long time scales.


I. INTRODUCTION
An efficient conversion of solar radiation into other usable forms of energy promises an unlimited clean energy source. Such conversion is achieved in living organisms by means of complex arrays of conjugated chromophores. 1,2 Synthesis of artificial organic materials [3][4][5] has shown a glimpse of the opportunity to mimic nature light-harvesting capabilities and high efficiency energy funnelling. Dendrimers are conjugated macromolecules with highly branched structures that can perform such a function. In particular, those based on poly(phenylene ethynylene) (PPE) have both the collection and transport features present in photosynthetic systems. [6][7][8][9][10] Linear PPE chromophore units with different conjugation lengths, linked by meta-substitution at the branching phenylene nodes, constitute the building blocks of the well-studied perylene-terminated dendrimer known as the nanostar. [11][12][13] These meta-branching vertices localize excitons within each linear PPE unit, allowing the study of PPE dendrimers as an ensemble of independent chromophore units with a weak coupling between them. 8,11,12 The dendrimer's well-defined tree-like architecture is crucial to achieve their highly efficient intramolecular energy funnelling. 7,11,13,14 On the one hand, compact PPE dendrimers composed of equivalent linear chromophore units do not experience an efficient periphery to core energy transfer. Their photoexcitation leads to an ultrafast spatial molecular scrambling of the excited state wave function with the exciton equally distributed among all linear chromophores. [15][16][17][18][19][20][21] On the other hand, extended PPE dendrimers, composed of PPE units of different conjugation lengths, undergo highly efficient and unidirectional energy transfer from higher energy (i.e., shorter) to lower energy (i.e., longer) units. [22][23][24][25] That is, the efficiency of the unidirectional energy transfer depends on how the different chromophoric units are assembled, into a supramolecular structure.
Subsequent to photoexcitation, PPE dendrimers display an ultrafast intramolecular energy relaxation and redistribution involving multiple possible pathways in which through-bond and throughspace energy transfer mechanisms can be distinguished. 24,26 Photophysical processes like exciton transient delocalization/localization, intramolecular migration and funnelling, hopping, and self-trapping can occur concurrently during energy relaxation. Such processes can be described by on-the-fly non-adiabatic molecular dynamic methods, 27,28 ranging in the approximation level and the computational expense. Surface Hopping (SH) algorithms [29][30][31][32] are among the less computationally expensive and have been used to study a wide variety of organic molecules including dendrimers. More recently, Quantum Direct Dynamics (QDD) 27 methods emerged as more accurate alternatives, however carrying extra computational expense. Compared to SH methods, the QDD techniques improve by including the effects of nuclear quantum dynamics. This is achieved by considering ensembles of coupled trajectory-guided Gaussian basis functions (TBFs) that attempt to cover the most important regions of configuration space.
Among different QDD methods, 33 the variational Multi-Configuration Gaussian (vMCG) approach, 34 which makes use of coupled non-classical variational trajectories, is computationally more expensive and has been applied only to relatively small molecules. The Ab Initio Multiple Spawning (AIMS) 35,36 makes use of TBFs that evolve in a simplified manner, allowing the expansion of the original TBFs according to transient coupling between states. Finally, the Multiconfigurational Ehrenfest (MCE), 37 which lies in between the vMCG and AIMS approaches, considers independent Ehrenfest mean field trajectories whose amplitude exchanges are obtained a posteriori, allowing its highly parallelized implementation. Within the MCE approach, its Time-Dependent Diabatic Basis (TDDB) implementation allows electronic excited states to naturally swap at the trivial unavoided crossings, 38 which is essential for simulations of photodynamics in multichromophoric molecular systems. Besides, different sampling techniques have been proposed to ensure the system to span through all relevant regions of configuration space. 39 In particular, the Ab Initio Multiple Cloning (AIMC) 40 sampling technique takes into account situations in which the Ehrenfest trajectories guided by a potential energy surface (PES) average become unphysical. More recently, a new AIMC-TDDB implementation of the MCE approach has been developed in order to deal with photoinduced dynamics simulations in large conjugated molecules. 40 It combines the MCE algorithms with onthe-fly calculations of excited-state energies, gradients, and nonadiabatic couplings using the Collective Electronic Oscillator (CEO) approach. 6,[41][42][43][44] In this paper, we compare photoinduced nonadiabatic molecular dynamics simulations performed on a set of different combinations of dendrimer building blocks composed of two-, three-, and four-ring linear PPE chromophoric units linked by meta-substitutions, which we denote as 234PPE, 243PPE, and 324PPE (see the insets of Fig. 1). By using the AIMC-TDDB 40 implementation of the MCE framework, we ensure that the excited states are suited to swap at trivial unavoided crossings and that the wave function is able to bifurcate in order to represent different relaxation pathways at once. Energy transfer is then monitored by following the changes in the spatial localization of the electronic transition density (TD). The energy flow between the individual chromophores is studied by means of the Statistical Minimum Flow (SMF), 45 developed originally for the analysis of the vibrational energy flow in polyatomic molecules, and recently adapted for the TD flux analysis. 46 By combining these techniques, we were able to obtain further insights into the role of the spatial organization of the building blocks of the nanostar and the energy transfer and the non-radiative relaxation pathways. In particular, we show that at the time scale of a few tens of femtoseconds, the energy transfer from the initially excited 2-ring unit is the fastest for 234PPE combination of building blocks. It reveals fast population exchange with the excited state localised on the 4-ring block. In the dendrimer nanostar, 2-ring segments are usually located at the ends of branches on the surface of the nanostar. Such fast energy transfer may increase the efficiency of light harvesting by funnelling the excitation toward the centre of the dendrimeric structure away from its surface, where the excitation can be quenched by the environment. At a longer time scale of several hundreds of femtoseconds, the energy relaxation in all three combinations (234PPE, 243PPE, and 324PPE) is similar. In all cases, we observe that correct description of wave function bifurcations is important and results in noticeable effects on the calculated energy relaxation rates at longer time scales. Thus, accurate representation of quantum nonadiabatic dynamics appears to be important for modeling of the light harvesting function in dendrimers.

ARTICLE scitation.org/journal/jcp
The paper is organized as follows: In Sec. II, we present an overview of the theory and working equations of Ehrenfest and MCE dynamics, together with numerical details related by our simulations. Our results are presented and discussed in Sec. III. Finally, our conclusions are summarized in Sec. IV.

A. Ehrenfest dynamics
For Ehrenfest dynamics (EHR), the configuration n of a given molecular system is described by means of a wave function ψn(t)⟩, where χn(t)⟩ is the nuclear part given by a Coherent State 47 (CS) moving along the classical mean field, while the electronic part is spanned in the adiabatic basis (n) The interaction between these two subspaces is given by a mean field having a weighted contribution from all excited states and the coupling between them. The force F driving the configuration n is such that the overall energy of the system is conserved, where V The first term Eq. (2) stands for the energy changes in the adiabatic states, while the second term stands for population exchange between the adiabatic states involved.
Along with the classical propagation of the nuclear coordinates, adiabatic amplitudes a

B. Multiconfigurational Ehrenfest (MCE) dynamics
Multiconfigurational Ehrenfest (MCE) is a generalization of EHR in which the wave function is written as a linear combination of Ehrenfest configurations (1), each of them moving along its own Ehrenfest (mean field) trajectory (2), The coherent state χn(t)⟩ in (1) is a Gaussian wave packet, which is delocalised. For large conjugated molecules, adiabatic states can change significantly within the CS width and the electronic overlaps must be taken into account. Two adiabatic states can also suddenly exchange their spatial localization at trivial unavoided crossings. Therefore, excited states must be suited to swap at these points in order to keep track of the right electronic overlaps. It is important to stress that although the electronic solutions are the adiabatic ones for the centre of the CS, they are suited to swap at trivial unavoided crossings following the diabatic path. Therefore, a more convenient representation of the electronic parts of the wave function is given by the Time Dependent Diabatic Basis (TDDB). 56 Each configuration evolves on its own mean field PES driven by force (2), allowing a highly efficient trivial parallelization in the numerical implementation. Configuration amplitudes and electronic overlaps can be calculated and postprocessed afterwards. However, for systems having different relaxation pathways, the average Ehrenfest PES does not represent its individual contributions to the overall relaxation process. A more adequate treatment of these individual contributions can be achieved by using the Ab Initio Multiple Cloning (AIMC), 40,57 which allows us to take into account situations in which the average Ehrenfest potential becomes unphysical, that is, regions of low non-adiabatic couplings presenting significant differences between the shapes of the potential energy surfaces of two or more substantially populated electronic states. This typically happens when a configuration passes through/near a conical intersection or region of strong nonadiabatic coupling. In this case, a sole nuclear wave packet splits into multiple parts, each dominated by a single electronic state.
The AIMC procedure quantifies how well the mean field represents the individual excited states dynamics. Whenever the mean field fails, the corresponding configuration is replaced by two new configurations, having the same nuclear wave function but different electronic populations and, therefore, their own distinct mean fields. Amplitudes corresponding to these new configurations are defined such that the wave function remains continuous at the cloning point. After the cloning procedure, each new configuration follows its own mean field and separates from the other, providing the desired bifurcation effect. More details and technical implementations of the AIMC method can be found elsewhere. 40,[57][58][59] C. Calculation of observables The MCE observables are obtained as expectation values of a given operator acting on the electronic subspace, the nuclear subspace, or the entire configuration. In the TDDB approach, the expectation value of an arbitrary operatorN acting over the electronic subspace can be calculated as 56 The population of the electronic excited state K can be represented using the population operatorPK, which leads to 56 The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp Within the CEO approach, each CIS eigenstate is represented by the respective TD matrix 44 written in the Atomic Orbital (AO) basis as ρ where (n) g ⟩ represents the electronic ground state wave function, whileĉ † i andĉj are the electron creation and annihilation operators, respectively, and indices i and j refer to the original atomic orbital (AO) basis functions. The diagonal elements ρ are related to the changes in the electronic density induced by photoexcitation to the Ith excited state. 44 In order to monitor the intramolecular energy redistribution that takes place during the non-radiative relaxation pathways, we need to monitor the time-evolution of the fraction of TD spatially localized over each chromophore unit. For this purpose, the operator ρX is defined such thatρ X (n) with eigenvalues ρ (n) I,X being the normalized fraction of the TD corresponding to excited state I localized in the chromophoric fragment X of the configuration n. That is, Therefore, the corresponding expectation value ⟨ρX⟩ is given by 40 I,X . (12) In order to calculate the fraction of TD corresponding to a single state K over a given chromophoric fragment X, we can define the operatorρ (K) X acting on the entire configuration space such that where the eigenvalues ρ (n) K,X are given by Eq. (11). The corresponding expectation value is then Finally, we calculate the chromophore-unit participation number 60,61 as where the subscript X on the summation represents the first and second units andρ is the operator associated with the TD corresponding to state K over the fragment X, which expectation value is calculated according to Eq. (12).

D. Statistical minimum flow (SMF)
Different energy transfer pathways can be identified using the statistical minimum flow (SMF) method. 45,46 The antisymmetric flow matrix f XY carries all the information related with the TD transfer between different chromophore units of the molecule. At each time interval ∆t during dynamics simulations, we classify the units according to whether their ⟨ρX⟩ increase or decrease. If ⟨ρX⟩ reduces, we say that the fragment X behaves as a Donor (D), while if ⟨ρX⟩ increases it behaves as an Acceptor (A). Assuming no transfer between donors or between acceptors and defining the transfer from D to A as positive, we have where ∆δX(t) is the change in the TD localized over the fragment X during ∆t, and ∆δ total is the total transfer defined as The accumulated flux FXY (t) from X to Y can be calculated as More details about the SMF method can be found elsewhere. 45,46 E. Computational details The photoinduced excited-state intramolecular energy transfer that takes place in a set of different combinations of dendrimer building blocks has been studied using the EHR and AIMC-TDDB simulations. Our molecular systems in question are composed of two-, three-, and four-ring linear PPE chromophore units linked by meta-substitutions. We consider 234PPE, 243PPE, and 324PPE molecules as depicted in the insets (a)-(c) of Fig. 1, respectively. For each molecular system, the initial conditions for exited state dynamics were generated from the equilibrated ground state molecular dynamics simulations carried out with the AMBER 12 software package 62-64 using the General Amber Force Field (GAFF). 65,66 During these simulations, a time step of 1 fs was used and temperature was equilibrated, employing a Langevin thermostat with a friction coefficient γ of 2.0 ps −1 . After minimization, each molecule was heated to 300 K during 1 ps. Then, 1 µs of NPT molecular dynamics simulations were performed. These simulations provide a good initial conformational sampling for the EHR and AIMC-TDDB simulations. The resulting configurations were used to compute the absorption spectra shown in Fig. 1. One hundred snapshots, equi-spaced on time, were then collected as initial geometries and momenta for the EHR and AIMC-TDDB simulations. While the ground state conformational sampling was obtained using GAFF, nonadiabatic excited-state EHR and AIMC-TDDB simulations were The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp performed using the semiempirical AM1 Hamiltonian, which has a slightly different ground-state potential energy surface compared to that in the GAFF approach. Therefore, the collected set of initial structures extracted from the GAFF needs to be further relaxed on the semiempirical AM1 Hamiltonian before photoexcitation.
In order to do that, a short time Langevin molecular dynamics (20 fs) using the ground state semiempirical AM1 Hamiltonian have been performed on each collected snapshot. The total energy and temperature have been monitored in order to confirm the final equilibration. Initial conditions for the electronic part of the wave function were selected according to the spatial localization of the TD for different adiabatic states. Figure 2 shows the TD localization in real space for each molecular system in its corresponding minimum energy configuration for the ground state. Since we are interested in the analysis of the different non-radiative relaxation pathways, our simulations start from the pure S 4 state, which is mainly localized in the 2-ring chromophore unit. This way our initial conditions are suited to follow the intramolecular energy fluxes from the 2-ring segment to other fragments of the molecule.
The EHR and AIMC-TDDB simulations were performed using a 0.05 fs time step (∆t). The six lowest electronic excited states (S 1 -S 6 ) have been included in our calculations. The same initial conditions have been used for both sets of simulations. For the AIMC-TDDB simulations, the CS widths [Eq. (1)] were chosen according to the average tested parameters reported by Thompson et al. 67 Each EHR trajectory was about 10 times more computationally expensive than the corresponding trajectory in the Surface Hopping (SH) simulations. Despite both types of simulations being performed using the same level of "on-the-fly" calculation for adiabatic energies, gradients, and NACRs, the EHR simulations require a shorter time step (i.e., ∆t = 0.05 fs compared to ∆t = 0.1 fs for SH simulations), and the simultaneous calculation of gradients for all the excited states was necessary for EHR. The shorter time step arises as a requirement for energy conservation due to the spiky behavior of the nonadiabatic contribution to the EHR force [second term of Eq. (2)], which is absent in the SH simulations. The AIMC-TDDB simulations require about four times more of computational cost than the EHR analog, depending on the frequency of cloning events. Notably, the additional EHR trajectories required for AIMC-TDDB simulations do not compromise parallelism of simulations.

It is well established that PPE dendrimers composed of linear polyphenylene ethynylene (PPE) units show negligible electron exchange (and charge transfer processes) across meta-substitutions.
Each chromophore unit retains its own electrons, and the molecular system may be analyzed in terms of the Frenkel exciton Hamiltonian. [68][69][70] Therefore, their optical excitations can be considered as the sum of contributions from separate 2-ring, three-ring, and fourring linear PPE units. 71 In the following discussion, "delocalization" of wave function refers to Frenkel-exciton-like delocalization appearing due to superposition of the electronic wavefunctions from the individual chromophores in the aggregate. Figure 1 shows the absorption spectra for our three dendrimer model building blocks, that is, 234PPE, 243PPE, and 324PPE. The three molecular systems are composed by two-ring, three-ring, and four-ring linear PPE chromophore units linked via meta substitution, differing only by the order they are combined. Therefore, equivalent contributions of each individual chromophore units lead to similar absorption spectra despite that every system underwent its own unique ground state sampling protocol. Figure 2 confirms that optical excitations in our three model systems are essentially localized on single linear PPE chromophore units. S 1 and S 2 states are mainly localized on the four-ring and three-ring units, respectively. In the case of 234PPE, higher energy S 3 and S 4 excited-states are also localized on the four-ring and two-ring units, respectively. S 3 has a node in the middle of the four-ring unit. We can rationalize this state using a simplified model considering an exciton as a "quasiparticle" in a box. Here S 3 is a second excited state of an exciton in a box with size corresponding to the four-ring linear PPE unit. However, S 3 and S 4 are delocalized between the tworing and four-ring units in 243PPE and 324PPE systems. In these two systems, both chromophore units linked by meta substitution are neighbouring each other. The electronic coupling between these linear units mixes near-degenerate electronic states leading to the observed wave function delocalization.
Our excited state dynamics simulations start by an instantaneous excitation of all systems to the S 4 state and use of the AIMC-TDDB method. On the time scale of less than 30 fs, this state undergoes an ultrafast electronic energy relaxation to the lower excited states, being more pronounced in 234PPE than 243PPE and 324PPE. This can be seen in Fig. 3(II) where the evolution in time of the average populations for the different excited states is displayed. In 234PPE, S 3 and S 4 states present coherent in-phase oscillations that are out-of-phase with respect to S 2 and S 1 oscillations. As it has been pointed out previously, all systems are initially excited to the S 4 state. All of them undergo an ultrafast S 4 → S 3 energy transfer within the first ∼3-4 fs after photoexcitation. After that, only 234PPE experiences a subsequent significant loss of the accumulated S 3 population from ∼0.32 to ∼0.2 that takes place concomitantly with a further energy transfer from S 4 to the other states. Both are not observed in the other two systems. Therefore, the 234PPE arrangement of chromophore units seems to favour the constructive electronic coherence effect compared to the other arrangements, inducing a slightly faster relaxation process at earlier times. After ∼30 fs, these oscillations are damped and both S 4 and S 3 experience slower relaxation while population on S 1 increases. S 2 states seem to participate more actively in the electronic relaxation of 234PPE than 243PPE and 324PPE, in agreement with the previously reported through-bond sequential transfer via S 4 /S 3 → S 2 → S 1 mechanism. 23 Despite the differences in the initial relaxation dynamics, at later times, the increase in the S 1 population with time and the relaxation rates to the S 1 state are remarkably similar in all three molecular systems. Fast initial relaxation from S 4 , localized at the two-ring unit, to S 1 localized further down the dendrimeric chain may be important for lossless energy funnelling in dendrimers as it can prevent the energy dissipation to the bath in contact with the photoexcited S 4 state.
In order to further analyze the impact that different combinations of chromophore units have on the exciton intramolecular redistribution after photoexcitation, Fig. 4 compares the evolution in time of the fraction of TD on each unit for our three model systems. As photoexcitation populates S 4 , the initial exciton spatial localization corresponds to the ones obtained by Frank-Condon excitation to that state using previously thermally equilibrated structures on the ground state. In the case of 234PPE, we have shown (see Fig. 2) that, at the minimum of ground state, S 4 is mainly localized on the tworing unit. Nevertheless, thermal fluctuations introduce additional delocalizations (see Fig. S1 in supplementary material) providing the initial averaged values shown in Fig. 4. Equivalent delocalizations are also experienced by 243PPE and 324PPE molecular systems (see Figs. S2 and S3). 234PPE arrangement has the most effective two-ring → four-ring energy transfer. Nevertheless, the final average exciton distribution among the units tends to be similar for all three systems.
At this point, it is interesting to analyze the evolution in time of the exciton localization within S 3 and S 4 states. In Fig. 2, we have seen that, at the minimum of 234PPE ground state, S 3 and S 4 are localized on the four-ring and two-ring segments, respectively. By contrast, both states are delocalized between the two-ring and fourring units for 243PPE and 324PPE. In order to monitor the extent of (de)localization of the excitation in these states, we further analyze the chromophore-unit participation number calculated for two-and four-ring units using Eq. (15).
Values of WK ≈ 1 indicate a complete localization of the TD corresponding to state K over a single chromophore unit, while values of WK ≈ 2 correspond to the TD fully delocalized between the two-ring and four-ring units. Figure 5 shows the variation of WK values for S 3 and S 4 states for all three systems. We observe that, in the case of 243PPE and 324PPE, a significant degree of exciton delocalization persists throughout the entire simulation. By contrast, the exciton is initially localized in both states for 234PPE. Its degree of delocalization gradually increases during the electronic relaxation process.
As we have previously mentioned, high energy states (S 3 /S 4 ) of 243PPE and 324PPE have a node in the middle of the four-ring unit (see Fig. 2). This is not the case for S 1 states. Therefore, the fraction of TD localized on the triple bond in the middle of the four-ring unit can elucidate the apparent lack of intramolecular exciton redistribution observed for 243PPE and 324PPE (see Fig. 4). Figure 6(a) shows the evolution in time of this value, revealing an equivalent behavior in the three molecular systems. This result, in addition to the evolution in time of the fraction of TD localized in the triple bond of the two-ring unit [ Fig. 6(b)], indicates a common intramolecular exciton redistribution for the three systems.
Previous studies, performed using the NAESMD method, 28 have shown that intramolecular vibrational energy redistribution of PPE dendrimers takes place concomitantly with the electronic energy redistribution. [22][23][24][25] Nuclear motions in the direction of the stretching of ethynylene bonds (C≡≡C) represent the main contributions to the non-adiabatic coupling between states (NACR). 72 Therefore, the ethynylene bond lengths spatially distributed across the system are a convenient descriptor for monitoring the intramolecular vibrational energy redistribution. Figures 6(c) and 6(d) show the evolution in time of the average over all simulations of C≡≡C bond lengths localized in the middle of the four-ring unit and in the two-ring unit, respectively. We observe that the C≡≡C bond of the two-ring unit in 234PPE is initially more excited than in the other systems. Gradually, owing to the exciton redistribution after photoexcitation, the C≡≡C bond of the four-ring unit becomes excited in all systems.
The different energy transfer pathways between units can be identified using the transition density flux analysis described in Sec. II D. Figure 7 shows the time evolution of the TD accumulated fluxes (FXY (t)) for each molecular system. In the case of 234PPE, an efficient through-space direct transfer two-ring → four-ring dominates the relaxation process. Besides, a through-bond sequential transfer two-ring → three-ring → four-ring is also observed. This agrees with previous surface hopping studies performed on similar PPE dendrimer fragments. 24 During the first ∼100 fs, the accumulated three-ring → four-ring flux is negative, indicating that actually  the transfer takes place in the opposite three-ring ← four-ring direction. During that time, an effective two-ring → three-ring energy transfer takes place. After 100 fs, a three-ring → four-ring energy transfer is achieved, while the two-ring → three-ring is interrupted. This is in agreement with a sequential mechanism in which an initial two-ring → three-ring energy transfer is followed by a subsequent three-ring → four-ring transfer. The inter-unit energy transfer pathways in 243PPE and 324PPE are different from the ones observed in 234PPE. In both cases, equivalent effective two-ring → three-ring and two-ring → fourring fluxes are observed. It is worth to note that, in previous articles, 72,73 we have shown that vibrations involving the common phenyl ring between two-ring and three-ring units contribute to the NACR vector and, therefore, actively participate in the energy transfer between these two chromophore units. That is, while the two-ring → three-ring transfer is a through-space energy transfer for 243PPE, it becomes a through-bond transfer for 324PPE. Besides, an initial effective flux three-ring ← four-ring is observed in the three systems, followed by a subsequent reverse flux three-ring → four-ring energy transfer.
Finally, it is interesting to analyze the modeling results calculated with the AIMC-TDDB used in this paper in comparison with that obtained with the EHR and NAESMD simulations used in the previous studies. Figure 8 compares the evolution in time of the average population on S 1 computed using both types of simulations. ∼13% faster relaxation time scales can be observed in the MCE simulations of the three molecular systems. As it can be seen in Fig. 8, despite the systematic faster relaxation times observed in the AIMC-TDDB simulations with respect to the EHR simulations of our three molecular systems, observed differences are within the standard deviation of our results. Nevertheless, Mann-Whitney statistical tests with p-values <10 −3 have been obtained between the EHR and AIMC-TDDB distributions of values at 300 fs, indicating the statistical significance of the differences between the AIMC-TDDB and EHR results for the three molecular systems considered. The AIMC-TDDB allows the system to account for more relaxation pathways compared to the EHR, enhancing those through which the relaxation process is faster. Cloning employed by the AIMC-TDDB allows us to evaluate the splitting of the wave packet between different potential energy surfaces better. Although in our previous work this effect was less important owing to a single-path relaxation process in a simple molecule, 40 here it results in noticeable change in the populations and faster relaxation rates. Previous NAESMD simulations, 25,72 based on the surface hopping approach, have shown an even faster relaxation for the 234PPE. While the nuclei move along an average mean field potential energy surface in the EHR modeling, they evolve on a single state during the NAESMD simulations. The AIMC-TDDB method results in between these extremes in the sense that nuclei evolve on a mean field potential like in the EHR method, but whenever the mean field potential energy surface is no longer a good representation of a particular excited state, the cloning procedure ensures that the dynamics recovers the right contribution of this state. As a result, the relaxation rates obtained with the AIMC-TDDB lie in between those obtained by EHR and NAESMD techniques. EHR and NAESMD are known to underestimate and overestimate the exact result, respectively, and this can be considered as a further argument in favour of the accuracy of AIMC calculations. Also it was shown previously 37,39 that in model systems a good convergence of multiconfigurational Ehrenfest calculations to the exact benchmark can be achieved, which also supports the results of this work.

IV. CONCLUSIONS
Here we study an intramolecular energy transfer in light harvesting dendrimers using the recently developed Ab Initio Multiple Cloning-Time Dependent Diabatic Basis (AIMC-TDDB) implementation of the Multiconfigurational Ehrenfest (MCE) approach. We compare the photoinduced nonadiabatic molecular dynamics simulations performed on a set of different combinations of dendrimer building blocks composed of two-, three-, and four-ring linear PPE oligomers linked by meta-substitutions, namely, 234PPE, 243PPE, and 324PPE molecular systems (shown at the insets in Fig. 1). These simulations address the effect of different spatial assemblies of chromophoric units on the energy transfer efficiency via investigating complex interplay of multiple bifurcating pathways involved in the non-radiative relaxation.
At the longer time scale of hundreds of femtoseconds [see Fig. 3(I)], similar relaxation rates are observed for all three arrangements of the units in the three systems studied. The final intramolecular exciton distribution is similar, and the average exciton delocalization patterns are the same. However, on the short time scale [ Fig. 3(II)], the dynamics of electronic energy relaxation depends on the molecular geometry. The nature of high energy excited states (S 3 /S 4 ) changes across the systems. In 234PPE, S 3 and S 4 are initially localized on the four-ring and two-ring units, respectively. By contrast, these states are delocalized between the two-ring The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp and four-ring segments in 243PPE and 324PPE, reflecting electronic coupling among the different fragments that compose these units. This difference in the exciton localization for S 3 and S 4 states affects the efficiency of two-ring → four-ring energy transfers. Sequential geometric arrangement of the linear units of different lengths in 234PPE leads to the most effective energy transfer network. 243PPE and 324PPE do not seem to experience significant intramolecular exciton redistributions after photoexcitation, gradually leading to the final excitonic state similar for the three systems considered. The electronic energy relaxation in 234PPE involves both a through-space direct transfer two-ring → four-ring and a throughbond sequential transfer two-ring → three-ring → four-ring with the major contribution coming from the former mechanism. Here both processes are spatially unidirectional, which contributes to the overall efficiency. Such sequential transfer two-ring → three-ring → four-ring is not observed neither in 243PPE nor in 324PPE, and the through-space direct transfer two-ring → four-ring seems to be less effective for the latter two systems.
Our simulations exemplify the effect of wave function bifurcations and the importance of cloning for accurate description of photoinduced dynamics of multichromophoric light harvesting dendrimers. The use of the AIMC-TDDB method enhances the electronic energy relaxation by ∼13% compared to the mean-field Ehrenfest results in all systems. Previous studies performed on 23PPE did not reveal such an effect. 40 Therefore, given a manageable rise in underlining computational cost, this work encourages the use of a more accurate AIMC-TDDB approach to other multichromophoric molecular systems featuring rich non-radiative dynamics toward quantitative modeling of experimental data.

SUPPLEMENTARY MATERIAL
See supplementary material for the distribution of the transition density fraction of each excited state localized on the three chromophore units evaluated for the ensemble of the initial ground state conformational sampling.