Modeling Ultrafast Exciton Migration within the ... - ACS Publications

Feb 21, 2017 - Department of Chemistry, DePaul University, 1110 West Belden Avenue, ... Franck Institute, Institute for Biophysical Dynamics, The Univ...
0 downloads 0 Views 4MB Size
Article pubs.acs.org/JPCC

Modeling Ultrafast Exciton Migration within the Electron Donor Domains of Bulk Heterojunction Organic Photovoltaics Mateusz Bednarz,† Joel Lapin,† Ryan McGillicuddy,‡ Kenley M. Pelzer,§ Gregory S. Engel,‡ and Graham B. Griffin*,† †

Department of Chemistry, DePaul University, 1110 West Belden Avenue, Chicago, Illinois 60614, United States Department of Chemistry, The James Franck Institute, Institute for Biophysical Dynamics, The University of Chicago, Chicago, Illinois 60637, United States § Materials Science Division, Argonne National Laboratory, 9700 Cass Avenue, Argonne, Illinois 60439, United States ‡

S Supporting Information *

ABSTRACT: Recent experimental studies revealed that charge carriers harvested by bulk heterojunction organic photovoltaics can be collected on ultrafast time scales. To investigate ultrafast exciton mobility, we construct simple, nonatomistic models of a common polymeric electron donor material. We first explore the relationship between the magnitude of energetic noise in the model Hamiltonian and the spatial extent of resulting eigenstates. We then employ a quantum master equation approach to simulate migration of chromophore-localized initial excited states. Excitons initially localized on a single chromophore at the center of the model delocalize down polymer chains and across pi-stacked chromophores through a coherent, wavelike mechanism during the first few tens of femtoseconds. We explore the dependence of this coherent delocalization on coupling strength and on the magnitude of energetic noise. At longer times we observe continued migration toward a uniform population distribution that proceeds through an incoherent, diffusive mechanism. A series of simulations modeling exciton harvesting in domains of varying size demonstrates that smaller domains enhance ultrafast exciton harvesting yield. Our nonatomistic model falls short of quantitative accuracy but demonstrates that excitons are mobile within electron donor domains on ultrafast time scales and that coherent exciton transport can enhance ultrafast exciton harvesting.

1. INTRODUCTION Bulk heterojunction organic photovoltaic devices (OPVs) have been widely studied since their introduction 20 years ago.1−4 The active layer of these devices contains a mixture of electrondonating molecules (frequently conjugated polymers) and electron-accepting molecules (usually fullerene derivatives) (Figure 1a). Photoexcitation of the electron donor creates bound electron−hole pairs that are conventionally described as diffusing to the donor−acceptor material interface through Förster resonance energy transfer.5 This incoherent, diffusive energy transport mechanism predicts exciton migration over tens of nanometers, on picosecond time scales (Figure 1b). Upon reaching the donor/acceptor interface, the excited electron is then transferred from the donor to the acceptor. This process is often described via Marcus’ nonadiabatic electron transfer theory, based on a Fermi’s Golden Rule calculation of the nonadiabatic rate for transfer of the excited electron from the donor to the acceptor.6 Due to the low dielectric constant of organic materials, the Coulombic binding energy between an excited electron in the acceptor domain and a hole in the donor domain is generally greater than the thermal energy at room temperature. This fact © XXXX American Chemical Society

inspired proposal of a mechanism in which exciton harvesting is mediated through formation of charge transfer excitons, wherein the excited electron is transferred into the acceptor but the exciton remains trapped at the material interface.7−10 Collection of fully separated charge carriers within this model requires further mechanisms to split the charge transfer exciton. Applying this conceptual framework, massive experimental efforts have optimized the molecular structure of donor and acceptor materials, device processing methods, and bulk heterojunction morphology,11 leading to a 10-fold increase in the power conversion efficiency of bulk heterojunction devices over the past decade.12,13 Recent ultrafast spectroscopy experiments have demonstrated that excitons created by photoexcitation of electron donor materials can be converted into free charge carriers within 100 fs.14−17 These studies build upon a substantial body of previous work on ultrafast polaron formation in heterojunction blends18,19 and polymer thin films.20−23 This ultrafast exciton Received: November 10, 2016 Revised: January 25, 2017 Published: February 21, 2017 A

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

interface, predicting that only a small percentage of excitons should be harvested on ultrafast time scales. Experimental results directly contradict this prediction, revealing that a substantial percentage of exciton harvesting occurs within 100 fs of photoexcitation.15,24 At such early times quantum mechanical effects can delocalize excitons, a mode of exciton migration not included in the incoherent, diffusive Förster mechanism. This delocalization determines both the spatial range over which exciton population can travel on ultrafast time scales and the shape of the exciton as it moves toward the donor/acceptor material interface. To design bulk heterojunction OPVs that optimize ultrafast exciton harvesting, a clearer picture of the interplay between the morphology of electron donor domains and the ultrafast behavior of excitons must be developed. To explore this relationship in an intuitive context, we created a nonatomistic model of electron donor domains, based on previous experimental measurements of the properties of crystalline poly-3-hexylthiophene (P3HT), a functional electron donor material. Our nonatomistic approach provides a simple and efficient way to model the domain, when contrasted with a fully quantum mechanical or QM/MM atomistic treatment, while still retaining meaningful connection to the parameters of the real system. This simplification permits us to build threedimensional domains comparable in size to those found in real bulk heterojunction devices and to simulate electronic dynamics within the model domains using desktop computers. We first use small ensembles of model electron donor domains to investigate the relationship between energetic disorder and the energetic and spatial properties of the eigenstates of the model domains. Next, we analyze how the spatial distribution of an initially localized exciton population evolves in the energy eigenbasis. We qualitatively examine the spatial pattern of exciton migration in the model donor domain as a function of the strength of the interchromophore couplings and the magnitude of energetic noise. We quantify the observed migration of exciton population by calculating the participation ratio of the density matrix as a function of time. We also calculate the quantum mechanical purity of the state of our model system as a function of time, employing this metric to characterize the time scales on which the observed exciton migration results from coherent, wavelike transport. Finally, we construct a series of model electron donor domains designed to interrogate the relationship between the size of an electron donor domain and the yield of ultrafast exciton harvesting. Our modeling strategy is inspired by methods used to explore ultrafast exciton migration in photosynthetic light harvesting systems,30−34 and the observed results are consistent with prior work on exciton delocalization in molecular aggregates.35,36 The patterns of migration observed in our time-domain simulations are also consistent with the results of previous frequency-domain modeling of solid-phase P3HT.37,38 The goals of this work are to demonstrate the importance of coherent exciton migration on ultrafast time scales and to stimulate the design of new devices in which bulk heterojunction composition and morphology has been optimized to enhance ultrafast exciton yield. While our simple model falls short of quantitative accuracy, our results provide a rational framework within which to consider the relationship between electron donor morphology and exciton transport as well as the potential impact of ultrafast exciton transport on exciton harvesting yields. The quantitative results serve to frame this consideration and may also provide a rough guide for

Figure 1. Emerging ultrafast channel for charge carrier collection in bulk heterojunction devices. (a) Bulk heterojunction organic photovoltaics consist of a blend of electron donor (blue) and electron acceptor (red) domains. Conjugated polymers or other conjugated organic molecules are commonly employed as electron donors, and fullerene derivatives are the most popular electron acceptors. (b) Excitation of conjugated organics leads to bound excitons formed from an electron (cyan) and a hole (yellow). Conventional exciton harvesting mechanisms are mediated by incoherent, diffusive exciton transport to the material interface. Formation and subsequent splitting of charge transfer excitons, with charge carriers separated but still Coulombically bound at the interface, may also play a role. This mechanism requires picosecond exciton harvesting time scales. (c) Ultrafast channel for exciton harvesting. Delocalization of the electron over multiple fullerene acceptors permits this faster separation of the nascent electron and hole. Coherent migration may permit excitons not initially located near the interface to participate in this process.

harvesting effect has also been demonstrated in bilayer devices using small conjugated organic molecules as electron donors.24 Friend et al.25 provided further insight into the mechanism of ultrafast exciton harvesting in OPV model systems, showing that photoexcitation of the electron donor can lead to an average electron−hole separation of at least 4 nm within 40 fs of photoexcitation. Delocalization of the transferred electron over several fullerene derivatives in the electron acceptor domain was proposed as a mechanism for this ultrafast charge separation (Figure 1c), a hypothesis that is supported by computational modeling of fullerene-based aggregates.26−29 This new mechanism for exciton harvesting necessitates reconsideration of exciton migration in bulk heterojunction OPVs. The incoherent, diffusive mechanism for exciton transport in electron donor domains does not support migration across nanometer length scales within 100 fs, such that only excitons generated within a monolayer or two of the interface could be harvested on ultrafast time scales. A back-ofthe-envelope consideration of the morphology and dimensions of a typical electron donor domain suggests that only a small percentage of photoexcitations should occur at the material B

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

along the lateral axis is caused by packing of the hexyl side chains. The non-90° crystal angle γ causes a small offset along the normal axis for sites with neighboring positions along the lateral axis and the same position along the longitudinal axis, as shown in Figure 2. A Hamiltonian representing the single-exciton manifold of each model domain was constructed according to eq 1.

design of optimal electron donor morphologies in new bulk heterojunction devices that optimize ultrafast charge separation.

2. THEORETICAL METHODS 2.1. Building Model Electron Donor Domains. We constructed model electron donor domains as nanometer-scale, three-dimensional arrays of chromophore sites, as shown in Figure 2. To keep the model as simple as possible, no

N

H=

∑ [ε + ϕm]|m⟩⟨m|+ ∑ [JTS + JTB]|m⟩⟨n| m=1

m≠n

(1)

The sums in eq 1 run over the set of chromophores in the model domain. The mean scalar site energy was chosen to be ε = 2.25 eV, based on the red-shifted absorption feature assigned to the crystalline region of regioregular P3HT thin films.45,46 In real polymers, many forms of noise modulate these site energies, including variations in the size, position, relative orientation, and specific molecular geometry of each polymer chromophore. We apply a simplified model to account for these effects, adding randomized Gaussian noise ϕm to the individual site energies. Noise distributions with standard deviations ranging from 0.01 to 1.0 eV were employed, and we refer to this value as the site energy noise magnitude. Our simplified noise model is intended to provide a reasonable estimate of the noise present in the system, without substantially increasing the complexity of our simple model. Couplings between the chromophores are a sum of two components: the through-space coupling (JTS) and the through-bond coupling (JTB). To calculate through-space coupling all chromophores were assigned identical transition dipole vectors. To vary the magnitude of the through-space coupling we varied the transition dipole strength from 5 to 15 D, in agreement with experimental determination of oligothiophene transition dipole magnitudes.47 Through-space couplings were first calculated as pairwise point dipole couplings between the chromophores. These initial calculations result in unphysically strong through-space couplings for chromophores adjacent along the normal axis of the model, due to the breakdown of the point dipole approximation at short interchromophore distances. To get more accurate through-space couplings, we scaled the point dipole couplings according to the work of Beenken and Pullerits48 to approximate the result of a distributed dipole calculation (see Supporting Information). After this adjustment, the largest through-space coupling for models employing a 15 D transition dipole vector was JTS = 69.4 meV and the largest through-space coupling for 5 D transition dipoles was 7.7 meV. This range of couplings represents a reasonable order of magnitude estimate of the appropriate nearest-neighbor couplings, with accuracy increasing at larger interchromophore distances. We evaluate our results in this context. An atomistic calculation of these couplings would undoubtedly improve the quantitative accuracy of our results but would be computationally intractable for the large number of chromphores needed to build and simulate dynamics in electron donor domains of realistic size. The through-bond coupling JTB was calculated by modeling the polymer strand as a donor−bridge−acceptor system employing a thiophene bridge.49,50 Chromophores that are not on the same polymer chain, i.e., not covalently bonded, are assigned a through-bond coupling of JTB = 0. Through-bond coupling for chromophores that are on the same polymer chain was calculated according to eq 2

Figure 2. Nonatomistic model of a nanometer-scale electron donor domain. We construct a nonatomistic model of an electron donor domain using a three-dimensional array of chromophore sites. Each polymer strand is represented by a chain of chromophores arranged along the longitudinal axis of the model. A cartoon of the atomic structure of two linked thiophene octamers, shown in the upper left corner, demonstrates the orientation of the chromophores along a polymer chain. Front layer of the model, outlined in black, consists of three parallel, evenly spaced polymer chains. Additional layers of these polymer sheets, shown as semitransparent gray rectangles, are stacked in lamellar fashion along the normal axis of the model. Model domain diagrammed in this figure contains 891 chromophores, arranged in an 11 × 3 × 27 array that spans 30.950 nm × 3.260 nm × 10.010 nm. Three axes of the model are not perfectly orthogonal, due to the crystal angle γ = 87° between the longitudinal and the normal axes. Crystal angles are diagrammed in the bottom left corner of the figure.

representations of individual atoms or molecular structure were included. Chromophores were assigned spatial positions according to experimental measurements of Phase I regioregular P3HT,39 yielding 1.630 and 0.385 nm interchromophore spacings along the lateral and normal axes, respectively. We based our model on this structure because highly ordered Phase I crystallinity is enhanced by the annealing process common during P3HT-based OPV construction,40 and enhanced crystallinity is known to improve charge carrier mobility41 and device performance.42 Site spacing along the longitudinal axis was 3.095 nm, chosen to be consistent with the physical dimensions of thiopheme octamers, as well as experimental measurements of hole−polaron magnitude in P3HT thin films.43,44 Figure 2 demonstrates the relative orientation of the chromophores along the nonorthogonal lateral, longitudinal, and normal axes of the model. Polymer strands run down the longitudinal axis. In real Phase I regioregular P3HT, polymer strands are closely spaced along the normal axis due to pi stacking. The larger spacing between adjacent polymer strands C

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C max JTB = JTB exp[−β(r − r0)]

migration, we built site-basis density matrices representing excitation of a single chromophore at the exact center of our model domains. We then employed an infinite-temperature Haken−Strobl−Reineker strategy to model the time evolution of these localized initial states.36,52 We simulated exciton migration in the single exciton manifold by using a quantum master equation to propagate the density matrix in time through numerical solution of the Liouville−von Neumann equation, as shown in eq 4

(2)

The maximum through-bond coupling Jmax TB was assigned values between 12 and 124 meV, in agreement with frequency domain modeling of P3HT thin films.38 The decay constant was fixed at a value of β = 3.7 nm−1, in agreement with experimental measurements of energy transfer in thiophenebridged macromolecules.50 The minimum distance between chromophores was set to r0 = 0.143 nm−1, representing the average bond length between two thiophenes in P3HT. The nearest edge to edge distance between the two chromophores is represented by r in eq 2. Using this approach, adjacent chromphores along a polymer chain have through-bond couplings equal to Jmax TB . The large interchromophore spacing along the polymer chain dictates that the through-bond couplings between nonadjacent chromophores are negligible, although they were still included in the model. 2.2. Delocalization in Eigenstates of the Hamiltonian. We calculated the eigenstates of a given model domain by diagonalizing its Hamiltonian. To examine the spatial extent of exciton delocalization for each eigenstate, an eigenbasis density matrix corresponding to excitation of that eigenstate was constructed and then transformed into the site basis using the transformation matrix from diagonalization. These site-basis representations of the eigenstates generally have population at multiple chromophores and also nonzero off-diagonal elements indicating superpositions of population at the two corresponding chromphores. These spatially delocalized population distributions are stationary with respect to time, due to the static nature of our model Hamiltonian. In contrast, dynamic fluctuation of the Hamiltonian can drive exciton migration through space in a real material excited to an eigenstate. This critical difference between the experimentally measured systems and our model deserves careful consideration, and our results are analyzed in this context. Previous simulations of coherent exciton migration in photosynthetic pigment−protein complexes have enacted dynamically fluctuating Hamiltonians through electron−phonon couplings, but these studies employ a much smaller number of chromophoric sites.32,33 Implementing a dynamic Hamiltonian for our large model domains would be computationally costly, beyond the scope of this work. As a quantitative metric for exciton delocalization, we calculated the participation ratio (PR) of the site-basis density matrix corresponding to each eigenstate of each model domain according to eq 3

PR =

d i ρ = − (Hρ − ρH ) − L(ρ) dt ℏ

where ρ is the site-basis density matrix and H is the site-basis Hamiltonian described in section 2.1. To simulate exciton migration, the site-basis density matrix was propagated in time according to eq 4 using the fourth-order Runge−Kutta method and a time step of 0.01 fs. To quantify exciton migration, we examine the time dependence of the PR calculated from the site-basis density matrix. The first term in eq 4 describes the evolution of a quantum system in the absence of dissipation such as damping of coherent processes through dynamic perturbations to the Hamiltonian, exciton recombination, or exciton harvesting. These dissipative processes are modeled in the second term of eq 4 through the Lindblad operator L.30,31 The effect the Lindblad operator on the elements of the derivative of the d density matrix, dt ρmn , is given in eq 5 [L(p)]mn = γ(1‐δmn)ρmn

(5)

Exponential damping of coherent processes was modeled through a single global dephasing parameter γ with values ranging from 300 to 1 × 106 cm−1. The dephasing parameter γ is an extremely simple model of dynamic variations in system− bath coupling, which cause coherent processes to dephase via exponential damping of the off-diagonal density matrix elements. The frequently employed value of 3000 cm−1 is consistent with system−bath coupling for a room-temperature bath with an Ohmic spectral density in the continuum limit and a system reorganization energy of a few hundred wavenumbers.30 Literature estimates of reorganization energy in P3HT range from hundreds to thousands of wavenumbers.53−55 Characterizing the mechanism for exciton migration required an additional metric beyond PR. The site-basis density matrix is inherently an ensemble-based description of excitation in the model system, and this fact can obscure the difference between the coherent migration we are interested in and the conventional diffusive transport mechanism. PR is sensitive only to the diagonal elements of the site-basis density matrix, which each represent exciton population at the corresponding chromophore in the model. Changes in these populations demonstrate exciton migration but do not report on the mechanism. The off-diagonal elements of the site-basis density matrix report on the behavior of superpositions of excitations at the two corresponding chromophores and are a better indicator of the nature of the observed transport. States supporting superpositions of excitation at multiple chromophores will have nonzero off-diagonal site-basis density matrix elements. These states support true quantum mechanical delocalization of population through space, facilitating coherent energy transport. Delocalized population distributions observed in the

1 N

∑i = 1 ρii2

(4)

(3)

where N is the number of sites in the model and ρii are the diagonal elements of the site basis density matrix. PR values can range from one, for complete localization on a single site, to a value of N, for excitation equally dispersed across all sites. Participation ratios and inverse participation ratios (calculated as the mathematical inverse of eq 1) are frequently employed to describe the spreading of exciton population, but the definitions are sometimes reversed in the literature. Scholes et al. recently presented a thorough analysis of various methods for quantifying delocalization and entanglement in quantum systems represented by density matrices, and we follow their lead by employing the original historical definition of a participation ratio here.51 2.3. Simulation of Ultrafast Exciton Migration. To characterize the mechanism and spatial pattern of exciton D

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

Figure 3. Properties of the eigenstates of model Hamiltonians with varying site energy noise magnitude. Examples of the participation ratio for each eigenstate (left) and the normalized, averaged absorption spectra (right) are shown for model Hamiltonians with site energy noise magnitudes of 10 (red), 30 (blue), 60 (green), 100 (black), 300 (turquoise), and 600 meV (brown). Participation ratios (left) are plotted as a function of an eigenstate index that arranges the 891 eigenstates in order of increasing energy. To simulate the absorption spectra (right), the transition dipoles for each eigenstate were squared to produce stick spectra, which were then convoluted with a Gaussian function employing a variance of 5 meV. Spectra shown in the figure were each summed over the five implementations of each site energy noise magnitude and then normalized to facilitate comparison.

interfacial sites. At these interfacial sites, the exciton harvesting rate ΓEH is set to a value of 2.5 × 1013 excitons harvested per second, or one exciton harvested every 40 fs, in agreement with Friend et al.25 2.5. Limitations of the Model. The model described here can characterize trends in ultrafast exciton delocalization but is not intended to produce quantitatively accurate predictions of exciton migration in a real OPV device. Real devices contain a mixture of crystalline and noncrystalline donor domains, mixed to varying degrees with aggregates of acceptor molecules that have their own variety of morphologies. Full modeling of this complex situation is beyond the scope of this work. To facilitate dynamic modeling in domains of realistic size, several additional simplifications were enacted. The model is not atomistic in nature, so it does not contain vibrational states and polaron formation is impossible. Polarons do form in OPVs on ultrafast time scales and may affect exciton transport and harvesting.18−23 The quantitative accuracies of our predictions of exciton migration and harvesting are dependent on the accuracy of the interchromophore couplings. The calculated through-space couplings we employ are an order of magnitude estimate of the appropriate coupling for nearest-neighbor chromophores. Quantitative descriptions of exciton migration are presented in this context. Dissipation of exciton population through recombination is also neglected in our model. Exciton−exciton annihilation events are impossible because we model only the single-exciton manifold. Single-exciton recombination occurs on much longer time scales than considered here56 and is excluded from the model. Exciton dissipation through exciton harvesting and damping of coherences via dephasing are included in the model but are implemented in the broadest possible sense using the global parameters γ and ΓEH. To preserve the simplicity and intuitive nature of our model, we also choose to inspect the behavior of straightforward states of our model system. Full modeling of the wavepacket of eigenstates excited by a broad-band ultrafast laser pulse is a complex task beyond the scope of this work. Time propagation of a single eigenstate is also computationally complex beyond the scope of this work, requiring implementation of a dynamic

absence of superposition character represent statistical mixtures of localized populations. Exciton migration in such a system proceeds through the incoherent Förster mechanism. Several related metrics for the “quantumness” of the state of a system have been applied in recent literature on exciton migration, such as von Neumann entropy, relative entropy, coherence length, and total tangle. Implementation and analysis of these metrics in systems described with density matrices were recently discussed in detail by Scholes et al.51 In this work we choose to employ purity as a metric because it has a wellaccepted mathematical definition and intuitive limits. More complex concurrence metrics for entanglement can provide additional information in more complex models,33 but for our simple model purity is sufficient. Purity is calculated as shown in eq 6

purity = Tr(ρ2 )

(6)

Purity is a basis set invariant quantity sensitive to the persistence of the off-diagonal elements of the density matrix. A value of one indicates pure quantum state, i.e., complete quantum mechanical delocalization of population. A value of zero indicates a completely mixed quantum state, i.e., a statistical mixture of localized populations with no superposition character. 2.4. Modeling Ultrafast Exciton Harvesting. Exciton harvesting is expected to occur at the material interface between the electron donor and the electron acceptor materials. To model this behavior, we designate selected faces of the model donor domains as donor−acceptor interfaces. The sites making up these faces form the monolayer of chromophores that are immediately adjacent to the electron acceptor domain, which is not modeled explicitly. An electron harvesting channel is then implemented only at these interfacial sites, as shown in eq 7 face [L(ρ)]mn = γ(1 − δmn)ρmn − iℏΓEHδmn

(7)

The first term on the right side of eq 7 is identical to the Lindblad operator shown in eq 5. The pure-imaginary second term removes population from the diagonal elements of the Hamiltonian at a rate determined by the parameter ΓEH. The ΓEH parameter is assigned a non-zero value only for the E

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

To further examine the degree of delocalization present in the eigenstates, we calculated the average PR observed for each noise magnitude. For a given value of the noise magnitude, we calculated a simple arithmetic mean and standard deviation of the PR over the eigenstates of all five model domains created. We also calculated Boltzmann-weighted average PRs for each set of model domains at 300 K as well as an absorptionweighted average PR that employs the square of the transition dipole for each eigenstate as a weighting factor. These data are shown in Figure 4.

Hamiltonian that would dramatically increase computation time. We choose instead to consider simpler states of the model system.

3. RESULTS AND DISCUSSION 3.1. Delocalization of Eigenstates. We first investigated the extent of population delocalization observed in the eigenstates of the model Hamiltonians as a function of the site energy noise magnitude. The chromophores of each model domain were arranged as shown in Figure 2. The strongest through-space interchromophore coupling was JTS = 69.4 meV, and the strongest through-bond coupling was JTB = 40 meV. Thirty model Hamiltonians were created, five each for site energy noise magnitudes 10, 30, 60, 100, 300, and 600 meV. Each Hamiltonian was diagonalized, producing the set of sitebasis eigenvectors for the model as a transformation matrix, along with their associated eigenenergies. To quantify the spatial delocalization of population observed in each eigenstate, PRs for each eigenstate of all 30 models were calculated according to eq 3. In Figure 3a we display the calculated PR values for six example model domains, each employing a different site energy noise magnitude. The PRs are plotted against an eigenstate index that arranges the eigenstates in order of increasing eigenenergy to best display trends. PRs calculated from the additional 24 model Hamiltonians follow similar patterns with minor quantitative variation due to the randomized noise. The spatial properties of the eigenstates derived from a given model Hamiltonian are highly varied but do exhibit some clear trends. For all site energy noise magnitudes, PR increases nonmonotonically with eigenenergy. Models with site energy noise magnitude greater than the largest interchromophore couplings exhibit PRs near one for low-energy states. In these cases the site energy noise strongly inhibits delocalization of the eigenstates. When site energy noise magnitude is much greater than coupling strength, even high-energy eigenstates are not substantially delocalized. Models with smaller noise magnitude have much larger PRs, even at low energies. In this case the site energy noise cannot inhibit delocalization, and even the lowest energy eigenstates are spread across several chromophores. Normalized, averaged absorption spectra for model domains employing each of the six site energy noise magnitudes employed are shown in Figure 3b. Calculation of these spectra is described in the Supporting Information. The spectra show the expected trend that a larger magnitude of energetic noise results in a broader absorption spectrum. When noise magnitude is greater than the largest interchromophore couplings, the noise strongly broadens the spectrum. When the noise magnitude is equivalent to or smaller than the couplings, broadening due to the couplings still enforces a minimum width for the spectrum. The observed minimum spectral widths of ∼400 meV (full width at half-maximum) are somewhat smaller than those observed experimentally,45,46 suggesting that site energy noise magnitude does play a role in broadening absorption in real P3HT thin films or that coupling is stronger in real materials than in the model. We speculate that the former is more likely. The fine structure observed in our simulated spectra is a trivial consequence of the site energy noise, because our model does not contain vibrational states. Absorption spectra of P3HT thin films have previously been modeled in the frequency domain,37,38 demonstrating that an information-rich vibrational structure emerges when vibrational levels are added to the model.

Figure 4. Mean participation ratio as a function of site energy noise magnitude. Arithmetic mean participation ratio (black squares), Boltzmann-weighted mean participation ratio (red circles), and absorption-weighted mean participation ratio (green triangles) are plotted as a function of the site energy noise magnitude. One standard deviation error bars are shown for the arithmetic mean. In each case an exponential dependence is observed, with color-coded best fit lines shown in the figure. Decay constants are 68 (±3) meV−1 for both the arithmetic and the absorption-weighted means. The Boltzmannweighted means yield a decay constant of 22.5 (±0.4) meV−1. The magnitude of the strongest interchromphore coupling, 69 meV, is marked with a dashed gray line. To demonstrate the spatial distribution of population at different noise levels, isosurface visualizations of the brightest eigenstate of two model domains are shown as insets. Isosurfaces were created using a population-weighted three-dimensional Gaussian basis, as described in the Supporting Information. The isosurface value employed here is 0.1. Semitransparent gray lines define the edges of the model domains for each isosurface as a guide to the eye. Each model domain contained 891 chromophores arranged in an 11 × 3 × 27 array (longitudinal × lateral × normal or vertical × horizontal × depth as shown).

All versions of the mean PR demonstrate that delocalization decreases exponentially with the magnitude of the site energy noise. All average PRs indicate that a very high degree of delocalization should be present in the eigenstates of a system with site energy noise smaller than the largest interchromophore coupling. As noise magnitude becomes greater than the coupling, delocalization is rapidly reduced. Isosurfaces of two eigenstates shown in the figure demonstrate the general trend that population distributions delocalize most strongly across chromophores adjacent along the normal and longitudinal axes of the model. Absorption of a single photon excites a system to an eigenstate of the Hamiltonian. The exciton populations observed in ultrafast spectroscopic experiments are created by ultrafast laser pulses, which have a finite bandwidth and create a wavepacket of eigenstates. In Figures 3 and 4, we analyze simpler, infinite time eigenstates of the Hamiltonian. While these states are not a perfect model of the exciton population observed in ultrafast spectroscopic excitations, they are a good F

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

real devices. The mechanism for transport will be determined by the dephasing of coherent processes, included in the model through the global dephasing parameter γ. The magnitude of exciton migration experienced by a localized exciton using a static Hamiltonian will certainly be quantitatively different from that resulting from the migration undergone by an eigenstate with a dynamically fluctuating Hamiltonian. We propagated the density matrix in time through numerical solution of the Liouville−von Neumann equation, as described in section 2.3. Simulations were performed for models employing a variety of coupling strengths (JTB and JTS), site energy noise magnitudes, and dephasing parameters (γ). Each simulation was repeated five times with five different realizations of the randomized site energy noise. Exciton harvesting through the ΓEH parameter was not included. Dynamic delocalization during these simulations is summarized in Figure 5a, which shows the average PR increasing as a function of simulation time for several different sets of parameters chosen to highlight key trends. For moderate values of the couplings and energetic noise (red lines, JTS = 69 meV, JTB = 40 meV, noise magnitude = 100 meV), the average PR is observed to rapidly increase during the first few tens of femtoseconds. An isosurface visualization of the average population distribution after 40 fs of simulated migration is shown in Figure 5b (red box). This isosurface visualization was created by modeling the spatial distribution of exciton population using a basis of three-dimensional Gaussian functions. Full mathematical details of this procedure are presented in the Supporting Information. The isosurface demonstrates that exciton population migrates primarily along the normal and longitudinal axes of the model, spreading across a greater number of chromophores in the normal dimension but traveling farther through space in the longitudinal dimension. Very little population migrates along the lateral dimension of the model. This qualitative pattern of exciton migration is consistent with the relative magnitude of the nearest-neighbor transition dipole couplings along each axis, such that stronger couplings drive more delocalization. The observed spatial migration pattern is not dependent on the quantitative accuracy of the couplings, although stronger couplings do generate more migration. If the magnitude of the through-bond coupling is increased (black lines, JTS = 69 meV, JTB = 124 meV, noise magnitude = 100 meV), the average PR increases by nearly a factor of 2 and the standard deviation is reduced. These metrics indicate stronger and more uniform exciton delocalization. The corresponding isosurface visualization demonstrates that exciton population spreads more strongly across the longitudinal dimension of the model in this case. Reducing the through-space coupling (blue lines, JTS = 31 meV, JTB = 40 meV, noise magnitude = 100 meV) is observed to greatly reduce the average PR, and the accompanying isosurface shows that this reduction occurs primarily along the normal axis of the model. Decreasing this coupling has a large effect on the PR because the exciton population accesses a large number of chromophores as it travels along the closely packed normal axis of the model. Increasing the site energy noise magnitude (green lines, JTS = 69 meV, JTB = 40 meV, noise magnitude = 300 meV) reduces the average PR, and the magnitude of this effect increases with simulation time. The accompanying isosurface demonstrates that the change in the spatial distribution is isotropic and relatively small when compared to changing the couplings by a similar factor.

representation of the spatial distribution of an exciton created via photon absorption. The PR of initially excited eigenstates in spectroscopic experiments is best modeled by the absorptionweighted metric, as this average PR accounts for the probability of absorption for each eigenstate. The Boltzmann-weighted average PR represents a population of excitons at thermal equilibrium, and the values observed are substantially lower than the arithmetic and absorption-weighted mean PRs. The lower values of the Boltzmann-weighted average PR indicate that in the “small noise” regime, excitons created via photon absorption can be expected to undergo substantial localization as they relax to thermal equilibrium. Such localization is often observed in experiments but is usually assigned to polaron formation. Our model does not include vibrational states, so true polaron formation is impossible. Figure 4 demonstrates that delocalization across nanometer length scales can be achieved in the eigenstates of model domains and that this effect is stronger if the site energy noise magnitude is less than the strongest interchromophore couplings. The amount of delocalization is observed to be highly dependent on the relative magnitudes of the noise and coupling. The wavepacket nature of the excited state in ultrafast experiments and dynamic behavior of the eigenstates, neither of which are considered here, should be expected to increase delocalization. The amount of delocalization observed in our model thus represents a lower bound to the amount of delocalization that should be expected during an ultrafast spectroscopic experiment. The positive correlation between PR and eigenenergy indicates that higher energy excitation could result in larger excitons. These larger excitons would have a greater likelihood of accessing the donor/acceptor material interface, leading to enhanced ultrafast exciton harvesting yield. However, the current data do not prove this concept as a general phenomenon, and full exploration must be left as a topic for future consideration. Note that this effect is distinct from the hot exciton effect frequently mentioned in OPV literature, which is mediated through vibrational states.57−59 Hot excitons cannot be created in our model, because the model does not contain vibrational states. 3.2. Ultrafast Exciton Migration. The stationary eigenstates discussed in section 3.1 provide a clear and intuitive starting point from which to consider delocalization of an exciton but do not account for the dynamic behavior of excitons. To characterize dynamic properties we examine the migration of excitons initially localized on a single chromophore in the center of our model domains. We choose this location for the initial excitation to provide the clearest view of the behavior of excitations traveling toward the donor− acceptor material interface, represented by the faces of our model domains. We are interested in the spatial patterns and mechanism of exciton migration as the excitation moves toward this interface. Our initially localized excitations are in some sense unphysical, because photoexcitation in real devices is unlikely to create excitons completely localized on a single chromophore. A single photon should excite an eigenstate of the Hamiltonian. As discussed in section 3.1, focusing on localized initial excitations allows us to access dynamic behavior while minimizing computation time. Given the simplicity of our nonatomistic modeling strategy, modeling more complex initial excitations may not yield more valid quantitative results. The spatial patterns of migration will be largely determined by the couplings between chromophores, both in the model and in G

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

red and green lines in Figure 5a demonstrates that increasing site energy noise by a factor of 3 results in a relatively small change in exciton migration. Our model does not provide a framework for consideration of the effects of the altered lattice packing due to a distribution of chromophore sizes, so this topic must be left as topic for future work. The general spatial pattern of exciton migration observed in Figure 5 and the dependence of migration on the couplings and noise magnitude are both consistent with conventional expectations for energy transport in molecular aggregates. The through-space coupling most strongly controls transport across the normal dimension of the model, while through-bond coupling drives transport along the longitudinal dimension. Transport is slowest across the lateral dimension of the model, where coupling is weakest. Site energy noise is observed to inhibit but not eliminate transport. Exciton population accesses more chromophores as it moves along the normal dimension than it does along the longitudinal dimension but travels farther in space along the longitudinal dimension due to the greater interchromophore spacing. To characterize the mechanism of ultrafast migration, we calculated the purity of the site basis density matrix as a function of simulation time. Average purity values for each model are plotted in Figure 5a. As discussed in section 2.3, purity is a basis set invariant quantity that is sensitive to the persistence of the off-diagonal elements of the density matrix. A larger value of the purity metric indicates that more of the population described by the density matrix is in delocalized superposition states. When superpositions are present, coherent energy transport can occur. As the purity of the system decreases, superpositions collapse into a statistical mixture of localized populations and incoherent, diffusive energy transport mechanisms take over. Figure 5a shows that for all models considered, the quantum state of the system retains substantial purity for a few tens of femtoseconds. In all cases, purity rapidly decreases by at least an order of magnitude during the first 40 fs of exciton migration. During this time period, exciton population rapidly spreads through the model domains, creating the population distributions visualized in Figure 5b. Purity exhibits trends opposite to those found for PR. Increasing couplings drive faster transport, but purity is reduced more quickly. Increasing noise inhibits transport, but purity decays more slowly. These results can be rationalized by considering the way that the dephasing parameter γ drives reduction of the purity metric. As population spreads, coherences are spawned in the density matrix. Exponential decay of each coherence as a function of time results in a faster decaying purity for systems that spread to more sites. Dependence of exciton migration on the value of the dephasing parameter is explored in the Supporting Information, demonstrating that delocalization is maximized in a “Goldilocks” regime around γ = 3000 cm−1. This type of behavior has been previously observed in models of photosynthetic pigment−protein complexes.30 When the dephasing parameter is increased to 1 × 106 cm−1, transport is almost nonexistent on ultrafast time scales, yielding a PR value of 1.1 after 95 fs for moderate coupling strengths and noise magnitude. In this case the large value of γ destroys coherences so rapidly that the purely incoherent limit of the model is achieved. The fact that almost no transport is observed in this case is further evidence that the exciton migration shown in Figure 5 is coherent.

Figure 5. Coherent exciton migration on ultrafast time scales. (a) Average participation ratio (solid lines, right vertical axis) and average quantum mechanical purity (dashed lines, left vertical axis) of four sets of model domains are shown as a function of simulation time. Colorcoded, shaded regions show one standard deviation error bars. Red lines show a model employing moderate couplings and site energy noise. Black lines show the effects of increasing through-bond coupling strength while holding other parameters constant. Green lines show the effects employing the moderate couplings and enhancing the noise magnitude. Blue lines show the effects of decreasing through-space coupling strength. (b) Isosurfaces visualizing the average spatial distribution of exciton population after 40 fs of simulated exciton migration are shown for each model described in part a. Isosurfaces were created using a population-weighted three-dimensional Gaussian basis, as described in the Supporting Information. Isosurface value employed here was 0.1. Semitransparent gray lines define the edges of the model domains as a guide to the eye. Chromophores were arranged within each model domain as diagrammed in Figure 2, and the value of the global dephasing parameter was γ = 3000 cm−1.

Quantitative measures of exciton migration in these simulations, such as the PR, are best interpreted as an upper bound for the amount of migration expected in a real electron donor domain. The moderate through-space coupling model employed in Figure 5 is expected to be an overestimate for short intersite distances, which speeds up exciton migration. As the figure shows, reducing through-space coupling has a strong effect on exciton migration along the normal axis of the model. Excitons created via photon absorption will also not be localized on a single chromophore, further slowing the type of migration demonstrated here. Real crystalline P3HT will also not consist of a perfect lattice of thiophene octamers but rather will contain chromophores of varying size and more irregular interchromophore spacing. Varying chromophore size can be directly pictured within the framework of the current model as an increase in the site energy noise magnitude. Comparing the H

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

at 100 fs to 12:1 by 200 fs. By 500 fs, this ratio has fallen to 1.4:1, indicating that the population distribution throughout the model domain is nearly uniform. By 950 fs this ratio falls further to a value of 1.02:1. In summary, employing moderate values for couplings, site energy noise magnitude, and the dephasing parameter results in localized initial excitons migrating over nanometer length scales on a time scale of tens of femtoseconds through a coherent mechanism. The average population distribution established by coherent migration during the first 100 fs extends primarily across a spheroid with a major semiaxis extending ∼12.4 nm along the longitudinal axis and a minor semiaxis extending ∼4.6 nm along the normal axis. Exciton migration across the lateral axes of the model is less strong on ultrafast time scales. On longer time scales population spreads roughly isotropically through the model domain, reaching a roughly uniform distribution within 0.5 ps. This slower migration is due to incoherent, diffusive transport. In real materials, delocalization across at least a few neighboring chromophores along the normal and longitudinal dimensions should be expected. Longer range exciton migration over nanometer length scales is supported by the calculations, but we caution that due to the simplicity of the model a quantitative estimate of the expected range, such as a diffusion constant, would be speculative. The spatial pattern of exciton migration is controlled by the relative magnitudes of the through-space and through-bond coupling and not affected by changes in the magnitude of the coupling, the site energy noise magnitude, or the global dephasing parameter. The spatial extent of the coherent delocalization observed over the first few tens of femtoseconds is strongly dependent on the coupling strength. The model clearly demonstrates that any ultrafast migration should be expected to occur both along the pistacked dimension of crystalline polymer aggregates and along the polymer chain itself and that this ultrafast migration will occur through coherent mechanisms. Enhanced crystallinity has been experimentally correlated with enhanced performance in P3HT-bsed OPV devices,60 and our results suggest that coherent exciton transport processes may play a role in this enhancement. 3.3. Ultrafast Exciton Harvesting. The exciton migration exhibited in Figure 5 demonstrates that, on ultrafast time scales, exciton population can migrate across nanometer lengths across both the normal and the longitudinal axes of the model. To investigate the dependence of ultrafast exciton harvesting yield on domain size, we created a series of models that vary the physical extent of the model domains in these dimensions. Models with varying length along the normal dimension were created by arranging five-by-five grids of sites along the lateral and longitudinal dimensions of the model, forming 25-site layers that can then be stacked along the normal axis as shown in Figure 7a. Exciton migration and harvesting was simulated in model electron donor domains extending nearly 10 nm across the normal axis dimension, a common size regime for bulk heterojunction OPVs. We also employed an analogous strategy to vary the physical extent of the model along the longitudinal axis. Due to larger interchromophore spacing along this axis, the largest of these models extended over 60 nm. Five instances of each size model domain were constructed, each with a unique site energy noise profile. For all simulations, the model parameters employed were JTS = 69 meV, JTB = 124 meV, noise magnitude = 100 meV, and γ = 3000 cm−1. Exciton migration and harvesting were simulated for 100 fs in each model domain.

Figure 6 shows results of continued exciton migration through the model donor domain beyond the first few tens of

Figure 6. Incoherent exciton migration on picosecond time scales. Average participation ratio (solid line, right vertical axis) and quantum mechanical purity (dashed line, left vertical axis) are plotted against a common simulation time. One standard deviation error bars are indistinguishable from the line width in this figure. Isosurfaces visualizing the spatial distribution of exciton population are shown for 100, 150, 200, and 500 fs of simulated exciton migration. PR increases rapidly toward the upper limit of 891 as exciton population spreads throughout the model. Coherent, ultrafast delocalization establishes a population distributed across a few neighboring sites along the normal and longitudinal axes, as shown by the 100 fs isosurface. Subsequently, the population spreads roughly isotropically in all directions. By 500 fs, population is nearly evenly distributed through all sites in the model. Purity decays from an initial value of 1 to a value near 0 by 50 fs, demonstrating the incoherent nature of the population migration shown. Isosurfaces were created using a population-weighted threedimensional Gaussian basis, as described in the Supporting Information. Isosurface value employed here was 0.3.

femtoseconds. The model parameters employed were JTS = 69 meV, JTB = 124 meV, noise magnitude = 100 meV, and γ = 3000 cm−1. The average PR and purity of the model system are plotted as a function of simulated time. Visualizations of the spatial distribution of exciton population at 100, 150, 200, and 500 fs are also shown as insets. By 100 fs, the population distribution established by coherent processes has completely developed. As demonstrated by the 100 fs isosurface shown in Figure 6, population is distributed most strongly across the sites adjacent to the initial excitation along the normal and longitudinal axes of the model. The initially excited state hosts the largest single site population, 1.3 (±0.4)%. The population falls monotonically with distance in all directions. Along the normal and longitudinal axes, population remains above 0.5% for sites within 2.310 and 6.190 nm of the initial excitation, respectively. Along the lateral axis the two sites immediately adjacent to the initial excitation host a population of 0.33 (±0.01)% and 0.327 (±0.004)%. The smallest populations in the model on are on the order of 10−3%. As simulation time increases, incoherent exciton migration drives the model system toward a uniform distribution, i.e., populations of 0.11% at all chromophores. The isosurfaces visualizing population distributions at 150 and 200 fs show that this migration occurs roughly isotropically through space. The elongated shape of the distribution created by coherent migration is preserved only because the model is smallest along the lateral axis. By 200 fs all site populations have fallen below 0.4%, and substantial populations have reached the lateral edges of the model. The ratio of the average populations at the most and least populated sites falls from a value of 223:1 I

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

separation rate in the first monolayer adjacent to the electron acceptor is an order of magnitude faster than the next closest monolayer.61 Computational modeling also indicates that exciton dissociation in OPVs should only be possible for exciton population near the interface.62 In each simulation, a localized excitation was placed in the chromophore at the center of the model. The average exciton populations harvested after 90 fs of simulation time are plotted in Figure 7b as a function of the size of the domains along the normal (green) and longitudinal (black) axes. The smallest domains were single five-by-five grids of chromophores, corresponding to zero length. In these models the exciton harvesting channel was active at all chromophores, and consequently, the harvesting yield was nearly 100%. In larger domains the initially excited chromophore did not have an active harvesting channel. As simulation time increased, population migrated to the edges of the model and was harvested at the donor/acceptor interfaces. Average ultrafast exciton harvesting yields decrease as model domains grow larger but remain above 10% for models extending 3.850 and 18.570 nm across the normal and longitudinal dimensions, respectively. Visualizations of the average exciton population distributions after 90 fs are shown as insets in Figure 7b. The isosurface visualizations shown are based on a normalized spatial function, as described in the Supporting Information. As a consequence, the isosurfaces show only the spatial distribution of the exciton population remaining in the system and do not report on the total exciton population remaining in the system. Creating the visualizations this way facilitates straightforward comparison of spatial distributions in both small and large model systems, which might otherwise be obscured by the larger fraction of exciton population that is harvested in smaller model domains. Isosurfaces of the population distribution remaining after 90 fs are shown as insets for several model domain shapes. Seven layers arranged along the normal axis and harvesting active at the extremes of this axis (tan boxes) results in an exciton harvesting yield of 40 (±1)% after 90 fs. The isosurface for this model shows that substantial population reaches many interfacial sites. A similar model extending across 21 layers (red boxes) harvests only 3.8 (±0.4)% of the exciton population within 90 fs. In this case the isosurface reveals that very little population can access the interfaces on ultrafast time scales. A similar effect is observed when the size of the model is varied along the longitudinal axis, with the smaller model (3 layers, purple boxes) showing substantial population at several interfacial sites, while in the larger model (11 layers, cyan boxes) the delocalized population does not access these sites. In this case the small and large models harvested 78.6 (±0.6)% and 0.79 (±0.04)% of the exciton population by 90 fs, respectively. These results demonstrate that smaller electron donor domains support substantially more ultrafast exciton harvesting than larger domains, and that ultrafast harvesting yield can be very high for domains with small spatial extent along the normal or longitudinal axes. In the smallest domains, excitons delocalize to the interface and ultrafast harvesting is limited by the harvesting rate at the interfaces rather than by migration. As domain size increases, exciton migration begins to limit harvesting. Ultimately the domains become so large that substantial exciton population cannot reach the interfacial layers within 90 fs, and ultrafast harvesting yields drop to nearly zero. In the size regime where substantial exciton harvesting

Figure 7. Ultrafast exciton harvesting yields in model electron donor domains of varying size. (a) Diagram of the model electron donor domains used to investigate the dependence of ultrafast exciton harvesting yield on domain size along the normal dimension. Exciton harvesting was enacted as a trapping channel active only at sites on the faces of the model at either end of the normal axis (highlighted in green). An analogous approach was used to find the dependence of harvesting yield on the size of the model along the longitudinal dimension. (b) Average harvested populations after 90 fs of simulated exciton migration and harvesting are plotted as a function of the size of the model domains along the normal axis (green) and longitudinal axis (black). One standard deviation error bars are indistinguishable from the line width. Isosurface visualizations of the average spatial distribution of exciton population remaining after 90 fs are shown as insets for 7 and 21 layers along the normal axis (tan and red) and for 3 and 11 layers along the longitudinal axis (purple and cyan). In all cases each layer of each model is a 5 × 5 grid of chromophores. Isosurfaces are created using a population-weighted three-dimensional Gaussian basis, as described in the Supporting Information. Note that the function used to create the isosurfaces is normalized prior to construction of the isosurfaces, so the visualized surfaces are not distorted by exciton harvesting. Isosurface value employed here was 0.3.

For this set of simulations, the two faces of the model domain at the extreme ends the axis being varied were designated as donor/acceptor interface layers. Exciton harvesting was enacted at the sites forming these faces. The computational mechanics of adding an exciton harvesting channel to the model are outlined in section 2.4. Because the exciton harvesting channel was only active at the chromophores designated as forming the donor/acceptor interface layers, only exciton population resident at these chromophores could be harvested. A recent ultrafast spectroscopic study of charge separation at OPV interfaces demonstrated that the charge J

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C



yield is observed, dependencies on parameters of the model such as coupling strength and noise magnitude will be the same as those observed in Figure 5. Enhancing through-space coupling will primarily drive more migration along the normal axis, enhancing harvesting in models that harvest at the extremes of this axis. Enhancing through-bond coupling will have a similar effect along the longitudinal axis. Larger site energy noise magnitude will reduce harvesting in both directions. Increasing or decreasing the dephasing parameter away from its “Goldilocks” value will also reduce harvesting. These results suggest that formation of smaller donor domains may in fact be quite important in promoting ultrafast exciton harvesting. Excitations in crystalline P3HT can be expected to spread across both pi-stacked chromophores (normal dimension) and down the polymer chains (longitudinal dimension). This transport will enhance access to the donor−acceptor material interface and increase ultrafast exciton harvesting yield. This model predicts that devices containing small electron donor domains dispersed within an electron acceptor substrate would maximize ultrafast harvesting.

Article

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcc.6b11332. Full mathematical details of model hamilonian, isosurface visualization, and simulated absorption spectra; additional data on ultrafast delocalization under varying conditions; table of simulations performed (PDF)



AUTHOR INFORMATION

Corresponding Author

*Phone (773)325-7099. Fax: (773)325-7421. E-mail: ggriffi6@ depaul.edu. ORCID

Kenley M. Pelzer: 0000-0001-8984-758X Gregory S. Engel: 0000-0002-6740-5243 Graham B. Griffin: 0000-0001-6586-6284 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS G.B.G. acknowledges support from the DePaul University CSH FSRG program. K.M.P. acknowledges support of the Aneesur Rahman Named Fellowship of Argonne National Laboratory. G.S.E. acknowledges funding research provided by Qatar National Research Foundation (QNRF), NPRP X-107-010027.

4. CONCLUSIONS Our simulations of exciton delocalization, migration, and harvesting show that excitons in highly ordered aggregates of conjugated organic electron donors such has crystalline P3HT may exhibit ultrafast mobility on nanometer length scales. The spatial distribution of exciton population evolves according to a coherent mechanism on ultrafast time scales, and the resulting distribution of population is determined by the morphology of the electron donor material. Eigenstates of the model Hamiltonian are frequently delocalized over many chromophores, and the degree of delocalization is highly dependent on the magnitude of energetic noise in the model. When exciton population is initially localized on a single-polymer chromophore, population coherently delocalizes over nanometer length scales within a few tens of femtoseconds. This effect is most prominent between pi-stacked chromophores and down polymer chains, with the pi-stacking direction accessing more chromophores but down-chain migration traveling farther through space. On longer time scales, exciton population spreads incoherently throughout the model domains in a roughly isotropic fashion. Ultrafast exciton harvesting yield is shown to be strongly dependent on the size and shape of electron donor domains, with smaller domains facilitating the largest exciton harvesting yields. These results can inform design of new bulk heterojunction OPV devices that seek to optimize ultrafast exciton harvesting. We have shown that domain size and morphology control ultrafast exciton harvesting and that very high ultrafast exciton harvesting yields may be possible for devices with electron donor domains on the scale of tens of nanometers. Conjugated organic nanoplates extending only a few nanometers across the pi-stacking and longitudinal dimension may be highly desirable. Further experimental studies could verify this result by measuring ultrafast exciton yields in devices where the morphology of the bulk heterojunction blend is characterized not just by the general size of the electron donor domains but also by the particular size distribution across the pi-stacking dimension and along polymer chains.



REFERENCES

(1) Yu, G.; Gao, J.; Hummelen, J. C.; Wudl, F.; Heeger, A. J. Polymer Photovoltaic Cells - Enhanced Efficiencies Via a Network of Internal Donor-Acceptor Heterojunctions. Science 1995, 270, 1789−1791. (2) Brédas, J.-L.; Norton, J. E.; Cornil, J.; Coropceanu, V. Molecular Understanding of Organic Solar Cells: The Challenges. Acc. Chem. Res. 2009, 42, 1691−1699. (3) Günes, S.; Neugebauer, H.; Sariciftci, N. S. Conjugated PolymerBased Organic Solar Cells. Chem. Rev. 2007, 107, 1324−1338. (4) Gao, F.; Inganas, O. Charge Generation in Polymer-Fullerene Bulk-Heterojunction Solar Cells. Phys. Chem. Chem. Phys. 2014, 16, 20291−20304. (5) Fö rster, T. Zwischenmolekulare Energiewanderung Und Fluoreszenz. Ann. Phys. 1948, 437, 55−75. (6) Marcus, R. A. On the Theory of Oxidation-Reduction Reactions Involving Electron Transfer. I. J. Chem. Phys. 1956, 24, 966−978. (7) Grancini, G.; Polli, D.; Fazzi, D.; Cabanillas-Gonzalez, J.; Cerullo, G.; Lanzani, G. Transient Absorption Imaging of P3ht:Pcbm Photovoltaic Blend: Evidence for Interfacial Charge Transfer State. J. Phys. Chem. Lett. 2011, 2, 1099−1105. (8) Deotare, P. B.; et al. Nanoscale Transport of Charge-Transfer States in Organic Donor-Acceptor Blends. Nat. Mater. 2015, 14, 1130−1134. (9) Müller, J. G.; Lupton, J. M.; Feldmann, J.; Lemmer, U.; Scharber, M. C.; Sariciftci, N. S.; Brabec, C. J.; Scherf, U. Ultrafast Dynamics of Charge Carrier Photogeneration and Geminate Recombination in Conjugated Polymer:Fullerene Solar Cells. Phys. Rev. B: Condens. Matter Mater. Phys. 2005, 72, 195208. (10) Ohkita, H.; et al. Charge Carrier Formation in Polythiophene/ Fullerene Blend Films Studied by Transient Absorption Spectroscopy. J. Am. Chem. Soc. 2008, 130, 3030−3042. (11) Clarke, T. M.; Durrant, J. R. Charge Photogeneration in Organic Solar Cells. Chem. Rev. 2010, 110, 6736−6767. (12) Yusoff, A. R. b. M.; Kim, D.; Kim, H. P.; Shneider, F. K.; da Silva, W. J.; Jang, J. A High Efficiency Solution Processed Polymer Inverted Triple-Junction Solar Cell Exhibiting a Power Conversion Efficiency of 11.83%. Energy Environ. Sci. 2015, 8, 303−316. K

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C (13) Chen, C.-C.; Chang, W.-H.; Yoshimura, K.; Ohya, K.; You, J.; Gao, J.; Hong, Z.; Yang, Y. An Efficient Triple-Junction Polymer Solar Cell Having a Power Conversion Efficiency Exceeding 11%. Adv. Mater. 2014, 26, 5670−5677. (14) Tong, M.; Coates, N. E.; Moses, D.; Heeger, A. J.; Beaupré, S.; Leclerc, M. Charge Carrier Photogeneration and Decay Dynamics in the Poly(2,7-Carbazole) Copolymer Pcdtbt and in Bulk Heterojunction Composites with Pcbm(70). Phys. Rev. B: Condens. Matter Mater. Phys. 2010, 81, 125210. (15) Kaake, L. G.; Moses, D.; Heeger, A. J. Coherence and Uncertainty in Nanostructured Organic Photovoltaics. J. Phys. Chem. Lett. 2013, 4, 2264−2268. (16) Provencher, F.; Bérubé, N.; Parker, A. W.; Greetham, G. M.; Towrie, M.; Hellmann, C.; Côté, M.; Stingelin, N.; Silva, C.; Hayes, S. C. Direct Observation of Ultrafast Long-Range Charge Separation at Polymer−Fullerene Heterojunctions. Nat. Commun. 2014, 5, 4288. (17) Causa, M.; De Jonghe-Risse, J.; Scarongella, M.; Brauer, J. C.; Buchaca-Domingo, E.; Moser, J.-E.; Stingelin, N.; Banerji, N. The Fate of Electron−Hole Pairs in Polymer:Fullerene Blends for Organic Photovoltaics. Nat. Commun. 2016, 7, 12556. (18) Brabec, C. J.; Zerza, G.; Cerullo, G.; De Silvestri, S.; Luzzati, S.; Hummelen, J. C.; Sariciftci, S. Tracing Photoinduced Electron Transfer Process in Conjugated Polymer/Fullerene Bulk Heterojunctions in Real Time. Chem. Phys. Lett. 2001, 340, 232−236. (19) Hwang, I.-W.; Moses, D.; Heeger, A. J. Photoinduced Carrier Generation in P3ht/Pcbm Bulk Heterojunction Materials. J. Phys. Chem. C 2008, 112, 4350−4354. (20) Miranda, P. B.; Moses, D.; Heeger, A. J. Ultrafast Photogeneration of Charged Polarons in Conjugated Polymers. Phys. Rev. B: Condens. Matter Mater. Phys. 2001, 64, 081201. (21) Grage, M. M. L.; Zaushitsyn, Y.; Yartsev, A.; Chachisvilis, M.; Sundstrom, V.; Pullerits, T. Ultrafast Excitation Transfer and Trapping in a Thin Polymer Film. Phys. Rev. B: Condens. Matter Mater. Phys. 2003, 67, 205207. (22) Banerji, N.; Cowan, S.; Leclerc, M.; Vauthey, E.; Heeger, A. J. Exciton Formation, Relaxation, and Decay in Pcdtbt. J. Am. Chem. Soc. 2010, 132, 17459−17470. (23) Guo, J.; Ohkita, H.; Benten, H.; Ito, S. Near-Ir Femtosecond Transient Absorption Spectroscopy of Ultrafast Polaron and Triplet Exciton Formation in Polythiophene Films with Different Regioregularities. J. Am. Chem. Soc. 2009, 131, 16869−16880. (24) Kaake, L. G.; Zhong, C.; Love, J. A.; Nagao, I.; Bazan, G. C.; Nguyen, T.-Q.; Huang, F.; Cao, Y.; Moses, D.; Heeger, A. J. Ultrafast Charge Generation in an Organic Bilayer Film. J. Phys. Chem. Lett. 2014, 5, 2000−2006. (25) Gelinas, S.; Rao, A.; Kumar, A.; Smith, S. L.; Chin, A. W.; Clark, J.; van der Poll, T. S.; Bazan, G. C.; Friend, R. H. Ultrafast Long-Range Charge Separation in Organic Semiconductor Photovoltaic Diodes. Science 2014, 343, 512−516. (26) Smith, S. L.; Chin, A. W. Ultrafast Charge Separation and Nongeminate Electron-Hole Recombination in Organic Photovoltaics. Phys. Chem. Chem. Phys. 2014, 16, 20305−20309. (27) Savoie, B. M.; Rao, A.; Bakulin, A. A.; Gelinas, S.; Movaghar, B.; Friend, R. H.; Marks, T. J.; Ratner, M. A. Unequal Partnership: Asymmetric Roles of Polymeric Donor and Fullerene Acceptor in Generating Free Charge. J. Am. Chem. Soc. 2014, 136, 2876−2884. (28) Cheung, D. L.; Troisi, A. Theoretical Study of the Organic Photovoltaic Electron Acceptor Pcbm: Morphology, Electronic Structure, and Charge Localization. J. Phys. Chem. C 2010, 114, 20479−20488. (29) Tamura, H.; Tsukada, M. Role of Intermolecular Charge Delocalization on Electron Transport in Fullerene Aggregates. Phys. Rev. B: Condens. Matter Mater. Phys. 2012, 85, 054301. (30) Rebentrost, P.; Mohseni, M.; Kassal, I.; Lloyd, S.; Aspuru-Guzik, A. Environment-Assisted Quantum Transport. New J. Phys. 2009, 11, 033003. (31) Plenio, M. B.; Huelga, S. F. Dephasing-Assisted Transport: Quantum Networks and Biomolecules. New J. Phys. 2008, 10, 113019.

(32) Ishizaki, A.; Fleming, G. R. Theoretical Examination of Quantum Coherence in a Photosynthetic System at Physiological Temperature. Proc. Natl. Acad. Sci. U. S. A. 2009, 106, 17255−17260. (33) Sarovar, M.; Ishizaki, A.; Fleming, G. R.; Whaley, K. B. Quantum Entanglement in Photosynthetic Light-Harvesting Complexes. Nat. Phys. 2010, 6, 462−467. (34) Pelzer, K. M.; Griffin, G. B.; Gray, S. K.; Engel, G. S. Inhomogeneous Dephasing Masks Coherence Lifetimes in Ensemble Measurements. J. Chem. Phys. 2012, 136, 164508. (35) Mukamel, S. Comment on “Coherence and Uncertainty in Nanostructured Organic Photovoltaics. J. Phys. Chem. A 2013, 117, 10563−10564. (36) Grover, M.; Silbey, R. Exciton Migration in Molecular Crystals. J. Chem. Phys. 1971, 54, 4843−4851. (37) Spano, F. C.; Silva, C. H- and J-Aggregate Behavior in Polymeric Semiconductors. Annu. Rev. Phys. Chem. 2014, 65, 477−500. (38) Yamagata, H.; Spano, F. C. Interplay between Intrachain and Interchain Interactions in Semiconducting Polymer Assemblies: The Hj-Aggregate Model. J. Chem. Phys. 2012, 136, 184901. (39) Dudenko, D.; Kiersnowski, A.; Shu, J.; Pisula, W.; Sebastiani, D.; Spiess, H. W.; Hansen, M. R. A Strategy for Revealing the Packing in Semicrystalline Pi-Conjugated Polymers: Crystal Structure of Bulk Poly-3-Hexyl-Thiophene (P3ht). Angew. Chem., Int. Ed. 2012, 51, 11068−11072. (40) Wu, Z.; Petzold, A.; Henze, T.; Thurn-Albrecht, T.; Lohwasser, R. H.; Sommer, M.; Thelakkat, M. Temperature and Molecular Weight Dependent Hierarchical Equilibrium Structures in Semiconducting Poly(3-Hexylthiophene). Macromolecules 2010, 43, 4646−4653. (41) Sirringhaus, H.; et al. Two-Dimensional Charge Transport in Self-Organized, High-Mobility Conjugated Polymers. Nature 1999, 401, 685−688. (42) Padinger, F.; Rittberger, R. S.; Sariciftci, N. S. Effects of Postproduction Treatment on Plastic Solar Cells. Adv. Funct. Mater. 2003, 13, 85−88. (43) Knupfer, M.; Fink, J.; Fichou, D. Strongly Confined Polaron Excitations in Charged Organic Semiconductors. Phys. Rev. B: Condens. Matter Mater. Phys. 2001, 63, 165203. (44) Holdcroft, S. A Photochemical Study of Poly(3-Hexylthiophene). Macromolecules 1991, 24, 4834−4838. (45) Zhokhavets, U.; Erb, T.; Gobsch, G.; Al-Ibrahim, M.; Ambacher, O. Relation between Absorption and Crystallinity of Poly(3Hexylthiophene)/Fullerene Films for Plastic Solar Cells. Chem. Phys. Lett. 2006, 418, 347−350. (46) Brown, P. J.; Thomas, D. S.; Köhler, A.; Wilson, J. S.; Kim, J.-S.; Ramsdale, C. M.; Sirringhaus, H.; Friend, R. H. Effect of Interchain Interactions on the Absorption and Emission of Poly(3-Hexylthiophene). Phys. Rev. B: Condens. Matter Mater. Phys. 2003, 67, 064203. (47) Kouki, F.; Spearman, P.; Valat, P.; Horowitz, G.; Garnier, F. Experimental Determination of Excitonic Levels in A-Oligothiophenes. J. Chem. Phys. 2000, 113, 385−391. (48) Beenken, W. J. D.; Pullerits, T. Excitonic Coupling in Polythiophenes: Comparison of Different Calculation Methods. J. Chem. Phys. 2004, 120, 2490−2495. (49) Kandada, A. R. S.; Grancini, G.; Petrozza, A.; Perissinotto, S.; Fazzi, D.; Raavi, S. S. K.; Lanzani, G. Ultrafast Energy Transfer in Ultrathin Organic Donor/Acceptor Blend. Sci. Rep. 2013, 3, 2073. (50) Albinsson, B.; Eng, M. P.; Pettersson, K.; Winters, M. U. Electron and Energy Transfer in Donor-Acceptor Systems with Conjugated Molecular Bridges. Phys. Chem. Chem. Phys. 2007, 9, 5847−5864. (51) Smyth, C.; Fassioli, F.; Scholes, G. D. Measures and Implications of Electronic Coherence in Photosynthetic LightHarvesting. Philos. Trans. R. Soc., A 2012, 370, 3728−3749. (52) Haken, H.; Strobl, G. An Exactly Solvable Model for Coherent and Incoherent Exciton Motion. Z. Phys. A: Hadrons Nucl. 1973, 262, 135−148. (53) Poelking, C.; Daoulas, K.; Troisi, A.; Andrienko, D. Morphology and Charge Transport in P3ht: A Theorist’s Perspective. In P3ht L

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C Revisited−from Molecular Scale to Solar Cell Devices; Ludwigs, S., Ed. Springer Berlin Heidelberg: Berlin, Heidelberg, 2014; pp 139−180. (54) Alberga, D.; Perrier, A.; Ciofini, I.; Mangiatordi, G. F.; Lattanzi, G.; Adamo, C. Morphological and Charge Transport Properties of Amorphous and Crystalline P3ht and Pbttt: Insights from Theory. Phys. Chem. Chem. Phys. 2015, 17, 18742−18750. (55) Liu, T.; Cheung, D. L.; Troisi, A. Structural Variability and Dynamics of the P3ht/Pcbm Interface and Its Effects on the Electronic Structure and the Charge-Transfer Rates in Solar Cells. Phys. Chem. Chem. Phys. 2011, 13, 21461−21470. (56) Liu, T.; Troisi, A. Absolute Rate of Charge Separation and Recombination in a Molecular Model of the P3ht/Pcbm Interface. J. Phys. Chem. C 2011, 115, 2406−2415. (57) Arkhipov, V. I.; Emelianova, E. V.; Bässler, H. Hot Exciton Dissociation in a Conjugated Polymer. Phys. Rev. Lett. 1999, 82, 1321− 1324. (58) Basko, D. M.; Conwell, E. M. Hot Exciton Dissociation in Conjugated Polymers. Phys. Rev. B: Condens. Matter Mater. Phys. 2002, 66, 155210. (59) Lee, J.; Vandewal, K.; Yost, S. R.; Bahlke, M. E.; Goris, L.; Baldo, M. A.; Manca, J. V.; Voorhis, T. V. Charge Transfer State Versus Hot Exciton Dissociation in Polymer−Fullerene Blended Solar Cells. J. Am. Chem. Soc. 2010, 132, 11878−11880. (60) Chu, C.-W.; Yang, H.; Hou, W.-J.; Huang, J.; Li, G.; Yang, Y. Control of the Nanoscale Crystallinity and Phase Separation in Polymer Solar Cells. Appl. Phys. Lett. 2008, 92, 103306. (61) Dutton, G. J.; Robey, S. W. Distance Dependence of Exciton Dissociation at a Phthalocyanine−C60 Interface. J. Phys. Chem. C 2013, 117, 25414−25423. (62) Caruso, D.; Troisi, A. Long-Range Exciton Dissociation in Organic Solar Cells. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 13498− 13502.

M

DOI: 10.1021/acs.jpcc.6b11332 J. Phys. Chem. C XXXX, XXX, XXX−XXX