Analytical expression for the current-voltage characteristics of organic bulk heterojunction solar cells

An expression to describe the current-voltage characteristics of organic bulk heterojunction (BHJ) solar cells is derived. The derivation is obtained by analytically solving the drift-diffusion model for organic BHJ solar cells with the assumption of uniform bimolecular recombination rate. The assumption of uniform bimolecular recombination rate leads to somewhat inaccurate, for example, carrier densities as functions of the position inside the device. However, we show that this assumption should still produce an expression for the current as a function of applied voltage as if the actual bimolecular recombination rate is considered in the derivation. Applying this analytical expression to experimental current-voltage data enable us to directly extract and analyze, for example, the recombination loss of an organic BHJ solar cell as a function of applied voltage.


I. INTRODUCTION
Organic bulk heterojunction (BHJ) solar cells have a huge potential as a source of clean energy in the future. 1,2In order to improve the performance, understanding of the current-voltage (J-V) characteristics is important.Organic solar cells differ from conventional inorganic solar cells.In organic solar cells, the absorbed light generates excitons first instead of free electrons and holes as in conventional inorganic solar cells.The excitons must be dissociated into free electrons and holes in order to produce useful electric current.Two materials (donor and acceptor) with different lowest unoccupied molecular orbitals (LUMOs) and different highest occupied molecular orbitals (HOMOs) are needed to assist the dissociation of excitons into free charge carriers.In order to generate and extract free charge carriers efficiently, the donor and the acceptor materials are blended together to create donor-acceptor interfaces throughout the device's active layer (called bulk heterojunction).
Drift-diffusion model with metal-insulator-metal picture 3 is a widely used model for organic BHJ solar cells. 4This model solves the continuity and the Poisson's equations, accounting for the drift and diffusion of the charge carriers, the effect of space charge, and the generation and loss of the charge carriers. 4,5A numerical procedure has been widely used in order to obtain the J-V characteristics from this model. 4,6][9][10] In previous analytical studies, [7][8][9] the electric field is assumed constant with negligible recombination.Hence, the J-V characteristics are obtained by adding the photogenerated and the dark current densities instead of solving the above mentioned drift-diffusion model.Recently, Chowdhury and Alam 10 proposed an analytical model based on the drift-diffusion model but did not incorporate the intrinsic bimolecular recombination.
The aim of this work is to derive an expression to describe the J-V characteristics of organic BHJ solar cells by solving the drift-diffusion model analytically, where the effect of bimolecular recombination is considered.Compared to previous analytical works, this expression should be more appropriate for analyzing the J-V characteristics of organic BHJ solar cells when the bimolecular recombination is not negligible.

II. DRIFT-DIFFUSION MODEL
First, let us present the operational mechanism for organic BHJ solar cells that is used in this study.After the light is absorbed, the generated excitons must reach the donor-acceptor interface before they can dissociate into free electrons and holes.At the interface, the exciton splits into an electron in the acceptor material and a hole in the donor material, where the electron and hole are presumed to form a bound geminate electron-hole pair (called polaron pair).The polaron pair needs to further dissociate into a free electron and a free hole in order to produce photocurrent.At the interface, the polaron pair can also decay, which is a geminate monomolecular (first order) recombination. 11Hence, the probability for a polaron pair to dissociate into free charge carriers can be defined as the ratio of the dissociation rate of a polaron pair to the sum of the dissociation rate and the decay rate of a polaron pair.The free electron in the acceptor and the free hole in the donor can also meet at the interface and recombine into a polaron pair again.This is a non-geminate bimolecular (second order) recombination where the recombination rate is the rate for free electrons and free holes finding each other, which can be predicted by the Langevin theory. 12here are also other possible recombination mechanisms (see.Ref. 11 for example) and their importance is still being studied and debated.It is generally believed that trap-assisted recombination can be significant if trap densities are high.For example, Cowan et al. 13 concluded that the recombination in organic BHJ solar cells evolves from being monomolecular at short circuit to being bimolecular at open circuit.This finding suggests that non-geminate monomolecular recombination (i.e.trap-assisted recombination) should be considered.However, Kniepert et al. 14 recently demonstrated that the recombination in their as-prepared and annealed poly(3-hexylthiophene) and [6,6]-phenyl-C 71 -butyric acid methyl ester (P3HT:PCBM) solar cells is strictly bimolecular and trap-assisted recombination cannot account for their results.In our study here, we assume the contribution of trap-assisted recombination in the studied organic BHJ solar cells is negligible.
Next, we present the basics of the drift-diffusion model for organic BHJ solar cells based on the operational mechanism described above.Details of the model can be referred to Refs. 4 and 5 for example.As usual, the transport, generation, and loss of the charge carriers are assumed to be 1-dimensional.The Poisson's equation is given by where F is the electric field, q is the elementary charge, ε is the effective permittivity of the organic active layer, p is the hole density, and n is the electron density.Based on the above operational mechanism, the continuity equations for electrons and holes in steady-state condition can be written as 4,15 where J n is the electron current density, J p is the hole current density, P is the dissociation probability of polaron pairs into free charge carriers, G is the generation rate of polaron pairs per unit volume, and R L is the non-geminate bimolecular recombination rate.The electron and the hole current densities consist of the drift and the diffusion components, and are given by 16 J n (x) = qn(x)µ n F(x) + qD n dn(x) dx , (3a) where µ n is the electron mobility, µ p is the hole mobility, D n is the electron diffusion coefficient, and D p is the hole diffusion coefficient.For simplicity, µ n and µ p are assumed to be constant, which is needed in order to solve the drift-diffusion model analytically.In a more realistic model, the carrier mobility is electric field, temperature, and carrier density dependent. 17The diffusion coefficients are given by where k is the Boltzmann constant, and T is the absolute temperature.Regarding the non-geminate bimolecular recombination in organic BHJ solar cells, it was found that the recombination rate determined experimentally is reduced compared to the Langevin's theory. 18Taking this into account, a reduction factor ς is often included, and the bimolecular recombination rate reads where γ = q µ n + µ p /ε is the Langevin recombination constant.
Regarding the polaron pair dissociation in organic semiconductor blends, the exact mechanism is still uncertain.The Onsager-Braun theory 19,20 is a commonly used quantitative model for considering the polaron pair dissociation at the interface (see Refs. 4, 5, and 21 for instance).In this study, we use the Onsager-Braun model where the dissociation probability of polaron pairs P is defined as where k d is the polaron pair dissociation rate, a is the polaron pair separation, and k f is the polaron pair decay rate.Besides the Onsager-Braun model, the role of charge delocalization in the dissociation of polaron pairs into free charge carriers has also been studied and proposed. 22n order to numerically solve the continuity and the Poisson's equations, boundary conditions at the contacts between the electrodes and the active layer are required.For the electrical potential ψ (where ψ/q = −  Fdx), the boundary condition is where V bi is the built-in voltage and V a is the applied voltage.Here, x = 0 and x = L are the contacts at the anode and the cathode, respectively.The built-in voltage V bi is defined as where χ a and χ c are the work functions of the electrodes at anode and at cathode, respectively.The carrier densities at the contacts (x = 0, L) can be simply given by the Boltzmann statistics with where N c and N v are the effective density of states for the LUMO of the acceptor and the HOMO of the donor, respectively, E g is the effective band gap (i.e. the difference between the LUMO of the acceptor and the HOMO of the donor), φ a is the hole injection barrier at the anode (i.e. the difference between the HOMO of the donor and the work function of the electrode at anode), and φ c is the electron injection barrier at the cathode (i.e. the difference between the LUMO of the acceptor and the work function of the electrode at cathode).For simplicity, the energetic disorder in the organic active layer is not considered in this study.
In order to obtain the J-V characteristic, the continuity and the Poison's equations are numerically solved using the Scharfetter-Gummel discretization 23 and an iteration approach based on the work of Gummel. 24

A. Analytical solution
In order to solve the drift-diffusion model analytically, several approximations and assumptions are required.The electric field F inside the active layer is assumed to be uniform.This should be valid since drift-diffusion simulations show that F is almost constant inside the active layer. 4The uniform electric field F is given by (see Ref. 7 for example) where L is the thickness of the active layer.The polaron pair dissociation probability P is also assumed to be uniform inside the device.For the Onsager-Braun model, a uniform F together with constant charge mobilities will automatically lead to a uniform P. Furthermore, we assume that the generation rate of polaron pairs per unit volume G is constant and uniform inside the device.
For organic BHJ solar cells, it has been shown that a constant generation profile gives identical electrical characteristics compared to the generation profile obtained from optical model as long as the thickness of the active layer is less than 250 nm. 25 Therefore, this assumption should be valid when the thickness of the active layer is less than 250 nm.The bimolecular recombination rate R L [Eq. ( 5)] is dependent on the position inside the active layer x because the product of n and p is different at different x.In order to obtain an analytical solution, let us also assume that the bimolecular recombination rate inside the device is uniform (denoted by R C ) and equal to the average R L .Therefore, R C is independent of x and its value is given by The recombination current density J R is given by where R is any recombination rate in general.From Eq. ( 12), for R = R C (where R C is independent of x) and integrating over x from 0 to L, we get Equation ( 13) gives the total recombination current density in the device for R = R C .Applying Eq. ( 11) to Eq. ( 13), we get Therefore, the total recombination current density using R = R C is the same as the one using R = R L .We can argue that for the current-voltage (J-V) characteristic of a solar cell, it is not important where inside the active layer the charge carriers are generated and lost, as long as the net charge carriers that are collected by the electrodes are the same.For example, by integrating Eq. (2a) over x from 0 to L (with P and G are independent of x), it can be seen that the net flow of the electron current density inside the device [i.e.J n (0) − J n (L)] depends on the total recombination current density.This means that we can replace R L with R C in Eq. (2a) and still be able to get the same net electron current density.Therefore, the basic idea of assuming uniform bimolecular FIG. 1. Ratio of the bimolecular recombination current density J RL to the polaron pair dissociation probability P as a function of applied voltage V a obtained using drift-diffusion numerical simulation.The values of the parameters used are shown in Table I.
recombination rate R C is to obtain an analytical solution that should be able to describe the J-V characteristics of organic BHJ solar cells.However, this assumption may not lead to very accurate descriptions of carrier densities, electron current density, and hole current density as functions the position inside the active layer x.
Figure 1 shows the ratio of the recombination current density for R = R L (denoted as J RL ) to the polaron pair dissociation probability P as a function of applied voltage V a obtained using drift-diffusion numerical simulation.In this study, all numerical simulations are performed using a commercial program called SETFOS, 26 which is specifically made for simulating organic BHJ solar cells based on the drift-diffusion model as described in Sec.II.The values of the parameters used in the simulation are shown in Table I.The simulation produces a short circuit current density J sc of -102.6 A/m 2 , an open circuit voltage V oc of 0.59 V, and a fill factor FF of 0.578, which are typically measured values for a P3HT:PCBM solar cell.
The behavior seen in Fig. 1 can be easily understood.As the applied voltage V a is increased, the electric field decreases.This makes it harder for the electrodes to extract the generated free charge carriers (given by the product of P, G, and L).As a result, this increases the carrier densities  per generated free charge carriers, thus increases the recombination per generated charge carriers.Therefore, it is expected that the ratio of the bimolecular recombination to the generated charge carriers (or to the dissociation probability P, since G and L are constants) increases with V a .As seen in Fig. 1, we find that the natural logarithmic of J RL /P is somewhat nonlinear as a function of applied voltage V a .Varying the values of various parameters, we find that the general behavior seen in Fig. 1 is the same at least between zero applied voltage and the open circuit voltage.From this, let us approximate and write ln (J RL /P) = ln (J 0 ) + gV s a , where J 0 , g, and s are (fitting) constants.Equation ( 15) can be rewritten as If we assume that the bimolecular recombination rate is uniform inside the device as given by Eq. ( 11), the total recombination current densities using R C and R L must be equal (i.e.J RC = J RL ), as shown before.Since Eq. ( 13) equal to Eq. ( 16), we get where R 0 = J 0 /qL.Inserting Eqs.(3a) and (3b) into Eqs.(2a) and (2b), and replacing R L [Eq. ( 5)] with R C [Eq. ( 17)], the continuity equations become Equations ( 18a) and (18b) can be solved using where A n , B n , A p , and B p are independent of x.Applying the boundary conditions at the contacts [i.e.Eqs.(9a) to (9d)] to Eqs. (19a) and (19b), we can solve for A n , B n , A p , and B p , which are given by It should be noticed that the values of A n , B n , A p , and B p are different at different applied voltages V a .Inserting Eqs.(19a) and (19b), and their differentiations with respect to x into Eqs.(3a) and (3b), we can get the expressions for J n and J p , where The total current density J is the sum of J n and J p , which is given by It can be seen that the total current density J is independent of x, which must be the case for a steady current density.Using Eq. ( 10), we can rewrite the term R C /F in Eq. ( 22) as Applying Eq. ( 11), this equal to Hence, it can be seen that although we assume uniform bimolecular recombination rate R C to obtain Eq. ( 22), it should still give a J-V characteristic as if we consider the actual bimolecular recombination rate R L .It is worth noting that Eq. ( 22) can be derived using any other form of bimolecular recombination rate [i.e.not limited to the form of Eq. ( 17)] as long as the recombination rate used is independent of x and equal to the average of the actual bimolecular recombination rates R L .Similarly, any other form of polaron pair dissociation (i.e.not limited to the Onsager-Braun model) can also be used as long as it is independent of x.

IV. DISCUSSION
In order to test this analytical expression, let us compare Eq. ( 22) with the J-V curve obtained using drift-diffusion simulation, which is widely used to analyze the experimental J-V characteristics of organic BHJ solar cells.Figure 2 compares the J-V characteristic obtained using drift-diffusion simulation with the one obtained using the analytical expression.The parameters used for the simulation and the analytical expression are given in Table I.It can be seen that the J-V characteristic obtained using the analytical expression is in good agreement with the J-V characteristic obtained using drift-diffusion simulation.In order to obtain the constants R 0 , g, and s [see Eq. ( 17)], let us rewrite Eq. ( 22) as  19b)] (symbols), and using drift-diffusion numerical simulation (lines) for an organic BHJ solar cell with properties as shown in Table I.

FIG. 2. Current-voltage characteristics of an organic BHJ solar cell with properties as shown in
where the values of m and c at any given V a can be determined from the values of the parameters in Table I.By inserting the simulated current density at short circuit into Eq.( 23), we can easily solve for R 0 .To obtain g and s, we can insert the simulated current densities at two other applied voltages to obtain two equations that enable us to solve for g and s.However, this requires numerical techniques where numerical solvers such as MATLAB can be helpful.Alternatively, we can insert the open circuit voltage and current density (i.e.J = 0 A/m 2 ) into Eq.( 23), and guess the value of s in order to get the value of g.Then, the obtained g and s, together with R 0 determined earlier are used to obtain the J-V curve.This process of guessing the values of s and g is repeated until a good fit to the simulated J-V curve is obtained.The dots in Fig. 2 are obtained with R 0 = 3.69 × 10 29 m −3 s −1 , g = 17.65V −1 , and s = 5. Figure 3 compares the carrier densities as functions of x obtained using drift-diffusion simulation with the ones obtained using Eqs.(19a) and (19b) with the values of R 0 , g, and s as mentioned above.As expected, the carrier densities obtained using the analytical approach are only in moderate agreement with the ones obtained using drift-diffusion simulation.As mentioned in Sec.III, the assumption of uniform bimolecular recombination rate may lead to somewhat inaccurate descriptions of carrier densities, electron current density, hole current density, and others as functions of x.
There are a few advantages of the analytical expression over numerical simulation of the drift-diffusion model.For example, it is practically impossible for the simulation to exactly reproduce the experimental J-V curve of an organic BHJ solar cell.Moreover, the quality of fittings to experimental J-V curves obtained using simulations may not be very good for some solar cells.This means that the recombination current densities obtained from simulations may differ quite significantly from the actual recombination current densities of the studied solar cells.Using the analytical expression, it is possible to directly extract the recombination current densities from experimental J-V characteristics.This can be done by inserting the experimental current densities and the associated applied voltages into Eq.( 22) [or Eq. ( 23)], and then directly determine the value of R C at each applied voltage.Then, the recombination current density can be obtained by multiplying the extracted R C with the thickness of the solar cell and the elementary charge.Therefore, this analytical expression enables us to directly extract the overall recombination loss of an organic BHJ solar cell from the experimental J-V curve.It is worth noting that series and shunt resistances in solar cells should be minimized and maximized, respectively.If the effects of these resistances are significant, then the experimental J-V curve should be adjusted accordingly before the analytical expression can be used to analyze it.Furthermore, we also find that our numerical program fails to produce the desired output when simulating an organic BHJ solar cell with a relatively large band gap and a high injection barrier.A large band gap organic BHJ solar cell has a useful application as a sub-cell of a tandem solar cell.Hadipour et al. 27 showed that the J-V characteristic of a tandem solar cell can be constructed from the J-V characteristics of individual sub-cells.Hence, it is important to analyze the J-V characteristic of each sub-cell.The values of the parameters used in the simulations are shown in Table II.The simulation with φ a = 0 eV produces J sc = −28.92A/m 2 , V oc = 1.06 V, and FF = 0.72.However, for a tandem solar cell with series configuration, it is unlikely that both φ a and φ c are zero since the anode (cathode) of a large band gap sub-cell is connected to the cathode (anode) of a smaller band gap sub-cell, and the effective band gaps of the sub-cells almost certainly overlap with each other.When φ a = 0.9 eV is used in the drift-diffusion simulation, we get J sc = −26.08A/m 2 , V oc = 0.367 V, and FF = 0.513.However, when φ a = 1.0 eV is used, the drift-diffusion simulation starts to become unstable, producing V oc = 0.267V but J sc = 0 A/m 2 .When φ a = 1.1 eV, we get J sc = 0 A/m 2 and V oc = 0 V.When φ a = 1.2 eV is used, the numerical simulation simply fails to produce any result.Since a commercial program is used for the simulations, we are only able to vary the numerical mesh size and the residual error.Increasing those iteration parameters decrease the accuracy, but increase the possibility of getting the desired output in general.The exact reason for the failure is uncertain since we are unable to get the desired output even if those two iteration properties are increased and decreased.It is possible that when certain values of certain parameters are used, the drift-diffusion simulation would fail to fully converge (which is unsurprising for a numerical approach).However, by using the analytical expression, we should be able to, for example, extract and analyze the recombination current density as a function of applied voltage if the experimental J-V curve of the large band gap solar cell with properties as shown Table II is available.
In order to clearly demonstrate the advantage of this analytical expression, let us try to estimate J-V curves of the solar cell with properties as shown in Table II when φ a = 1.0eV and φ a = 1.1eV.Since we do not have the experimental J-V data, let us use the simulated J-V curve of the solar cell when φ a = 0.9 eV in order to estimate the R C [refer Eq. ( 17)] when φ a = 1.0eV and φ a = 1.1 eV.In order to satisfactorily reproduce the simulated J-V curve when φ a = 0.9 eV by using the analytical expression, we set V oc = 0.35V instead of V oc = 0.367V.Therefore, it is worth noting that in order to satisfactorily reproduce a given J-V curve by using the analytical expression with R C as given by Eq. ( 17), we may need to set the V oc and/or J sc to be slightly different compared to the V oc and J sc of the given J-V curve.We find that R 0 = 2 × 10 28 m −3 s −1 , g = 3.423V −1 , and s = 1 can satisfactorily reproduce the simulated J-V curve for the solar cell when φ a = 0.9eV.From this, let us write R C [refer Eq. ( 17)] for the solar cell when φ a = 1.0eV and φ a = 1.1eV as Since R 0 , g, and s are obtained from the simulated J-V curve when φ a = 0.9eV, V diff is introduced to account for the difference between the built-in voltage when φ a = 0.9eV and the built-in voltages when φ a = 1.0eV and φ a = 1.1eV.Here, we simply take V diff = 0.1V when φ a = 1.0eV and V diff = 0.2V when φ a = 1.1eV.Figure 4 shows the J-V curves estimated using the analytical expression with R C given by Eq. ( 24) when φ a = 1.0eV and φ a = 1.1eV.When φ a = 1.0eV, we get J sc = −24.42A/m 2 , V oc = 0.245V, and FF = 0.53.When φ a = 1.1eV, we get J sc = −22.12A/m 2 , V oc = 0.145 V, and FF = 0.482.Therefore, we demonstrate that when the numerical simulation fails  II and with R C given by Eq. ( 24) (see text for details).The solid line shows when φ a = 1.0eV and the dashed line shows when φ a = 1.1eV.
to produce the J-V curve of a specific organic BHJ solar cell, the analytical expression is useful in order to estimate the J-V curve and analyze the organic BHJ solar cell.

V. CONCLUSIONS
In conclusion, we have proposed an analytical approach to describe the current-voltage (J-V) characteristics of organic BHJ solar cells.The derivation is obtained by analytically solving the drift-diffusion model for organic BHJ solar cells with the assumption of uniform bimolecular recombination rate that equal to the average of the actual bimolecular recombination rates.We show that when the uniform bimolecular recombination rate is assumed, the derivation should still produce an expression for the current as a function of applied voltage as if the actual bimolecular recombination rate is used.Compared to previous analytical expressions, this analytical expression should be more suitable for analyzing the J-V characteristics of organic BHJ solar cells when the bimolecular recombination is not negligible.There are a few advantages of this analytical expression compared to numerical simulation of the drift-diffusion model.As discussed in Sec.IV, applying the analytical expression to experimental J-V data enable us to directly extract and analyze the overall recombination loss of an organic BHJ solar cell.Furthermore, as shown in Sec.IV, we find that when the injection barrier of a large band gap organic BHJ solar cell is high, our numerical program fails to give the desired output.We then demonstrated that this analytical approach can be used to analyze the large band gap organic BHJ solar cell.Therefore, when the numerical simulation fails to produce the desired results for a specific organic BHJ solar cell, the analytical expression is useful in order to analyze the organic BHJ solar cell.

FIG. 4 .
FIG. 4. Current-voltage characteristics of an organic BHJ solar cell with properties as shown in TableIIand with R C given by Eq. (24) (see text for details).The solid line shows when φ a = 1.0eV and the dashed line shows when φ a = 1.1eV.

TABLE I .
Values of the parameters used in numerical simulation.

TABLE II .
Values of the parameters used in numerical simulation.Other values that are not shown here are the same as the values in TableI.