Analysis of the dependence of critical electric field on semiconductor bandgap

Understanding of semiconductor breakdown under high electric fields is an important aspect of materials’ properties, particularly for the design of power devices. For decades, a power-law has been used to describe the dependence of material-specific critical electrical field (Ecrit\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mathcal{E}}_{\text{crit}}$$\end{document}) at which the material breaks down and bandgap (Eg). The relationship is often used to gauge tradeoffs of emerging materials whose properties haven’t yet been determined. Unfortunately, the reported dependencies of Ecrit\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mathcal{E}}_{\text{crit}}$$\end{document} on Eg cover a surprisingly wide range in the literature. Moreover, Ecrit\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mathcal{E}}_{\text{crit}}$$\end{document} is a function of material doping. Further, discrepancies arise in Ecrit\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mathcal{E}}_{\text{crit}}$$\end{document} values owing to differences between punch-through and non-punch-through device structures. We report a new normalization procedure that enables comparison of critical electric field values across materials, doping, and different device types. An extensive examination of numerous references reveals that the dependence Ecrit\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mathcal{E}}_{\text{crit}}$$\end{document} ∝ Eg1.83 best fits the most reliable and newest data for both direct and indirect semiconductors.


Introduction
Avalanche breakdown due to carrier-induced impact ionization is an important phenomenon as it governs the maximum sustainable voltage in semiconductor devices. The maximum electric field at breakdown (also called the critical electric field, E crit ) is known to increase with the semiconductor bandgap E g [1,2]. This relationship between E g and E crit is typically fitted to a power-law: E crit ∝ E γ g . This dependence is what has ultimately driven the adoption of wide bandgap (WBG) and, more recently, ultrawide bandgap (UWBG) semiconductors in power electronics. In such applications, system performance dramatically benefits from an increase in the breakdown voltage (V BD ) of active devices used, and V BD is directly related to critical electric field.
This paper aims to refine the power-law dependence of the critical field on bandgap and to provide a physics-based explanation for this dependence. This will give a benchmark for future first-principles models and enable a more reasonable prediction of critical electric field for WBG and UWBG semiconductors whose ionization properties have not yet been directly measured. This work pays particular attention to the selection of the material dataset-discussing why the materials were selected or excluded-and to the assumptions necessary to construct and validate the power-law fit of E crit vs. E g .
The discussion in this work strictly applies to p + -n-n + or n + -p-p + diodes for which the heavily-doped regions and the middle drift region are plane-parallel to each other. For vertical diodes, the doping in the drift layer should be less than approximately 10 18 cm −3 to avoid tunneling and degeneracy effects [3,4]. This is the case for which the vast majority of E crit versus E g data has been obtained or calculated. However, the identified dependence of E crit on E g should approximately hold for other geometries and doping profiles as well.
We focus on devices exhibiting avalanche breakdown from impact ionization, where breakdown is uniform throughout the channel and characterizes bulk material properties. Narrow bandgap materials can also break down due to carrier tunneling through the junction, which makes it hard to isolate the tunneling mechanism from impact ionization [4]. A related Invited Feature Paper-Review breakdown mechanism exists where carriers tunnel through the contact-this is of particular concern in Schottky barrier devices, or WBG or UWBG materials operating under high voltages [4]. Material defects can cause localized elevated electric fields, leading to early on-set of breakdown. Device design can also be a source of high electric fields, due to sharp metal lines concentrating field potentials or inadequate edge termination to control fields at junction edge [4]. Surface state-assisted breakdown can create issues in two-dimensional devices. This work limits its scope to bulk devices and so surface state-assisted breakdown is not considered. Lastly, devices can break down due to thermal instability if they operate in a high-temperature regime.
Discussion in this work does not apply to high electron mobility transistor (HEMT) and other lateral structures. These structures generally have very different and much more complicated field profiles compared to vertical diodes [4]. Of note is the spike in electric field at the gate edge in HEMTs which can lead to localized breakdown. To complicate the physics further, the presence of surface states affects breakdown in lateral structures, and this effect is difficult to separate experimentally from avalanche breakdown.
A normalization technique is introduced to establish equivalency between measurements under different conditions (punch-through vs. non-punch-through, doping level, and material types) in §2. It is meant as a simple and intuitive tool for comparing the wide range of structures reported in the literature. However, this technique does not address underlying errors that arise due to poor field management or material imperfections. Some of the shortcomings of the data used for previously derived conclusions on the dependence of E crit on E g are briefly discussed. This paper then explores the historical dependency of reported E crit values on E g ( §3). The normalization technique is applied to historical experimental data to derive a revised power-law dependence of E crit on E g . Finally, a first-order theoretical explanation for why the identified powerlaw is the appropriate one is covered in §4.

Normalization Impact ionization theory
Impact ionization events inside the depletion region of a reverse biased diode increase the current inside the device. An expression for indicating this increase in current is given by the socalled multiplication factor equation. Following Willardson and Beer, for a p + n junction as shown in Fig. 1 under reverse bias so the electric field points left with holes traveling in negative x and electrons traveling in positive x directions, multiplication factor for holes is given by [3]: (1) where the variable x is the linear position within the semiconductor, x D is the depletion depth, α n (x) and α p (x) are the ionization coefficients (cm −1 ) for electrons and holes respectively. For electrons the multiplication factor is given by a complimentary expression: Both expressions approach infinity when the integral in the denominator approaches unity, indicating avalanche breakdown. For holes, the breakdown condition is given by: For electrons the equivalent expression is: Different, yet equivalent, forms of the multiplication factor equations exist. Willardson and Beer state alternatives to Eqs. 1 and 2 in their textbook [3]. Or, as in Sze and Ng, derivation of the multiplication factor in a differing coordinate orientation for the p + n junction yields another variant [4].
Following Selberherr's impact ionization model [5]: the impact ionization coefficients, α n and α p , are a function of the electric field E at position x in the depletion region, with α n0,p0 , β n0,p0 , and E n0,p0 fitting parameters for the model. This model is considered accurate for the scope of this study, but it assumes that the electric field E is constant along the lateral dimensions of the device. It thereby ignores more localized means of breakdown and only considers breakdown that occurs across the entire device. Generally, α n and α p are different and their magnitude is a function of the carrier scattering mechanisms of the material in question. Such carrier scattering mechanisms, including momentum scattering processes, also influence other material properties, such as minority carrier lifetimes and momentum relaxation times. The functional dependence of the electric field on position depends on the device dimensions, the doping level, and on the permittivity of the semiconductor. Complete description of avalanche breakdown conditions requires an iterative solution to the ionization integrals described by (3) and (4), as for example described by Cooper and Morisette [6]. However unless full accuracy is needed, 1-dimensional simplified analysis is sufficient to describe well-behaved devices [4]. (2)

Invited Feature Paper-Review
Two general profiles are of interest: the so-called punchthrough (PT) and non-punch-through (NPT) configurations, as shown in Fig. 1. Solving for both (3) and (4) provides a value for the depletion depth, x D = W NPT or x D = W PT , depending on the device design. This will correspond to a peak electric field at x = 0. This maximum electric field value is referred to as the critical electric field ( E crit,NPT or E crit,PT ) of the device.
Using boundary condition E = 0 at x D = W NPT for a simplified 1D analysis, the breakdown voltage in NPT diodes is related to critical electric field by [7]: where q is the fundamental charge, ǫ is the dielectric permittivity, and N D the doping of the drift layer. The breakdown voltage in a punch-through diode is [7]: Assuming abrupt junction profiles the extension of the depletion layer into the p + and n + layers can be neglected. For devices where drift layer doping is on the same order of magnitude as the p + or n + layers or for extremely thin drift layer thicknesses, the field penetration into these layers must be considered as part of the drift region.
Due to the relation between the critical field and breakdown shown in Eqs. (6) and (7), the critical electric field is considered an important material parameter of interest. An ideal approach to derive this value is to use the impact ionization coefficients of a material and solve Eqs. (1) and (2). Unfortunately, these coefficients are often unknown, so a more common approach is to fabricate a device, measure a breakdown voltage, and use either Eq. (6) or (7) to predict the critical field of the material. As will be shown in the next section, this can lead to incorrect or misleading values for the critical electric fields of semiconductor materials.

Normalization theory
There are two normalizations that typically need to be made to fairly compare E crit data from the literature: (1) for doping and (2) for electric field profile. Researchers often assume that E crit is invariant with doping due to the relatively weak dependence of the critical field on doping (described theoretically by Baliga and Ozbek as E crit ∝ N D 1/8 ) [7,8]. This assumption can skew the accuracy of power-law fits to E crit vs. E g data if it is not accounted for. The effect of this correction will be shown in the next section. If the reported critical electric field is extracted directly from measured impact ionization coefficients for a given material, then one must simply compare the results of Eqs. (3) and (4) using the same value of the drift layer doping N D .
The second necessary correction accounts for the differing electric field profiles of non-punch-through (NPT) and punchthrough (PT) cases. As Fig. 1 illustrates, the field profiles of the two structures are not directly comparable. For equal doping levels, the critical field in a PT device will always be greater than that of an NPT device. This subtle increase can be seen by solving (3) and (4) using x D = W PT for a PT design and x D = W NPT for an NPT design. Additionally, Eq. (5) must be used with the appropriate electric field profile, for NPT and for PT. Essentially the reduction in the drift region thickness suppresses impact ionization events and leads to an increase in the critical electric field needed for (3) and (4) to remain true. The magnitude of E crit,PT /E crit,NPT above unity increases as W PT / W NPT decreases or E n0,p0 from Eq. (5) decreases. The extent of the change in E crit,PT can be negligible in some cases, but nevertheless must be calculated before it is ignored.
For these reasons, E crit data reported in literature may need to be corrected to account for structure (PT vs. NPT) as well as doping. The basic principle used to correct for both cases is to equate the one-dimensional ionization integrals, Eqs. (1) and (2), for the case of different doping levels and PT or NPT. This correction only considers the lower-doped device drift region. As the doping of the adjacent layers is typically much higher than the drift region, depletion into those regions is negligible and the analysis does not account for it.
We consider first the correction for doping in an NPT configuration, where two materials with doping levels N D1 and N D2 are compared. In order to proceed, a critical assumption that the electron and hole ionization coefficients are equal (α n = α p ) must be made. For the vast majority of materials this is not the case. Nevertheless, it has often been used as a starting point in previous analyses, allows an intuitive understanding of general trends, and enables us to develop a simple analytical solution. Furthermore, it can be shown that this assumption gives the lower bound for the critical electric field in a material. Using this assumption makes Eq. (3) equal to Eq. (4), allowing for the following relation: where x D1 and x D2 are the depletion depths at breakdown for doping levels N D1 and N D2 respectively and α is assumed to be invariant with doping. For the non-punch-through case, Eq. (8) can be used to translate the integration variable from position to electric field, and one obtains: where E crit1 and E crit2 are the critical electric fields for doping levels N D1 and N D2 in the non-punch-through case. At this point, a second critical assumption is made, which is that the ionization rate α follows a power-law in electric field [9]: This work sets δ = 7 which has been shown to be consistent with ionization data for a wide range of semiconductors, including SiC [7]. However, the specific value of δ used, and indeed whether Eq. (12) is truly valid for all semiconductors, is open to debate. A comparison of Fulop's expression and ionization rates derived directly from impact ionization parameters is the subject of future work. Nevertheless, as was the case with the assumption α n = α p , this is necessary to obtain an analytic expression by which to equate critical electric field, which is the goal of the derivation. The validity of both assumptions, examined using numerical evaluation of the ionization integral, is currently underway in our group. Inserting this into the integral and performing the integration, one obtains: Next, we account for the second correction, i.e., NPT compared to PT. For the purposes of this derivation, it is assumed that the doping in both cases is equal (the doping correction will be brought in at the end). The equality of the ionization integral for the two cases, Eq. (10), is again the starting point, with equal electron and hole ionization rates assumed. This may be written as: where x NPT and x PT are the depletion widths for the non-punchthrough and punch-through cases, respectively. In the latter case, this is approximately equal to the physical thickness of the drift layer. In both cases, the electric field profile is a linear function of position, although for the PT case the field does not go to zero, so the field distribution is trapezoidal rather than triangular. If the critical fields for the NPT and PT cases are denoted as E crit,NPT and E crit,PT , respectively, transforming the integration variable from space to electric field yields: which is the analog of Eq. (11) for the PT-to-NPT transformation. Again, assuming the Fulop power-law dependence holds [8], and performing the integration, one obtains: The expressions for the two corrections may be multiplied together to obtain the expression where both corrections are considered: In this general expression, the critical field is normalized to an NPT condition with doping N D,NPT , from a PT condition with doping N D,PT . For the calculations described below, δ = 7 was used, consistent with the value used by Fulop.
The same derivation for normalization between NPT and PT devices, but with integration over x rather than the electric field E , is covered in Appendix A.

Example of normalization correction
Our normalization procedure requires the device to be a single region of uniform doping, which is not always the case for reported devices. The device may need to be reduced to an equivalent drift region as illustrated by normalization of the devices fabricated by Ohta et al. [10,11] and Volpe et al. [12] Ohta et al. reported GaN devices with multiple regions of different doping levels, with the most advanced device using three drift layers [10,11]. Numerical simulation to replicate the results and detailed comparison of the two reports strongly suggests that the doping profile in Fig. 1 of the 2015 paper [10] is non-compensated. Figure 1 of the 2019 papers shows lower, presumably compensated, doping of the three drift layers: 5.5 μm at N D = 5 × 10 14 cm −3 , 22 μm at N D = 4 × 10 15 cm −3 , and 5.5 μm at N D = 1.5 × 10 16 cm −3 , To obtain an equivalent drift region, we calculate layer charge multiplying each drift layer thickness by its doping: 2.75 × 10 11 cm −2 , 8.8 × 10 12 cm −2 , 8.25 × 10 12 cm −2 . These are added for a total charge of 1.73 × 10 13 cm −2 over a 33 μm thick region. Assuming uniform doping over this thickness results in N D = 5.25 × 10 15 cm −3 . With a 4.9 kV breakdown in such a region, Eq. (6) predicts a depletion depth of ~ 31 μm making this a non-punch-through device.
In some instances, grown epitaxial layers are highly nonuniformly doped as Volpe et al. [12] reported for their growth of C (diamond). As shown in their Fig. 1, C-V measurements were used to obtain the compensated doping profile. Using numerical integration on the sampled data gives a total charge of ~ 2.2 × 10 13 cm −2 over the 13.6 μm thick region. Dividing the total charge by the drift layer thickness gives a uniform N D ~ 1.65 × 10 16 cm −3 for a 13.6 μm drift layer. With a 9.8 kV breakdown in such a region, Eq. (6) predicts depletion of 19.6 μm for a NPT device, indicating that this device is well into the punch-through regime.
Together these examples illustrate the necessity of reporting complete compensated doping profiles, as these values affect many extracted devices parameters including E crit .
Without proper normalization correction, reported critical field values can be unrealistic. In 2006 and 2007, researchers published data on sets of nearly identical p + -n-n + Al x Ga 1−x N diodes were grown wherein the Al composition of the middle layer was varied from 0 to 57% across multiple wafers, with E g of 3.4-4.6 eV, as measured by photoluminescence [13][14][15]. These diodes were strongly punch-through as the middle layer of each was only 225 nm thin with a nominal background doping level of 2 × 10 16 cm −3 . The measured E crit values reported in Nishikawa et al. [15] did not correct for the PT case, resulting in overestimations of the critical field values in Al 0.57 Ga 0.43 N by as much as 2.5 MV/cm. Using these corrected values, the dependence of the critical field on bandgap is closer to E crit ∝ E g 2.0 with an R 2 value of 0.46. This indicates a need to consider all reported dependencies by first normalizing the data to the same drift region doping level and accounting for PT vs. NPT conditions. Table 1 compares calculated E crit values using the normalization procedure outlined above for the work by Nishikawa et al., a few other GaN and Al x Ga 1−x N devices, a β-Ga 2 O 3 Schottky diode by Yang et al., and the C (diamond) Schottky diode by Volpe et al. [10,12,[15][16][17][18][19] Drift layer width, doping (compensated if so stated in the paper), and breakdown voltage are used to calculate E crit,PT . The corrected E crit,NPT was found for the same junction along with the extrapolated depletion depth W NPT , and then normalized to N D = 10 16 cm −3 to allow comparison between devices with different drift layer width and doping values. In some devices the NPT correction is indeed minoronly noticeable after several digits past the decimal point, but it is not the case for all and is in fact quite dramatic in several devices with extremely thin drift layers. Further, note that for the work by Armstrong et al. [17] the value of the GaN critical field listed here was extracted directly from the reported breakdown voltage and then normalized, while the value reported in the paper was determined from the Baliga Figure of Merit, which involves both breakdown voltage and specific on-resistance. The latter approach requires knowledge of quantities such as carrier mobility and effective device area, while the former does not, so these two approaches may not necessarily lead to the same critical field value even after normalization.

Historical fitting of E crit vs. E g dependencies
This work carefully examined numerous relevant publications in the avalanche breakdown literature. Historically, several power-law fits were used to describe the relationship between critical electric field and semiconductor bandgap: The relationship arose from initial characterization of material impact ionization coefficients, in a way simplifying interplay of those parameters and device breakdown to a single value: critical electric field. This kind of fitting has been successful, but that does not mean it can be applied limitlessly, nor that it indicates a physical basis for the power-law relationship between critical electric field and bandgap [20]. between closely spaced p + and n + layers. Note that E crit, PT > E crit, NPT as explained in the text.
The first established the power-law fit is traceable to a 1966 paper by Sze and Gibbons [21]. The authors used the measured ionization rates for electrons and holes in Ge, Si, GaAs, and GaP to calculate the breakdown voltage as a function of doping N D . Their use of the relation . This equation has been presented unaltered in all editions of Sze's textbooks [4].
In 1989 Kyuregyan et al. updated the ionization coefficients for a number of materials by compiling and aggregating earlier published results [22]. As Maes et al. discussed in their definitive 1990 paper on impact ionization in Si [23], this approach is better than utilizing multiple reports each relying on a single device, since it allows for exclusion of outliers and covers a wider range of electric fields. Kyuregyan et al. also reported R 2 values of their fitting. For a number of materials, these ionization coefficients were the most recent we could identify and are the basis of many of the E crit values utilized here.
In A summary of many semiconductor parameters, including critical field as noted above, impact ionization coefficients, bandgaps, and dielectric constants can be found in Levinshtein et al. [26][27][28], though WBGs and especially UWBGs are still being actively investigated. We have surveyed the literature to update and provide the best critical electric field values for a range of materials. The materials included here will be discussed in order of ascending bandgap. A comparison of the normalized critical fields verses the values used by Sze et al. [21], Hudgins et al. [2], and Wang [29] are shown in Table 2. It should be noted that Table 2 lists normalized critical field values that are now felt to be more accurate than those published by some of the present authors in 2018 [30].

Updated E crit values
Demonstrated avalanche breakdown due to impact ionization is necessary for consideration of a semiconductor material to be used in the power-law fit of E crit vs. E g . For many materials these parameters have been measured and can be used to make E crit vs. doping plots. Where ionization coefficients were not available, we looked for a positive correlation between breakdown voltage and temperature. This is a reliable indicator that breakdown is caused by avalanche from impact ionization events, rather than  16  some other mechanism. Materials that have not yet shown this behavior will be noted and are excluded. In the end, the majority of the data is from established materials. Unless otherwise mentioned, normalized data were obtained by analyzing the more descriptive E crit vs. doping plots found at the end of chapters in Semiconductor Parameters Vols. 1 and 2, or from V BD vs. doping data in the same volumes [26,27]. In many cases, the original data are no longer available, and only figures can be found. In such cases, E crit values were obtained by line-tracing the digitized plots.
We surveyed the literature for reports on WBG and UWBG devices and materials. For many of these, only individual devices are reported and we used our normalization correction to extract and compare E crit . Due to incomplete understanding of impact ionization in GaN and C (diamond), critical electric field values for these materials were obtained from individually reported devices [10,12].
Since the 1980's newer techniques have been developed to grow cleaner semiconductor materials and to more accurately measure breakdown parameters and ionization coefficients. Availability of these new methods is especially important for materials that may still contain significant levels of defects, like SiC, GaN and emerging UWBG semiconductors [7,8,31,32]. One of the newer techniques is electron-beam-inducedcurrent (EBIC) wherein a scanning-electron microscope (SEM) is used both to apply electron pulses to stimulate direct carrier generation as well as to image test diodes. EBIC has become a key tool for characterizing materials with potentially high defectivity. It is used to reject diodes with 'hotspots' that result in premature breakdown [33], and it may be used to obtain ionization coefficients. Though valuable, these data are not necessary for building our dataset. In addition, although some devices may be engineered to manipulate material properties such as minority carrier lifetimes to artificially increase the critical field, these data would only create outliers in our data set and have therefore been excluded wherever possible. Although it is recognized that such material properties as minority carrier lifetime and momentum scattering have some relationship with device breakdown, these properties are not directly relevant to our normalization procedure and are therefore left for further consideration.
Note that in a few instances even reference plots are not up-to-date. For Si, most recent coefficients are from Maes [22,26]. And Sze and Ng still report the original 1966 values as far as authors can tell [4]. 4H/6H-SiC E crit values were taken from Raghunathan et al. [32] There are noted issues with that experiment, and more recent measurements of impact ionization coefficients are found in literature [34]. This necessarily introduces error in fitting the exponent γ, no more than 0.1 in authors' estimate. However, a thorough review and update spanning all materials would constitute its own work and is left for the future. Finally, keeping outlook toward development in WBG & UWBG materials, critical electric field values are discussed assuming doping of 10 16 cm −3 . Such doping is fairly typical for these semiconductors as lower levels can be hard to achieve.
Each of the semiconductors is discussed in order of ascending bandgap below.

InSb (0.17 eV) and InAs (0.354 eV)
In narrow-bandgap semiconductors, injection of carriers from tunneling cannot be decisively excluded as a contributing factor in carrier multiplication under high reverse bias [35,36]. Due to this, we believe that extremely narrow-bandgap materials, such as InSb and InAs, cannot be fairly compared to other semiconductors. Therefore, while these are included in Hudgins et al. analysis [2], we excluded them in this work.

Ge (0.661 eV)
The data in Sze & Gibbons [21] for Ge were obtained from devices made in the mid-1950s. The value used in this study for Ge was obtained from the impact ionization parameters utilized by Kyuregyan et al. [22] GaSb (0.726 eV) GaSb impact ionization coefficients have been measured, but it seems E crit vs. doping not calculated from them and was shown in Semiconductor Parameters Vol. 1 [26]. As no other reliable data in the literature was found, this material is also excluded from analysis. Semiconductor Parameters Vol. 2 reports both ionization coefficients and measurement of breakdown voltage vs. doping for In 0.53 Ga 0.47 As [27]. E crit for In 0.53 Ga 0.47 As was found using Eq. (6) from the breakdown voltage in Vol. 2 [27].

Si (1.12 eV)
As can be expected, many researchers have made measurements to determine ionization coefficients in Si over the decades; an excellent review of this work is provided in Maes et al. [23] As it stands as common reference, the E crit vs doping plot of Sze and Ng [4] was used in this work.

InP (1.344 eV)
The value for InP was obtained from the impact ionization parameters utilized by Kyuregyan et al. [22] These data were based on the results of four device papers published from 1979 to 1982.

GaAs (1.424 eV)
The data in Sze & Gibbons [21] for GaAs were for devices fabricated in the mid-1960s and as such the defect density was likely high. The value for GaAs was obtained from the impact ionization parameters utilized by Kyuregyan et al. [22]

GaP (2.26 eV)
The reported values for the critical field of GaP published in Vol. 1 [26] are from work by Sze and Gibbons in 1966 [21]. In turn, Sze and Gibbons relied on Logan and White [37], which was the sole source for Kyuregyan et al. as well [22]. No other power devices with reported values for the critical electric field or impact ionization parameters were found in the literature. Normalized E crit of 0.69 MV/cm looks like a large outlier compared to rest of the data. For these reasons GaP was excluded from our analysis.

3C-SiC (2.36 eV) 6H-SiC (3 eV) 4H-SiC (3.23 eV)
SiC is a WBG semiconductor that exists in several polytypes, but the primary focus has been on 3C, 4H, and 6H. The critical field values used by Hudgins et al. [2] for SiC appear traceable to Baliga [24]. Neudeck et al. published data for 3C-SiC p + -n-n + diodes at different drift layer doping levels and listed derived E crit values [38]. We were unable to determine the width of the depletion region of this device and thus cannot conclude whether it is PT or NPT. Without this information, E crit cannot be accurately determined, and only an estimate can be made, up to ~ 0.98 MV/cm. Without further information the 3C polytype of SiC is excluded from the dataset. However, Raghunathan et al. performed extensive studies of impact ionization in 4H-and 6H-SiC devices using pulsed EBIC and corrected for PT structure [32]. Several characterizations of impact ionization coefficients of 4H-SiC using photomultiplication [5], most recently in 2014 by Niwa et al. are also reported [34]. As other studies do not report the derived E crit values the results of Raghunathan et al. are used here.

GaN (3.45 eV)
GaN epitaxial growth and device fabrication have undergone significant development in recent years. Positive dependence of breakdown voltage on temperature was reported by a few groups [11,39], and impact ionization coefficients for GaN were recently obtained by Ji et al. [40] Still, uncertainty remains regarding the critical electric field of GaN. Hudgins et al. cites 3 MV/cm, and 3.3 MV/cm is often quoted in the literature. Work by Avogy indicates that E crit is higher than this, at least 3.5 MV/ cm [39]. Surveying the literature, the highest E crit found is 3.4 MV/cm, after PT correction and normalization, as reported in a device by Ohta et al. [11] It must be noted that critical electric field predicted by Ji et al. is lower, < 3 MV/cm [40], than this champion device. As most of literature suggests E crit > 3 MV/ cm, we kept the higher value. Further work needs to be done to reconcile these results, but this is beyond the scope of this paper.  N (3.45-6.1 eV) The authors have not been able to find reports of impact ionization measurements on the Al x Ga 1−x N system, including AlN, although a variety of breakdown measurements have been reported for different structures and doping profiles [13-15, 41, 42]. Unfortunately, none of these works show breakdown vs. temperature data to indicate true avalanche behavior so we exclude this material system from our analysis.

β-Ga 2 O 3 (4.9 eV)
β-Ga 2 O 3 is an emerging UWBG material with significant potential. Impact ionization measurements on the β-Ga 2 O 3 system are unavailable, nor could we locate any publications experimentally confirming temperature-dependent behavior indicative of true avalanche breakdown in β-Ga 2 O 3 . Its critical electric field value is reported to be 8 MV/cm, but this is only an estimate [43], most likely using predictions by Hudgins et al. Reports of Schottky barrier diodes [19,[44][45][46] yield much lower E crit values, suggesting that defects limit performance. The best value reported in the literature from Yang et al. [19] is 3.0 MV/cm after normalization, and is likewise below theoretical value. Thus, we excluded this material from our analysis.

C (diamond) (5.5 eV)
The high bandgap and path to dopability of C (diamond) make this UWBG semiconductor an attractive candidate for future use in power electronics [30]. There is evidence of breakdown by impact ionization via a positive trend of breakdown voltage with increasing temperature in diamond diodes reported by Suzuki et al. [47]. While shallow acceptor/donor doping hasn't yet been achieved, 2DHG looks a viable approach to doping the material.
Critical electric field for diamond is not yet settled with a range of values reported. Unfortunately, the second diamond E crit value of 7 MV/cm used by Hudgins et al. [2] appears incorrectly referenced as we cannot find any experimental basis in Ref. [17] to get this value. Landstrass et al. [48] refers in passing to a breakdown field of 20 MV/cm in diamond diodes, but we have not been able to confirm a critical field this high via the information given, so these data are excluded from further consideration as well. Similarly, Liu et al. [49] report a breakdown field in diamond of 21.5 MV/cm, but this is from a laser measurement where breakdown was detected via a flash of light observed by the naked eye. As undetected avalanching could have been occurring at lower fields, this result is also excluded. In 2010 Volpe et al. published a diamond Schottky barrier diode with a breakdown voltage of 9.8 kV [12]. PT correction of this device is discussed in §2.C. Correcting for PT and normalizing gives an E crit = 10 MV/cm. This is the best value we could confirm, and as such use it in our power-law fit. In 2014 Traore et al. published results on diamond Schottky barrier devices with even better characteristics [42]. Unfortunately, they do not report a breakdown voltage, due to power supply limitations. PT correction and normalization gives only a lower bound of E crit = 7 MV/cm.

Power-law fitting of updated E crit values
Data points in Table 2 were fit by first taking the logarithm of the data, and then performing a linear least-squares fit to the resulting log 10 (E crit ) versus log 10 (E g ) data points. Assuming the postulated relationship E crit ∝ E γ g , the slope of this fit yields the exponent γ. The bulk of the data (7 of 9 materials) are wellestablished and from reference tables. As shown in Fig. 2 the full corrected dataset is reproduced by a power-law fit with exponent of ~ 1.83. The fit limited to only reference data Ge through SiC gives γ ~ 1.74.
Appendix Fig. B1 shows comparison of normalized data and fitting to doping level 10 16 , 10 15 and 5 × 10 14 cm −3 covering a range more typical of non-WBG semiconductor devices. No doping dependence is seen in experimental data.
As can be seen in Fig. 2, the updated exponent γ provides a better fit to the normalized data than the previous literature estimates. Appendix Fig. B2 shows the separated direct-and indirect-gap data, along with their respective fits. Exponents are 1.66 with R 2 = 0.954 and 1.89 R 2 = 0.992. However with reduced sample sizes and small difference between the direct-and indirect-gap cases, a single fit works best to explain the E crit vs. E g relationship for all semiconductors.
Based on the physics of impact ionization, the direct vs. indirect nature of the bandgap is not expected to affect critical electric field values for a semiconductor material. Carriers undergo impact ionization at energies much higher than the bandgap energy (E > 1.5E g ) [23] where the direct vs. indirect nature of the bandgap is no longer decisive. Hudgins et al. concluded on using different fits for the direct vs. indirect cases solely from R 2 values of their fits. Our updated fitting shows that this is no longer warranted. The critical electric field of GaN is still debated, and many UWBG materials are not near their predicted potential [30], so it is realistic to expect that further clarification on the relationship between bandgap and critical electric field will be obtained in the future. Furthermore, the power-law fit is not expected to hold indefinitely for arbitrarily high bandgaps.
To our knowledge a theoretical explanation of the observed power-law dependence of the critical electric field with bandgap has not been reported in the literature. In the following section, a theory of avalanche ionization first presented by Ridley [50] is generalized to a dependence on semiconductor bandgap and is used to calculate the normalized critical electric fields as shown in Table 2. These data are shown in Fig. 3, with a slightly lower power-law fit of 1.82 with R 2 = 0.999.

First-order model to explain E crit vs. E g 1.83 dependence
To determine the theoretical dependence of E crit on E g , it is necessary to explore the avalanche breakdown phenomenon and to parameterize the expression as a function of bandgap. Avalanche is assumed to occur when the multiplication factor for either holes or electrons approaches infinity as discussed in Sect. 2 and described by Eqs. (1) and (2).
Specifically, for the NPT case of a p + -n − one-sided step junction, the electric field can be approximated by a linear field profile described by Eq. (8) in Sect. 2. Using this to transform Eq. (3) from an integration over position x into an integration over electric field E , one obtains: A complimentary equation describes the case for the avalanche multiplication of electrons. Equation (18) is simply the generalization of Eq. (1), without the assumption that α n = α p and for a single material only. Also, the permittivity ǫ has been explicitly written as the product of the relative permittivity of the semiconductor ǫ r and the permittivity of free space ǫ 0 . Equation (18) indicates that the ionization multiplication in a p + -n − junction is dependent only on the dielectric constant of the material ( ǫ r ), doping (N D ), the ionization rate (α n , α p ), and the electric field (which depends on the applied voltage). For a (18) given device, the critical field can be found by evaluating the integral in Eq. (18) until the expression approaches unity.
In developing the theoretical model, as with the derivation of normalization in §2.2, the ionization rates of electrons and holes are assumed to be equal. This assumption gives a minimum E crit value for the material in question. For moderate-towide bandgap semiconductors, the value of E crit can increase slightly ( 10% for up to a magnitude of difference between α n and α p ) if α n = α p . The effects of this α n = α p assumption on critical electric field is an ongoing area of research. In the case of α n = α p , the avalanche breakdown condition simplifies to: Of these parameters, all but the ionization rates are well characterized for many materials. In order to model the electron and hole ionization rates α n and α p for different materials, we utilize the lucky-drift model of the electron as reported by Ridley [50]. This is seen as a more accurate representation than either the lucky electron theory of Shockley [51] or the thermalized electron model of Wolff [52]. In fact, both the Shockley and Wolff theories can be recreated by limiting approximations of the lucky-drift theory.
The lucky-drift theory describes the ionization rate as [50]: where x = E th eE , ζ = P T 2rx 2 , P T = 1 − e −2rx(x−3) x≥3 , where E th is the threshold energy required for ionization, λ is the mean relaxation length for the carrier, and r is the ratio of average energy loss per collision to the threshold energy, E th . Therefore, Eq. (20) depends only on the applied electric field ( E and three material parameters: the threshold energy required for ionization (E th ), the mean carrier relaxation length (λ), and the ratio of average energy loss per collision to the threshold energy (r)). In order to understand the relationship of the ionization rate to material bandgap, we transformed E th , λ, and r into functions that depend on bandgap.
The threshold energy E th is the energy that a hot carrier must possess to create an electron-hole pair. While the assumption by many is that ionization can be initiated by any electron with energy > E g , carrier energy must actually generally be 1.5E g or more as explained in Maes et al. [23] This relationship is derived from conservation of energy and momentum, and assumes that the effective masses of electrons and holes are equal (relaxing the latter assumption results in a constant factor between 1.0 and 2.0 multiplying E g , rather than 1.5). While it generally holds for many materials, it is not always true, but is nevertheless the best starting point for a general, intuitive theory. Fully relaxing this relationship requires more complicated treatment of band structure and is beyond this paper.
The ratio of average energy loss per collision to the threshold energy (r) is described by: where ℏ is the reduced Planck constant, ω is phonon angular frequency, and n(ω) is the quantization number.
Of the variables in Eq. (20), only E th is a direct function of E g . Since E th is proportional to E g , Eq. (21) reduces to: The parameter λ describes the mean free path length of a carrier before thermal relaxation, as a hot carrier must interact with an electron-hole pair before it is thermalized to an energy below E th . The mean free path is actually energy dependent [53] and is equal to the product of the group velocity of the carrier (v g ) and the scattering time (τ).
In general, scattering time (τ) depends on the particular scattering mechanism with several competing mechanisms (ionized impurity, dislocation, acoustic phonon, optical phonon, etc.) occurring simultaneously [54]. The intervalley and interband scattering processes, which result in the absorption or emission of optical phonons, are the dominant mechanisms at high temperatures and electric fields [55]. For non-polar optical phonon scattering, which is important for the majority of semiconductor materials, the scattering time is energy independent and depends inversely on the effective mass (m * ) [55]: The group velocity of carriers away from the band edge is dependent on energy (E) and is given by [56]: If we consider a carrier at the threshold energy, then Eq. (23) becomes: The mean free path (λ) is then a function of E th (which has already been shown to be proportional to bandgap) and effective mass (m * ). To assess the approximate dependence of effective mass on energy gap for a wide range of materials, we have surveyed the reported effective masses for the primary bands for electrons, light holes, and split-off holes for the materials shown in Table 3 and have computed the exponent k using the equation: for each material pair, where m * x represents the effective mass and the subscript '2' denotes the wider bandgap material. The data for this survey comes from Refs. [21,22,50], and it should be noted that fitted results for heavy holes were excluded from Table 3 as that trend is quite sublinear; note that heavy holes are not expected to participate as strongly in avalanche ionization due to higher mass, lower group velocity, and higher scattering rate. Striking in Table 4 is that the dependence of effective mass for both electrons and the lighter holes from InSb through the mid-and wide bandgap materials is nearly linear, with the numerical average for k for all cases shown being k = 0.93. A linear relationship between m * and E g can also be derived via Bloch Theory for a periodic potential, although this treatment is not shown here [57].
With E th and m * both being linearly related to E g , Eq. (26) can be directly related to E g : With E th , r, and λ all related to E g , it is possible to solve Ridley's avalanche equation for a variety of material systems. Toward this end, to determine the predicted E crit for a variety of materials, a computation program using Python was developed that incorporates the equations presented here and iterates voltage to find critical field. A brief flow diagram for the program is shown in Fig. 4. Although the three parameters of the lucky-drift model (E th , r, and λ) are functions of material parameters, for simplicity these parameters were taken to be values having a dependency on bandgap with proportionality constants based on those measured for Si (as shown in Table 4). The threshold energy E th was assumed to be proportional to bandgap with a proportionality constant of 1.5 [23]. The ratio of average energy loss per collision to the threshold energy r was taken to be proportional to E g −1 .
Ridley lists room temperature values of r for Si ranging from 0.049 to 0.063 [53]. Assuming a middle value of r = 0.056 and correcting for the Si bandgap (E g = 1.12), we assume a proportionality constant of 0.063. Last, the mean relaxation length λ is assumed to vary as E g −3/2 . For Si (E g = 1.12), the carrier mean free path has been calculated by multiple authors and been found to be in the range of 7.1 to 11.94 nm [58][59][60][61][62]. We assume that the mean free path is approximately 10 nm, which would give a proportionality constant of 12 × 10 -7 cm. We assume this mean free path holds for materials of other bandgaps. Deviation from the assumption that materials constants scale directly with bandgap due to particular materials properties will lead to deviations off the trend and data scatter dependent on the particular material of interest. Results from the simulation algorithm are shown in Table 2. These results are plotted (Fig. 3, solid trace) against the experimental results and best fit (dashed trace) as described in §3. The simulated results give excellent agreement with the experimental data for materials with E g > 0.5 eV. The simulated critical field varies as a power-law in bandgap with a slope of 1.67. This shows reasonable agreement with the slope of 1.83 derived in Sect. 3.3.

Asymmetries in electron-hole ionization coefficients
In the 1st order approximation of the theoretical model, as with the derivation of normalization in §2.2, the ionization rates of electrons and holes are assumed to be equal: α n = α p . In reality, for almost all materials studied, there is an asymmetry in the ionization coefficients between electrons and holes. In some instances, this asymmetry can be many orders of magnitude. In this section, we consider relaxation of this assumption with a preliminary evaluation of the effect of carrier ionization asymmetry on dependence E crit ∝ E γ g . A full treatment of other cases of electron-hole ionization coefficient asymmetry (including   asymmetries that are non-constant with electric field) will be fully described in future work.
To evaluate the case of ionization coefficient asymmetry, we consider α p to be linearly related to α n by a constant C.
This linear relationship can be used to explore the effect of ionization asymmetries and the previous assumption can be reduced to a sub-case where C is unity. For these simulations, α n was calculated as described previously in Eq. (20). α p was then derived using Eq. (23) utilizing values of C ranging from 10 -4 to 10 4 . The results of these simulations are shown in Fig. 5. Ridley states that the lucky-drift model loses accuracy in low-bandgap materials, though omits the precise cutoff [50]. Thus, as with the experimental data analysis in §3, we chose to only include materials with bandgaps equal or greater than Ge (E g = 0.661 eV). Appendix Fig. B3 includes materials below Ge and shows these being far off trend, validating their exclusion.
In the chosen bandgap span, 0.661 to 5.5 eV, all fits to simulated data points are excellent with R 2 > 0.999. We take the agreeableness of R 2 to mean that physics are well-behaved, and the impact ionization described by Ridley's lucky-drift model is self-consistent. As can be seen in supplementary Fig. B3, data points below 0.661 begin to deviate from the fit, supporting the loss in accuracy of Ridley's model at low bandgaps.
Over the range of C = 10 -4 to 10 4 , the power exponent γ goes from 1.66 to 1.99. While there is a dependence of γ on C, that it only varies 0.33 over a range of 10 8 shows remarkable capability of Ridley's luck-drift model of describing impact ionization physics and avalanche breakdown. Moreover, as in actual semiconductors C is material-dependent and does not follow a neat relationship with bandgap-instead varying all over; the close match of normalized empirical data with simulated data is strong support for the validity of both approaches.

Conclusion
This work has carefully examined the relevant literature on avalanche breakdown. The analysis has shown that many of the previous reports of the behavior of E crit vs. E g have been influenced by non-optimal experimental data, or non-optimal fits, or uneven comparisons between different fabricated devices. In this work, we have introduced a normalization technique that can be used to correct for differences in doping and device design to develop a fair comparison between breakdown measurements.
By normalizing the available data, the best relationship between E crit ∝ E γ g was found to be a power-law with γ = 1.83. Contrary to long-standing projections by Hudgins et al., we see no need for separate power-law fits for direct-and indirectbandgap semiconductors. Further, the dependence of critical electric field on bandgap is slightly weaker that reported by Hudgins et al., with γ < 2.0.
The relationship between E crit and E g was then derived via a first-principles calculation of the avalanche mechanism using the expression for the ionization coefficient derived by Ridley and applied to materials with different bandgaps. Simulations of E crit vs. E g using these equations can re-create the relationship shown by the normalized experimental data.
This new relationship has implications for the usage of WBG devices for power electronics as well as RF applications. For example, based on the E crit ∝ E g 1.83 dependence, the relationship between specific on-resistance (R ON,sp ), breakdown voltage (V BD ), and E g for power switches over this bandgap range is best described by R ON,sp ∝ V BD 2 E g −5.49 for both direct and indirect semiconductors. A literature review of the latest data shows that the historical relationship (R ON,sp ∝ V BD 2 E g −7.5 ) for direct semiconductors may be overly optimistic [2], placing too much emphasis on the breakdown performance of devices. This new analysis and theory also has application for emerging ultrawide bandgap semiconductors, for which accurate measurements of the impact ionization coefficients and critical electric field have yet to be made. Finally, as best as possible we stated our assumptions and reasoning for building the power-law dataset, deriving the normalization procedure, and simulating the semiconductor critical electric field from first-principles. These assumptions are meant to be relaxed and challenged in subsequent work. We hope this initial analysis enables discussion and opens avenues of future research, driving further scientific progress.

Acknowledgments
Mark Hollis would like to thank William Loh and Erik Duerr at MIT Lincoln Laboratory for their helpful discussions on this work. Oleksiy Slobodyan, Jack Flicker, Jeramy Dickerson, Andrew Binder, Trevor Smith, and Robert Kaplar would like to thank Albert Baca and Mary Crawford of Sandia

Data availability
Data that support findings of this study are available from the corresponding author upon reasonable request. These data come from published work, including materials handbooks, which are referenced herein. The parameters and equations necessary to reproduce the findings of this study can too be found or derived from published work or textbooks. However, no single source exists for all the material, and the authors kept their own aggregate record during writing of this article.

Conflict of interest
On behalf of all authors, the corresponding author states that there is no conflict of interest.

Open Access
This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.