Materials selection rules for optimum power factor in two-dimensional thermoelectrics

Increase in demand for energy sources inspired a vigorous search for new and efficient technologies that can maximize the utilization of existing resources. Recovering waste heat to generate electricity using thermoelectric (TE) materials is a viable path to increase efficiency in existing power generation methods, power wearable sensors, and even supply some of the energy needs in buildings [1]. However, doings so will require materials with high TE power factor that will maximize TE energy conversion efficiency. The efficiency of TE conversion is measured in terms of the dimensionless figure-of-merit [2] ZT = S²σ T/(κ_e + κ_p), where σ is electrical conductivity, S is Seebeck coefficient, κ_e is electronic thermal conductivity, and κ_p is lattice (phonon) thermal conductivity. Controlled engineering of materials by introducing grain boundaries [3, 4], surface roughness [5–7], interface scattering [8–10], and alloy mass disorder [11–13] reduces κ_p without significantly changing the electronic transport properties. For maximum enhancement, all of these techniques can be combined in hierarchical nanostructuring [14, 15]. Many TE materials studied thus far rely on some variant of this strategy for reducing κ_p to improve ZT.

On the other hand, the product S²σ is referred to as power factor and it is mostly dependent on the electronic band structure of material; increasing the power factor or discovering new materials with intrinsically high power factors also improves ZT. To achieve that, a broad swath of approaches, ranging from band-engineering techniques like moving to low-dimensional materials [16–18], to modulation doping [19, 20], electron filtering [21], and electrostatic gating [22] have been proposed. One of the factors that limits our ability to improve power factor by band engineering is the strong coupling between S and σ. Especially in semiconductors, they have opposite behavior with respect to the carrier density and cannot be tuned independently: $S\propto \langle E-{E}_{F}\rangle$ so increasing carrier concentration to boost σ also brings the Fermi level closer to the band edge, lowering S. Hence, careful selection of materials with desired properties is required in order to discover new materials with an intrinsically higher power factor. Once such materials are discovered, their ZT could be further improved via aforementioned extrinsic techniques to lower the lattice thermal conductivity; materials with both a large power factor and a large lattice thermal conductivity can be often be further engineered with minor impact on the electronic properties. Towards this end of identifying new materials with high power factors, previous studies proposed simple material property criteria to narrow the search time for identifying the bulk materials or compounds. Among them, the widely accepted criteria for bulk materials are low effective mass [23, 24], a band-gap ≈10 k_BT [24–26], high band degeneracy [24, 27], and a carefully tuned composition [28] or a complex unit cell [29]. Mahan and Sofo derived that a delta-function TDF is ideal [17]. Subsequently, Zhou et al showed that, for low dimensional materials, a window-like TDF with a narrow transport bandwidth [30] ≈5 k_BT is favorable; however, the optimal bandwidth [31] and bandgap [32] depend on how scattering processes were modeled.

More sophisticated material descriptors or parameters have been proposed to combine different material properties that can be used to evaluate a material's TE capabilities. One of the first is the B-factor introduced by Chasmar and Stratton [33], that can be used as a material descriptor independent of the carrier density (therefore independent of the Fermi energy). B-factor, sometimes also denoted by β, is a material dependent parameter that combines the intrinsic electronic and phononic properties, including mobility, thermal conductivity, and effective mass, into a single parameter that can be used to evaluate a material's TE capability. It is related to ZT via ZT = uβ/(vβ + 1), where u and v depend on the carrier concentration of the material, which is normally optimized via doping to maximize ZT. A material with a larger β and optimal carrier density was shown to provide higher ZT. Subsequently, to understand the anisotropy contribution to TE performance, the b-factor [34, 35] was defined by expanding the B-factor (β) to include any form of band structure and band degeneracy. The b-factor is defined as $b={{BF}}_{1/2}\left[({E}_{F}-{E}_{C})/({k}_{B}T)\right]$, where F_1/2 is a Fermi–Dirac integral of order 1/2, E_F is the Fermi energy and E_C is the conduction band edge. The higher band degeneracy requirement for high performing TEs [24, 28] is often proved using a simple constant relaxation time approximation or constant mean-free-path approximation. First principles band structure calculations, along with comprehensive electron–phonon scattering rates, improved our understanding of the effect complex band structures and band engineering on TE properties [36, 37]. In a material with multiple anisotropic valleys, the relative strength of intra- versus inter-valley electron–phonon scattering was shown to control the TE performance more than the degeneracy alone [34]. Parameters such as the Fermi surface complexity factor was then proposed to quantify the effect of anisotropy and band degeneracy on TE performance in complex band structure materials [38].

In an effort to formulate a metric that can help to identify materials with single cut-off energy, as originally proposed by Mahan [39], the shape factor (SF) was proposed as a measure of the asymmetry in the density of states (DOS) [40]. SF, for which we give the explicit formula later in the text, is closest to unity when the DOS has abrupt changes and lower when the DOS is smooth and gradual, which follows the conventional wisdom that abrupt changes near the Fermi level lead to better TE properties. A high throughput search among the known compounds confirmed that the SF should approach unity to maximize ZT [40]. The required shape for the DOS to achieve it resembles a step function, which is also the shape of the DOS in 2D materials with a parabolic bandstructure. With the advent of new 2D materials like graphene and single-layer transition metal dichalcogenides (TMDCs), the TE efficiency of these materials has garnered increasing interest. Even though reports suggest an improvement in the power factor of these materials [41–43], with some of the highest power factors ever reported being in graphene encapsulated with hexagonal boron nitride (hBN) [44], high lattice thermal conductivity restricts materials like graphene from achieving higher ZT. Apart from graphene, hBN, and TMDCs, there are many other 2D materials that have exceptional electronic and thermal properties, with the total number of potentially stable 2D van der Waals materials approaching 2000, and a few exhibit good power factor as well. But it is neither obvious which among them are good starting points for TE studies nor is it feasible to synthesize and test every 2D material. We need more insight into the relation between material parameters and TE efficiency, which can help in formulating the material selection rules for new efficient TE materials. Using publicly available material databases, studies [45–47] have predicted numerous new materials with weak van der Waals forces that can be exfoliated into single-layer forms With this extensive and growing knowledge, there is a need to explore the possibility of their use in future TEs.

In our recent study [18], we simulated transport in the inversion layer of a gated silicon-on-insulator device, where the DOS was two-dimensional due to confinement of carriers. We found that inelastic scattering from optical phonons and transitions between subbands are the dominant mechanisms and operating the device in a regime where both equal 5 k_BT significantly enhances the TE power factor. There is a need for further research aimed at formulating more complete material descriptors in 2D materials to identify potential candidates with high power factors. The main aim of this study is to use the distinct features of 2D materials and formulate simple yet powerful material selection rules that can provide guidance and identify candidates with high intrinsic power factor for further computational study, synthesis, and measurements. We focus on those material parameters, such as effective mass and deformation potentials, that are easily available in existing 2D materials databases or that require relatively computationally inexpensive first principles calculations in order to reduce the time required for initial evaluation and screening of materials for higher TE performance. Our focus is on semiconducting parabolic-band 2D materials with a step-like DOS and strong inelastic and intravalley scattering. We identify combinations of materials properties that lead to improvements in the power factor by restricting the energy 'bandwidth' of those electrons that contribute to transport but do not have enough energy to participate in inelastic and intervalley scattering. Our approach is one avenue through the vast field of new 2D materials toward more efficient TEs with higher intrinsic power factors and does not eliminate other complementary avenues via extrinsic modification, such as reducing lattice thermal conductivity or forming vertical and lateral heterostructures [48].

Improving the power factor requires identifying the material parameters that simultaneously effect the Seebeck coefficient and electrical conductivity. Starting from the Boltzmann transport equation under the relaxation time approximation and assuming a small deviation from equilibrium, S and σ are be related to the transport integrals as [22]

$$\begin{eqnarray}&&\sigma ={e}^{2}{L}^{(0)},\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&S=\left(-{L}^{\left(1\right)}/{L}^{\left(0\right)}\right)/{eT},\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&{L}^{(\alpha )}=\sum _{n}\mathop{\mathop{\int }\limits^{\infty }}\limits_{0}\left(-\partial {f}_{0}/\partial E\right){(E+{E}_{n}-{E}_{F})}^{\alpha }\sigma (E){dE},\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&{\sigma }_{n}(E)={\tau }_{n}\left(E\right){D}_{n}\left(E\right){v}_{n}^{2}(E),\end{eqnarray} \tag{ 4 }$$

where L^(α) is the transport integral and σ_n(E) is the transport distribution function (TDF), sometimes referred to as the differential conductivity; furthermore, E_F is the Fermi energy, E_n is the energy of the bottom of valley n, D_n(E) is DOS of that valley, τ_n(E) is the energy-dependent relaxation time, and v_n(E) is the group velocity of the carriers. The Fermi window ∂f₀/∂E is symmetric about E_F and has a spread of about 10 k_BT. The relaxation time of the carriers depends on dominant scattering mechanisms in the material and has varying dependency on energy (E). Assuming that the scattering completely randomizes the transport, the relaxation time is calculated by taking an inverse of the sum of all the scattering rates.

In the effective mass approximation, both carrier velocity and the DOS can be obtained from the effective mass of a particular (main or satellite) valley and the carrier energy. Now observing the transport integral, it is apparent that TDF has to be maximized to improve σ. On the other hand, it was shown that a δ-shaped [17] or narrow-bandwidth [30] TDF is needed for maximum achievable S. But in realistic conditions, a delta-shaped or even narrow TDF is not feasible in 3D materials due to the square-root dependence of both the relaxation time (τ_n) and DOS (D_n) on energy, stemming from the parabolic bands, which cancel each other's gains (as τ_n ∝ 1/D_n). In our previous work, we proposed an alternative approach [18] using an extrinsic 2D material, in which quantum confinement in the material plane causes a step-like DOS with a separation of E_cs between subbands. We further found that this effect is enhanced when inelastic scattering is the dominant mechanism. The delay in the onset of inelastic scattering, where ${{\rm{\Gamma }}}^{{inel}}\propto {D}_{n}(E+{E}_{{ph}})$, creates a window of energies from the bottom of the conduction band up to E_ph where only elastic scattering occurs. The sharp onset of inelastic processes above E_ph means that scattering, and with it the TDF is much higher for carriers below E_cs = E_ph and lower for energies above this cut-off, which results in a narrow TDF function with a limited energy bandwidth (E < E_cs = E_ph) over which most of the transport occurs. This, in addition to ${v}_{n}^{2}\propto {E}$ in the parabolic bands, gives a narrow window-like shape to the TDF, where most of the transport occurs in a restricted range of carrier energies, as illustrated schematically in figure 1. To simultaneously improve S and σ, the Fermi window should coincide with the TDF. Maximizing this overlap of the half-window width of 5 k_BT around E_F requires TDF to be 5 k_BT wide and ${E}_{F}\cong 0$. Using this approach, we showed a peak enhancement [18] in S at low temperatures in an inversion layer created by confinement that has a step-like DOS.

Zoom In Zoom Out Reset image size

Here we extend our approach to intrinsic 2D materials and formulate selection rules that can be used to shortlist candidate 2D materials with potentially higher power factors. As discussed above, the important components to be optimized are TDF (which include ${\tau }_{n}$(E), D_n(E), and v_n(E)) and the Fermi window. Table 1 summarizes the contributing factors to TDF for single valley in two different band profiles possible in 2D regime: a Dirac material (graphene) and a semiconducting material (such as MoS₂). The electronic band structure of a Dirac material can be defined using massless Dirac equation with a Fermi velocity of v_F, while for the semiconducting 2D material, we use the effective mass approximation (effective mass of ${{m}}_{n}^{* }$) with elliptical valleys in the conduction band. This results in massless DOS for Dirac material, and an effective-mass-dependent DOS for semiconducting material with a band degeneracy of ${{g}}_{v}^{n}$. Also, ${v}_{n}^{2}$(E) is constant for Dirac materials and linearly dependent on energy (${v}_{n}^{2}(E)\propto E$) for the semiconducting material with parabolic bands.

Table 1.  Comparison of important material parameters between a Dirac material and a semiconducting 2D material.

The scattering of carriers especially from phonons determines the relaxation time (τ_n) used to calculate TDF. Carrier scattering from phonons can occur in two forms: elastic and inelastic scattering. Elastic scattering, where there is little or no change in the energy of carriers, is generally from acoustic phonons due to their low energy. The elastic scattering rate (Γ^el) is proportional to the density of final states, which depends on D_n(E) of the band, and a pre-factor Γ₀, which is proportional to deformation potentials (Ξ₀). On the other hand, inelastic scattering (Γ^inel) involves emission or absorption of phonons by carriers, resulting in a change of their energy. Controlling the relative strength of these two processes by targeting materials with a small acoustic and large optical phonon deformation potentials helps in achieving a narrow TDF. Δ₁ (for emission) and Δ₂ (for absorption) are used to make sure that the energy conservation is achieved in inelastic scattering. In this work, we seek to identify the properties of 2D materials that are intrinsically good candidates for TE applications. We include acoustic phonons as the source of quasi-elastic scattering; however, in realistic 2D materials, the presence of charged impurities can be another source of elastic scattering. Its impact on TE performance is detrimental as it reduces mobility without enhancing the Seebeck coefficient, somewhat analogous to our acoustic phonons. Ongoing developments of fabrication techniques, better substrates, and encapsulation by hexagonal boron nitride (hBN) [44] have been shown to reduce the effect of impurities and could be applied to future 2D materials as well. As shown in table 1, inelastic scattering in 2D materials can occur in the form of intervalley and intravalley transitions. Intervalley transitions require carriers from one valley (with DOS of D_i(E)) to other valley (with DOS of D_f(E) and Z represents number of available valleys).

In Dirac materials, inelastic scattering from both intravalley and intervalley is proportional to the ratio of initial and final D_n(E), while in a semiconducting material, intravalley scattering is independent of D_n(E), but intervalley scattering has the same dependence as a Dirac material. The linear DOS and a constant v_n(E) in Dirac materials result in a linear TDF. Graphene, when encapsulated by hBN to improve mobility, has been shown to have a high TE power factor [44], about four times higher than MoS₂. While both Dirac and semiconducting 2D materials can exhibit good TE performance, our focus here is on the further enhancements that could be observed in 2D materials in which inelastic phonon scattering is strong and leads to a restricted transport bandwidth. Inelastic scattering, primarily from electron–electron interactions, has been shown to play a role in enhancing the thermopower of graphene [49] beyond values expected based on the Mott formula, valid for transport dominated by elastic processes. The bipolar nature of the transport near the Dirac point, stemming from the absence of a band-gap, somewhat limits their Seebeck coefficient; the record-high power factor of graphene stems largely from high carrier mobility after encapsulation by hBN to mitigate impurity scattering. Since Dirac materials lack a bandgap and typically do not have satellite valleys near the Dirac point, optical phonon emission is always possible while intervalley scattering rarely occurs. Consequently, the type of enhancements that we study in this work are less effective in identifying Dirac materials with even higher power factors. Hence, the rest of our analysis is focused on semiconducting 2D materials with parabolic bands.

Consolidating the variables to be analyzed, we divide them into two categories: parameters that control (a) the width and (b) the height of the TDF. The width of the TDF depends on the onset of inelastic scattering: optical phonon emission is significantly more likely than absorption by a ratio $({N}_{{ph}}+1)/{N}_{{ph}}\approx \exp ({E}_{{ph}}/{kT})$, where N_ph is the Bose–Einstein distribution function of phonons. An electron must have at least E > E_ph kinetic energy to emit the phonon, the onset of emission causes a sharp drop in the TDF which can be tuned by phonon energy. In addition to step-like DOS and E_ph, the separation between main and satellite valleys (E_cs) also determines the width of TDF through scattering into the satellite valley, as shown schematically in figure 1. Along with controlling the width of TDF, E_ph and E_cs help to maximize the Fermi window utilization. D_n(E), which depends on the effective mass ${{m}}_{n}^{* }$ and band degeneracy ${g}_{v}^{n}$, controls the height of TDF. The impact of intervalley scattering of carriers depends on the ratio of initial and final state DOS (D_i/D_f). Here, to effectively understand the effect of D_i/D_f the total DOS ($\sum _{n}$ D_n) is kept constant using m_tot = $\sum _{n}$ ${{m}}_{n}^{* }$. Now using these parameters to tune the power factor, our model uses MoS₂ [50] as prototype material, where E_cs = 70 meV ≈ 2.7 k_BT, E_ph = 47 meV ≈ 1.8 k_BT, ${{g}}_{v}^{1}$ = 6, ${{g}}_{v}^{2}$ = 6, m_tot = 1.4, D₁/D₂ = 0.8, and isotropic circular valleys are assumed for conduction band. We chose MoS₂ far our starting point as it has been shown to have a high TE power factor [51] and vary the material parameters and find the optimum combination of properties for a 2D material that can provide further increases in the power factor at room temperature (T = 300 K).

To obtain the window-shaped TDF by delaying the onset of scattering in a material, inelastic scattering should be dominant over the elastic processes. In general, we would seek materials with small acoustic phonon deformation potentials; here we achieve this favorable feature indirectly by reducing the elastic scattering by 50% through lowering the Γ₀. This results in an improvement in the power factor, as shown in figure 2(a), mainly due to the increase in σ. We plot our resulting Seebeck coefficients S on the vertical axis and conductivity σ on the horizontal axis, while varying the carrier concentration from 10¹¹ to 3 × 10¹³ cm⁻². Doing so emphasizes the fundamental trade-off between S and σ and allows us to visualize the peak power factor: further toward the top-right corner the S versus σ curve extends, higher the power factor. Next, the optical phonon energy is tuned to shift the onset of inelastic scattering to create a window-shaped TDF and to maximize the Fermi window utilization. By increasing the energy (E_ph) of phonons that participate in inelastic scattering a maximum in power factor is observed when E_ph = 5 k_BT. Any increase in E_ph further resulted in little or no improvement in power factor. This is a consequence of better overlap of window-shaped TDF and Fermi window (which has a half width of 5 k_BT). Once a window shape for TDF is achieved in the material with E_ph = 5 k_BT, the DOS is tuned to further improve the height of TDF. By keeping the m_tot constant, the ratio of D₁/D₂ is varied. This is possible due to the linear dependence of DOS on effective mass in 2D materials. As seen in figure 2(b), a lower D₁/D₂ results in higher power factor. This shows that the higher separation between D₁ and D₂ is desirable to maximize power factor. Next, while keeping D₁/D₂ constant, the m_tot is varied and we found that lower m_tot gives higher power factor. Lower effective mass reflects higher velocity of carriers that improves the height of TDF function (TDF $\propto {v}^{2}$). This requirement has already been established in 3D materials [27]; we show here that it applies to 2D materials as well.

Zoom In Zoom Out Reset image size

Our approach tries to create a rectangular shape for the product of relaxation time (τ_n) and DOS (D_n) using the step-like DOS. The valley splitting (E_cs) controls the width of first step, which in turn decides the width of TDF. We already know that the ideal shift in the onset of scattering occurs when E_ph = 5 k_BT. Now by further matching this with the width of the first step (when E_cs = E_ph), transport can be effectively restricted to first step or valley. By calculating the power factor in a 2D material in which inelastic scattering is dominant with D₁/D₂ = 0.2 and m_tot = 0.5 when varying the valley splitting, a maximum is observed when E_cs = E_ph = 5 k_BT as shown in figure 2(c). Further increase in E_cs beyond 5 k_BT while limiting E_ph to 5 k_BT, does not help since the width is limited by E_ph. Till now, we studied the effect of different material parameters on power factor assuming the degeneracy of the valleys to be 1 (circular valleys). Now to study the effect of valley degeneracy, ${{g}}_{v}^{n}$ is varied but D₁/D₂ is still maintained at 0.2 by adjusting the ${{m}}_{n}^{* }$. As shown in figure 2(d), high valley degeneracy for D₁ (the valley that contributes to transport the most) improves the power factor. However, a higher degeneracy for the secondary valley deteriorates σ, thereby reducing the power factor. This is due to the availability of more states in the second (satellite) valley to scatter into from the first valley. A related strategy of increasing band degeneracy by valley alignment was studied to demonstrate an improvement in the power factor of half-Heusler alloys, showing that including intervalley processes is critical in obtaining accurate results [37]. While our approach also relies on inelastic scattering, we observe largest improvements when both intervalley and optical phonon emission are aligned (E_ph = E_cs).

Figure 3(a) shows the step-by-step contribution towards Seebeck coefficient and electrical conductivity from each material selection criteria proposed here that can be achieved or maximum achievable power factor enhancement. The power factor gains seen by introducing these selection criteria and the influence of carrier concentration is shown in figure 3(b). The dominance of inelastic processes along with the optimal optical phonon energies provide at least ten-fold increase in power factor. A low D₁/D₂ along with low effective mass gives a further six-fold increase in power factor. Implementing the rest of the proposed selection criteria can boost the overall power factor by two orders in magnitude. Finally, to understand the correlation between TDF and the power factor achieved, TDF is plotted in figure 3(c) for each of the identified selection criteria. The progression in figure 3(c) shows that implementing these selection criteria results in narrowing transport bandwidth, with increasing proportion of transport contributed by states between the bottom of the conduction band and 5 k_BT (0.13 eV at room temperature), corresponding to increases in power factor. To compare the material selection criteria with already defined metrics, SF is calculated while introducing the features identified in this study (as shown in figure 3(d)). The value of the SF is calculated by finding a cut-off energy E* that maximized the ${\rm{SF}}=1-\left({\int }_{0}^{{E}^{* }}D(E){dE}\right)/{E}^{* }D({E}^{* });$ consequently, the SF is independent of carrier concentration. Comparing the SF values (y-axis in figure 3(d)) with our calculated power factors (x-axis), shows the inability of SF to predict the power factor enhancement seen here. In 2D systems, the SF can be approximately reduced to D₂/(D₁ + D₂) due to their step-like DOS. Since SF is only controlled by the asymmetry in the DOS, it is not able to completely capture the changes that occur in the TDF, primarily those related to inelastic scattering. Instead, SF is only affected by D₁/D₂ and increases when D₁/D₂ is reduced. Using only SF as a selection criterion will result in a broad spectrum of 2D materials due to their intrinsic step-like DOS. Hence there is a need to look beyond SF for the features like the ones identified in this study, that can be used as material selection rules for 2D materials for better power factor prediction.

Zoom In Zoom Out Reset image size

Using the simple transport model as described here along with the recent development of 2D databases [47] that are comprehensive, identifying higher power factor materials could be done at a faster pace, especially if materials databases add transport properties such as deformation potentials, which can be obtained from first principles calculations. To corroborate the trends we predict, we contrast two popular TMDCs, MoS₂ and WSe₂. Experimental observations showed a larger power factor MoS₂ (around 8.5 mW m⁻¹ K⁻²) [51] than in WSe₂ (3.2–3.7 mW m⁻¹ K⁻²) [52]. Using the selection criteria proposed here and properties calculated from first principles, we can predict this decrease in WSe₂, even though the two have similar electronic bandstructure [53]. With small valley separation of 16 meV in WSe₂ versus 70 meV in MoS₂, coupled with smaller optical phonon energies (28.5 versus 47 meV) hurts the power factor despite low elastic scattering (due to lower acoustic deformation potential) and smaller effective mass of WSe₂ [53]. While our model produces higher power factors (44 and 6 mW m⁻¹ K⁻² for MoS₂ and WSe₂, respectively) than what was measured, mainly due the presence of impurities in the experiments, the relative ordering is reproduced correctly.

In conclusion, the search for new highly efficient TE materials requires that we identify candidates that have a high power factor. Here, a novel approach is proposed to ease the process of this material search by identifying those materials properties that ensure high power factor. A simple numerical model based on the widely-used effective mass approximation is implemented; the model captures inelastic scattering in 2D semiconducting materials with step-like DOS. We find that inelastic processes delay the onset of scattering, thereby maximizing the utilization of the Fermi window for transport. Our simulations show that materials where inelastic processes are dominant (or elastic processes are suppressed), specifically when optical phonon energy E_ph equals 5 k_BT, improve the power factor. A further enhancement is possible if the separation between the DOS of valleys is large (lower D₁/D₂), the effective mass is low, and the conduction-band valley that participate in transport have high degeneracy. By introducing this combination of material features in our prototype material, based initially on MoS₂ for comparison, the power factor is boosted from 0.044 Wm⁻¹ K⁻² in MoS₂ to 22 Wm⁻¹ K⁻² in a hypothetical material that has all of the aforementioned features, which is an increase of more than two orders in magnitude. We believe the material selection rules formulated here can help researchers to reduce the time and resources in the process of predicting and preselecting 2D materials that have high power factor.