Skip to content

Original Research

Frontiers in Nanotechnology ·

Volume 8 - 2026 | https://doi.org/10.3389/fnano.2026.1885901

Statistics and physics of area-dependent dielectric breakdown: an overlooked yet inevitable natural phenomenon

  • Ernest Y. Wu1*corresponding author
  • Takashi Ando2
  1. 1University of Vermont, Burlington, VT, United States
  2. 2IBM Research Division, Yorktown Heights, NY, United States
Article metrics0 viewsView details

Abstract

The traditional inverse area-scaling law for dielectric breakdown (BD) has been recognized for more than a century since, but its physical basis is rarely discussed, aside from an earlier pioneering study. This may be due to the long-standing misconception that area-dependent BD is purely statistical, rather than acknowledging the transformation of a dielectric layer under electric stress from pristine to degraded and ultimately to fragmented because of BD. This unique area dependency distinguishes the BD mechanism from other failure mechanisms such as bias-temperature stability (BTI) and hot-carrier degradation (HCD). Recent experimental observations of reverse and diminished area-scaling phenomena in hydrogen-doped HfO 2 have prompted investigations into these effects presented in this article. Through probability density function (PDF) analysis, we examined the interplay between defect generation and annihilation from both statistical and physical viewpoints. Our findings show that time-dependent decay in annihilation centers is consistent with maximum-value statistics and reduces the likelihood of annihilation over time, thus supporting defect survival and altered BD area-scaling behaviors. Both defect creation and annihilation scale with area, making their combined dynamics compatible with a thermodynamic description of defect-generation–annihilation and leading to diverse area-scaling phenomena in conjunction with the joint minimum-value and maximum-value statistics. We show the traditional inverse area-scaling law is merely a special case assuming infinite annihilation times. Based on the agreement between our BD data using hydrogen-doped HfO 2 stacks and our defect-generation–annihilation model, we propose that intrinsic oxygen interstitials and extrinsic species such as hydrogen can act as annihilation centers. Notably, the estimated time span of the defect annihilation PDF aligns well with the estimated times from ∼1.5 ns to 1.5 μ s for hydrogen substitution into oxygen vacancies based on density functional theory (DFT) calculations. Overall, contrary to conventional wisdom, the root cause of BD area dependency is the most critical and inevitable for a comprehensive physical understanding of dielectric BD.

Keywords

1 Introduction

Dielectric BD is a failure mechanism that can lead to malfunction of devices andcircuits in a wide range of applications such as transistors, electronic memories, and high-energy capacitors. The time to formation of a conducting filament between two electrodes across a dielectric, often called time-to-breakdown (t<sub>BD</sub>), is statistically distributed. They are traditionally characterized by the Weibull model, F t B D = 1 − exp − t B D / t 63 % β , where F is the cumulative distribution function (CDF). t63% is the t<sub>BD</sub> at 63.2%, and β is the shape of the distribution or Weibull slope (Gumbel, 1958; Wolters, 1958; Degraeve et al., 1998; Dumin et al., 1994; Suñé, 2001; Chen, 2014; Suñé et al., 2007; Vici et al., 2023; Wu et al., 2019; Wu, 2019b; Vollertsen and Kleppmann, 1990; Wu and Vollertsen, 2002). A critical characteristic of dielectric BD is its dependence on dielectric area, which was known a century ago (Peirce and Tensile Tests for Cotton Yarns, 1926). Nevertheless, physical understanding of area-dependent dielectric BD is a subject that remains nearly untouched or trivial to some. The view of the triviality of area-dependent BD is so permeative partly because most consider this phenomenon as a simple random process. This perspective fails to recognize that competition among BD events drives area-dependent dielectric BD (Suñé et al., 1990; Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024; Wu, 2026; Wu et al., 2026; Wolters, 1958; Degraeve et al., 1998; Dumin et al., 1994; Suñé, 2001; Chen, 2014; Suñé et al., 2007; Vici et al., 2023; Wu et al., 2019; Wu, 2019a; Vollertsen and Kleppmann, 1990; Wu and Vollertsen, 2002). This dependence traditionally follows an inverse area-scaling law as t B D ∼ A − η , where η = 1/β if and only if the Weibull model is applicable to describe the BD data.

A possible cause of this disregard is the lack of understanding that the origin of the traditional area-scaling law arises from the weakest link of the BD process, which is sometimes misinterpreted as “the winner takes all.” It is often not appreciated or not widely understood that, in fact, BD statistics involve a distribution within many underlying distributions; that is, many underlying near-BD events in competition ultimately control the outcome of the resultant extreme-value distribution (Gumbel, 1958; Wolters, 1958; Degraeve et al., 1998; Dumin et al., 1994; Suñé, 2001; Chen, 2014; Suñé et al., 2007; Vici et al., 2023; Wu et al., 2019; Wu, 2019a; Vollertsen and Kleppmann, 1990; Wu and Vollertsen, 2002). The weakest-link nature of dielectric BD has been directly demonstrated by an experimental grouping method of random BD events (Wolters, 1958; Vollertsen and Kleppmann, 1990; Wu and Vollertsen, 2002). Contrary to BTI and HCD, competition-driven area dependence is unique to dielectric BD, and its atomistic microscopic mechanisms remain unknown.

Furthermore, BD spots shown in Figure 1a appear and disappear during stress, eventually leading to a hard BD (HBD) event defined at a high BD current of 10–100 μA (Tsang and Linder, 2004). These results conclusively demonstrate the competing nature of the filament-forming BD process in HfO2. Recently, the reconstructed images from conductive atomic force microscopy (c-AFM) (Wei et al., 2023) vividly reveal multiple completed (89%) and incomplete (11%) filaments with various sizes in a 25 μm2 ReRAM device that directly confirm the notion of a fragmental state because of dielectric BD. Furthermore, after BD occurrence, electrical measurements in Figure 1b have revealed large current fluctuations by orders of magnitude, with current recovering to near-original levels prior to the final HBD. Until now, there has been no clear account for these recovery events in the post-BD stage, whereas the defect annihilation process proposed in our model can provide a physically plausible explanation for these effects in Figure 1b as discussed throughout this work.

**(a)** Hot-carrier emission images from a pFET device of SiO<sub>2</sub> at different times with t<sub>0</sub><t<sub>1</sub><t<sub>2</sub><t<sub>3</sub> (Tsang and Linder, 2004). Note that a BD spot in t<sub>1</sub> disappears later at t<sub>2</sub>. **(b)** Time-evolution of stress currents of stressed HfO<sub>2</sub>/SiO<sub>2</sub> samples from Wu et al. (2011). This adapted figure with permission from (Tsang and Linder, 2004).
Figure 1. (a) Hot-carrier emission images from a pFET device of SiO2 at different times with t0<t1<t2<t3 (Tsang and Linder, 2004). Note that a BD spot in t1 disappears later at t2. (b) Time-evolution of stress currents of stressed HfO2/SiO2 samples from Wu et al. (2011). This adapted figure with permission from (Tsang and Linder, 2004).

1.1 Overview of contemporary BD modeling

It might be surprising for many to appreciate the fact that over the century since the first realization of the area-dependent BD effect (Peirce and Tensile Tests for Cotton Yarns, 1926), there was only one key study (Suñé et al., 1990) using a physics-based percolation model that provides an explicit explanation for Weibull statistics and the resulting area-dependent dielectric BD. While the percolation approach (Suñé et al., 1990; Degraeve et al., 1998; Dumin et al., 1994; Suñé, 2001; Chen, 2014; Suñé et al., 2007) identifies defects as sole participants in BD filament formation, it was considered for lacking atomic insights (Vandelli et al., 2013; Larcher and Padovani, 2017), while we have not been able to find any other microscopic modeling work in literature to address this critical area-dependent BD process. On the contrary, to our best knowledge, all contemporary modeling work (Vandelli et al., 2013; Larcher and Padovani, 2017; Strand et al., 2022; Pavodomi, 2017; Padovani et al., 2023; Padovani et al., 2024; Guy et al., 2015; Li et al., 2012; Qin et al., 2015; Xu et al., 2020; Huang et al., 2013; Pan et al., 2011; Lv et al., 2014; Dirkmann et al., 2018; Sadi et al., 2018; Aldana et al., 2020; Guseinov et al., 2020; Raffone and Cicero, 2018; Min et al., 2021; Zeumault et al., 2022; Fyrigos et al., 2021) sidesteps or omits this important issue of area-dependent BD, as we will briefly discuss below.

A large body of theoretical work using different formulation schemes aims at microscopic insight into BD filament formation. A common scheme is the kinetic Monte Carlo (kMC) simulation of redox reactions and various conduction mechanisms, essentially defect-assisted tunneling, across a dielectric layer (Guy et al., 2015; Li et al., 2012; Qin et al., 2015; Xu et al., 2020; Huang et al., 2013; Pan et al., 2011; Lv et al., 2014; Dirkmann et al., 2018; Sadi et al., 2018; Aldana et al., 2020; Guseinov et al., 2020; Raffone and Cicero, 2018; Min et al., 2021; Zeumault et al., 2022; Fyrigos et al., 2021). Unfortunately, these modeling approaches mainly rely on the comparison with experimental IV measurements to justify their modeling approaches, rather than focusing on statistical distributions of BD filament formation, essentially avoiding the most critical issue of area-dependent dielectric breakdown.

Density functional theory (DFT) calculations have been widely used to investigate defect properties of dielectrics at atomic levels (Shluger and McKenna, 2013; Freysoldt et al., 2014; Gunkel et al., 2020; Chen and Kuo, 2011; Gao et al., 2016; Strand et al., 2017; Clima et al., 2014; Bradley et al., 2015; Traore et al., 2016; Guo and Robertson, 2014). Oxygen vacancies (VO) have been found both theoretically (Strand et al., 2022; Gunkel et al., 2020; Chen and Kuo, 2011; Gao et al., 2016; Strand et al., 2017; Clima et al., 2014; Bradley et al., 2015; Traore et al., 2016; Guo and Robertson, 2014) and experimentally (Li et al., 2008; Kumar et al., 2016) as a result of BD and are widely identified as defects responsible for BD filament formation and switching operations in ReRAM devices. Particularly, it was shown that two injected electrons are trapped in SiO2 (Gao et al., 2016), whereas four electrons are required to stabilize neutral oxygen vacancies in HfO2 (Strand et al., 2022; Bradley et al., 2015). Nonetheless, these proposed mechanisms based on DFT calculations must be thoroughly compared with experimental BD data because no other experimental means of verification exist, as suggested by Strand et al. (2022).

To our knowledge, the most sophisticated modeling work involving multiscale atomistic simulation (Vandelli et al., 2013; Larcher and Padovani, 2017; Strand et al., 2022; Pavodomi, 2017; Padovani et al., 2023; Padovani et al., 2024) using Ginestra software combines DFT-based defect formation energies with the kMC technique to model carrier transport and then to compare with experimental CDF data by extrapolating the theoretical HBD CDFs from small theoretical cell of 49 nm2 to the device area of 1 μm2 using t B D ∼ A − 1 / β. Despite the appearance of agreement with experimental BD data, the simulated progressive BD (PBD) time, aka residual time, is extremely long, on the order of ∼108 s (Strand et al., 2022). Based on successive BD statistics (Suñé and Wu, 2003) (Supplementary Appendix A), such a long PBD time would yield a BD spot density of ∼3400 nm−2, which is unphysically large.

In typical BD measurements, a high compliance current of 10–100 μA at an HBD event is commonly defined and used to limit further damage to a device. It is well established (Lombardo et al., 2005; Monsieur et al., 2002; Linder et al., 2002; Miranda et al., 1998; Suñé et al., 2006; Wu et al., 2019; Wu and Sune, 2006; Roussel et al., 2007; Raghavan et al., 2013; Kauerauf et al., 2011; Pompl et al., 2006; Kaczer et al., 2004; Kerber et al., 2007; Tous et al., 2008; Wu et al., 2007; Wu et al., 2012; Joshi et al., 2019; Vici et al., 2025; Sahcchi et al., 2026) that the entire BD process leading to the final HBD is a two-stage process: an initial stage of defect generation, and a post-BD stage. The initial stage occurs as the so-called first soft BD event (t1BD), marked by a discontinuity in the current time evolution. The post-BD stage can occur in two different modes: 1) a continuous growth of the same BD filament, so-called single spot (filament) progressive BD; 2) opening multiple filaments and simultaneously growing into a final HBD event. Over last 25 years, extensive work has been carried out to develop statistical methodologies (Lombardo et al., 2005; Monsieur et al., 2002; Linder et al., 2002; Miranda et al., 1998; Suñé et al., 2006; Wu et al., 2019; Wu and Sune, 2006; Roussel et al., 2007; Raghavan et al., 2013; Kauerauf et al., 2011; Pompl et al., 2006; Kaczer et al., 2004; Kerber et al., 2007; Tous et al., 2008; Wu et al., 2007; Wu et al., 2012; Joshi et al., 2019; Vici et al., 2025; Sahcchi et al., 2026) to characterize these two post-BD modes that have also been confirmed by electron scanning microscopy (Tsang and Linder, 2004; Monsieur et al., 2002) whereas HBD distributions have been shown to be non-Weibull (Suñé et al., 2006; Wu et al., 2019; Wu and Sune, 2006; Roussel et al., 2007; Raghavan et al., 2013; Kauerauf et al., 2011; Pompl et al., 2006; Kaczer et al., 2004; Kerber et al., 2007; Tous et al., 2008; Wu et al., 2007; Wu et al., 2012; Joshi et al., 2019; Vici et al., 2025; Sahcchi et al., 2026). This contradicts the simulated HBD distributions being Weibull-derived from multiscale modeling work (Padovani et al., 2023; Strand et al., 2022).

1.2 Thermodynamics of defect-generation–annihilation

The considerations in the previous section reveal a large gap between experimental BD measurements and ongoing BD modeling work that ignores the critical issue of area-dependent dielectric BD. A realistic microscopic understanding of dielectric BD can only be achieved by incorporating an area-dependent defect-generation–annihilation competition process. It must be pointed out that a thermodynamic framework for the generation and recombination of oxygen vacancies (VO) and interstitial oxygens (Oi), often called a Frenkel (or anti-Frenkel for ionic crystals) pair (schematically shown in Figure 2 is well established and discussed in a textbook by Maier (2004) and thoroughly reviewed by Menzel and Waser (2019), Dittmann et al. (2022), and Schie et al. (2018).

Schematic view of oxygen-vacancy–interstitial (V<sub>O</sub>–O<sub>i</sub>) pair generation in a dielectric. A hexagonal lattice is chosen for ease of visualization. The solid arrow is for the generation of the oxygen vacancy, whereas the dashed arrow represents a recombination of the vacancy and interstitial oxygen. The grey balls are for either hafnium or oxygen atoms.
Figure 2. Schematic view of oxygen-vacancy–interstitial (VO–Oi) pair generation in a dielectric. A hexagonal lattice is chosen for ease of visualization. The solid arrow is for the generation of the oxygen vacancy, whereas the dashed arrow represents a recombination of the vacancy and interstitial oxygen. The grey balls are for either hafnium or oxygen atoms.

The thermodynamics of VO–Oi pair generation in solids such as HfO2 closely resemble that of electron–hole pair generation in semiconductors, as shown in Figure 3a without a bias and in (b) with an applied bias. In the case of VO–Oi pair generation, the separation between two chemical potentials, μ O i and μ V O, is the Gibbs energy (ΔG) that plays the role of the energy gap (EG) in electron–hole pair generation in semiconductors. Interstitial oxygen can act as an annihilation center as a vacancy is readily recombined with it. DFT calculations have shown a recombination time of ∼500 fs for positively charged oxygen vacancies due to Coulombic attraction in intrinsic undoped HfO2 (Clima et al., 2014; Bradley et al., 2015). Therefore, annihilation centers can also affect the net defect-generation process in a competing manner during BD filament formation.

Chemical potential representation of oxygen vacancy, Vo, and oxygen interstitials, Oi, **(a)** without and **(b)** with a bias (Maier, 2004; Menzel and Waser, 2019; Dittmann et al., 2022; Schie et al., 2018). ΔG is the Gibbs energy for V<sub>O</sub>–O<sub>i</sub> pair generation.
Figure 3. Chemical potential representation of oxygen vacancy, Vo, and oxygen interstitials, Oi, (a) without and (b) with a bias (Maier, 2004; Menzel and Waser, 2019; Dittmann et al., 2022; Schie et al., 2018). ΔG is the Gibbs energy for VO–Oi pair generation.

It is worthwhile noting that in electroforming processes or switching operations, an oxygen-exchange reaction as an “extended Frenkel pair” can take place first and be followed by subsequent oxygen-vacancy migration via its recombination (annihilation) with neighboring oxygens (Menzel and Waser, 2019; Dittmann et al., 2022; Cooper et al., 2017). Extensive experimental evidence, including bubble formation and Ti oxidation, has been reported for SrTiO3-based TiO2 and SiO2-based and Ti/HfO2/TiN VCM ReRAM cells (Cooper et al., 2017; Mehonic et al., 2018; Szot et al., 2006; Yang et al., 2009; Sowinska et al., 2012) in support of this oxygen-exchange and oxygen-vacancy migration model (Menzel and Waser, 2019; Dittmann et al., 2022; Cooper et al., 2017; Mehonic et al., 2018; Szot et al., 2006; Yang et al., 2009; Sowinska et al., 2012; Traore et al., 2016). Thus, the defect annihilation process is not totally a foreign concept.

1.3 Key experimental observations of reversing and diminishing area dependence

Unfortunately, from a technological viewpoint, the traditional inverse area-scaling law is becoming one of the most critical bottlenecks for technology scaling using ReRAM chips and metal/insulator/metal capacitors (MIMCAP) applications in CMOS technology (Chen, 2014; Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024; Wu, 2026; Wu et al., 2026; Ando et al., 2025). Recently, we experimentally demonstrated that this area-scaling law can be reversed and diminished by process innovation through hydrogen plasma-treated HfO2. Figure 4a displays the flow for this new process (Consiglio et al., 2021; Kim et al., 2021a; Kim et al., 2021b; Clark et al., 2007). Compared to the traditional reservoir-driven process, this new process is intended to introduce oxygen vacancies and interstitials in single and bilayer stacks. ReRAM devices have been fabricated using this new process under various hydrogen exposure times. The forming voltages, that is, BD voltages (VBD) obtained from ramped voltage stress (RVS), are shown in Figure 4b over a variety of process conditions and single/bilayer stacks with two areas. Note the reverse and diminishing area dependence for samples A to D. In addition, crossbar arrays with smaller device areas of a 150 nm diameter have been fabricated to demonstrate for AI deep learning demonstration (Ando et al., 2025; Consiglio et al., 2021; Kim et al., 2021a; Kim et al., 2021b; Clark et al., 2007). In Figure 5, we show the TBD data obtained from constant stress voltage (CVS) tests for three typical area-scaling schemes in comparison with our physics-based statistical model by using the CDF analysis (Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024; Wu, 2026; Wu et al., 2026). However, due to the monotonic increasing nature of the CDF time-dependence, we were limited to examining the detailed physics of competition between different mechanisms such as defect generation and annihilation processes. Because area-dependent dielectric BD also involves various aspects, previous studies (Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024; Wu, 2026; Wu et al., 2026) concentrated primarily on the statistical features of our physics-based statistical model. The physical components that we aim to address in this article were not thoroughly discussed in our earlier reports (Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024; Wu, 2026; Wu et al., 2026).

**(a)** Process flow and schematics for the reservoir-driven approach using a reactive metal for the reference process and the hydrogen-plasma-driven approach in single-layer and bilayer structures. **(b)** Average V<sub>BD</sub> values for a variety of process splits (Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024; Wu, 2026; Wu et al., 2026) for 25 μm<sup>2</sup> and 10<sup>4</sup> μm<sup>2</sup>. The splits are marked with A, B, C, D, E, and F. BD is defined at 100 μA (Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024).
Figure 4. (a) Process flow and schematics for the reservoir-driven approach using a reactive metal for the reference process and the hydrogen-plasma-driven approach in single-layer and bilayer structures. (b) Average VBD values for a variety of process splits (Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024; Wu, 2026; Wu et al., 2026) for 25 μm2 and 104 μm2. The splits are marked with A, B, C, D, E, and F. BD is defined at 100 μA (Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024).
Comparison of the experimental BD data with a physics-based statistical model in the form of cumulative distribution fractions (CDF) obtained from Equation 7. **(a)** Reverse area scaling, **(b)** diminishing area scaling, and **(c)** forward area scaling.
Figure 5. Comparison of the experimental BD data with a physics-based statistical model in the form of cumulative distribution fractions (CDF) obtained from Equation 7. (a) Reverse area scaling, (b) diminishing area scaling, and (c) forward area scaling.

Instead, we will use a probability density function (PDF) analysis, that is, probability density per unit time, to quantitatively unravel the competition between defect generation and annihilation in order to develop a deep physical understanding of these new findings in a holistic and comprehensive picture that can unify various area-dependent BD schemes. This article is structured as follows. In Section 2, we will give a brief overview of our physics-based statistical model with an extensive discussion of interplay between defect generation and its annihilation process compared to the classic defect-generation-only scheme (Suñé et al., 1990; Degraeve et al., 1998; Dumin et al., 1994; Suñé, 2001; Chen, 2014; Suñé et al., 2007; Vici et al., 2023; Wu et al., 2019; Wu, 2019a). In Section 3, the PDF analysis will be presented with various features of defect-generation–annihilation interactions from reverse to diminishing BD area-scaling schemes. This section will also cover the PDF analysis from filament areas at nanoscale to device area at micrometer scale, eventually linking to different variabilities of area-dependent BD phenomena. In Section 4, we will discuss the physical roles of defect generation and annihilation involved in the BD process by either an intrinsic recombination mechanism or an extrinsic mechanism, such as hydrogen-related vacancy passivation. Finally, we will summarize our key findings, highlighting that dielectric BD leaves pristine dielectrics in a fragmented state that has no resemblance to any other degradation phenomena, such as BTI and HCD.

2 Defect-generation–annihilation statistic model

Figure 6 displays a schematic of defect generation and defect annihilation toward a BD event. Figures 6a,f represent the traditional BD picture, in which a continuous defect generation directly leads to the formation of a conducting path without any intermediate steps (c)–(e). In contrast to this conventional view, in competition with defect generation, the VO annihilation process can also occur by either intrinsic processes such as recombination of VO–Oi pairs or extrinsic means such as substitutional hydrogens into VO, as shown in Figures 6b,c. The defect annihilation can passivate some VOs, yielding incomplete filaments as in Figures 6d,e. This constant competition continues until a filament forms, as seen in Figure 6f. This entire physical process described in Figure 6 can be mathematically captured by extending the cell-based percolation model (Suñé et al., 1990) to include both defect-generation and annihilation processes, as shown in Figure 7. In this model, defects and annihilation centers are randomly distributed across a device area.

Schematic drawings of interplay of V<sub>O</sub> generation and annihilation in evolution into partial filaments, eventually evolving into a final BD filament. **(a)** Represents V<sub>O</sub> generation directly to BD **(f)** without annihilation. **(b–f)** represent the interplay and competition between V<sub>O</sub> generation **(a)** and its annihilation, leading to partial and complete BD events. See text for discussion. The dashed circles are for a V<sub>O</sub>, whereas the grey balls with solid lines are for O<sub>i</sub>, and the small red ball stands for a hydrogen atom or molecule. The grey balls are for either Hf or oxygen atoms.
Figure 6. Schematic drawings of interplay of VO generation and annihilation in evolution into partial filaments, eventually evolving into a final BD filament. (a) Represents VO generation directly to BD (f) without annihilation. (b–f) represent the interplay and competition between VO generation (a) and its annihilation, leading to partial and complete BD events. See text for discussion. The dashed circles are for a VO, whereas the grey balls with solid lines are for Oi, and the small red ball stands for a hydrogen atom or molecule. The grey balls are for either Hf or oxygen atoms.
Cell-based percolation model incorporating both defect-generation and annihilation processes. The red box is for a randomly formed filament at BD. *a* <sub>0</sub> is the area of the unit cell for any one column, and t<sub>diel</sub> is the dielectric thickness.
Figure 7. Cell-based percolation model incorporating both defect-generation and annihilation processes. The red box is for a randomly formed filament at BD. a 0 is the area of the unit cell for any one column, and tdiel is the dielectric thickness.

The time evolution of the average generated defect density ρ G t and annihilation-center density ρ A t are modeled by a power law (Degraeve et al., 1998; Vici et al., 2023) as ρ G = N G t G α G and ρ A = N A t A − α A as plotted in Figure 8. Then, the time-dependence of the average numbers of defect and annihilation center can be readily obtained as λ G = tdiel a 0ρG and λ A = tdiel a 0ρA. For the case of reverse area scaling, Figure 8 shows that defect density increases with time, while annihilation-center density decreases with time as more annihilation centers diminish due to recombination. Eventually, the magnitude of defect generation exceeds that of the annihilation center at later times, so net defect generation dominates BD filament formation.

Time-dependent average density of generated defects, ρ G t = N G t α G , and annihilation centers, ρ A t = N A t − α A, used in the defect-generation–annihilation statistical model. The parameters are given in Supplementary Table C1 of Supplementary Appendix C.
Figure 8. Time-dependent average density of generated defects, ρ G t = N G t α G , and annihilation centers, ρ A t = N A t − α A, used in the defect-generation–annihilation statistical model. The parameters are given in Supplementary Table C1 of Supplementary Appendix C.

The time-dependent PDFs, f u t G f u t A, for generated defects and annihilation centers of the unit cells (Figure 7) can be obtained by taking the derivatives, f t = d F t / d t, of their CDFs, F u t G and F u t A . Their expressions are given below: where N CG and N CA are the critical numbers of generated defects and annihilation centers, respectively. Any one of these columns represents a BD event as generated defects and annihilation centers in a column reach a critical number.

F u t G = 1 ‐ ∑ i = 0 N CG ‐ 1 λ G i exp ‐ λ G i ! , (1)

F u t A = ∑ i = 0 N CA ‐ 1 λ A i exp ‐ λ A i ! , (2)

To transform f u t G of a unit cell to the PDF, f dev t G, of a device area, we apply minimum-value statistics in Equation 3 for defect generation as done by Suñé et al. (1990). For annihilation processes, maximum-value statistics are applied for transformation of f u t A to obtain f Dev t A, as given in Equation 4, where N = A dev / a 0 and A dev is the device area.

f dev t G = N f u t G 1 − F u t G N − 1 , (3)

f dev t A = N f u t A F u t A N ‐ 1 , (4)

The PDF for the net defect generation can be readily obtained as given below:by integrating over the annihilation time into the joint PDF of Equation 6, given as

f Net def t G = ∫ 0 ∞ f t G , t A d t A (5)

f t G , t A = f d e v t G f d e v t A . (6)

Then, the overall CDF of a device area for the comparison with experimental data is given in Equation 8,where Z is the normalization constant and can be obtained by matching the T BD data. For simplicity, the PDFs of defect generation and defect annihilation are assumed to be independent as discussed below. Figure 5 shows good agreement between experimental data and the modeling results of Equation 8 for all different areas using a single set of parameters (Supplementary Table C1 of Supplementary Appendix C) for various splits, demonstrating both the validity and predictability of this model because each parameter set allows us to simultaneously fit multiple data sets of different areas, showing good agreement.

F t BD = Z ∫ 0 t B D f Net def t G d t G , (7)

2.1 Assumptions and validity of the parameters used in the defect-generation–annihilation model

It is worthwhile to discuss the assumptions and validity of the parameters used in modeling the defect annihilation process. First, the power-law time-dependence of ρ G = N G t G α G has been well established for the defect-generation process in corroboration with experimental measurements of stress-induced leakage currents (SILC) (Degraeve et al., 1998; Vici et al., 2023). Due to a lack of knowledge of a functional form of the time-dependence of the defect annihilation process in the literature, we also adopted this power-law dependence for defect annihilation as ρ A = N A t A − α A with a negative exponent of − α A. The time-decaying function for defect annihilation is physically plausible, as the number of annihilation centers is expected to decrease over time because they are eliminated during the ongoing annihilation process. This further corroborates the fact that this process, along with Equation 2 for f u t A, is naturally compatible with the maximum-value statistics we introduce in Equation 4 for modeling area-transformation from f u t A to f d e v t A. In other words, both the time-decaying function and the maximum-value statistics would yield a physically meaningful positive f d e v t A, rather than a negative f d e v t A that is physically meaningless. This demonstrates that the time-decaying function for defect annihilation is completely compatible with maximum-value statistics. Thus, they provide a natural framework for experimental observations of reverse or diminishing area scaling of dielectric BD shown in Figures 4, 5. Furthermore, this framework is fully consistent with the widely accepted thermodynamic description of defect-generation–annihilation (Maier, 2004; Menzel and Waser, 2019; Dittmann et al., 2022; Schie et al., 2018). Thus, our model represents a unique advance by incorporating BD area dependency into a fundamental understanding of dielectric BD since the work of Suñé et al. (1990).

Third, to obtain the net defect-generation PDF of Equation 5, the two processes of defect generation and defect annihilation are assumed to be independent in Equation 6. This assumption is justified because we consider at least two different sources for annihilation of oxygen vacancies: intrinsic oxygen interstitials and extrinsic species like hydrogens. Because both reverse and diminishing area-scaling dependences are found on the HfO2 bilayer stacks with hydrogen plasma treatment shown in Figure 4b, it is physically reasonable to expect that interstitial oxygen and hydrogen as extrinsic species act independently from each other. This consideration is consistent with the DFT results that show hydrogen can annihilate oxygen vacancies as substitute hydrogen (Sklénard et al., 2023; Magyari-Köpe et al., 2016; Duncan et al., 2016; Kaviani et al., 2017; Jung et al., 2017). To our best knowledge, existing literature lacks quantitative descriptions of defect annihilation. Thus, assuming independence between f D e v t G and f D e v t A is a reasonable first-order approximation adopted in our model.

As for the parameters used in our model, NCG is the critical defect number used previously by Suñé et al. (1990). However, neither NCG nor NCA is an experimentally accessible quantity. Here, we assume NCA to be equal to NCG for simplicity, so that we retain a single parameter NC. This should not be mistaken as meaning that there are equal concentrations of defects and annihilation centers. Note that the introduction of defect annihilation is not intended to replace the generation process; rather, it is intended to modulate the overall defect generation in a manner leading to a variety of experimental BD area-dependences in Figure 5. Overall, our model contains six parameters: a 0, NC, NG, αG, NA, and αA.

Using the reliable experimental technique of isotope exchange annealing, the range of vacancy density is determined to be between 2.8 × 1017 cm−3 and 40 × 1017 cm−3 by a secondary ion mass spectrometry (SIMS) study (Mueller and De Souza, 2018) and further supported by the results of both DFT calculations (Clima et al., 2014) and molecular dynamics (MD) (Schie et al., 2017). This is in good agreement with the range of NG values in Supplementary Table C1 that we used to match the experimental data in Figure 5. The value of αG from our work for defect generation is smaller than the previously estimated 0.4 based on SILC measurements (Degraeve et al., 1998; Vici et al., 2023) for the traditional BD area-scaling law. Thus, this value corresponds to a situation in which the annihilation process is absent. These discrepancies can contribute to the smaller value obtained here that involves defect generation in bilayer HfO2 stacks. In addition, the values of NA and αA are obtained by the best fitting of the experimental data. In general, we found that the NA values are larger than those of NG. This is reasonable because we consider two sources of annihilation centers: intrinsic Oi and extrinsic hydrogen, as we will discuss further. Finally. 25 nm2 for a 0 is chosen for the column area in Figure 7 based on a range of 9–100 nm2 based on c-AFM images of BD spots (Hou et al., 2016).

3 Analysis of PDF results and comparison with experiments

Now we will unravel new understandings for the underlying physics and statistics between defect-generation and annihilation mechanisms by examining various PDFs of defect generation of Equations 1, 3 and defect annihilation of Equations 2, 4 in conjunction with the net defect-generation PDF (Equation 5) to compare with CDF data in Figure 4.

3.1 Reverse area scaling at fixed stress voltages

Figure 9 exhibits the relative PDFs of defect generation and annihilation using the derivatives of Equations 1, 2 for a unit cell and using Equations 3, 4 for a device area of 100 μm2. While the PDF of defect generation shows a monotonic decreasing trend with time, like the one in Supplementary Figure B1a for β = 0.7 in Supplementary Appendix B, the defect annihilation PDF increases and then decreases with increasing time. This means that at shorter times, generated defects are quickly annihilated so that they cannot contribute to BD. This process reaches a peak at ∼0.1 μs and is longer than the intrinsic recombination time (∼500 fs) in HfO2 (Clima et al., 2014; Bradley et al., 2015).

Competition-driven PDFs for sample B in Figure 4. **(a)** PDFs of a unit cell and 100 μm<sup>2</sup> using Equations 3, 4 for the PDFs of defect generation and annihilation, respectively. **(b)** The PDFs of three areas to compare with the net defect generation of Equation 5. The grey area in Figure 9a represents an experimental time window. Note that the decline of the PDF of defect annihilation corresponds to the peak of net defect generation.
Figure 9. Competition-driven PDFs for sample B in Figure 4. (a) PDFs of a unit cell and 100 μm2 using Equations 3, 4 for the PDFs of defect generation and annihilation, respectively. (b) The PDFs of three areas to compare with the net defect generation of Equation 5. The grey area in Figure 9a represents an experimental time window. Note that the decline of the PDF of defect annihilation corresponds to the peak of net defect generation.

On the other hand, this time coincides with the estimated times from 1.6 ns to 1.6 μs for substituting hydrogen (HO) into VO based on the formation energy of a HO–VO swap using the DFT calculation (Sklénard et al., 2023) as discussed in Section 4. Beyond the peak, the PDFs for defect annihilation gradually decrease with time. The rate of defect annihilation is much steeper than that of defect generation. This suggests that generated defects start to survive or are not recombined by annihilation centers so that they can contribute to BD filament formation. Note also that the crossovers between the PDFs of defect generation and annihilation are independent of area, consistent with the crossover in ρG(t) and ρA(t) in Figure 8 for sample B.

In Figure 9a, the PDF of defect annihilation at a device area drops drastically in magnitude from the PDF at a unit cell and shifts to longer times, as expected from maximum-value theory (Gumbel, 1958) because more annihilation centers are available in a larger area, causing more defects to be annihilated. In contrast to the small difference in the PDFs for defect generation, a drop of approximately seven orders of magnitude in the defect annihilation PDF is much larger due to the large difference in N A and N G values shown in Figure 8 and in Supplementary Table C1. This suggests that a large imbalance in densities of generated defects and annihilation centers is a dominant factor responsible for reverse area scaling.

In contrast, competition among various generated defects also occurs and causes an increase in both time and magnitude for the PDFs in a large device area, as seen in Figure 9a. This process also competes with their own survival times that are affected by the presence of annihilation centers. As stress continues to late times, the PDF of defect generation eventually exceeds that of defect annihilation so that generated defects can survive and contribute to filament formation. Although this crossover occurs with a small probability of ∼10−15 at ∼103 s at the unit cell, the use of a large device area increases its probability by five orders of magnitude, so that it is experimentally observable in Figure 5a.

In Figure 9b, we compare the PDFs of the net defect generation with the individual PDFs of defect generation and annihilation. First, it is seen that the peaks of the PDFs of the net defect generation coincide with the peaks of the PDFs of defect annihilation for all three areas. Second, despite the slow increase in defect generation, it ultimately exceeds the fast decline of annihilation. Third, these peaks also coincide with the onset of BD data found experimentally, as shown in Figure 5a. This suggests that the decline in the PDFs of the annihilation process is directly responsible for the net defect generation. In other words, the annihilation process becomes less effective so that more defects can survive and thus participate in BD. Hence, the PDFs of defect generation exceed those of annihilation, resulting in a reverse BD area dependence. Note that both slopes (or exponents) of the PDFs of defect generation and annihilation are negative and much larger, while a positive slope in the PDF vs. time relation for the annihilation process would cause defects to disappear so that BD can be slowed or even prevented.

Figure 10a shows another example (sample C in Figure 4) of the reverse area-scaling trend. In contrast to Figure 9, the defect-generation PDF in Figure 10b shows a positive slope, whereas the defect annihilation PDF is like those of sample B. This positive defect generation rate increases the net defect-generation rate and results in a weaker reverse area dependence than in other cases (Figure 11). Note that the peaks of PDFs due to the annihilation process coincide with the onset of increase in the CDFs of three areas in Figure 10a, further facilitating that the decline in annihilation PDFs is largely responsible for the effectiveness of defect generation and causing reverse BD area scaling. It must be emphasized that the above discussion is based on the consideration of separate contributions of various PDFs. In other words, the net PDF of defect generation must consider the contributions from all the annihilation times with the proper normalization in comparison with the experimental CDF results because effects of defect generation and annihilation on the net defect generation are inherently intertwined.

Area-dependent dielectric BD of sample C. **(a)** Comparison of the model is specifically defined by Equation 7 and CDF T<sub>BD</sub> data. **(b)** Competition-driven PDFs for defect generation and defect annihilation of Equations 3, 4 compared to the net defect generation.
Figure 10. Area-dependent dielectric BD of sample C. (a) Comparison of the model is specifically defined by Equation 7 and CDF TBD data. (b) Competition-driven PDFs for defect generation and defect annihilation of Equations 3, 4 compared to the net defect generation.
Joint PDFs of Equation 6 vs. generation and annihilation times for three areas of two splits with 2 nm/3 nm stacks showing reverse area scaling. Hydrogen exposure time of 50 s **(a)** and 40 s **(b)**. Note that the PDF changes slope along the direction of generation time (see the text for discussion).
Figure 11. Joint PDFs of Equation 6 vs. generation and annihilation times for three areas of two splits with 2 nm/3 nm stacks showing reverse area scaling. Hydrogen exposure time of 50 s (a) and 40 s (b). Note that the PDF changes slope along the direction of generation time (see the text for discussion).

To gain more insight into the overall interaction and competition between defect-generation and annihilation processes, we examine their overall effects on the joint probability density function (joint PDF) of Equation 5 as seen in Figure 12. In Figure 12a, the magnitude in the PDFs decreases in the time-axis of defect generation, whereas the magnitude in the PDFs in Figure 12b increases with time, due to the negative and positive rates of the net defect generation shown in Figures 9, 10, respectively. Nevertheless, in both cases, the peaks of the joint PDF decrease in magnitude but shift to longer times with increasing areas in the time-axis of annihilation, leading to remarkable reverse area-scaling effects in the CDF TBD data found experimentally.

Voltage–dependent effect on BD CDFs and PDFs for sample B. **(a)** Comparison of CDFs between the model of Equation 7 and TBD data. **(b)** Comparison of individual PDFs of defect generation and annihilation with those of the net defect generation. Note that the peaks in net defect generation increase with voltage.
Figure 12. Voltage–dependent effect on BD CDFs and PDFs for sample B. (a) Comparison of CDFs between the model of Equation 7 and TBD data. (b) Comparison of individual PDFs of defect generation and annihilation with those of the net defect generation. Note that the peaks in net defect generation increase with voltage.

3.2 Voltage effects for reverse area scaling

Figure 13a plots the CDFs at three stress voltages for the device area of 104 μm2 showing the expected trend of voltage dependence as in the case of traditional BD area scaling. The same voltage dependence was also found using a smaller area of 900 μm2 (Wu et al., 2025). These voltage effects are reproducible for a different bilayer stack of sample A, which also shows reverse area-scaling behavior (Wu et al., 2025). The modeling results for CDFs also nicely capture the experimental TBD data, as seen in Figure 13a.

Competition-driven PDFs for sample E showing the forward area scaling in Figure 5c. **(a)** PDFs of a unit cell using the derivative of Equations 1, 2 and of 100 μm2 using Equations 3, 4 for the PDFs of defect generation and annihilation, respectively. **(b)** The PDFs of three areas in comparison with that of the net defect-generation Equation 5. Note that the PDFs of defect generation are generally higher than those of defect annihilation for all three areas in the experimental window, as shown in (b).
Figure 13. Competition-driven PDFs for sample E showing the forward area scaling in Figure 5c. (a) PDFs of a unit cell using the derivative of Equations 1, 2 and of 100 μm2 using Equations 3, 4 for the PDFs of defect generation and annihilation, respectively. (b) The PDFs of three areas in comparison with that of the net defect-generation Equation 5. Note that the PDFs of defect generation are generally higher than those of defect annihilation for all three areas in the experimental window, as shown in (b).

To match experimental CDF data, we found that the annihilation-center density, NA, drops with increasing voltage, as seen in Supplementary Table C1 in Supplementary Appendix C. In addition, we also found that a change in defect-generation density, NG, does not significantly alter the overall fitting results, suggesting that annihilation centers play a dominant role in reversing traditional BD area scaling. For intrinsic and undoped HfO2, oxygen vacancies and interstitial oxygens are generated in pairs so that their concentrations should always be equal. In our bilayer stacks, one of the HfO2 layers is exposed to H2 plasma treatment. It is known that hydrogen can easily reside in Vo sites to form a substitutional hydrogen defect, HO (Sklénard et al., 2023; Magyari-Köpe et al., 2016; Duncan et al., 2016; Kaviani et al., 2017; Jung et al., 2017). These considerations may suggest the imbalance between concentrations of generated defects and annihilation centers is so large that more annihilation centers are pushed out by a higher applied bias, as revealed above.

3.3 Forward area scaling

The TBD area dependence of a single-layer stack with an H2 time of 30 s of sample E in Figure 4c shows an appearance of the inverse area-scaling behavior, t B D ∼ A − η, similar to those seen in the traditional area scaling. However, the reciprocal value of the exponent does not yield the value of the Weibull slope, that is, η ≠ 1 / β, as required by the Poisson area-scaling rule (Degraeve et al., 1998; Dumin et al., 1994; Suñé, 2001; Chen, 2014; Suñé et al., 2007; Vici et al., 2023; Wu et al., 2019; Wu, 2019a; Vollertsen and Kleppmann, 1990; Wu and Vollertsen, 2002). To examine this effect further, we compare the individual PDFs of defect generation and annihilation for both the unit cell and the device area in Figure 14a, while the corresponding comparison to the net defect generation is given in Figure 14b. The PDFs of defect generation in Figure 14a increase and then decrease with time, like those for β > 1 in Supplementary Figure B1a of Supplementary Appendix B.

Competition-driven PDFs for sample E showing the diminishing area scaling in Figure 5B. **(a)** PDFs of a unit cell using the derivative of Equations 1, 2 and 100 μm2 using Equations 3, 4 for the PDFs of defect generation and annihilation, respectively. **(b)** PDFs of three areas compared with that of the net defect generation. Note that the slight decrease in the PDF of defect annihilation combined with the sharp drop in the PDF of defect generation gives rise to diminishing areadependent PDFs for the net defect generation.
Figure 14. Competition-driven PDFs for sample E showing the diminishing area scaling in Figure 5B. (a) PDFs of a unit cell using the derivative of Equations 1, 2 and 100 μm2 using Equations 3, 4 for the PDFs of defect generation and annihilation, respectively. (b) PDFs of three areas compared with that of the net defect generation. Note that the slight decrease in the PDF of defect annihilation combined with the sharp drop in the PDF of defect generation gives rise to diminishing areadependent PDFs for the net defect generation.

Figure 14b displays the voltage dependence of PDFs of net defect generation compared with the PDF of annihilation centers, showing the coincidence of their peaks with the onset of the experimental CDFs, following a consistent picture with other cases of the CDF TBD data presented earlier.

Moreover, the crossovers in the individual PDF of defect generation and the net defect-generation PDF can be seen with a change in device areas. For only defect generation, a change in PDFs is of more than two orders of magnitude at shorter times for a change in area ratio of ∼200× in Figure 14b. In contrast, a change in the PDF of the net defect generation is only 10× at shorter times. This difference is caused by interference with the defect annihilation process. Note that in short times, the slopes of net defect-generation PDFs are smaller and negative than the positive slopes of the defect-generation PDF. In general, the magnitude of defect-generation PDFs is much higher than that of defect annihilation PDFs. Their differences become larger with increasing time, suggesting that the annihilation process plays a less important role than in the case of reverse area scaling. Nonetheless, the presence of defect annihilation slows the defect-generation process so that the shape of the PDF of the net defect generation follows the defect annihilation PDF rather than the defect-generation PDF.

It is worth mentioning that VO is often more intentionally generated via an oxygen-exchange layer (Menzel and Waser, 2019; Dittmann et al., 2022; Bradley et al., 2015) for ReRAM, resulting in V O • • without a paired interstitial oxygen, and both cases may coexist. A very high concentration of V O • • without interstitial pairs is expected (TiO2 is formed instead), and a concentration of annihilation centers created by hydrogen plasma may be too low to reverse the trend, resulting in very low forming voltages with forward area scaling.

3.4 Diminishing area scaling

In the above sections, we have observed that the defect-annihilation process has a significant impact on BD area-scaling behavior to the extent of reversing the decades-long inverse area-scaling law. A diminishing BD area scaling as shown in Figure 5b would suggest that devices with different areas can have the same BD reliability, thus greatly providing a pathway to overcome the technology scaling limit for the semiconductor industry (Wu and Sũné, 2005; Stathis and DiMaria, 1998). According to traditional Weibull statistics for the inverse area scaling, t B D ∼ A − 1 / β, the diminishing BD area dependence would suggest that BD distributions become extremely steeper, that is, approaching an infinitely large β value. In practice, for a β value of 8, a change in CDF from 5% to 95% only leads to a small change of ∼1.5× in time, whereas the diminishing BD area dependence spans more than five orders of magnitude in time, as seen in Figure 5b. This implies other mechanisms such as annihilation also play a role in this case. To understand the interaction of defect generation and annihilation for the diminishing area dependence, we plot various PDF contributions from defect generation and annihilation for the unit cell and a device area of 100 μm2 in Figure 15a. Despite a huge change of 4 × 106 in the area ratio between these two areas, the two PDFs for defect generation only shift by ∼10×, indicating a much weaker area dependence. On the other hand, the shift in the two corresponding PDFs of defect annihilation is approximately six orders of magnitude over time. As a result, the effect of the annihilation process becomes more pronounced at the device level.

Normalized TBD values at 10 s and 104 μm2 for samples showing reverse, diminishing, and forward area scaling. Note the distinction between traditional and forward area scaling.
Figure 15. Normalized TBD values at 10 s and 104 μm2 for samples showing reverse, diminishing, and forward area scaling. Note the distinction between traditional and forward area scaling.

Figure 15b shows the individual PDFs (defect generation vs. annihilation) of different areas compared with the PDFs of their net defect generation. In general, at short times prior to these peaks, the rates of increases in the individual PDFs of defect generation are comparable to those of the net defect generation, whereas beyond these peaks, the rates in the PDFs of individual defect annihilation control the net defect generation. This means that the higher relative PDFs of defect annihilation prolong the defect-generation process, thus leading to a wide spread of BD times for filament formation at the same time without apparent area dependence, as observed in Figure 5B. In intrinsic and undoped HfO2, there is only a unique lifetime of an oxygen vacancy controlled by the recombination time of its counterpart, oxygen interstitial, for a single charge state (Clima et al., 2014; Bradley et al., 2015). In other words, there is a single quantity of recombination lifetime shared by both VO and Oi. However, in doped HfO2, hydrogens can also participate in the process of VO being annihilated; thus, more channels are available for defect annihilation, and most likely also for defect generation. Thus, our current understanding of defect generation or its recombination must be modified to reflect this reality.

3.5 Relationship between the variability of BD area dependence and the net defect-generation rate

In Figure 11, we show the normalized TBD values as a function of area for a wide range of BD area-scaling schemes from reverse area-scaling of three different samples (A, B, and C) to the diminishing area dependence (sample D) and then to forward area scaling (sample E) and finally traditional inverse area-scaling law (sample F). These results clearly show that the hydrogen exposure time and bilayer stacks influence the sensitivity of BD area dependence. To establish an overall connection with the defect-generation–annihilation model, we define the following quantity called the variability of the TBD-vs-area relationship, that is, T B D ∼ A κ,with κ = − 1 / β for the traditional inverse area scaling using the Weibull model.

κ = ∂ Ln t B D ∂ Ln A , (8)

Due to the physical connection between defect generation and annihilation centers, either intrinsic interstitial oxygen or extrinsic species like hydrogens, a more natural or complete picture of dielectric BD is that defects can physically be annihilated rather than being permanently created as previously understood. We have shown that the very nature of defect annihilation described by maximum-value statistics is compatible with the well-established minimum-value statistics for defect generation. Thus, it is natural to examine any correlation between the variability of the TBD-vs-area relation as defined in Equation 8 and the net defect generation as defined in Equation 5. As discussed above, the net defect-generation PDFs over the experimental time window in Figures 9b, 10b, 14b can be well fitted by a power-law function. Thus, we define a quantity called the net defect-generation rate, that is, the exponent of a power-law fit of each case as

χ = ∂ Ln f N e t _ ⁡ def ∂ Ln t . (9)

In the case of defect generation only using Weibull statistics, it can be shown that χ = β − 1 for CDF well below 63.2%. Figure 16 plots the experimental BD variability (κ) for five different sets of the t<sub>BD</sub>-vs-area data in Figure 11 versus their corresponding net defect-generation rates, χ. Here, the χ values are the averages over various areas for each split.

Variability (κ) defined in Equation 8 of BD area dependence plotted against the net defect-generation rate (χ) defined in Equation 9, which includes both defect generation and annihilation.
Figure 16. Variability (κ) defined in Equation 8 of BD area dependence plotted against the net defect-generation rate (χ) defined in Equation 9, which includes both defect generation and annihilation.

For the case of forward area scaling in Figure 13b, the PDF of net defect generation shows two χ values of −0.495 and −1.129 at short and long times, respectively; an average value of −0.812 is used in Figure 16. While the sign of χ separates the reverse area scaling from forward or traditional inverse area scaling, the magnitude of the net defect-generation rate (χ) offers a measure of the effectiveness of defect annihilation in affecting the net defect generation. A remarkable correlation can be seen in Figure 16, demonstrating the capability of our defect-generation–annihilation statistical model to capture various experimental BD area-scaling schemes in Figure 11.

4 Discussion

In Section 3, using PDF analysis, we investigated statistics and physics of area-dependent dielectric BD. To our knowledge, in nearly all literature on dielectric BD, defects such as VO are the sole contributors in consideration of breakdown mechanisms. It is important to recognize that an oxygen vacancy cannot exist by itself, while the role of its counterpart, oxygen interstitial, is often ignored. In this regard, the dynamic nature of the creation of VO–Oi pairs should be properly considered, as a thermodynamic description has been developed for some time (Maier, 2004; Menzel and Waser, 2019; Dittmann et al., 2022; Schie et al., 2018). For intrinsic undoped HfO2 with nearly spontaneous recombination (Clima et al., 2014; Bradley et al., 2015), one may conclude oxygen-vacancy generation would be almost impossible. On the other hand, doping is known to significantly impact forming voltages and ReRAM switching (Sklénard et al., 2023; Magyari-Köpe et al., 2016; Duncan et al., 2016; Kaviani et al., 2017; Jung et al., 2017). In our ReRAM devices with HfO2 bilayer stacks treated with hydrogen plasma, hydrogen is expected to play a dominant role. In amorphous solids, there is an almost infinite number of mechanisms for drift of hydrogen through a HfO2 network.

One critical mechanism revealed by the DFT calculation (Sklénard et al., 2023) is the migration of a hydrogen atom in a HO defect into an oxygen vacancy, that is, the hydrogen and oxygen vacancies switch their positions. It was found that the activation energy for the swap between HO and VO requires a lower barrier of ∼0.5 eV, implying that this process can occur at a rapid scale (Sklénard et al., 2023). Thus, using this result (Sklénard et al., 2023), we can estimate the hydrogen hopping time by calculating the attempt frequency via ν = ν 0 ⁡ exp ⁡ − E A / k B T with ν 0 = 10 TeraHz , where ν 0 is the typical phonon frequency for HfO2 (Zhao and Vanderbilt, 2002). For estimation of filament temperature, a range of values from 350 K to 600 K have been reported in the literature for HfO2 and Ta2O5-x (Rodriguez-Fernandez et al., 2018; Sharma et al., 2014). Using these two temperature values, hydrogen hopping times can vary from 1.6 μs to 1.6 ns, respectively. This range reasonably corresponds to the annihilation times in the PDF with a peak time at ∼0.1 μs, as shown in Figure 9a used in our model. Moreover, in SiO2, the passivation of oxygen vacancies involving hydrogen, such as H2–VO and H2O–VO defects revealed by the DFT calculations, can take place as hole trapping (Freire et al., 2024). These results suggest there might be additional mechanisms as annihilation centers for passivation of oxygen vacancies.

Regarding the multiscale simulation work using Ginestra software (Strand et al., 2022), showing extremely long residual time (or progressive BD time) of 108 s, this may be caused by the fact the formation time of the first filament (t1st_BD) is extremely short, at ∼100 s, because the residual time (tRes) is the difference between and the hard breakdown time, (tHBD), that is, tRes = tHBD − t1st_BD. As we pointed out earlier, contemporary BD models only consider the generation of immobile defects and the defect-assisted tunneling mechanism without any competition effects. In other words, the formation time of the first filament, t1st_BD, can be delayed by the defect-annihilation process if it were included in the multiscale simulation, thus resulting in a more reasonable residual time. This illustrates that the defect annihilation process can still play a role in the case of traditional inverse area-scaling BD.

A main challenge in dielectric BD detection is that it almost exclusively relies on measurements of leakage currents passing through a dielectric layer. Thus, understanding of leakage current conduction through defect states is very important. Figure 17a shows the formation energies of various oxygen vacancies of an intrinsic HfO2 dielectric, showing that V O • • has the lowest formation energies over a ∼3.5 eV span from the valence-band edge. On the other hand, Figure 17b displays the HfO2 band diagram with the TiN electrode material under a bias, along with various charge states of various oxygen vacancies. It is evident that the closest defect levels to the cathode interface are V O • •, whereas the defect levels of oxygen vacancies of other charge states across the bandgap of HfO2 are located further away from the interface.

**(a)** Formation energies of various oxygen vacancies of an intrinsic HfO<sub>2</sub> (Strand et al., 2017). This adapted figure with permission from (Strand et al., 2017). **(b)** Energy band diagram of intrinsic HfO<sub>2</sub> of different charge states with a TiN cathode at the interface (Strand et al., 2022). This adapted figure with permission from (Strand et al. 2022).
Figure 17. (a) Formation energies of various oxygen vacancies of an intrinsic HfO2 (Strand et al., 2017). This adapted figure with permission from (Strand et al., 2017). (b) Energy band diagram of intrinsic HfO2 of different charge states with a TiN cathode at the interface (Strand et al., 2022). This adapted figure with permission from (Strand et al. 2022).

The DFT calculation reveals a spontaneous recombination time of ∼500 fs for V O • •, whereas the DFT calculation assuming four electrons injected into HfO2 can generate neutral oxygen vacancies. Unfortunately, no additional DFT results are available for other charge states. These two results probably represent two extreme cases for the dynamics of these different oxygen vacancies. Nevertheless, they clearly demonstrate the dynamical nature of generated defect states in competition, leading to competition-driven dynamical conduction channels. It is unclear whether and how this dynamical nature of these conduction mechanisms is included in ongoing BD modeling approaches either by an exclusive kMC approach (Guy et al., 2015; Li et al., 2012; Qin et al., 2015; Xu et al., 2020; Huang et al., 2013; Pan et al., 2011; Lv et al., 2014; Dirkmann et al., 2018; Sadi et al., 2018; Aldana et al., 2020; Guseinov et al., 2020; Raffone and Cicero, 2018; Min et al., 2021; Zeumault et al., 2022; Fyrigos et al., 2021) or multiscale BD modeling (Strand et al., 2022; Pavodomi, 2017; Padovani et al., 2023; Padovani et al., 2024). Nevertheless, these considerations are consistent with the competition-based, area-dependent BD discussed so far.

Moreover, it is important that the DFT calculations discussed above in Figure 17 represent a system of HfO2 remaining in a pristine state. It remains to be observed whether such a DFT calculation is truly representative of a dielectric experiencing a BD event in a fragmental state. Therefore, our current understanding of dielectric BD is still in an infant stage, probably due to neglecting the importance of area-dependent dielectric BD.

5 Conclusion

In this work, we highlight that defect annihilation is a natural process compatible with the thermodynamics of VO–Oi pair formation because generated defects and annihilation centers coexist. Despite defect generation being the triggering mechanism for dielectric BD, defect annihilation can also participate in filament formation or the BD process in competition and should not be ignored in pursuit of understanding dielectric BD. We show that a time-decaying function for the annihilation process compatible with the maximum-value statistics in competition with the minimum-value statistics for defect generation is the ultimate cause leading to BD. The defect annihilation process can take place in the forms of intrinsic defect-pair recombination and extrinsic defect-passivation by doped species such as hydrogen. Furthermore, an imbalance in generated defects and annihilation centers due to pre-existing species can further complicate the entire breakdown process.

Using PDF analysis, we show that the traditional inverse area-scaling law based solely on minimum-value statistics is merely a special case involving an infinitely long annihilation time. The peaks of the annihilation PDFs were found to generally coincide with the onset of an increase in experimental CDFs. For various area-scaling trends reported (Wu et al., 2023a; Wu et al., 2023b; Wu et al., 2025; Wu et al., 2024; Wu, 2026; Wu et al., 2026), we found a strong correlation between the variability (κ) of the TBD-vs-area relation and the net defect-generation rate (χ), which is strongly altered or modulated by the presence of defect annihilation. This implies that in the traditional inverse area-scaling case, defect annihilation not only exists but also plays an equally critical role in BD filament formation. Area-dependent defect-generation–annihilation competition is an essential and dynamic process ultimately leading to dielectric breakdown, and its role should not be neglected even in the traditional area-scaling case. Our work can stimulate needed debate and interest in bridging the large gaps between experimental reality and the theoretical modeling framework, particularly for large-scale and high-throughput modeling approaches.

In summary, there are several setbacks of ongoing theoretical BD modeling efforts. We think that a multiscale modeling approach can be viable. Our PDF analysis highlights that the primary gap is the lack of a microscopic or atomistic theoretical framework addressing area-dependent dynamics between defect generation and annihilation. Thus, the ongoing BD models address the issue of defect-assisted conduction, rather than true BD behavior.

Generative AI statement

The author(s) declared that generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fnano.2026.1885901/full#supplementary-material

Supplementary material Supplementary material at the source

Footnotes

Edited by: Lijun Liu, Interuniversity Microelectronics Centre (IMEC), Belgium

Reviewed by: Marco A. Villena, University of Granada, Spain

Andrea Vici, Interuniversity Microelectronics Centre (IMEC), Belgium

Statements

Conflict of interest

Author TA was employed by IBM Research Division.

The remaining author(s) declared that this work was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Funding

The author(s) declared that financial support was not received for this work and/or its publication.

Data availability

The data analyzed in this study are subject to the following licenses/restrictions: Data are available from authors upon request. Requests to access these datasets should be directed to ewu@uvm.edu.

Author contributions

EW: Conceptualization, Data curation, Formal Analysis, Investigation, Methodology, Project administration, Software, Supervision, Validation, Writing – original draft, Writing – review and editing. TA: Conceptualization, Investigation, Project administration, Writing – review and editing.

Acknowledgements

The authors are grateful to Jing Li and Magyari-Köpe for a stimulating discussion regarding the work of Sklénard et al. (2023) and Magyari-Köpe et al. (2016) on oxygen vacancy and hydrogen in amorphous HfO2.

References

  1. Aldana S., García-Fernández P., Romero-Zaliz R., González M. B., Jiménez-Molinos F., Gómez-Campos F. (2020). Resistive switching in HfO2-based valence change memories, a comprehensive 3D kinetic monte carlo approach. J. Phys. D. Appl. Phys. 53, 2250106. doi:10.1088/1361-6463/ab7bb6
  2. Ando T., Fagbohungbe O., Imakita K., Miyazoe H., Aoki S., Pujari R. (2025). “ Transfer learning on edge using 14nm CMOS-compatible ReRAM array and analog In-memory training algorithm,” in 2025 IEEE International Electron Devices Meeting (IEDM). San Francisco, CA, USA, 1 – 4. doi:10.1109/IEDM50572.2025.11353900
  3. Bradley S., Shluger A., Bersuker G. (2015). Electron-injection-assisted generation of oxygen vacancies in monoclinic HfO2. Phys. Rev. Appl. 4, 064008. doi:10.1103/PhysRevApplied.4.064008
  4. Chen A. (2014). “ Area and thickness scaling of forming voltage of resistive switching memories,” in IEEE Electron Device Letters 35 (1), 57 – 59. doi:10.1109/LED.2013.2288262
  5. Chen T., Kuo C. (2011). First prinsicples study of the structural, electronic, and dielectric properties of amorphous HfO 2. J. Appl. Phys. 110, 064105. doi:10.1063/1.3636362
  6. Clark R. D., Wajda C. S., Leusink G. J., Edge l. F., Faltermeier J., Jamison P. (2007). Process and electrical characteristics of MO-ALD HfO2 films for high-k gate applications grown in a production worthy 300 mm deposition system. ECS Trans. 11 (3), 55 – 69. doi:10.1149/1.2778650
  7. Clima S., Govoreanu B., Jurczak M., Pourtois G. (2014). HfO X as RRAM material-first principal insights on the working principles. Microelectron. Eng. 120, 13 – 18. doi:10.1016/j.mee.2013.08.002
  8. Consiglio S., Higuchi H., Ando T., Jamison P., Seo S. C., Kong D. (2021). Process-induced ReRAM performance improvement of atomic layer deposited HfO 2 for analog In-Memory computing applications. ECS Trans. 102 (2), 1928 – 2028. doi:10.1149/10202.0019ecst
  9. Cooper D., Baeumer C., Bernier N., Marchewka A., La Torre C., Dunin-Borkowski R. E. (2017). Anomalous resistance hysteresis in oxide ReRAM: oxygen evolution and reincorporation revealed by in-situ TEM. Adv. Mater. 29, 1700212. doi:10.1002/adma.201700212
  10. Degraeve R., Groeseneken G., Bellens R., Ogier J., Depas M., Roussel P. (1998). New insights in the relation between electron trap generation and the statistical properties of oxide breakdown. IEEE Trans. Electron Devices 45 (4), 904 – 911. doi:10.1109/16.662800
  11. Dirkmann S., Kaiser J., Wenger C., Mussenbrock T. (2018). Filament growth and resistive switching in hafnium oxide memristive devices. ACS Appl. Mater. Interfaces 10, 14857 – 14868. doi:10.1021/acsami.7b19836
  12. Dittmann R., Menzel S., Waser R. (2022). Nanoionic memristive phenomena in metal oxides: the valence change mechanism. Adv. Phys. 70, 155 – 349. doi:10.1080/00018732.2022.2084006
  13. Dumin D. J., Maddux R. J., Scott R. S., Subramaniam R. (1994). A model relating wearout to breakdown in thin oxides. IEEE Trans. Electron Devices 41 (9), 1570 – 1580. doi:10.1109/16.310108
  14. Duncan D., Magyari-Köpe B., Nishi Y. (2016). Hydrogen doping in HfO 2 resistance change random access memory. Appl. Phys. Lett. 108, 043501. doi:10.1063/1.4940369
  15. Freire T. C., Mahapatra S., Shluger A. L. (2024). “ Atomistic simulation of defect processes involved in the reaction-diffusion-drift model of NBTI,” in 2024 IEEE International Integrated Reliability Workshop (IIRW), 1 – 4. doi:10.1109/IIRW62856.2024.10947139
  16. Freysoldt C., Grabowski B., Hickel T., Neugebauer J., Kresse G., Janotti A. (2014). First-principles calculations for point defects in solids. Rev. Mod. Phys. 86, 253 – 305. doi:10.1103/RevModPhys.86.253
  17. Fyrigos I.-A., Ntinas V., Sirakoulis G. C., Dimitrakis P., Karafyllidis I. G. (2021). Quantum mechanical model for filament formation in metal-insulator-metal memristors. IEEE Trans. Nanotechnol. 20, 113 – 122. doi:10.1109/TNANO.2021.3049632
  18. Gao D., El-Sayed A., Shluger A. (2016). A mechanism for frenkel defect creation in amorphous SiO2 facilitated by electron injection. Nanotechnology 27 (50), 505207. doi:10.1088/0957-4484/27/50/505207
  19. Gumbel E. J. (1958). Statistics of Extremes. New York, NY, USA: Columbia Univ. Press.
  20. Gunkel F., Christensen D., Chen Y., Pryds N. (2020). Oxygen vacancies: the (in)visible friend of oxide electronics. Appl. Phys. Lett. 116, 120505. doi:10.1063/1.5143309
  21. Guo Y., Robertson J. (2014). Materials selection for oxide-based resistive random-access memories. Appl. Phys. Lett. 105, 223516. doi:10.1063/1.4903470
  22. Guseinov D., Korolev D. S., Belov A. I., Okulich E. V., Okulich V. I., Tetelbaum D. I. (2020). Flexible Monte-Carlo approach to simulate electroforming and resistive switching in filamentary metal-oxide memristive devices, modelling simul. Mater. Sci. Eng. 28, 015007. doi:10.1088/1361-651X/ab580e
  23. Guy J., Molas G., Blaise P., Bernard M., Roule A., Le Carval G. (2015). Investigation of forming, SET, and data retention of conductive-bridge random-access memory for stack optimization. IEEE Trans. Electron Devices 62 (11), 3482 – 3489. doi:10.1109/TED.2015.2476825
  24. Hou Y., Celano U., Goux L., Liu L., Degraeve R., Cheng Y. (2016). Multimode resistive switching in nanoscale hafnium oxide stack as studied by atomic force microscopy. Appl. Phys. Lett. 109, 023508. doi:10.1063/1.4954258
  25. Huang P., Liu X. Y., Chen B., Li H. T., Wang Y. J., Deng Y. X. (2013). A physics-based compact model of metal-oxide-based RRAM DC and AC operations. IEEE Trans. Electron Devices 60 (12), 4090 – 4097. doi:10.1109/IMW.2011.5873222
  26. Joshi K., Lee Y. H., Yao Y. C., Chang S. W., Bian S. S., Liao P. J. (2019). A statistical learning model for accurate prediction of time-dependent dielectric degradation for low failure rates. Proc. IEEE Int. Rel. Phys. Symp. (IRPS), Monterey, CA, U. S. A., 1 – 6. doi:10.1109/IRPS.2019.8720424
  27. Jung K., Magyari-Köpe B., Nishi Y. (2017). Hydrogen-induced oxygen vacancy bistability and its impact on RRAM device operation. IEEE Electron Dev. Lett. 38, 728 – 731. doi:10.1109/LED.2017.2693368
  28. Kaczer B., Degraeve R., O'Connor R., Roussel P., Groeseneken G. (2004). Implications of Progressive wear-out for Lifetime Extrapolation of Ultra-thin (EOT/spl sim/1 Nm) Sion Films, IEDM Technical Digest. San Francisco, CA, USA: IEEE International Electron Devices Meeting, 2004, 713 – 716. doi:10.1109/IEDM.2004.1419270
  29. Kauerauf T., Degraeve R., Ragnarsson L. A., Roussel P., Sahhaf S., Groeseneken G. (2011). “ Methodologies for sub-1nm EOT TDDB evaluation,” in 2011 International Reliability Physics Symposium. Monterey, CA, USA, 2A.2.1 – 2A.2.10. doi:10.1109/IRPS.2011.5784444
  30. Kaviani M., Afanas’ev V., Shluger A. (2017). Interaction of hydrogen with amorphous hafnium oxide. Phys. Rev. B 95, 075117. doi:10.1103/PhysRevB.95.075117
  31. Kerber A., Rohner M., Pompl T., Duschl R., Kerber M. (2007). “ Lifetime prediction for CMOS devices with ultra-thin gate oxides based on progressive breakdown,” in Proc. 2007 IEEE International Reliability Physics Symposium, 217. doi:10.1109/RELPHY.2007.369894
  32. Kim Y., Seo S. C., Consiglio S., Jamison P., Higuchi H., Rasch M. (2021a). Resistive memory process optimization for high resistance switching toward scalable analog compute technology for deep learning. IEEE Electron Device Lett. 42 (5), 759 – 762. doi:10.1109/LED.2021.3066181
  33. Kim Y., Seo S.-C., Jamison P., Consiglio S., Higuchi H., Miyazoe H. (2021b). “ Hardware algorithm Co-optimization on resistive RAM cross-bar array for scalable analog compute technology for AI application, (invited),”. Kyoto, Japan, 1 – 2. 2021 Silicon Nanoelectron. Workshop (SNW).
  34. Kumar S., Wang Z., Huang X., Kumari N., Davila N., Strachan J. (2016). Conduction channel formation and dissolution due to oxygen thermophoresis/diffusion in hafnium oxide memristors. ACS Nano 10, 11205 – 11210. doi:10.1021/acsnano.6b06275
  35. Larcher L., Padovani A. (2017). Multiscale modeling of oxide RRAM devices for memory applications: from material properties to device performance. J. Comput. Electron 16, 1077 – 1084. doi:10.1007/s10825-017-1095-3
  36. Li X., Tung C., Pey K. (2008). The nature of dielectric breakdown. Appl. Phys. Lett. 93, 072903. doi:10.1063/1.2974792
  37. Li D., Li M., Zahid F., Wang J., Guo H. (2012). Oxygen vacancy filament formation in TiO 2: a kinetic monte carlo study. J. Appl. Phys. 112, 073512. doi:10.1063/1.4757584
  38. Linder B., Lombardo S., Stathis J., Vayshenker A., Frank D. (2002). Voltage dependence of hard breakdown growth and the reliability implication in thin dielectrics. IEEE Elec. Dev. Lett. 23, 661 – 663. doi:10.1109/LED.2002.805010
  39. Lombardo S., Stathis J., Linder B., Pey K. L., Palumbo F., Tung C. H. (2005). Dielectric breakdown mechanisms in gate oxides. Appl. Phys. Rev. 98, 121301. doi:10.1063/1.2147714
  40. Lv S., Wang H., Zhang J., Liu J., Sun L., Yu Z. (2014). An analytical model for the forming process of conductive-bridge resistive-switching random-access memory. IEEE Trans. Electron Devices 61 (9), 3166 – 3171. doi:10.1109/TED.2014.2341274
  41. Magyari-Köpe B., Duncan D., Zhao L., Nishi Y. (2016). “ Doping technology for RRAM — opportunities and challenges,” in International Symposium on VLSI Technology, Systems and Application (VLSI-TSA), Hsinchu, Taiwan, 1. doi:10.1109/VLSI-TSA.2016.7480494
  42. Maier J. (2004). Physical Chemistry of Ionic Materials: Ions and Electrons in Solids. Chichester: Wiley.
  43. Mehonic A., Shluger A. L., Goa D., Valov I., Minranda E., Ielmini D. (2018). Silicon oxide (SiO X): a promising material for resistance switching? Adv. Mater. 28, 1801187. doi:10.1002/adma.201801187
  44. Menzel S., Waser R. (2019). “ Chap. 4, mechanism of memristive switching in OxRAM,” in Advances in Non-volatile Memory. Woodhead Publishing Series in Electronic and Optical Materials (Cambridge, England: Elsevier), 137 – 170. doi:10.1016/B978-0-08-102584-0.00005-X
  45. Min K., Jung D., Kwon Y. (2021). Investigation of switching uniformity in resistive memory via finite element simulation of conductive-filament formation. Sci. Rep. 11, 2447. doi:10.1038/s41598-021-81896-z
  46. Miranda E., Sũné J., Rodríguez R., Nafría M., Martín F., Aymerich X. (1998). Soft breakdown in ultrathin SiO 2 layer: the conduction problem from a new point of view. Jpn. J. Appl. Phys. 38, 2223. doi:10.1143/JJAP.38.2223
  47. Monsieur F., Vincent E., Roy D., Bruyère S., Pananakakis G., Ghibaudo G. (2002). A thorough investigation of progressive breakdown in ultra-thin oxides, physical understanding and application for industrial reliability assessment. Proc. IRPS, 45 – 54. doi:10.1109/RELPHY.2002.996609
  48. Mueller M., De Souza R. A. (2018). SIMS study of oxygen diffusion in monoclinic HfO 2. Appl. Phys. Lett. 112, 051908. doi:10.1063/1.5020370
  49. Padovani A., Torraca P. L., Strand J., Shluger A., Milo V., Larcher L. (2023). “ Towards a universal model of dielectric breakdown,” in 2023 IEEE International Reliability Physics Symposium (IRPS). Monterey, CA, USA, 1 – 8. doi:10.1109/IRPS48203.2023.10117846
  50. Padovani A., La Torraca P., Strand J., Larcher L., Shluger A. L. (2024). Dielectric breakdown of oxide films in electronic devices. Nat. Rev. Mater 9, 607. doi:10.1038/s41578-024-00702-0
  51. Pan F., Yin S., Subramanian V. (2011). “ A comprehensive simulation study on metal conducting filament formation in resistive switching memories,” in Proc. 2011 International Memory Workshop. doi:10.1109/IMW.2011.5873221
  52. Pavodomi A. (2017). A microscopic mechanism of dielectric breakdown in SiO 2 films: an insight from multi-scale modeling. J. Appl. Phys. 121, 155101 2017. doi:10.1063/1.4979915
  53. Peirce, Tensile Tests for Cotton Yarns (1926). “The Weakest Link” theorems on the strength of long and of composite specimens. J. Text. Inst. 17, 355. doi:10.1080/19447027.1926.10599953
  54. Pompl T., Kerber A., Röhner M., Kerber M. (2006). Gate voltage and oxide thickness dependence of progressive wear-out of ultra-thin gate oxides. Microelec. Relia. 46, 1603 – 1607. doi:10.1016/j.microrel.2006.07.042
  55. Qin S., Liu Z., Zhang G., Zhang J., Sun Y., Wu H. (2015). Atomistic study of dynamics for metallic filament growth in conductive-bridge random access memory, Phys. Chem. Chem. Phys. 17, 8627--8632 8632. doi:10.1039/c4cp04903a
  56. Raffone F., Cicero G. (2018). Unveiling the fundamental role of temperature in RRAM switching mechanism by multiscale simulations. ACS Appl. Mater. and Interfaces 10, 7512. doi:10.1021/acsami.8b00443
  57. Raghavan N., Fantini A., Degraeve R., Roussel P., Goux L., Govoreanu B. (2013). Statistical insight into controlled forming and forming free stacks for HfO x RRAM. Microelec. Eng. 109, 177 – 181. doi:10.1016/j.mee.2013.03.065
  58. Rodriguez-Fernandez A., Muñoz-Gorriz J., Suñé J., Miranda E. (2018). A new method for estimating the conductive filament temperature in OxRAM devices based on escape rate theory. Micro. Relia. 88–90, 142 – 146. doi:10.1016/j.microrel.2018.06.120
  59. Roussel Ph.J., Degraeve R., Sahhaf S., Groeseneken G. (2007). A consistent model for the hard breakdown distribution including digital soft breakdown: the noble art of area scaling. Microelectron. Eng. 84, 1925 – 1928. doi:10.1016/j.mee.2007.04.067
  60. Sadi T., Mehonic A., Montesi L., Buckwell M., Kenyon A., Asenov A. (2018). Investigation of resistance switching in SiO x RRAM cells using a 3D multi-scale kinetic monte carlo simulator. J. Phys. Condens. Matter 30, 084005. doi:10.1088/1361-648X/aaa7c1
  61. Sahcchi S., Varanasi A., Degraeve R., Vici A., Molinaro G., Franco J. (2026). Markov model describing progressive degradation of local percolation path in thin oxides. Solid-State Electron. 231, 109266. doi:10.1016/j.sse.2025.109266
  62. Schie M., Müller M. P., Salinga M., Waser R., De Souza R. A. (2017). Ion migration in crystalline and amorphous HfO x. J. Chem. Phys. 146, 094508. doi:10.1063/1.4977453
  63. Schie M., Menzel S., Robertson J., Waser R., De Souza R. (2018). Field-enhanced route to generating anti-Frenkel pairs in HfO2. Phys. Rev. Mater. 2, 035002. doi:10.1103/PhysRevMaterials.2.035002
  64. Sharma A. A., Noman M., Abdelmoula M., Skowroński M., Bain J. A. (2014). Electronic instabilities leading to electroformation of binary metal oxide‐based resistive switches. Adv. Funct. Mater. 24 (35), 5522 – 5529. doi:10.1002/adfm.201400461
  65. Shluger A. L., McKenna K. P. (2013). “ Models of oxygen vacancy defects involved in degradation of gate dielectrics,” in 2013 IEEE International Reliability Physics Symposium (IRPS). Monterey, CA, USA, 5A.1.1 – 5A.1.9. doi:10.1109/IRPS.2013.6532018
  66. Sklénard B., Cvitkovich L., Waldhoer D., Li J. (2023). Oxygen vacancy and hydrogen in amorphous HfO 2. J. Phys. D. Appl. Phys. 56, 245301. doi:10.1088/1361-6463/acc878
  67. Sowinska M., Bertaud T., Walczyk D., Thiess S., Schubert M. A., Lukosius M. (2012). Hard x-ray photoelectron spectroscopy study of the electro-forming in Ti/HfO2-based resistive switching structures. Appl. Phys. Lett. 100, 233509. doi:10.1063/1.4728118
  68. Stathis J. H., DiMaria D. J. (1998). “ Reliability projection for ultra-thin oxides at low voltage,” in International Electron Devices Meeting 1998. San Francisco, CA, USA: Technical Digest Cat. No.98CH36217, 167 – 170. doi:10.1109/IEDM.1998.746309
  69. Strand J., Kaviani M., Shluger A. L. (2017). Defect creation in amorphous HfO 2 facilitated by hole and electron injection. Microelectron. Eng. 178, 279 – 283. doi:10.1016/j.mee.2017.05.005
  70. Strand J., La Torraca P., Padovani A., Larcher L., Shluger A. (2022). Dielectric breakdown in HfO 2 dielectrics: using multiscale modeling to identify the critical physical processes involved in oxide degradation. J. Appl. Phys. 131, 234501. doi:10.1063/5.0083189
  71. Suñé J. (2001). New physics-based analytic approach to the thin oxide breakdown statistics. IEEE Electron Dev. Lett. 22, 296 – 298. doi:10.1109/55.924847
  72. Suñé J., Wu E. (2003). Statistics of successive breakdown events in gate oxides. IEEE Elec. Dev. Lett. 24, 272 – 274. doi:10.1109/LED.2003.812139
  73. Suñé J., Placencia I., Barniol N., Farrés E., Martín F., Aymerich X. (1990). On the breakdown statistics of very thin SiO films. Thin Solid Films 185, 347 – 362. doi:10.1016/0040-6090(90)90098-X
  74. Suñé J., Wu E., Lai W. (2006). Statistics of competing post-breakdown failure modes in ultrathin MOS devices. IEEE Trans. Electron Devices 53 (2), 224 – 234. doi:10.1109/TED.2005.861597
  75. Suñé J., Wu E., Tous S. (2007). A physics-based deconstruction of the percolation model of oxide breakdown. Microelec. Eng. 84, 1917 – 2007. doi:10.1016/j.mee.2007.04.008
  76. Szot K., Speier W., Bihlmayer G., Waser R. (2006). Switching the electrical resistance of individual dislocations in single-crystalline SrTiO 3. Nat. Mater. 5, 312. doi:10.1038/nmat1614
  77. Tous S., Wu E. Y., Sune J. (2008). A compact model for oxide breakdown failure distribution in ultrathin oxides showing progressive breakdown. IEEE Electron Device Lett. 29, 949 – 951. doi:10.1109/LED.2008.2001178
  78. Traore B., Blaise P., Vianello E., Perniola L., De Salvo B., Nishi Y. (2016). HfO 2 -based RRAM: electrode effects, Ti/HfO 2 interface, charge injection, and oxygen (O) defect diffusion through experiment and Ab Initio calculations. IEEE Trans. Electron Device 63, 360 – 368. doi:10.1109/TED.2015.2503145
  79. Tsang J., Linder B. P. (2004). Characterization of breakdown in ultrathin oxides by hot electron emission. Appl. Phys. Lett. 84, 4641 – 4643. doi:10.1063/1.1759383
  80. Vandelli L., Padovani A., Larcher L., Bersuker G. (2013). Microscopic modeling of electrical stress-induced breakdown in poly-crystalline hafnium oxide dielectrics. IEEE Trans. Electron Devices 60 (5), 1754 – 1762. doi:10.1109/TED.2013.2255104
  81. Vici A., Degraeve R., Franco J., Kaczer B., Roussel P. J., De Wolf I. (2023). Analytical markov model to calculate TDDB at any voltage and temperature stress condition. IEEE Trans. Electron Devices 70 (12), 6512 – 6519. doi:10.1109/TED.2023.3326430
  82. Vici A., Pantisano L., Roussel P. J., Chasin A., Franco J., Kaczer B. (2025). Unified analysis of time-dependent dielectric breakdown to extend reliability margins. IEEE Trans. Electron Devices 72 (11), 5780 – 5792. doi:10.1109/TED.2025.3613295
  83. Vollertsen R.-P., Kleppmann W. G. (1990). “ Dependence of dielectric time to breakdown distributions on test structure area,” in Proceedings of the 1991 International Conference on Microelectronic Test Structures, Kyoto, Japan, 75 – 79. doi:10.1109/ICMTS.1990.161716
  84. Wei T., Lu Y., Zhang F., Tang J., Gao B., Yu P. (2023). Three-dimensional reconstruction of conductive filaments in HfO 2 -based memristor. Adv. Mater. 35, 2209925. doi:10.1002/adma.202209925
  85. Wolters D. (1958). Breakdown and Wearout Phenomena in SiO 2, Insulating Films on Semiconductors, 180, 1981. Gumbel, Statistics of Extremes. New York, NY, USA: Columbia Univ. Press.
  86. Wu E. Y. (2019a). Facts and myths of dielectric breakdown processes – part I: statistics, experimental, and physical acceleration models. IEEE Tran. Elec. Dev. 66, 4523. doi:10.1109/TED.2019.2933612
  87. Wu E. (2019b). Facts and myths of dielectric breakdown processes – part II: post-breakdown and variability. IEEE Tran. Elec. Dev. 66, 4535 – 4546. doi:10.1109/TED.2019.2933615
  88. Wu E., Sune J. (2006). “ Statistical and voltage scaling properties of post-breakdown for ultra-thin-oxide PFETs in inversion mode ” in 2006 Proc. IEEE International Reliability Phys. Symp, 54 – 62. doi:10.1109/RELPHY.2006.251191
  89. Wu E., Sũné J. (2005). Power-law voltage acceleration: a key element for ultra-thin gate oxide reliability. Microelectron. Rel. 45, 1809 – 1834. doi:10.1016/j.microrel.2005.04.004
  90. Wu E. Y., Vollertsen R.-P. (2002). On the weibull shape factor of intrinsic breakdown of dielectric films and its accurate experimental determination. Part I: theory, methodology, experimental techniques. IEEE Trans. Electron Devices 49 (12), 2131 – 2140. doi:10.1109/TED.2002.805612
  91. Wu E., Jamison P., Consiglio S., Tsunomura T., Ando T. (2025). “ Modulating competition between defect generation and annihilation in dielectric breakdown showing a full range of reverse, diminishing, and forward area scaling trends,” in 2025 IEEE International Reliability Physics Symposium, 1 – 8. doi:10.1109/IRPS48204.2025.10983848
  92. Wu E., Tous S., Suñé J. (2007). “ On the progressive breakdown statistical distribution and its voltage acceleration,” in Digest of the 2007 International Electron Devices Meeting, 493 – 496. doi:10.1109/IEDM.2007.4418982
  93. Wu E., Suñé J., Linder B., Achanta R., Li B., Mittl S. (2011). “ Post-breakdown statistics and acceleration characteristics in high-K dielectric stacks,” in 2011 International Reliability Physics Symposium. Monterey, CA, USA, 2A.1.1 – 2A.1.6. doi:10.1109/IRPS.2011.5784443
  94. Wu E., Suñé J., LaRow C. (2012). Generalized successive failure methodology for non-weibull distributions and its applications to SiO 2 or high-k/SiO 2 bilayer dielectrics and extrinsic failure mode. Proc. IRPS 6A (1), 20012 – 6A.2.7. doi:10.1109/IRPS.2012.6241863
  95. Wu E., Ando T., Li B., Southwick R., Stathis J. (2019). Fundamental roles of extreme-value distributions in dielectric breakdown and memory applications (minimum-value versus maximum-value statistics). Jpn. J. Appl. Phys. 59, SM083. doi:10.35848/1347-4065/ab918b
  96. Wu E., Jamison P., Ando T. (2023a). Reversing A decades-long scaling law of dielectric breakdown using hydrogen-plasma-treated HfO 2 ReRAM devices. Adv. Electron Mater. 9, 2300296. doi:10.1002/aelm.202300296
  97. Wu E., Ando T., Jamison P. (2023b). “ Reversing a decades-long scaling law of dielectric breakdown for ReRAM forming voltage reduction - modeling competition among defect generation and annihilation processes (invited),” in 2023 International Electron Devices Meeting (IEDM), 1 – 4. doi:10.1109/IEDM45741.2023.10413665
  98. Wu E., Ando T., Jamison P. (2024). “ General statistical model for dielectric breakdown including reverse area scaling - the role of area-dependent dynamic competition,” in 2024 IEEE International Reliability Physics Symposium (IRPS), 1 – 8. doi:10.1109/IRPS48228.2024.10529419
  99. Wu E. Y. (2026). Area-dependent dielectric breakdown: an overlooked yet inevitable physical phenomenon (invited). HsinChu, Taiwan 1, 4. doi:10.1109/VLSITSA69131.2026.11527700
  100. Wu E. Y., Menzel S., Ando T. (2026). “ Is Dielectric (Breakdown) Research at a Crossroad? Its Impact on Technology Innovation,” in IEEE Transactions on Electron Devices. doi:10.1109/VLSITSA69131.2026.11527700
  101. Xu X., Rajendran B., Anantram M. P. (2020). “ Kinetic Monte Carlo simulation of interface-controlled hafnia-based resistive memory,” in IEEE Trans. on Electron Devices, 67, 118. doi:10.1109/TED.2019.2953917
  102. Yang J. J., Miao F., Pickett M., Ohlberg D., Stewart D., Lau C. (2009). The mechanism of electroforming of metal oxide memristive switches. Nanotechnology 20, 215201. doi:10.1088/0957-4484/20/21/215201
  103. Zeumault A., Alam S., Omar Faruk M., Aziz A. (2022). Memristor compact model with oxygen vacancy concentrations as state variables. J. Appl. Phys. 131, 124502. doi:10.1063/5.0087038
  104. Zhao X., Vanderbilt D. (2002). First-pinciples study of structural, vibrational, and lattice dielectric properties of Hafnium oxide. Phys. Rev. B 65, 233106. doi:10.1103/PhysRevB.65.233106

Citation

Wu, E. Y. & Ando, T. (2026). Statistics and physics of area-dependent dielectric breakdown: an overlooked yet inevitable natural phenomenon. Frontiers in Nanotechnology, Volume 8. https://doi.org/
Received
Accepted
Published
Volume
8 - 2026

Copyright

Copyright © 2026 Wu and Ando.

License: https://creativecommons.org/licenses/by/4.0/

*Correspondence: Ernest Y. Wu

Article metrics

Views
0
Downloads
0
Citations
0
Shares
0