A Fractional Diffusion Model for Dye-Sensitized Solar Cells.

Dye-sensitized solar cells have continued to receive much attention since their introduction by O`Regan and Grätzel in 1991. Modelling charge transfer during the sensitization process is one of several active research areas for the development of dye-sensitized solar cells in order to control and improve their performance and efficiency. Mathematical models for transport of electron density inside nanoporous semiconductors based on diffusion equations have been shown to give good agreement with results observed experimentally. However, the process of charge transfer in dye-sensitized solar cells is complicated and many issues are in need of further investigation, such as the effect of the porous structure of the semiconductor and the recombination of electrons at the interfaces between the semiconductor and electrolyte couple. This paper proposes a new model for electron transport inside the conduction band of a dye-sensitized solar cell comprising of TiO 2 as its nanoporous semiconductor. This model is based on fractional diffusion equations, taking into consideration the random walk network of TiO 2 . Finally, the paper presents numerical solutions of the fractional diffusion model to demonstrate the effect of the fractal geometry of TiO 2 on the fundamental performance parameters of dye-sensitized solar cells, such as the short-circuit current density, open-circuit voltage and efficiency.

1. Introduction

Dye-sensitized solar cells (DSSCs) were first introduced by O’Regan and Grätzel in their fundamental 1991 paper [1], providing a viable low-cost alternative for renewable solar energy. Functionally, DSSCs operate by a photosensitive dye using absorbed sunlight to inject excited electrons into a nanoporous semiconductor. This approach relives DSSCs of the need for a costly and high purity semiconductor as opposed to the predominant silicon solar cells that have been at the forefront of solar energy since 1954 [2].

Mathematical modelling for DSSCs began shortly afterwards and continues to offer unique insight into the current–voltage relationship and ultimately the efficiency of the DSSC. While some early mathematical models of DSSCs used Maxwell’s equations and electric potentials [3], Gregg [4] outlined the prevalence of the photoinduced potential as more influential than electric fields in DSSCs. Södergren et al. [5] were the first to model DSSCs by a simple diffusion equation together with analytical expressions for the short-circuit current density and the open-circuit voltage . Cao et al. [6] extended the model to include time-dependence, leading to the partial differential equation where is the conduction band electron density at depth and time , is the diffusion coefficient, is the incident photon flux, is the absorption coefficient, is the recombination constant and is the equilibrium electron density.

Anta et al. [7] considered a fully nonlinear diffusion equation for DSSCs and investigated the effect of different power-law diffusivities, which was analysed using Lie symmetry by Maldon et al. 2020 [8]. Andrade et al. 2011 [9] and Gacemi et al. 2013 [10] proposed diffusion equations for modelling the electrolyte concentrations, which were solved analytically by Maldon and Thamwattana in 2019 [11].

Fractional calculus has existed conceptually as long as calculus itself. For close to three centuries, fractional calculus had only one known application in Abel’s 1823 tautochrone problem [12] until Nigmatullin [13] suggested a fractional diffusion equation for modelling media exhibiting a fractal geometry. In 2000, Henry and Wearne [14] developed the standard fractional diffusion equation based on continuous-time random-walk (CTRW) models. Henry and Wearne’s model [14], together with the CTRW simulation of by Nelson [15], provides motivation for this paper to model DSSCs using fractional partial differential equations.

Modelling DSSCs with fractional calculus is a relatively unexplored research direction. So far, there is only one paper by Sibatov et al. in 2014 [16]. In Sibatov et al. [16], they consider the role of trap states within the network and its effect on the electron hole density. To address the fact that previous mathematical models for DSSCs do not generally consider the effect of the porous network in [6,17], we develop our model based on generalised fractional diffusion-reaction equations, taking this effect into account.

For this study, we use the Caputo fractional derivative, which is given by [18] where is the order of the fractional derivative, , denotes the classical derivative of f with respect to its variable of order and is the usual Gamma function

The literature employs several definitions for the fractional derivative in modelling diffusion in media exhibiting a fractal geometry. Henry and Wearne [14] suggest the Riemann–Liouville definition, but also remark that the Caputo definition has also seen some use [14,19]. However, Baeumer et al. [20] comment on the use of Caputo derivatives for space-fractional diffusion models that positivity is not preserved under vanishing Neumann boundary conditions. We find that problem is alleviated by the use of a time-fractional derivative on the diffusion term of our equation. Also, the presence of spatially dependent source terms and the combination of Dirichlet and Neumann boundary conditions enjoy a greater level of compatibility with the Caputo fractional derivative [21].

2. Mathematical Model

In this paper, we adopt a Caputo fractional derivative in time on the diffusion term as shown by Henry and Wearne [14], resulting in the fractional partial differential equation (FPDE) where is the order of the Caputo fractional derivative in time and all other parameters retain their values as in Equation (1). Physically, the parameter is the exponent in the mean square displacement of the CTRW simulation [14]. The special case recovers the standard diffusion equation. Lower values for correspond to longer path lengths for electron transport through the nanoporous semiconductor [22]. The parameter is also strongly influenced by the porosity of the semiconductor. Benkstein et al. [22] found that increasing porosity led to a decrease in . We note that Equation (3) does not feature the term where denotes the inverse Laplace transformation. Though this term is critical for physically meaningful fractional diffusion equations [14], it vanishes in Equation (3) under the Caputo fractional derivative.

We prescribe boundary conditions as found in [5], with a Dirichlet boundary condition at and a Neumann boundary condition at together with a prescribed initial condition, namely where q is the standard electron charge, V is the bias voltage, m is the diode ideality factor, is Boltzmann’s constant and T is the temperature of the DSSC.

The diode equation is commonly used to compute the current–density relationship for solar cells, in which the current J as a function of bias voltage V is given by where is the dark saturation current density, given by Södergren et al. [5] in the form

To compute the short-circuit current density we use noting that the standard flux is recovered in the special case of linear diffusion ().

Given that the open-circuit voltage satisfies , we may compute the open-circuit voltage from

Maximising the power output over V, we obtain the maximum power point by where W is the Lambert-W function and . With , we compute the efficiency of the DSSC by where is the power of incident light.

3. Finite Difference Method

Finite difference methods (FDMs) have been used by Hu et al. [23] to solve parabolic FPDEs under the Caputo fractional derivative for fractional time derivatives, and Takeuchi et al. [24] have used finite difference methods to solve FPDEs under fractional spatial derivatives. We refer the reader to Li and Zeng [25] for a finite difference scheme for solving fractional ordinary differential equations over a finite interval under boundary conditions defined at the left boundary.

In this paper, we solve Equation (3) under a FDM scheme using expressions given by Oldham and Spanier [26]. All numerical computations are performed using numerical values of constants provided in Table 1 except for , which is given several values in Table 2.

Table 1

Parameter values for the Dye-Sensitized Solar Cell (DSSC) model.

Parameter	Value	Unit	Reference
D_0	10^−11	m^2s^−1	[7]
α	10^5	m^−1	[10]
d	5×10^−5	m	[7]
k_R	4×10^−8	s^−1	[7]
m	1	-	[5]
n_eq	10^22	m^−3	[27]
P_i	10	Wm^−2	[10]
φ	10^21	m^−2s^−1	[28]

Table 2

Values for , and for several values of at time s.

γ	J_sc(Am^−2)	V_oc(V)	η(%)
0.25	56.0248	0.6056	2.8144
0.39	70.7900	0.6117	3.5967
0.5	82.6370	0.6156	4.23
0.612	94.6085	0.6191	4.8741
0.75	108.9019	0.6227	5.6481
1	134.2799	0.6281	7.03

To numerically solve Equation (3) under boundary conditions (4)–(6) with a finite difference scheme, we use the L1 approximation for the fractional derivative given by Oldham and Spanier [26] in which where is given by

Let be the final simulation time. Discretise into nodes and into nodes, and let estimate the solution to Equation (3) under the boundary conditions (4)–(6) at the point . That is,

3.1. Nodes Determined by Boundary Conditions

To satisfy the initial condition (6), we set for all . For the boundary condition (4) at , we set for all . Finally, for the Neumann boundary condition (5) at we employ a ‘ghost node’ at and a central difference approximation for the first derivative at to set for all .

3.2. Iteration Algorithm

For the second row , we have

Given and , the finite difference iteration for numerically solving Equation (3) is given by

3.3. Estimate for Short-Circuit Current Density

To compute the short-circuit current density, we must estimate the electron flux at . We achieve this by a 10 point estimate of the form where is a constant for each . We determine these constants so that polynomials are perfectly estimated up to degree 9. This leads to the constants

4. Results and Discussion

Using the finite difference method, we numerically solve Equation (3) for several values of using 10 spatial nodes and 5000 temporal nodes over and s. To investigate the effect of the parameter , we plot numerical solutions for the special cases in Figure 1.

Figure 1

Plots of the numerical solution to Equation (3) against x and t for (top-left), (top-right), (bottom-left), (bottom-right).

Using the numerical solutions to Equation (3) and the flux estimate given by Equation (12), we are able to compute the short-circuit current density and the open-circuit voltage , leading to the overall efficiency . Table 2 shows the effect of the parameter on these DSSC performance parameters.

From Table 2, we see that the short-circuit current density increases when increases. Though the open-circuit voltage is not affected to the same extent, the efficiency is notably lower for decreased values of . We note the special case is equivalent to the standard diffusion equation without fractional derivatives, and the efficiency is in agreement with expected efficiencies for DSSCs [8].

As decreases, Figure 1 presents two primary trends to the numerical solution to Equation (3). Firstly, the time required to reach steady-state increases as decreases by comparison to the numerical solutions for the cases . This result is consistent with the observation that lower values for imply slower diffusion, based on the CTRW simulations. Benkstein et al. [22] show that decreases when porosity increases from to .

Secondly, the overall electron density is remarkably higher for the cases compared to . This suggests that the electron density is highly sensitive to the order of the fractional derivative. The standard gradient for flux would consequently produce extreme results, requiring the fractional derivative to rectify this issue. Numerically, the finite difference scheme presents stability problems and is computationally expensive for lower values of .

The parameter denotes the exponent for the power-law in the mean-square displacement [14]. We conclude from Figure 1 that lower values of lead to progressively less realistic behaviour for nanoporous semiconductors used in DSSCs, as the electron density dramatically increases when decreases. This observation is consistent with the longer path lengths associated with low values for in CTRW simulation. In particular, Ni et al. [27] show that extremely low porosities (such as ) show a significant reduction in efficiency (below ). From the numerical solution to Equation (3) with using 25 spatial nodes, the flux estimate given by Equation (12) and 50,000 temporal nodes to s, we find , and .

5. Conclusions

We propose a new mathematical model for evaluating the efficiency of dye-sensitized solar cells by using fractional diffusion to incorporate the fractal geometry of the semiconductor. Our results show that lower values of the mean square exponent lead to lower efficiencies, a result that is consistent with the literature [11,27]. In particular, Figure 11 of Ni et al. [27] shows that efficiency decreases when porosity increases above or decreases below , which suggests the relationship between porosity and efficiency is nonlinear. We note that the solution profile of the electron density presented here is similar to those obtained from nonlinear diffusion modelling [8], though the orders of magnitude differ significantly. This is due to the longer waiting times associated with lower values for , which slows down the diffusion process.

We also develop a finite difference scheme to numerically solve the fractional diffusion equation and provided a tenth-order estimate for obtaining the short-circuit current density. Together, this provides a comprehensive model for incorporating the effect of the random-walk behaviour of the nanoporous semiconductor on the performance of dye-sensitized solar cells.

Future consideration includes incorporating the role of the electrolyte couple by a pair of standard diffusion equations as seen in Maldon and Thamwattana [11].