Quantum dynamic simulation of triplet formation in an effective model of Y6 (BTP-4F)

Abstract: We construct a five-state model for photoexcitation in Y6 (BTP-4F) dimers, and then solve the non-adiabtic dynamics using the Hierarchical Equations of Motion (HEOM) method. We find that triplets are populated mainly via a transiently excited \textit{intermolecular} charge-transfer singlet to triplet Frenkel exciton route; this route is not available to the monomer. Analysis of one-particle transition density matrices suggests that the charge-transfer states are spatially distinct to the Frenkel exciton states, indicating that the large spin-orbit-coupling for this transition is due to it being permitted by an associated change in orbital character. Aggregation in Y6 therefore directly enables fast and high-yield intersystem crossing. We selenise our model dimers, significantly enhancing spin-orbit-coupling, which then accelerates this charge-transfer mediated route. Looking forwards to simulations on larger aggregates, we show that, though Marcus theory gives qualitatively correct dynamics, the long-time yields are incorrect due to it missing quantum recurrences. Instead, we show that the recently developed memory-kernel projector\cite{Gestsson2025-ez} method can produce semi-classical rates directly from the HEOM equations which lead to quantitatively correct dynamics and yields.

• The paper uses a five-state effective Hamiltonian and non-Markovian HEOM dynamics to model six Y6 dimers, finding ultrafast singlet Frenkel-exciton to charge-transfer conversion within about 100 fs.
• The paper shows that triplets form predominantly through charge-transfer singlet to triplet Frenkel-exciton intersystem crossing, with rates of 10⁵–10⁷ s⁻¹ in Y6 that increase nearly tenfold after selenisation.
• The paper finds that packing strongly controls triplet generation, while Nakajima–Zwanzig-projected HEOM rates reproduce full quantum dynamics more accurately than Marcus theory and offer a lower-cost route to larger simulations.

Overview and motivation

This paper constructs a five-state effective Hamiltonian for Y6 (BTP-4F) dimers—extracted from the OHUBUR crystal structure (2606.12221)—and propagates the non-adiabatic excited-state dynamics with the Hierarchical Equations of Motion (HEOM) method. The central question is the origin of triplet excitations in Y6, which bear directly on two open problems: the debated role of triplets as loss channels in organic photovoltaics (OPVs), and the mechanism behind the record solid-state triplet–triplet annihilation upconversion observed by Izawa and Hiramoto in Y6/rubrene bilayers. Two competing hypotheses exist: triplets formed via charge separation and recombination at the interface (spin-statistical 3/4 yield), or directly via intersystem crossing (ISC) from the photoexcited singlet. The distinction matters practically: if ISC dominates, a simple design rule follows—heavy-atom substitution increases spin-orbit coupling (SOC) as $\mathcal{O}(Z^2)$, accelerating ISC as $\mathcal{O}(Z^4)$.

The paper also takes an explicit position against an alternative account: Souza et al. attribute monomer ISC to extreme excited-state torsion of the key intramolecular torsion to 90°, generating large SOC. The authors argue such distortion is sterically implausible in the close $\pi$-$\pi$ stacked solid state, and note that Souza et al.'s Franck–Condon analysis (with Herzberg–Teller corrections) yields rates insufficient to explain experiment.

Effective Hamiltonian construction

Six "contact pair" dimers (D1–D6) from Giannini et al., with intermolecular distances below 3.5 Å, are projected onto five effective states: $S_0$, singlet/triplet Frenkel exciton ($S_{FE}$, $T_{FE}$), and singlet/triplet charge transfer ($S_{CT}$, $T_{CT}$). Electronic structure is B3LYP/6-31G(d,p) TDDFT with the Tamm–Dancoff approximation, a combination the authors previously validated against optimally tuned range-separated hybrid calculations for these dimers; they emphasize precision and consistent approximation across matrix elements over absolute accuracy. Effective site energies are arithmetic means of contributing TDDFT states; couplings are pooled in quadrature to respect Fermi's golden rule. FE–CT couplings use counterpoise-corrected HOMO/LUMO projection integrals ($V_h$, $\mathcal{O}(Z^4)$0); SOC matrix elements (SOCMEs) come from ORCA's SHARK analytic $\mathcal{O}(Z^4)$1 integrals under spin-orbit mean-field theory.

Two features of the resulting parameter sets stand out:

1. State ordering: singlet CT states lie below singlet FE states (~1.66–1.82 eV vs ~1.87–1.98 eV), while triplet FE states lie below triplet CT states (~1.46–1.51 eV vs ~1.67–1.85 eV). This inversion is essential to the dynamics.
2. SOC selection rule: SOCMEs between states of different character (e.g., CT(S)↔FE(T), 0.12–0.36 cm⁻¹ in Y6) exceed those within the same character (0.01–0.14 cm⁻¹). The authors interpret this as a generalization of El-Sayed's rule: the change in orbital location between CT and FE states changes orbital symmetry, permitting spin-orbit coupling.

By contrast, the Y6 monomer has negligible $\mathcal{O}(Z^4)$2–$\mathcal{O}(Z^4)$3 SOCMEs (0.01 and 0.07 cm⁻¹). To rationalize this, the paper introduces a scalar spatial-overlap metric on fragment-resolved one-particle transition density matrices,

$\mathcal{O}(Z^4)$4

which equals 1 for identical electron-hole distributions and 0 for spatially orthogonal ones. Monomer $\mathcal{O}(Z^4)$5 overlaps its triplets strongly ($\mathcal{O}(Z^4)$6), while dimer CT↔FE pairs have near-zero overlap ($\mathcal{O}(Z^4)$7) yet carry the largest SOCMEs. This correlation supports the claim that aggregation-induced spatial distinctness of CT and FE states enables fast ISC—a route structurally unavailable to the monomer.

For selenised analogues (Y6Se, outer sulfurs replaced by selenium with no structural relaxation, mirroring T9SBN-F chemistry), relativistic ZORA tests showed that apparent energy shifts stem mainly from disabling TDA rather than scalar relativity, so non-relativistic B3LYP/TDA was retained for all matrix elements.

HEOM dynamics

Dynamics are propagated with HierarchicalEOM.jl at hierarchy depth $\mathcal{O}(Z^4)$8 (convergence verified), coupling each state to two baths: a Drude-Lorentz bath ($\mathcal{O}(Z^4)$9 eV, $\pi$0 eV ≈ $\pi$1) for slow intermolecular modes, and an under-damped Brownian oscillator ($\pi$2 eV, empirically chosen $\pi$3 eV, peak at 0.160 eV) for high-frequency intramolecular modes. Phenomenological Lindblad terms model fluorescence decay to $\pi$4 at $\pi$5 s⁻¹ (FE) and $\pi$6 s⁻¹ (CT), taken from transient absorption data; supplementary results show the 1 ns triplet yield is insensitive to their inclusion.

Three robust findings emerge across all six dimers:

• Ultrafast internal conversion from FE(S) to CT(S) within ~100 fs.
• Triplet formation predominantly via CT(S) → FE(T) ISC, with a minor direct FE(S) → CT(T) secondary channel. This pathway exists only because intermolecular CT states exist—i.e., it is a property of the aggregate, not the molecule.
• Nanosecond-scale triplet populations consistent with transient absorption and electroabsorption measurements on films.

Triplet formation rates extracted from steady-state fluxes span roughly $\pi$7–$\pi$8 s⁻¹ for Y6 and rise by nearly an order of magnitude upon selenisation (e.g., D1: $\pi$9 from $\pi$0 to $\pi$1 s⁻¹), exactly the scaling expected if ISC controls the rate. The rates depend strongly on packing: the D1 H-aggregate, with the largest transfer integrals, forms triplets about an order of magnitude faster than other motifs. The implication drawn is that micromorphology—crystallinity, deposition protocol—should strongly modulate both upconversion efficiency and triplet losses in OPV devices, which may connect to Y6's known processing sensitivity.

Rate models: Marcus theory versus projected HEOM rates

Because full HEOM scales factorially and cannot extend to device-scale aggregates, the paper benchmarks two routes to classical master-equation rates. Marcus theory with an empirical reorganisation energy of 0.5 eV reproduces the initial ~50 fs FE(S)→CT(S) transfer and long-time detailed balance, but—as a weak-coupling, Markovian, perturbative theory—it misses the early-time FE/CT hybridisation around 100 fs and therefore fails to capture the minor FE(S)→CT(T) channel, yielding incorrect long-time yields.

The alternative applies the recent Nakajima–Zwanzig projection method of Gestsson et al., extracting time-independent effective rates $\pi$2 directly from the HEOM memory kernel, where $\pi$3 is the resolvent of the irrelevant-space Liouvillian evaluated at $\pi$4. These rates reproduce the full HEOM dynamics quantitatively—including coherent beating integrated out but quantum recurrences retained—at a fraction of the computational cost. Notably, convergence of the projection required increasing the Drude-Lorentz Matsubara truncation from $\pi$5 to $\pi$6, and the authors recommend retaining Marcus rates as cross-validation where HEOM itself is infeasible.

Limitations and open questions

Several caveats qualify the conclusions. The five-state coarse-graining pools multiple microscopic states into single effective FE and CT manifolds; the authors concede that the lack of clear dimer-to-dimer trends in the CT(S)–CT(T) and FE(T)–CT(T) exchange fluxes likely reflects this conflation, motivating larger Hamiltonians requiring GPU-accelerated or tensor-network HEOM. The reorganisation energies and the under-damped bath parameters are partly empirical, and Marcus modelling uses a 0.5 eV $\pi$7 acknowledged to exceed what the electronic structure supports—an implicit correction for operating outside the perturbative regime. The El-Sayed generalization rests on the assumption that macroscopic spatial orthogonality implies orbital-symmetry permission for SOC, which is asserted rather than derived. Finally, the mechanism makes a falsifiable prediction left to experiment: triplet formation should be much faster in the solid state than in solution (where the intermolecular CT route is absent), whereas the vibrational-torsion mechanism of Souza et al. predicts the opposite trend; comparative TAS/ESA measurements across phases would discriminate the two accounts.

Conclusion

This work establishes, through exact non-Markovian quantum dynamics on a first-principles-parametrised dimer model, that Y6 forms triplets primarily by ISC from a transiently populated intermolecular singlet charge-transfer state to triplet Frenkel excitons—a route enabled by aggregation, governed by a spatial-overlap-based selection rule analogous to El-Sayed's rule, and accelerated an order of magnitude by heavy-atom substitution. It further demonstrates that Nakajima–Zwanzig-projected HEOM rates deliver quantitatively correct semi-classical dynamics where Marcus theory fails, providing a practical bridge toward simulating exciton diffusion and full Y6–rubrene upconversion interfaces.