Charge and spin delocalisation in photoinduced polarons of polymer donors and non-fullerene acceptors for organic photovoltaics – a multi-frequency pulse EPR study

Abstract

Organic photovoltaics rely on efficient charge separation and transport, processes facilitated by charge delocalisation across the π-conjugated backbone of donor and acceptor molecules. By probing the interactions between the unpaired electron spin associated with photoinduced charged states and magnetic nuclei in their molecular environment, electron paramagnetic resonance (EPR) spectroscopy enables precise experimental quantification of spin and charge delocalisation. We first characterise the EPR spectral signatures of the positive and negative polarons on polymer donors and non-fullerene acceptors in the PBDB-T:ITIC and PM6:Y6 blends, as well as the corresponding fullerene-based blends PBDB-T:PC₆₁BM and PM6:PC₆₁BM, by multi-frequency EPR spectroscopy. Reliable separation of overlapping donor and acceptor signatures is enabled by EDNMR-induced EPR spectroscopy exploiting unique nuclear hyperfine couplings in the non-fullerene acceptors ITIC and Y6. Then, by combining the measurement of electron–nuclear hyperfine couplings by ENDOR with DFT modelling and a regularised least-squares fitting approach, we quantify the extent of spin and charge delocalisation. The experimental results reveal intramolecular delocalisation of the positive polarons on the PBDB-T and PM6 donors across approximately 6 nm. Delocalisation of the negative polarons, on the other hand, depends on the nature of the acceptor: for ITIC, the electron spin is found to be localised on a single molecule, whereas for Y6, contributions from spins localised on a single molecule, as well as spins delocalised over two adjacent molecules in different π–π stacked configurations, are required to explain the experimental ENDOR and HYSCORE data. Validation of computational predictions by experimental results is shown to be crucial for the accurate estimation of charge delocalisation, and therefore conclusions on its relevance in determining device efficiency in organic photovoltaics.

---

1 Introduction

Organic electronics exploit the high degree of synthetic tuneability of organic molecules and polymers for the design of materials tailored to applications ranging from light-emitting diodes and field-effect transistors to photodiodes and photovoltaic cells. In the field of organic photovoltaics, research efforts are focused on the design of matched pairs of organic donor polymers and small molecule acceptors enabling efficient light absorption, charge separation and charge transport. Significant advances over the last decade focused on the design of acceptor molecules, replacing fullerenes with acceptors based on fused-ring structures with alternating electron-donating and electron-accepting moieties. These non-fullerene acceptors have now led to power conversion efficiencies of around 20% for research cells^1–4 and to commercial applications in building-integrated photovoltaics and solar-powered consumer electronics.⁵ Despite this remarkable progress, an in-depth understanding of fundamental processes at the basis of energy conversion and of correlations between molecular (and supramolecular) structure and device performance is still missing, with non-fullerene acceptor blends challenging previously established theories on the photovoltaic mechanism.^6–8

State-of-the-art organic photovoltaics are based on an intermixed blend of a donor polymer and a small molecule acceptor, forming a so-called bulk heterojunction. Light absorption leads to the formation of a singlet exciton, followed by charge transfer at the donor:acceptor interface and separation of the hole and electron into free charge carriers in the donor and acceptor domains, respectively. These charges can be extracted at the electrodes to generate electricity. A complete understanding of the photovoltaic mechanism crucially relies on information on the local electronic and molecular structure of the states involved in the photovoltaic energy conversion. Charge delocalisation, in particular, has been suggested as a crucial factor determining efficient charge separation and charge transport based on several theoretical studies,^9–15 but direct experimental evidence correlating increased delocalisation with improved performance remains limited. By exploiting selective detection of the electron spin associated with the charge on the donor and acceptor, electron paramagnetic resonance (EPR) can provide valuable insights into the different charged states involved in the photovoltaic mechanism – from the charge-transfer state at the donor:acceptor interface to the separated charges.^16–18 Characterisation of hyperfine couplings between the electron spin and nuclear spins in its molecular environment provides a way to measure the spatial distribution of the electron spin density, allowing experimental quantification of the extent of spin and charge delocalisation and validation of quantum chemical approaches used for its estimation. Experimental EPR studies of spin delocalisation in materials for organic electronics, including photovoltaics, have so far focused on investigating delocalisation along the conjugated backbone of organic semiconducting donor polymers,^19–25 however, delocalisation across multiple acceptor molecules is now being proposed as one of the reasons for the exceptional performance of the most promising non-fullerene acceptors.^26–28

In this work, we use pulse EPR spectroscopy across X-, Q-, and W-band frequencies to characterise photogenerated charged states (polarons) on donor and acceptor molecules in two photovoltaic blends based on non-fullerene acceptors, PBDB-T:ITIC and PM6:Y6 (Fig. 1), and the two corresponding blends based on fullerene acceptors, PBDB-T:PC₆₁BM and PM6:PC₆₁BM. ITIC was the first non-fullerene acceptor achieving a comparable performance to fullerene acceptors, which were previously almost exclusively used.²⁹ This inspired further research into non-fullerene acceptors, with Y6 and related Y-series acceptors now serving as key components in the currently best-performing blends.³⁰

Fig. 1 Molecular structures of the investigated donor (PBDB-T and PM6) and acceptor (PC₆₁BM, ITIC, Y6) molecules. The orientations of the principal g-axes with respect to the molecular structures are also shown (see SI Section S9 for details).

We first identify the EPR spectral signatures of the photogenerated charged states on PBDB-T, PM6, ITIC and Y6 (Fig. 1). In contrast to photovoltaic blends based on the fullerene PC₆₁BM acceptor, where distinct differences in g-values allow clear separation of contributions from polarons on donor and acceptor molecules,^20,21,23,31 in the case of non-fullerene acceptors, the similarities in molecular building blocks, and therefore similar g-values for donors and acceptors, lead to extensive spectral overlap. To disentangle these overlapping contributions, we exploit the presence of ¹⁴N nuclear spins on the ITIC and Y6 acceptor molecules, but not on the donor polymers, to extract the acceptor EPR spectrum through ELDOR-detected NMR (EDNMR) experiments performed as a function of magnetic field, following the approach previously described by Van Landeghem et al.^31,32

We then probe the extent of charge delocalisation in donor and acceptor polarons through a combination of EPR hyperfine spectroscopy and density functional theory (DFT) calculations. The electron–nuclear hyperfine couplings are directly related to the spin density distribution, and therefore allow experimental quantification of spin delocalisation.^19–25 DFT, though commonly used to estimate charge delocalisation in extended π-conjugated systems, such as donors and acceptors for organic photovoltaics, tends to overestimate delocalisation due to self-interaction errors.³³ This has prompted the development of range-separated hybrid functionals, designed to more accurately model long-range interactions, which are now widely used for computational modelling of organic semiconductors.^34–39 We compare the hyperfine couplings predicted by DFT for a range of molecular models and both standard and range-separated hybrid functionals with experimental orientation-selective ¹H and ¹⁹F Electron Nuclear DOuble Resonance (ENDOR) spectra recorded for both donor and acceptor molecules and ¹⁴N HYperfine Sublevel CORrElation (HYSCORE) spectra recorded for ITIC˙⁻ and Y6˙⁻. This allows us to validate computational models and correlate hyperfine couplings with spin density distributions, mapping out the extent of spin and charge delocalisation for PBDB-T˙⁺, PM6˙⁺, ITIC˙⁻ and Y6˙⁻. In order to account for the intrinsic heterogeneity of donor:acceptor blends, which leads to a distribution of spin delocalisation lengths not reflected in a single DFT-optimised geometry, we employ a regularised least-squares approach guided by the results of DFT calculations, as introduced by Pribitzer et al.⁴⁰ This combined experimental and theoretical approach reveals delocalisation across approximately 6 nm (three repeat units) on the π-conjugated backbone of the donor polymer for PBDB-T˙⁺ and PM6˙⁺, while the extent of delocalisation varies for polarons on the acceptors: in ITIC, the electron spin is localised on a single molecule; in Y6, both electron spins localised on single molecules and spins delocalised across π–π stacked dimers are present. Our findings offer detailed insights into the local electronic structure of polarons, enabling further investigations of the role of delocalisation in determining charge separation and charge transport efficiencies, ultimately affecting the power conversion efficiency of a given donor:acceptor blend.

2 Results and discussion

2.1 EPR spectra of donor and acceptor polarons

The EPR spectral signatures of positive polarons on the donor polymers PBDB-T and PM6 and of negative polarons on the non-fullerene acceptors ITIC and Y6 are first identified based on EPR experiments performed across multiple frequencies, providing initial information on electronic structure through the characteristic g-values.

2.1.1 Multi-frequency EPR characterisation. The EPR spectra recorded for the photoinduced charged states in the donor:acceptor blends PBDB-T:PC₆₁BM, PBDB-T:ITIC, PM6:PC₆₁BM and PM6:Y6 at X-, Q- and W-band under steady-state illumination are shown in Fig. 2.

Fig. 2 Light-induced echo-detected X-, Q- and W-band EPR spectra recorded at 20 K for the PBDB-T:PC₆₁BM, PBDB-T:ITIC, PM6:PC₆₁BM and PM6:Y6 blends under steady-state illumination. The experimental spectra shown were obtained after subtraction of a spectrum recorded prior to illumination (see Section S5 in the SI) and are compared to simulations with the parameters listed in Table 1 (dashed black lines). The individual contributions of polarons on the donor and acceptor molecules are shown as shaded light/dark blue and red/orange areas, respectively, and the field positions corresponding to the principal g-values are indicated. Further experimental details are described in Section S2.1 in the SI.

The fullerene-containing blends clearly show two contributions already at X-band, which are fully separated at higher frequencies. The high-field signal is assigned to the negatively-charged polaron on the PC₆₁BM molecule based on extensive literature on the characterisation of PC₆₁BM˙⁻.^20,41–44 The low-field contribution is therefore assigned to the positively-charged polaron on the donor molecules, PBDB-T˙⁺ and PM6˙⁺, respectively. The shape of the PC₆₁BM˙⁻ contribution is consistent across the PBDB-T:PC₆₁BM and PM6:PC₆₁BM blends and reflects an axial g-anisotropy with g_x ≈ g_y > g_z and significant g_z-strain, as evident in particular from the Q- and W-band measurements. The shape of the PBDB-T˙⁺ and PM6˙⁺ spectra changes considerably with microwave frequency, suggesting that at X-band the lineshape is mostly determined by unresolved electron–nuclear hyperfine couplings between the electron spin and the large number of protons on the polymer backbone, while at higher frequencies the spectra are dominated by the g-anisotropy characteristic of a near-axial system with g_x > g_y, g_z.

The clear separation of the donor and acceptor spectral signatures for PBDB-T:PC₆₁BM and PM6:PC₆₁BM at Q- and W-band allows straightforward determination of the g-values by global fitting of the spectra recorded at the three microwave frequencies to ensure an accurate and consistent set of EPR parameters for the donor and acceptor polarons (Fig. 2). Additionally, the PC₆₁BM˙⁻ contribution in both blends was fit using the same set of parameters. The extracted g-values are listed in Table 1, where the principal g-value pointing along the out-of-plane direction of the aromatic π-system is denoted as g_z, and g_y and g_x are ordered based on magnitude. The experimentally observed line broadening is due to a combination of unresolved hyperfine couplings and distributions of g-values arising from variations in the molecular environment, which were accounted for by an isotropic convolutional broadening with a Gaussian of indicated peak-to-peak linewidth and anisotropic g-strain as implemented in EasySpin,⁴⁵ respectively. The substantial broadening of the g_z feature in the EPR spectrum of PC₆₁BM˙⁻ and, to a lesser extent, of the g_x feature in the spectra of the polarons on the donor polymers, PBDB-T˙⁺ and PM6˙⁺ (particularly prominent at W-band, see Fig. 2, right panel) indicates the presence of a distribution of g-values. A common g_z- or g_x-strain was therefore included in the fitting of the spectra across all microwave frequencies, whereas strains along the remaining two principal g-axes were accounted for within the convolutional broadening, and fit individually for each microwave frequency.

Table 1 Simulation parameters for the light-induced echo-detected EPR spectra for PBDB-T˙⁺, PM6˙⁺, PC₆₁BM˙⁻, ITIC˙⁻ and Y6˙⁻ obtained through global fitting of the experimental X-, Q- and W-band echo-detected EPR spectra and the W-band EDNMR-induced EPR spectra for PBDB-T:ITIC and PM6:Y6. The magnetic field was calibrated using an N@C₆₀g-standard⁴⁷ (see Section S2 in the SI)

Material	g-values^a			g i -strain^a		linewidths^b (mT)
	g x	g y	g z			X	Q	W
The relative donor : acceptor weights used in the simulations are 1 : 0.86 (PBDB-T : PC61BM), 1 : 0.92 (PM6 : PC61BM), 1 : 1.09 (PBDB-T : ITIC) and 1 : 1.07 (PM6 : Y6).a Estimated uncertainties of ±0.0001 for g-values and g-strains.b FWHM for convolutional Gaussian broadening, uncertainties of ±0.02 mT.c Obtained for the blends with PC61BM.d Obtained for the blends with the non-fullerene acceptors ITIC and Y6.e Hyperfine couplings included for X-band simulations, taken from DFT calculations (nuclei with maximum coupling ≥10 MHz, four ^1H nuclei for ITIC˙^− and two ^1H and two ^19F nuclei for Y6˙^−).
PBDB-T˙^+^c	2.0033	2.0023	2.0020	0.0009 (x)		0.26	0.32	0.74
PBDB-T˙^+^d						0.23	0.29	0.63
PM6˙^+^c	2.0033	2.0022	2.0019	0.0008 (x)		0.22	0.24	0.48
PM6˙^+^d						0.23	0.22	0.46
PC61BM˙^−	2.0002	2.0001	1.9990	0.0019 (z)		0.10	0.16	0.32
ITIC˙^−^e	2.0039	2.0031	2.0021	0.0013 (x)	0.0007 (y)	0.19	0.41	0.58
Y6˙^−^e	2.0038	2.0031	2.0021	0.0013 (x)	0.0004 (y)	0.22	0.63	0.90

---

Integration of the separated spectral contributions of donor and acceptor reveal small deviations from a 1 : 1 ratio expected based on formation of a pair of oppositely charged states following each photoinduced charge separation event, with a slight excess of the donor polymer contribution (footnote to Table 1). This has been observed in other donor:acceptor blends previously.^31,46 Differences in relaxation times have been explored as a possible origin, however differences in phase memory time are too small to significantly impact echo signal intensities and effects of different T₁ were excluded through careful selection of the measurement conditions. Therefore, the differing signal intensities of donor and acceptor signals are attributed to differences in decay rates and mechanisms for the two types of polarons.

For the donor:acceptor blends based on the non-fullerene acceptors ITIC and Y6, the contributions of the positive and negative polarons are not spectrally resolved even at W-band. At X-band, the EPR spectra for both blends consist of a single peak with broad low-intensity wings extending, in particular, towards lower fields. At Q- and W-band, the spectral shape evolves to reveal a shoulder on the low-field side of the spectrum. Knowledge of the donor polaron signature from the analysis of the fullerene-containing blends suggests the contribution of the acceptor polaron extends to lower magnetic fields, but with considerable spectral overlap. Assuming a similar molecular environment of the donor polymer in the different blends, an initial estimate of the acceptor EPR spectra is possible, however, precise determination of the principal g-values from this data alone is challenging. To accurately characterise the ITIC˙⁻ and Y6˙⁻ spectral contributions, more sophisticated pulse EPR experiments aimed at separating the signals of the donor and acceptor polarons are required. We adopt a similar approach to Van Landeghem et al.,^31,32 who successfully disentangled the contributions of ITIC˙⁻ in PBDB-T:ITIC and 2,4-diCN-Ph-DTTzTz˙⁻ in MDMO-PPV:2,4-diCN-Ph-DTTzTz by leveraging the presence of unique electron–nuclear hyperfine couplings in one of the components.

2.1.2 Disentangling spectral contributions by EDNMR. The presence of ¹⁴N nuclei in ITIC and Y6, but not PBDB-T and PM6, allows extraction of the acceptor contributions to the EPR spectra of the PBDB-T:ITIC and PM6:Y6 blends by ELDOR-detected NMR-induced EPR spectroscopy.³¹ ELDOR-detected NMR probes the transitions of nuclei coupled to the electron spin, and can be used to encode information on the presence of coupled nuclei into the EPR spectral dimension.

In a typical ELDOR-detected NMR experiment, a high-turning angle (HTA) pulse is used to selectively excite ‘forbidden’ electron–nuclear transitions, while monitoring the spin polarisation across an ‘allowed’ electron spin transition with a two-pulse echo sequence (Fig. S1).^48,49 When the HTA pulse is on-resonance, spin polarisation is transferred from the electron to the nuclear transition, resulting in a reduction in echo intensity. Detection of the echo intensity as a function of the offset between the frequencies of the HTA pulse and the observer sequence (Δν = ν_HTA − ν_obs) yields a nuclear frequency spectrum. Additionally, saturation of the allowed EPR transition by the HTA pulse at Δν = 0 leads to a central blind spot.^48,49

Fig. 3 displays the W-band ELDOR-detected NMR spectra for PBDB-T:ITIC and PM6:Y6, measured as a function of magnetic field. The right panels display the EDNMR spectra extracted at a single field position, which are characterised by the strong central blind spot and by smaller peaks centred at Δν ≈ ±10 MHz, corresponding to the ¹⁴N Larmor frequency at the selected field position. In theory, EDNMR spectra of nuclei weakly coupled to the electron spin exhibit signals centred at the nuclear Larmor frequency and split by the electron–nuclear hyperfine interaction, and by the nuclear quadrupole interaction in the case of I ≥ 1 nuclei such as ¹⁴N.⁴⁸ The observation of a single peak at ±10 MHz without resolved splitting and a full-width half maximum of ca. 4.5 MHz for ITIC˙⁻ and ca. 4.0 MHz for Y6˙⁻ indicates the presence of a distribution of relatively small hyperfine couplings for the four equivalent ¹⁴N nuclei in ITIC³¹ and the three groups of equivalent ¹⁴N nuclei in Y6.

Fig. 3 W-band ELDOR-detected NMR spectra of PBDB-T:ITIC and PM6:Y6 recorded as a function of magnetic field and microwave frequency offset (Δν = ν_HTA − ν_obs) at 20 K. The central region including the contribution from the central blind spot is masked to highlight the contributions centred at the ¹⁴N Larmor frequency. The right panel shows the EDNMR spectrum at a fixed magnetic field (3353.3 mT), while the top panel shows the EPR spectrum measured as the difference between experiments with off- (20 MHz) and on-resonant (10 MHz ≈ ν_¹⁴N) HTA pulses as well as the corresponding simulation with the parameters reported in Table 1 (dashed black lines). The simulations of the donor and acceptor contributions as well as of the full EPR spectrum are also shown as shaded areas.

The ¹⁴N signal in the EDNMR spectra for both PBDB-T:ITIC and PM6:Y6 extends over most of the magnetic field range of the corresponding EPR spectrum without change in spectral shape, indicating the absence of orientation-selection effects. This enables the use of an EDNMR-induced EPR experiment performed with a fixed microwave frequency offset Δν for the extraction of the spectral signature of ITIC˙⁻ and Y6˙⁻, since the signal intensity of the ¹⁴N peak is proportional to the contribution of the nitrogen-nuclei containing species.³¹ The spectra displayed in the top panels of Fig. 3 are obtained as the difference of EPR signals recorded with an off-resonant (Δν = 20 MHz) and a resonant (Δν = 10 MHz) HTA pulse. A similar separation of donor and acceptor contributions can also be achieved by exploiting differences in T₁ relaxation times in relaxation-filtered EPR experiments (SI Section S6).^50,51

Comparison of the extracted acceptor spectral signatures and the simulations of the donor polymer contributions, obtained from the analysis of the EPR spectra of the fullerene-based blends, with the echo-detected EPR spectra of PBDB-T:ITIC and PM6:Y6 reveal extensive overlap of the donor and acceptor signatures across most of the donor spectrum, with the acceptor contribution extending further to lower magnetic fields (Fig. 3, top panels). The g-values of ITIC˙⁻ and Y6˙⁻ obtained from simulations of the EDNMR-induced EPR spectra are reported in Table 1 and reflect orthorhombic symmetry. The g-values extracted for ITIC˙⁻ are in agreement with literature values,³¹ and the g-values determined for Y6˙⁻ are very similar to ITIC˙⁻, as expected based on the similarities in molecular structure and predictions from DFT calculations.

Knowledge of the g-values of ITIC˙⁻ and Y6˙⁻ now allows simulation of the X-, Q- and W-band EPR spectra obtained for the PBDB-T:ITIC and PM6:Y6 blends, as shown in Fig. 2. Accurate simulation of the Q- and W-band acceptor EPR spectra required g_x- and g_y-strain in addition to convolutional broadening accounting for unresolved hyperfine couplings in the acceptor molecules. The PBDB-T˙⁺ and PM6˙⁺ contributions to the donor:non-fullerene acceptor blends can be simulated with the same g-values and g_x-strain as for the fullerene acceptor blends and with only minor changes in line broadening, suggesting a similar molecular environment of the donor molecules in both types of blends. The broad wings observed in the X-band EPR spectra cannot be reproduced by considering g-anisotropy and isotropic line broadening alone, but require the inclusion of hyperfine couplings. Simulations were therefore performed based on information on hyperfine couplings and principal axis orientations obtained from DFT calculations, validated by ENDOR experiments below. For ITIC˙⁻, the simulations shown in Fig. 2 include coupling to two protons on the vinyl linkers and two protons on the outer edges of the central indacenodithienothiophene (IDTT) unit (largest hyperfine couplings of −16 MHz and −10 MHz for H_vinyl and H_IDTT,o, respectively) and, for Y6˙⁻, to the two protons on the vinyl linkers and two fluorines on the terminal 1,1-dicyanomethylene-3-indanone (INCN) groups (largest hyperfine couplings of −15 MHz and −16 MHz for H_vinyl and F_INCN).

2.1.3 g-Values of polarons on donors and acceptors. The g-values determined by simulation provide initial information on the electronic structure of the polarons, if analysed in combination with the results of DFT calculations. The orientations of the principal g-axes with respect to the molecular structures are shown in Fig. 1 (see also Section S9 in the SI). The similarities in the molecular building blocks of the donor polymers, PBDB-T and PM6, and the non-fullerene acceptors, ITIC and Y6, explain the similarities in g-values and therefore the extensive spectral overlap. In all cases, the smallest g-value, g_z, is associated with the out-of-plane axis, in agreement with the expectation that the out-of-plane g-value in aromatic hydrocarbons shows the least deviation from the free electron g-value.^52,53 The g-anisotropy in organic semiconducting molecules including sulfur-containing thiophene units is attributed to spin–orbit interactions involving sulfur and shows the largest deviations from g_e in the plane of the molecule.^54,55 For organic molecules with spin density localised on a sulfur atom, g-values in the range from 2.014 to 2.020 have been reported.⁵⁶ Therefore, the much smaller deviations from g_e observed here for PBDB-T˙⁺ and PM6˙⁺ (g_y ≈ 2.0022, g_x ≈ 2.0033), suggest extensive spin delocalisation and the absence of significant spin density on the sulfur atoms, in agreement with DFT predictions (Fig. S8). Since PBDB-T and PM6 have identical molecular backbones and differ only by replacement of two protons on the thiophene side chains with fluorines, the g-values do not differ significantly within the uncertainty. The larger deviations from g_e and the increased orthorhombic g-anisotropy observed for ITIC˙⁻ and Y6˙⁻ (g_x ≈ 2.0038–2.0039) can be attributed to increased localisation of spin density on the sulfur atoms of the backbone (Fig. S8).

The g-anisotropy and distinct g-values in PC₆₁BM˙⁻ have previously been investigated extensively and are mostly determined by the presence of near-degenerate orbitals.^43,44 The spin density in PC₆₁BM˙⁻ is distributed in a belt around the fullerene cage perpendicular to the side chain (Fig. S8), with the smallest g-value, g_z = 1.9990, in the plane of maximum spin density, and the two near-identical g_x and g_y values along the axis of the phenyl-butyric acid methyl ester side chain and along the second in-plane axis. The significant g_z-strain has been explained in terms of variation of mainly the g_z-value due to deformation of the fullerene cage in low-energy breathing vibrational modes.⁴⁴

2.2 Spin delocalisation probed by ENDOR

The assignment of the EPR spectral signatures for the donor and acceptor molecules now allows us to investigate spin delocalisation on the polymer donors PBDB-T and PM6 and the non-fullerene acceptors ITIC and Y6 using ¹H Davies ENDOR measurements combined with DFT calculations.

Davies ENDOR detects changes in spin polarisation across an EPR transition following excitation of nuclear transitions with a radiofrequency pulse.⁵⁷ Detection of the intensity of an inverted echo as a function of radiofrequency (Fig. S1) yields high-resolution nuclear frequency spectra, where the contribution of each coupled nucleus results in peaks centred at the nuclear Larmor frequency and split by the electron–nuclear hyperfine coupling. Comparison of the experimental ENDOR data with simulations based on DFT calculations provides information on the spin density distribution in the molecule and therefore on the extent of spin and charge delocalisation.

2.2.1 Davies ENDOR on donors and acceptors. In the first instance, we recorded Davies ENDOR spectra at W-band at a series of field positions across the donor EPR spectrum for the PBDB-T:PC₆₁BM and PM6:PC₆₁BM blends (Fig. 4, left). Due to the limited excitation bandwidth of pulses compared to the width of the W-band EPR spectra, at each field position, a subset of molecules with a specific orientation with respect to the static magnetic field vector is excited, allowing determination of the electron–nuclear hyperfine coupling along this direction. These orientation selection effects provide more detailed information on the hyperfine interactions and on the relative orientation of the corresponding principal axes with respect to the g-frame.

Fig. 4 W-band ¹H Davies ENDOR spectra for the fullerene-based donor:acceptor blends PBDB-T:PC₆₁BM and PM6:PC₆₁BM (left) and for the donor:non-fullerene acceptor blends PBDB-T:ITIC and PM6:Y6 (right) recorded at different magnetic field positions at a temperature of 20 K (see section S2.4 in the SI for additional experimental details). The bottom panels show ENDOR spectra recorded at different magnetic field positions, overlaid on the top to highlight field-dependent differences and offset below to more clearly visualise the individual spectral shapes. The top panels show the field positions with respect to the experimental W-band EPR spectra and the simulations of the donor and acceptor contributions.

The clear separation of the EPR spectra of PBDB-T˙⁺ or PM6˙⁺ and PC₆₁BM˙⁻ at W-band allows characterisation of the ¹H hyperfine interactions in the donor polymers. The ¹H ENDOR spectra of PBDB-T˙⁺ cover a range of about 5 MHz and are relatively featureless except for two shoulders centred around the ¹H Larmor frequency (ν_¹H ≈ 142.7 MHz, Fig. 4). A small degree of orientation selectivity is evident from the distance between maxima of the shoulders decreasing from about 1.8 MHz on the high-field side of the EPR spectrum (g_z, g_y) to about 1.1 MHz at lower magnetic fields (g_x). The ENDOR spectra for PM6˙⁺, differing from PBDB-T˙⁺ only by the presence of a fluorine nucleus on the thiophene side chains of the benzodithiophene unit, exhibit very similar shapes, and only a slightly decreased overall width (see Fig. S4 in the SI), indicating a similar extent of spin delocalisation in the two donor polymers. The Davies ENDOR spectra recorded for PM6˙⁺ show hints of signals at the fluorine Larmor frequency (ca. 8.4 MHz lower than ν_¹H), which can be observed more clearly by exploiting the higher sensitivity for small hyperfine couplings of Mims ENDOR (Fig. S6 in the SI). Only narrow ¹⁹F ENDOR peaks with a width <0.5 MHz were detected, indicating weak long-range couplings and therefore the absence of significant spin density on the thiophene side chains.

ENDOR spectra recorded at the field position corresponding to PC₆₁BM˙⁻ only show a narrow peak (full-width half maximum of ca. 0.5 MHz, Fig. 4, grey) attributed to very weakly coupled protons on the side chain and in the surrounding molecular environment, as expected for a spin density distribution localised almost exclusively on the fullerene cage (Fig. S12).^43,44

The ¹H Davies ENDOR spectra recorded for the blends with the non-fullerene acceptors ITIC and Y6 exhibit clear differences (Fig. 4, right). For PBDB-T:ITIC and PM6:Y6, the significant spectral overlap of the donor and acceptor EPR signatures leads to overlapping contributions of both donor and acceptor polarons to the ENDOR spectra at most field positions. The ENDOR spectra recorded at the high-field end of the EPR spectra, where contributions of the donor polaron are dominant, are almost identical to the ENDOR spectra recorded for the corresponding blends with the PC₆₁BM acceptor (see Fig. S5 in the SI), confirming the absence of significant changes in spin delocalisation for donor polarons in different blends. Clear differences can be observed in the ENDOR spectra for lower magnetic field positions and are attributed to the polarons on the non-fullerene acceptors. At the lowest field positions, the contributions of the non-fullerene acceptors are almost exclusively selected. For both PBDB-T:ITIC and PM6:Y6, the acceptor ENDOR spectra are characterised by broader features, starting to cover a range of at least 20 MHz for PBDB-T:ITIC and 18 MHz for PM6:Y6. The shape of the spectra close to the Larmor frequency also starts to change, with a sharp peak centred at the Larmor frequency for both acceptors. Additionally, for ITIC˙⁻, pronounced peaks spaced apart by ca. 3 MHz clearly contribute to the spectra at the three lowest field positions. For Y6, on the other hand, the spectra at the lowest field positions consist of less pronounced shoulders centred around the ¹H Larmor frequency, with a gradual decrease in intensity towards the edges of the ENDOR spectrum. Additionally, for the PM6:Y6 blend, broad contributions to the ENDOR spectra centred around the fluorine Larmor frequency (ν_¹⁹F ≈ 134.3 MHz, ca. 8.4 MHz lower than ν_¹H) are also detected and assigned to the fluorine nuclei of the Y6 backbone.

The significantly increased width of the ¹H ENDOR spectra of ITIC˙⁻ and Y6˙⁻ compared to PBDB-T˙⁺ and PM6˙⁺ indicates the presence of relatively strongly coupled protons for the acceptor polarons, while the hyperfine couplings of protons on the polymer donors do not exceed 5 MHz, indicating more extensive spin delocalisation. The relatively featureless ENDOR lineshapes prevent clear determination of individual proton hyperfine couplings without further simulations, and are indicative of contributions of protons across a range of coupling strengths.

2.2.2 Simulation of ENDOR spectra based on DFT results. In order to gain further insights into the electron spin density distribution in the investigated donor and acceptor molecules, the experimental ENDOR results were compared to simulations based on hyperfine couplings predicted by DFT calculations.

Spin delocalisation for polarons on the donor polymers. A series of molecular models with increasing chain lengths were constructed for the polymer donors PBDB-T and PM6 and the spin density distributions and ¹H hyperfine couplings were calculated using two different functionals: the hybrid PBE0 functional and the range-separated hybrid ωB97X-D4 functional with tuned range-separation parameter ω (see Section S3 in the SI for details). While the PBE0 functional predicts complete delocalisation of the electron spin across the full backbone for all of the considered models (n = 1–4), the ωB97X-D4 functional consistently predicts localisation of the electron spin mainly on one electron-donating 2-alkylthiophene-substituted benzodithiophene (BDT) and the two adjacent electron-accepting bis(thiophenyl)-benzodithiophenedione (BDD) units for models with two or more repeat units (Fig. S9 and S10 in the SI).

Fig. 5 compares the ENDOR spectra recorded for PBDB-T˙⁺ and PM6˙⁺ at the field position corresponding to the spectral maximum with simulations based on PBE0/EPRII-predicted hyperfine couplings for models with spin densities delocalised over two to four repeat units. The simulations only include contributions from protons on the polymer backbone (H_BDT and H_thio in Fig. 5) and on the thiophene side groups, since the contributions of weakly coupled protons on the alkyl side chains are mostly suppressed by the central blind spot.^58–60 The overall width of the ENDOR spectrum is determined by the largest hyperfine couplings, and therefore a measure of the extent of spin delocalisation. Since hyperfine couplings are directly proportional to the spin density, increased delocalisation of the electron spin would lead to decreased hyperfine couplings and therefore ENDOR spectra with reduced overall width. The largest coupling predicted for spin delocalisation over two repeat units exceeds the largest coupling measured experimentally, whereas delocalisation over four repeat units would lead to a narrower ENDOR spectrum than observed experimentally. The overall shape and width of simulated spectra for models of intermediate length is in reasonable agreement with experiment, with the best match for PBDB-T_2.5A and PM6₃, respectively. Simulations of the full set of ENDOR spectra recorded across the EPR spectrum for these models reveal similar orientation selection effects as observed experimentally, with a narrowing of the spectrum towards lower fields (g_x), confirming accurate modelling of the hyperfine anisotropies and principal g-axis orientations by DFT (Fig. S13). Remaining discrepancies between experiment and simulations based on DFT predictions are likely due to distributions in hyperfine couplings expected based on the intrinsic heterogeneity of molecular environments in bulk-heterojunction blends, which is not reflected in a single DFT-optimised structure and will be addressed below.

Fig. 5 Comparison of the experimental W-band ¹H ENDOR spectra recorded for PBDB-T˙⁺ and PM6˙⁺ at the field position corresponding to the maximum of the EPR spectrum (g_y, g_z) with simulations based on hyperfine couplings predicted by DFT for models of increasing chain length (number of repeat units = 2, 2.5, 3 and 4) at the PBE0/EPRII level (see Section S3 in the SI for details and Fig. S14 and S15 for a comparison including additional models and results for both the PBE0 and ωB97X-D4 functional). The corresponding molecular structures and spin density distributions are also shown, as well as the orientation of the hyperfine tensors with respect to the molecular structure and g-frame for a selected model.

In contrast to the results of DFT calculations with the PBE0 functional, predictions based on the range-separated ωB97X-D4 functional widely used to model organic semiconductors^34,38,39 fail to match the experimental results in this instance, predicting more localised spin density distributions that results in ENDOR spectra with different widths and lineshapes compared to those observed experimentally (Fig. S14 and S15 in the SI).

Spin delocalisation for polarons on non-fullerene acceptors. The ENDOR spectra of the acceptor molecules ITIC and Y6 contain information on the spatial distribution of the electron spin in the photoinduced negative polarons and are compared to simulations based on DFT-predicted hyperfine couplings in Fig. 6.

Fig. 6 Comparison of the experimental W-band ENDOR spectra recorded for the PBDB-T:ITIC and PM6:Y6 blends at the field position corresponding to the maximum of the EPR spectrum of ITIC˙⁻ and Y6˙⁻ (g_y) with simulations based on hyperfine couplings predicted by DFT at the PBE0/EPRII level for a single ITIC or Y6 molecule and for different types of dimers extracted from the ITIC⁶¹ and Y6²⁶ crystal structures. The contribution of the donor ¹H ENDOR spectrum at the selected field position is shown as a blue shaded area and the simulations of the acceptor contribution are shown as red/orange shaded areas. For Y6, simulations of the ¹⁹F contributions to the ENDOR spectrum are also shown in green. The overall simulations, including the donor contribution, are displayed as dashed lines. The molecular structures and spin density distributions corresponding to the different models are also shown, as well as the orientation of the hyperfine tensors with respect to the molecular structure and g-frame for the isolated molecules.

DFT calculations performed for a radical anion on the ITIC molecule predict complete delocalisation across the conjugated backbone, independent of functional used, and the corresponding largest hyperfine couplings for the protons on the vinyl linker and on the outer edges of the central indacenodithienothiophene (IDTT) unit (H_vinyl and H_IDTT,o in Fig. 6, [−6, −13.4, −18]  MHz and [−4, −8.4, −10.8]  MHz, respectively) agree reasonably well with the width of the broad wings in the experimental ¹H ENDOR spectrum. Contributions of the remaining more weakly coupled protons on the central indacene and on the terminal 1,1-dicyanomethylene-3-indanone (INCN) groups (H_IDTT,i and H_INCN, −4 to 1 MHz) give rise to the pronounced narrower shoulders in the centre of the ENDOR spectrum. The DFT-predicted dependence of the ENDOR spectra on the field position is in qualitative agreement with experiment (Fig. S16 in the SI), however, while DFT predicts clearly resolved contributions for each of the more strongly coupled protons in the wings of the ENDOR spectrum, the experimental spectra suggest the presence of a distribution of couplings, implying some degree of heterogeneity. Possibilities are contributions of spins either localised on part of the ITIC molecule or delocalised over two molecules in crystalline regions of the blend, for example across ITIC dimers with π–π stacked terminal groups identified in the crystal structure of the pristine material.⁶¹ DFT calculations performed on an ITIC dimer modelled based on the crystal structure predict delocalisation over both molecules with increased spin density in the regions of overlap (Fig. 6). While the largest ¹H hyperfine couplings in the region of terminal–terminal interaction of two adjacent ITIC molecules are of similar magnitude to the ones in the isolated molecule, the additional delocalisation leads to the increased contribution of smaller couplings and therefore an ENDOR spectrum dominated by peaks at ca. ±0.5 MHz and shoulders at ±3 MHz (Fig. 6 and Fig. S16 in the SI). This simulated spectrum clearly deviates from the experimental ENDOR spectrum (Fig. 6), where the pronounced shoulders extending up to ±8 MHz imply dominant contribution of spins localised on a single ITIC molecule. The distribution of hyperfine couplings required to fully explain the experimental data is therefore likely due to the influence of the molecular environment on the nature of the spin density distribution within single ITIC molecules.

In the case of a radical anion on the Y6 molecule, DFT again predicts complete delocalisation of the spin density (Fig. 6). In this case, only three types of protons are expected to contribute to the ENDOR spectrum: two strongly coupled protons on the vinyl linkers (H_vinyl in Fig. 6, [−6, −13.5, −18.8]  MHz) and two groups of two protons on the terminal INCN groups (H_INCN, [−1.4, −3.1, −3.6]  MHz and [0.1, 0.4, 0.9]  MHz), in addition to two types of fluorine nuclei on the terminal INCN groups (F_INCN, [−1.4, −3.5, 20.5]  MHz and [−0.4, 0.9, −5.9]  MHz). Since both ITIC and Y6 contain terminal INCN acceptor groups, the types of protons contributing to the Y6 ENDOR spectra are a subset of the protons determining the appearance of the ITIC ENDOR spectra, with similar hyperfine coupling strengths predicted by DFT. This is at odds with the broad and relatively featureless ¹H ENDOR spectra observed experimentally for Y6, and, while the largest predicted ¹H couplings match the width of the wings and the smaller couplings reproduce the two narrower central shoulders, the experimental data suggests the presence of additional intermediate couplings between ca. 4 MHz and 10 MHz, which are absent in the simulation. Since morphological and spectroscopic studies of Y6 in pristine and blend films have demonstrated increased conformational rigidity and uniformity compared to other non-fullerene acceptors,^26,62 the presence of a large distribution of hyperfine couplings implied by the experimental ENDOR spectra is unlikely to be due to structural heterogeneity of the Y6 molecule alone. In order to explain the presence in the experimental spectrum of a range of intermediate hyperfine couplings, of smaller magnitude compared to the H_vinyl couplings predicted for an isolated Y6 molecule, contributions with either asymmetric delocalisation over a single Y6 molecule or with delocalisation extending beyond a single Y6 molecule may be considered. In the first case, a spin partially localised on one half of the molecule would lead to both a larger hyperfine coupling for one of the vinyl linker protons and a smaller hyperfine coupling for the other one. While this could account for contributions at intermediate coupling strengths, the experimental data does not contain any clear evidence for the presence of the accompanying more strongly coupled protons. Spin delocalisation extending across a pair of Y6 molecules, on the other hand, would lead to an overall reduction of all hyperfine couplings, by an extent determined by the degree of delocalisation.

Experimental and theoretical studies of molecular packing of Y6 have revealed an extended regular 3D network with increased π–π interactions for this non-fullerene acceptor, which is maintained even in blends, and has been proposed to determine the increased efficiency of Y6-based organic photovoltaics by improved charge separation and charge transport due to increased charge delocalisation.^26–28,63 Analysis of the Y6 crystal structure and molecular dynamics calculations have identified three types of π–π stacking motifs: pairs with overlapping terminal groups, Y6_2,TT, pairs where the terminal group of one molecule overlaps with the core of the other, Y6_2,CT–CT, and pairs with more extensive overlap of both core and terminal groups, Y6_2,CC–TT (see Fig. 6 and Fig. S11 in the SI).^26,27 The spin density distributions and ENDOR spectra calculated for these three types of Y6 pairs are shown in Fig. 6. DFT predicts spin density distributions spanning both molecules in the dimer, but localised mainly on the overlapping parts of the molecules, therefore the largest predicted hyperfine couplings are only slightly reduced compared to the monomer. However, asymmetric delocalisation across the two molecules of the pair leads to a larger spread towards smaller hyperfine couplings (Fig. 6 and Fig. S17 in the SI). A detailed comparison of the experimental results with simulations based on DFT predictions for the different models suggests the presence of spins localised on a single Y6 molecule as well as spins delocalised across different types of Y6 pairs, with contributions from Y6_2,TT and Y6_2,CC–TT required to account for the presence of both strong and intermediate hyperfine couplings across the full range of probed field positions (Fig. S17). The absence of distinct features in the experimental ENDOR spectra is again an indication for the presence of a distribution of hyperfine couplings due to different local spin environments in individual Y6 molecules and π–π stacked pairs influencing the extent of spin delocalisation. While an independent experimental characterisation of the spin density distribution in isolated Y6 molecules could further support this interpretation, such measurements were prevented so far by the strong propensity of this acceptor molecule to aggregate under the conditions required for ENDOR spectroscopy.^63–65 To validate our conclusion, we have instead extended our modelling approach to reduce reliance on the results of DFT calculations in interpreting the ENDOR data and additionally characterised the ¹⁴N hyperfine couplings using HYSCORE spectroscopy, as discussed in the following sections.

2.2.3 Modelling distributions by regularised fitting. In order to fully explain the experimental ENDOR data for positive and negative polarons in organic semiconducting materials, the complex structural and electronic heterogeneity at the molecular level, leading to varying extents of spin delocalisation and therefore distributions of hyperfine couplings, needs to be taken into account in the analysis. We employ a regularised least-squares fitting approach, introduced by Pribitzer et al. for the analysis of ENDOR data of systems with several nuclei and distributions of hyperfine couplings,⁴⁰ that exploits prior knowledge on the hyperfine interactions derived from quantum chemical modelling to guide fitting of the experimental results.

The ENDOR data is analysed in terms of a multidimensional probability distribution over the parameters determining the hyperfine interaction: the isotropic hyperfine coupling a_iso, the axial dipolar hyperfine coupling T, the rhombicity of the dipolar hyperfine interaction η and, in order to model orientation selection effects, the three Euler angles defining the orientation of the principal hyperfine axes with respect to the g-frame. The probability distribution over this set of hyperfine parameters for the spin system under investigation is obtained by combining least-squares fitting with a penalty term determined by the deviation between the best-fit distribution and a Bayesian prior distribution derived from DFT calculations, weighted by a regularisation parameter λ (see Section 4.3 in the SI for details). The prior probability distribution was constructed based on the DFT results for the different model systems of the donor polymers and acceptor molecules introduced above (Fig. S19 in the SI). Significant broadening was included to account for the uncertainty of the DFT predictions and additional structural variability.

The results of simultaneous fits of the full set of field-dependent ENDOR spectra for PBDB-T˙⁺, PM6˙⁺, ITIC˙⁻ and Y6˙⁻ are shown in Fig. 7, along with projections of the final probability distributions along a_iso and T, compared to DFT predictions for individual protons. An excellent fit of the experimental ENDOR spectra is obtained for probability distributions that are mostly centred around the DFT predictions for the considered molecular models, but include relatively broad distributions of hyperfine couplings around these individual solutions. In all cases, a clear correlation between the isotropic hyperfine coupling a_iso, determined by the spin density on the nucleus, and the strength of the anisotropic hyperfine coupling T, determined by the distribution of spin density on the molecular structure surrounding the nucleus, is maintained, as expected for protons directly bonded to an aromatic carbon in a π-conjugated system.^53,66

Fig. 7 Comparison of the experimental W-band ENDOR spectra recorded for the PBDB-T:PC₆₁BM, PBDB-T:ITIC, PM6:PC₆₁BM and PM6:Y6 blends with simulations resulting from regularised least-squares fitting with a penalty term based on the agreement of the probability distribution P(a_iso,T,η,[α,β,γ]) with a prior constructed based on DFT calculations (see Section S4.3 and Fig. S19 in the SI for details). The ENDOR contributions of the donor and acceptor polarons are shown as shaded areas and the overall simulation is displayed as a dashed black line. In the bottom panels, a projection of the final probability distribution onto the a_iso, T space is compared to the results of DFT calculations for individual models (colour-coded dots represent different types of protons, and regions included in the prior for the protons with the largest hyperfine couplings are highlighted). The prior was based on polymer models with n = 1–4 repeat units for PBDB-T and PM6, on a single molecule for ITIC and on the single molecule as well as three pair configurations for Y6. For additional details, see sections S4.3 and S12 in the SI.

In the case of the polarons on the donor polymers, PBDB-T˙⁺ and PM6˙⁺, the obtained probability distributions encompass predictions of models where the electron spin is delocalised mostly over three (6 nm) and up to four (8 nm) repeat units, with almost negligible contributions for delocalisation over two repeat units. In contrast to DFT models for a single static model structure, this distribution of polarons characterised by different extents of spin delocalisation is able to fully account for the experimental observations, including the observed orientation-selection effects. The relatively broad distribution of hyperfine parameters is likely determined by differences in local chain conformations and the surrounding molecular environment, in agreement with previous findings for other organic semiconducting polymers.²⁵

For the non-fullerene acceptor ITIC, the probability distribution determined by regularised least-squares fitting shows that the broad shoulders extending from ±2.5 MHz to ±10 MHz can be explained by distributions of hyperfine couplings centred around a_iso = −13 MHz, T = 4 MHz and a_iso = −8 MHz, T = 2.5 MHz, both slightly larger than DFT predictions for the protons on the vinyl linker and the outer edges of the IDTT core, respectively. This may suggest a slight redistribution of spin density compared to DFT predictions, either towards the core or towards one of the end groups, the latter resulting in an asymmetric delocalisation with increased spin density on one half of the molecule. Including the ITIC dimer with π–π stacked terminal groups in the prior for the regularised least-squares fitting does not significantly impact the resulting probability distribution (see Fig. S21 in the SI). There is, therefore, no evidence to suggest significant contributions from spins delocalised over pairs of ITIC molecules.

In the case of the non-fullerene acceptor Y6, on the other hand, comparison of the experimental ENDOR results with DFT predictions already suggested the likely contribution from spins delocalised over more than a single molecule. The broad shoulders in the ENDOR spectra, extending from ±3 MHz to ±9 MHz can only be due to strongly coupled protons on the vinyl linker (highlighted in orange in Fig. 6). Regularised least-squares fitting of the ENDOR data with a penalty term considering a prior based on DFT results for a single Y6 molecule results in a hyperfine distribution deviating significantly from DFT predictions, in particular extending towards decreased coupling strengths (Fig. S21 in the SI). A much closer agreement between the fitted probability distribution and the prior is obtained after inclusion of the three pair models constructed based on the crystal structure, Y6_2,TT, Y6_2,CT–CT and Y6_2,CC–TT (Fig. 7). The more extensive, and partly asymmetric, delocalisation across the two molecules of the pair leads to decreased hyperfine couplings of the vinyl linker protons (Fig. S19), which are required in addition to contributions from spins localised on single Y6 molecules to explain the broad, extended shoulders in the ENDOR spectrum. By considering both single molecules and different π–π stacked dimers, an excellent agreement between DFT predictions and the hyperfine distribution determined for the region of strong isotropic and dipolar coupling is obtained (Fig. 7).

Our results contribute to increasing experimental evidence for the ability of charge to delocalise across multiple Y6 molecules in donor:acceptor blends, following the identification of Y6_2,TT dimers by solid-state NMR⁶⁷ on PM6:Y6 blends, of J-aggregates (Y6_2,CT–CT) in UV-vis spectra of blend films⁶³ and evidence of exciton delocalisation in different Y6 derivatives.⁶⁸ By mapping out the spin density distributions based on electron–nuclear hyperfine couplings, ENDOR directly probes delocalisation in the photoinduced charged states. The ENDOR spectra clearly suggest contributions not only of spins localised on a single Y6 molecule, but also spins delocalised over two adjacent Y6 molecules in different possible configurations, including Y6_2,TT, Y6_2,CT–CT and Y6_2,CC–TT. The significant overlap of the corresponding ENDOR spectra prevents reliable identification of contributions from individual Y6 pair models or determination of relative weights. However, the high probability of relatively strong hyperfine couplings in the a_iso = −8 to −14 MHz, T = 2 to 4 MHz region indicates that, in most cases, even in Y6 pairs, the charge is not completely evenly delocalised over both molecules. Nevertheless, the ability for the charge to delocalise across different Y6 molecules, in contrast to observations for ITIC that suggest localisation of the charge on a single molecule, can explain the high charge separation and charge transport efficiencies determined for organic photovoltaic blends based on Y-series acceptors.

The analysis of the ENDOR data of polarons on donor and acceptor molecules for organic photovoltaics by regularised least-squares fitting provides a general means to account for the heterogeneity of molecular environments in donor:acceptor blends. While the results should not be over-interpreted, the probability distributions across the hyperfine parameter space provide insights beyond those available from comparison of experimental ENDOR spectra with simulations for single static DFT models, yet still benefit from the connection between spin density distributions predicted by DFT and experimentally observed hyperfine couplings.

2.3 HYSCORE on non-fullerene acceptors

For the negative polarons on the non-fullerene acceptors ITIC and Y6, analysis of the ¹⁴N hyperfine couplings can provide further evidence for the conclusions drawn from the analysis of the ¹H ENDOR data. While the W-band ELDOR-detected NMR spectra only showed broad unresolved peaks, which are challenging to interpret in terms of the strength of hyperfine couplings due to the additional contribution of the nuclear quadrupole interaction for this I = 1 nucleus, HYSCORE experiments offer increased resolution by correlating nuclear frequencies in different electron spin manifolds in a two-dimensional experiment, allowing easier determination of both nuclear quadrupole and hyperfine interaction parameters from the observed correlation patterns.

The results of Q-band HYSCORE experiments for PBDB-T:ITIC and PM6:Y6 are shown in Fig. 8. Due to the large excitation bandwidth of the pulses used for the HYSCORE experiment, orientation selection effects are negligible for experiments performed at the maximum of the ITIC˙⁻ or Y6˙⁻ EPR signal. The HYSCORE spectra for the two blends show strong similarities, as expected given ITIC and Y6 share the same nitrogen-containing INCN terminal groups. The main signals are confined to the (+, +) quadrant of the HYSCORE spectrum, indicating the ¹⁴N hyperfine interactions are in the weak coupling limit at Q-band (|A| < 2 ν_¹⁴N, ν_¹⁴N = 3.7 MHz). For an I = 1 nucleus such as ¹⁴N, each electron spin sublevel is split into three nuclear sublevels connected by three nuclear transitions (two single-quantum Δm_I = ±1 transitions, and one double-quantum Δm_I = ±2 transition). The HYSCORE experiment correlates the nuclear transition frequencies across the two electron spin manifolds, providing information on the electron–nuclear hyperfine and nuclear quadrupole interaction causing the sublevel splitting. In the spectra shown in Fig. 8, features around 2.1 MHz correspond to correlations between single quantum transitions, those around 5.7 MHz to single-double quantum correlations, and the ridges perpendicular to the diagonal at around 8 MHz are assigned to the double quantum transitions (dq,dq). Signals due to ¹³C hyperfine couplings centred at the corresponding Larmor frequency (ν_¹³C = 13 MHz, ) are also visible. The most precise information on the ¹⁴N hyperfine couplings can be extracted from the extent of the double quantum ridges, highlighted in Fig. 8, and suggests larger couplings for ITIC compared to Y6.

Fig. 8 Experimental Q-band ¹⁴N HYSCORE spectra (top) for the PBDB-T:ITIC and PM6:Y6 blends recorded at a field position corresponding to the maximum of the ITIC˙⁻ and Y6˙⁻ EPR spectrum and corresponding simulations (bottom) based on ¹⁴N hyperfine couplings predicted by DFT calculations for a negative charge on an isolated ITIC molecule and a dimer of Y6 molecules with overlapping terminal groups (Y6_2,TT). Only the (+, +) quadrant is shown, as no significant contributions were observed in the (−, +) quadrant. The molecular structures and calculated ¹⁴N hyperfine tensors are shown at the bottom of the figure. Simulations for the other Y6 models are shown in Fig. S23 in the SI and additional experimental details can be found in Section S2.5.

The simulation of the HYSCORE spectrum based on the hyperfine couplings predicted by DFT for the single ITIC molecule reproduces the main features observed in the experimental spectrum well, including the approximate width of the double quantum ridges. This suggests the spin density distribution predicted by DFT is reflective of the negative ITIC˙⁻ polaron in the blend film. In agreement with the ENDOR results, the experimental HYSCORE data shows increased broadening of all features indicating a distribution of nuclear quadrupole and hyperfine parameters centred around the values predicted by DFT.

In the case of Y6, which, in addition to the nitrogen nuclei on the terminal INCN groups, also contains nitrogens on the central core, larger deviations between simulation and experiment are observed for the DFT model of the single Y6 molecule (Fig. S23). The extent of the double-quantum ridge observed experimentally is reduced compared to that observed for ITIC and that predicted for spin density distributed across a single Y6 molecule, indicating a reduction in hyperfine coupling strength and therefore suggesting contributions with increased spin delocalisation. The best agreement with experiment is obtained for the Y6 dimer model with overlapping terminal groups (Y6_2,TT), where the more extensive delocalisation across both molecules leads to a reduction of the largest ¹⁴N hyperfine couplings and therefore a reduced width of the double-quantum feature (Fig. 8). The HYSCORE simulations for the remaining Y6 models mainly differ in the extent of the double-quantum ridge, which is still always centred around 2ν_¹⁴N and therefore overlaps at the centre (see Fig. S23 in the SI). A linear combination of contributions from the single Y6 molecule as well as the different pair configurations, with an increased weight for Y6_2,TT, provide a good agreement with the experimental spectrum. This further confirms the conclusions drawn based on the ENDOR results, suggesting the presence of Y6 polarons characterised by different extents of delocalisation across one or two molecules.

3 Conclusions

In this study, we have used a series of pulse EPR techniques for the selective detection and in-depth characterisation of photoinduced charged states in the donor:non-fullerene acceptor bulk heterojunction blends PBDB-T:ITIC and PM6:Y6, as well as the corresponding fullerene-based blends PBDB-T:PC₆₁BM and PM6:PC₆₁BM. The electron spin provides a sensitive probe for the local molecular and electronic structure of positive and negative polarons involved in photovoltaic energy conversion.

Initially, we have presented a full characterisation of the EPR spectra and principal g-values of PBDB-T˙⁺, PM6˙⁺, ITIC˙⁻ and Y6˙⁻, exploiting EDNMR-induced EPR³¹ for the accurate determination of the ITIC and Y6 spectra based on the hyperfine couplings to ¹⁴N nuclei present on the non-fullerene acceptors but absent on the donor polymers. While the W- and Q-band EPR spectra of the non-fullerene acceptors are dominated by the g-anisotropy, significant hyperfine couplings in ITIC˙⁻ and Y6˙⁻ contribute to the EPR lineshape at X-band. The accurate identification of the donor and acceptor spectral signatures across all three frequency bands is essential for studies of charge-transfer states by transient EPR^17,21,69,70 and of spin-dependent charge transport and recombination processes by optically or electrically detected magnetic resonance.^71–73

The measurement of electron–nuclear hyperfine couplings to ¹H and ¹⁹F by W-band Davies ENDOR and to ¹⁴N by Q-band HYSCORE, combined with DFT modelling and regularised least-squares fitting, has then allowed us to quantify the extent of spin and charge delocalisation. For the positive polarons on the donor polymers PBDB-T and PM6, we have identified a distribution of slightly different extents of intramolecular delocalisation along the π-conjugated polymer backbone with the highest probability for charge delocalisation across about three repeat units (ca. 6 nm). Our experimental ENDOR results have allowed us to validate different DFT approaches for the estimation of polaron delocalisation in donor polymers, concluding that, in the case PBDB-T˙⁺ and PM6˙⁺, the range-separated hybrid ωB97X-D4 functional underestimates the extent of spin delocalisation, whereas a close agreement with experiment is obtained for calculations with the standard hybrid PBE0 functional on truncated polymer models.

The ENDOR and HYSCORE results for the negative polarons on the non-fullerene acceptors ITIC and Y6 have revealed dominant contributions from spins localised on a single molecule for ITIC˙⁻, whereas a broader distribution of hyperfine couplings for Y6˙⁻ indicates additional contributions from spins partially delocalised over two π–π stacked Y6 molecules, in agreement with different dimer configurations previously identified in pristine films and bulk heterojunction blends.^26,27,63,67 The measurement of hyperfine interactions provides direct evidence of intermolecular charge delocalisation for photoinduced polarons on Y6, with the electron spin density extending over one to two molecules and localising predominantly in regions of overlap of the conjugated backbones in the case of Y6 pairs. The complex morphology of Y6 in the solid state and in blends thus results in a distribution of local polaron environments, motivating future studies to further disentangle the interplay between local morphology, density of states and the extent of charge delocalisation. Since charge delocalisation is believed to play an essential role in promoting efficient charge separation and charge transport, additional experimental information on the charge distribution in photoinduced polarons can support multiscale modelling aimed at unravelling the mechanistic details at the basis of photovoltaic energy conversion.

Conflicts of interest

There are no conflicts to declare.

Data availability

Supplementary information (SI): all experimental details, description of analysis, simulation and computational methods, as well as additional supporting experimental results. See DOI: https://doi.org/10.1039/d5cp03276h.

The experimental data files, results of DFT calculations and the main simulation and fitting code are available on the Oxford Research Archive (ORA) at https://doi.org/10.5287/ora-obj9g5nxo.

Acknowledgements

The authors gratefully acknowledge the support of the Oxford-Berlin Research Partnership for funding enabling the W-band EPR measurements in the laboratory of Prof. Robert Bittl at the Freie Universität Berlin and would like to thank Prof. Bittl and Dr Christian Teutloff for access to the instrumentation, support and insightful discussions. The authors would also like to thank Prof. Stefan Stoll and an anonymous reviewer for valuable comments on the manuscript. The authors acknowledge the use of the University of Oxford Advanced Research Computing (ARC) facility in carrying out this work (https://doi.org/10.5281/zenodo.22558). C. E. T. is thankful to the Royal Society for a University Research Fellowship (URF\R1\201071) and to Balliol College, Oxford, for an Early Career Fellowship.

References

1. Y. Jiang, S. Sun, R. Xu, F. Liu, X. Miao, G. Ran, K. Liu, Y. Yi, W. Zhang and X. Zhu, Nat. Energy, 2024, 9, 975–986.
2. W. Su, X. Zhou, Q. Wu, Y. Wu, H. Qin, Z. Liang, H. Li, H. Bai, J. Guo, L. Jiang, Y. Liu, R. Ma, Y. Li, W. Zhu and Q. Fan, Adv. Funct. Mater., 2024, 2415090.
3. W. Liang, S. Zhu, K. Sun, J. Hai, Y. Cui, C. Gao, W. Li, Z. Wu, G. Zhang and H. Hu, Adv. Funct. Mater., 2024, 2415499.
4. L. Zhu, M. Zhang, G. Zhou, Z. Wang, W. Zhong, J. Zhuang, Z. Zhou, X. Gao, L. Kan, B. Hao, F. Han, R. Zeng, X. Xue, S. Xu, H. Jing, B. Xiao, H. Zhu, Y. Zhang and F. Liu, Joule, 2024, 8, 3153–3168.
5. E. M. Herzig, F. Gao, J. Bergqvist, M. A. Loi and S. B. Meier, Joule, 2024, 8, 2171–2178.
6. A. Armin, W. Li, O. J. Sandberg, Z. Xiao, L. Ding, J. Nelson, D. Neher, K. Vandewal, S. Shoaee, T. Wang, H. Ade, T. Heumüller, C. J. Brabec and P. Meredith, Adv. Energy Mater., 2021, 11, 20003570.
7. S. Shoaee, H. M. Luong, J. Song, Y. Zou, T.-Q. Nguyen and D. Neher, Adv. Mater., 2023, 2302005.
8. J. Wang, Y. Xie, K. Chen, H. Wu, J. M. Hodgkiss and X. Zhan, Nat. Rev. Phys., 2024, 6, 365–381.
9. A. A. Bakulin, A. Rao, V. G. Pavelyev, P. H. M. van Loosdrecht, M. S. Pshenichnikov, D. Niedzialek, J. Cornil, D. Beljonne and R. H. Friend, Science, 2012, 335, 1340–1345.
10. M. Casalegno, R. Pastore, J. Idé, R. Po and G. Raos, J. Phys. Chem. C, 2017, 121, 16693–16701.
11. S. Athanasopoulos, H. Bässler and A. Köhler, J. Phys. Chem. Lett., 2019, 10, 7107–7112.
12. D. Balzer, T. J. Smolders, D. Blyth, S. N. Hood and I. Kassal, Chem. Sci., 2021, 12, 2276–2285.
13. D. Balzer and I. Kassal, Sci. Adv., 2022, 8, eabl9692.
14. D. Balzer and I. Kassal, Chem. Sci., 2024, 15, 4779–4789.
15. X. Jia, L. Soprani, G. Londi, S. M. Hosseini, F. Talnack, S. Mannsfeld, S. Shoaee, D. Neher, S. Reineke, L. Muccioli, G. D'Avino, K. Vandewal, D. Beljonne and D. Spoltore, Mater. Horiz., 2024, 11, 173–183.
16. J. Niklas and O. G. Poluektov, Adv. Energy Mater., 2017, 7, 1602226.
17. F. Kraffert and J. Behrends, Mol. Phys., 2017, 115, 2373–2386.
18. S. Weber, eMagRes, 2017, 6, 255–270.
19. A. Aguirre, P. Gast, S. Orlinskii, I. Akimoto, E. J. J. Groenen, H. El Mkami, E. Goovaerts and S. Van Doorslaer, Phys. Chem. Chem. Phys., 2008, 10, 7129–7138.
20. J. Niklas, K. L. Mardis, B. P. Banks, G. M. Grooms, A. Sperlich, V. Dyakonov, S. Beaupré, M. Leclerc, T. Xu, L. Yu and O. G. Poluektov, Phys. Chem. Chem. Phys., 2013, 15, 9562–9574.
21. S. A. J. Thomson, J. Niklas, K. L. Mardis, C. Mallares, I. D. W. Samuel and O. G. Poluektov, J. Phys. Chem. C, 2017, 121, 22707–22719.
22. R. Steyrleuthner, Y. Zhang, L. Zhang, F. Kraffert, B. P. Cherniawski, R. Bittl, A. L. Briseno, J.-L. Brédas and J. Behrends, Phys. Chem. Chem. Phys., 2017, 19, 3627–3639.
23. J. Niklas, T. Zheng, A. Neshchadin, K. L. Mardis, L. Yu and O. G. Poluektov, J. Am. Chem. Soc., 2020, 142, 1359–1366.
24. M. Arvind, C. E. Tait, M. Guerrini, J. Krumland, A. M. Valencia, C. Cocchi, A. E. Mansour, N. Koch, S. Barlow, S. R. Marder, J. Behrends and D. Neher, J. Phys. Chem. B, 2020, 124, 7694–7708.
25. C. E. Tait, A. Reckwitz, M. Arvind, D. Neher, R. Bittl and J. Behrends, Phys. Chem. Chem. Phys., 2021, 23, 13827–13841.
26. G. Zhang, X. K. Chen, J. Xiao, P. C. Y. Chow, M. Ren, G. Kupgan, X. Jiao, C. C. S. Chan, X. Du, R. Xia, Z. Chen, J. Yuan, Y. Zhang, S. Zhang, Y. Liu, Y. Zou, H. Yan, K. S. Wong, V. Coropceanu, N. Li, C. J. Brabec, J.-L. Brédas, H. L. Yip and Y. Cao, Nat. Commun., 2020, 11, 3943.
27. G. Kupgan, X. K. Chen and J.-L. Brédas, Mater. Today Adv., 2021, 11, 100154.
28. D. Kroh, S. Athanasopoulos, V. Nádaždy, F.-J. Kahle, H. Bässler and A. Köhler, Adv. Funct. Mater., 2023, 2302520.
29. Y. Lin, J. Wang, Z. G. Zhang, H. Bai, Y. Li, D. Zhu and X. Zhan, Adv. Mater., 2015, 27, 1170–1174.
30. J. Yuan, Y. Zhang, J. Yuan, Y. Zhang, L. Zhou, G. Zhang, H.-L. Yip, T.-K. Lau, X. Lu, C. Zhu, H. Peng, P. A. Johnson, M. Leclerc, Y. Cao, J. Ulanski, Y. Li and Y. Zou, Joule, 2019, 3, 1140–1151.
31. M. Van Landeghem, W. Maes, E. Goovaerts and S. Van Doorslaer, J. Magn. Reson., 2018, 288, 1–10.
32. M. Van Landeghem, J. Kudrjasova, W. Maes, E. Goovaerts and S. Van Doorslaer, Appl. Magn. Reson., 2019, 50, 1253–1265.
33. K. R. Bryenton, A. A. Adeleke, S. G. Dale and E. R. Johnson, Wiley Interdiscip. Rev.:Comput. Mol. Sci., 2023, 13, e1631.
34. T. Körzdörfer, J. S. Sears, C. Sutton and J.-L. Brédas, J. Chem. Phys., 2011, 135, 204107.
35. T. Körzdörfer and J.-L. Brédas, Acc. Chem. Res., 2014, 47, 3284–3291.
36. K. Do, M. K. Ravva, T. Wang and J.-L. Brédas, Chem. Mater., 2017, 29, 346–354.
37. V. Coropceanu, X.-K. Chen, T. Wang, Z. Zheng and J.-L. Brédas, Nat. Rev. Mater., 2019, 4, 689–707.
38. T. Wang, G. Kupgan and J.-L. Brédas, Trends Chem., 2020, 2, 535–554.
39. K. J. Thorley, J. Phys. Chem. B, 2023, 127, 5102–5114.
40. S. Pribitzer, D. Mannikko and S. Stoll, Phys. Chem. Chem. Phys., 2021, 23, 8326–8335.
41. J. De Ceuster, E. Goovaerts, A. Bouwen, J. C. Hummelen and V. Dyakonov, Phys. Rev. B:Condens. Matter Mater. Phys., 2001, 64, 195206.
42. O. G. Poluektov, S. Filippone, N. Martín, A. Sperlich, C. Deibel and V. Dyakonov, J. Phys. Chem. B, 2010, 114, 14426–14429.
43. K. L. Mardis, J. N. Webb, T. Holloway, J. Niklas and O. G. Poluektov, J. Phys. Chem. Lett., 2015, 6, 4730–4735.
44. J. Niklas, K. L. Mardis and O. G. Poluektov, J. Phys. Chem. Lett., 2018, 9, 3915–3921.
45. S. Stoll and A. Schweiger, J. Magn. Reson., 2006, 178, 42–55.
46. F. Kraffert, R. Steyrleuthner, S. Albrecht, D. Neher, M. C. Scharber, R. Bittl and J. Behrends, J. Phys. Chem. C, 2014, 118, 28482–28493.
47. J. J. Wittmann, T. V. Can, M. Eckardt, W. Harneit, R. G. Griffin and B. Corzilius, J. Magn. Reson., 2018, 290, 12–17.
48. D. Goldfarb, eMagRes, 2017, 6, 101–114.
49. N. Cox, A. I. Nalepa, W. Lubitz and A. Savitsky, J. Magn. Reson., 2017, 280, 63–78.
50. T. Maly and T. F. Prisner, J. Magn. Reson., 2004, 170, 88–96.
51. T. Maly, F. MacMillan, K. Zwicker, N. Kashani-Poor, U. Brandt and T. F. Prisner, Biochemistry, 2004, 43, 3969–3978.
52. A. J. Stone, Mol. Phys., 1964, 7, 311–316.
53. N. M. Atherton, Principles of Electron Spin Resonance, PTR Prentice Hall, New York, 1993.
54. V. I. Krinichnyi, H.-K. Roth and M. Schrödner, Appl. Magn. Reson., 2002, 23, 1–17.
55. V. I. Krinichnyi, P. A. Troshin and N. N. Denisov, J. Chem. Phys., 2008, 128, 164715.
56. T. S. Cameron, R. C. Haddon, S. M. Mattar, S. Parsons, J. Passmore and A. P. Ramirez, J. Chem. Soc., Dalton Trans., 1992, 1563–1572.
57. A. Schweiger and G. Jeschke, Principles of pulse electron paramagnetic resonance, Oxford University Press, 2001.
58. A. V. Astashkin and A. Kawamori, J. Magn. Reson., 1998, 135, 406–417.
59. C. Fan, P. E. Doan, C. E. Davoust and B. M. Hoffman, J. Magn. Reson., 1992, 98, 62–72.
60. P. E. Doan, N. S. Lees, M. Shanmugam and B. M. Hoffman, Appl. Magn. Reson., 2010, 37, 763–779.
61. T. J. Aldrich, M. Matta, W. Zhu, S. M. Swick, C. L. Stern, G. C. Schatz, A. Facchetti, F. S. Melkonyan and T. J. Marks, J. Am. Chem. Soc., 2019, 141, 3274–3287.
62. J. Wu, J. Lee, Y. C. Chin, H. Yao, H. Cha, J. Luke, J. Hou, J. S. Kim and J. R. Durrant, Energy Environ. Sci., 2020, 13, 2422–2430.
63. D. Kroh, F. Eller, K. Schötz, S. Wedler, L. Perdigón-Toro, G. Freychet, Q. Wei, M. Dörr, D. Jones, Y. Zou, E. M. Herzig, D. Neher and A. Köhler, Adv. Funct. Mater., 2022, 32, 2205711.
64. L. Xu, S. Li, W. Zhao, Y. Xiong, J. Yu, J. Qin, G. Wang, R. Zhang, T. Zhang, Z. Mu, J. Zhao, Y. Zhang, S. Zhang, V. Kuvondikov, E. Zakhidov, Q. Peng, N. Wang, G. Xing, F. Gao, J. Hou, W. Huang and J. Wang, Adv. Mater., 2024, 36, 2403476.
65. G. Cai, Y. Li, Y. Fu, H. Yang, L. Mei, J.-L. Brédas, M.-C. Tang and X. Chen, Nat. Commun., 2024, 15, 2784.
66. J. R. Morton, Chem. Rev., 1964, 64, 453–471.
67. B. R. Luginbuhl, P. Raval, T. Pawlak, Z. Du, T. Wang, G. Kupgan, N. Schopp, S. Chae, S. Yoon, A. Yi, H. Jung Kim, V. Coropceanu, J.-L. Brédas, T.-Q. Nguyen and G. N. M. Reddy, Adv. Mater., 2022, 34, 2105943.
68. R. J. E. Westbrook, A. J. Levin, W. Gao, U. Bothra, S. M. Pratik, B. Fan, F. R. Lin, Q.-Q. Zhang, K. Ngo, W. Kaminsky, J.-L. Brédas, V. Coropceanu, A. K.-Y. Jen, M. F. Toney and D. S. Ginger, J. Am. Chem. Soc., 2025, 147, 30199–30209.
69. J. Behrends, A. Sperlich, A. Schnegg, T. Biskup, C. Teutloff, K. Lips, V. Dyakonov and R. Bittl, Phys. Rev. B:Condens. Matter Mater. Phys., 2012, 85, 125206.
70. J. Niklas, S. Beaupré, M. Leclerc, T. Xu, L. Yu, A. Sperlich, V. Dyakonov and O. G. Poluektov, J. Phys. Chem. B, 2015, 119, 7407–7416.
71. D. R. McCamey, S.-Y. Lee, S.-Y. Paik, J. M. Lupton and C. Boehme, Phys. Rev. B:Condens. Matter Mater. Phys., 2010, 82, 125206.
72. A. Schnegg, J. Behrends, M. Fehr and K. Lips, Phys. Chem. Chem. Phys., 2012, 14, 14418.
73. F. Kraffert, R. Steyrleuthner, C. Meier, R. Bittl and J. Behrends, Appl. Phys. Lett., 2015, 107, 043302.

---