[picture] | Computer Modeling in Engineering & Sciences | [picture]

DOI: 10.32604/cmes.2023.028992

ARTICLE

Nonlinear Analysis of Organic Polymer Solar Cells Using Differential Quadrature Technique with Distinct and Unique Shape Function

Ola Ragb¹, Mokhtar Mohamed², Mohamed S. Matbuly¹ and Omer Civalek³,*

1 Department of Engineering Mathematics and Physics, Faculty of Engineering, Zagazig University, Zagazig, 44519, Egypt
2 Basic Science Department, Faculty of Engineering, Delta University for Science and Technology, Gamasa, 11152, Egypt
3 Department of Medical Research, China Medical University Hospital, China Medical University, Taichung, 404327, Taiwan
Corresponding Author: Omer Civalek. Email: civalek@yahoo.com
Received: 22 January 2023; Accepted: 23 April 2023; Published: 03 August 2023

Abstract: Four numerical schemes are introduced for the analysis of photocurrent transients in organic photovoltaic devices. The mathematical model for organic polymer solar cells contains a nonlinear diffusion–reaction partial differential equation system with electrostatic convection attached to a kinetic ordinary differential equation. To solve the problem, Polynomial-based differential quadrature, Sinc, and Discrete singular convolution are combined with block marching techniques. These schemes are employed to reduce the problem to a nonlinear algebraic system. The iterative quadrature technique is used to solve the reduced problem. The obtained results agreed with the previous exact one and the finite element method. Further, the effects of different times, different mobilities, different densities, different geminate pair distances, different geminate recombination rate constants, different generation efficiencies, and supporting conditions on photocurrent have been analyzed. The novelty of this paper is that these schemes for photocurrent transients in organic polymer solar cells have never been presented before, so the results may be useful for improving the performance of solar cells.

Keywords: Sinc; differential quadrature; block marching; organic solar cells; discrete singular convolution

Nomenclature

Jn | Flux density of electron
Jp | Flux density of proton
J | Total current density
μh,μp | Carrier mobilities for electron and hole
Gn,Gp | Carrier generation rates for electron and hole
Rnn,Rpp | Recombination rates for electron and hole
Dn,Dp | Charge carrier diffusion coefficients for electron and hole
q | Elementary charge (q>0)
t | Time
kdiss | Dissociation of charge pairs rate
G | Charge pairs generation rate
EB | Electron–hole pair binding energy
a | Pair distance
εr | Relative permittivity of polymer
ε0 | Permittivity of free space
KB | Boltzmann’s constant
T | Absolute temperature (300 k)
γ | Bimolecular recombination rate
krec | Monomolecular recombination rate
NC,Nv | State effective density
Bn,Bp | Electron and hole energy barriers
λ | Wave length
Egap | Band gap
Va | Applied voltage

1  Introduction

Energy disasters, the release of greenhouse gases, and global warming are known to be some of the most extreme threats to life in the future. Therefore, renewable energy has found great attention from consumer countries to reduce carbon-based energy supplies [1]. One of the renewable energy sources that has attracted great attention from research communities is organic solar cells (OSCs) [2–4].

Third-generation photovoltaic devices [5] divide into two branches: organic polymer cells [6–8] and electrochemical cells [9–11]. We will concern with organic polymer cells in this paper. This type of renewable energy source has numerous advantages, such as low cost of manufacturing [12,13] and easy processing on flexible substrates [14–16]. Bulk-heterojunction (BHJ) structure improves the performance of organic polymer cells where electron donor and acceptor materials are mixed in a solution and cast into a thin film sandwiched between two electrodes [7,8], as shown in Fig. 1. The dynamics of charge behaviors in photovoltaic cells are displayed by measurements of transient photocurrent [17]. The classic time-of-flight technique measures the behavior of charge carriers generated according to a light, short pulse. Carriers are generated beside the electrode to allow extracting mobility and transit time [18,19].

Figure 1: Schematic of a heterojunction organic solar cell with its energy band diagram

Many researchers have examined the design of solar cells experimentally [20], but only a limited number of cases have been solved analytically and numerically. de Falco et al. [21] numerically studied photocurrent transients in organic polymer solar cells using the finite element and Newton-Raphson techniques, and they also produced an analytical solution for the problem. Buxton et al. [22] used the finite difference method to present computer simulations of polymer solar cells. Hwang et al. [23] investigated the photocurrent transients of organic solar cells using numerical simulations of the drift-diffusion formulas. van Mensfoort et al. [24] described two alternative iterative methods for obtaining solutions to the drift-diffusion equation. Koster et al. [25] presented the theory of Braun for showing the temperature and field dependence of the photocurrent in PPV: PCBM blends. However, for such problems, computational ill-conditioning is to be expected [26,27]. El Karkri et al. [28] used one dimension program [AMPS-1D] to optimize the performance of organic solar cells based on (carbazole-methylthiophene), benzothiadiazole and thiophene [(Cbz-Mth)-B-T]2 as electron donors, and [6,6]-phenyl-C61-butyric acid methyl ester (PCBM) as an electron acceptor. Yousuf et al. [29] analyzed the p-type Cu2O layer as a hole transport layer on Cu (In, GA) Se2 (CIGS) absorber layer. The analysis was performed by utilizing the solar cell capacitance simulator, SCAPS-1D. Overall, the solar cell conversion efficiency improved from 18.72% to 26.62% by adding a Cu2O hole transport layer and substituting CdS with ZnSe. Islam [30] investigated the thin-film organic solar cell (OSC) performances in detail by improved analytical computation in this work. The generation of excitons inside OSC has been estimated by using the optical transfer matrix method (OTMM) to include the optical phenomena of the incident light. The effect of mobility and active layer thickness on the properties of bulk heterojunction solar cells was investigated using the drift-diffusion model by Sadoogi et al. [31].

The differential quadrature method (DQM) is able to achieve accurate results without exerting high effort [28]. DQM comes in a variety of versions, each with its own shape function and influence domain for each point. The most used versions are Polynomial (PDQM) [32–34], Sinc (SDQM) [35], and Discrete singular convolution (DSCDQM) [36–44]. Ragb et al. [45–46] studied the performance of composite solar cell model consisting of perovskite absorber layer (CH3NH3PbI3), electron transport layer (TiO2 or PCBM), and hole transport layer (Spiro-OMeTAD or CuI) via sinc and discrete singular convolution quadrature techniques. Also, they presented different numerical methods that are used for the first time in solving Perovskite solar cells (PSCs). Ragb et al. [47] offered Improved electromagnetism-like algorithm and Differential quadrature technique to estimate of photovoltaic system and tracking power peaks of PV array under partial shading for three tests such as Schutten STP6-120/36 (Polycrystalline), Shell SM-55 (Monocrystalline) and PVM 752 GaAs (Thin film). Further, Ragb et al. [48] demonstrated improved electromagnetism-like algorithm and differential quadrature approach to evaluate the parameters of photovoltaic single, double- and three-diode model such as Kyocera polycrystalline (KC200GT), polycrystalline (Solarex MSX-60) and monocrystalline (R.T.C France).

According to the knowledge of the authors, the previous methods are not presented for photocurrent transients in organic polymer solar cells. A numerical scheme based on these methods is introduced to solve systems of nonlinear differential equations. A block marching technique is presented to approximate the time derivative. These methods are used to reduce the problem to a nonlinear algebraic system, which is then solved iteratively [32]. For each scheme, a MATLAB program is used to solve the problem, followed by a comparison of previous exact and finite element method results. The efficiency and convergence have been verified for each scheme (Polynomial-based differential quadrature, Sinc, and Discrete singular convolution) via presented results. Furthermore, the effects of different times, different mobilities, different densities, different geminate pair distances, different geminate recombination rate constants, different generation efficiencies, and different supporting conditions on photocurrent are investigated in detail.

2  Formulation of the Problem

Consider that bulk heterojunction solar cells [20,21,23,49] have been modeled using one-dimensional drift-diffusion equations. These equations involve the solution of the continuity equations for hole (p), electron (n), and charge pair (X) densities.

Based on equation of Poisson the electrostatic potential (ϕ) can be written as [20,21,23,49]:

−div(ε∇ϕ)=q(p−n),(1)

The continuity equations that are governing the charge transport in the device [20,21,23,49]:

∂n∂t=divJn+Gn−Rnn,(2)

∂p∂t=divJp+Gp−Rpp,(3)

where Gn=Gp=kdissX, Rnn=Rpp=γpn

The volume density of geminate pairs is [20,21,23,49]:

∂X∂t=G(x,t)+γpn−(kdiss+krec)X,(4)

The relation between electron (hole) flux densities and electrostatic terms (drift and diffusion) [20,21,23,49]:

Jn=Dn∇n−µhn∇ϕ,(5)

Jp=Dp∇p+μpp∇ϕ,Dn=μhkBTq,Dp=μpkBTq(6)

J=q(Jp−Jn),q>0(7)

We consider the coefficient krec as a given constant, but kdiss related to the electric field magnitude according to the following equations [6]:

E=−∂ϕ∂x,(8)

kdiss=3kR4πa3e−EB/kBT(1+b+b23+b318+b4180+⋯),EB=q24πεrε0a,KR=q(μh+μp)εrε0(9)

where all parameters used in the simulation of a polymer solar cells are defined in Table 1 [20,21,23,49].

According to Braun’s theory [25], the charge transfer state with binding energy EB is considered an intermediate state, leading to modification of the free charge carrier generation G(x,t)=I0G, where G is the exciton generation rate and I0=kdisskdiss+krec is the probability of the dissociation of the charge transfer (charge generation efficiency). The effect of a charge transfer state can be observed in experiments with a strongly pronounced electric field and a temperature-dependent photocurrent [50].

Substituting from Eqs. (5) and (6) into (1)–(4), the governing equations can be written as:

∂n∂t=∂∂x[µhkBTq∂n∂x−µhn∂ϕ∂x]+kdissX−γnp,(10)

∂p∂t=∂∂x[µpkBTq∂p∂x+µpp∂ϕ∂x]+kdissX−γnp,(11)

∂X∂t=I0G+γnp−(kdiss+krec)X,(12)

−ε∂2ϕ∂x2=q(p−n),(13)

To get the initial conditions, we will solve the system of equations at steady state by setting ∂X∂t=∂n∂t=∂p∂t=0 in Eqs. (10)–(12), and then the governing equations can be written as:

∂∂x[µhkBTq∂n∂x−µhn∂ϕ∂x]+kdissX−γnp=0,(14)

∂∂x[µpkBTq∂p∂x+µpp∂ϕ∂x]+kdissX−γnp=0,(15)

(kdiss+krec)X−γnp=I0G,(16)

The boundary conditions can be described as [20,21,23,49]:

[n(0)p(0)n(L)p(L)]=[Nce(−Bn/kBT)Nve−(Egap−Bn/kBT)Nce−(Egap−Bp/kBT)Nve(−Bp/kBT)],(17)

The energy barrier for carriers is zero when the contact for carrier was ohmic. The potential boundary condition is:

ϕ(L)−ϕ(0)=Egap−Bn−Bpq−Va,(18)

3  Method of Solution

3.1 Polynomial Based Differential Quadrature Method (PDQM)

The shape function is used in this method is Lagrange interpolation polynomial such that the unknown u and its nth derivatives can be approached as a weighted linear sum of nodal values, ui, as follows [32–34]:

u(xi)=∑j=1N∏k=1N⁡(xi−xk)(xi−xj)∏j=1,j≠kN⁡(xj−xk)u(xj),(i=1:N),(19)

∂nu∂xn|x=xi=∑j=1NCij(n)u(xj)(i=1:N)(20)

where u expresses to (n,p,Xandϕ). N is mesh size numbers. Cij(1) is the 1st derivative weighting coefficient, which can be calculated as [34]:

Cij(1)={1(xi−xj)∏k=1,k≠i,jN⁡(xi−xk)(xj−xk)i≠j−∑j=1,j≠iNCij(1)i=j(21)

Higher order derivatives weighting coefficients, can be determined as:

[Cij(n)]=[Cij(1)][Cij(n−1)],(n=2,3,…)(22)

3.2 Sinc Differential Quadrature Method (SDQM)

The shape function in this method is cardinal sine function. The unknown u and its nth derivatives can be approached as a weighted linear sum of nodal values, ui as follows [35]:

Sj(xi,hx)=sin⁡[π(xi−xj)/hx]π(xi−xj)/hx,where hx>0 is the step size(23)

u(xi)=∑j=−NNsin⁡[π(xi−xj)/hx]π(xi−xj)/hxu(xj),(i=−N:N),(24)

∂u∂x|x=xi=∑j=−NNCij(1)u(xj),∂2u∂x2|x=xi=∑j=−NNCij(2)u(xj),(i=−N:N),(25)

where u expresses to (n,p,Xandϕ). N is mesh size numbers. hx is grid size. Cij(1)andCij(2) can be calculated by differentiating Eqs. (23) and (24) as:

Cij(1)={(−1)i−jhx(i−j)i≠j0i=j,Cij(2)={(−1)i−j+12[hx(i−j)]2i≠j−13[πhx]2i=j(26)

3.3 Discrete Singular Convolution Differential Quadrature Method (DSCDQM)

Based on singular convolution defined as:

F(t)=(T∗η)(t)=∫−∞∞T(t−x)η(x)dx(27)

where T(t−x)and η(t) are the singular kernel and the test function’s space element, respectively.

The discrete singular convolution (DSC) technique depends on several kernels. The shape functions used in this method depend on these kernels so that the unknown u and its derivatives are approached as a weighted linear sum of ui [51–54]. Two types of kernels will be used as:

3.3.1 Delta Lagrange Kernel (DLK)

To approximate the derivative of a given function with respect to a space variable at a discrete point, the DSC is usually using a weighted linear combination of the function values at 2M+1 points in the direction of the space variable [32,33,53]. Delta Lagrange Kernel (DLK) can be used as a shape function such that the unknown (u) and its derivatives can be approximated as follows:

u(xi)=∑j=−MM∏k=−MM⁡(xi−xk)(xi−xj)∏j=−M,j1kM⁡(xj−xk)u(xj),(i=−N:N),M≥1(28)

∂u∂x|x=xi=∑j=−MMCij(1)u(xj),∂2u∂x2|x=xi=∑j=−MMCij(2)u(xj),(29)

where 2M+1 is the effective computational band width which is centered around x and is usually smaller than the whole computational domain.

Cij(1)andCij(2) are defined as [32,33,53]:

Cij(1)={1(xi−xj)∏M(xi−xk)(xj−xk)k=−M,k≠i,ji≠j−∑j=−M,j≠iMCij(1)i=j,Cij(2)={2(Cij(1)Cii(1)−Cij(1)(xi−xj))i≠j−∑j=−M,j≠iMCij(2)i=j,(30)

3.3.2 Regularized Shannon Kernel (RSK) [53]

To comparisons and illustrations, the regularized Shannon’s kernel, is used to solve this problem. DSCDQM-RSK is assumed that the unknown u with its derivatives is the approximation weighted linear sum of nodal values. Hence, the regularized Shannon’s delta kernel is discretized by [53]:

u(xi)=∑j=−MM⟨sin⁡[π(xi−xj)/hx]π(xi−xj)/hxe−((xi−xj)22σ2)⟩u(xj),(i=−N:N),σ=(r∗hx)>0(31)

∂u∂x|x=xi=∑j=−MMCij(1)u(xj),∂2u∂x2|x=xi=∑j=−MMCij(2)u(xj),(32)

where σ is regularization parameter and (r) is a computational parameter. It is also known that the truncation error is very small due to the use of the Gaussian regularizer, the above formulation given by Eq. (31) is practical and has an essentially compact support for numerical interpolation.

Cij(1)andCij(2) can be defined as [32,33]:

Cij(1)={(−1)i−jhx(i−j)e−hx2((i−j)22σ2),i≠j0i=j,Cij(2)={(2(−1)i−j+1hx2(i−j)2+1σ2)e−hx2((i−j)22σ2),i≠j−1σ2−π23hx2i=j(33)

Discrete form of Eqs. (13)–(16), can be write as below for each methods using the related weighting coefficients defined by Eqs. (19)–(33):

µhkBTq∑j=1NCij(2)nj−µhni∑j=1NCij(2)ϕj−µh∑k=1NCik(1)ϕk∑j=1NCij(1)nj+kdiss∑j=1NδijXj−γ∑k=1Nδiknk∑j=1Nδijpj=0,(34)

µpkBTq∑j=1NCij(2)pj+µppi∑j=1NCij(2)ϕj+µp∑k=1NCik(1)ϕk∑j=1NCij(1)pj+kdiss∑j=1NδijXj−γ∑k=1Nδiknk∑j=1Nδijpj=0,(35)

−ε∑j=1NCij(2)ϕj=q(∑j=1Nδijpj−∑j=1Nδijnj),(36)

(kdiss+krec)∑j=1NδijXj−γ∑k=1Nδiknk∑j=1Nδijpj=I0G,(37)

Therefore, we can get the initial values for (n,p,X,ϕ) at t=0.

3.4 Block-Marching Technique with Differential Quadrature Discretization

Eqs. (10)–(12) are considered one-dimensional, time-dependent equations. To solve these equations, the block-marching technique is used [34]. In time direction, the block-marching method divides the semi-infinite domain into many time intervals, δt1,δt2,δt3,.... Each block consists of one time interval δt and x-direction domain between zero and Lx. At block 1, the initial condition is considered the functional values at the bottom boundary (t=0). The top boundary (t=δt1) functional values are unknowns obtained by solving Eqs. (10)–(13). Block 2 considers the functional values at the bottom boundary (t=δt1) to be the top boundary of block 1. The numerical solution at the top of block 2 is obtained by repeating the steps from block 1. To meet the deadline, we repeat the preceding procedures several times. The differential quadrature (DQ) weighting coefficients are based on grid point coordinates. The grid point distribution for all blocks is equal, and these weighting coefficients can be used at all blocks by setting δt=δt1=δt2=δt3.... In the x and t-directions of Kth block, the following grid sizes are presented [54]:

xi=Lx2[1−cos⁡(i−1N−1π)],(i=1:N)(38)

tk=(K−1)δt+δt2[1−cos⁡(k−1L−1π)],(k=1:L)(39)

where N is the grid size number in x direction and L is the time level number in the block.

On suitable substitution from Eqs. (19)–(33) into (10)–(13), the problem can be reduced to the following:

∑j=1LC¯ij(1)nj=µhkBTq∑j=1NCij(2)nj−µhni∑j=1NCij(2)ϕj−µh∑k=1NCik(1)ϕk∑j=1NCij(1)nj +kdiss∑j=1NδijXj−γ∑k=1Nδiknk∑j=1Nδijpj,(40)

∑j=1LC¯ij(1)pj=µpkBTq∑j=1NCij(2)pj+µppi∑j=1NCij(2)ϕj+µp∑k=1NCik(1)ϕk∑j=1NCij(1)pj +kdiss∑j=1NδijXj−γ∑k=1Nδiknk∑j=1Nδijpj,(41)

−ε∑j=1NCij(2)ϕj=q(∑j=1Nδijpj−∑j=1Nδijnj),(42)

∑j=1LC¯ij(1)Xj=I0G−(kdiss+krec)∑j=1NδijXj−γ∑k=1Nδiknk∑j=1Nδijpj,(43)

where C¯ij(1) is the weighting coefficient of 1st derivative with respect to time. For all results, the boundary conditions Eqs. (17) and (18) are augmented in the governing Eqs. (40)–(43). Then, using the iterative quadrature technique [32,33] to obtain linear algebraic problem as:

Firstly, solving the Eqs. (40)–(43) as linear system:

∑j=1LC¯ij(1)nj=µhkBTq∑j=1NCij(2)nj+kdiss∑j=1NδijXj,(44)

∑j=1LC¯ij(1)pj=µpkBTq∑j=1NCij(2)pj+kdiss∑j=1NδijXj,(45)

−ε∑j=1NCij(2)ϕj=q(∑j=1Nδijpj−∑j=1Nδijnj),(46)

∑j=1LC¯ij(1)Xj=I0G−(kdiss+krec)∑j=1NδijXj,(47)

Then, the following iterative system is performed until reaching the required convergence [32,33]:

|ns+1ns|<1|ps+1ps|<1 where s = 0, 1, 2, ...

∑j=1LC¯ij(1)ns+1,j=µhkBTq∑j=1NCij(2)ns+1,j−µhns,i∑j=1NCij(2)ϕj−µh∑k=1NCik(1)ϕk∑j=1NCij(1)ns,j+kdiss∑j=1NδijXj−γ∑k=1Nδikns,k∑j=1Nδijps+1,j,(48)

∑j=1LC¯ij(1)ps+1,j=µpkBTq∑j=1NCij(2)ps+1,j+µpps,i∑j=1NCij(2)ϕj+µp∑k=1NCik(1)ϕk∑j=1NCij(1)ps,j+kdiss∑j=1NδijXj−γ∑k=1Nδikns+1,k∑j=1Nδijps,j,(49)

−ε∑j=1NCij(2)ϕj=q(∑j=1Nδijps+1,j−∑j=1Nδijns+1,j),(50)

∑j=1LC¯ij(1)Xj=I0G−(kdiss+krec)∑j=1NδijXj−γ∑k=1Nδikns,k∑j=1Nδijps,j(51)

4  Numerical Results

The following numerical results of each scheme explain convergence and efficiency for the analysis of transient photocurrent to improve power efficiency. The computational characteristics of each scheme are modified to obtain accurate results with error of order ≤10−8. For the present results, material parameters are taken from previous studies [21,23]. For the PDQ scheme [32–34], the non-regular grid, with Gauss-Chebyshev-Lobatto discretization as Eq. (38) is used to solve this problem. The grid sizes (N) are between 3 and 11. Table 1 shows that increasing the number of grid points leads to a very good convergence rate for electron density (n). The obtained results agreed with the previous exact one [23] and finite element method [21] over grid size ≥9 and execution time of 6.35 s, as shown in Table 1. This table also demonstrates that the current method produces accurate results.

Table 2 for the SincDQ scheme shows convergence of the obtained results with uniform grid sizes (N) ranging from 3 to 11. Furthermore, it demonstrates that the convergence percentage for electron density (n) increases with the mesh size numbers. Over grid size ≥9 and execution time of 4.2 s, these results agreed with the analytical one [23] and the finite element method [21]. We notice that the SincDQ scheme takes less time to execute than the PDQ scheme. As a result, it is more efficient than PDQM for analyzing organic polymer solar cells.

For DSCDQ scheme depended on kernel of delta Lagrange Table 3 demonstrates convergence of the results, which have been obtained with uniform grids (N) between 3 and 11, and bandwidth 2M+1 ranges between 3 and 9. These findings agreed with those of the analytical one [23] and finite element approach [21] over N≥7, time for performance 2.7 s, and 2M+1≥3. We notice that the DSCDQ-DLK scheme takes less time to execute than the SincDQ scheme. As a result, it is more efficient than PDQM and SincDQM for analyzing organic polymer solar cells. The DSCDQ scheme is based on the kernel of regularized Shannon (RSK) Tables 4 and 5, which explain convergence of the obtained results with regular grids (N) ranging from 3 and 9. With regularization parameter 1.0hx≤σ≤1.75hx the bandwidth is (3≤2M+1≤7). Over mesh size ≥7, bandwidth ≥3, execution time 2.1 s, and regularization parameter σ=1.75hx, these results agreed with analytical one [23] and finite element approach [21]. We notice that DSCDQM-RSK scheme takes less time to execute than the DSCDQM-DLK, SincDQM and PDQM schemes. Therefore, DSCDQ-RSK scheme has the best efficiency comparing the examined schemes for organic polymer solar cells analysis.

Each scheme combined with block marching technique to ensure the accuracy of the obtained results. Table 6 demonstrates transient currents due to DSCDQM-RSK with block marching at high and low intensities (G=4.3×1026m−3s−1, G=4.3×1029m−3s−1) with different times at δt=0.1 µs,number of timelevels (L=4). The results agreed with analytical one [23] over number of blocks ≥10 and execution time 3.63 s. Table 7 shows transient currents at low (G = 4.3×1026m−3s−1) and high (G = 4.3×1029m−3s−1) intensities with different times at δt=0.1µs,L=8. The results agreed with analytical one [23] over number of blocks ≥6 and execution time 3.15 s. Therefore, increasing number of time levels in the block (L) is more efficient for photocurrent transient analysis. Table 8 shows transient currents at low (G = 4.3×1026m−3s−1) and high (G = 4.3×1029m−3s−1) intensities with different times at δt=0.05 µs, L=12. The results agreed with analytical one [23] over number of blocks ≥3 and execution time 2.93 s. Therefore, decreasing time interval (δt) is more efficient for photocurrent transient analysis. Tables 6–8 show that the value of transient currents is larger at high intensity. Further, the accuracy and efficiency of the result depend on choice of time interval and number of blocks. From Tables 3–8, DSCDQM-RSK with block marching is best choice for efficient for photocurrent transients analysis at (grid size ≥7, bandwidth ≥3, regularization parameter σ=1.75hx, δt=0.05 µs, L=12 and Execution time is 2.93 s). Then, effect of different times, mobilities, densities, geminate pair distances, geminate recombination rate constants, generation efficiency and supporting conditions on photocurrent are investigated in a detailed parametric study.

Firstly, Fig. 2 is introduced to show the best wavelength (λ) at different distances from the cathode for photocurrent transients analysis. This figure shows that the electron density increased with increasing distances from the cathode. Its mean that distances from cathode of 220 nm, the values of active layer optical absorptivity is higher. This ensures that more photons are absorbed by thicker layers than by thinner ones. Increasing wavelength (λ) makes the curves flatter. In addition, decreasing distance from the cathode leads to faster stability for curves, which explains that the internal quantum efficiency (IQE) is sensitive to the absorption layer thickness. Furthermore, IQE is simply the quotient of absorptivity, and a thicker active layer means the photogenerated carrier must travel a longer average distance out of the device. Therefore, the best value for wavelength is above (λ≥600 nm) at a distance from the cathode = (70 nm). As well as, we will introduce the parametric study based on these values where high electron density is obtained at low and high intensity as shown in Fig. 2.

Figure 2: Variation of the free electron density with wavelength (λ) at different distances from cathode at (a) low intensity G=4.3×1026m−3s−1 (b) high intensity G=4.3×1029m−3s−1 with krec=105s−1 and µh=µp=2×10−4 cm2V−1s−1

Fig. 3 explains the normalized photocurrent transients at low density G = 4.3×1026m−3s−1 for different monomolecular recombination rate (krec), different geminate pair distances (a) and times. The normalized photocurrent transients increase with increasing monomolecular recombination rate (krec), different geminate pair distances (a) and times. Further, curves achieve a faster stability with increasing geminate pair distances. So, we will choose (a = 1.5 nm) in the parametric study. Figs. 4 and 5 demonstrate the normalized photocurrent transients at low and high intensities for different monomolecular recombination rate (krec) and mobilities µh,µp. In Figs. 4 and 5, at high intensity, we notice that the time increases. The increasing in monomolecular recombination rate (krec) and mobilities µh,µp leads to increase in time because the high rate of monomolecular recombination controls the dynamics of charge pair. From Fig. 5, at high intensities, we notice that bimolecular recombination has a great effect by the difference between efficiency (I0) and net efficiency (Pnet) of charge generation. According to the increasing in mobility, the normalized photocurrent transients increase continuously because increasing mobility leads to increase recombination.

Figure 3: Variation of transient currents at G=4.3×1026m−3s−1 (low) with time (µs), different geminate pair distances (nm) and different monomolecular recombination rate constants (a) krec=106s−1 and (b) krec=105s−1. µh=µp=2×10−4cm2s−1

Figure 4: Variation of transient currents with times, t, generation rates, G, and geminate recombination rate constants, krec, and different mobilities (a) μh=μp=2 ×10−4cm2v−1s−1 (b) μh=2×10−5,μp=2×10−4cm2v−1s−1

Figure 5: Transient currents at low intensity (G = 4.3 × 1026 m−3s−1) and high intensity (G = 4.3 × 1029 m−3s−1) for (a) and (G = 4.3 × 1030 m−3s−1) for (b) with different parameters (µh,µp,krec,G,I0 and Pnet)

The velocity of infinite surface recombination causes more losses in surface recombination. This shows that mobility increasing will not effect on power conversion efficiency (PCE). So, Figs. 3–5 show the active-importance of mobility and layer thickness in calculating the PSC efficiency. Therefore, we will choose the value of mobility is µh=µp=2×10−4cm2s−1 in the parametric study. Fig. 6 demonstrates the charge pair densities at high intensity (G = 4.3 × 1029 m−3s−1) and electric field with different distances from cathode. The electric field has a maximum value in the middle. The rate of charge pair dissociation depends on the electric field. Then, the charge pair increases in the middle. Further, this figure shows that the value of the device electric field for high light intensity its deviation around its mean value (E = −ΔV/distance from cathode). ΔV is the total voltage drop through the device.

Figure 6: Variation of charge pair and electric field at high intensity (G = 4.3 × 1029 m−3s−1) with distance from cathode. µh=µp=2×10−4cm2s−1

Fig. 7 depicts the effect of various gap energies (Egap) and the recombination prefactor (ξ) on free electron density at low intensity and high charge generation efficiency (G=4.3×1026m−3s−1,krec=105s−1). The electron density increases at the lowest gap energy (Egap = 1) and the recombination prefactor (ξ=0). The effective band gap of the absorption layer plays a pivotal role in improving PSC performance. Across a large band gap, the electrons are thermally excited, and therefore the electron density in a higher band gap is less than that in a lower band gap. Otherwise, increasing the temperature leads to the electron being excited into the conduction band, so the electron density will increase. This leads to higher efficiency of the solar cell. Also, all of the parameters examined increase as recombination strength decreases. Fig. 8 explains the variation of transient currents with times (t), energy barriers (Bn, Bp), and different intensities. Increased energy barriers result in more normalized photocurrent in less time. The normalized photocurrent at high intensity requires less time than at low intensity. Figs. 2–8 demonstrate that the interface between the active region and the contact is critical in determining PSC performance. Fig. 9 shows the variation of the electron density at high intensity with distance from the cathode in case of krec=105s−1 and krec=107s−1. The electron density is consistent at low times. As the density of electron approaches the steady state, the electron density is peaked in the center at krec=105s−1 but more uniform at krec=107s−1. Electron density increases with time in both cases, as shown in Fig. 9.

Figure 7: Variation of the free electron density with distance from cathode at (a) different gap energies (Egap) (b) different Recombination prefactor (ξ) with krec=105s−1 and µh=µp=2×10−4 cm2V−1s−1 at G=4.3×1026m−3s−1

Figure 8: Variation of transient currents with times, t, energy barrier, (Bn, Bp), and different intensities (a) Low intensity, G = 4.3 × 1026 (b) High intensity, G = 4.3 × 1029. μh=μp=2×10−4cm2v−1s−1,krec=105s−1,a=1.5nm

Figure 9: The electron distribution at high intensity with distance from cathode (a) krec=105s−1, and (b) krec=107s−1, G=4.3×1029m−3s−1, µh=µp=2×10−4 cm2V−1s−1

Fig. 10 explains the electron distribution with high and low charge generation efficiency at low intensity, different distances from cathode and times. The shape of electron density at high charge generation efficiency is similar to the shape at low efficiency.

Figure 10: Variation of the free electron density with distance from cathode, different times (t), at (a) krec=105s−1 (b) krec=107s−1, µh=µp=2×10−4 cm2V−1s−1,G=4.3×1026m−3s−1

To exit the device, all electrons move to the cathode. As a result, the electron has a high density towards the cathode. It is worth noting that equal mobilities result in hole density figures that are an exact mirror reflection of the electron density. The free electron generation rate kdissX with krec=105s−1 and krec=107s−1 charge generation efficiency at high intensity, different distances from cathode and times are shown in Fig. 11. For high dissociation probabilities, the values are sickly-looking in the midpoint of the device, but for low dissociation probabilities, the values are regular during the central percentage of the device. At low and high charge generation efficiency, the free electron generation rate increases with time. Fig. 12 shows the variation of the electron effective bimolecular recombination (1/γp) with time, distance from cathode and both of high and low efficiency of charge generation at high intensity. As time increases, (1/γp) decreases at low and high efficiency of charge generation.

Figure 11: The free electron generation rate, kdissX, at high intensity with distance from cathode (a) krec=105s−1, and (b) krec=107s−1. µh=µp=2×10−4 cm2V−1s−1, G=4.3×1029m−3s−1

Figure 12: The variation of the electrons bimolecular recombination (1/γp) with distance from cathode, time at (G = 4.3 × 1029 m−3s−1) with (a) krec = 105 s−1 and (b) krec = 107s−1. µh=µp=2×10−4 cm2V−1s−1

5  Conclusions

Polynomial-based differential quadrature, Sinc and Discrete Singular Convolution techniques have been successfully employed for the analysis of photocurrent transients in organic photovoltaic devices. To solve one-dimensional, time-dependent equations, the block-marching technique is used. Further, the iterative quadrature technique is used to solve the reduced problem. A MATLAB program is performed to get a solution for this problem at each scheme. Then, block marching technique accuracy is validated by comparing the obtained results with an analytical one [23] over a number of blocks ≥3 at δt=0.05 µs, L=12, and execution time 2.93 s. From computational analysis it found that DSCDQM-RSK with (N  ≥  7,2M  +  1  ≥  3  and  σ  =  1.75  *  hx) is the most accuracy and efficiency values for solving this problem. The maximum error for these schemes is ≤10−8. Also, the minimum execution time for DSCDQM-RSK is 2.1 s. Effect of different times, different mobilities, different densities, different geminate pair distances, different geminate recombination rate constants, different generation efficiency, and supporting conditions on photocurrent are investigated in a detailed parametric study by using the efficient technique DSCDQM-RSK. For all results, it is found that:

•   The best value for wave length is above (λ≥600 nm) at distance from cathode = (70 nm).

•   The normalized photocurrent transients increase with increasing monomolecular recombination rate (krec), different geminate pair distances (a) and times (t).

•   In the central of the device, the electric field has a maximum value. The charge pair dissociation rate depends on the electric field. Then, the charge pair population increases in the central of the device. Across large band gap, the electrons are thermally excited, then the electron density increases at the lowest gap energy (Egap = 1) and recombination prefactor (ξ=0).

•   The increasing in energy barrier leads to more normalized photocurrent at less time.

•   On the way to the steady state, the electron density increases; at more efficiency of charge generation, the electron density reaches its peak at the midpoint, whereas at less efficiency of generation, the electron spreading becomes more regular.

•   For high generation efficiency the values are peaked in the middle, but for low generation efficiency the values are regular across the middle of the device.

•   As time increases, the effective bimolecular recombination time for electrons decreases at low and high charge generation efficiency.

•   Equal mobilities make hole density figures being the exact mirror image of the electron density. The active-layer thickness, mobility, energy barrier and band gap are vital components in defining the significant efficiency and performance of the PSC device.

These techniques, which aided in the creation of diverse and novel parametric studies with less effort and time, will be useful for improving the performance of solar cells. This type of renewable electrical energy has numerous advantages for use in many fields of the industrial revolution.

Funding Statement: The authors received no specific funding for this study.

Conflicts of Interest: The authors declare that they have no conflicts of interest to report regarding the present study.

References

1. Waldau, A. J. (2005). PV status: Research, solar cell production and market implementation of photovoltaics. Refocus, 6(3), 20–23.

2. MacKenzie, R. C. I., Balderrama, V. S., Schmeisser, S., Stoof, R., Greedy, S. et al. (2016). Loss mechanisms in high efficiency polymer solar cells. Advanced Energy Materials, 6(4), 1501742.

3. Jørgensen, M., Norrman, K., Krebs, F. C. (2008). Stability/degradation of polymer solar cells. Solar Energy Materials and Solar Cells, 92(7), 686–714.

4. Deibel, C., Dyakonov, V. (2010). Polymer-fullerene bulk heterojunction solar cells. Reports on Progress in Physics, 73(9), 96401.

5. Pillai, S. A., Green, M. A. (2010). Plasmonics for photovoltaic applications. Solar Energy Materials and Solar Cells, 94(9), 1481–1486.

6. Mihailetchi, V. D., Koster, L. J. A., Hummelen, J. C., Blom, P. W. M. (2004). Photocurrent generation in polymer-fullerene bulk heterojunctions. Physical Review Letters, 93(21), 216601.

7. Günes, S., Neugebauer, H., Sariciftci, N. S. (2007). Conjugated polymer-based organic solar cells. Chemical Reviews, 107(4), 1324–1338.

8. Mayer, A. C., Scully, S. R., Hardin, B. E., Rowell, M. W., McGehee, M. D. (2007). Polymer-based solar cells. Materials Today, 10(11), 28–33.

9. Bai, Y., Cao, Y., Zhang, J., Wang, M., Li, R. et al. (2008). High-performance dye-sensitized solar cells based on solvent-free electrolytes produced from eutectic melts. Nature Materials, 7(8), 626–630.

10. Grätzel, M. (2005). Solar energy conversion by dye-sensitized photovoltaic cells. Inorganic Chemistry, 44(20), 6841–6851.

11. Grätzel, M. (2001). Photoelectrochemical cells. Nature, 414(6861), 338–344.

12. Krebs, F. C., Senkovskyy, V., Kiriy, A. (2010). Preorganization of nanostructured inks for roll-to-roll-coated polymer solar cells. IEEE Journal of Selected Topics in Quantum Electronics, 16(6), 1821–1826.

13. Krebs, F. C. (2009). All solution roll-to-roll processed polymer solar cells free from indium-tin-oxide and vacuum coating steps. Organic Electronics, 10(5), 761–768.

14. Wagner, J., Fritz, T., Böttcher, H. (1993). Computer modelling of organic thin film solar cells. I. Exciton model of photocurrent generation. Physica Status Solidi, 136(2), 423–432.

15. Hu, T., Li, F., Yuan, K., Chen, Y. (2013). Efficiency and air-stability improvement of flexible inverted polymer solar cells using ZnO/poly(ethylene glycol) hybrids as cathode buffer layers. ACS Applied Materials & Interfaces, 5(12), 5763–5770.

16. Liu, Y. F., Feng, J., Cui, H. F., Yin, D., Song, J. F. et al. (2012). Highly flexible inverted organic solar cells with improved performance by using an ultrasmooth Ag cathode. Applied Physics Letters, 101(13), 133303.

17. Redecker, M., Bradley, D. D. C., Inbasekaran, M., Woo, E. P. (1998). Nondispersive hole transport in an electroluminescent polyfluorene. Applied Physics Letters, 73(11), 1565–1567.

18. Coehoorn, R., Pasveer, W. F., Bobbert, P. A., Michels, M. A. J. (2005). Charge-carrier concentration dependence of the hopping mobility in organic materials with Gaussian disorder. Physical Review B, 72(15), 155206.

19. Campbell, A. J., Bradley, D. D. C., Antoniadis, H. (2001). Dispersive electron transport in an electroluminescent polyfluorene copolymer measured by the current integration time-of-flight method. Applied Physics Letters, 79(14), 2133–2135.

20. Liang, C., Wang, Y., Li, D., Ji, X., Zhang, F. et al. (2014). Modeling and simulation of bulk heterojunction polymer solar cells. Solar Energy Materials and Solar Cells, 127, 67–86.

21. de Falco, C., Sacco, R., Verri, M. (2010). Analytical and numerical study of photocurrent transients in nanoscale organic solar cells. Computer Methods in Applied Mechanics and Engineering, 199(25–28), 1722–1732. https://doi.org/10.1016/j.cma.2010.01.018

22. Buxton, G. A., Clarke, N. (2006). Computer simulation of polymer solar cells. Modelling and Simulation in Materials Science and Engineering, 15(2), 13–26.

23. Hwang, I., Greenham, N. C. (2008). Modeling photocurrent transients in organic solar cells. Nanotechnology, 19(42), 424012.

24. van Mensfoort, S. L. M., Coehoorn, R. (2008). Effect of Gaussian disorder on the voltage dependence of the current density in sandwich-type devices based on organic semiconductors. Physical Review B, 78(8), 85207.

25. Koster, L. J. A., Smits, E. C. P., Mihailetchi, V. D., Blom, P. W. M. (2005). Device model for the operation of polymer/fullerene bulk heterojunction solar cells. Physical Review B, 72(8), 85205.

26. Lazarov, R. D., Zikatanov, L. T. (2012). An exponential fitting scheme for general convection-diffusion equations on tetrahedral meshes. arXiv Preprint arXiv:1211.0869.

27. Gatti, E., Micheletti, S., Sacco, R. (1998). A new Galerkin framework for the drift-diffusion equation in semiconductors. East-West Journal of Numerical Mathematics, 6(2), 101–135.

28. El Karkri, A., El Malki, Z., Bouachrine, M., Serein-Spirau, F., Sotiropoulos, J. M. (2020). Characterization and simulation study of organic solar cells based on donor-acceptor (D–π–A) molecular materials. RSC Advances, 10(32), 18816–18823.

29. Yousuf, M. H., Saeed, F., Tauqeer, H. A. (2022). Numerical investigation of Cu2O as hole transport layer for high-efficiency CIGS solar cell. Preprints.org, 2021100326.

30. Islam, M. S. (2021). In-depth analysis of organic solar cells using transport equation and optical transfer matrix method with detailed analytical derivations. Energies, 14(3), 735.

31. Sadoogi, N., Rostami, A., Faridpak, B., Farrokhifar, M. (2021). Performance analysis of organic solar cells: Opto-electrical modeling and simulation. Engineering Science and Technology, an International Journal, 24(1), 229–235.

32. Ragb, O., Mohamed, M., Matbuly, M. S. (2019a). Free vibration of a piezoelectric nanobeam resting on nonlinear Winkler-Pasternak foundation by quadrature methods. Heliyon, 5(6), e01856. https://doi.org/10.1016/j.heliyon.2019.e01856

33. Ragb, O., Mohamed, M., Matbuly, M. S. (2019b). Vibration analysis of magneto-electro-thermo NanoBeam resting on nonlinear elastic foundation using sinc and discrete singular convolution differential quadrature method. Modern Applied Science, 13(7), 49–79.

34. Shu, C. (2012). Differential quadrature and its application in engineering. Berlin: Springer Science & Business Media.

35. Korkmaz, A., Dağ, İ. (2011). Shock wave simulations using sinc differential quadrature method. Engineering Computations, 28(6), 654–674.

36. Allam, O., Draiche, K., Bousahla, A. A., Bourada, F., Tounsi, A. et al. (2020). A generalized 4-unknown refined theory for bending and free vibration analysis of laminated composite and sandwich plates and shells. Computers and Concrete, 26(2), 185–201.

37. Ouakad, H. M., Sedighi, H. M., Al-Qahtani, H. M. (2020). Forward and backward whirling of a spinning nanotube nano-rotor assuming gyroscopic effects. Advances in Nano Research, 8(3), 245–254.

38. Bakoura, A., Bourada, F., Bousahla, A. A., Tounsi, A., Benrahou, K. H. et al. (2021). Buckling analysis of functionally graded plates using HSDT in conjunction with the stress function method. Computers and Concrete, 27(1), 73–83.

39. Bendenia, N., Zidour, M., Bousahla, A. A., Bourada, F., Tounsi, A. et al. (2020). Deflections, stresses and free vibration studies of FG-CNT reinforced sandwich plates resting on Pasternak elastic foundation. Computers and Concrete, 26(3), 213–226.

40. Bekkaye, T. H. L., Fahsi, B., Bousahla, A. A., Bourada, F., Tounsi, A. et al. (2020). Porosity-dependent mechanical behaviors of FG plate using refined trigonometric shear deformation theory. Computers and Concrete, 26(5), 439–450.

41. Chuan, M. W., Wong, K. L., Hamzah, A., Rusli, S., Alias, N. E. et al. (2021). Device modelling and performance analysis of two-dimensional AlSi3 ballistic nanotransistor. Advances in Nano Research, 10(1), 91–99.

42. Ebrahimi, F., Kokaba, M., Shaghaghi, G., Selvamani, R. (2020). Dynamic characteristics of hygro-magneto-thermo-electrical nanobeam with non-ideal boundary conditions. Advances in Nano Research, 8(2), 169–182.

43. Low, S., Shon, Y. S. (2018). Molecular interactions between pre-formed metal nanoparticles and graphene families. Advances in Nano Research, 6(4), 357–375.

44. Singh, P. P., Azam, M. S. (2021). Size dependent vibration of embedded functionally graded nanoplate in hygrothermal environment by Rayleigh-Ritz method. Advances in Nano Research, 10(1), 25–42.

45. Ragb, O., Mohamed, M., Matbuly, M. S., Civalekm, O. (2021). Sinc and discrete singular convolution for analysis of three-layer composite of perovskite solar cell. International Journal of Energy Research, 46(4), 1–22.

46. Ragb, O., Mohamed, M., Matbuly, M. S., Civalekm, O. (2021). An accurate numerical approach for studying perovskite solar cells. International Journal of Energy Research, 45(11), 16456–16477.

47. Ragb, O., Bakr, H. (2023). A new technique for estimation of photovoltaic system and tracking power peaks of PV array under partial shading. Energy, 268, 126680.

48. Ragb, O., Bakr, H., Civalek, O. (2022). Parameters identifcation for photovoltaic system via improved electromagnetism like approach and quadrature technique. International Journal of Energy and Environmental Engineering, 1913, 1–25.

49. Blom, P. W. M., Mihailetchi, V. D., Koster, L. J. A., Markov, D. E. (2007). Device physics of polymer: Fullerene bulk heterojunction solar cells. Advanced Materials, 19(12), 1551–1566.

50. Braun, C. L. (1984). Electric field assisted dissociation of charge transfer states as a mechanism of photocarrier production. The Journal of Chemical Physics, 80(9), 4157–4161.

51. Glatthaar, M., Riede, M., Keegan, N., Sylvester-Hvid, K., Zimmermann, B. et al. (2007). Efficiency limiting factors of organic bulk heterojunction solar cells identified by electrical impedance spectroscopy. Solar Energy Materials and Solar Cells, 91(5), 390–393.

52. Tress, W., Leo, K., Riede, M. (2012). Optimum mobility, contact properties, and open-circuit voltage of organic solar cells: A drift-diffusion simulation study. Physical Review B, 85(15), 155201.

53. Civalek, Ö. (2017). Free vibration of carboforced (CNTR) and functionally graded shells and plates based on FSDT via discrete singular convolution method. Composites Part B: Engineering, 111(3), 45–59.

54. Shu, C., Richards, B. E. (1992). Application of generalized differential quadrature to solve two-dimensional incompressible Navier-Stokes equations. International Journal for Numerical Methods in Fluids, 15(7), 791–798.

[picture] | This work is licensed under a Creative Commons Attribution 4.0 International License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.