Sebastian
Fernandez-Alberti
a,
Dmitry V.
Makhov
*b,
Sergei
Tretiak
c and
Dmitrii V.
Shalashilin
b
aUniversidad Nacional de Quilmes, Roque Saenz Peña 352, B1876BXD Bernal, Argentina
bSchool of Chemistry, University of Leeds, Leeds LS2 9JT, UK. E-mail: d.makhov@leeds.ac.uk
cCenter for Nonlinear Studies (CNLS), and Center for Integrated Nanotechnologies (CINT), Los Alamos National Laboratory, Los Alamos, NM 87545, USA
First published on 10th March 2016
Photoinduced dynamics of electronic and vibrational unidirectional energy transfer between meta-linked building blocks in a phenylene ethynylene dendrimer is simulated using a multiconfigurational Ehrenfest in time-dependent diabatic basis (MCE-TDDB) method, a new variant of the MCE approach developed by us for dynamics involving multiple electronic states with numerous abrupt crossings. Excited-state energies, gradients and non-adiabatic coupling terms needed for dynamics simulation are calculated on-the-fly using the Collective Electron Oscillator (CEO) approach. A comparative analysis of our results obtained using MCE-TDDB, the conventional Ehrenfest method and the surface-hopping approach with and without decoherence corrections is presented.
The time evolution of electronic excitation in organic conjugated materials, such as dendrimers, is determined by non-adiabatic dynamics involving multiple coupled electronic excited states. Following photoexcitation, the multiple photoinduced pathways to electronic and vibrational relaxation involve energy and/or charge transfer, internal conversion, and transition density localization/delocalization. Different theoretical approaches can deal with such photophysical processes.
A number of models, such as Pariser–Parr–Pople,8 Frenkel exciton,9 and Frenkel–Holstein,10 have been used to describe quasi one-dimensional conjugated systems. The real-time Density Matrix Renormalization Group (DMRG)11 technique is a powerful numerical method for studying these model systems. In particular, PPE dendritic models can be described using Frenkel exciton Hamiltonians12 with Coulombic coupling between excitations localized on the linear units.13
In atomistic molecular dynamics simulations, surface hopping has been proved to be an efficient method to accurately deal with non-adiabatic processes in a large variety of systems. Nevertheless, it presents certain shortcomings. First, electronic decoherence effects are not accounted for in the standard surface hopping algorithm resulting in a number of computational inconsistencies that have been addressed by implementing a large variety of alternative decoherence corrections.14 The inadequate treatment of decoherence can, in many cases, influence the accuracy of the simulated results. Moreover, there is no clear and direct strategy to improve surface hopping methods in order to introduce nuclear quantum effects such as delocalization, interferences and tunneling.
Among the most accurate fully quantum methods is time-dependent Hartree (MCTDH),15 which has been recently applied to simulate exciton dynamics in extended conjugated polymers.16 However, this approach, in its original version, requires the parameterization of potential energy surfaces. A useful alternative is represented by a variety of methods relying on trajectory-guided Gaussian basis functions (TBF). The main advantage of these methods is that they can deal with direct or ab initio molecular dynamics17 where excited state energies, gradients, and non-adiabatic coupling terms are evaluated on the fly, simultaneously with the nuclear propagation. In a variational multi-configuration Gaussian (vMCG)18 approach, the trajectories for all TBFs are calculated together from the variational principle. vMCG is similar in spirit to MCTDH and it is among the most elegant methods of direct dynamics. However, because of its complexity, on the fly vMCG so far has been implemented only for a small number of basis Gaussians.
Another group of methods17d,19 applies TBFs moving along classical trajectories. Despite that these methods are still fully quantum because these trajectories are employed only for propagating the basis, while the evolution of their amplitudes and, thus, of the total wave-function is determined by the time-dependent Schrödinger equation. The trajectories here are calculated independently allowing for a large and well-sampled basis. The methods involving branching, such as multiple spawning17d,19a or multiple cloning19g algorithms can hardly be used for dendrimers because the large number of electronic states and intersections between them would lead to the uncontrolled growth of the number of trajectories. But, on the other hand, the same makes dendrimers and other large conjugated molecules a nearly ideal object for Ehrenfest based approaches,19b–f,20 where nuclear motion is defined by the average field, and amplitudes for all electronic states are propagated together with nuclear coordinates. In such systems, a wave packet undergoes frequent transition between many coupled electronic states, and as a result its motions can be frequently well described by the Ehrenfest trajectory even on long time scales.
In recent years we have been developing a group of methods based on multicofigurational Ehrenfest ansatz for the total wave function.19d–g The central idea is to use Gaussian wavepackets guided by Ehrenfest average of the electronic states. Several variations of this idea have been presented. Our original version of ab initio MCE was using adiabatic electronic states generated by the electronic structure. Here we suggest a modification of the approach, which now has new features: each Gaussian trajectory carries its own electronic basis, which is adiabatic only in the center of the Gaussian. This is similar to the approach described in the appendix in ref. 19g, except that now we accurately take into account the overlaps between electronic states belonging to different trajectories, instead of assuming them to be Kronecker's delta as in ref. 19g. This modification makes the method applicable to large conjugated molecules where electronic states can change significantly on the length-scale of the Gaussian width. In particular, the wave function can change instantly at trivial unavoided crossings between two electronic states localized at remote parts of a large molecule.21 To distinguish the new approach from the previous work we term it as multiconfigurational Ehrenfest in time-dependent diabatic basis (MCE-TDDB). We would like to emphasize that “time-dependent diabatic basis” used here should not be confused with “diabatic basis” commonly used in many studies: the trajectories here are still calculated in adiabatic basis, and the new representation affects only matrix elements between different TBFs.
We use the proposed MCE-TDDB method to simulate the excited state dynamics of the system composed of two- and three-ring linear polyphenylene ethynylene (PPE) chromophore units linked through meta-substitution as shown in Fig. 1(a). These chromophore units correspond to the building blocks of more complex phenylethynylperylene-terminated dendrimers, such as the nanostar.5b,22 The meta-branching localizes the electronic excited states within each linear fragment, and the non-adiabatic coupling between linear fragments are responsible for quantum transitions resulting in the two-ring → three ring unidirectional electronic energy flow. Thus, this system is a good model for analyzing intramolecular electronic and vibrational energy transfer between chromophore units subsequent to the photoexcitation of dendrimer macromolecules.
We run the MCE-TDDB calculations in the basis of a swarm of interacting coherent state trains,19g,23 which is the most advanced basis for a given number of calculated trajectories. The excited state energies,13b,24 gradients25 and non-adiabatic coupling terms25a,26 needed to run the dynamics are calculated on the fly using the Collective Electron Oscillator (CEO) method.5a,27 The results are compared with those given by the canonical Ehrenfest approach (EHR) and by the NA-ESMD (non-adiabatic excited-state molecular dynamics)25a,28 method. The latter technique utilizes surface hopping guided by Tully's fewest-switches (FSSH) algorithm29 and was specifically developed to describe photoinduced dynamics in large organic conjugated molecules. In particular, the NA-ESMD framework has been extensively applied before to study the intramolecular energy transfer in PPE dendrimers.30
The paper is organized as follows. In Sections II and III we describe the MCE-TDDB method and computational details of simulations performed on our model system. In Section IV we present and discuss our results. Conclusions are given in Section V.
![]() | (1) |
![]() | (2) |
![]() | (3) |
The basis |ϕ(n)I〉for the electronic part of wave-function can be specified in different ways. In the ab initio multiple spawning method and in the original version of the ab initio MCE approach, the electronic part is represented in the basis of adiabatic eigenstates:
|ϕ(n)I〉 = |ϕI(r;R)〉. | (4) |
|ϕ(n)I〉 = |ϕI(r;![]() | (5) |
However, this is not a case for the system under consideration, and a more rigorous approach must be developed. We will use electronic wave functions (5) and will take into account that the overlaps between the electronic eigenstates belonging to different Gaussians 〈ϕ(n)I|ϕ(m)J〉 can be very far from Kronecker's δIJ, even when these Gaussians are sufficiently close to each other and nuclear parts have a significant overlap.
The electronic overlaps 〈ϕ(n)I|ϕ(m)J〉 can be either calculated directly or propagated together with the basis. The latter approach is more convenient in this case, and the following equation for the time-dependence of the overlap integrals is used:
![]() | (6) |
As before, the motion of the centers of Gaussians is determined by the usual set of Hamilton's equations
![]() | (7) |
![]() | (8) |
![]() | (9) |
Also similar to our previous studies,19f,g the evolution of the Ehrenfest amplitudes a(n)I is determined by the equation
![]() | (10) |
The above eqn (7)–(10) describe the evolution of the basis functions (2). Let us now consider the evolution of the whole wave-function (1) in the above-defined basis (2). Substituting eqn (1) into the time-dependent Schrodinger equation and taking into account non-orthogonality of the basis, we get
![]() | (11) |
![]() | (12) |
![]() | (13) |
![]() | (14) |
![]() | (15) |
![]() | (16) |
The term 〈Ψm(t)|d/dt|Ψn(t)〉 in eqn (11), which originates from the time dependence of the basis, can be written as
![]() | (17) |
![]() | (18) |
![]() | (19) |
![]() | (20) |
![]() | (21) |
![]() | (22) |
![]() | (23) |
![]() | (24) |
![]() | (25) |
![]() | (26) |
A similar approach can be applied for calculations of any other electronic properties, such as electronic transition density matrices reflecting spatial extent of the respective many-body electronic wave-functions. If the quantity is described by a quantum operator , then
![]() | (27) |
![]() | (28) |
N(m)IK = 〈ϕ(m)I|![]() | (29) |
![]() | ||
Fig. 2 The basis sampling. (A) – a single Ehrenfest configuration; (B) – a swarm of trajectories; (C) a swarm of coherent state trains. The principle basis function is shown in orange. |
![]() | (30) |
The FSSH algorithm implemented in NA-ESMD simulations is designed to ensure the internal consistency between the so-called classical populations and the quantum populations. The classical populations are defined as a fraction of trajectories evolving on each potential energy surface:
PClI(t) = NI(t)/Ntraj, | (31) |
PQI(t) = 〈|aI(t)|2〉Ntraj, | (32) |
Empirical corrections can be introduced to account for electronic decoherence.14 Following extensive numerical tests,14a we adopt the instantaneous decoherence (ID) approach, performed by resetting the quantum amplitude of the current state to unity after every attempted hop (regardless of whether hops are allowed or forbidden), thus removing the coherence of the quantum coefficients. This approach is based on the assumption that wavepackets traveling on different surfaces should immediately separate in the phase space and evolve independently. It has been previously shown14a that ID improves FSSH internal consistency between classical and quantum systems while at the same time provides results that do not depend on external parameters and maintains physical relevance.
The accuracy of excited-state energies in the AM1/CIS approximation has been benchmarked for a large variety of conjugated molecules,32 such as PPV33 and polyfluorenes.34 In the particular case of our model PPE dendritic molecule (shown at the Fig. 1(a)), vertical excitation energies for the two lowest excited states calculated at the AM1/CIS level for an AM1-optimized structure have been compared with the results obtained using TDDFT (B3LYP/6-31G* model) as well as with the experimental results,30,35 and a reasonable agreement has been achieved. More details about the benchmark calculations by SCS-AD(2), TDHF, TDDFT methods can be found in ref. 36, where the mechanism of intramolecular energy transfer suggested by NA-ESMD simulations at the AM1/CIS level30 has been confirmed.
The evaluation of the accuracy of NACME calculations is a more difficult task for any electronic structure method, because direct experimental measurements of NACME are absent. NACME can only be evaluated indirectly from electron–phonon coupling elements and Frank–Condon overlaps. Such an evaluation suggests reasonable accuracy of NACME calculated by AM1/CIS for excited state potential energy surfaces in conjugated hydrocarbons.33 Moreover, semi-quantitative direct comparisons of non-radiative relaxation rates calculated using our NA-ESMD approach and experimental time-resolved transient spectroscopy data37 supports the conclusions of ref. 33. Also, both AM1/CIS calculations38 and TDDFT36 show that the major contribution to d12(n) comes from the ethynylene stretches, which represent the reactive coordinates of the energy transfer process.
Six singlet electronic states (S1–S6) and their corresponding non-adiabatic couplings have been included in the simulations. The excited-state trajectories of 150 fs duration are propagated at constant energy in order to obtain a reasonable statistics. Time steps of 0.5 fs and 0.02 fs have been used in the ground and excited state dynamics, respectively. The Verlet finite difference algorithm39 is used to integrate the nuclear equations of motion, while the electronic coefficients are propagated using the Runge–Kutta–Verner fifth- and sixth-order method based on the code designed by Hull, Enright, and Jackson.40
In order to generate the initial conditions for excited-state dynamics, 1 ns of ground state molecular dynamics at 300 K was first performed using a Langevin friction coefficient γ of 2.0 ps−1. Snapshots of nuclei positions and momenta (conformational phase space) have been collected and used as initial conditions for subsequent photoexcitation dynamics modeling. The excited-state trajectories have been started from these initial configurations by instantaneously promoting the system to the state I selected according to a Frank–Condon window:
gI(R) = fI![]() | (33) |
In order to make the comparison consistent, the same set of initial positions, momenta and excited-state populations have been used both for principle trajectories in Ehrenfest dynamics and for NA-ESMD surface hopping. Equations of motions were integrated with the same time step, and at the same level of electronic structure theory.
We performed the MCE-TDDB simulations for 100 trajectory swarms. Each swarm consists of 10 trajectories: the principle one and 9 satellites. This procedure is similar to the one used previously in model calculation,19d where first the initial wave packet created in the excited state after absorbing a photon was split into “bits” and then each “bit” was propagated separately on using a swarm of Ehrenfest configurations.
The initial conditions for the principle trajectories were generated as described above. The initial conditions for the satellites were taken23 very close to those for the principle trajectory with a random shifts ΔR and ΔP generated according to the Gaussian distribution
![]() | (34) |
![]() | ||
Fig. 3 The decay of the absolute value of the average overlap between the principle and satellite trajectories in a swarm for MCE-TDDB calculations. |
For each trajectory in a swarm, we use a coherent state train basis of 11 basis functions. Thus each “bit” was propagated on the basis of 10 × 11 = 110 TBFs in total. A half of additional TBFs in the train are placed before and a half – after the main TBF with a time-shift of 0.6 fs, which provide the overlap of ∼0.7 between the nearest basis functions in the train. In order to avoid the necessity of backward propagation, the train basis was applied only after 3 fs of dynamics, when there was enough trajectory length to accommodate the tail half of the train; for first 3 fs, calculations were run in the basis of just 10 TBFs.
All the trajectories were calculated independently. Subsequently, the amplitudes cn(t) for all configurations in the basis were propagated using eqn (11) in a “post-processing” procedure.19g This makes the numerical procedure particularly efficient because a very large number of trajectories can be calculated in parallel. Although the total CPU time can reach years, most of it is spent on calculating trajectories using our ab initio on the fly procedure described by eqn (7)–(10) above. Quantum equations for the amplitudes use approximations (16) and (20), which allow to calculate coupling matrix elements in eqn (11) on the basis of the information already available from the calculation of the Ehrenfest trajectories. This “post-processing” technique is the central advantage of the group of methods based on our MCE methodology developed by us. The sampling techniques used in this work were tested on model systems, and it was demonstrated that the MCE can reproduce accurately fully quantum benchmark calculations.42 Thus, MCE-TDDB and related methods19f,g represent a very efficient ab initio direct dynamics methodology.
The results of ab initio MCE-TDDB calculations were compared with those given by the standard Ehrenfest (EHR) and the NA-ESMD surface hopping methods. The EHR results were obtained from the same trajectories as in the MCE-TDDB, but only the principle trajectory in each swarm was taken into account, and the populations and transition densities for these principle trajectories were simply averaged over the ensemble. The NA-ESMD surface hopping calculations28,29b,c were performed for a set of 100 excited-state trajectories with the same initial conditions as for principle trajectories in MCE-TDDB calculations. Simulations were also performed using the instantaneous decoherence corrections. More details about these kinds of simulations can be found elsewhere.14a
The initial excited states are populated using a Frank–Condon window defined by a Gaussian shaped laser centered at 348 nm, which corresponds to the maximum of the contribution of the S2 state. Fig. 4 presents the average populations of the lowest four electronic states as a function of time calculated using the NA-ESMD (a, b), MCE-TDDB (c), EHR (d) methods. For NA-ESMD, both classical PClI(t) and quantum PQI(t) populations are presented. One can see that there is a significant difference between these two populations, which indicates that the so-called internal consistency requirement is violated here.
The dependencies for MCE-TDDB and EHR are very close to each other, except that MCE-TDDB, which relies on swarms and trains basis sets, produces slightly smother curves. On the one hand, this may be an indication that the Ehrenfest approach works really well for large complex systems with a singular photoproduct, such as PPE dendrimers. On the other hand, we know from our previous experience with model systems that using swarm/train basis in combination with trajectory branching (e.g., via cloning procedure19g) can further improve the convergence of the results and provide even quantitative agreement with multidimensional benchmark calculations. Cloning could introduce additional decoherence allowing basis functions for different electronic states to run away from each other. Such approaches, however, are very expensive for on the fly use, especially in the case of large molecules with many coupled electronic states, as a rapidly growing number of trajectories would require a huge amount of electronic structure calculations. We do not expect branching to change the result dramatically, but working on implementing even larger basis sets to prove that. Presented here MCE-TDDB technique provides a promising starting point for future developments of methods incorporating fully-quantum vibronic dynamics.
All population dependencies presented in Fig. 4 are similar qualitatively. In particular, one can see a noticeable oscillatory interchange of electronic populations between S1 and S2 states during the first 20 fs of dynamics after photoexcitation. The amplitude of the interchange between these states is revealed to be larger for Ehrenfest approaches with respect to surface hopping simulations. In all cases the ultrafast non-radiative decay of the electronic population on S2 takes place, and by the time of 150 fs the S1 population becomes significantly higher compared to its initial value. One can see that at the end of our excited-state dynamics, the S1 population given by the MCE-TDDB and EHR methods lies in between the classical and quantum populations PClI and PQI given by the NA-ESMD.
In order to gain insight into the origin of the main differences between surface hopping and Ehrenfest approaches, we have performed the NA-ESMD simulations using decoherence corrections. The importance of incorporating decoherence schemes in the NA-ESMD simulations applied to large extended conjugated molecules have been previously demonstrated.14a Decoherence has a large effect on the accuracy of the simulated dynamics both in terms of time scales for nonradiative relaxation and in promoting internal consistency of the FSSH algorithm. Fig. 5 shows the time-dependence of the average quantum populations of different electronic states for these NA-ESMD simulations using decoherence corrections. We observe that the oscillatory interchange of electronic populations between S1 and S2 (see Fig. 4) is damped in these simulations. This is expected because periodic collapses of the electronic wave-function to a single state eliminate the interference between these two states. After that, the population of the S1 state increases faster compared to quantum populations given by other approaches.
![]() | ||
Fig. 5 Average populations of different electronic states as a function of time obtained from NA-ESMD (quantum populations) using decoherence corrections. The colours refer to the same states as in Fig. 4. |
On the one hand, decoherence is not accounted for in the standard FSSH simulations, leading to internal inconsistencies of the method.14a,c,n On the other hand, instantaneous decoherence assumes that the divergent wavepackets are immediately decoupled. According to the results shown in Fig. 4 and 5, the Ehrenfest approach data fall in between quantum populations given by the NA-ESMD simulations with and without decoherence corrections. The possible explanation is that the Ehrenfest force depends on the amplitudes of all electronic states taken into account, including sufficiently low populated higher-energy states. Then, continuous population changes due to numerous intersections between a fairly large number of electronic states lead to weak quasi-random fluctuations of Ehrenfest force, which provide natural decoherence. This does not happen in the surface hopping approach where the force is usually determined by a single occupied adiabatic state, and the hops are not frequent and mostly reflect only significant population changes.
The CEO approach43 calculates the many-body eigenstates by diagonalizing the respective tetradic CIS operator in the molecular orbital (MO) basis. These CIS eigenstates are represented in the form of Nocc × Nvirt matrices (where occ and virt refer to occupied and virtual MOs), which can be transferred to the original atomic orbital (AO) basis. The resulting quantities are frequently denoted as transition density matrices (or electronic normal modes) and can be formally written as (ρI)mn ≡ 〈ϕg|cm+cn|ϕI〉, where |ϕg〉 and |ϕI〉 are the many-electron wave-functions for the ground and Ith excited states, and cm+ and cn are the electron creation and annihilation operators with indices n and m referring to the AO basis functions. Notably, owing to a unitary transformation between AO and MO basis and an orthogonality of all AO in the AM1/CIS approach, the normalization condition 24b,25a is maintained.44 Matrix elements of the transition density matrices are subject to a straightforward interpretation: the diagonal elements of (ρI)nn reflect45 light-induced changes in an orbital n due to promotion of the molecule to the Ith excited state defining its transition dipole moment (the expectation value of the dipole operator on the transition density matrix) and, subsequently, an oscillator strength fI. In contrast, the off-diagonal elements of (ρI)nm reflect probability amplitudes of electron transfer between orbitals n and m, describing the optically induced electronic coherences and charge-transfer phenomena. Overall, the transition density matrices directly and compactly reflect the properties of many-body wave-functions even beyond the CIS framework. Subsequently, these are convenient quantities to analyze photoinduced spatial electronic/excitonic delocalization, charge-transfer, etc.13b
In order to obtain the fraction of the transition density localized either on the two-ring or the three-ring linear PPE units, we sum up the atomic contributions belonging to each fragment as
![]() | (35) |
Fig. 6 shows localization of the electronic transition density of the four lowest electronic states at the minimum of the ground state potential energy surface. One can see that S1 and S2 states are mainly localized on the three-ring and two-ring units, respectively. This means that the dendritic PPE molecule can be represented as an ensemble of weakly coupled linear chromophore units,30c and the absorption spectra shown in Fig. 1(b) are approximately a sum of the absorption of three-ring and two-ring units:2,27a,b,46 the peak at 376 nm corresponds to the absorption of three-ring unit (S1 state), and the shoulder around 348 nm is due to the additional contribution of the two-ring unit (S2 state). Therefore, the excitation at our laser wavelength of 348 nm corresponds to an initial localization of the electronic transition density almost equally distributed between the two-ring and three-ring units. The electronic transition density of higher energy excited-states (S3 and S4) is qualitatively different from that of S1 and S2 states: it is delocalized over both two- and three-ring fragments for S3, and strongly localized for the S4 state.
![]() | ||
Fig. 6 Initial localization of the electronic transition densities for the four lowest excited states. |
Nuclear motion mixes electronic states providing ultrafast energy transfer between the photoactive units. Fig. 7 compares the time-dependencies of the average transition density fractions localized on the two- and three-ring linear PPE units ((ρI)X-ring2) given by MCE-TDDB and NA-ESMD approaches. Note that the NA-ESMD density is calculated by averaging those of the “current” electronic state for all trajectories in the ensemble, similar to eqn (31) for classical electronic populations. The time-dependencies of transition densities are consistent with our results for the electronic populations (see Fig. 4(a) and (c)): the NA-ESMD approach shows more efficient energy migration toward the three-ring unit, as well as more efficient relaxation from S2 to the lowest state S1 associated with this transfer.
![]() | ||
Fig. 7 Fraction of the transition density localized in the two- and three-ring linear PPE units as a function of time obtained from NA-ESMD (solid) and MCE-TDDB (dashed) simulations. |
It was shown previously that the electronic energy transfer is concomitant to the intramolecular vibrational energy redistribution, and the nuclear motions in the direction of the ethynylene bonds is playing a critical role in the process.30bFig. 8 compares the distributions of the two-ring ethynylene bond length in the ground state and at S2–S1 crossing, defined as hop points for surface-hopping trajectories and as points with high NACMEs for Ehrenfest trajectories. Both types of calculations show that, in agreement with previous results,30b S2–S1 crossing is associated with ethynylene bond stretching in the two-ring part of the molecule.
In Fig. 9, we plot the ethynylene bond lengths for both two- and three-ring units as a function of time. Accompanying the initial excitation of S2, this bond of the two-ring unit is initially excited. After that, its excitation relaxes together with the ultrafast decay of the S2 electronic population (Fig. 4). In contrast, the ethynylene bond of the three-ring unit becomes excited increasing its average amplitude of motion.
![]() | ||
Fig. 9 Length of particular ethynylene bonds (labeled X in the insets) as a function of time obtained from averaged NA-ESMD (solid) and Ehrenfest (dashed) simulations. |
The difference in the amplitude of motion for the ethynylene bonds given by different approaches is consistent with our other results. At early times of the dynamics (∼20 fs), the amplitude of motion of the ethynylene bond of the two-ring unit is larger for Ehrenfest than that for surface hopping simulation. This time corresponds to the time at which the oscillatory interchange of electronic populations between S1 and S2 states takes place (see Fig. 4). At longer times, we observe the reduction of the average vibrational amplitudes. This decay is faster for the Ehrenfest trajectories where the nuclear motion is coupled with amplitudes a(n)I for all electronic states providing better vibrational energy redistribution through continuous partial electronic transitions between many states.
Fig. 8a indicates that the ethynylene bond of the two ring unit behaves differently for NA-ESMD and the Ehrenfest approach. The non-adiabatic contribution to the Ehrenfest force n (eqn (9)), not present in NA-ESMD simulations, introduces additional vibrational coherence at earlier times of simulations in the regions with large d12(
n) values. At later times, the Ehrenfest force dependence on the continuous low population changes among higher-energy states seems to result in quasi-random fluctuations that lead to vibrational decoherence. The same additional vibrational coherence is responsible for higher amplitude of the oscillation of the electronic population, which both Ehrenfest and MCE-TDDB reveal as shown in Fig. 4. This is consistent with the results of ref. 30, where a strong coupling between S1 and S2 states was shown to exist only during first 15–20 fs after photoexcitation.
The proposed ab initio MCE-TDDB approach is very efficient for large conjugated molecules, such as dendrimers, where the dynamics involving many continuously interacting electronic states, is well described by mean-field Ehrenfest trajectories. Although even larger basis sets with branching trajectories can be introduced, such calculations would be extremely expensive and unlikely to produce significantly different results.
As an application, here we simulate the exited state dynamics for a dendrimer molecule composed of two- and three-ring PPE chromophore units following laser pulse excitation at 348 nm. Our ab initio MCE-TDDB results are compared with those of surface hoping. Our calculations show that for dynamics involving a large number of electronic states, the approaches based on Ehrenfest trajectories naturally account for decoherence phenomena. Thus, our MCE-TDDB calculations validate decoherence corrections, previously introduced to the surface hoping technique in an ad hoc fashion.
Altogether, this makes Ehrenfest-based techniques, such as presented here MCE-TDDB, a method of choice for on the fly simulations of the excitation dynamics in large photoactive molecules.
This journal is © the Owner Societies 2016 |