Doping induced enhanced density of states in bismuth telluride

Considerable efforts have been devoted in recent years to improve the thermoelectric efficiency. The figure of merit zT commonly used to characterize the performance of a given material or structure is defined as $S2σT/κ$⁠, where S is the Seebeck coefficient, σ the electronic conductivity, κ the thermal conductivity, and T the absolute temperature.¹ The strategy to maximize the efficiency can evidently take two paths, viz., optimize κ and/or the power factor term, $S2σ$⁠. The primary focus thus far has been on reducing the lattice contribution to κ through nanostructuring.^2,3 However, the reduction in κ is constrained both by the electronic thermal conductivity and the fact that the lower bound on the phonon mean free path is the interatomic distance. The alternative route via the power factor is similarly fraught with challenges. In particular, the Seebeck coefficient experiences degradation when enhancement in the conductivity is sought by increasing the carrier concentration as characterized by the Pisarenko relation.⁴ The reduction in S must be avoided since the power factor is proportional to its square.

The potential approaches to circumvent this limiting inter-dependency include band structure engineering through low-dimensional structures,^5,6 carrier pocket engineering,⁷ hot carrier filtering,⁸ and introduction of resonant levels.^9,10 Here, the investigation is focused on the improvement through resonant doping. The Mahan-Sofo theory illustrates that the Seebeck coefficient in a material can be maximized when the density of states (DOS) has a Dirac delta-like peak with the Fermi level aligned 2.4k_BT away from the peak and the contributions from the rest of the electronic bands are held at the minimum.¹¹ While it is impossible to have a Dirac delta function for the DOS, a close approximation is realizable by introducing resonant impurity states through doping of appropriate atoms. In addition, for the resonant level to be conductive as necessary, it is preferable for the impurity state to have an s- or a p-like character rather than a highly localized f-like state. The increase in the Seebeck coefficient due to the resonant doping can be quantified by the Mott expression as¹² 

The first term simply accounts for the increased DOS g(E), whereas the second refers to the increased energy dependence of the mobility [$μ(E)$] due to the resonant impurity scattering which is relevant only at low temperatures. Here, k_B is the Boltzmann constant, n the carrier density, and E_f the Fermi energy. It should be noted that the Mott relation is derived assuming that the transport equation is a smoothly varying function of energy near the Fermi level and breaks down when there are sharp features in the DOS. The predicted improvement in zT as a result of the enhanced local DOS has been experimentally demonstrated in PbTe with 2% doping of Tl.⁹ This observation naturally raises a similar prospect for Bi₂Te₃, a material with the highest zT among room-temperature bulk thermoelectrics¹³ and thus the focus of significant research interests.

In this paper, a comprehensive theoretical investigation is undertaken to identify the dopants in bulk Bi₂Te₃ that are likely to realize the resonant states and thus the desired increase in zT. Experimental studies with regard to resonant doping in Bi₂Te₃ have so far shown that doping with Sn increases S which was attributed to the formation of a resonant level in the valence band (VB).¹⁰ On the other hand, no theoretical reports are available in the literature. As part of this work, a fully first-principles approach is adopted to study the doping induced modifications to the DOS for different candidate dopant atoms. Following the initial screening, partial DOS and charge density analyses are carried out for those atoms with promising results. Transport properties are calculated to show the effect of the resonant states on the power factor for doped Bi₂Te₃.

The ab initio calculations are performed in the density functional theory (DFT) framework of the plane wave based QUANTUM ESPRESSO.¹⁴ Fully relativistic norm conserving pseudopotentials are used that account for the spin-orbit interaction. Spin polarization, on the other hand, is not considered. The adopted exchange-correlation function is based on the generalized gradient approximation as parametrized by Perdew, Burke, and Ernzerhof.¹⁵ The long-ranged van der Waals (vdW) interaction relevant to the material of interest (i.e., Bi₂Te₃) is accounted for by including the semiempirical DFT-D2 correction.¹⁶ The kinetic energy cutoff of 60 Ry along with the Methfessel-Paxton smearing of 0.01 Ry is employed. Doped Bi₂Te₃ calculations are carried out by constructing a 60-atom supercell and replacing one of the host atoms with the dopant. The supercell is fully relaxed until the total energy is minimized, and the total force on the atoms is less than 10⁻⁴Ry. The Brillouin zone is sampled with a (12 × 12 × 6) grid. While the cases considered correspond to relatively high doping concentrations (⁠$∼4$% and $∼2.7$% for Bi and Te replacement, respectively), a larger supercell is not examined due to the computational constraints. As Bi₂Te₃ has a layered crystal structure with each quintuple layer consisting of two Bi atoms and three Te atoms, three substitutional sites are possible (see Fig. 1), namely, Bi atom, ionic Te atom (sandwiched between two Bi atoms in the quintuple layer), and vdW Te atom (weakly bonded to the Te atom of the next quintuple). Formation energy for each possible combination is evaluated to determine the most energetically favorable substitution site.

The dopant atoms considered are In, Pb, As, I, Cs, Po, Na, Tl, and S. These candidates are selected based on their outer shell properties. More specifically, the In, Pb, As, I, Po, and Tl atoms possess either similar outer shell characteristics or atomic radii to Bi and Te. On the other hand, Cs and Na are expected to be the representatives of the atoms with the s-orbital outer shell. All three substitution sites are considered for each of these dopant atoms. Only the impurity states present within 500 meV of the band edge are of interest since those located outside that energy range are practically irrelevant to the enhancement of thermoelectric power factor. Of the screened dopant atoms, In, Po, and Na have resonant peaks in the energy range of interest and are selected for further analysis. Figure 2(a) provides the DOS obtained for the In atom substituting the Bi atom in Bi₂Te₃. As shown, a resonant peak is clearly visible approximately 180 meV away from the VB edge (set at E = 0 before doping). The case of In substituting either kind of Te atom does not exhibit any resonance-like features. Doping with Po also shows favorable results as illustrated in Fig. 2(b). More interestingly, this case introduces, through substituting the Bi atom, twin resonant peaks: one that is 200 meV below the VB edge and the other 80 meV above the conduction band (CB) edge. As with In, no enhancement in the DOS is observed when Po substitutes Te atoms. Finally, the third dopant atom of interest, Na, has a resonant peak 200 meV away from the VB edge [Fig. 2(c)]. This peak is formed when Na substitutes either kind of Te atom, ionic or vdW, but absent for Na substituting the Bi atom. Note that Bi₂Te₃ no longer exhibits a bandgap after doping with the resonant atoms in all three cases, while a value of 80 meV is obtained for the pristine case. The disappearance of the bandgap is attributed to the large doping density assumed in the calculation (see above). Once the dopant level is lowered sufficiently, the bulk bandgap is expected to be restored and the resonant peaks will become narrower and sharper.¹⁷ 

A partial DOS analysis is carried out to probe the nature or the origin of the resonant peaks which result from the interaction between the wavefunctions of the dopant atoms and those of the host atoms. At a resonant level, the impurity state coincides with the band energy of the host atoms, unlike a normal dopant state where the impurity energy falls within the bandgap. For the case of In, the examination reveals that the resonant peak in the VB can be attributed to the s orbital of the dopant [Fig. 2(d)]. A closer look further suggests sp³-like hybridization with the neighboring out-of-plane Te atoms, as can be seen from Fig. 2(g). Partial DOS for the Po atom is displayed in Fig. 2(e). While the VB resonant peak is composed mainly of the p_z orbital of Po, there is an equal contribution from both p_z and p_xy orbitals toward the CB resonant peak. Unlike In, none of the resonant states show hybridized characteristics; instead, they are purely p-like states. The analysis also reveals that the Po p-like impurity peaks simply cause an increase in the DOS in the neighboring Te atoms in the same energy range [Fig. 2(h)], which is associated with a p-like character likewise. An additional point to note is the presence of a very pronounced peak in the impurity p_xy orbital centered around 600 meV away from the CB edge [not fully shown in Fig. 2(e)]. This, however, does not make a significant change for its contribution nearly shadowed by the large DOS of the host atoms already present in the corresponding energy range [see also Fig. 2(b)].

As for Na, the origin of the peak differs when Na substitutes either an ionic Te atom or a vdW Te atom. In the case of ionic Te substitution, the VB resonant peak is entirely p-like composed of the Na p orbital and Bi p orbital as shown in Figs. 2(f) and 2(i). Here, the Na atom modifies the DOS in the neighboring Bi atoms, similar to the Po doping described above. Na has a couple of impurity peaks close to the CB edge arising from p_z and p_xy orbitals, which does not lead to a resonant state due to an overlap between the two peaks resulting in a broadened DOS instead. When Na substitutes a vdW Te atom, on the other hand, the enhanced DOS is due to hybridization between the Na s orbital and the Bi p_z orbital. Since this hybridization is weak, the resonant peak formed is smaller. Note that the inclusion of the spin-orbit interaction is essential for the accurate description of the resonant states. It is because both Bi and Te p orbitals (which are involved in the hybridization with impurity atoms) experience spin-orbit induced energy shift. However, a universal trend is not anticipated in the impact of this interaction on the properties of the resonant states with different dopant atoms.

To estimate the ease of doping as well as to determine the most stable configuration, the formation energy is calculated for the different substitution sites. Since resonant states are formed only for certain substitution sites, it is essential to determine the most favorable growth condition resulting in the required configuration. The formation energy can be obtained from the following expression:¹⁸ 

where “d” and “ud” denote the doped and undoped cases, respectively. In addition, $ΔnA$ refers to the number of atoms of species A removed or added to the system with chemical potential μ_A. The chemical potentials for In, Po, and Na atoms [denoted as “X” in Eq. (2)] are readily obtained from the bulk unit cells (with bct, sc, and bcc symmetry, respectively), as μ in each case corresponds to the total energy per atom of the elemental ground state in the DFT calculation. The values for Bi and Te depend on the growth conditions and must follow the constraint $2μBi+3μTe=μBi2Te3$⁠. For the Bi-rich condition, $μBi$ is set to its bulk value (R$3¯$m symmetry). Correspondingly, $μTe$ in the Te-rich case is determined likewise (⁠$P$3₁21 symmetry). The calculated formation energies for different substitution sites are shown in Table I. When a dopant atom substitutes a Bi atom, it is found that E_fr is lower for the Te-rich condition than for the Bi-rich counterpart. The opposite also appears to hold for the case of dopant atom substituting the Te atom. Further, E_fr is always lower while substituting a Te vdW site compared to the Te ionic site in both Bi-rich and Te-rich conditions. It is due to the fact that the vdW Te is more loosely bound than the ionic Te and can thus be more easily substituted. For In and Po atoms, the lowest E_fr occurs when substituting the Bi atom in both Bi-rich and Te-rich conditions. Incidentally, this is the configuration that is predicted with the resonant peaks for both dopants. In the case of Na, replacing a vdW Te atom yields the lowest E_fr for the Bi-rich condition, while Bi substitution has the lowest E_fr for the Te-rich condition. Since the resonance peak is observed only for Te substitution, Bi-rich growth is necessary.

The bonding between the dopants and the host atoms is also examined. Figure 3 shows the electron charge density along the [010] plane of the Bi₂Te₃ cell for the dopant atoms (In, Po, Na) substituting the host atom. In an ideal Bi₂Te₃ crystal, the bonds between the neighboring intralayer Bi and Te atoms are of a mixed covalent-ionic nature. When In replaces a Bi atom [Fig. 3(a)], a charge density shift toward the Te atoms (away from the In site) is observed due to the comparatively lower electronegativity of In (vs. Bi). In contrast, Bi substitution by Po depletes the charge density from the neighboring Te atoms resulting in a large electron accumulation around the Po atom as shown in Fig. 3(b). For both In and Po doping, an enhancement of the ionic nature is observed in the covalent-ionic dopant-Te bond. On the other hand, Na substituting an ionic Te site has a very weak ionic bond with the neighboring Bi atoms. As seen in Fig. 3(c), the large electronegativity difference between Bi and Na causes complete electron depletion around the Na atom. The weak bonding of Na with its neighboring atoms in Bi₂Te₃ is not surprising, given the small size of Na in comparison to Bi and Te.

It is expected that the resonant peaks formed by proper doping as shown above will increase the power factor and in turn zT of bulk Bi₂Te₃. The effect of these peaks on the transport properties is quantitatively evaluated by calculating S and σ based on the solution to Boltzmann transport equation in the relaxation time approximation¹⁹ 

where q is the unit charge and the differential conductivity $σ(E)$ is given by

In the above expression, $τ−1(E)$ is the scattering rate, $vξ(E)$ the carrier velocity along the transport direction ξ (in this study, x), and f(E) the Fermi-Dirac distribution function. For g(E), the DOS obtained from the DFT is used as shown in Figs. 2(a)–2(c). The total scattering rate is estimated according to Matthiessen's rule accounting for the interactions with acoustic phonon, optical phonon, and ionized impurities, where the directional dependence in the momentum space is ignored. Further details of the scattering calculation and the adopted parameters (i.e., deformation potential constants and effective mass values) can be found in Ref. 20. For the doped Bi₂Te₃, the rate is adjusted to reflect the larger DOS near the band edge (i.e., more scattering). This is achieved by replacing the effective mass in the scattering rate calculation with a value that is extracted from the DOS curve. It is assumed that the increase in the effective mass due to the resonant peaks does not impact the mobility significantly in Bi₂Te₃.¹⁰ A detailed account of the scattering rate with dopant induced states requires first-principles calculation of the carrier-phonon interaction as a function of the doping concentration; this is beyond the scope of the current investigation. As for $vξ(E)$⁠, it is estimated from the derivative of the energy-momentum dispersion based on a non-parabolic Kane model with the effective mass parameters discussed above.²⁰ The obtained values show good agreement with those obtained directly from the DFT band calculation. Once $τ(E), vξ(E)$⁠, and g(E) are obtained as described, Eqs. (3) and (4) can be readily evaluated.

The power factor evaluated at 300 K is shown in Fig. 4(a) as a function of the Fermi energy in the VB. The direction of carrier transport is assumed to be along the basal plane. Note that the minority carriers are also included in the calculation since Bi₂Te₃ is a low bandgap material. The region of interest for E_f is near the resonant peak [0.18 eV away from E_v for In and 0.2 eV for Po and Na; see Figs. 2(a)–2(c)]. The obtained power factor clearly exhibits a maximum about 2.4k_BT below the resonant peak position in accordance with the Mahan-Sofo theory (i.e., $Ev−Ef≈0.24–0.26$ eV).¹¹ A similar observation is made for the resonant peak in the CB with the Po atom. The observed upturn in the power factor is expected to become much larger and sharper for lower doping concentrations as the DOS enhancement will tend to a more Dirac delta-like form; the energy range showing an increase (in the power factor) directly correlates with the broadening of the impurity peak in the DOS. One uncertainty, though, may be with the exact position of the resonant states. While a limited number of simulations conducted in the present study show no evidence of significant shift in energy, a more comprehensive investigation is needed for a definitive answer. For instance, a case with the doping concentration of $∼17$% In substituting Bi atoms results in a resonant peak 185 meV away from the VB edge which is actually very close to the $∼4$% case (180 meV). In addition, the width of the resonant level becomes 60% narrower with the corresponding change in the doping concentration (from 17% to 4%). For a more detailed examination of transport properties, Fig. 4(b) shows the constituents of the calculated power factor (S and σ) for Bi₂Te₃ doped with In. As illustrated, the enhancement in the power factor near the resonant DOS peak (i.e., $Ev−Ef≈0.2–0.3$ eV) is predominantly due to the increase in S. In contrast, the observed gain close to the VB edge in Fig. 4(a) is due to the higher σ that is a consequence of the elevated doping concentration adopted to enable the DFT computation as described earlier. Accordingly, this peak in the power factor is not expected to survive at low doping unlike that related to the resonant states.

In summary, our DFT study identifies Po as the dopant that can potentially realize the desired resonant peak in the DOS of n-type Bi₂Te₃ while pinpointing In, Po, and Na for the p-type. Transport calculations indeed predict the anticipated increase in the power factor that is directly attributed to the resonant features in the DOS. An investigation with a much lower doping concentration for In, Po, and Na will help to quantify the precise enhancement in the power factor.

The authors would like to thank the Tandy supercomputer center for the many hours of computing time. This work was supported in part by U.S. Air Force Office of Scientific Research (FA9550-12-1-0225) and the National Science Foundation (EEC-1160483, ECCS-1351533, and CMMI-1363485).