Non-Wiedemann-Franz Behavior of the Thermal Conductivity of Organic Semiconductors

Organic semiconductors have attracted increasing interest as thermoelectric converters in recent years due to their intrinsically low thermal conductivity compared to inorganic materials. This boom has led to encouraging practical results, in which the thermal conductivity has predominantly been treated as an empirical number. However, in an optimized thermoelectric material, the electronic component can dominate the thermal conductivity in which case the figure of merit $ZT$ becomes a function of thermopower and Lorentz factor only. Hence design of effective organic thermoelectric materials requires understanding the Lorenz number. Here, analytical modeling and kinetic Monte Carlo simulations are combined to study the effect of energetic disorder and length scales on the correlation of electrical and thermal conductivity in organic semiconductor thermoelectrics. We show that a Lorenz factor up to a factor $\sim 5$ below the Sommerfeld value can be obtained for weakly disordered systems, in contrast with what has been observed for materials with band transport. Although the electronic contribution dominates the thermal conductivity within the application-relevant parameter space, reaching ZT>1 would require to minimize both the energetic disorder but also the lattice thermal conductivity to values below $\kappa_\text{lat}<0.2$W/mK.


Introduction
Semiconducting organic materials have attracted increasing interest as thermoelectric (TE) converters in the recent years due to their potentially low material and fabrication costs 1 and non-toxicity. 2,3 Another key benefit of organic semiconductors compared to typical inorganic materials is their intrinsically low thermal conductivity that is due to structural disorder and weak intermolecular coupling. The latter can remarkably increase the thermoelectric conversion efficiency, governed by the dimensionless figure of merit where S is the Seebeck coefficient, σ is the electrical conductivity, T is the absolute temperature, and is the total thermal conductivity, composed of electronic (κ el ) and lattice (κ lat ) contributions.
For inorganic materials, it is well known that the parameters S, σ and κ are coupled, so that optimizing one parameter tends to compromise the others. This makes the optimization of the TE efficiency an extremely challenging task, requiring strategies such as engineering of the electronic band structure and lattice thermal conductivity.
Such correlations are less established for organic semiconductors. In the recent years, lots of effort in the field of organic TEs went into studying the relation between σ and S, both experimentally and theoretically. [3][4][5][6] However, little is known about the correlation between electrical and thermal conductivity in these materials. In inorganic materials, the electronic contribution of the thermal conductivity is related to the electrical conductivity by the Wiedemann-Franz law, κ el /σ = LT , with a Lorenz number L that is typically close to the Sommerfeld value L 0 = π 2 3 k B q 2 . 7 In contrast, for organic materials there are indications that the Wiedemann-Franz law does not hold, and a large range of Lorenz numbers has been reported. While some experimental studies suggest L to be equal to the Sommerfeld value for a free electron gas, 8,9 others show large deviations from this value in both directions. [10][11][12] This experimental divergence is accompanied by a lack of theoretical understanding. Lu et al. 13 as well as Upadhyaya et al. 14 reported large deviations from the Sommerfeld value, which was mainly attributed to the effect of disorder but the results are contradictory. The low attention paid to the electronic contribution of the thermal conductivity and Lorenz factor is surprising. First, in the application relevant regime the electronic contribution can dominate the thermal conductivity. Second, in the limit when electronic thermal conductivity dominates κ lat κ el , the maximum ZT for any material is given by ZT = S 2 /L. Thus, reducing the Lorenz number will be crucial for maximum ZT .
Here, we combine kinetic Monte Carlo (kMC) simulations and analytical modeling to study the electronic thermal conductivity κ el of disordered organic semiconductors. All simulations are based on a hopping formalism that takes into account the specific shape of the density of states (DOS). In particular, we also investigate a Gaussian DOS modified by Coulomb trapping of mobile charge carriers by ionized dopants. We show that the presence of energetic disorder as well as the magnitude of the localization length can lead to significant derivations from the Wiedemann-Franz law, corresponding to effective Lorenz numbers that can either be smaller or larger than L 0 .

Theoretical Framework
The macroscopic observables σ, S and κ el can be determined from: 15 with the electronic thermal conductivity defined as Here f FD is the Fermi-Dirac distribution function, E F the Fermi level, σ (E) the energydependent differential conductivity and κ 0 the electronic thermal conductivity when the electrochemical potential gradient inside the sample is zero. is here used as an ansatz and its applicability is discussed below.

Numerical Model
The kMC model has been extensively described before. 18,19 In brief, the kMC simulations account for variable-range hopping on a random lattice with a mean intersite distance of a NN = N The thermal conductivity κ el was calculated numerically exact using the definitions in Eq. (3)-(6). Using σ (E) = j(E)/F , where F is the constant electric field and j(E) the current density at an energy E, we arrive at: To handle integrals of the form dE E x j(E) with x = [0 − 2], we use the fact that under steady-state conditions all (differential) currents are constant and integrate over time. The resulting double integrals reflect the total amount of charge that passes through a unit cross section in a given time t, weighted by E x . Numerically, they can be evaluated as sums over all hopping events i in the time t, where A and L are the cross section and length of the simulation box, respectively, and ∆z i the displacement of the ith hopping event in the direction of the electric field applied along z. The energy of a hop between sites i and j is calculated as (E i + E j )/2 The thermopower S is determined analogously, further details can be found in Ref. 19 Analytical Model Our analytical model extends the work of Schmechel 20 and Ihnatsenka et al., 21 that expresses the Seebeck coefficient for a hopping system in a differential conductivity σ (E) to also describe the heat conductivity of the hole or electron gas. For completeness, the key expressions are given in the following. To calculate the total conductivity σ, thermopower S and the electronic contribution to the thermal conductivity κ el , the resulting differential conductivity is inserted in expressions for the transport coefficients stemming from the Boltzmann transport equation (see Eq. (3)-(6)).
The differential escape rate distribution ν esc from a state with energy E 0 via a state with E 0 + W is given by the product of the thermal activation rate from the initial state at energy Here, ν 0 is the carrier attempt-to-hop frequency, k B is the Boltzmann constant, T is the temperature and α describes the decay length of the localized wavefunction. It was shown in Ref. 20 that Eq. (9) is equivalent to the common Miller-Abrahams expression if the total transfer rate from an initial state with energy E 0 to a specific target state at energy E target , separated by a distance R is calculated. The parameterR(E) is the mean tunneling distance for an electron at energy E to reach a target state at lower energy: where f FD is the Fermi-Dirac distribution function and g( ) the DOS distribution. Eq. (10) is commonly employed to identify the critical hop in the infinite percolating network, 22,23 while it is used here as a measure of the range of all hops starting at energy E. In the former case, the parameter B = 2.8 reflects the critical number of bonds on the percolating network.
The total escape rate can then be determined from the differential escape rate [Eq. (9)] via Moreover, the mean level at which the carrier is released from its initial state, can be calculated with the help of the differential escape rate distribution.
From the differential form of the generalized Einstein relation one obtains the energydependent carrier mobility, where η is a fitting constant 24  is determined by Eq. (10). As a side note, we would like to stress that in principle the expectation value for the squared hopping distance < R 2 > should be used instead of < R > 2 . 26 However, the difference is negligible and did not affect the results here.
With the energy-dependent differential conductivity

Density of States (DOS)
The energy distribution of the localized sites through which hopping takes place is typically assumed to be Gaussian in shape, where E i and σ DOS are the central energy and the width of the Gaussian DOS, respectively, and N i is the total site density. Ionized dopants are known to act as Coulomb traps and therefore lead to a perturbation of the DOS in the form of exponential tail states and broadening of the main DOS peak. Such modifications of the DOS shape can have an enormous impact on the thermoelectric properties, e.g. result in a qualitative change of the S vs. σ curve. 5, 27 Arkhipov et al., 28 developed an approximation for the ion-perturbed DOS of a doped semiconductor, is the Coulomb energy, and ε r is the relative dielectric constant of the semiconductor. Here, we use an extension to Eq. (15) developed by Zuo et al., 19,29 to account for energy level differences between the dopant and the semiconductor ∆E = E d − E i , with E d being the relevant energy level of the dopant:

Results
In Fig. 1a), the electronic contribution of the thermal conductivity calculated from the analytical model for various values of the localization length α is compared to the results from the kMC simulation. We find that both the absolute value as well as the concentration dependence of the kMC data are accurately reproduced by the analytical model. For increasing carrier concentration, the thermal conductivity κ el increases, in contrast to the trends observed by Lu et al. 13 but consistent with the general trends which could be expected from the Wiedemann-Franz law and are observed in experiments. 9,12 As shown in Fig. 1b), changing the localization length to larger values increases κ el and σ simultaneously, as it implies stronger wavefunction overlap between adjacent sites and is therefore favorable for the hopping process. However, the scaling is not entirely linear, which results in a conductivity, and with this also carrier concentration dependent Lorenz factor, as shown in Fig. 1c). In organic semiconductors, the disorder of the DOS is known to have a significant influence on the electrical conductivity 29 and it is therefore likely that it also impacts the thermal conductivity κ el . Fig. 2a) shows that the absolute value of κ el as well as its dependence on the carrier concentration is strongly influenced by the disorder parameter σ DOS . Moreover, the latter affects the slope of κ el vs. σ, leading to a Lorenz factor that strongly depends on the disorder. As shown in Fig. 2c), weakly disordered systems and high carrier concentrations would allow a Lorenz factor of that is substantially below the Sommerfeld value. This is a promising result for the development of organic TE materials, as the condition L L 0 occurs in the high-conductivity part of the parameter space that is most relevant for TE where the effect of disorder on the Lorenz factor in organic semiconductors is studied and with findings from the inorganic community regarding the trend of the Lorenz factor with the bandwidth of the charge carrier dispersion. Therefore, we suggest that Eqs. (5) and (6) provide a reasonable approximation of the electronic part of the thermal conductivity also in organic semiconductors. Having shown in the previous section that the analytical model can well describe the results of the kMC simulation for a purely Gaussian DOS, we used the former model to study the Lorenz factor and the figure of merit ZT as a function of the energetic disorder and the localization length. We focused here on the behavior of ZT at T = 300 K, as the performance around room temperature is probably the most relevant case for organic TE materials in contrast to the high-temperature regime in which traditional inorganic crystalline TE materials are applied. 32,33 Moreover, polaronic effects in terms of Marcus theory where not explicitly taken into account. Nevertheless, we have thus far been able to very well describe experimental data without explicitly accounting for polaronic effects, using the simpler Miller-Abrahams rates. 5 Similar conclusions regarding the minor differences between the two rates were also drawn by Cottaar et al. 34 and Mendels and Tessler. 35 Fig. 4a) shows that the localization length has no significant influence on L for low energetic disorder, but gains importance for σ DOS > 3k B T . Meanwhile, the energetic disorder impacts the Lorenz factor in all parameter combinations considered within this study.
Independent on the exact value of the lattice contribution to the thermal conductivity κ lat , a low energetic disorder is decisive for an optimized figure of merit ZT , as depicted in Fig. 4b) and c). However, whether a low or high localization length is favorable depends on κ lat , as the absolute value of κ el found with this parameter set is approximately two orders of magnitude lower than κ lat = 0.2 W/mK, which is typically assumed for organic TE materials to date. 36 In the limit κ lat κ el , ZT is proportional to the power factor P F = S 2 σ (see Supplemental Material, 25 Fig. S2a) and maximized for large localization length, as this leads to (exponentially) higher σ, while S is hardly affected (see Supplemental Material, 25 Fig. S1). In contrast, for κ lat κ el , ZT is given by by P F/κ el T = S 2 /L.
Thus, in this limit the thermopower is decisive as the decrease in S 2 with increasing α is much greater than the one in L (see Supplemental Material, 25 Fig. S2b). Note, that the lattice thermal conductivity κ lat and the localization length α could be coupled via the strength of the inter-and intra-molecular bond but are here considered as independent parameters. An important implication of Figs. 4b) and c) is that it appears unlikely that organic semiconductors will reach ZT values approaching unity while maintaining a lattice thermal conductivity of ∼ 0.2 W/mK. Although the latter is already a rather low value, achieving application-relevant ZT will require even lower values.
Chemical or electrochemical doping of organic semiconductors is an effective way to increase the electrical conductivity by increasing the charge carrier concentration and is increasingly used in applications. As discussed above, adding ionized dopants to the intrinsic semiconductor will change the shape of the DOS due to long-range Coulombic interactions.  Fig. S3a), this is not the case for the thermal conductivity, as shown in Fig. 5a) and the thermopower (see Supplemental Material, 25 Fig. S3b). However, the general trends, i.e. a Lorenz number of L < L 0 for σ DOS < 3k B T (see Fig. 5b) and c)) are reproduced by the analytical model, which results from a cancellation of errors in κ el and σ.
Note that a first reaction scheme was used for the kMC simulations in Fig. 5, meaning that rates and energies were only recalculated for the moving particle, which significantly reduces calculation time. This approximation is accurate up till relative doping levels of 10 −2 , as shown by Zuo et al. 19 Extrapolating the electrical conductivity in Fig to illustrate the influence of the electronic contribution to the thermal conductivity. Charge carriers are introduced by doping and parameters are set according to the standard parameter set but with ν 0 = 3 × 10 16 s −1 and κ lat = 0.2 W/mK. The inset shows a zoom-in to the simulated data and experimental data for (n-type) poly(nickel-ethylenetetrathiolate) (black circle, ref. 38 ) and for (p-type) PEDOT-Tos (black square, ref. 37 ).

Conclusion
The electrical contribution to the thermal conductivity of disordered organic semiconductors has been systematically analyzed by means of analytical modeling and numerical kinetic Monte Carlo simulations. We show that optimizing organic thermoelectrics with respect to thermal properties requires to consider their energetic disorder as well as the relevant length scales of the system, i.e. the localization length relative to the typical inter-site distance. Large derivations from the Wiedemann-Franz law can be observed that can be characterized by an effective Lorenz number L that can be larger as well as substantially smaller than the Sommerfeld value L 0 . In consequence, our results are in contradiction with the universal Lorenz number found for materials with band transport but can give an explanation for the large range of Lorenz numbers reported from experimental studies.
Minimizing energetic disorder can be a viable strategy to reduce the Lorenz factor and obtain higher thermoelectric performance. Moreover, we show that for most cases reported to date κ el < κ lat but the electronic contribution to the thermal conductivity will unavoidably dominate the thermal conductivity for higher electrical conductivities that are needed to reach application-relevant figures of merit. At the same time, reaching beyond ZT = 1 will likely require suppressing the lattice thermal conductivity to values below κ lat = 0.2 W/mK, a situation that is conceptually similar to that for inorganic thermoelectrics.