Random walk numerical simulation for hopping transport at finite carrier concentrations : diffusion coefficient and transport energy concept w

The random walk numerical simulation (RWNS) method is used to compute diffusion coefficients for hopping transport in a fully disordered medium at finite carrier concentrations. We use Miller–Abrahams jumping rates and an exponential distribution of energies to compute the hopping times in the random walk simulation. The computed diffusion coefficient shows an exponential dependence with respect to Fermi-level and Arrhenius behavior with respect to temperature. This result indicates that there is a well-defined transport level implicit to the system dynamics. To establish the origin of this transport level we construct histograms to monitor the energies of the most visited sites. In addition, we construct ‘‘corrected’’ histograms where backward moves are removed. Since these moves do not contribute to transport, these histograms provide a better estimation of the effective transport level energy. The analysis of this concept in connection with the Fermi-level dependence of the diffusion coefficient and the regime of interest for the functioning of dye-sensitised solar cells is thoroughly discussed.


I. Introduction
The theoretical description of electron transport in disordered materials is a challenging issue with implications in the fields of dye-sensitised solar cells (DSC), 1 plastic solar cells, 2 organic light emitting diodes 3 and organic electronics. 4In these materials, transport of charge occurs by jumps of electrons between localized states, although extended states may also play a role.The transport rates are determined by two kinds of microscopic disorder: (1) energetic disorder characterized by a broad distribution of localized states 5 and (2) spatial disorder, related to the morphological features of the material. 6,7The correct description of the influence of these two kinds of disorder and their microscopic parameters on the transport features of the material is crucial to design of better performing devices.
Two main approaches have been applied so far to describe electron transport in these materials.The first is the classical multiple-trapping model, [8][9][10] in which transport occurs via extended states along a mobility edge (or conduction band) but it is slowed down by a succession of trapping-detrapping events in localized states.In this model, only energetic disorder is taken explicitly into account by means of the distribution of energies (relative to the mobility edge) characteristic of the ensemble of localized states.2][13][14] In the hopping model electron transport occurs by direct jumps between localized states and the hopping rates depend explicitly on both energy difference and spatial distance. 15 obtain usable analytical expressions for electron mobilities and diffusion coefficients requires making averages over spatial and energy disorder.This analysis is especially cumbersome in the context of the hopping model since both energetic and spatial disorder must be taken into account.However the analysis can be simplified if the distribution of energies for the localized states is very steep.7][18][19][20][21] The existence of an effective transport level reduces the hopping transport to multiple trapping, with the transport energy playing the role of a mobility edge.The transport energy concept has been utilized to derive a theoretical expression for the diffusion coefficient of electrons hopping in an exponential distribution of localized states. 16he transport energy has been shown to be affected by the fact that the system is not ideal, that is, correlations between carriers may play an important role.These correlations can be due to exclusion effects, which makes the transport energy depend on Fermi level position, 17,22 or due to energetic correlations between charges and dipoles. 23,24opping transport in amorphous semiconductors has been amply studied over the last decades mainly in relation to disordered inorganic semiconductors such as amorphous silicon, and in recent years also for organic conductors. 25ecently, the interest in electronic transport in the presence of an exponential distribution of states has increased with the advent of nanostructured wide bandgap semiconductors applied in DSC. 26,27Indeed, electron transport in DSC has been largely described using multiple trapping, arguments. 28or DSC using relatively thick TiO 2 porous nanocrystalline layers, electron transport may impose limitations to charge extraction. 29Since DSC operate at large electron densities, it is crucial to further determine the transport mechanism in these systems as a function of charge density and, especially at high Fermi levels, beyond the analytical approximations adopted previously. 16,28A recent report shows experimental results of diffusion coefficients in TiO 2 at large electron densities. 30n this paper we apply the random walk numerical simulation (RWNS) method [31][32][33][34][35][36] to obtain the jump diffusion coefficient in a hopping system with an exponential distribution of localized states and at finite carrier concentration.We use our calculations to cast light on the foundations of the transport energy approximation in this case.The RWNS method is a stochastic technique that permits us to analyse the transport mechanism for a particular transport model from first principles and with no approximations (as those sometimes applied to compute magnitudes such as the transport energy and the mobility 18,19,37,38 ).The density of localized states (energy distribution) is used as an input to construct a three-dimensional network of sites whose energies are allocated according to this distribution.The simulation is performed by implementing jumping rates characteristic of the selected transport mechanism.In this case we implement the hopping mechanism via the Miller-Abrahams jumping rates. 15The RWNS calculations yield the jump diffusion coefficient as a function of Fermi level and temperature. 31On the other hand, we have carried out our simulations on a network of randomly distributed sites instead of a simple cubic lattice.Placing the sites on an ordered spatial arrangement has been shown to affect the results for the carrier mobility. 39To work with a fully disordered system permits us to eliminate the effect of introducing an artificial spatial order on the simulation results.
In this work we have used the simulations to construct histograms of most visited energies so that the probability for the electrons to jump to target sites of specific energy can be calculated.The form of this histogram for jumps upward in energy will allow us to identify the existence of a well-defined maximum and how it depends on carrier concentration and Fermi level.As noted by Arkhipov et al., 19 the transport energy can differ noticeably from the energy of the most probable jump due to the influence of neighbored sites close in energy.These sites make carriers hop back and forth many times so that those moves do not contribute to transport and hence to the computation of the diffusion coefficient.The RWNS method makes it possible to remove those jumps from the calculation so that a better approximation to the ''effective transport energy'' can be obtained for the studied cases.
This numerical work is aimed at understanding recent experimental studies of electron transport in nanostructured semiconductors where exponential distributions and strong Fermi level variations are very common. 40,41We believe that the results here presented will be quite useful to guide and interpret future experimental work in DSC and related systems.

II. Random walk numerical simulation for hopping transport A Method and simulation details
The RWNS is a stochastic computational procedure that allows for a flexible description of transport of charge carriers in a network of traps without huge computational demands.This is especially useful in the context of nanostructured materials since the existence of spatial disorder coupled with a broad distribution of trap energies is characteristic of these systems.5,42 In this work we run the random walk simulation on a threedimensional network of traps distributed randomly and homogeneously in space.As mentioned above in this way we avoid the undesired influence of an artificial ordering in the system. 39However in a fully disordered network there exist traps that happen to be very close to each other, which is not likely to occur in real materials.Nevertheless, as we will show below, the effect of these very close pairs do not have a strong effect on the diffusion coefficient because they produce back-and-forth moves which do not contribute to transport.
In the RWNS calculation a certain number of carriers (which can represent either electrons or holes) are allowed to jump between neighboring traps.The formula here used to compute the hopping times for carriers jumping from a trap i to a trap j is derived from the well-known Miller-Abrahams hopping rates. 15In this work we use the following formulation of the Miller-Abrahams formula based on times (inverse of rates) where R is a random number distributed uniformly between 0 and 1, t 0 is the inverse of the attempt-to-jump frequency, r is the distance between the traps, a is the localization radius, and E j , E i are the energies of the target and starting traps, respectively.A random walk simulation based on times rather than on rates or probabilities 11 leads to the same results because a jump between two traps i and j for which eqn (1) predicts a long hopping time is equivalent to consider a very small jumping rate between the two traps (and vice versa).
Trap energies are extracted from the usual exponential distribution where N L is the total trap density, k B T 0 the width of the distribution, E is the trap energy (negative) and E 0 is the lower (higher) edge of the conduction (valence) band, for electrons and holes respectively.E 0 indicates the energy of extended states, in case they exist (such states are not a necessary assumption in our model).Hereafter we take E 0 = 0.
In this work we aim to obtain the properties of a hopping system as a function of the carrier concentration or Fermi level.To achieve this, the simulations are run under the condition that no more than one carrier is allowed per site at the same time.As we will see below this makes the carrier occupancy function to follow Fermi-Dirac statistics. 31,33In contrast to previous studies 17,22 the Fermi-Dirac function is not imposed a priori, but it arises naturally from the calculation instead.
The simulations were performed as follows (see Fig. 1).A random network of traps is generated and energies are allocated to the traps according to eqn (2).Carriers are then initially placed at random on the network of traps.For each carrier hopping times to neighbored traps are computed via eqn (1).This calculation is restricted to non-occupied traps within a certain cut-off radius r cut .The minimum of these hopping times and its corresponding target trap is identified and stored.The procedure is repeated for all carriers so that the jump with the minimum hopping time, called t min , can be executed.The hopping times of the rest of the carriers are then reduced by t min and the process is repeated in such a way that, for each simulation step, the carrier that happens to have the minimum hopping time is moving along the network and the simulation is advanced by time intervals of variable size t min .
Calculations were carried out with 1-100 carriers and the size of the simulation box ranged between 10 and 65 nm.A density of traps of N L = 10 27 m À3 was used in all cases.This corresponds to an average distance between traps of 1 nm.It must be stressed that, as traps are distributed randomly, hops can be executed for distances either longer or shorter than this averaged distance.Hereafter, the simulations are described by a label N/a L 3 where N is the number of carriers and a L the size of the simulation box in nm's.
The jump diffusion coefficient 5 D J is obtained from the mean square displacement according to 31 The mean square displacements are observed to be linear at longer times (normal diffusion).This allows extracting the diffusion coefficient from the slope of the curve in the time plot.

B Convergence tests
As mentioned above, to save computing time a certain cut-off distance r cut is introduced.Neighbors located beyond this distance are not considered as target sites.Since the hopping times in eqn (1) do depend on distance between traps, the cut-off distance should be large enough to ensure that the results are not significantly affected.In Fig. S1 in the ESIw the diffusion coefficient as a function of r cut is plotted for two values of the localization radius.
As it could be expected, a larger localization radius requires a larger cut-off radius to ensure convergence.Hence, for a = 0.5 nm and 2.5 nm a cut-off radius of 2.5 nm and 4.5 nm were found to be sufficient respectively.These are the parameters used henceforth.

C Energy level populations and one particle approximation
By running a long enough RWNS calculation it is possible to construct a histogram of the number of carriers that occupy levels of energy E. From this the corresponding occupancy probabilities can be extracted.As stated above, in our simulations it is observed that this probability resembles a Fermi-Dirac distribution with a well-defined Fermi level (see Fig. 2).The Fermi level is a monotonic function of the carrier density.
In previous work for RWNS with multiple-trapping it was found that it is possible to reproduce the diffusion coefficient of multi-carrier calculations by running a random walk simulation with just a single carrier and a modified trap energy distribution in which all traps with energies below E F are ignored. 31This approximation is found to work well for the hopping model here considered, although a constant shift is observed in the one-particle calculations (see Fig. S2 in the ESIw).This shift is not surprising if we take into account that the zero temperature approximation neglects the influence of unoccupied traps in the vicinity of the Fermi level that might contribute to transport with a constant weight that would depend on temperature but not on the position of the Fermi level.In any case, it must be noted that we are mainly interested in the behavior of the diffusion coefficient with respect to the Fermi level and temperature rather than in absolute values.

A Energy of the most probable jump and transport energy concept
][21] In this approach it is assumed that in equilibrium the transport is governed by a single energy level related to the fastest hop of a charge carrier.The most probable upward jump corresponds to an optimized combination of the distance and energy difference.For an exponential distribution of localized levels, the result 18,43 is that the fastest hops occur in the vicinity of the so-called transport energy, given by For a certain carrier (black circle), hopping times to neighbored traps are computed according to eqn (1).This computation is restricted to those traps within the cut-off r cut that are not occupied.Once all release times are computed the minimum hopping time and its corresponding target trap are looked for and stored.The process is repeated for all carriers in the simulation box so that the carrier having the minimum hopping time is identified.This carrier is then moved to its corresponding target trap.See text for more details.
independently of the energy of the starting site.This expression is obtained by maximizing the upward hopping rate for an average hopping distance.Alternatively the transport energy can be obtained by averaging the hopping rate below a certain energy value as reported by Arkhivov. 19,37his latter procedure has been put into question 44 due to the difficulty of considering the effect on transport of all relevant hops.In any case, the existence of a transport energy implies that the hopping model should behave in a very similar way to the multiple trapping model, where there is a transport level by definition.
In connection with the transport energy approximation, we have monitored the energies of the target sites for jumps upward in energy in the RWNS calculations.These values were used to construct a histogram of energies.Results can be found in Fig. 2 and Fig. S3 (ESIw) for two test cases (a = 0.5 nm and a = 2.0 nm, T 0 = 800 K, T = 275 K and densities corresponding to labels 100/12 3 and 100/15 3 ).The results reveal that most carrier moves take place in the vicinity of a certain energy that always lies (as expected) above the Fermi energy for each particular case.
In this work we make a critical analysis of the following assumption: the maximum of the energy histogram, E max , can be assimilated to the value of the transport energy.We must note that the former is just a simulation result whereas the latter is a theoretical concept obtained under certain approximations whose origin we want to test in this work using numerical simulation.Monte Carlo simulation has been used by Cleve et al. 38 and Novikov and Malliaras with similar purposes. 23However, Cleve et al. 38 investigate an empty system with no influence of the concentration of carriers.The work in ref. 23 investigates a Gaussian distribution that applies in organic conductors.
The most relevant feature of the present calculations is that E max is found to move upwards in the energy scale when the Fermi level is raised.A similar effect has been described recently for the transport energy with a Gaussian distribution of states. 22The variation of E max with density and Fermi level is shown in Fig. 3 for two characteristic temperatures (T 0 = 600 K and T 0 = 800 K).The calculations have been extended to the regime of very low densities, with Fermi levels between À0.17 and À0.61 eV and densities up to 7 Â 10 16 cm À3 .It must be noted that at low densities the statistics of the simulation is very poor, which increases the uncertainty of E max .This is extracted when the population distribution if found to relax to a Fermi-Dirac distribution with a well-defined Fermi level as explained in Section II.
The poor statistics in the low density limit are related to the occurrence of spurious peaks in the energy histograms.These are due to carriers jumping many times back and forth between sites that happen to be close in distance and in energy and tend to disappear when the simulation is very long.As a matter of fact the RWNS predictions at low densities do not converge to the classical value of eqn ( 6) as it could be expected.The reasons for this disagreement, in connection with the concept of effective transport energy of Arkhipov et al. 19 will be discussed in Section IIIC below.
In any case, if we assume that E max can be assimilated to the transport energy, the same behavior is found by Arkhipov and coworkers 17 and Li and coworkers. 22The carrier density dependence of E max is a result of the progressive filling of the localized states, that prevent carriers from hopping to neighbored sites for which the Miller-Abrahams formula yields high probability.The carriers are then forced to jump to levels of higher energies, hence producing a larger value of Fig. 2 Occupation probabilities (f(E), full symbols) and histograms of the energies of target sites (N(E), open symbols) in RWNS calculations for jumps upward in energy.The latter have been normalized with respect to the maxima.Simulations were carried out for a = 0.5 nm, T 0 = 800 K, T = 275 K and densities corresponding to labels 100/12 3 (circles) and 100/15 3 (squares).The following values are obtained from the simulations for both densities: E F = À0.22 eV (E max = À0.12) and E F = À0.26eV (E maz = À0.18eV), respectively.The solid line stands for the exponential trap distribution of T 0 = 800 K. the transport energy.At very low concentrations this filling effect is negligible and E max remains constant.Nevertheless, the real connection between E max and the transport energy is subtle and requires further analysis, as discussed below.

B Fermi level dependence of the diffusion coefficient
As explained before, the jump diffusion coefficient for carriers can be computed from the RWNS calculations as a function of Fermi level.Results in reduced units for two test cases (a = 0.5 nm, T 0 = 800 K, T = 275 K and a = 2.0 nm, T 0 = 800 K, T = 275 K) are presented in Fig. 4.
The simulation data show that the logarithm of diffusion coefficient scales almost linearly with Fermi level.Diffusion coefficients are found to be higher for large localization radius.This is explained by the fact that delocalization favors jumps to traps further apart and produces shorter average hopping times.
The exponential dependence of the diffusion coefficient with respect to the position of the Fermi level is analogous to the typical behavior of the multiple trapping model.This result indicates that there should exist a well-defined transport level that controls the transport of carriers under equilibrium conditions.However, the results presented in the previous subsection reveal that the energy of the most probable jump does move to higher energies when the trap distribution becomes progressively filled.This appears to be contradictory to the fact that there is a fixed transport energy.In the next subsection this issue is discussed and clarified.

C Diffusion coefficient and transport energy
The theory of diffusion [45][46][47] often allows to separate the kinetic or jump diffusion coefficient in two factors: where hr 2 i is an average hopping distance and hni is an average hopping frequency.In hopping transport, there is not a well defined separation between hopping at different distances and hopping at different energy levels.However, the rationale for the transport energy approximation is that the relevant jumps occur to a well defined level, and in this case eqn (7) may provide a useful approach to obtain analytical expressions for hopping transport as a function of Fermi level.The numerical simulations performed in this work constitute an excellent tool to check the validity of such approximations.Therefore, following the work from previous authors, 16,17 we compute the jump diffusion coefficient using eqn (7).According to the transport energy concept both quantities can be calculated from where n(E,E tr ) is the frequency for an upward hop from the energy E to the transport energy E tr (inverse of eqn ( 1)) at fixed distance r = hri.By applying the zero-temperature limit of the Fermi-Dirac distribution in eqn (8b) and introducing the classical value of eqn (6) for the transport energy, Bisquert found the following expression for the diffusion coefficient: 16 This theoretical expression predicts an exponential behavior with respect to the Fermi energy, in analogy with the multipletrapping result and in accordance with the simulation (see Fig. 4).However, the theoretical slope (27.71 eV À1 for T 0 = 800 K and a = 0.5 nm) is slightly larger than the simulation result.
In spite of this encouraging result, the exponential behavior of the diffusion coefficient is not consistent with the upward shift of the average hopping energies when the Fermi level is increased.As it can be observed in Fig. 2 and 3, the maximum of the energy histogram E max lies always above and approximately at a constant distance with respect to the Fermi level.If we would assume that E max can be assimilated to the transport energy, this behavior would lead to a constant diffusion coefficient according to eqn (9).

D Effective transport energy
To disentangle from the paradox posed in the previous subsection, the concept of effective transport energy of Arkhipov and coworkers 19 is especially useful.These authors make a distinction between the energy that controls transport at equilibrium conditions and the energy of the most probable View Online jumps.That these two are different has been already observed in Monte Carlo simulations for hopping systems in a Gaussian density of states. 43s mentioned above, RWNS calculations at low densities produce energy histograms with spurious peaks in the low energy region.These peaks arise from carriers jumping back and forth between neighboring sites.The consequence in the numerical simulation is that these ''oscillatory'' moves do not contribute to the diffusion of the carriers and therefore should be excluded in the estimation of the transport energy.
Bearing this is mind, we have extended the computation of the histograms of hopping energies to the situation in which backward jumps are ignored.To achieve that, the coordinates of the starting site are stored for every move so that when the carrier returns to its original position, the target energy is not used to compute the energy histogram, since these jumps do not produce the spatial displacement of the electrons.
Results for both types of energy histograms are presented in Fig. 5 for calculations with a single carrier in an empty exponential trap distribution and for a finite density corresponding to label 10/15 3 .The most visible feature is that the spurious peaks tend to disappear when backward jumps are ignored.However, sharp peaks are not completely removed.This is due to the fact that oscillatory moves between pairs of sites are not the only moves that do not contribute to transport.Carriers can get ''trapped'' between small groups of sites and follow circular trajectories before escaping, especially at lower energies.Nevertheless to remove these ''second-order'' moves is much more difficult in the numerical computation and goes beyond the scope of the present work.The occurrence of spurious peaks is magnified in the present calculations by the fact that we perform our simulations on a random network of traps.As mentioned above, this leads to the possibility of traps that happen to be very close to each other.This problem does not appear in the simulations of Ba¨ssler and coworkers, 11,38,43 which are executed on a cubic lattice.Simulations on-lattice reduce the numerical demands and produces results more in accordance to the assumptions of the theory (see eqn (8a) for instance) but at the cost of losing the subtleties of the positional disorder implicit to these kind of systems. 39 second feature of the corrected histograms is that the maximum, that we call E 0 max , lies at higher energies than in the original histogram.That the effective transport energy lies above the energy of the most probable jump is the main conclusion of the work of Arkhipov et al. 19 and it is confirmed in the present calculations.The simulations of Hartenstein and Ba¨ssler 43 and Cleve et al. 38 also predict energies for the most probable jump below the classical value of eqn (6).On the contrary, the computation of the histogram without backward jumps for a single carrier leads to a maximum much closer to the theoretical value of À0.26 eV predicted by eqn (6) (see Fig. 5).It must be born in mind that eqn (6) is obtained under the assumption that all hops occur at a constant average distance whereas in the simulation traps can be occasionally very close to each other and this induces the appearance of the oscillatory moves mentioned above.
The energy of the maximum of the corrected histograms, E 0 max , allows us to propose a better estimate for the transport energy that is implicit to the diffusion coefficient dependence on the Fermi level.Results for this are collected in Fig. 6 together with the values of the most probable jump as computed in subsection IIIA.Here it is observed that E 0 max , lies always above E max and that it converges to the classical value of eqn ( 6) at low densities.

E Simulated diffusion coefficient versus theoretical predictions
The concepts introduced in the previous subsections allow us to use the E 0 max values from the simulated histograms to produce theoretical values of the diffusion coefficient according to eqn (7) and (8).The results, together with the simulated data and the predictions of the approximate formulas ( 6) and ( 9) can be found in Fig. 4.
We observe that eqn (7) and (8) with the transport energy assimilated to E 0 max reproduce Bisquert's formula at low Fermi levels.This is not surprising if we take into account that the simulation reproduces the classical value of eqn (6) in this regime as explained above.The agreement between the theories and the simulation is also good in the low Fermi level region for the localized case.However, as we move towards large carrier densities the theoretical values separate from Bisquert's formula although they tend to remain close to the simulated data.This effect is basically a consequence that eqn (9) is derived under the assumption that the Fermi level is well below the transport level.By introducing the proper Fermi-Dirac function in eqn (8) the match with respect to the simulation is improved.This effect is more visible in the delocalized case (a = 2 nm) for which the classical transport energy is À0.55 eV, than in the localized case (a = 0.5 nm) for which the classical value equals À0.26 eV.
Due to this saturation effect, we find that eqn ( 7)-( 8) in combination with the transport energy values obtained from the simulated histograms do predict a linear dependence at low values of the Fermi level only.Nevertheless the simulation predicts an almost linear dependence at all regimes.To understand this we have to take into account that at high occupations a substantial amount of the upward hopping moves go to levels close to the conduction band level (see Fig. 2).This introduces a distortion in the average implicit to eqn (9) because no hops above E = 0 are allowed.To ascertain the magnitude of this distortion we have performed calculations with eqn ( 7)- (8) assuming that the transport level coincides with the conduction band level, i.e., E tr = 0.This calculation renders a linear dependence in the full density range.The agreement with the simulation data is good at high Fermi levels (where upwards hopping moves are controlled by the upper limit of E = 0) but poor at low Fermi levels, where transport is controlled by jumps to the transport energy level.
The results shown in Fig. 4 indicate that the real transport energy should lie between the classical value of eqn (6) and the conduction band level E = 0.The values of E 0 max obtained from our corrected histograms are close but not the same as E tr .To obtain this we should distinguish moves that contribute effectively to transport from those that do not.This calculation would require to remove also the ''second-order'' moves discussed in subsection IIID.
F Temperature dependence of the diffusion coefficient RWNS calculations were performed to obtain the effect of ambient temperature on the diffusion coefficient.Arrhenius plots for these calculations are shown in Fig. 7 and 8 in the temperature range 260-340 K. Nearly linear plots are obtained, with an activation energy that is larger for deeper Fermi levels, as it could be expected.The Arrhenius behavior is characteristic of the multiple-trapping transport. 31,48This is an indication, as discussed above, that at a fixed Fermi level, there is a well-defined transport energy that makes transport to occur effectively via thermal activation to a transport level.A similar result has been obtained by Vissenberg and Matters using percolation theory. 49 must be noted that the theoretical framework contained in eqn ( 7)-( 9) is shown to predict a quasi Arrhenius behavior as well.This is due to the fact that the temperature dependences of the prefactors and the transport energy are much weaker than the energetic exponential factor.Furthermore, the transport energy is either a constant (at low occupations) or it moves towards higher values (at high occupations).In both cases an Arrhenius behavior with respect to temperature is expected.
The Arrhenius behavior is maintained if the characteristic temperature of the distribution is lower.Another important feature is that the activation energy is smaller for the delocalized case.This indicates that carrier percolation becomes facilitated when the range of the mean jump is larger, so that sites of similar energies are available for carriers.

E Implications in DSC functioning
In this work we want to make a connection with the relevant regime in DSC and related devices.It is known that at 1 sun illumination the electron density inside the semiconductor oxide is approximately equal to 10 17 cm À3 = 10 À4 nm À3 (1 electron per nanoparticle 50 ).For a characteristic temperature of T 0 = 600-800 K and a trap density of 10 21 cm À3 , which are realistic values 36,41 for nanocrystalline TiO 2 , this density corresponds to Fermi energies below À0.60 eV.As it can be observed in Fig. 3 and 6, this value corresponds to the regime for which the effective transport energy converges with the classical value given by eqn (6).Hence the predicted behavior for the diffusion coefficient is close to that yielded by the approximate formula (9) and thus indistinguishable from that predicted by the multiple-trapping model.
Furthermore, Arrhenius behavior with typical activation energies of 0.10-0.15eV are commonly found in the experiments 48 for nanocrystalline TiO 2 .Best agreement with the simulation data is found for T 0 = 800 K and a = 0.5 nm.Again the fact that there exist a well-defined transport level in the Fig. 6 Energy of the most probable jump (triangles), E max , and estimation of the effective transport energy, E 0 max , (circles) as a function of Fermi level.The first are extracted from the maxima of the energy histograms whereas the latter are extracted from the maxima of the ''corrected'' histograms with backward jumps between pair of sites removed.The horizontal line represents the classical value predicted by eqn (6).The parameters used were T = 275 K, T 0 = 800 K and a = 0.5 nm.

IV. Conclusions
The RWNS method with Miller-Abrahams hopping rates and exponential distribution of energies on a random network of traps has been utilized to test the transport characteristics in random media and to obtain the jump diffusion coefficient versus Fermi level and temperature.An approximate exponential dependence is found for the former and Arrhenius behavior for the latter.
The simulation helps to distinguish between the energy of the most probable jump and an estimation of the effective transport energy that determines the transport properties of the system.This latter value is found to move upward as the carrier density is increased except at low occupations where it converges with the classical value predicted in the literature.We found that in numerical modeling aiming to detect the transport energy at high densities, it is essential to remove from the computation back-and-forth jumps between near sites, otherwise the more probable target site displays a large distortion with respect to the sites contributing to diffusive transport.
Comparison of the present results with the conditions of interest in the functioning of photovoltaic devices based on nanocrystalline TiO 2 reveal that in this case the effective transport energy is approximately independent of the Fermi level.Hence the observed behavior is similar to that found with the multiple-trapping model and demonstrates that a hopping mechanism can also explain the experimental behaviour of the diffusion coefficient.

Fig. 1
Fig. 1 Illustration of the random walk method employed in this work.Traps (open circles) are distributed on a simulation box of size a L .Some of these traps are occupied by charge carriers (grey circles).For a certain carrier (black circle), hopping times to neighbored traps are computed according to eqn (1).This computation is restricted to those traps within the cut-off r cut that are not occupied.Once all release times are computed the minimum hopping time and its corresponding target trap are looked for and stored.The process is repeated for all carriers in the simulation box so that the carrier having the minimum hopping time is identified.This carrier is then moved to its corresponding target trap.See text for more details.

Fig. 3
Fig. 3 Energy of the most probable jump versus Fermi level (upper panel) and carrier density (lower panel) as obtained from RWNS calculations with Miller-Abrahams hopping rates with a = 0.5 nm.Results shown correspond to T 0 = 600 K (circles) and T 0 = 800 K (triangles).The dashed and dotted lines represent the classical values as obtained from eqn (6).

Fig. 4
Fig. 4 Jump diffusion coefficient vs. Fermi level as obtained from RWNS calculations with Miller-Abrahams hopping rates (full circles) and several theoretical predictions (see text for details): eqn (9) (solid line), eqn (8) with E tr taken from the classical value of eqn (6) (open circles), eqn (8) with E tr = E 0 max (squares), eqn (8) with E tr = 0 (times).The dashed line is a linear fit of the simulation data.Results shown correspond to T = 275 K and T 0 = 800 K and localization radii of a = 2 nm and 0.5 nm.

Fig. 5
Fig. 5 Histograms of the energies of the target sites N(E), (squares) and the same without considering backward jumps between pair of sites N 0 (E) (triangles, see text for details).Results for simulations at a finite carrier density (10/15 3 ) (upper panel) and for a single carrier (lower panel) are shown.The parameters used were T = 275 K, T 0 = 800 K and a = 0.5 nm.

Fig. 7
Fig. 7 Jump diffusion coefficient vs. of ambient temperature as obtained from one-particle RWNS calculations with Miller-Abrahams hopping rates at T 0 = 800 K and a = 0.5 nm.Results shown correspond to E F = À0.3 eV (circles) and E F = À0.4 eV (squares).The activation energies derived from both set of data are 0.15 and 0.24 eV, respectively.