Exciton transport in molecular organic semiconductors boosted by transient quantum delocalization

Abstract

Designing molecular materials with very large exciton diffusion lengths would remove some of the intrinsic limitations of present-day organic optoelectronic devices. Yet, the nature of excitons in these materials is still not sufficiently well understood. Here we present Frenkel exciton surface hopping, an efficient method to propagate excitons through truly nano-scale materials by solving the time-dependent Schrödinger equation coupled to nuclear motion. We find a clear correlation between diffusion constant and quantum delocalization of the exciton. In materials featuring some of the highest diffusion lengths to date, e.g. the non-fullerene acceptor Y6, the exciton propagates via a transient delocalization mechanism, reminiscent to what was recently proposed for charge transport. Yet, the extent of delocalization is rather modest, even in Y6, and found to be limited by the relatively large exciton reorganization energy. On this basis we chart out a path for rationally improving exciton transport in organic optoelectronic materials.

Fig. 1. Organic molecular crystals investigated in this work.

The chemical structures of the single molecules (left panels), natural transition orbitals (NTOs) for the first singlet excited state (middle panels), and the experimental crystal structure (right panels) are shown for anthracene (ANT, a), α-sexithiophene (a6T, b), perylenetetracarboxylic diimide (PDI, c), dicyanovinyl-capped S,N-heteropentacene (DCVSN5, d), and (Y6, e). Directions of strongest excitonic couplings are indicated in the crystal structures and labeled according to their crystallographic direction (e.g., P_a, P_b, T, see also Table 1). In Y6, P_a2 is eclipsed by P_a1. The unit cell axes a, b, c are shown in red, green, and blue. The percentages indicated denote the contribution to the S1 singlet excited states due to the NTOs shown. The NTOs labeled “Hole” and “Electron” resemble very closely the HOMO and LUMO of each molecule.

Fig. 2. Accuracy of fast excitonic coupling estimators.

a Correlation between the reference total excitonic coupling V (in absolute value), obtained using MS-FED-FCD (Eq. (7) in Supplementary Method 1), and V^Coulomb (Eq. (11)) for neighboring molecules along the π-stacking direction (P_b of ANT and a6T and P_a, P_a1, P_a2 for PDI, DCVSN5, and Y6 respectively). Excitonic couplings were calculated for configurations taken from MD trajectories generated for supercells of the crystals at 300 K. V, couplings correspond to the off-diagonal block ${\mathbb{H}}_{\mathrm{XT1,XT2}}^{″\prime }$ in MS-FED-FCD related to the lowest localized excited states in molecules k and l, respectively and the diabatization was performed using the first 20 adiabatic states, which proved to be sufficient for convergence (see Supplementary Note 3). b Correlation between V^Coulomb (Eq. (11)) and V_TrESP (Eq. (12)). Mean relative unsigned error (%) is defined as MRUE = (∑_n(∣y_calc − y_ref∣/y_ref))/n. In both panels the MRUE is ca. 7%. This small deviation found for the systems investigated justifies the use of V_TrESP in FE-SH simulations.

Fig. 3. Hopping and transient delocalization mechanisms for exciton transport.

Panels a, b show time series along a single representative FE-SH trajectory for ANT and DCVSN5, respectively. Top panels: energy of the active potential energy surface on which the nuclear dynamics is run in FE-SH simulations, E_a (red lines), referred to the ground state. Middle panels: index of the active exciton band state, a (blue lines, index a = 0 corresponds to the ground state which is the bottom of the excitonic band). Bottom panels: IPR(t), Eq. (8), of the excitonic wavefunction Ψ(t), Eq. (3), (black lines) and average IPR of Ψ(t) over the swarm of FE-SH trajectories (horizontal dashed gray lines). Note the correlation between the intra-band excitations, i.e., eigenstate index a, excitation energy, and IPR(t) of Ψ(t). Panel c depicts a typical nearest-neighbor hopping event of the exciton in ANT and panel d a representative transient delocalization event of the exciton in DCVSN5. Notice the much larger spatial displacement from the initial position (indicated with orange circles) in the latter case. The excitonic wavefunctions in panels (c, d), Ψ(t), Eq. (3), are represented by isosurfaces of the transition density on each molecule scaled by the expansion coefficient u_l. The transition density is approximated, for visualization purposes, by the conjugated product of HOMO and LUMO orbitals.

Fig. 4. Exciton diffusion constant and impact of transient delocalization.

a Exciton diffusion constant (D) as a function of the average IPR, (〈IPR〉), for the systems investigated (see Fig. 1 for their chemical and crystal structures). The diffusion constants along the crystallographic directions a and b obtained from FE-SH simulations are represented by filled circles and filled triangles pointing upwards, respectively. The average of D and 〈IPR〉 values over different system sizes is plotted. Numerical values for specific system sizes are shown in Supplementary Fig. 9. Error bars indicate the corresponding standard deviations. Diffusion constants obtained from experimental data, as described in Table 2, are indicated by empty symbols and are placed at the computed IPR value since experimental estimates for exciton delocalization are not available. The crystallographic direction for the experimental diffusion constant for Y6 and a6T was not reported. b IPR distribution obtained from 600 FE-SH trajectories for a rod-like DCVSN5 nano-crystal comprised of 300 molecules. c The full MSD of the exciton wavefunction for DCVSN5 is shown with a blue line and is calculated according to Eq. (7). The MSD obtained when transient delocalization events are excluded is shown with a thin red line. In the latter case, the IPR threshold was set to IPR_thr = 〈IPR〉 + 1. Notably, a threefold increase in the diffusion constant is obtained when transient delocalization events are retained. d Accumulated percentage contribution of wavefunction delocalization events to D (i.e., (D_thr/D) × 100), plotted as a function of the IPR threshold. The green region indicates the contribution of nearest-neighbor hopping to D, whereas the yellow region indicates the contribution due to transient delocalization.

Fig. 5. Comparison charge vs exciton transport.

Correlation between a the diffusion constant (D, Eq. (6)) and b the time average inverse participation ratio (IPR, Eq. (8)) of charges and excitons against the activation barrier (Eq. (1)) for different systems, respectively. Data from present exciton transport simulations are depicted in blue (ANT (1), a6T (2), PDI (3), DCVSN5 (4), and Y6 (5)). Data for charge transport simulations for the 2D conductive layers of OSs are taken from ref. and depicted in green (ANT (1), naphthalene (6), perylene (7), pMSB (8), rubrene (RUB), and pentacene (PEN)). For 2D simulations, two data points are shown for diffusion constant and activation barrier, one for the direction of highest coupling (unprimed number) and one for the direction of second highest coupling (primed number). PEN-h⁺-T₁ and RUB-h⁺-a denote hole transport along a crystallographic direction indicated. For comparison, 1D charge transport simulations (taken from ref. ) are indicated by empty black circles and the values reported along the chain direction. The best linear fits to the data are indicated by dashed gray lines.