Effect of spatial distribution of generation rate on bulk heterojunction organic solar cell performance: A novel semi-analytical approach

Organic Electronics 46 (2017) 226e241

Contents lists available at ScienceDirect

Organic Electronics journal homepage: www.elsevier.com/locate/orgel

Effect of spatial distribution of generation rate on bulk heterojunction organic solar cell performance: A novel semi-analytical approach Mahnaz Islam a, Sumaiya Wahid a, Mokter Mahmud Chowdhury b, Faysal Hakim a, Md. Kawsar Alam a, * a b

Department of Electrical and Electronic Engineering, Bangladesh University of Engineering and Technology, Dhaka, 1205, Bangladesh Department of Electrical and Computer Engineering, University of British Columbia, Vancouver, BC, V6T1Z4, Canada

a r t i c l e i n f o

a b s t r a c t

Article history: Received 7 January 2017 Received in revised form 16 March 2017 Accepted 18 April 2017 Available online 19 April 2017

We present a physics-based semi-analytical model of bulk heterojunction (BHJ) organic solar cell (OSC) for predicting the electrical characteristics of the device, taking into account the space dependency of generation rate profiles. The model enables us to derive the J-V characteristics of BHJ OSC without the need of a closed form expression of arbitrary carrier generation rate (which may not exist), hence avoiding the cumbersome numerical fitting method employed in literature previously. Using the proposed model, we perform an extensive analysis to study the effect of spatial distribution of generation rate profiles on the device performance. For this purpose, we use Gaussian shaped profiles that have a common average value thus retaining the total number of generated carriers same. We vary the position of the generation peak and its sharpness (width of the Gaussian peak) as well as number of peaks to analyze their effects on device efficiency. For the considered profiles, results show that the optimized profile has a peak carrier generation rate exactly halfway through the active layer and falls off sharply on either side. In the end, we propose methods of controlling the generation profiles by modifying the device structure and perform optical simulation to show the corresponding generation profiles. Thus, we show the usefulness of our derived model in finding the spatial distribution of a given number of carriers along the active layer that yields the best device performance. © 2017 Elsevier B.V. All rights reserved.

Keywords: Organic solar cell Bulk heterojunction Current-voltage characteristics Spatial carrier generation Semi-analytical model

1. Introduction Over the past two decades, increased research interest has been observed in the field of organic photovoltaics owing to features such as inexpensive mass production techniques, thin film structures, and room temperature processing leading to the production of low-cost, flexible and light weight solar cells [1e4]. The initial hindrance of low power conversion efficiency (PCE) of organic solar cells (OSC) (due to low dielectric constants [5]) has been overcome with the introduction of bulk heterojunction (BHJ) configuration [6e9]. BHJ structures are fabricated by blending two organic materials in solution, one electron donating (usually conjugated polymer) and the other electron accepting (fullerene or fullerene derivatives), followed by spin coating them as a single layer [1,2]. The increased donor-acceptor (D-A) interfacial area throughout the

* Corresponding author. E-mail addresses: [email protected], (Md.K. Alam). http://dx.doi.org/10.1016/j.orgel.2017.04.021 1566-1199/© 2017 Elsevier B.V. All rights reserved.

[email protected]

active layer leads to efficient dissociation of the photo generated excitons, since the length scale of the two phases is comparable to the exciton diffusion length (~10 nm), thus enhancing the PCE [6,7,9]. The reported efficiency for thin film organic solar cells has approached 11.5% in 2015 [10e12]. Also, 10.6% efficiency has been achieved for tandem polymer BHJ solar cells [13]. However, Park et al. [8] have reported the internal quantum efficiency of BHJ composites reaching 100%, implying the potential of polymer BHJ solar cells that remains unexploited. In relation to optical engineering of OSCs, numerous past attempts have been made in efforts to maximize the light harvesting capability of the active layer in BHJ structures. While increasing the thickness of the active layer appears to be a natural solution in order to enhance the light absorption in the active layer, it may deteriorate the efficiency due to increased electrical losses in the active layer, owing to enhanced charge carrier recombination and poor charge extraction resulting from the low carrier mobility and small exciton diffusion length in typical polymer blends [1,14e16]. Thus, optical engineering techniques such as light trapping

M. Islam et al. / Organic Electronics 46 (2017) 226e241

mechanisms [17,18], use of optical spacers such as ZnO2 [19,20] or TiO2 [15,21,22], light manipulation using micro/nanostructures to invoke surface plasmonic resonance effect [23e27] etc. are employed to enhance the photon absorption in the active layer and improve the light harvesting efficiencies in OSCs. The physics behind the role of these techniques in enhancing the device performance is still a matter of further research. However one of the causes most commonly cited is the spatial redistribution of the optical electric field energy dissipation (measured by jEj2, where E is the optical electric field) that causes the maximum light intensity to fall within the active layer and thus increase the photon absorption [15,16,19,22,28]. The free carrier generation rate is a direct consequence of the photo generation of excitons, which in turn is proportional to the optical electric field energy dissipation (jEj2 ). Thus, the carrier generation rate has an irrefutable effect on the overall device performance and consequently is an important parameter to study and exploit for maximizing the device properties. Despite the integral effect of the carrier generation rate on the device performance parameters, the number of theoretical models eligible to study the effect of a change in the position distribution of the carrier generation rate available in literature fall short. While numerical modeling of BHJ solar cells have been immensely studied upon [1,9,29,30], simple and usable analytical models correctly representing the operation of such devices are still scarce in literature [14,31e36]. Their limitations include ignoring the percolated paths in BHJ structures [31], the drift and diffusion transport of charges and optical propagation [32] and/or neglecting the spatial dependency of the carrier generation rate and using an average value that remains constant throughout the active layer [14]. To assist more correct and realistic studies of BHJ OSCs' electrical characteristics, in this paper, we present a semi-analytical model that reflects the drift-diffusion physics in BHJ structures and is capable of incorporating a position dependent arbitrary carrier generation rate. In a previous work, we have analyzed the importance of considering the spatial dependency of generation rate in the active layer and reported that its negligence may result in up to 10% error in calculated device efficiency of BHJ solar cells for some particular cases [33]. Following our initial discussion on the optical manipulation of light to increase the carrier generation rate in the active layer, it is therefore of utmost importance that one considers the distribution of the generation rate along the active layer to conduct a theoretical study of the effect of employing light harvesting techniques on the device performance parameters. In doing so, previous models such as that proposed by Arnab and Kabir [34] have used an exponentially decreasing generation rate (BeerLambert's law) in the direction of propagation of light, thus neglecting effects such as reflection and interference. However, since the thickness of active layer in thin film OSCs is usually less than the wavelength of incident light, these optical phenomena are significant and cannot be ignored in evaluation of the generation rate [28]. Accordingly, previous works [33,36] have used the transfer matrix formalism (TMF) to calculate the optical absorption profile. However, closed form expression obtained from this method involves several assumptions, which is not satisfied for all cases thus invalidating the optical model [6]. Furthermore, in cases where a closed form expression of the generation rate is unavailable, it may be necessary to incorporate an arbitrary spacedependent generation profile found from complex numerical simulations such as finite difference time domain (FDTD) or from experimental data in the electrical charge transport analytical model. In this regard, the only work that is available in literature is our previously proposed empirical formula to model the given generation rate using a sum-series of sine functions [35]. The lengthy procedure of numerical fitting in order to determine the

227

empirical series coefficients of the generation rate profile used in this method is inevitable to apply the analytical formulae presented. This approach is quite laborious and involves numerical methods for finding a suitable series approximation of the generation rate. Also the number of coefficients required to fit a generation rate will vary with the profile, and may fall short of correctly representing it. In this work, we have eliminated this numerical complexity and limitation using a novel method and derived an improved semi-analytical expression for J-V characteristics that can incorporate any arbitrary shaped generation profile. In addition to that, we use this proposed model to conduct a study of the effect of a change in the distribution of the carrier generation rate on the device efficiency. As mentioned previously, employing light harvesting techniques in BHJ structures results in enhanced carrier generation in the active layer, depicted by the various experimental works in this field [15,19,22]. Following our previous discussion, it can be inferred that the carrier generation rate and its spatial distribution in the active layer significantly affect the device efficiency [33]. However, a comprehensive theoretical study of this effect and its extent of altering the device performance has not yet been reported in literature. Additionally, there is a lack of an appropriate, simple and usable physics-based analytical model that incorporates the spatial distribution of any random generation rate distribution used, for example, to test proposed theories involving new and advanced device architectures intended to optimize the carrier generation rate. To overcome this hindrance, at first we derive a semi-analytical expression to model the electrical characteristics of BHJ OSCs, capable of incorporating the spatial dependency of any generation rate. With the help of this model, we then study the effect of a change in the position of the peak generation rate along the active layer on the device efficiency. In this work, we consider a Gaussian type generation profile and analyze the dependency of BHJ characteristics on spatial distribution of the generation rate. Thus, we find the position of the peak carrier generation rate and the number of peaks that yield the best device performance and show that an optimum generation rate profile maybe devised using our proposed model. In the end, we also show optical simulation results along with electrical characteristics and suggest possible methods of controlling the spatial distribution of generation profiles in OSCs using different light management techniques. 2. Model description and theory In this section, we present our model starting with an explanation of our derived technique in order to incorporate the spatial distribution of the available generation rate followed by the resulting analytical equations. We have studied the BHJ configuration used by Monestier et al. [30] in Sections 3.1 and 3.2 [Fig. 1(a)]. A similar structure has been used extensively in literature for numerical and analytical studies [14,30e33,35]. The active layer consists of a blend of poly (3-hexylthiophene) and [6,6]-phenylC61-butyric acid methyl ester (P3HT:PCBM) in a uniform weight ratio of 1:1, which has been modeled as a single meta-material consisting of the highest occupied molecular orbital (HOMO) of the donor (P3HT) and the lowest unoccupied molecular orbital (LUMO) of the acceptor (PCBM), functioning as the conduction and valence band, respectively [9,30,31]. The constituents of the structure include ITO (180 nm), PEDOT:PSS (50 nm), P3HT:PCBM (varying thickness), Al (100 nm) where poly(3,4ethylenedioxythiophene):poly(styrene sulfonate) (PEDOT:PSS) has been used as the injecting layer, and Indium tin oxide (ITO) and aluminum (Al) have been used as anode and cathode, respectively. The contacts have been assumed to be perfectly ohmic [9,32] implying that the boundary conditions

228

M. Islam et al. / Organic Electronics 46 (2017) 226e241

Fig. 1. (a) Configuration of the BHJ OSC studied in this paper. (b) Schematic energy band diagram of the structure before the materials are brought into contact.

reduce to thermal concentration [37]. The energy levels of the materials in the device before they have been brought into contact are shown in Fig. 1(b) and have been collected from literature [1,6,7]. In Section 3.3, we have modified the above structure by using a graded blend layer instead of an uniform one (by either varying the blend ratio or annealing time) (Section 3.3.1) or by inserting an optical spacer in between the active layer and the metal electrode (Al) and varying the lengths of the transport layers (TL) (Section 3.3.2). In order to derive the J-V characteristics of BHJ OSCs, the driftdiffusion formalism has been used to model the charge carrier transport in the semiconductor device [6,9,14,30] resulting in the following current density equations:

Jp ¼ qmp pE  qDp

dp ; dx

(1)

dn ; dx

(2)

Jn ¼ qmn nE þ qDn

where mpðnÞ and DpðnÞ are the mobility and diffusion coefficient of D holes (electrons), that obey the Einstein's relation, m pðnÞ ¼ kT q ¼ VT , k pðnÞ being the Boltzmann constant, T is the temperature in Kelvin and VT is the thermal voltage [6,9]. Further, pðnÞ is the concentration of holes (electrons) and JpðnÞ is the hole (electron) current density, all defined in the active layer. The charge carriers have been assumed to be transported along only one spatial dimension (x, distance along the depth of the active layer defined from the PEDOT:PSS/ active layer interface) since the ratio of the thickness to the lateral dimension (of the planar device) is very small [9,14,30,33]. The charge carriers in the active layer are drifted by the electric field, E which can be assumed to be uniform in OSCs that have thin active layers (up to 250 nm) [14] with similar electron and hole mobility values. Since in our analysis the largest active layer thickness used is 200 nm, the uniform electric field approximation has been used where E is given by Eq. (3).

E¼

Va  Vbi ; d

(3)

Here Va is the terminal voltage (conventionally termed as applied bias in literature) and Vbi is the built-in potential of the device, given by the difference in the electrodes' work functions, as shown in Fig. 1(b) [14,30,31]. Similar assumption has also been made in previously published works [14,30e33,35]. In the next section, we have analyzed the effect of this approximation on the

PCE. The current density equations are solved self consistently with the current continuity equations [14,34,35], to explain the electrical behavior of OSCs. Under steady state, the current continuity equations maybe written as:

1 dJp ¼ G  Rp ; q dx

(4)

1 dJn ¼ G  Rn ;  q dx

(5)

where G and R denote generation rate and recombination rate, respectively. In this analysis, we have neglected the recombination of charge carriers in the active layer (the terms Rp and Rn in the continuity equations) which is the case for OSCs where the carrier mobility values do not differ by more than a factor of ~50 [14], thus maintained in our analysis by using mn =mp ¼ 15. The negligible recombination effect has been illustrated in details by Altazin et al. [14] and employed in many previous analytical [14,33,35], as well as numerical/experimental works, such as Monestier et al. [30] who explained their experimental results by ignoring the recombination loss. However, it is to be noted that when the mobility values are too different, the carriers have very different transport properties resulting in charge pile-up in the active layer, leading to the eventual formation of space charge. In the latter case, the recombination of charge carriers is more pronounced and must be considered in the continuity equations [14]. Also, the uniform electric field approximation is then no longer valid [9,14], which necessitates solving the Poisson equation self-consistently in addition to the current density and current continuity equations to determine the electric field profile in the active layer. Using Eqs. (1) and (2) in Eqs. (4) and (5), we arrive at the following:

mp VT

d2 p dp ¼ GðxÞ þ mp E dx dx2

(6)

mn VT

d2 n dn ¼ GðxÞ  mn E dx dx2

(7)

In order to solve the second order differential Eqs. (6) and (7), we have discretized the given generation rate profile by assuming it to be constant in each region i of the active layer over a very small distance Dx, where i varies from 1 (at x ¼ 0 nm) to l ¼ Ddx þ 1(at x ¼ d nm) and l is the total number of regions, as demonstrated in Fig. 2. This results in a series of ordinary

M. Islam et al. / Organic Electronics 46 (2017) 226e241

229

and Gm is the value of the discretized generation rate in region m. Following through, we define a summation relation for the coefficients in region ðm þ rÞ, m  r  l  m, given by:

Amþr ¼ Am 

! VT iDx ðGi  Giþ1 Þ  mp E2 mp E

(14)


X VT mþr1 EðiDxÞ ðG  G Þexp  i iþ1 VT mp E2 i¼m

(15)

mþr1 X i¼m

A0mþr ¼ A0m þ

Equations (14) and (15) imply that if the value of the coefficients in one region is known, then all other coefficients maybe found out. If m ¼ 1 and m þ r ¼ l and we define quantities D and D0 equal to the summation terms in Eqs. (14) and (15), respectively defined over i ¼ 1 to i ¼ l  1 then:

Fig. 2. Discrete approximation of arbitrary generation rate profile for demonstration purposes only.

differential equations, in each region i that are solved using appropriate boundary conditions to find the carrier concentration profiles. As we will show in the next section, the least number of regions required to yield an acceptable percentage of error is small enough to compute the results by hand. Naturally, a better approximation may be obtained by using a very small value of Dx with the help of a simple computer program; thus the model may be deemed “semi-analytical”. In each region i, we can solve Eqs. (6) and (7) to find the electron and hole concentrations in that region as a function of x and Gi , given by:

Ad ¼ A1  D

(16)

A0d ¼ A01 þ D0

(17)

Using the standard boundary conditions at the metal contacts in
 Eq. (8) for the hole concentration [38] p1 ð0Þ ¼ Nv exp  fV3T and
 pl ðdÞ ¼ Nv exp  VfT4 , where f3 and f4 are as defined in Fig. 1(b) and Nv is effective density of states for holes, and using Eqs. (16) and (17) we find the value of the coefficients for hole concentration in region 1 to be:


f A1 ¼ Nv exp  3  A01 VT

(18)

(19)

(20)


 Ex Gx þ i ; pi ðxÞ ¼ Ai þ A0i exp VT mp E

(8)

  L  Nv exp  Vf3T   A01 ¼ exp Ed VT  1


0 Ex Gx ni ðxÞ ¼ Bi þ Bi exp   i ; VT mn E

(9)



Ed f Gd þ Nv exp  4  l where L ¼ D  D0 exp VT VT mp E

where Ai ; A0i ; Bi and B0i are the coefficients of the i-th region needed to evaluate electron and hole concentrations. In order to find them, we have considered the boundary conditions, namely the continuity of the concentration profile and its derivative. Then, for the carrier concentration, c (where c ¼ n or p) at the interface of the mth and ðm þ 1Þth regions, at a distance mDx from the origin (of the discretized profile) we have:

cm ðmDxÞ ¼ cmþ1 ðmDxÞ

(10)

   dcm  dc ¼ mþ1   dx mDx dx mDx

(11)

Solving Eqs. (10) and (11) for hole type carrier, we arrive at the following recurrence relations:

!

Am  Amþ1 ¼

A0mþ1



A0m

VT mDx  ðGm  Gmþ1 Þ; mp E2 mp E


V EðmDxÞ ðGm  Gmþ1 Þ; ¼ T 2 exp  VT mp E

(12)

(13)

where the subscripts are such that Am and A0m are the coefficients

Having found the value of the coefficients in region 1, we are now able to find the coefficients of every other region using the recurrence relation in Eqs. (12) and (13), thus determining the hole concentration profile in the active region as depicted by Eq. (8). Using the same approach for the electron concentration and stan
 dard boundary conditions, n1 ð0Þ ¼ Nc exp  VfT2 and
 nl ðdÞ ¼ Nc exp  VfT1 where f1 and f2 are as defined in Fig. 1(b) and Nc is the effective density of states for electrons, we may find the coefficients in Eq. (9) in every discrete region. The corresponding derived expressions have been provided in Table 1. To find the expressions for Jn and Jp in each region i, we differentiate Eqs. (8) and (9) and replace them in Eqs. (1) and (2) to derive:


   Ex qV G  T i Jpi ðxÞ ¼ qmp pi E  qmp A0i E exp VT E

(21)


  Ex qV G Jni ðxÞ ¼ qmn ni E  qmn B0i E exp   T i VT E

(22)

The total current density is then the summation of the hole and electron current density and remains constant throughout the device. Thus, in any region i,

230

M. Islam et al. / Organic Electronics 46 (2017) 226e241

Table 1 Derived expressions to find the coefficients required to determine the electron concentration in Eq. (9). Eq. No. for hole type carrier

Corresponding equation for electron type carrier



 Gm  Gmþ1

 DxÞ Gm  Gmþ1 B0mþ1  B0m ¼ mVET 2 exp Eðm VT n

 P VT iDx Bmþr ¼ Bm  mþr1 Gi  Giþ1 i¼m mn E2 þ mn E
 P ðGi  Giþ1 Þexp EðiVDT xÞ B0mþr ¼ B0m þ mVET 2 mþr1 i¼m

(12)

VT m Dx mn E2 þ mn E

Bm  Bmþ1 ¼

(13) (14) (15)

n

Bd ¼ B1  K B0d ¼ B01 þ K 0
 B1 ¼ Nc exp  fV2T  B0 1  

(16) (17) (18) (19)

f

B01 ¼

MNc exp V2



exp



VEd T

T

1

Fig. 3. J-V characteristics of the OSC for d ¼ 100 nm using our analytical model (solid line) and numerical simulation (circular marker).



  M ¼ K  K 0 exp VEdT þ Nc exp  fV1T þ mGl dE

(20)

n




Ex JT ¼ JTi ðxÞ ¼qE mn ni þ mp pi  qE mn B i exp  VT


Ex 2qVT Gi 0  þ mp Ai exp VT E 



0

(23)

where the coefficients A0i and B0i have already been derived. Thus far, we have assumed that the rate of free charge carrier generation rate is equal to the rate of exciton generation, similar to previous reports [6,30]. However, the photo generated excitons that dissociate into bound e-h pair at a D-A interface, driven by the difference in the LUMO levels of the donor and acceptor materials, is only metastable and may decay to the ground state before they can dissociate in to free charge carriers [6,9]. According to the Onsager-Braun theory [38,39], if the probability of this electric field dependent [6] dissociation of e-h pair is given by P(E) then we can include this phenomenon in our analysis by rewriting Eqs. (4) and (5) with P(E)G(x) in place of G(x). Comprehensive analysis of P(E) and its effects on device parameters may be found in our previous work [33]. As for this model, one may simply replace all terms containing Gi with P(E)Gi in the above derived expressions to account for the dissociation probability. 3. Results and discussion 3.1. Model verification using conventional generation profile To generate the J-V characteristics of the BHJ structure, as shown in Fig. 3 for d ¼ 100 nm, we have used the generation rate profiles produced by Monestier et al. [30] who have performed complex numerical simulation to calculate them from the electromagnetic field distribution. In this analysis, we have considered P(E) ¼ 1, i.e. the rate of exciton generation is equal to the free charge carrier generation rate, without any loss through decay [30,33,35,40]. The parameters used in this analysis are [32,41]: ε (dielectric constant) ¼ 3  1011 F/m, Eg (effective bandgap) ¼ 1.0 eV, Nc ; Nv (effective density of states) ¼ 1  1026 m3, mn (electron mobility) ¼ 3  107 m2/V-s, mp (hole mobility) ¼ 2  108 m2/V-s, and T (temperature) ¼ 300 K. The results obtained using our model are then compared with numerical solution computed using the Finite Element Method (FEM) in Comsol Multiphysics [42]. Agreement between the two verifies the mathematical accuracy of our model, as depicted in Fig. 3, establishing our technique used to

derive the semi-analytical equations to be correct and viable. Furthermore, we have also produced the J-V characteristics for different generation rate profiles (data also taken from Monestier et al. [30]) by varying the active layer thickness (60 nm, 140 nm and 200 nm) which have shown the same matching of results with that generated numerically. The number of regions that need to be considered to correctly model the space dependency of G(x) using a staircase approximation depends on the particular profile. For example, for the generation profile used in our analysis for d ¼ 100 nm (corresponds to the J-V curve shown in Fig. 3) using Dx ¼ 11 nm ðl ¼ 10Þ results in an error of ~2.9% in the solar efficiency. Increasing the number of regions increases the accuracy, as expected. For very fine discretization step, a simple computer program may be used to evaluate the derived equations. As for our analysis, using Dx ¼ 0:01 nm ðl ¼ 1000Þ gives an error of ~0.3%. We have used the above analysis to study the performance of OSC for different active layer thicknesses such as d ¼ 60 nm, 100 nm, 140 nm and 200 nm. The results are listed in Table 2(a). The Jsc Voc efficiency (h) of the device is given by h ¼ FF AM 1:5G  100%, where Jsc and Voc are the short circuit current density and open circuit Jsc Voc voltage, respectively and FF is the fill factor, (FF ¼ Jmp Vmp , where Jmp

and Vmp are the current density and voltage at maximum power respectively), while AM1.5G ¼ 1000 W/m2. The short-circuit current density, Jsc and device efficiency, h as a function of thickness have also been plotted in Fig. 4 as shown by the blue lines and conform to previously established findings [33,35,41]. The behavior of these

Table 2 Effect of varying active layer thickness on performance parameters of BHJ OSC. d (nm)

Jsc (A/m2)

Voc (V)

hð%Þ

FF

(a) Data obtained using our analytical model assuming electric field (E) to be uniform. 60 100 140 200

67.78 69.30 69.32 85.02

0.571 0.586 0.595 0.610

2.84 3.02 2.97 3.57

0.736 0.743 0.721 0.699

(b) Data obtained numerically considering actual electric field (E). 60 100 140 200

68.40 70.15 70.98 87.10

0.542 0.543 0.543 0.545

2.82 2.91 2.86 3.12

0.761 0.764 0.743 0.658

M. Islam et al. / Organic Electronics 46 (2017) 226e241

231

Fig. 4. Effect of varying active layer thickness, d on (a) Short circuit current density (Jsc) and (b) Efficiency (h). The blue line represents the results obtained using our model while the circular markers (cyan) depict the values obtained numerically assuming uniform electric field, the matching between the two is apparent. The square marked red line is obtained numerically by calculating actual electric field. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

parameters have also been explained in Ref. [33] assuming uniform electric field approximation. Here, the main purpose of presenting these results along with rigorous numerical simulation results (i.e. without the uniform field approximation) is to investigate the range of thickness over which our derived expressions are valid. To include the actual electric field, Poisson equation was solved self consistently with the current density equations and the current continuity equations using Comsol Multiphysics [42]. The results are given in Table 2(b) and Fig. 4 (red lines, square marker). The comparison illustrates that the device performance parameters show insignificant difference due to the assumption of uniform electric field for relatively thin active layer thicknesses, as expected. For instance, the efficiency of the device shows approximately 1% to 3% discrepancy for thicknesses between 60 nm and 140 nm. As mentioned earlier, many previous studies have used this approximation and it has been argued that the condition may be applied in BHJ structures of active layer thickness up to ~250 nm [14]. We notice that the actual electric field (Fig. 5) within the bulk of the device remains fairly constant for smaller

thicknesses. However, for thicker active layers (~200 nm), the electric field shows slight increase within the bulk so that the uniform electric field approximation may no longer correctly represent the device physics; thus resulting in unacceptable discrepancies. For larger values of thickness, although the photoabsorption is enhanced, the low carrier mobility in the active layer causes a rise in the internal resistance thereby leading to reduced device efficiency [1,14]. On the other hand, a higher carrier collection efficiency is achieved for very small values of thickness at the expense of lower light absorption capability. Thus, to optimize the device performance yielding increased photo-absorption and lower electrical losses [15,43], the commonly used active layer thickness in literature for the study of BHJ solar cell lies roughly within 70 nm to 160 nm [1,9,44e46]. Since our model produces acceptable results in this range, the validity of our model is not hampered by the uniform electric field approximation. 3.2. Application example: study using Gaussian generation profile With a view to finding an optimized generation profile, in the following analysis we use a Gaussian type generation profile in our derived model to study the effect of a change in the carrier generation rate distribution on the device performance for an active layer thickness of 100 nm. In order to compare the results obtained, we keep the total number of carriers generated constant by retaining the average value of the different generation profiles same. Further, it is to be noted that our results are not affected by the consideration of a standard Gaussian function in place of an arbitrary one. This is because any standard generation rate that has an analytical form or any random distribution that doesn't have a closed form expression is treated in the same manner by our proposed technique. The versatility of our model has been further established in Section 3.3 where we have used realistic generation profiles obtained from FDTD optical simulation and compared the resulting cell performance.

Fig. 5. Actual electric field in the blend layer for varying values of the blend thickness, d for Va ¼ 0.2 V found by numerical simulation.

3.2.1. Effect of the position of generation peak (Gpeak) At first, we analyze the optimum position of the peak carrier generation rate in the active layer. For this purpose, the Gaussian shaped generation profiles used are shown in Fig. 6, generated using Eq. (24).

232

M. Islam et al. / Organic Electronics 46 (2017) 226e241

Fig. 6. Gaussian generation profiles used to study the effect of position of peak on the device performance. All profiles have the same standard deviation (c ¼ 10 nm) and different value of position of peak (b). From the left, G1(b ¼ 5 nm), G2 (b ¼ 25 nm), G3 (b ¼ 35 nm), G4 (b ¼ 50 nm), G5 (b ¼ 65 nm), G6(b ¼ 75 nm) and G7(b ¼ 95 nm), and all profile have the same average value (Gavg). The green straight line corresponds to a constant generation profile, Gc throughout the blend whose value is equal to Gavg. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

2 ðxbÞ 2c2

GðxÞ ¼ ae

;

(24)

where a is an arbitrary constant, b is the position of the Gaussian peak and c is the standard deviation. The position of the peak (b) has been varied to appear at different distances along the active layer to produce the generation profiles. For example, the peak of generation profile 1 (G1-red line with solid circular markers) occurs at b ¼ 5 nm, while that of G4 (blue dashed line) occurs at b ¼ 50 nm in the active layer where, d ¼ 100 nm. An increase in the total carrier generation in the blend will inevitably provoke enhanced device performance. Thus, to find the sole effect of the position of the peak carrier generation rate, the total carrier generation in the active layer is maintained constant, by keeping the average generation rate in the active layer same for all the profiles in Fig. 6 (G1 to G7), using the parameter a in Eq. (24). As a point of reference, we have chosen this average value to be the same as that for the generation profile in 100 nm of active layer thickness taken from Monestier et al. [30] which was also used to produce the J-V characteristics in Fig. 3. Another property of Gaussian function is standard deviation (c) which is a measure of the sharpness of its peak. A higher standard deviation implies that carriers have been generated throughout the active layer in comparable proportion. In contrast, a low standard deviation means that most carriers have been generated near the Gaussian peak and the rest of the active region has very few or no generated carriers. All the generation profiles shown in Fig. 6 have the same standard deviation (c ¼ 10 nm). Table 3 shows the results achieved using the generation profiles as shown in Fig. 6. The different columns correspond to the different generation profiles, while the rows give Jsc, Voc, Jmp and Vmp along with FF and h. The results show that the device efficiency increases as the peak of the generation rate profile moves towards

the mid-position of the active layer (from G1 to G4) and decreases as it moves away from the mid-position. Thus, maximum efficiency (~3.42%) occurs at b ¼ 50 nm (G4), i.e. exactly halfway through the active layer for the given profiles. Further, Voc and Vmp show slight changes with variation in the peak position while Jsc and Jmp show similar trend as h and FF are maximum for G4. However, comparing the results for G1 or G7, where the peak is close to an edge of the active layer and that for G4 where the peak is exactly at the middle of the active layer, the increase in Jmp (~55.7%) is more than that in Jsc (~21.6%). As a result, fill factor follows the trend in Jmp and is also maximum (~80.05%) for the generation profile, G4. Thus, we attribute the increase in device efficiency due to a change in peak position towards the middle of the active layer due to a similar increase in Jmp. To further investigate the physics behind the increase in device efficiency as Gpeak moves towards the middle of the active layer, we study the drift (Jdr,mp) and diffusion (Jdff,mp) components of the current at maximum power point (Jmp) for each generation profile (G1 to G7, Fig. 6). To present the change in Jdr,mp and Jdff,mp with a change of peak position, we calculate the difference in each of these components relative to those for G1, DJdr;mp and DJdff ;mp . We consider the components for generation profile G1 as the reference since G1 has peak generation rate closest to the edge, resulting in the least device efficiency. Thus DJdr;mpðGnÞ ðxÞ ¼

Jdr;mpðGnÞ ðxÞJdr;mpðG1Þ ðxÞ , Jdr;mpðG1Þ ðxÞ

where n represents the generation profile number in Fig. 6 (except 1). The resulting plots are shown in Fig. 7 (a) and (b) where DJdff ;mpðGnÞ ðxÞ has been defined similarly. The plots reveal that generation profiles that yield higher efficiencies than G1 also have greater drift and diffusion components of current compared to that of G1, yielding positive differences. The difference in the drift and diffusion components for G7 is ~0% all throughout the active layer since it gives a similar efficiency. We note here that considering the difference with current components

Table 3 Results obtained using the proposed model and the generation profiles shown in Fig. 6. Position of peak in the active layer, b

5 nm (G1)

25 nm (G2)

35 nm (G3)

50 nm (G4)

65 nm (G5)

75 nm (G6)

95 nm (G7)

Constant Gc ¼ Gavg

Jsc (A/m2) Voc (V) Jmp(A/m2) Vmp (V) FF(%) hð%Þ

59.55 0.573 43.94 0.469 60.85 2.07

71.39 0.587 63.64 0.495 74.04 3.12

72.21 0.590 66.92 0.502 78.53 3.35

72.41 0.591 68.40 0.505 80.05 3.42

72.21 0.590 66.92 0.502 78.53 3.35

71.38 0.587 63.63 0.495 74.03 3.12

59.14 0.573 43.63 0.468 60.83 2.05

67.49 0.584 58.08 0.492 72.70 2.85

M. Islam et al. / Organic Electronics 46 (2017) 226e241

233

Fig. 7. Difference in (a) drift current component, DJdr;mp and (b) diffusion current component, DJdff;mp of the current at maximum power point, Jmp for each of the generation profiles shown in Fig. 6 with respect to those for G1.

under G7 as reference results in the same plots, also predicted by the symmetry seen in the values of Table 3. The two sets of plots in Fig. 7 (a) and 7 (b) are identical, implying that drift and diffusion components are affected in a similar manner and to the same extent. The drift component is proportional to the charge concentration while the diffusion component is proportional to the derivative of the charge concentration. Since the charge concentration profile shows an exponential decay (results not shown), its derivative is exponential as well. Thus, both the current components show similar behavior and are equally contributing to the increase in efficiency due to a shift of peak generation rate to the middle of the active layer. Since the efficiency achieved was highest under G4, the difference in drift and diffusion components are highest for G4 (blue dashed line). It is also seen that the difference in both drift and diffusion components themselves show spatial distribution: the difference is largest at the mid-position of the active layer and decreases sharply as we move towards either edge. The last column of Table 3 shows the results obtained from considering a constant generation rate throughout the active layer, Gc equal to the average value (Gavg) common to all the other generation profiles (G1 to G7). Gc is represented by the green line in Fig. 6 and the resulting h ¼ 2:85% and FF ¼ 72.7%, both of which are close to the average h (2.9%) and FF (72.4%) of the other values shown in the table. This illustrates the significance of considering the spatial distribution of the generation profiles, as we have done

in our proposed model, instead of approximating it as remaining constant, as employed in literature previously [14]. The latter approximation would overstate the efficiency by ~38% for G1 or G7 and understate that for G4 by ~17%. The difference in drift and diffusion components of Jmp for each generation profile (G1 to G7) shown in Fig. 6 relative to those obtained from Gc are plotted in Fig. 8. In contrast to Fig. 7, we notice DJdr;mpðGnÞ ðxÞ and DJdff ;mpðGnÞ ðxÞ for some generation profiles (G1 and G7) show negative values. The device efficiencies under G1 and G7 are less than that under Gc implying comparatively smaller drift and diffusion components, hence leading to negative differences. 3.2.2. Effect of uniformity of the spatial carrier distribution Next, to further optimize the carrier generation rate, we vary the standard deviation (c) of the Gaussian profiles while the peak generation rate is made to occur at the mid-position of the active layer for all the profiles (Fig. 9). Again, the total carrier generation in the active layer is kept constant by maintaining the same average value of all the profiles, since a higher total number of carriers will inevitably lead to higher efficiency. While the analysis so far shows that the peak carrier generation rate should occur halfway through the active layer for maximizing the device efficiency, it is still unclear how strong the peak generation rate should be relative to the generation rate in the rest of the active layer. For a Gaussian function, the standard deviation is a measure of the width of its peak. A

Fig. 8. Difference in (a) drift current component, DJdr;mp and (b) diffusion current component, DJdff;mp of the current at maximum power point, Jmp for each of the generation profiles shown in Fig. 6 with respect to those for Gc.

234

M. Islam et al. / Organic Electronics 46 (2017) 226e241

Fig. 9. Gaussian generation profiles used to study the effect of the distribution of generated carriers about the mid-position of the active layer. All the profiles have the same position of peak generation rate (b ¼ 50 nm) and different values of standard deviation (c). From the top, G10 (c ¼ 7.5 nm), G20 (c ¼ 10 nm), G30 (c ¼ 20 nm), G40 (c ¼ 30 nm) and G50 (c ¼ 40 nm) and all profiles have the same average value (Gavg).

smaller deviation implies a stronger and narrower peak, suggesting most carriers are generated close to the peak and the curve falls off sharply on either side, implying almost no carriers are generated in regions away from the peak (dead regions). A higher standard deviation translates to higher generation rates away from the peak and thus a broader peak. The generation profile G10 in Fig. 9 has the lowest standard deviation (c ¼ 7.5 nm), which is increasing onwards to G50 which has the highest standard deviation (c ¼ 40 nm) while they have the same average value equal to Gavg. One may notice that the generation profile G20 of Fig. 9 is the same as G4 in Fig. 6. The profiles show that a portion of the active region has almost no carrier generation (dead region) under G10 at each edge

with a very sharp peak at b ¼ 50 nm while G50 contributes carriers throughout the active layer. It has been suggested previously in literature that the efficiency is expected to be higher if the proportion of dead region in the active layer may be minimized [22]. However, the results achieved using our model for the generation profiles in Fig. 9, shown in Table 4, depict that the efficiency achieved under G10 is higher by ~14% compared to that under G50. This implies that a very sharp peak generation rate at the optimum mid-position, even with dead edge regions will result in a higher efficiency compared to a more spread out distribution of the same total number of carriers. The rows in Table 4 depict the same quantities as in Table 3 for the different generation profiles shown in Fig. 9. This change in the device efficiency due to a change in the standard deviation is also traced back to a change in Jmp as before and its drift and diffusion components are analyzed as shown in Fig. 10. The difference in the drift (D Jdr,mp) and diffusion (D Jdff,mp) components under each generation profile in Fig. 9 are now calculated with respect to those under G50, since it yields the lowest efficiency. The results show that the variation in these relative differences takes the same trend as the generation profiles themselves.

3.2.3. Effect of the number of generation peaks in the active layer Lastly, we investigate the effect of number of peaks in the generation profile that fall within the active layer on the device performance. Fig. 11 shows the corresponding generation profiles used for this purpose, where Profile 1 has one peak at b ¼ 50 nm, with c ¼ 10 nm, Profile 2 is the sum of two Gaussian functions which have peaks at 30 nm and 70 nm and Profile 3 is the sum of three Gaussian functions with peaks at 15 nm, 50 nm and 85 nm, each having c ¼ 9.5 nm. The three profiles under investigation have the same average value, Gavg. The resulting performance parameters are shown in Table 5. The values depict that an increase in the number of peaks deteriorate the device efficiency as the FF suffers.

Table 4 Results obtained using the proposed model and the generation profiles shown in Fig. 9 showing the effect of uniformity of carrier distribution. Standard deviation, c 2

Jsc (A/m ) Voc (V) Jmp(A/m2) Vmp (V) FF(%) hð%Þ

7.5 nm (G10 )

10 nm (G20 )

20 nm (G30 )

30 nm (G40 )

40 nm (G50 )

72.30 0.591 68.75 0.505 80.52 3.44

72.16 0.591 68.40 0.505 80.24 3.42

71.61 0.589 65.62 0.501 77.79 3.28

70.42 0.587 62.36 0.497 75.44 3.12

69.36 0.586 60.26 0.495 74.18 3.02

Fig. 10. Difference in (a) drift current component, DJdr;mp and (b) diffusion current component, DJdff;mp of the current at maximum power point, Jmp for each of the generation profiles shown in Fig. 9 with respect to those for G50.

M. Islam et al. / Organic Electronics 46 (2017) 226e241

235

3.3. Management of optical generation

Fig. 11. Gaussian generation profiles used to study the effect of the number of peaks in the generation rate that fall within the active layer. Profile 1 has one, Profile 2 has two, and Profile 3 has three peak generation rates while all of them have a common average value equal to Gavg.

Table 5 Results obtained using the proposed model and the generation profiles shown in Fig. 11. Generation profile used

Profile 1

Profile 2

Profile 3

Jsc (A/m2) Voc (V) Jmp(A/m2) Vmp (V) FF(%) hð%Þ

72.41 0.591 68.40 0.505 80.32 3.42

72.06 0.589 65.39 0.499 76.91 3.27

69.30 0.585 59.29 0.492 72.28 2.91

So far, we have considered standard Gaussian profiles to show the efficacy of our model and emphasize the role of spatial distribution of carriers in the active layer. In the next section, we suggest modifications in the device structure described in section 2 in order to control the photo-carrier generation rate and obtain an optimum profile. Furthermore, we perform optical simulation using FDTD on the modified device structure and show the generation profiles obtained for two active layer thicknesses, to ¼ 100 and 140 nm.

As suggested in literature, the maximum carrier generation in thin film solar cells may occur deeper in the active layer instead of at the edge due to various optical phenomenon such as reflection and interference [22,28]. These optical phenomena may also give rise to multiple peaks in the carrier generation rate distribution [30]. At first, we perform optical simulation (using FDTD in Lumerical Solutions) on the device structure described in section 2, which has been used by Monestier et al. [30] to calculate the electromagnetic field and the carrier generation rate profiles using numerical simulation. In absence of any experimental data on generation rate profiles in literature, we have used published simulation data to benchmark our FDTD simulator settings. It was found that the spatially distributed generation profiles obtained using our optical simulation is able to follow those produced by Monestier et al. [30]. The oscillating nature in these profiles also demonstrate the interference effect seen in thin film OSCs as mentioned previously. While researchers have performed similar FDTD simulation on BHJ structures [47,48], they didn't explicitly study the carrier generation rate or the means of controlling it. We propose plausible solutions to manipulate the generation rate profile by either (i) selectively varying the carrier generation rate in different regions of the active layer and/or (ii) shifting the peak position within the active layer. 3.3.1. Selective variation of generation rate We consider attaining control of the generation rate profile as selectively increasing/decreasing the optical absorption in respective (desired) areas of the blend. This translates to using variable absorption coefficients in different regions of the blend layer. The morphology of polymer thin films significantly affect their optoelectronic properties and the eventual absorption characteristics when used as an active layer in solar cells [49e52]. It has been reported that factors such as the donor: acceptor ratio (blend ratio) [53,54], the solvent used [55,56], processing conditions of blend mixture from solution [55,57,58] and annealing processes [54,59e61] have a large influence on the active layer morphology. In particular, for P3HT: PCBM blend, absorption characteristics and the corresponding device performance have been extensively studied for different blend ratios [62e67] and annealing conditions

Fig. 12. Generation profiles obtained by optical simulation by varying the blend ratio for (a) to ¼ 100 nm and (b) to ¼ 140 nm. In both cases, uniform blend ratios of 1:1 and 1:4 throughout the active layer as well as graded blend structures with 1:4 ratio at the edge region and 1:1 ratio in the middle region have been employed. The length of each region is shown in the legend for different cases considered in the form 1:4/1:1/1:4.

236

M. Islam et al. / Organic Electronics 46 (2017) 226e241

[48,60,62,63,68e70]. For the former case, it has been reported that the value of the extinction coefficient (k) decreases significantly in the visible range with increasing concentration of PCBM in the blend, particularly between 450 and 600 nm which is the peak absorption range for P3HT [64,65,71,72]. Since the absorption coefficient is proportional to the extinction coefficient, the absorption spectrum shows a decrease in this region as well. The value of the extinction coefficient represents the degree of crystallization in active layer, and inclusion of PCBM disrupts the crystallinity of P3HT, resulting in low absorption [71,72]. Thus, we propose using different concentration of PCBM in the active layer to selectively vary the optical absorption and hence the carrier generation rate; i.e. a graded blend ratio structure. For example, we use a lower concentration of PCBM in the middle of the active layer where the peak generation rate is required. Fig. 12 shows the generation profiles obtained from FDTD simulations of such a graded blend structure with 1:4 blend ratio at the edge region and 1:1 ratio in the

middle region of the active layer (1:4/1:1/1:4), along with uniform blend structures of ratios 1:1 and 1:4. Optical constants for 1:1 and 1:4 P3HT:PCBM blend ratio have been taken from Ref. [30] and Ref. [73], respectively. It is seen that by controlling the thickness of each region, it is possible to sharpen the peak and shift its position towards the middle of the active layer. Graded blend ratio structures have already been proposed for many polymer-based solar cells, namely, P3HT:PCBM [74], zinc phthalocyanine: fullerene (ZnPc: C60) [75e77], copper phthalocyanine: fullerene (CuPc: C60) [78,79], boron subphthalocyanine chloride:fullerene (SubPc:C60) [80,81] and poly (3-octylthiophene) [6,6]:-phenyl-C61 butyric acid methyl ester ((P3OT:PCBM) [82]. The graded composition of the active layer may be attained by annealing [70,83], thermal induction of interdiffusion between two layers [82,84,85] or by adequate solvent adjustment [84,86]. Selective absorption in the blend can also be obtained by using different annealing conditions. With increasing annealing

Fig. 13. Generation profiles obtained by optical simulation by varying the annealing duration in the blend for (a) to ¼ 100 nm and (b) to ¼ 140 nm. The uniform profiles represent full blend layers annealed for constant duration of 0 min (Uniform 1, Uniform 10 ), 30 min (Uniform 2, Uniform 20 ) and 60 min (Uniform 3, Uniform 30 ). The graded profiles represent blend layers where different regions have been annealed for different durations: Graded 1: 60(0)/40(30) nm (mins), Graded 2: 50(60)/25(30)/25(0) nm (mins), Graded 1': 20 (0)/ 30(30)/90(0) nm (mins) and Graded 2’: 50 (60)/40(30)/50(0) nm (mins). The annealing temperature is kept constant at 140  C.

Fig. 14. Optically simulated generation rate profiles for (a) to ¼ 100 nm and (b) to ¼ 140 nm. The solid curve is for uniform blend ratio structures (1:1) and the dashed curve represents profile obtained for a graded blend ratio structure 1:4/1:1/1:4 where the lengths shown in legend (in that order) have been chosen so that the graded generation profiles have the same average as the uniform one.

M. Islam et al. / Organic Electronics 46 (2017) 226e241

237

Fig. 15. Optically simulated generation rate profiles for (a) to ¼ 100 nm and (b) to ¼ 140 nm. The red short dashed curve is for blend structures annealed for constant durations of 60 min (Uniform 3, Uniform 30 ). The green long dashed curve represents profiles obtained for blend structures where different regions of varying length have been annealed for different durations such that the resulting generation profile has the same average value as its corresponding uniform one: Graded 3: 30 (0)/40(10)/30(0) nm (mins) and Graded 3': 35(10)/35(60)/35(30)/35 (0) nm (mins). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

temperature, the refractive index decreases and the extinction coefficient increases in the absorption range 450e600 nm [72,87]. The highest absorption coefficient for P3HT:PCBM is obtained at an annealing temperature of 140  C [48,60,70] at which the optimum phase separation for solar cell application [60] and highest crystallinity is attained [48]. Also, the absorption coefficient varies with the annealing duration and is found to increase with increasing

annealing time [60,63]. At 140  C, the highest absorption coefficient for 1:0.8 P3HT:PCBM blend has been reported for an annealing duration of 30 min, which reduces for an extensive annealing duration of 60 min due to disrupted crystallinity [60]. Therefore, we propose selective annealing of layered regions of the blend for different durations. Fig. 13 shows the generation rate profiles obtained for blend layers uniformly annealed for 0, 30 and 60 min,

Fig. 16. Effect of inserting an optical spacer that also behaves as an electron transport layer (ETL) of varying lengths while the hole transport layer (HTL) (Pedot:PSS) is kept at a constant length (50 nm) for to ¼ 100 nm.(a) TiO2 spacer of lengths ¼ 0, 10 and 20 nm (b) TiO2 spacer of lengths ¼ 30, 40 and 50 nm (c) TiO2 spacer of lengths ¼ 60, 70, 80, 90 and 100 nm and (d) Position of peak generation rate in the active layer (100 nm) versus thickness of spacer inserted.

238

M. Islam et al. / Organic Electronics 46 (2017) 226e241

Fig. 17. Effect of inserting an optical spacer that also behaves as an electron transport layer (ETL) of varying lengths while the hole transport layer (HTL) (Pedot:PSS) is kept at a constant length (50 nm) for to ¼ 140 nm. (a) TiO2 spacer of lengths ¼ 0, 10, 20, 30, and 40 nm (b) TiO2 spacer of lengths ¼ 50, 60, 70, 80, 90 and 100 nm (d) Position of peak generation rate in the active layer (140 nm) versus thickness of spacer inserted.

Fig. 18. Effect of varying the length of HTL in (a) to ¼ 100 nm and (b) to ¼ 140 nm for a fixed electron transport layer lengths (7 nm in (a) and 100 nm in (b)).

and for annealing times in different regions of different durations. Optical data for P3HT: PCBM (1:0.8) annealed for different durations at 140  C are taken from Ref. [60]. As suggested in graded blend ratio structure, the sharpness and position of the peak can be controlled by choosing the best absorption criteria and varying the length of the regions. In order to confirm that optically simulated generation rates show the same trends as those found in Section 3.2, we use our proposed model to calculate the device parameters in this case too. For instance, Fig. 14 shows the generation rate profiles in both uniform (1:1) and graded blend ratio structures (1:4/1:1/1:4), where the graded profile has the same average generation rate value as the uniform blend ratio profile. The use of different blend ratio in the active layer allows us to obtain an optimum-like profile, where the peak generation occurs close to the middle of the active layer and falls off sharply on both sides. It has been found that there is an increase in efficiency from 4.59% to 4.70% for 100 nm and 4.35% to 4.51% for 140 nm. Similarly, in Fig. 15 we show generation rate profiles obtained in blend structures prepared under uniform annealing time (60 min) and graded annealing time, where the annealing time and thickness of each region has been chosen such that the resulting profile has the same average as that obtained from the uniform blend, while retaining an optimum-like shape. Applying the two profiles to our model, we find an efficiency increase of 6% in case of

to ¼ 100 nm and slight increase for to ¼ 140 nm. Since the two profiles have the same overall absorption, we trace this improvement to the appropriate distribution of the generated carriers throughout the active layer. 3.3.2. Shifting the peak generation rate As mentioned previously, optical spacers like TiO2 and ZnO have been used extensively in order to shift the peak of the optical electric field, and hence the peak carrier generation rate, inside the active layer, resulting in enhanced device performance [15,19e22]. Here, we show the generation rate profiles obtained by performing FDTD simulation on the device structure in Fig. 1 after inserting TiO2 spacers of variable length in between the active layer and the metal cathode (Al). It is seen in Fig. 16 (to ¼ 100 nm) and Fig. 17 (to ¼ 140 nm) that insertion of a spacer may have either or both of two effects on the generation rate profiles depending on its length. Certain ranges of spacer length may cause the existing generation rate peak to shift towards right i.e. the active layer/ spacer interface. For example, for to ¼ 100 nm, this shift is seen in Fig. 16 (a) for spacer lengths of 10e20 nm, and Fig. 16(c) for spacer lengths of 60e100 nm. As for to ¼ 140 nm, spacer lengths of 50e100 nm cause this right shift of the generation rate peak. [Fig. 17 (b)]. The other effect observed is a simultaneous decrease in the value of the existing peak and the emergence of a steadily increasing peak at x ¼ 0 nm (PEDOT:PSS/active layer interface)

M. Islam et al. / Organic Electronics 46 (2017) 226e241

[Figs. 16 (b) and 17(a)]. Thus, we conclude that, depending on the active layer thickness, certain thicknesses of TiO2 spacer enable us to achieve a shift in the generation peak away from light. Oxides like TiO2 and ZnO are also used as electron transport layers [ETL] with P3HT:PCBM as the active layer in a BHJ OSC. Such structures are called p-i-n structures, in contrast to the p-i-metal structure described in Fig. 1, where the active layer is modelled as an intrinsic layer, with hole (HTL) and electron transport layers (ETL) on either sides [88]. It has been reported that, by manipulating the lengths of the two transport layers, it is possible to shift the peak of the generation rate profile [88]. We have already shown the effect of changing the ETL length for a fixed HTL length in Figs. 16 and 17. In Fig. 18, the TiO2 ETL length is kept constant at a value that yields a generation rate peak very close to the middle of the active layer. These values are 7 nm for to ¼ 100 nm and 100 nm for to ¼ 140 nm. After that, the PEDOT:PSS HTL length is varied over a range of 10 nm to 100 nm, few of which are shown in Fig. 18 for both active layer lengths. It is seen that the shape of the curve for both cases remains same and no shift in the position of the peak values are seen. However, there is a decrease in the value of the generation rate peak with an increase in the HTL length for

239

to ¼ 100 nm. As for, to ¼ 140 nm, there are essentially two peaks in the active layer, one at the edge (x ¼ 0 nm) and another at the middle (x ¼ 70 nm). For small lengths of PEDOT:PSS (up to 20 nm), the peak at the edge is greater than the middle peak, however for further increase in the HTL length, the peak at the middle increases over that at the edge. In Fig. 19, we manipulate the lengths of both transport layers and take advantage of their peak shifting and decrease of peak value effects, in order to keep the average generation rate same in each of to ¼ 100 nm [Fig. 19(a)] and to ¼ 140 nm [Fig. 19(b)]. The profiles show that by controlling both TL lengths, it is possible to retain the shape of the profile while attaining the desired shift towards the middle of the curve [in contrast to Figs. 16 and 17 where the curve shapes changes]. Lastly, we combine the two methods of control we have shown here and in Section 3.3.1 and present the generation profiles in Fig. 20 where the device structure has been modified with a graded blend ratio structure (1:4/1:1/1:4) for selectively increasing the generation rate in the middle region and using a spacer as well in order to show shifting effect in a graded structure. The solid curves in Fig. 20(a) and (b) are found for uniform blend ratio structures of

Fig. 19. Varying lengths of both transport layers to retain the shape of the generation profile and obtain a shift in the peak generation rate towards the middle of the active layer. (a) to ¼ 100 nm. Profile 1': ETL ¼ 0 nm; HTL ¼ 95 nm; Profile 2': ETL ¼ 8 nm; HTL ¼ 47 nm; Profile 3': ETL ¼ 10 nm; HTL ¼ 40 nm; Profile 4': ETL ¼ 15 nm; HTL ¼ 18 nm. (b) to ¼ 140 nm. Profile 1'': ETL ¼ 0 nm; HTL ¼ 95 nm; Profile 2'': ETL ¼ 70 nm; HTL ¼ 60 nm; Profile 3'': ETL ¼ 75 nm; HTL ¼ 50 nm; Profile 4'': ETL ¼ 100 nm; HTL ¼ 60 nm.

Fig. 20. Generation profiles obtained using a graded blend ratio structure (1:4/1:1/1:4) along with an optical spacer. (a) to ¼ 100 nm. S.G.1: 18/52/30 nm (spacer ¼ 20 nm); S.G.2: 30/ 40/30 nm (spacer ¼ 30 nm); S.G.3: 28/32/40 nm (spacer ¼ 66 nm) (b) to ¼ 140 nm. S.G.1': 20/100/20 nm (spacer ¼ 42 nm) and S.G.2': 50/60/30 nm (spacer ¼ 70 nm).

240

M. Islam et al. / Organic Electronics 46 (2017) 226e241

to ¼ 100 nm and 140 nm respectively with no spacer. The graded structures along with optical spacers in each figure have the same average as their corresponding uniform profile. The two control mechanisms allow us to obtain profiles like S.G. 3 [Fig. 20(a)] and S.G. 2' [Fig. 20(b)] where a very sharp peak occurs close to the middle of the active layer and consequently result in better device efficiencies of 3.58% and 4.45% respectively. 4. Conclusion In this paper, we proposed a recursive solution to the driftdiffusion model that enables us to find the electrical characteristics of BHJ OSC. Our method of deriving the semi-analytical expression involved discretizing a given generation profile with assumed discretization step to generate a piecewise constant generation profile, thus retaining the spatial distribution of the profile. Then we solved a series of ordinary differential equations with appropriate continuity boundary conditions. Using the characteristic of our proposed model to account for the spatial distribution of arbitrary carrier generation rate, we performed an elaborate study to analyze the effect of spatial distribution of generation profile on the device efficiency considering Gaussian shaped generation profiles that have the same average value. Our results show that a given total number of carriers should be distributed such that the resulting generation profile has a single narrow peak occurring around the middle of the active layer and then falls off sharply on both sides to attain the maximum device performance. Lastly, we also perform optical simulation to suggest modifications in the device structure in order to control the spatial distribution of the generation rate profiles. Thus, our results emphasize that manipulating the local value of the carrier generation rate along the active layer is of utmost importance to realize the true potential of BHJ OSCs. References €usermann, E. Knapp, M. Moos, N.A. Reinke, T. Flatz, B. Ruhstaller, J. Appl. [1] R. Ha Phys. 106 (2009) 104507. das, Energy Environ. Sci. 2 (2009) 251. [2] B. Kippelen, J.-L. Bre [3] Y. Yang, F. Wudl, Adv. Mater. 21 (2009) 1401e1403. [4] F.C. Krebs, Sol. Energ. Mat. Sol. Cells 93 (2009) 394e412. [5] B.A. Gregg, M.C. Hanna, J. Appl. Phys. 93 (2003) 3605. [6] G. Li, L. Liu, F. Wei, S. Xia, X. Qian, IEEE J. Photovoltaics 2 (2012) 320e340. [7] G. Li, L. Liu, IEEE Nanotechnol. Mag. 5 (2011) 18e24. , S. Cho, N. Coates, J.S. Moon, et al., Nat. Photonics 3 [8] S.H. Park, A. Roy, S. Beaupre (2009) 297e302. [9] L.J.A. Koster, E.C.P. Smits, V.D. Mihailetchi, P.W.M. Blom, Phys. Rev. B 72 (2005) 085205. [10] H. Hu, K. Jiang, G. Yang, J. Liu, Z. Li, H. Lin, et al., J. Am. Chem. Soc. 137 (2015) 14149e14157. [11] M.A. Green, K. Emery, Y. Hishikawa, W. Warta, E.D. Dunlop, Prog. Photovoltaics Res. Appl. 24 (2016) 905e913. [12] M.C. Scharber, Adv. Mater. 28 (2016) 1994e2001. [13] J. You, L. Dou, K. Yoshimura, T. Kato, K. Ohya, T. Moriarty, et al., Nat. Commun. 4 (2013) 1446. [14] S. Altazin, R. Clerc, R. Gwoziecki, G. Pananakakis, G. Ghibaudo, C. Serbutoviez, Appl. Phys. Lett. 99 (2011) 143301. [15] A. Roy, S.H. Park, S. Cowan, M.H. Tong, S. Cho, K. Lee, et al., Appl. Phys. Lett. 95 (2009) 013302. [16] Q.-D. Ou, Y.-Q. Li, J.-X. Tang, Adv. Sci. 3 (2016) 1600123. [17] L. Chen, W.E.I. Sha, W.C.H. Choy, Opt. Express 20 (2012) 8175. [18] K. Tvingstedt, V. Andersson, F. Zhang, O. Ingan€ as, Appl. Phys. Lett. 91 (2007) 123514. [19] J. Gilot, I. Barbu, M.M. Wienk, R.A.J. Janssen, Appl. Phys. Lett. 91 (2007) 113520. [20] A.K.K. Kyaw, D.H. Wang, D. Wynands, J. Zhang, T.-Q. Nguyen, G.C. Bazan, et al., Nano Lett. 13 (2013) 3796e3801. €nsel, H. Zettl, G. Krausch, R. Kisselev, M. Thelakkat, H.W. Schmidt, Adv. [21] H. Ha Mater. 15 (2003) 2056e2060. [22] J.Y. Kim, S.H. Kim, H.H. Lee, K. Lee, W. Ma, X. Gong, et al., Adv. Mater. 18 (2006) 572e576. [23] P.-P. Cheng, L. Zhou, J.-A. Li, Y.-Q. Li, S.-T. Lee, J.-X. Tang, Org. Electron. 14 (2013) 2158e2163. [24] J.-D. Chen, L. Zhou, Q.-D. Ou, Y.-Q. Li, S. Shen, S.-T. Lee, et al., Adv. Energy

Mater. 4 (2014) 1301777. [25] R. Wang, L.-H. Xu, Y.-Q. Li, L. Zhou, C. Li, Q.-D. Ou, et al., Adv. Opt. Mater. 3 (2015) 139, 139. [26] J.-D. Chen, C. Cui, Y.-Q. Li, L. Zhou, Q.-D. Ou, C. Li, et al., Adv. Mater. 27 (2015) 1035e1041. [27] R.F. Service, Science 332 (2011) 293, 293. [28] L.A.A. Pettersson, L.S. Roman, O. Ingan€ as, J. Appl. Phys. 86 (1999) 487e496. [29] D.W. Sievers, V. Shrotriya, Y. Yang, J. Appl. Phys. 100 (2006) 114509. [30] F. Monestier, J.-J. Simon, P. Torchio, L. Escoubas, F. Flory, S. Bailly, et al., Sol. Energ. Mat. Sol. Cells 91 (2007) 405e410. [31] P. Kumar, S.C. Jain, V. Kumar, S. Chand, R.P. Tandon, J. Appl. Phys. 105 (2009) 104507. [32] R.A. Marsh, C.R. McNeill, A. Abrusci, A.R. Campbell, R.H. Friend, Nano Lett. 8 (2008) 1393e1398. [33] M.M. Chowdhury, M.K. Alam, Sol. Energ. Mat. Sol. Cells 132 (2015) 107e117. [34] S.M. Arnab, M.Z. Kabir, J. Appl. Phys. 115 (2014) 034504. [35] M.M. Chowdhury, M.K. Alam, Sol. Energ. 126 (2016) 64e72. [36] M.S. Islam, Org. Electron. 41 (2017) 143e156. [37] D.A. Neamen, Semiconductor Physics and Devices: Basic Principles, fourth ed., Irwin, Chicago, 1997. [38] C.L. Braun, J. Chem. Phys. 80 (1984) 4157. [39] L. Onsager, Phys. Rev. 54 (1938) 554e557. [40] V.D. Mihailetchi, L.J.A. Koster, J.C. Hummelen, P.W.M. Blom, Phys. Rev. Lett. 93 (2004) 216601. [41] G. Namkoong, J. Kong, M. Samson, I.-W. Hwang, K. Lee, Org. Electron. 14 (2013) 74e79. [42] Comsol Inc., Comsol Multiphysics 3.5a, 2008. https://www.Comsol.Com. [43] H. Shen, P. Bienstman, B. Maes, J. Appl. Phys. 106 (2009) 073109. [44] C.H. Chou, W.L. Kwan, Z. Hong, L.-M. Chen, Y. Yang, Adv. Mater. 23 (2011) 1282e1286. [45] P.P. Boix, G. Garcia-Belmonte, U. Muneecas, M. Neophytou, C. Waldauf, R. Pacios, Appl. Phys. Lett. 95 (2009) 233302. [46] C.W. Liang, W.F. Su, L. Wang, Appl. Phys. Lett. 95 (2009) 133303. [47] D. Duche, P. Torchio, L. Escoubas, F. Monestier, J.-J. Simon, F. Flory, et al., Sol. Energ. Mat. Sol. Cells 93 (2009) 1377e1382. [48] W.H. Lee, S.Y. Chuang, H.L. Chen, W.F. Su, C.H. Lin, Thin Solid Films 518 (2010) 7450e7454. [49] H. Hoppe, N.S. Sariciftci, J. Mater. Res. 19 (2004) 1924e1945. [50] H. Hoppe, N.S. Sariciftci, J. Mater. Chem. 16 (2006) 45e61. [51] M. Campoy-Quiles, T. Ferenczi, T. Agostinelli, P.G. Etchegoin, Y. Kim, T.D. Anthopoulos, et al., Nat. Mater. 7 (2008) 158e164. [52] H. Hoppe, T. Glatzel, M. Niggemann, W. Schwinger, F. Schaeffler, A. Hinsch, et al., Thin Solid Films 511e512 (2006) 587e592. [53] J.K.J. van Duren, X. Yang, J. Loos, C.W.T. Bulle-Lieuwma, A.B. Sieval, J.C. Hummelen, et al., Adv. Funct. Mater. 14 (2004) 425e434. [54] J.A. Renz, T. Keller, M. Schneider, S. Shokhovets, K.D. Jandt, G. Gobsch, et al., Sol. Energ. Mat. Sol. Cells 93 (2009) 508e513. [55] G. Li, V. Shrotriya, J. Huang, Y. Yao, T. Moriarty, K. Emery, et al., Nat. Mater. 4 (2005) 864e868. [56] Q. Liu, Z. Liu, X. Zhang, N. Zhang, L. Yang, S. Yin, et al., Appl. Phys. Lett. 92 (2008) 223303. [57] X. Yang, J. Loos, Macromolecules 40 (2007) 1353e1362. [58] J. Peet, J.Y. Kim, N.E. Coates, W.L. Ma, D. Moses, A.J. Heeger, et al., Nat. Mater. 6 (2007) 497e500. [59] F. Padinger, R.S. Rittberger, N.S. Sariciftci, Adv. Funct. Mater. 13 (2003) 85e88. [60] Y.-C. Huang, S.-Y. Chuang, M.-C. Wu, H.-L. Chen, C.-W. Chen, W.-F. Su, J. Appl. Phys. 106 (2009) 034506. [61] L.H. Nguyen, H. Hoppe, T. Erb, S. Günes, G. Gobsch, N.S. Sariciftci, Adv. Funct. Mater. 17 (2007) 1071e1078. €renklau, G. de With, H. Hoppe, J. Loos, Adv. Funct. Mater. [62] S.S. van Bavel, M. Ba 20 (2010) 1458e1463. [63] D. Chirvase, J. Parisi, J.C. Hummelen, V. Dyakonov, Nanotechnology 15 (2004) 1317e1323. [64] B. Kadem, A. Hassan, Energy Procedia 74 (2015) 439e445. [65] A. Supriyanto, A. Mustaqim, M. Agustin, A.H. Ramelan, Suyitno, E.S. Rosa, et al., in: IOP Conference Series: Materials Science and Engineering, vol. 107, 2016, p. 012050. [66] E. Lafalce, P. Toglia, C. Zhang, X. Jiang, Appl. Phys. Lett. 100 (2012) 213306. [67] O. Ourida, B.M. Said, in: Science Academy Transactions on Renewable Energy SystemsEngineering and Technology, vol. 1, 2011, pp. 90e92. e, P. Lu, M. Piche , P. Mascher, Q. Chen, et al., 7099 (2008) [68] Y. Ding, R. Valle 709919. [69] B. Kadem, A. Hassan, W. Cranton, J. Mater. Sci. Mater. Electron. 27 (2016) 7038e7048. [70] Y. Kim, S.A. Choulis, J. Nelson, D.D.C. Bradley, S. Cook, J.R. Durrant, Appl. Phys. Lett. 86 (2005) 063502. [71] G. Li, V. Shrotriya, Y. Yao, J. Huang, Y. Yang, J. Mater. Chem. 17 (2007) 3126. [72] S.-Y. Chuang, H.-L. Chen, W.-H. Lee, Y.-C. Huang, W.-F. Su, W.-M. Jen, et al., J. Mater. Chem. 19 (2009) 5554. [73] W.-Y. Wong, X.-Z. Wang, Z. He, A.B. Djurisi c, C.-T. Yip, K.-Y. Cheung, et al., Nat. Mater. 6 (2007) 521e527. [74] T. Fukuda, K. Suzuki, N. Yoshimoto, Y. Liao, Org. Electron. 33 (2016) 32e39. [75] B. Beyer, R. Pfeifer, J.K. Zettler, O.R. Hild, K. Leo, J. Phys. Chem. C 117 (2013) 9537e9542. [76] W. Tress, J. Photonics Energy 1 (2011) 011114.

M. Islam et al. / Organic Electronics 46 (2017) 226e241 [77] W. Tress, K. Leo, M. Riede, Sol. Energ. Mat. Sol. Cells 95 (2011) 2981e2986. [78] L. Chen, Y. Tang, X. Fan, C. Zhang, Z. Chu, D. Wang, et al., Org. Electron. 10 (2009) 724e728. [79] P. Sullivan, S. Heutz, S.M. Schultes, T.S. Jones, Appl. Phys. Lett. 84 (2004) 1210e1212. [80] F. Jin, B. Chu, W. Li, Z. Su, X. Yan, J. Wang, et al., Org. Electron. 15 (2014) 3756e3760. [81] R. Pandey, R.J. Holmes, Adv. Mater. 22 (2010) 5301e5305. [82] M. Kaur, A. Gopal, R.M. Davis, J.R. Heflin, Sol. Energ. Mat. Sol. Cells 93 (2009) 1779e1784. [83] S. Pfuetzner, J. Meiss, A. Petrich, M. Riede, K. Leo, Appl. Phys. Lett. 94 (2009)

241

253303. [84] A. Loiudice, A. Rizzo, G. Latini, C. Nobile, M. de Giorgi, G. Gigli, Sol. Energ. Mat. Sol. Cells 100 (2012) 147e152. [85] M. Drees, R.M. Davis, J.R. Heflin, J. Appl. Phys. 97 (2005) 036103. [86] D.H. Wang, H.K. Lee, D.-G. Choi, J.H. Park, O.O. Park, Appl. Phys. Lett. 95 (2009) 043505. [87] G.R. Fowles, Introduction to Modern Optics, second ed., Holt, Rinehart and Winston Press, New York, 1975. [88] B. Maennig, J. Drechsel, D. Gebeyehu, P. Simon, F. Kozlowski, A. Werner, et al., Appl. Phys. A 79 (2004) 1e14.