Organic photovoltaics (OPVs) can potentially provide a cost-efficient means of harnessing solar energy. However, optimum OPV performance depends on understanding the process–structure–property (PSP) correlation in organic semiconductors. In the working of bulk-heterojunction OPVs, the morphology plays a crucial role in device performance. In order to understand PSP linkage, a theoretical framework has been developed. We first established process–structure correlations by generating a range of morphologies with various blend ratios of donor and acceptor organic semiconductors for various annealing periods. Second, we calculated the effective electronic properties corresponding to the simulated structures using a diffuse interface approach that is numerically more robust and straightforward than the classical sharp interface method. This novel framework, wherein both the process–structure and the structure–property relationship have been established using the diffuse interface approach, completes the theoretical PSP linkage, allowing the optimization of process parameters for device applications. The theoretical PSP linkage is then benchmarked qualitatively with experimental results on a model P3HT:PCBM system. We have been able to identify the morphological characteristics that maximize device performance. This work is carried out in the broad overview of the integrated computational materials engineering framework wherein the processing parameters are optimized by determining the process–structure–property relationships.

The fabrication of organic solar cells is mainly carried out at low temperature using cost-efficient solvent-based techniques.1 Solution-processed bulk-heterojunction (BHJ) organic photovoltaics (OPVs) came into being in 1995, offering the advantages of low-cost processing, a higher interfacial area that aids in exciton dissociation, and a thicker photoactive layer leading to greater absorption of the incident light.2,3 Fabrication of OPVs is in contrast to the energy-intensive production of their inorganic counterparts. Therefore, OPVs have the potential for high energy return on investment. Another advantage of organic materials is that they have a high optical absorption coefficient. Therefore, thin films can capture sufficient solar energy, thus reducing material cost. However, these advantages are offset by the challenges faced due to low stability, low efficiency, and process sensitivity. In one of the popular methods of device fabrication, the donor (D) and the acceptor (A) are dissolved in a suitable solvent and then spin-coated onto the substrate. During spin-coating, the solvent evaporates, resulting in the D–A materials spontaneously phase-transforming to D and A domains within the active organic-semiconductor layer. The principal process parameters that influence OPV performance include the annealing time and temperature, blend ratio of the donor to the acceptor material, and nature of the solvent.4,5

Until now, the optimization of OPVs has mostly been done on a trial and error basis. Another way of improving device performance is by a bottom-up approach, wherein a theoretical framework is developed to understand the process–structure–property relationship that is then utilized for optimizing device performance. Negi et al.6 have developed a phase-field model to simulate the microstructure evolution during spin-coating. Several theoretical structure–property models at the microscopic7–17 and continuum scale18–26 have been developed so far. In some of the models, the active layer or the BHJ has been treated as a material having effective electron and hole properties such as mobility and carrier density. The kinetic Monte Carlo (kMC) method has been implemented by various groups7,11–13,15,16,27,28 for deriving structure–property correlations. With respect to continuum models, Barker et al.18 developed a numerical model to calculate the current density–voltage characteristics for a bilayer polymer solar cell and were able to explain specific experimental trends and identify a few design rules for OPVs. A quantitative model to simulate the current density–voltage characteristics for polymer/fullerene BHJ solar cells, considering an active layer comprising a semiconductor having the LUMO of the acceptor as the valence band edge and the HOMO of the donor as the conduction band edge, was developed by Koster et al.19 Buxton and Clarke20 predicted structure–property correlation for a diblock copolymer morphology using a two-dimensional model. Furthermore, the optimization of device thickness and interfacial area for optimal device performance was carried out by Martin et al.21 by using the drift-diffusion model for a nanostructured heterojunction. A two-dimensional optical and electrical model was formulated by Williams and Walker22 for simulating the electrical characteristics of an interdigitated morphology. Kirchartz et al.23 developed a one-dimensional electro-optical model for BHJ solar cells, considering the active layer to be a homogeneous material. The structure–property relationship for two-dimensional rod-coil block copolymers was determined by Shah and Ganesan.24 Van Der Holst et al.29 modeled the single-carrier current density in three-dimensional sandwich-type organic semiconductor devices. The effect of blend stoichiometry, carrier density, electric field, and disorder on charge transportation in organic blends was explored by Koster30 by solving the Pauli master equation.31 Ray et al.26 applied the model developed by Buxton and Clarke20 to simulate the effect of annealing on device performance for one particular blend ratio of the donor and acceptor. Kodali and Ganapathysubramanian32 extended the formulation developed by Buxton and Clarke33 and implemented a new formulation involving a transformation of variables for computing the net current density.

In the following, we briefly describe the working mechanism within the active organic layer in order to appreciate the importance of morphology and, consequently, its effect on electronic properties.

1. Charge conduction mechanism

The charge conduction mechanism in the active layer of BHJ-OPV (Fig. 1) is as follows:

  1. As light shines on the OPV, an exciton (a hole and electron pair bound together by coulombic force) is created predominantly in the donor material.

  2. The exciton diffuses at the D–A interface where it experiences an energy level offset.

  3. The exciton dissociates into a charge transfer (CT) complex which then further dissociates into an electron and a hole. The dissociation is due to the electric field created by the energy difference between the lowest unoccupied molecular orbital (LUMO) of the p-type and the n-type material. The diffusion length of the exciton is about 10 nm–20 nm.

  4. The charge carriers are then carried away to their respective electrodes.

FIG. 1.

Schematic depicting the working mechanism of the BHJ OPV.

FIG. 1.

Schematic depicting the working mechanism of the BHJ OPV.

Close modal

In order to optimize the steps as mentioned above, strategies such as tailor-making the molecules and copolymers to get desired bandgaps and energy levels have been adopted. However, in this work, the focus is on the role of BHJ morphology on OPV performance.

2. Morphology engineering

Having understood the working mechanism in the active layer, the importance of morphology can be better appreciated. As is evident from the charge conduction mechanism of OPVs, the key to an efficient solar cell is the optimization of morphology (also referred to as morphology engineering).34 The point to note here is that there are two important processes involved in the working of OPVs, i.e., exciton dissociation and carrier conduction through the morphology to the respective electrodes. Therefore, optimizing the morphology in order to maximize exciton dissociation and transportation of the resulting electrons and holes to their respective electrode via a well-connected network of phases with minimal leakage of current and recombination of charge carriers is key to an efficient solar cell.

In this work, therefore, we have a twofold objective. The first objective is to derive a process–structure correlation that describes the variation in the microstructure with change in the blend ratios. Here, we use a diffuse interface theoretical framework for simulating microstructures that would occur during spinodal decomposition of blends with different compositions. The second objective is to derive the structure–property correlations corresponding to the simulated microstructures. Here, again, we have adopted a diffuse interface model to incorporate heterogeneities in electronic properties for determining the structure–property correlations, which, in this case, is the current density–voltage relationship. In addition, the dynamic variable in our model is electrochemical potential instead of concentration, and this leads to a more straightforward numerical treatment since, as the steady state approaches, at the interface, the gradient in the electrochemical potential diminishes. In contrast, the concentration sees variations differing by orders of magnitude. Another advantage of the diffuse interface approach is that the boundary conditions at the interface are imposed implicitly in the model, thus allowing easier computational treatment of complex geometries, as in the case of BHJ OPVs. In this model, the electron transport layer (ETL) and the hole transport layer (HTL) have also been considered in the device architecture since they are crucial for deriving optimum OPV performance. Using the combined diffuse interface framework, we derive PSP linkages that correlate the process of device fabrication to the BHJ morphology and its resultant properties. The soundness of the framework is then verified by performing experiments on a model P3HT:PCBM system. The organization of this paper is as follows: Sec. II details the modeling framework comprising the process–structure and structure–property models used in this work. In Sec. III, we have benchmarked the simulation results from our diffuse interface electronic model against the analytical solution for a bilayer device, consisting of p-type and n-type semiconductors sandwiched between electrodes. Thereafter, current density–voltage (J–V) characteristics for different BHJ morphologies in 2D and 3D are simulated. We then derive the influence of blend ratio and annealing period on device performance. Section IV describes the experimental verification of the theoretical PSP correlation. The morphological characteristics that optimize device performance are analyzed in Sec. V. Finally, we have listed the assumptions made in our model, their consequence on the results obtained, and the future direction of work.

In order to understand and quantify the morphological features that influence device performance, we modeled the process–structure–property relationship in organic semiconductors. For investigating the effect of processing parameters (blend ratio and annealing time) on the resulting microstructure, the phase-field model has been used since it is most appropriate for modeling microstructure evolution. It has been assumed that the polymer–fullerene blend undergoes spinodal decomposition.6,20,26,35 In the phase-field model, the free energy of the polymer blend has been modeled using the Flory–Huggins function (2),36 and the Cahn–Hilliard formulation37 is used to model the free-energy functional (1). The Cahn–Hilliard free energy functional is as follows:

(1)

where ϕi represents the volume fraction of the polymer (p), fullerene (f), or solvent (s), and κi is the gradient energy coefficient of p, f, or s. The Flory–Huggins bulk free energy density f is

(2)

where Vm represents the reference molar volume, T is the temperature, R is the universal gas constant, Ni is the relative size of the polymer, fullerene, or solvent, and χpf is the Flory interaction parameter between the polymer and the fullerene. Upon extending the formulation implemented by Bhattacharyya38 to the ternary polymer blend, the evolution equations for the volume fraction of the polymer and fullerene are

(3)
(4)

where the mobilities are derived using the following relations:

(5)

Λp, Λf, and Λs represent the chemical mobility of the polymer, fullerene, and solvent, respectively.

1. Non-dimensionalization

We utilize the following scales for energy f*, length l*, and time t* for non-dimensionalizing the free-energy density and the kinetic parameters: f*=RT*Vm, l*=σ*f*, and t*=l*2D*. In our system, the annealing temperature is 423 KT*, the molar volume Vm of chlorobenzene is 101 cm3/mol, the surface energy scale σ* is 25 mN/m, and the diffusivity scale D* is 5.0 × 10−9 m2/s. With this scale, the following set of non-dimensional parameters (see Table I) has been utilized for simulating the microstructures during spinodal decomposition, while the phase diagram is depicted in Fig. 2(a). In order to generate a range of morphologies corresponding to compositions that range from P3HT rich to PCBM rich, one particular tie-line in the phase diagram of P3HT, the solvent, and the PCBM ternary blend has been chosen [see Fig. 2(a)], and morphologies corresponding to various blend ratios of P3HT and PCBM have been generated [see Figs. 2(b)–2(f)].

TABLE I.

Process–structure model parameter table.

Parameterχpfκp,f,sΛaaΛbbΛabNpNfNs
Value 0.15 0.1 1.0 1.0 0.3 100 10 
Parameterχpfκp,f,sΛaaΛbbΛabNpNfNs
Value 0.15 0.1 1.0 1.0 0.3 100 10 
FIG. 2.

(a) Isothermal section of the phase diagram of P3HT, PCBM, and CB. BHJ morphology with a blend ratio of (b) 1:1.25, (c) 1:1.5, (d) 1:2.0, (e) 1:2.2, and (f) 1:5.4. On the color bar, 0 corresponds to the PCBM rich phase, whereas 1 corresponds to the P3HT rich phase.

FIG. 2.

(a) Isothermal section of the phase diagram of P3HT, PCBM, and CB. BHJ morphology with a blend ratio of (b) 1:1.25, (c) 1:1.5, (d) 1:2.0, (e) 1:2.2, and (f) 1:5.4. On the color bar, 0 corresponds to the PCBM rich phase, whereas 1 corresponds to the P3HT rich phase.

Close modal

In order to model the microstructure–property relationship, we begin with the discussion of basic semiconductor equations.

1. Semiconductor equations

The equations governing carrier transportation are

(6)
(7)

where (6) represents Poisson’s equation, (7) represents the net current density of carriers due to drift and diffusion, i represents the charged species that is either an electron or a hole, j represents the phase that is p-type or n-type, ϵj is the dielectric constant of the phase, ξ is the electrostatic potential, qi is the charge density, given by zF, ci corresponds to the carrier density, NA and ND represent the acceptor and donor densities, respectively, Ji represents the net current density of carriers, Di represents the diffusivity of the charged species, and λi is the electronic mobility of ith species, where the mass conservation equation is written in terms of the composition field of the charged species. This is the classical sharp interface model. However, the shortcoming of this model is that it works well for simple geometries, but its application to complex morphologies is difficult. In order to deal with complex morphologies and low concentration of carriers, as is typical in organic semiconductors (OS), we adopt the diffuse interface model. Since the thickness of the active layer is typically 100 nm and the entire morphology is depleted, the depletion width is 60 nm,39 and the concentration of carriers is low; it becomes easier to deal with the equations numerically if the effective chemical potential (diffusion potential) is used instead of concentration. Therefore, if the mass conservation equation is written in terms of electrochemical potential (which is the sum of the diffusion potential and qiξ), then, as the steady state approaches, at the interface, the gradient in the electrochemical potential diminishes. In contrast, the concentration sees variations differing by orders of magnitude. Therefore, the numerical treatment of the model is much easier if electrochemical potential is used instead of composition. In order to write the evolution equations in terms of electrochemical potential, we define the electronic free energy of the system in Sec. II B 2.40 

2. Free-energy

The electronic free energy of the system, considering that electrons and holes are much lower in concentration than the host atoms, is given by

(8)

where Ei is the Fermi level. Therefore, the electrochemical diffusion potential of ith species (μ̃i), which is the sum of the chemical diffusion potential μi and the electrostatic potential qiξ, turns out to be μĩ=RTlnci+Ei+qiξ. Thus, the composition can be calculated as ci=expμiEiRT.

3. Model formulation

In this section, we will describe the drift-diffusion model governing the transportation of electrons and holes. The mass conservation equation given by the following equation and Poisson’s equation given by (6) are solved in a coupled manner in order to derive device characteristics. The mass conservation equation is

(9)

where μ̃i is the electrochemical diffusion potential, μi is the diffusion potential, Mi is the mobility (for relaxation of electrochemical potential) of the ith species, which is interpolated in the diffuse interface model as Mi=ψMiα+(1ψ)Miβ, where Mij=Dijcijμij, where μij is the chemical diffusion potential of the ith species in the jth phase, and Ġ and represent the generation and recombination rate of carriers, respectively. Ġ is given by G0(6ψ(1 − ψ)), where G0 is the coefficient of the generation rate of carriers, and ψ represents the phase field order parameter that varies smoothly between 0 ≤ ψ ≤ 1. Here, it is assumed that the dissociation probability of the excitons is one at the interface, demarcated by the region where [ψ(1 − ψ)] is non-zero; thereby, the exciton generation rate is directly the carrier generation rate. In addition, in this particular formulation, the diffusion of excitons is not considered. This is because the microstructural length scale is comparable to the typical exciton diffusion length. is modeled using the modified Langevin recombination model in low mobility materials, given by

(10)

where ce and ch represent the electron and hole concentration, respectively, and ni denotes the intrinsic carrier concentration. The recombination coefficient γ is given by the following equation 41:

(11)

4. Non-dimensionalization for the structure–property model

The carrier densities are non-dimensionalized by using a density scale of N* = 3.99 × 10−2 mol/m3; thereby, the charge density ρ is scaled by zFN*, where z is the valency of the atom, and F is Faraday’s constant. The electrostatic potential is normalized by [RT/(zF)]. With this, the natural length scale of the problem can be derived from the Poisson equation, which is l*=(RT)1/2(ϵ)1/2N*1/2(zF), where ϵ is the permittivity of one of the semiconducting phases in the vacuum. Therefore, the scale for the current density becomes (zF)l*(Me+Mh)(RT), where T is the temperature, R is the universal gas constant, and Mi represents the mobility (with units mol2 kg−1 m−3 s). Dimensionalized quantities can be retrieved by multiplying the normalized quantities, as mentioned in Table II, with their respective dimensional counterparts, where the typical numbers for a P3HT:PCBM solar cell are as follows: Vbi—0.7 V, M—1.1 × 10−11 mol2 kg−1 m−3 s, W—60 × 10−9 m, ϵ—3.1 × 10−11 F m−1, and T—300 K.

TABLE II.

Structure–property model parameter table.

ParameterValue
DeP3HT 10 
DePCBM 15 
DhP3HT 15 
DhPCBM 10 
ceP3HT 5 × 10−5 
cePCBM 5 × 10−3 
chP3HT 5 × 10−3 
chPCBM 5 × 10−5 
G0 5 × 10−6 
R0 
Vbi 
ϵ 
ParameterValue
DeP3HT 10 
DePCBM 15 
DhP3HT 15 
DhPCBM 10 
ceP3HT 5 × 10−5 
cePCBM 5 × 10−3 
chP3HT 5 × 10−3 
chPCBM 5 × 10−5 
G0 5 × 10−6 
R0 
Vbi 
ϵ 

In this section, the structure–property diffuse interface model has been benchmarked against analytical results obtained from the classical sharp interface model. Thereafter, BHJ device characteristics have been simulated using the model described in Sec. II B 3.

We simulated a bilayer device, i.e., one consisting of p-type and n-type semiconductors sandwiched between electrodes (see Fig. 3), to verify whether the device physics are being captured by the model.

FIG. 3.

Morphology of the bilayer device. 1 represents p-type OS, and 0 represents n-type OS with the junction having a diffuse profile.

FIG. 3.

Morphology of the bilayer device. 1 represents p-type OS, and 0 represents n-type OS with the junction having a diffuse profile.

Close modal

In order to benchmark the model, the electrostatic potential and charge density profile as obtained from the numerical model and analytically were compared against each other in the open-circuit condition. Excellent agreement between the two is observed [Figs. 4(a) and 4(b)]. As the p- and n-type semiconductors are brought in contact, the holes flow from the p-type to the n-type semiconductor, as governed by the chemical potential gradient, whereas the electrons flow in the opposite direction. The holes from the p-type semiconductor recombine with the electrons in the n-type semiconductor, resulting in exposed positive charge at the n-side of the device due to donor atoms (the mobile carriers were previously shielding this charge). Simultaneously, electrons from the n-type semiconductor recombine with the holes in the p-type semiconductor, resulting in exposed negative charge at the p side of the device due to acceptor atoms. This diffusion of carriers leads to an electrostatic potential build-up (known as the built-in potential Vbi) across the PN junction such that it opposes the driving force due to the chemical potential gradient. Equilibrium is reached when the electrostatic potential gradient opposes the chemical potential gradient. In the case of charge density, as depicted in Fig. 4(b), the analytical results and those of the model match well, except at the PN junction; this is because the analytical equation considers a sharp interface where there is naturally a discontinuity in the charge density, whereas the variation is continuous in the diffuse interface model. Thereafter, the device is simulated in forward and reverse bias in order to observe the variation in the depletion width with bias [Fig. 4(c)]. In the case of forward bias, the electrostatic potential gradient aids the chemical potential gradient, resulting in reduction in the depletion width. In contrast, in reverse bias, the electrostatic potential gradient opposes the chemical potential gradient, resulting in a majority of carriers flowing away from the PN junction, thereby increasing the width of the depletion region. It should be noted that all the above profiles have a diffuse interface. The advantage of adopting the diffuse interface approach is that the boundary conditions at the interface are imposed implicitly by the model, thus allowing an easier computational treatment of complex geometries, as in the case of BHJ OPVs. Thereafter, the J–V characteristics for the bilayer device under dark and light conditions were simulated. Non-linear J–V characteristics were observed, as is typical of semiconductor devices. The nonlinearity arises from the fact that the number of mobile carriers is a function of the applied bias that increases with the magnitude of the forward bias and vice-versa with that of the reverse bias. In addition, the effect of light intensity on this bilayer device was simulated [see Fig. 5(b)]. As expected, with an increase in light intensity, the maximum current density that can be drawn from the device and the open circuit voltage increase.

FIG. 4.

(a) Electrostatic potential profile, (b) charge density profile, and (c) variation in the depletion width with bias.

FIG. 4.

(a) Electrostatic potential profile, (b) charge density profile, and (c) variation in the depletion width with bias.

Close modal
FIG. 5.

(a) Dark J–V characteristics and (b) effect of light intensity on device characteristics. Go is the coefficient of the generation rate of carriers.

FIG. 5.

(a) Dark J–V characteristics and (b) effect of light intensity on device characteristics. Go is the coefficient of the generation rate of carriers.

Close modal

Having benchmarked the model for a simple bilayer device, we further went on to simulate the J–V characteristics for a complex BHJ morphology. In the present morphology [see Fig. 6(a)], the blend ratio is 1:2.5 (p-type:n-type OS). In order to model the device characteristics, the active layer (polymer–fullerene blend) is placed between the buffer layers, as is typical in OPVs. The buffer layer aids in allowing the transportation of only one charge carrier (i.e., electron or hole), blocking the other one. Hence, the p-type semiconducting buffer layer is placed near the electrode responsible for extracting holes, and the n-type semiconducting buffer layer is placed near the electrode responsible for collecting electrons. J–V characteristics of this morphology in the dark and light are shown in Fig. 6(b).

FIG. 6.

(a) BHJ morphology with a blend ratio of 1:2.5 and (b) J–V characteristics of the device with a BHJ morphology of the active layer.

FIG. 6.

(a) BHJ morphology with a blend ratio of 1:2.5 and (b) J–V characteristics of the device with a BHJ morphology of the active layer.

Close modal

The simulations mentioned in Secs. III A and III B are performed for 2D microstructures using an in-house C code. Since 3D microstructures give a more realistic understanding of the structure–property correlation, the device characteristics for the same were simulated. Given the computationally intensive nature of the problem, we have utilized the OpenFOAM software that is preloaded with efficient multigrid solvers.42 The 3D microstructures are shown in Fig. 7. In order to derive the effect of morphology on device characteristics, the J–V characteristics for the morphologies shown in Fig. 7 were simulated [see Fig. 8(a)]. Using these characteristics, the power conversion efficiency (PCE) of each device was calculated, and the variation in efficiency with the volume fraction of the donor phase can be seen in Fig. 8(b). In order to verify these theoretical results, experiments are performed, and the current density–voltage characteristics are measured for devices with various blend ratios and annealing periods.

FIG. 7.

BHJ morphology with a blend ratio of (a) 1:5.4, (b) 1:3.5, (c) 1:2.0, (d) 1:1.5, (e) 1:1.25, and (f) 1:0.8.

FIG. 7.

BHJ morphology with a blend ratio of (a) 1:5.4, (b) 1:3.5, (c) 1:2.0, (d) 1:1.5, (e) 1:1.25, and (f) 1:0.8.

Close modal
FIG. 8.

(a) J–V characteristics of morphologies with different volume fractions of P3HT and (b) efficiency vs volume fraction.

FIG. 8.

(a) J–V characteristics of morphologies with different volume fractions of P3HT and (b) efficiency vs volume fraction.

Close modal

We have chosen a model P3HT:PCBM system as the active layer to verify the theoretically derived PSP correlations. The fabricated device architecture is ITO/PEDOT:PSS/P3HT:PCBM/Al [see Fig. 9(a)]. The devices have an active area of 6 mm2. The device fabrication procedure is as follows: ITO-coated glass substrates are sonicated in soap water, DI water, acetone and then IPA each for 10 min. They are dried using compressed air and then treated by UV-ozone for 20 min. Immediately after the UV-ozone treatment, PEDOT:PSS is spin-coated at 4500 rpm for 1 min and then annealed at 150 °C for 30 min. After that, the substrates are transferred to a JACOMEX nitrogen glovebox, and the active layer is spin-coated at 1100 rpm for 1 min. P3HT obtained from Rieke metals (4002-EE) has a molecular weight of 52 kg/mol, a polydispersity of 2.4, and a regio-regularity of 91%. It is dissolved with PCBM (from Nano-C) at 60 °C overnight in chlorobenzene (CB). The polymer concentration is 10 mg/ml. The substrates are then thermally annealed at 150 °C for 10 min. 100 nm of aluminum is evaporated at approx 1.5 Å/s using an Angstrom engineering thermal evaporator. The devices are then characterized using a solar simulator. The effect of blend ratio on device efficiency is depicted in Fig. 9(b), where we find that the trend with the change in blend ratios is in good qualitative agreement with the simulated results. Upon annealing the device at 150 °C, it is observed that Jsc first increases and then decreases with time. In order to gain insights into the nature of these trends, the simulated data are analyzed in Sec. V.

FIG. 9.

(a) Device architecture, (b) efficiency vs weight % of PCBM, and (c) Jsc vs annealing time. Note: The error bars represent the standard deviation.

FIG. 9.

(a) Device architecture, (b) efficiency vs weight % of PCBM, and (c) Jsc vs annealing time. Note: The error bars represent the standard deviation.

Close modal

As was explained in the working of OPVs, the key to efficient solar cell performance is to maximize exciton dissociation and carrier conduction through the morphology to the respective electrodes. In other words, we would want to maximize the interfacial area (since excitons dissociate at the interface) while having percolating channels to the respective electrodes. In order to investigate the variation of Jsc and derive any possible correlations with changes in the morphology, the interfacial area and percolation of the various simulated microstructures were calculated and compared with their device performance. The percolation of the morphology is derived using the Hoshen–Kopelman algorithm.43 We define a parameter called percolation fraction to quantify percolation of p-type and n-type semiconductor clusters. The percolation fraction is calculated by summing the fraction of p-type and n-type percolating clusters; for example, if either the p-type or n-type semiconductor is completely percolating, the percolation fraction would be 0.5, whereas if both the semiconductors are completely percolating, the percolation fraction would be 1. The device performance of semi-percolating morphologies, i.e., either P3HT rich (volume fraction—0.75) or PCBM rich (volume fraction—0.20), is found to be poor (volume fractions are of P3HT). All other morphologies are bi-continuous, thereby implying that the device performance is controlled by exciton dissociation. This is precisely what is observed, i.e., Jsc and Voc for different morphologies follow the same trend as that of the interfacial area (see Fig. 10). It can be noted here that the variation in Voc is relatively small, as is observed experimentally in organic photovoltaics. The correlation in the trend can be explained by the fact that since the morphologies corresponding to a volume fraction of 0.30, 0.41, 0.44, 0.49, 0.54, and 0.59 are completely percolating, charge carrier conduction to electrodes will not be a limiting factor; therefore, if the interfacial area increases, the generation of carriers becomes larger, resulting in a higher Jsc and Voc.

FIG. 10.

(a) Percolation fraction vs volume fraction, (b) Jsc and interfacial area vs volume fraction, and (c) Voc and interfacial area vs volume fraction. Note: The lines have been drawn just to follow the points and have no physical significance.

FIG. 10.

(a) Percolation fraction vs volume fraction, (b) Jsc and interfacial area vs volume fraction, and (c) Voc and interfacial area vs volume fraction. Note: The lines have been drawn just to follow the points and have no physical significance.

Close modal

It was experimentally observed that initially during annealing, Jsc increases and then decreases, thus resulting in an optimal annealing time [see Fig. 9(c)]. The initial increase in Jsc corresponds to the existence of a better interface between the electrode and the active layer due to post-production annealing of the device.44 The improvement in the interface is evident from the fact that there is an enhancement in the fill factor. A better interface naturally leads to an enhancement in carrier collection at the electrodes. The reduction in Jsc upon further annealing stems from coarsening of the morphology, resulting in lower exciton dissociation. The simulation results indicate that as the morphology is annealed, due to coarsening (resulting in a lower interfacial area), Jsc decreases (see Fig. 11).

FIG. 11.

Effect of annealing on Jsc.

FIG. 11.

Effect of annealing on Jsc.

Close modal

The results from the theoretical framework seem to follow trends similar to those from the experimental measurements. Additionally, the simulations allow us to derive correlations between the changes in the microstructural features and changes in the composition and annealing to the effective efficiency of the devices. However, some assumptions presently limit a quantitative match between experiments and simulations.

First, in the process–structure modeling framework, it has been assumed that spinodal decomposition is the mechanism of phase-transformation, given the theoretical and experimental works by several groups.6,20,26,35,45,46 However, some groups have carried out in situ studies on drop-cast samples and contend that the mechanism is nucleation and growth.47,48 It is important to note that the kinetics of phase-transformation change drastically upon spin-coating as the thickness of the film reduces to 100 nm from 1 mm (as obtained by drop-casting); moreover, this confinement of the morphology inhibits PCBM crystallization. Therefore, the conclusions drawn from in situ drop-casting studies about the mechanism of phase-transformation would differ if they were carried out while spin-coating. Instead of spinodal decomposition, if crystallization-induced phase transformation is assumed, the morphology picture would be different. At low PCBM loadings, P3HT crystallites would form in an amorphous matrix of P3HT and PCBM. Owing to the lack of crystalline phases (and, therefore, paths for carriers to travel to electrodes), such devices would exhibit lower PCE. Intermediate PCBM loading would result in the formation of three phases: crystalline P3HT, crystalline PCBM, and an amorphous phase consisting of intermixed P3HT and PCBM. These devices would exhibit the highest PCE. A higher PCBM loading tends to suppress P3HT crystallinity, and the PCE would see a setback. These aspects are certainly exciting and form the motivation for future work, wherein apart from morphology engineering, an additional complexity is the number of different phases with different electronic properties that may arise out of the fabrication process.

Second, we would like to mention here that exciton diffusion has not been considered in our model and the exciton density has been treated to be constant throughout the morphology. This assumption is justified since all morphologies discussed here have domains within 20 nm (exciton diffusion length), and therefore, it is only the interfacial area and percolation of the phases that affect device performance. Therefore, for the present situation, this is a reasonable approximation but can easily be incorporated in the model for a future theoretical exercise wherein the microstructural scales are larger.

Third, while we have simulated different microstructures derived from the spinodal decomposition of blends with varying compositions, in reality, the composition shift into the spinodal region of the phase diagram is triggered due to enhanced evaporation of the solvent. Evaporation would naturally bring in a variation in the microstructural scales in the vertical direction in the context of thin-film device fabrication, wherein only the top surface is exposed to the atmosphere. Negi et al.6 have carried out comprehensive work in this direction, and the incorporation of this additional complexity into our model is part of our future work.

Fourth, all domains (or phases) considered in the model are either pure P3HT or PCBM. Chemical inhomogeneity that results in domains having semiconducting characteristics intermediate to pure p-type and n-type has been ignored. In our model, we have considered semiconductor characteristics of only pure phases and not mixed ones. The incorporation of this chemical inhomogeneity in the model requires additional information on the relation between the electronic properties and the composition of the blend. This could bring in additional complexity, and the theoretical framework will be useful in deriving critical insights into the blend compositions with optimal electronic properties.

In this work, we have, for the first time, developed a diffuse interface framework for deriving both the structure–property and process–structure relationships that provide for an elegant methodology for deriving PSP linkages. The theoretical PSP linkage in OPVs was established, and we were able to get guidelines on the processing parameters for optimizing the BHJ morphology and, subsequently, electronic properties. This linkage is found to be in qualitative agreement with our experimental results. We find that the necessary condition for fabricating an efficient OPV device is the existence of a bi-continuous network of donor and acceptor phases as this leads to percolating channels for electrons and holes to their respective electrodes. Among the percolating morphologies, the ones with a higher interfacial area lead to higher short circuit current density (Jsc) and eventually to a higher efficiency since the open-circuit potential (Voc) is a weak function of the morphology. Therefore, the key to an efficient solar cell is the optimization of the two critical processes involved in the working mechanism of OPV, i.e., exciton dissociation and carrier conduction through the morphology to the respective electrodes. The interfacial area governs exciton dissociation. Efficient carrier conduction results from the transportation of the electron and the hole to their respective electrode via a well-connected network of phases with minimal leakage of current and recombination of charge carriers. Hence, efficient carrier conduction is governed by a higher percolation fraction of the morphology. Another essential step during device fabrication is annealing the device after production. It is experimentally observed that initially during annealing, Jsc increases and then decreases, thus resulting in an optimal anneal time. The initial increase in Jsc corresponds to the existence of a better interface between the electrode and the active layer. The reduction in Jsc stems from the coarsening of morphology, resulting in lower exciton dissociation that has been verified numerically by the simulations.

The raw/processed data required to reproduce these findings cannot be shared at this time as the data are part of an ongoing study.

The authors would like to acknowledge the support from the Ministry of Communication and Information Technology under a grant from the Centre of Excellence in Nanoelectronics, Phase II, Grant No. MITO/EEC/NKB/0096, and Grant No. DST-SERB-O1679. The authors are grateful to the TUE-CMS and the SERC, IISc, for the provision of high-performance computing resources. F.K. is grateful to his labmates Sumeet, Bhalchandra, Sandeep, Jagdish, Arul, and Arun at the IISc, as well as Upendra Pandey, Shiv Nadar University, for insightful and critical discussions.

1.
F. C.
Krebs
, “
Fabrication and processing of polymer solar cells: A review of printing and coating techniques
,”
Sol. Energy Mater. Sol. Cells
93
(
4
),
394
412
(
2009
).
2.
G.
Yu
,
J.
Gao
,
J. C.
Hummelen
,
F.
Wudl
, and
A. J.
Heeger
, “
Polymer photovoltaic cells: Enhanced efficiencies via a network of internal donor-acceptor heterojunctions
,”
Science
270
(
5243
),
1789
1791
(
1995
).
3.
X.
Yang
,
J.
Loos
,
S. C.
Veenstra
,
W. J. H.
Verhees
,
M. M.
Wienk
,
J. M.
Kroon
,
M. A. J.
Michels
, and
R. A. J.
Janssen
, “
Nanoscale morphology of high-performance polymer solar cells
,”
Nano Lett.
5
(
4
),
579
583
(
2005
).
4.
Y.
Yao
,
J.
Hou
,
Z.
Xu
,
G.
Li
, and
Y.
Yang
, “
Effects of solvent mixtures on the nanoscale phase separation in polymer solar cells
,”
Adv. Funct. Mater.
18
(
12
),
1783
1789
(
2008
).
5.
J.
Alstrup
,
M.
Jørgensen
,
A. J.
Medford
, and
F. C.
Krebs
, “
Ultra fast and parsimonious materials screening for polymer solar cells using differentially pumped slot-die coating
,”
ACS Appl. Mater. Interfaces
2
(
10
),
2819
2827
(
2010
).
6.
V.
Negi
,
O.
Wodo
,
J. J.
van Franeker
,
R. A. J.
Janssen
, and
P. A.
Bobbert
, “
Simulating phase separation during spin coating of a polymer-fullerene blend: A joint computational and experimental investigation
,”
ACS Appl. Mater.
1
(
2
),
725
735
(
2018
).
7.
P. K.
Watkins
,
A. B.
Walker
, and
G. L. B.
Verschoor
, “
Dynamical Monte Carlo modelling of organic solar cells: The dependence of internal quantum efficiency on morphology
,”
Nano Lett.
5
(
9
),
1814
1818
(
2005
).
8.
R. A.
Marsh
,
C.
Groves
, and
N. C.
Greenham
, “
A microscopic model for the behavior of nanostructured organic photovoltaic devices
,”
J. Appl. Phys.
101
(
8
),
083509
(
2007
).
9.
C. R.
McNeill
,
S.
Westenhoff
,
C.
Groves
,
R. H.
Friend
, and
N. C.
Greenham
, “
Influence of nanoscale phase separation on the charge generation dynamics and photovoltaic performance of conjugated polymer blends: Balancing charge generation and separation
,”
J. Phys. Chem. C
111
(
51
),
19153
19160
(
2007
).
10.
C.
Groves
,
R. A.
Marsh
, and
N. C.
Greenham
, “
Monte Carlo modeling of geminate recombination in polymer-polymer photovoltaic devices
,”
J. Chem. Phys.
129
(
11
),
114903
(
2008
).
11.
L.
Meng
,
Y.
Shang
,
Q.
Li
,
Y.
Li
,
X.
Zhan
,
Z.
Shuai
,
R. G.
Kimber
, and
A. B.
Walker
, “
Dynamic Monte Carlo simulation for highly efficient polymer blend photovoltaics
,”
J. Phys. Chem. B
114
(
1
),
36
41
(
2009
).
12.
C.
Groves
,
R. G. E.
Kimber
, and
A. B.
Walker
, “
Simulation of loss mechanisms in organic solar cells: A description of the mesoscopic Monte Carlo technique and an evaluation of the first reaction method
,”
J. Chem. Phys.
133
(
14
),
144110
(
2010
).
13.
R. G. E.
Kimber
,
A. B.
Walker
,
G. E.
Schröder-Turk
, and
D. J.
Cleaver
, “
Bicontinuous minimal surface nanostructures for polymer blend solar cells
,”
Phys. Chem. Chem. Phys.
12
(
4
),
844
851
(
2010
).
14.
C.-K.
Lee
,
C.-W.
Pao
, and
C.-W.
Chu
, “
Multiscale molecular simulations of the nanoscale morphologies of P3HT: PCBM blends for bulk heterojunction organic photovoltaic cells
,”
Energy Environ. Sci.
4
(
10
),
4124
4132
(
2011
).
15.
A.
Pershin
,
S.
Donets
, and
S. A.
Baeurle
, “
A new multiscale modeling method for simulating the loss processes in polymer solar cell nanodevices
,”
J. Chem. Phys.
136
(
19
),
194102
(
2012
).
16.
F.
Wei
,
L.
Liu
,
L.
Liu
, and
G.
Li
, “
Multiscale modeling and simulation for optimizing polymer bulk heterojunction solar cells
,” in
2012 IEEE 38th Photovoltaic Specialists Conference (PVSC) PART 2
(
IEEE
,
2012
), pp.
1
10
.
17.
X.
Meng
,
L.
Zhang
,
Y.
Xie
,
X.
Hu
,
Z.
Xing
,
Z.
Huang
,
C.
Liu
,
L.
Tan
,
W.
Zhou
,
Y.
Sun
 et al., “
A general approach for lab-to-manufacturing translation on flexible organic solar cells
,”
Adv. Mater.
31
(
41
),
1903649
(
2019
).
18.
J. A.
Barker
,
C. M.
Ramsdale
, and
N. C.
Greenham
, “
Modeling the current-voltage characteristics of bilayer polymer photovoltaic devices
,”
Phys. Rev. B
67
(
7
),
075205
(
2003
).
19.
L. J.
Koster
,
E.
Smits
,
V.
Mihailetchi
, and
P.
Blom
, “
Device model for the operation of polymer/fullerene bulk heterojunction solar cells
,”
Phys. Rev. B
72
(
8
),
085205
(
2005
).
20.
G. A.
Buxton
and
N.
Clarke
, “
Predicting structure and property relations in polymeric photovoltaic devices
,”
Phys. Rev. B
74
(
8
),
085207
(
2006
).
21.
C. M.
Martin
,
V. M.
Burlakov
,
H. E.
Assender
, and
D. A. R.
Barkhouse
, “
A numerical model for explaining the role of the interface morphology in composite solar cells
,”
J. Appl. Phys.
102
(
10
),
104506
(
2007
).
22.
J.
Williams
and
A. B.
Walker
, “
Two-dimensional simulations of bulk heterojunction solar cell characteristics
,”
Nanotechnology
19
(
42
),
424011
(
2008
).
23.
T.
Kirchartz
,
B. E.
Pieters
,
K.
Taretto
, and
U.
Rau
, “
Electro-optical modeling of bulk heterojunction solar cells
,”
J. Appl. Phys.
104
(
9
),
094513
(
2008
).
24.
M.
Shah
and
V.
Ganesan
, “
Correlations between morphologies and photovoltaic properties of rod- coil block copolymers
,”
Macromolecules
43
(
1
),
543
552
(
2009
).
25.
W.
Tress
,
M.
Furno
,
K.
Leo
, and
M.
Riede
, “
Numerical drift-diffusion modeling of organic solar cells in comparison with experimental data series
,” in
Numerical Simulation of Optoelectronic Devices
(
IEEE
,
2010
), pp.
33
34
.
26.
B.
Ray
,
P. R.
Nair
, and
M. A.
Alam
, “
Annealing dependent performance of organic bulk-heterojunction solar cells: A theoretical perspective
,”
Sol. Energy Mater. Sol. Cells
95
(
12
),
3287
3294
(
2011
).
27.
M.
Casalegno
,
G.
Raos
, and
R.
Po
, “
Methodological assessment of kinetic Monte Carlo simulations of organic photovoltaic devices: The treatment of electrostatic interactions
,”
J. Chem. Phys.
132
(
9
),
094705
(
2010
).
28.
V.
Negi
, “
Exploring the process-structure-property relationship of organic solar cells by high performance computer simulations
,” Ph.D. thesis (
Technische Universiteit Eindhoven
,
Eindhoven
,
2018
).
29.
J.
Van Der Holst
,
M.
Uijttewaal
,
B.
Ramachandhran
,
R.
Coehoorn
,
P.
Bobbert
,
G.
De Wijs
, and
R.
De Groot
, “
Modeling and analysis of the three-dimensional current density in sandwich-type single-carrier devices of disordered organic semiconductors
,”
Phys. Rev. B
79
(
8
),
085203
(
2009
).
30.
L.
Koster
, “
Charge carrier mobility in disordered organic blends for photovoltaics
,”
Phys. Rev. B
81
(
20
),
205318
(
2010
).
31.
Z.
Yu
,
D.
Smith
,
A.
Saxena
,
R.
Martin
, and
A.
Bishop
, “
Molecular geometry fluctuations and field-dependent mobility in conjugated polymers
,”
Phys. Rev. B
63
(
8
),
085202
(
2001
).
32.
H. K.
Kodali
and
B.
Ganapathysubramanian
, “
Computer simulation of heterogeneous polymer photovoltaic devices
,”
Modell. Simul. Mater. Sci. Eng.
20
(
3
),
035015
(
2012
).
33.
G. A.
Buxton
and
N.
Clarke
, “
Computer simulation of polymer solar cells
,”
Modell. Simul. Mater. Sci. Eng.
15
(
2
),
13
(
2006
).
34.
N. E.
Jackson
,
B. M.
Savoie
,
T. J.
Marks
,
L. X.
Chen
, and
M. A.
Ratner
, “
The next breakthrough for organic photovoltaics?
,”
J. Phys. Chem. Lett.
6
(
1
),
77
84
(
2014
).
35.
O.
Wodo
and
B.
Ganapathysubramanian
, “
Modeling morphology evolution during solvent-based fabrication of organic solar cells
,”
Comput. Mater. Sci.
55
,
113
126
(
2012
).
36.
R. A. L.
Jones
,
Polymers at Surfaces and Interfaces
(
Cambridge University Press
,
1999
).
37.
J. W.
Cahn
, “
Phase separation by spinodal decomposition in isotropic systems
,”
J. Chem. Phys.
42
(
1
),
93
99
(
1965
).
38.
S.
Bhattacharyya
, “
Evolution of multivariant microstuctures with anisotropic misfit
,” Ph.D. thesis,
Indian Institute of Science
,
Bangalore
,
2009
.
39.
S. R.
Gollu
,
R.
Sharma
,
G.
Srinivas
,
S.
Kundu
, and
D.
Gupta
, “
Incorporation of SiO2 dielectric nanoparticles for performance enhancement in P3HT:PCBM inverted organic solar cells
,”
Org. Electron.
24
,
43
50
(
2015
); URL http://www.sciencedirect.com/science/article/pii/S1566119915002141.
40.
A.
Slimani
,
A.
Iratni
,
H.
Henry
,
M.
Plapp
,
J.-N.
Chazalviel
,
F.
Ozanam
, and
N.
Gabouze
, “
Macropore formation in p-type silicon: Toward the modeling of morphology
,”
Nanoscale Res. Lett.
9
(
1
),
585
(
2014
).
41.
L. J. A.
Koster
,
V. D.
Mihailetchi
, and
P. W. M.
Blom
, “
Bimolecular recombination in polymer/fullerene bulk heterojunction solar cells
,”
Appl. Phys. Lett.
88
(
5
),
052104
(
2006
).
42.
Openfoam user guide,
O Limited
,
2018
.
43.
J.
Hoshen
and
R.
Kopelman
, “
Percolation and cluster distribution. i. cluster multiple labeling technique and critical concentration algorithm
,”
Phys. Rev. B
14
,
3438
3445
(
1976
).
44.
D.
Chen
,
A.
Nakahara
,
D.
Wei
,
D.
Nordlund
, and
T. P.
Russell
, “
P3HT/PCBM bulk heterojunction organic photovoltaics: Correlating efficiency and morphology
,”
Nano Lett.
11
(
2
),
561
567
(
2010
).
45.
O.
Wodo
and
B.
Ganapathysubramanian
, “
How do evaporating thin films evolve? Unravelling phase-separation mechanisms during solvent-based fabrication of polymer blends
,”
Appl. Phys. Lett.
105
(
15
),
153104
(
2014
).
46.
O.
Wodo
,
S.
Tirthapura
,
S.
Chaudhary
, and
B.
Ganapathysubramanian
, “
A graph-based formulation for computational characterization of bulk heterojunction morphology
,”
Org. Electron.
13
(
6
),
1105
1113
(
2012
).
47.
W.
Yin
and
M.
Dadmun
, “
A new model for the morphology of P3HT/PCBM organic photovoltaics from small-angle neutron scattering: Rivers and streams
,”
ACS Nano
5
(
6
),
4756
4768
(
2011
).
48.
P.
Kohn
,
Z.
Rong
,
K. H.
Scherer
,
A.
Sepe
,
M.
Sommer
,
P.
Müller-Buschbaum
,
R. H.
Friend
,
U.
Steiner
, and
S.
Hüttner
, “
Crystallization-induced 10-nm structure formation in P3HT/PCBM blends
,”
Macromolecules
46
(
10
),
4002
4013
(
2013
).