Simplified Charge Separation Energetics in a Two-Dimensional Model

Dec 7, 2004 - Kristian O. Sylvester-Hvid*. Department of Chemistry, H. C. Ørsted Institute, UniVersity of Copenhagen,. DK-2100 Copenhagen Ø, Denmark...
0 downloads 0 Views 316KB Size
200

J. Phys. Chem. B 2005, 109, 200-208

Simplified Charge Separation Energetics in a Two-Dimensional Model for Polymer-Based Photovoltaic Cells Kristian O. Sylvester-Hvid* Department of Chemistry, H. C. Ørsted Institute, UniVersity of Copenhagen, DK-2100 Copenhagen Ø, Denmark

Mark A. Ratner Department of Chemistry, Center for Nanofabrication and Molecular Self-Assembly, and Materials Research Center, Northwestern UniVersity, 2145 Sheridan Road, EVanston, Illinois 60208-3113 ReceiVed: August 11, 2004; In Final Form: October 12, 2004

An extension of our two-dimensional working model for photovoltaic behavior in binary polymer and/or molecular photoactive blends is presented. The objective is to provide a more-realistic description of the charge generation and charge separation processes in the blend system. This is achieved by assigning an energy to each of the possible occupation states, describing the system according to a simple energy model for exciton and geminate electron-hole pair configurations. The energy model takes as primary input the ionization potential, electron affinity and optical gap of the components of the blend. The underlying photovoltaic model considers a nanoscopic subvolume of a photoactive blend and represents its p- and n-type domain morphology, in terms of a two-dimensional network of donor and acceptor sites. The nearest-neighbor hopping of charge carriers in the illuminated system is described in terms of transitions between different occupation states. The equations governing the dynamics of these states are cast into a linear master equation, which can be solved for arbitrary two-dimensional donor-acceptor networks, assuming stationary conditions. The implications of incorporating the energy model into the photovoltaic model are illustrated by simulations of the short circuit current versus thickness of the photoactive blend layer for different choices of energy parameters and donor-acceptor topology. The results suggest the existence of an optimal thickness of the photoactive film in bulk heterojunctions, based on kinetic considerations alone, and that this optimal thickness is very sensitive to the choice of energy parameters. The results also indicate space-charge limiting effects for interpenetrating donor-acceptor networks with characteristic domain sizes in the nanometer range and high driving force for the photoinduced electron transfer across the donor-acceptor internal interface.

I. Introduction This contribution is concerned with the theoretical description and modeling of thin-film organic solar cells, currently the subject of numerous experimental investigations (for reviews and overviews, cf. refs 1-3). We report improvements in our working model and demonstrate how these are manifested in simulations of the short-circuit current (Isc) for the ITO|MDMO-PPV:PCBM|Al system.4 Our model for the photovoltaic (PV) effect in disordered polymer (and/or molecular)-based thin-film solar cells has an entirely microscopic origin. Based on molecular parameters, it allows the stationary kinetics of the photogenerated charge carriers to be solved for in arbitrary two-dimensional (2D) donor-acceptor (DA) networks. The model is particularly aimed at analyzing how microscopic polymer morphology influences the PV performance of thin-film based solar cells. In our initial contribution, we were qualitatively able to account for the performance improvement changing from a heterojunction to a bulk heterojunction device structure.5 However, in that study, we arrived at an unphysical limiting behavior for the scaling of Isc with the size of the polymer volume considered. * Author to whom correspondence should be addressed. E-mail: ksh@ theory.ki.ku.dk, www.sylvesterhvid.dk/kristian.

We believe that this is attributable to the improper assumption made there that all states of the system, a priori, were considered equally likely. Obviously, states describing short and long geminate electron-hole pairs (e-/h+), respectively, cannot be equally likely. Failing to include this leads to exaggeration of the bulk charge carrier generation. To improve the model, we must consider more carefully charge separation (CS) and recombination in disordered molecular and polymeric materials. Recent interest in this area stems (in part) from efforts made to understand and optimize the performance of organic lightemitting diodes (OLEDs) and PV devices, based on the experimentally accessible energy parameters of the constituent components. Among tunable material parameters, the optical (∆Eo), single particle (∆Es), and transport gaps (∆Et),6,7 along with the ionization potential (Ip) and the electron affinity (Ea), are of importance. The issue has spawned much controversy as epitomized by the ongoing debate as to the proper magnitude of the exciton binding energy (b) for PPV and its derivatives.7-9 In molecular van der Waals bound systems (typically obtained by vacuum deposition), excitons are considered as strongly bound and localized to within a single molecule or oligomer segment (Frenkel excitons), with b values as large as ∼1.5 eV.10 These conditions favor (radiative) charge recombination in OLEDs. Conversely, in ordered polymer systems, such as

10.1021/jp0463767 CCC: $30.25 © 2005 American Chemical Society Published on Web 12/07/2004

Charge Separation in Polymer-Based Photovoltaics carefully prepared poly(p-phenylene-vinylene) (PPV) samples, excitons are considered as loosely bound delocalized entities that can be described within the one-electron approximation, with b values in the meV range.8,11,12 For these systems, it is argued that the band picture can be applied for transport along the polymer strands, which have the characteristics of a onedimensional crystal lattice. Such conditions are interesting from a PV point of view, favoring photoinduced CS and high mobility in low-dielectric systems. Between these extreme cases lies the reality of most polymer-based PV devices that use, e.g., derivatives of PPV with b values in the range of 0.1-0.9 eV,6,7,9,13,14 where charge transfer (CT) excitons are also argued to be part of the CS process.15,16 Fundamental to our model is the assumption of disorder and localization of charge carriers on sites. Charge carrier transport is consequently treated as activated hopping among sites. We assume localization of excitons and thus their energetic binding to within individual hopping sites. Thus, the transition from an exciton via a geminate e-/h+ pair to an unbound e-/h+ pair is an overall endothermic process within a single material. The nature of the particular e-/h+ interaction along the entire CS coordinate has received much attention, since the original models for geminate e-/h+ recombination as covered in the monograph by Pope and Swenberg.17 In the simplest approximation of a screened Coulomb potential for the e-/h+ pair, even on including the polarization stabilization of the charges (i.e., polarons), complete CS will require energies orders of magnitude larger than kBT. This is at variance with the occurrence of photocurrents seen already at the absorption edge, e.g., in PPV derivatives.18 Hence, mechanisms involving hot excitons,19-22 tandem photogeneration,18 and Gaussian distributions of b 8 have been proposed to account for this and the subsequent increase in photocurrent ∼1 eV above the absorption edge. However, in discussing these mechanisms, one should recall that the systems considered are pristine materials with low extrinsic doping levels and thus few dissociation sites. In the present study, we consider the e-/h+ dynamics in DA blends, i.e., systems with extremely high concentration of dissociation sites (the A sites), which provide the required driving force for CS. Consequently, the relative importance of CS within either D or A domains is reduced compared to CS mediated by the DA interface. In the context of DA composites, where the focus is morphology effects, approximating the e-/h+ pair interaction with a screened Coulomb potential is justified, provided we account for the energy change as carriers cross domain DA boundaries. This approximate approach involves only a few input parameters, a condition that is not met by more elaborate electrostatic/quantum mechanical models for photoinduced charge separation and recombination.13,23-25 By adopting a simple Coulomb potential for CS, our model becomes biased, favoring charge carrier generation at DA interfaces over carrier generation within domains. Combined with the neglect of exciton diffusion (energy migration) by which excitons generated only at the DA interface are likely to dissociate, thus makes the model quite sensitive to variations in the DA morphology. To illustrate the implications of adopting this simple e-/h+ interaction model, we simulate how Isc varies with the thickness (d) of the photoactive film for different choices of Ea and nearest-neighbor interaction energies (vide infra). With respect to device fabrication and optimization of PV efficiency, the ability to predict the optimal thickness of the photoactive film is an important asset. From optical device modeling, it is wellknown that interference effects that are due to the device

J. Phys. Chem. B, Vol. 109, No. 1, 2005 201

Figure 1. Schematic representation of a PV device composed of a ∼100-nm-thick photoactive polymer blend layer sandwiched between front and back electrodes. The nanosized volume is approximated as a distribution of noninteracting 2D patches, each of which subsequently is treated as an individual microscopic PV device.

interfaces lead to (non-Lambert-Beer like) variations in the absorption profile within the photoactive film.26-28 In such studies, the contribution (if any) to the absorption profile due to variations in the microscopic charge carrier dynamics with d are not accounted for. This effect, which may well be morphology-dependent, is investigated here in terms of Isc and its dependence on d. Simulations were performed for two distinctly different topologies that we previously used to model heterojunction and bulk-heterojunction device structures, respectively. Finally, we demonstrate how device performance can be improved by optimization of a given subset of energy parameters. In section II, we briefly outline the essentials of our simulation model and describe the e-/h+ interaction model. In section III, we then present numerical simulations for Isc(d), the results of which are discussed in section IV, with concluding remarks given in section V. II. The Model The geometry of the PV device that we seek to model is shown schematically in Figure 1, where the photoactive film is sandwiched between a transparent indium tin oxide (ITO) electrode and a metal back contact. In the most general case, this film may be a blend of p- (dark) and n-type (white) materials, typically with a thickness in the range of 100 nm. Focus is on some nanoscopic representative volume of this film, which is approximated as a distribution of noninteracting twodimensional PV patches, as also shown in Figure 1. It is assumed that each of these 2D patches are individual, active PV units, contributing to the total PV current from the nanoscopic volume in Figure 1. Thus, the blend morphology of this volume is represented collectively, in terms of the 2D morphology for each of the individual patches. In the implementation presented here, the model is used to perform PV simulations for single patches. This is done by imposing a grid representation of the blend morphology as p- (dark) and n-type (white) sites arranged between virtual electrodes, as illustrated in Figure 2. The sites are enumerated either according to a coordinate representation, with the photoactive layer being A sites long and B sites wide, or by consecutive numbering with k ) 1, 2, ..., AB, where the total number of sites is K ) AB. We further assume a stratified blend morphology by which all a ) 1 sites are p-type and all a ) A are n-type as shown in Figure 2c. The grid in panels b and c in Figure 2 has molecular resolution such that dark sites refer to donor (D) molecules or polymer segments and white sites correspond to acceptor (A) molecules or polymer segments.

202 J. Phys. Chem. B, Vol. 109, No. 1, 2005

Figure 2. Transition from the 2D blend morphology for a PV patch to a virtual microscopic PV device: (a) 2D representation the p- (dark) and n-type (white) blend morphology; (b) discretized representation of this morphology by a grid A sites long and B sites wide. Sites are either p- (dark) or n-type (white) and are numbered consecutively by k ) 1, 2, 3, ..., K with K ) AB or in coordinate representation as (a, b); (c) the virtual 2D PV device made by interfacing the PV patch with structure-less electrodes and external circuitry. Sites in the first (a ) 1) and last (a ) A) columns are entirely p- and n-type making for hole and electron membranes, respectively. The nearest-neighbor hopping sphere around site k is encircled.

Sylvester-Hvid and Ratner two-level systems) on each site are responsible for exciton formation, and are indicated with broken arrows in Figure 3. The kinetics is formulated in terms of electron transfer (ET) processes. Hence, for each type of transition (arrow) in Figure 3, we associate a generic rate, and, similarly, for each type of allowed transition among sites in a general DA network. Determination of the actual rates for these transitions has been described previously,5 and here we only consider the necessary modifications that are due to the implementation of the energy model. B. Occupation State Space. To cast the kinetic equations governing the e-/h+ dynamics of the illuminated PV patch into a linearized master equation, for the particular DA network in question, we need to build the occupation state space, according to the procedure described previously.5 It consists of the reference state |1〉 and all possible single excitations |n〉 ) a+ q ap|1〉 out of this state, where index p refers to HOMO levels and index q to LUMO levels. The total number of occupation states thus generated is 1 + K2. We use a reduced representation, accounting only for the K active electrons distributed among 2K levels. Hence, for a given occupation state

we specify the occupations according to

Onp )

{

for 1 e p e K and

Onq ) Figure 3. Schematic representation of the simple M|DDA|M structure, in terms of HOMO (hA,D) and LUMO (lA,D) energy levels and the Fermi levels of the electrodes. Energies are given relative to the vacuum in electron volts. The different hopping processes involved are shown by arrows. Filled arrows are for all nearest-neighbor (horizontal) electron-hopping processes and dotted arrows are for the injection/ ejection processes. The optically induced (vertical) excitation/deexcitation processes are shown with broken arrows. The numbering of levels is according to the consecutive numbering scheme.

A. Charge Carrier Dynamics. Upon illumination, charge carriers are generated and allowed to perform nearest-neighbor hops among the D and A sites. In Figure 2c, the nearest-neighbor hopping sphere is shown for site k. For sites that are the nearest neighbor to either the anode or cathode, charge carriers may be exchanged with the electrodes under the condition that the PV patch, as a whole, remains uncharged. To illustrate the hopping dynamics in closer detail, consider the energy diagram for the hypothetical one-dimensional (1D) PV system M|DDA|M as depicted in Figure 3. As shown, D and A sites are represented in terms of their highest occupied molecular orbital (HOMO) (hA,D) and lowest unoccupied molecular orbital (LUMO) levels (lA,D) and the electrodes in terms of the Fermi level (0). Hopping of charge carriers between neighboring sites occurs either via HOMO or LUMO levels, as indicated by filled arrows in Figure 3. Exchange of charge carriers with the electrodes only occurs via HOMO levels on a ) 1 donor sites and via LUMO levels on a ) A acceptor sites, as also shown with dotted arrows in Figure 3. Finally, optically induced (vertical) transitions between HOMO and LUMO levels (optical pumping of

}

(2a)

}

(2b)

1 (for doubly occupied HOMO) 0 (for singly occupied HOMO)

{

1 (for singly occupied LOMO) 0 (for vacant LOMO)

for K < q e 2K. Hence, the reference occupation state takes the form

with K doubly occupied HOMO levels and K vacant LUMO levels. In this representation, no level can be more than doubly occupied and each of the K sites (X ) D or A) can be either neutral (X), cationic (X+), anionic (X-), or excited (X*). Different spin states are not distinguished. Limiting the occupation state space to singly excited configurations is a good approximation for solar cell applications and the patch sizes as typically considered here.5 Also, in this approximation, each occupation state represents, at most, one geminate e-/h+ pair configuration, significantly reducing the complexity of the energy model. C. The Energy Model. To each particular occupation state vector |n〉, we assign an electrostatic energy, taken relative to the reference state |1〉, the energy of which we arbitrarily set to zero. Occupation states are categorized as excitonic, chargetransfer, charge-separated states, and the reference state |1〉. Excitonic occupation states (EX states) result from (vertical) excitation of |1〉 and have an energy of ∆Eo[k], where k is the index of the site at which the photoexcitation occurs. Charge-separated occupation states (CS states) involve an e-/h+ pair on sites X and Y separated by distance R, thus

Charge Separation in Polymer-Based Photovoltaics

J. Phys. Chem. B, Vol. 109, No. 1, 2005 203

forming a cation-anion pair that we generally represent as X+-//-Y-. With {X,Y} ) {A or D} we have four different types of such pairs, which, in the simple Coulomb picture, have the electrostatic energy10 g g EXY CS (R) ) Ip(X) - Ea(Y) - (P+ + P_) -

e2 +∆P(R) 4πχ0R (4)

where Igp(X) is the gas-phase ionization potential of X, Ega(Y) the gas-phase electron affinity of Y, and (P+ + P-) the polarization energy of a fully relaxed h+ and e-, at infinite separation in the bulk materials of X and Y, respectively. The fourth term in eq 4 is the Coulomb attraction of the unscreened e-/h+ pair and ∆P(R) is the subsequent reduction in polarization stabilization as the h+ and e- are brought together at distance R in the bulk. We make the approximation that the charge carriers experience only the bulk of either the D or A material, and, thus, the solid-state ionization potential and electron affinity are Ip(X) ) Igp(X) - P+ and Ea(Y) ) Ega(Y) + P-.17 If the blend has a dielectric constant χ, the two last terms in eq 4 can be approximated as -e2/(4πχR),29 and

EXY CS (R) ) Ip(X) - Ea(Y) -

e2 4πχR

(5)

At this level of approximation, the exciton binding energy in either material (X ) D or A) becomes Xb ) Ip(X) - Ea(X) ∆EXo , i.e., the difference between the single particle gap and the optical gap.6 By adopting eq 5 in the PV context, we implicitly assume that charge carriers are injected directly into polaron levels from the electrodes. The solid-state values Ip(X) and Ea(Y) in eq 5, in principle, include all electronic, vibrational, and lattice polarization contributions. However, using experimental values for Ip(X) and Ea(Y) requires awareness, because not all polarization contributions are necessarily captured by a given technique,10,30 making the simplicity of eq 5 somewhat deceptive. For weakly interacting materials, the approximation leading to eq 4 holds almost down to nearest-neighbor distances,29 where partial or full CT then may occur.10,15 Hence, for occupation states describing nearest-neighbor e-/h+ configurations (CT states), we generally add a correction term, ∆EXY CT , in eq 5 that accounts for interaction contributions beyond the electrostatic picture. Therefore, for CT states, the analogue of eq 5 becomes

EXY CT (R)) Ip(X) - Ea(Y) -

e2 + ∆EXY CT 4πχR

Figure 4. Energy ordering of the 10 possible occupation states {|1〉,|2〉,..., |10〉} resulting from single excitation out of |1〉 for the DDA structure. Also shown are the reduced occupations of the HOMO and LUMO levels on D and A sites for each occupation state. The energy parameters are as follows: Ea[D] ) 2.6 eV, Ea[A] ) 2.8 eV, Ip[D] ) 5.0 eV, Ea[A] ) 6.1 eV, ∆EXY CT ) 0.0 eV, χ ) 3.0, and a site separation of 1 nm.

D. The Master Equation. The kinetics of the illuminated PV system, under stationary conditions, may be described by the master equation

d P(t) ) MP(t) ) 0 dt

(7)

where Pn(t) is the average probability for the system to be in occupation state |n〉 at time t. M is the transition matrix, and element Mmn then is the average rate for the system to make the transition |m〉 r |n〉. For a closed system where the number of electrons, 2K, is conserved, the diagonal elements of Mmn take the form

Mmm ) -

∑ Mmn

(8)

n*m

To solve the set of homogeneous linear equations in eq 7 for Pn, we consider the supplementary condition of probability conservation N

(6)

+ where ∆EXY CT is specific to the cation-anion pair (X -//-Y ) in question and must be inferred from either quantum chemical methods25 or spectroscopic data.31 As a simple example, Figure 4 shows the energy ranking of the possible occupation states {|1〉, |2〉, |3〉, ..., |10〉} for the M|DDA|M structure derived with ∆EXY CT ) 0.0 eV and Ip and Ea as specified in section III. In Figure 4, we see that |4〉, |6〉, and |8〉 are EX type, |2〉 and |10〉 are CS type, and |3〉, |5〉, |7〉, and |9〉 are CT type. This ranking immediately illustrates the importance of ∆EXY CT : i.e., in the purely electrostatic picture, ) 0.0 eV, the D*DA configuration is higher in where ∆EDD CT energy than the D+D-A and D-D+A configurations, at variance with the general case for π-conjugated polymers.25

∑ Pn ) 1

(for 0 e Pn e 1)

(9)

n)1

and then replace one of the N equations in eq 7 with eq 9, where N is the total number of occupation state vectors. Solving for Pn then becomes a matter of solving a set of inhomogeneous equations, which is done by standard procedures. To obtain the transition rates that constitute M, we need identify the various types of site couplings, as, e.g., indicated by arrows in Figure 3. Generally, we identify three types of couplings. Horizontally coupled occupation states include pairs of states (|n〉, |m〉) that can be connected by a nearest-neighbor move of an electron. Rates for these transitions are specified by the initial and final location of the electron transferred. Hence, if the nearestneighbor coupling moves an electron from a donor LUMO level

204 J. Phys. Chem. B, Vol. 109, No. 1, 2005

Sylvester-Hvid and Ratner

in state |n〉 to an acceptor LUMO level in state |m〉, the generic ET rate lWAD is assigned to it. As |n〉 and |m〉 generally differ in energy by an amount ∆E[n,m], the endothermic transition occurs at a rate of lWAD exp{-∆E[n,m]/(kBT)}, where T is the device temperature. Vertically coupled occupation states designate pairs of states connected by on-site HOMO-LUMO transitions and, here, always occur between |1〉 and some singly excited occupation state |m〉. The rates for these optically induced excitation (absorption) and de-excitation (stimulated and spontaneous emission) processes at site k are5 ex

W(k) )

4π(∆Eo[k])2

dex

c2hme

[ ( ) ] ( )

Rfk exp

W(k)) exW(k) exp

∆Eo[k] -1 kBTrad

∆Eo[k] kBTrad

-1

(10)

(11)

where it is assumed that the intensity of the incident solar radiation can be approximated as Planckian blackbody radiation at temperature Trad. In eqs 10 and 11, R is the fine structure constant, me the electron mass, fk the oscillator strength for the transition at site k, and ∆Eo[k] the corresponding HOMOLUMO gap. It is further assumed that the intensity remains constant throughout the photoactive film and that light absorption occurs only at the monochromatic frequency corresponding to ∆Eo[k] at site k. Hence, hot excitons are not taken into consideration here. Current generating transitions are due to pairs of occupation states coupled through simultaneous injection or ejection of an e-/h+ pair. Injection refers to the process where an electron is transferred from a LUMO level on an a ) A acceptor site to the anode, simultaneously with an electron is transferred from the cathode to a HOMO level on an a ) 1 donor site. The net process, i.e., injection of an e-/h+ pair into the respective electrodes, is assigned the single generic rate inW. The reverse process is termed ejection and occurs with the generic rate ejW. The present study is limited to simulation of Isc, and, under these conditions, we assume that there is no built-in field, i.e., flat-band conditions. Considered together with the energy levels as specified in Figure 3, this implies that injection of an e-/h+ pair generally will be energetically favorable. Conversely, the ejection process is endothermic with a barrier of (0 - lA) and occurs with the rate ejW exp[-(0 - lA/(kBT)]. III. Numerical Simulations Simulations were performed for the double-cable and laminar structures depicted in Figure 5 for different thicknesses (d) of the photoactive layer, using the current software implementation of the model (PhotoSim 1.12). With regard to the grid length in Figure 2, the film thickness is given as d ) A × 1 nm. The use of structures only two sites wide does not notably effect our conclusions: as confirmed elsewhere, Isc scales linearly with the width.32 Also, the validity of eq 5 is still justified as we always imagine the PV patch as just a representative subunit of the bulk. External and internal model parameters were chosen as done previously,5 specifically to mimic the ITO|MDMO-PPV: PCBM|Al4 polymer blend system. The parameters required by the energy model are the solid-state ionization potentials Ip[MDMO-PPV] and Ip[PCBM] and electron affinities Ea[MDMO-PPV] and Ea[PCBM] for the pristine components along with ∆EXY CT for the nearest-neighbor non-Coulombic energy correction terms. Ionization potentials were kept constant

Figure 5. Laminar and double-cable DA structures investigated in this study and their orientation relative to the fictitious electrodes; the PV patches generally are two sites wide and A sites long, the latter of which relates directly to the thickness d of the photoactive layer.

Figure 6. Electrostatic energy of CS and CT occupation states versus cation-anion separation as obtained in simulations for 2 × 30 doublecable structures. Energy parameters were chosen as Ea[MDMO-PPV] ) 2.60 eV, Ea[PCBM] ) 2.80 eV, and ∆EXY CT ) 0.00 eV. Also shown are the optical band gaps ∆Eo for MDMO-PPV and PCBM.

at Ip[MDMO-PPV] ) 5.0 eV and Ip[PCBM] ) 6.1 eV, and as a starting point ∆EXY CT was chosen to 0.0 eV for all cationanion combinations. For MDMO-PPV, we used Ea ) 2.60 eV as estimated from measurements on MEH-PPV.33 Ea for PCBM was estimated by adding the difference of the half-wave potentials for C60 (E1 ) -1.11 V) and the methyl ester of PCBM (E1 ) -1.19 V) measured in an identical solvent system34 to the value of Ea for solid C60.35 Hence, Ea[PCBM] ) 2.8 eV was used. In Figure 6, we show the electrostatic energy of the possible CS and CT occupation states, relative to the reference state |1〉, plotted versus the cation-anion separation. As this separation increases, the individual potential curves converge to their limiting value Ip(X) - Ea(Y). For the largest separations shown in Figure 6, the energies are less than a few percent off this limit. Note that for X ) Y, Ip - Ea corresponds to the single particle gap (∆Es) for component X. Also shown in Figure 6 with broken and dotted lines are the ∆Eo values for MDMOPPV and PCBM, respectively. To assess the importance of ∆EXY CT , simulations of Isc for double-cable and laminar structures were performed for ∆EXY CT ) 0.00, (0.05, (0.10, and (0.20 eV, for d up to 30 nm. Results obtained on changing ∆EDA CT are shown in Figure 7, whereas AA , ∆E , and ∆EAD changing ∆EDD CT CT CT showed little or no effect

Charge Separation in Polymer-Based Photovoltaics

Figure 7. Isc versus thickness of the photoactive layer, as simulated for (a) double-cable and (b) laminar structures for different choices of ∆EDA CT and with Ea[MDMO-PPV] ) 2.60 eV and Ea[PCBM] ) 2.80 eV.

Figure 8. Isc versus thickness of the photoactive layer, as simulated for (a) double-cable and (b) laminar structures for different choices of Ea for PCBM and with ∆EXY CT ) 0.00 eV and Ea[MDMO-PPV] ) 2.60 eV.

and, thus, those results are omitted here. Comparing panels a and b in Figure 7, it is obvious that the dependence of Isc on the thickness of the photoactive layershenceforth referenced as Isc(d)sdiffers distinctly for the two topologies investigated. For double-cable structures, there are indications of an optimal thickness, i.e., the thickness dopt, where Isc(d) is at maximum, which, in turn, is dependent on the particular choice of ∆EDA CT . For laminar structures, Isc simply decreases with d at a rate also dependent on ∆EDA CT . For both topologies in Figure 7, larger Isc(d) result for positive ∆EDA CT , and it seems that double-cable structures are more sensitive to variations in ∆EDA CT than laminar structures. Electron affinities are frequently argued to be determining for the CS dynamics, and, therefore, simulations of Isc for different choices of Ea[MDMO-PPV] and Ea[PCBM] were performed. The simulations were extended to include structures with d up to 100 nm, to demonstrate further that the model can address layer thicknesses of experimental relevance. In Figure 8, we show results for simulations of Isc(d) performed for Ea[MDMO-PPV] ) 2.60 eV, ∆EXY CT ) 0.00 eV, and different choices of Ea for PCBM. Similarly, simulations for Ea[PCBM] ) 2.80 eV, ∆EXY CT ) 0.00 eV, and different choices of Ea[MDMO-PPV] are shown in Figure 9. From both Figures 8

J. Phys. Chem. B, Vol. 109, No. 1, 2005 205

Figure 9. Isc versus thickness of the photoactive layer as simulated for (a) double-cable and (b) laminar structures for different choices of Ea for MDMO-PPV and with ∆EXY CT ) 0.00 eV and Ea[PCBM] ) 2.80 eV.

and 9, we again see markedly different behavior of Isc for the two topologies: Isc generally decreases with d starting from ∼2 × 10-19 A for laminar structures. For double-cable structures, Isc generally increases at first, also starting from ∼2 × 10-19 A, to reach a maximum after which its value decreases with further increases in d. Increasing Ea[PCBM] has different effects on Isc(d) for double-cable and laminar structures, as observed when comparing panels a and b in Figure 8. In the laminar case, small increases in Isc(d) result, saturating for Ea > 2.85 eV, whereas for double-cable structures, increasing Ea maximizes Isc(d) for Ea ) 2.75 eV, after which Isc(d) decreases with further increases in Ea. Increasing Ea[MDMO-PPV] has almost no effect on Isc(d) for laminar structures, as seen from Figure 9a, whereas for the double-cable results in Figure 9b, Isc(d) increases up to Ea ) 2.75 eV, after which Isc(d) also decreases with further increases in Ea. The results in Figures 7-9 indicate that for double-cable structures, one can optimize Isc(d), with respect to the energy parameters. Therefore, we performed simulations within the parameter space spanned by Ea[PCBM] ) {2.60, 2.65, 2.70, 2.75, 2.80, 2.85, 2.90} eV, Ea[MDMO-PPV] ) {2.50, 2.55, 2.60, 2.65, 2.70, 2.75, 2.80} eV, and ∆EDA CT ) {0.00, 0.05, 0.10, 0.15, 0.20} eV, and found that values of Ea[MDMO-PPV] ) 2.75 eV, Ea[PCBM] ) 2.85 eV, and EDA CT ) 0.15 eV optimize Isc(d). In Figure 10, we show Isc(d) simulated using the optimized parameters (denoted by 0) along with Isc(d) simulated using the initial parameters (denoted by ]) for double-cable structures: Optimization leads to a 3-fold increase in the maximum for Isc and a shift of dopt from ∼10 nm to ∼40 nm. To illustrate further the consequences of the energy model, in Figure 10 we also show Isc(d) as simulated using the previous model5 for laminar structures (denoted by f) and double-cable structures (denoted by O), along with Isc(d) simulated for laminar structures in the present model (denoted by +). Apart from energy parameters investigated here, numerous simulations were performed where each of the other internal parameters were varied separately. These simulations showed that dopt is influenced additionally by variations in the electron and hole mobilities for the pristine components, which is addressed elsewhere.32

206 J. Phys. Chem. B, Vol. 109, No. 1, 2005

Figure 10. Isc versus thickness of the photoactive layer, as simulated for double-cable and laminar device structures with and without the energy model. For simulations where the energy model was used, energy parameters were chosen as Ea[PCBM] ) 2.80 eV, Ea[MDMOPPV] ) 2.60 eV, and ∆EXY CT ) 0.00 eV, except for the optimized case (denoted by 0) where these parameters were Ea[PCBM] ) 2.85 eV, Ea[MDMO-PPV] ) 2.75 eV, and ∆EDA CT ) 0.15 eV.

IV. Discussion Limiting the separation of e-/h+ pairs and, in particular, bulk charge carrier generation by imposing an energy ordering of occupation states clearly alters the behavior of Isc(d) for both the laminar and double-cable structures, as evident from Figure 10 (denoted by + and ], respectively). For double-cable structures, at small d, Isc almost follows the length dependence of the old model (denoted by O) but becomes drastically attenuated at larger d. In the laminar case, the almost constant value of Isc obtained with the old model (denoted by /) is replaced by a rapid decline with increasing d (denoted by +), which eventually approaches the double-cable result (denoted by ]) at ∼100 nm. We generally observe that double-cable structures always yield larger Isc values than laminar structures of similar size and composition, with the difference in Isc, however, decreasing with increasing d. The screened Coulomb potential of the energy model clearly removes the unphysical limiting behavior of Isc with d, as observed in our initial model, where all states were equally likely. The attenuation of Isc with increasing d can be viewed as a gradual increase in series resistance of the device.36 The distinguishing feature for all double-cable simulations is the maximum in Isc(d) which, qualitatiVely, is in agreement with observations made for polyfluorene/PCBM-based bulkheterojunction devices, where, however, a maximum was observed at d ≈ 230 nm.37 It must be stressed that dopt, as determined in this study, is based solely on kinetic considerations, and we thus determine the thickness at which we obtain the best compromise between charge generation capacity and recombination losses. Our approach does not account for an absorption profile through the blend reflecting variable absorptivity and interference effects due to interfaces, both of which may influence Isc(d). The influence of the interfaces on the intensity profile have been addressed using macroscopic optical modeling for heterojunction26 and bulk-heterojunction27,28 devices, and eventually should be superposed onto our microscopic approach. In passing, we note that the characteristic initial increase in Isc(d) for double-cable structures is not due to an initial lack of occupation states connected by current yielding transitions for the very thin structures. This is clear from simulations on wider (B > 2), but otherwise similar, double-

Sylvester-Hvid and Ratner cable structures (not shown here), which showed similar features.32 The dependence of Isc on the thickness of the photoactive bilayer for laminar devices has not been addressed experimentally for the specific MDMO-PPV|PCBM system. Findings for the related PPV|BBL heterojunction device indicates decreasing efficiency with the thickness of n-type material.38 Generally, one would expect laminar device structures to be most affected by limiting bulk charge carrier generation, because they have a smaller “DA interface area” to volume ratio than double-cable structures. This is in accord with our simulations, as observed in Figure 10. With respect to device fabrication, the ability to optimize d is desirable; there are both mechanical and processing reasons for this inclination. In particular, films should not be too thin, simply to reduce the risk of short circuits due to rough substrate surfaces. Note that, in discussing dopt, our concern should be the thickness at which the most power is produced. Here, we refer to Isc instead, because this is more straightforwardly defined in our model. In Figure 10, the results of double-cable simulations, using the generic energy parameters (denoted by ]) indicate dopt to be in the range of 10 nm, which clearly is unrealistic. Subsequent optimization of the energy parameters shows that we may improve Isc(dopt) substantially. Hence, optimization such that Ea[PCBM] increases from 2.80 eV to 2.85 eV, Ea[MDMOPPV] increases from 2.60 eV to 2.75 eV, and ∆EDA CT increases from 0.00 eV to 0.15 eV increases the value of dopt to ∼40 nm and triples the value of Isc, as shown in Figure 10 (denoted by 0). Note that only the general trends are useful here, because characteristic lengths determined in 2D models may not necessarily be carried into three dimensions. As shown by the optimization of the energy parameters, the CS dynamics is critically dependent on the potential curves in Figure 6. Because of the simplicity of the energy model, these individual potential curves are modulated only by Ea and + ∆EXY CT . Hence, increasing Ea[PCBM] makes the A -//-A + and D -//-A curves shift (globally) downward in energy. Similarly, increasing Ea[MDMO-PPV] leads to a downward shift for the D+-//-D- and A+-//-D- curves in Figure 6. The shallowness of the X+-//-Y- potential curve at the nearestneighbor distance may be adjusted by ∆EXY CT , and, as noted, only variations in ∆EDA CT had a notable effect on Isc. This is in accord with the fact that the driving force for CS arises primarily at the DA interface. On analyzing the density of states (not shown here) for laminar and double-cable structures of similar size and composition, it appears that double-cable structures generally have more CT states that involve the -D+A-sequence. Consequently, double-cable structures are expected to be more susceptible to variations in ∆EDA CT than laminar structures, as understood from Figure 5 and confirmed by Figure 7. For both topologies, better PV performance is achieved for + positive ∆EDA CT , because, in this case, the D -//-A potential is more shallow and, thus, counteracts trapping of CT excitons. On increasing ∆EDA CT beyond 0.2 eV (not shown here), Isc(d) is not increased further, because the A+-//-D- and D+-//-Dpotentials then intercept. This implies an increase in the population of the latter states, which do not serve as a direct pathway to CS. In more elaborate investigations, one would therefore need to consider the simultaneous variation of all the ∆EXY CT parameters, to avoid such artifacts. As Ea[PCBM] is increased for laminar structures, the simulations generally lead to larger Isc(d) values, as observed from Figure 8b, although this effect ceases for Ea[PCBM] > 2.85

Charge Separation in Polymer-Based Photovoltaics eV. Qualitatively, this is in accord with expectations on increasing Ea for the electron-accepting (n-type) material. For double-cable structures, however, the value of Isc(d) increases only with Ea[PCBM], up to ∼2.75 eV, after which it decreases. At first, this seems contra-intuitive, because, for Ea[PCBM] > 2.8 eV, the process D*-/∞/-A f D+-/∞/-A- (where ∞ refers infinite separation of D and A), becomes overall energetically favorable. However, as the driving force for CS is increased, so is the likelihood of finding only partially separated D+-//-A- configurations (i.e., short geminate e-/h+ pairs). In the double-cable structures with domains only 1 nm wide as considered here, the subsequent polarization of the DA interface due to short e-/h+ pairs, for large Ea[PCBM], begins to block the transport of charge carriers to the electrode regions. This limitation on the efficiency, which is due to a buildup of space charge, is not observed to the same extent in laminar structures, since there the charge transporting regions are much larger compared to the charge generating DA interface. The problem of charge congestion, which is observed for the double-cable structures, illustrates a general problem of interpenetrating DA networks with nanoscopic domain sizes, in particular when combined with a DA couple with strong driving force for photoinduced CS.39 To improve the situation, the hole and electron mobilities of the p- and n-type materials must be increased and/or domains should become increasingly larger when approaching the electrodes. In the latter case, however, we compromise, by decreasing the internal DA area of the device. Increasing Ea for the p-type (electron-donating) material is expected to facilitate exciton dissociation in the D material and, thus, bulk charge carrier generation. However, in binary DA systems, this process is not primarily responsible for e-/h+ generation, and furthermore, because we neglect the migration of excitons, this mechanism can be considered less important here. As suggested by the potential curves in Figure 6, increasing Ea[MDMO-PPV] will shift the D+-//-D- curve closer to the D+-//-A- potential curve and consequently increase the likelihood of (endothermic) transitions from D+-//-A- to noncharge-generating D+-//-D- occupation states. Upon increasing Ea[MDMO-PPV], we therefore expect a decrease in PV efficiency; however, from Figure 9b, we see that, for laminar structures, this has no effect al all. For the double-cable results in Figure 9a, Isc(d) increases with Ea[MDMO-PPV], up to 2.75 eV, after which point it falls off again. We can understand this only in terms of an indirect increase of the electron mobility in the DA interface region: i.e., as transitions from D+-//-A- to D+-//-D--type configurations become more likely with increasing Ea[MDMO-PPV], the indirect coupling of different, initially uncoupled D+-//-A- configurations is enhanced via connecting D+-//-D--type configurations. Effectively, this facilitates electron transport along the DA interface and, for double-cable structures, this therefore improves the PV performance. However, as Ea[MDMO-PPV] > 2.75 eV, for doublecable devices, we again encounter the situation that has also been observed for large values of Ea[PCBM] in Figure 8b, namely, that the p- and n-type domains become hole- and electron-congested, inhibiting further increases in Isc(d). It is noteworthy that changing Ea by only 50 meV, in concert with adjusting ∆EDA CT by 150 meV, results in significant changes in the PV performance. This clearly underscores the importance of a proper (and possibly more-detailed) description of the energetics of the CS process. The optimized value for Ea[MDMO-PPV] implies b ≈ 50 meV for MDMO-PPV, which is qualitatively in agreement with recent assessments

J. Phys. Chem. B, Vol. 109, No. 1, 2005 207 made for PPV.8 However, the implications of b values this small are in direct contradiction with our initial assumption of binding of excitons in a disordered system. Thus, within our simple energy model, we are not able to maintain a fully consistent description of the energetics. Our initial assumption of a disordered molecular film rigorously is not an argument to exclude energy transfer by exciton diffusion. However, we choose to neglect exciton diffusion, to emphasize the importance of the internal DA interface: this then becomes the primary locus for charge carrier generation. In MDMO-PPV:PCBM systems, exciton diffusion lengths may vary, depending on sample preparation; however, often, ∼10 nm is mentioned.40 Such a length implies that neglecting exciton diffusion will affect simulations of doublecable structures for each additional layer, as the thickness of the film is increased. Hence, this (lacking) contribution to Isc increases as the size of the device. In contrast, in the laminar case, this (lacking) contribution remains constant as the MDMO-PPV layer exceeds the exciton diffusion length. In any case, we conjecture that including exciton diffusion into our model, as well as a possible built-in field, will increase the PV efficiency. Isc, as simulated with the present model, should thus be taken as a lower limit. V. Conclusion By introducing an ordering of the occupation states, according to a simple energy model, the excessive bulk-charge-carrier generation of our original model5 has been corrected and photoinduced e-/h+ pair generation now occurs primarily at the internal DA interface. These modifications cause a general attenuation of Isc, laminar device structures being most affected, and more importantly, give a physically correct limiting behavior of Isc with the thickness of the photoactive film. It was found that increasing Ea[PCBM] for laminar structures increased Isc, whereas increasing Ea[MDMO-PPV] had little effect. For both topologies, the variation of ∆EDA CT had significant impact on Isc and generally better PV performance resulted from positive corrections. When the electron affinities for double-cable structures were increased, the results were more complex. However, it was generally observed that Isc(d) displayed an optimum at a film thickness of dopt, and that this optimum is strongly dependent on the particular choice of Ea[PCBM], Ea[MDMO-PPV], and ∆EDA CT . Varying these parameters has shown that Isc(dopt) could be maximized, and furthermore, that this optimization was accompanied by an increase in dopt. The double-cable results indicate that the common perception that increasing Ea for the n-type (electron accepting) medium leads to better PV performance is an oversimplification. In fact, for the MDMO-PPV:PCBM system, with its extreme driving force for the initial charge carrier generation, effects such as charge carrier congestion and excessive interface polarization become limiting factors when considering DA networks with nanoscopic domain sizes. Also, the double-cable simulations reveal the inevitable tradeoff between increasing charge carrier generation capacity and the series resistance, as the film thickness is increased. More generally, this underscores the necessity of balancing the internal donor-acceptor (DA) interface area (where e-/h+ pair generation occurs) and the volume of p- and n-type domains required for efficient charge carrier extraction. Therefore, our results suggest that increasing Ea for the n-type material in bulk-heterojunction devices must be done in a balanced way, taking mobilities and/or characteristic domain sizes into account. A quantitative understanding

208 J. Phys. Chem. B, Vol. 109, No. 1, 2005 of this interplay has wide implications for the efficiency optimization in these devices, in particular, for bulk-heterojunction devices that are based on dyadic double-cable polymers.41-45 Clearly, there is need to perform full 3D simulations of morphology effects in photoactive composites. Acknowledgment. K.O.S. is indebted to Carlsbergfondet and the Danish Technical Research Council for financial support and the Danish Center for Scientific Computing for providing computational resources. K.O.S. also thanks S. Rettrup and J. C. Hummelen for helpful discussions during the preparation of this work. M.R. thanks the Chemistry Divisions of the NSF, and the DoE. References and Notes (1) Nelson, J. Curr. Opin. Solid State Mater. Sci. 2002, 6, 87. (2) Spanggaard, H.; Krebs, F. C. Sol. Energy Mater. Sol. Cells 2004, 83, 125. (3) Brabec, C. J. Sol. Energy Mater. Sol. Cells 2004, 83, 273. (4) Brabec, C. J.; Sariciftci, N. S.; Hummelen, J. C. AdV. Funct. Mater. 2001, 11, 15. (5) Sylvester-Hvid, K. O.; Rettrup, S.; Ratner, M. A. J. Phys. Chem. B 2004, 108, 4296. (6) Campbell, I. H.; Hagler, T. W.; Smith, D. L.; Ferraris, J. P. Phys. ReV. Lett. 1996, 76, 1900. (7) Conwell, E. M. Synth. Met. 1996, 83, 101. (8) Moses, D.; Schmechel, R.; Heeger, A. J. Synth. Met. 2003, 139, 807. (9) Ba¨ssler, H.; Arkhipov, V. I.; Emelianova, E. V.; Gerhard, A.; Hayer, A.; Im, C.; Rissler, J. Synth. Met. 2003, 135-136, 377. (10) Hill, I. G.; Kahn, A.; Soos, Z. G.; Pascal, R. A., Jr. Chem. Phys. 2000, 327, 181. (11) Moses, D.; Wang, J.; Heeger, A. J.; Kirova, N.; Brazovski, S. Synth. Met. 2002, 125, 93. (12) Moses, D.; Wang, J.; Heeger, A. J.; Kirova, N.; Brazovski, S. Proc. Natl. Acad. Sci. U.S.A. 2001, 98 (24), 13496. (13) Chandross, M.; Mazumdar, S.; Jeglinski, S.; Wei, X.; Vardeny, Z. V.; Kwock, E. W.; Miller, T. M. Phys. ReV. B 1994, 50, 50. (14) Rossi, L.; Alvarado, S. F.; Riess, W.; Schrader, S.; Lidzey, D. G.; Bradley, D. D. C. Synth. Met. 2000, 111-112, 527. (15) Schweitzer, B.; Ba¨ssler, H. Synth. Met. 2000, 109, 1. (16) Hoffmann, M.; Schmidt, K.; Fritz, T.; Hasche, T.; Agranovich, V. M.; Leo, K. Chem. Phys. 2000, 258, 73. (17) Pope, M.; Swenberg, C. E. Electronic Processes in Organic Crystals and Polymers; Oxford University Press: New York, 1999.

Sylvester-Hvid and Ratner (18) Arkhipov, V. I.; Emelianova, E. V.; Ba¨ssler, H. Chem. Phys. Lett. 2001, 340, 517. (19) Basko, D. M.; Conwell, E. M. Synth. Met. 2003, 139, 819. (20) Basko, D. M.; Conwell, E. M. Phys. ReV. B 2002, 66, 155210. (21) Arkhipov, V. I.; Emelianova, E. V.; Barth, S.; Ba¨ssler, H. Phys. ReV. B 2000, 61, 8207. (22) Arkhipov, V. I.; Emelianova, E. V.; Ba¨ssler, H. Phys. ReV. Lett. 1999, 82, 1321. (23) Moore, E.; Yaron, D. Synth. Met. 1997, 85, 1023. (24) Wu, M. W.; Conwell, E. M. Chem. Phys. 1998, 227, 11. (25) Soos, Z. G.; Hennessy, M. H.; Wen, G. Chem. Phys. 1998, 227, 19. (26) Pettersson, L. A.; Roman, L. S.; Ingana¨s, O. J. Appl. Phys. 1999, 86, 487. (27) Hoppe, H.; Arnold, N.; Sariciftci, N. S.; Meissner, D. Sol. Energy Mater. Sol. Cells 2003, 80, 105. (28) Persson, N.-K.; Schubert, M.; Ingana¨s, O. Sol. Energy Mater. Sol. Cells 2004, 83, 169. (29) Bounds, P. J.; Siebrand, W. Chem. Phys. Lett. 1980, 75, 414. (30) Cahen, D.; Kahn, A. AdV. Mater. 2003, 15, 271. (31) Bulovic´, V.; Burrows, P. E.; Forrest, S. R.; Cronin, J. A.; Thompson, M. E. Chem. Phys. 1996, 210, 1. (32) Sylvester-Hvid, K. O.; Rettrup, S. 2004, manuscript to be published. (33) Cervini, R.; Li, X.-C.; Spencer, G. W. C.; Holmes, A. B.; Moratti, S. C.; Friend, R. H. Synth. Met. 1997, 84, 359. (34) Ravaine, S.; Pecq, F. L.; Mingotaud, C.; Delhaes, P.; Hummelen, J. C.; Wudl, F.; Patterson, L. K. J. Phys. Chem. 1995, 99, 9551. (35) Lezius, M. Int. J. Mass Spectrosc. 2003, 223-224, 447. (36) Aernouts, T.; Geens, W.; Poortmans, J.; Heremans, P.; Borghs, S.; Mertens, R. Thin Solid Films 2002, 403-404, 297. (37) Yohannes, T.; Zhang, F.; Svensson, M.; Hummelen, J. C.; Andersson, M. R.; Ingana¨s, O. Thin Solid Films 2004, 449, 152. (38) Jenekhe, S. A.; Yi, S. Appl. Phys. Lett. 2000, 77, 2635. (39) Brabec, C. J.; Zerza, G.; Cerullo, G.; De Silvestri, S.; Luzzati, S.; Hummelen, J. C.; Sariciftci, S. Chem. Phys. Lett. 2001, 340, 232. (40) Haugeneder, A.; Neges, M.; Kallinger, C.; Spirkl, W.; Lemmer, U.; Feldmann, J. Phys. ReV. B 1999, 59, 15346. (41) Nierengarten, J.-F.; Eckert, J.-F.; Nicoud, J.-F.; Ouali, L.; Krasnikov, V.; Hadziioannou, G. Chem. Commun. 1999, 617. (42) Peeters, E.; van Hal, P. A.; Knol, J.; Brabec, C. J.; Sariciftci, N. S.; Hummelen, J. C.; Janssen, R. A. J. J. Phys. Chem. B 2000, 104, 10174. (43) Eckert, J.-F.; Nicoud, J.-F.; Nierengarten, J.-F.; Liu, S.-G.; Echegoyen, L.; Barigelletti, F.; Armaroli, N.; Ouali, L.; Krasnikov, V.; Hadziioannou, G. J. Am. Chem. Soc. 2000, 122, 7467. (44) Ramos, A. M.; Rispens, M. T.; van Duren, J. K. J.; Hummelen, J. C.; Janssen, R. A. J. J. Am. Chem. Soc. 2001, 123, 6714. (45) Cravino, A.; Zerza, G.; Neugebauer, H.; Maggini, M.; Bucella, S.; Menna, E.; Svensson, M.; Andersson, M. R.; Brabec, C. J.; Sariciftci, N. S. J. Phys. Chem. B 2002, 106, 70.