Correlated insulator collapse due to quantum avalanche via in-gap ladder states

The significant discrepancy observed between the predicted and experimental switching fields in correlated insulators under a DC electric field far-from-equilibrium necessitates a reevaluation of current microscopic understanding. Here we show that an electron avalanche can occur in the bulk limit of such insulators at arbitrarily small electric field by introducing a generic model of electrons coupled to an inelastic medium of phonons. The quantum avalanche arises by the generation of a ladder of in-gap states, created by a multi-phonon emission process. Hot-phonons in the avalanche trigger a premature and partial collapse of the correlated gap. The phonon spectrum dictates the existence of two-stage versus single-stage switching events which we associate with charge-density-wave and Mott resistive phase transitions, respectively. The behavior of electron and phonon temperatures, as well as the temperature dependence of the threshold fields, demonstrates how a crossover between the thermal and quantum switching scenarios emerges within a unified framework of the quantum avalanche.

The significant discrepancy observed between the predicted and experimental switching fields in correlated insulators under a DC electric field far-from-equilibrium necessitates a reevaluation of current microscopic understanding. Here we show that an electron avalanche can occur in the bulk limit of such insulators at arbitrarily small electric field by introducing a generic model of electrons coupled to an inelastic medium of phonons. The quantum avalanche arises by the generation of a ladder of in-gap states, created by a multi-phonon emission process. Hot-phonons in the avalanche trigger a premature and partial collapse of the correlated gap. The phonon spectrum dictates the existence of two-stage versus single-stage switching events which we associate with chargedensity-wave and Mott resistive phase transitions, respectively. The behavior of electron and phonon temperatures, as well as the temperature dependence of the threshold fields, demonstrates how a crossover between the thermal and quantum switching scenarios emerges within a unified framework of the quantum avalanche.
In spite of recent progress in our understanding of the nonequilibrium many-body state of matter, one of the long-standing problems that has remained unresolved concerns the microscopic mechanism behind the insulator-to-metal transition (IMT) of stronglycorrelated electronic systems, driven by a DC electric field. For more than five decades, the community has fiercely debated the origin of the dielectric breakdown in charge-density-wave (CDW) systems [1][2][3][4][5][6], for which the reported threshold electric fields are orders of magnitude smaller than theoretical estimates based on the Landau-Zener [7] mechanism. In the late seventies, this problem led to the formulation of the classical theory of depinning, in which the CDW order parameter is understood as being pinned by the presence of disorder and can be abruptly unpinned under the action of a static electric field [8][9][10]. More recently, a similar mismatch between theory and experiments has also been found in studies of transition-metal compounds such as Mott insulators [11][12][13][14][15][16]. In spite of the emerging potential of these materials for applications such as non-volatile neuromorphic computing [17], the lack of understanding of the microscopic origins of their resistive transitions is a bottleneck to the development of such technologies.
From the early days of CDW research, the theoretical paradigm for the resistive transition in correlated systems has been the Landau-Zener mechanism [3,7]. This model predicts switching fields on the order of with the (zero-field) insulating gap ∆ and the Fermivelocity v F . Using typical electronic energy scales in the above expression, we obtain a rough estimate of E LZ ∼ 10 −2 V/Å = 10 6 V/cm, many orders of magnitude larger than the switching fields of 10 kV/cm found in Mott insulators [11] and 1 − 10 3 V/cm in CDW materials [1][2][3][4]. Over the years, various theoretical attempts have been made to improve the description of the resistive switching transition in a many-body context. Existing theories include the explicit time-evolution of the Hubbard chain [18], the multi-band Hubbard model [19], disorder-driven [20] and nonequilibrium phase transitions [21], spatial inhomogeneity [22], and Coulomb blockade of multi-solitons [5,23]. While bandedge broadening [21,24,25] and metal-insulator filament dynamics [22] by a uniform electric field are known to reduce the switching fields, these effects are quite modest. Recently, a nucleation mechanism for metallic domains, assisted by electron-phonon coupling, has been reported [26]. These various microscopic models have been unable to settle the aforementioned energy-scale problem, or to provide clarity to the age-old debate over the role of thermal [27] versus electronic [28,29] mechanisms in the resistive transition. It has therefore long been speculated that an important ingredient must be missing.
With a wide class of switching materials [11] the mechanisms could be just as diverse. The impact ionization in semiconductors and the electro-migration in nanoscale oxide devices have been well studied. In correlated oxides much of the phenomenology has been understood in terms of the filamentary dynamics [11,12,14,27] The recent report on the nano-scale resolution of nonfilamentary patterns in Ca 2 RuO 4 [16], however, suggests quite a different non-thermal mechanism, reminiscent of the switching in organic solids [30] in which perpendicular patterns to bias have been observed. The importance of the pattern formation along the bias direction in device geometry has been pointed out by recent theories [31][32][33]. Given the diverse switching phenomena, we limit the scope of this paper to a new switching mechanism in the bulk correlated insulators where the switching is understood as phase transition controlled by the electric field. The colossal mismatch of the switching fields in the field-driven Mott and CDW transitions motivates us to look for a common underlying mechanism that is shared between them. Despite the differences in CDW and Mott phenomenology, it is often hard to disentangle one mechanism from the other [34]. Here we propose that the resistive transitions in CDW and Mott insulators can be coherently explained in terms of the role of inelastic many-body processes. We show how these phononemitting transitions can lead to the generation of in-gap states, resulting in a quantum avalanche and a destruction of the correlated gap, at much smaller electric field scales than previously predicted.
For a conceptual understanding of the avalanche mechanism, we start by considering a one-dimensional conduction band of mass m separated by a gap ∆ from the Fermi energy, see Fig. 1 a. The electrons are subject to three processes. First, they are accelerated by a DC electric field E. Second, they can emit optical phonons of energy ω ph by means of an electron-phonon scattering process controlled by the coupling parameter g ep . Finally, dephasing introduces a finite lifetime /Γ to the electrons, setting a cutoff on the time scale over which they can form the multi-phonon state. The dephasing may arise from electron tunneling to a substrate, virtual transitions to higher electronic states, and other manybody mechanisms separate from the electron-phonon interaction. Hereafter, we adopt a unit system in which = 1, the electric charge e = 1, the Boltzmann constant k B = 1, and the lattice constant a = 1.
Let us now give a back-of-the-envelope criterion for an avalanche that can be triggered by a combination of the three aforementioned electronic processes. We note that, in an insulator subject to an external bias, there will exist dilute, yet nonetheless finite, charge excitations. Due to the Franz-Keldysh effect [24], the bandedge smears into the gap and the transition of energetic electrons into those states by phonon emission is substantially enhanced over that in the equilibrium limit. Importantly, electron-phonon scattering processes create a ladder of intermediate replica bands that are equally spaced by ω ph . The minimal number of successive phononscattering events to bridge the gap is N ladder ∼ ∆/ω ph , see Fig. 1b. The timescale it takes to accelerate an electron in a dispersive band, to gain the energy of a single phonon, can be estimated as τ ∼ (m ω ph ) 1/2 /E with the electron band mass m. The electron-phonon scattering rate is roughly proportional to the dimensionless electron-phonon mass-renormalization factor λ ≡ (g ep /ω ph ) 2 ν 0 , where the typical density of states ν 0 ∼ m. The number of phonons that are emitted during the electronic lifetime 1/Γ is then N ep ∼ λ/Γτ . Consequently, the criterion for the avalanche becomes N ep = N ladder , yielding the following electric field required to trigger the avalanche This outcome is highly non-trivial in that the dissipative electron decay rate Γ plays a crucial role in the onset of this nonequilibrium quantum phase transition. The Γ → 0 limit is singular in nonequilibrium. By taking Γ → 0 before or after the E → 0 limit, we arrive at an effectively infinite [35] or a zero temperature (equilibrium) limit, respectively. This singular nature of the Γ = 0 limit plays a fundamental role in the avalanche, as will be confirmed below. The switching field is very different in nature from that of Eq. (1), and we note that the avalanche via inelastic scattering is also distinct from the conventional mechanism that has been discussed for high-field transport in semiconductors [36]. Most importantly, as we shall see below, the avalanche fields predicted by Eq. (2) are much smaller than those expected from Eq. (1).

Results
Quantum avalanche in a rigid-band model To demonstrate the existence of the proposed avalanche mechanism, we investigate a model of a one-band tightbinding chain under a DC electric-field E, in the Coulomb gauge [37,38] where d † i /d i is the creation/annihilation operator for an electron at site i, for which the site position x i = ia. ∆ 0 is the bare gap and t is the tight-binding parameter. The electrons are locally coupled to phonons, which are modeled by a collection of harmonic oscillators given by the Hamiltonian with ϕ k the amplitude, p k the momentum, k the continuum index, and ω k the frequency of the phonon. In this section we consider the case of optical phonons, with ω k = ω ph , deferring a discussion of acoustic phonons until later. The on-site electron-phonon coupling is given by with the coupling constant g ep . We use the Schwinger-Keldysh formulation of the dynamical mean-field theory (DMFT [37,39,40]) which bypasses the transient dynamics and directly yields the homogeneous nonequilibrium steady-state of the manybody dynamics [41,42]. The fermion baths enter the computation of the electronic Green's function via local retarded and lesser self-energies at site i as while the Ohmic baths [43] enter the phonon Green's function via local self-energies [38,44] as In the above expressions, f 0 (ω) = (e ω/T + 1) −1 and n 0 (ω) = (e ω/T − 1) −1 are the Fermi-Dirac and Bose-Einstein distributions at the bath temperature T , respectively, and τ P [45] is the phonon decay time. We compute the second-order self-energy to electrons and phonons on the same footing. We refer the reader to the Methods Section and to the Supplementary Information for further details on the electron-phonon calculations.
In the calculations that follow, we choose the model parameters ∆ 0 = 1, t = 1, g ep = 0.25, ω ph = 0.3, τ −1 P = 0.001, and T = 0.001 in units of eV. We use the lattice constant a = 5Å to compute the electricfield. We verify the schematic estimate of the avalanche b, Post-avalanche spectra of electrons and phonons, plotted using semi-log scales. The electron spectral function A el (ω) with its band edge at the gap ∆ ≈ 0.9 (slightly reduced from ∆0 = 1) extends down to the reservoir level ω = 0, and the distribution function f el (ω) reveals in-gap state occupation with multi-phonon edges (marked by red arrows) separated by the phonon energy ω ph . The distribution function n ph (ω) for hot-phonons is equilibrium-like.
field, Eq. (2), via systematic comparison with full manybody calculations over a wide range of parameters, as detailed in the Supplementary Information. There, we demonstrate the conceptual validity of the avalanche as arising from a competition between inelastic transport and dephasing and justify the limitations of Eq. (2) in the high-field limit. Fig. 2a shows the electronic excitations n ex per site (0 ≤ n ex ≤ 1) when ramping up the electric field starting from the insulator state. Filled (empty) symbols denote results obtained by discarding (including) the nonequilibrium self-energy of the phonons. At the onset of the avalanche at the threshold field E av , n ex increases abruptly, yet continuously, in a similar fashion as in critical phenomena. Notably, E av is found to be roughly proportional to Γ, consistent with the back-ofthe-envelope estimate made in Eq. (2). Furthermore, the pre-avalanche excitations are inversely proportional to E av [see inset to (a)], which demonstrates that the avalanche is not initiated by thermal excitation but is of quantum origin. Remarkably, the onset of the avalanche is not affected by hot-phonon effects since they are only infinitesimally excited at the continuous avalanche transition, while the incoherence of phonons impedes the avalanche resulting in a linear increase of n ex well after the avalanche.
The nontrivial behavior of the switching field on Γ (see Fig. 2a) should be distinguished from phonon-assisted tunneling [46,47], which is commonly manifested in resonant tunneling in heterostructures. With the electron oc-cupation quickly reaching beyond 0.1 after the avalanche, we are in a metallic bulk limit through a phase transition, instead of in the perturbative tunneling regime. The avalanche shown here only emerges after a fully selfconsistent solution is reached after hundreds of iterations, unlike what is expected from a low-order perturbative theory [46]. Figure 2b elucidates nature of the nonequilibrium state after the avalanche. The energy distribution f el (ω) reveals a sizable nonequilibrium occupation of the in-gap states around energies at multiples of ω ph away from the band edge. This confirms the involvement of multiphonon emission in the avalanche mechanism. Moreover, the strong deviation of f el (ω) from a Fermi-Dirac distribution points to the non-thermal nature of the avalanche. In contrast, the phonon distribution n ph (ω) appears mostly thermal.

Insulator-to-metal transition induced by avalanche
Having established the existence of the avalanche, we now turn our attention to the implications of the avalanche in the context of the resistive transition in correlated insulators. The strong charge fluctuations initiated by the avalanche are expected to generate phonon excitations, which then profoundly perturb the interorbital mixing and destabilize the charge gap. To address this point, we extend our model to a two-band correlated insulator where the gap ∆ between the bands is generated by electronic interactions. (See Fig. 3a.) Specifically, we consider the following Hamiltonian which allows us to address both the CDW transitions and resistive switching on an equal footing, Here, α is the band index and i is the site index on a square lattice with ij denoting nearest neighbors. The last term in Eq. (8) is the decoupling of an inter-orbital Hubbard-like interaction via an auxiliary quantum field ξ. The magnitude of ξ respresents the order parameter that opens a gap, i.e. the strength of the density modulation in CDW or the charge fluctuations in correlated insulators. The phase of ξ captures the phase slip in CDW materials or the Goldstone excitations in correlated insulators. In this work, we investigate the instability of a uniform symmetry-broken state, and ignore the phase fluctuations of ξ for simplicity. Consideration of ξ-fluctuations would only weaken the order parameter and lead to an even earlier collapse of the insulating state, and thus further supporting our conclusions. The static mean-field condition becomes The electron-phonon coupling is given by H ep = g ep i,α ϕ αi d † αi d αi , with the independent phonons coupling to each orbital. The fermion/phonon baths are set up as described previously with the Fermi energy at the band crossing.
The spectrum of the phonons has important consequences for the IMT. Let us first discuss the electricfield driven IMT in the presence of optical phonons with energy ω ph , which we will associate with the multi-stage transitions that are commonly observed in CDW systems. Later, we will turn to the case of acoustic phonons which we will relate to Mott systems. With increasing electric field, the insulating system undergoes a two-stage transition to a metal as manifested by the spectral function in Fig. 3b-d. Here, we adjust the electron-phonon coupling g ep such that, given U = 2, the initial gap ∆ 0 is tuned to 1.0 at E = 0 and at the lowest temperature. The strength of this coupling corresponds to a moderate (dimensionless) mass-renormalization factor of λ ≈ (g ep /ω ph ) 2 ν 0 = 0.32 with the phonon energy ω ph = 0.3 eV and the 2D density of states ν 0 ≈ (8t) −1 . (The damping is set to Γ = τ −1 P = 0.001 eV [48], and the chemical potential to µ = t.) In the low-field limit, the spectrum features a well-defined charge gap. As E increases beyond the first threshold field E T , in-gap states develop via the avalanche mechanism. This results in a pseudogap phase where the system becomes metallic while the charge gap mostly remains intact. The energy-resolved electron occupation n ex (ω) = (2π) −1 ImG < (x i = 0, ω) (the red shaded area shown in the inset to panel c) indicates that the electric current is carried mainly by the states of the conduction band above the gap while the in-gap states provide the pathway for this population inversion. The pseudo-gap regime persists until E = E ∆ , when the gap ∆ collapses to zero in a strongly discontinuous transition.
The current-voltage (I − V ) relation in Fig. 3e shows evidence of a two-stage IMT. First, the system continuously becomes metallic at the threshold field E = E T . Later, at the higher threshold E ∆ , the current rises discontinuously. The avalanche behavior discussed in Fig. 2 is responsible for the threshold behavior at E T . This behavior is similar to that in CDW systems, in which it has been widely attributed to the depinning transition [8,9]. Here, we propose an alternative mechanism of electron avalanche via inelastic scattering that does not require any disorder or reduced dimensionality. With coupling to optical phonons, the avalanche is not sufficiently disruptive to the gap, and the intermediate gapped state is instead sustained over a wide range of electric field (E T < E < E ∆ ). We note that E T is around two orders of magnitude smaller than the switching field expected for zero electron-phonon coupling strength E ∆ (λ = 0, ∆ 0 = 1) ≈ 1.6 MV/cm.
The non-thermal nature of the avalanche transition is illustrated by the electric-field dependence of the effective temperatures for electrons and phonons, as shown in Fig. 3f. (See Methods for the definition of effective temperatures.) As the electric field is increased beyond E = E T , the electrons heat up immediately while the phonons stay cold. This clearly demonstrates that heat exchange does not trigger the avalanche. On the other hand, the full IMT at E ∆ is initiated once the electron and phonon temperatures equilibrate, suggesting that this second transition can be described in terms of a thermal mechanism. It is remarkable to demonstrate that the mechanisms for electronic and thermal switching, the topic of intense discussions in the literature, are not exclusive of each other but can arise simultaneously from a single microscopic model. The phase diagram defined by the two switching fields (E T and E ∆ ) is plotted in Fig. 3g. The most notable observation here is that E T remains constant near zero temperature, but increases at higher temperatures until it merges with E ∆ . This may seem counter-intuitive, since the gap is expected to decrease with increasing T . The gap, however, remains nearly constant until close to a critical temperature T c , meaning that a thermal argu-ment is not applicable. We find that, as suggested by the hot-phonon effects apparent in Fig. 2, thermal decoherence may be detrimental to the avalanche, meaning that a stronger electric field is required to induce an avalanche. This is a further evidence that the threshold behavior is of quantum-mechanical origin. In contrast to E T , E ∆ decreases with T , and is conventional and thermally driven.
The existence of bias-driven multiple-stage transitions in CDW systems has been intensely debated [1,[49][50][51][52]. For instance, recent ARPES studies on the NbSe 3 system have shown unambiguously that the various CDW gap energies are constant (marginally increasing) with increasing temperature and observable even beyond their respective T c values [53,54]. This suggests a nonthermal mechanism of gap formation below T c . At temperatures beyond T c , the gap gradually closes suggesting a thermal mechanism, where these predictions show an interesting parallel to our results. What our model as a conceptual framework showed are that the nature of the initial threshold [51] in the low field limit is quantum in the sub-gap energy scale, and that, after an intermediate pseudo-gapped phase, there are subsequent discontinuous resistive transitions that thermally destroy the order parameter.
Finally, let us discuss the case of electron coupling to acoustic phonons, with a continuous spectrum of the form ω k ∝ k and cutoff Debye energy chosen as ω D = 0.6 eV.  [55], y-axis arbitrarily scaled). c, I − V relation in the strong-coupling limit showing non-monotonic dependence of E∆ on T . d, Non-monotonic E∆ versus T . The initial increase of E∆ at low T demonstrates the non-thermal nature of the IMT in the strong-coupling limit. kBTc = 0.025. The blue circles are switching fields of TiS3 with Tc = 290 K (adapted from Ref [56,57], y-axis arbitrarily scaled). The red squres are the same as (b).
Compared to optical phonons, the influence of acoustic phonons on the nonequilibrium dynamics is more dramatic. We discuss this situation for both the weak and strong limits of electron-phonon coupling, as shown in Fig. 4. In the weak-coupling limit (λ = 0.113), corresponding to panels (a) and (b), the IMT is dominated by a strongly discontinuous collapse of the gap. While the signature of the avalanche is still present as a precursor to the IMT, the range of the avalanche region is much narrower than that found for coupling to optical phonons. The avalanche current is also very small so that its effect is insignificant. E ∆ decreases monotonically with increasing temperature and the E−T phase diagram is conventional and fully consistent with the thermal switching scenario. We tested this result against experimental data by overlaying the experimental switching fields in the Vanadium oxides W x V 1−x O 2 [55] at temperatures close to the transition temperature. In addition to an overall agreement of the trends, the lightly concave curve shape [55,58] found in the experiment is reproduced by the theory.
In the strong electron-phonon coupling limit (λ = 0.602), on the other hand, E ∆ is strongly nonmonotonic, increasing with T in the low-T limit. As is especially apparent from panel (c), the intermediate field region between E T and E ∆ has diminished dramatically so that the IMT appears as a fully single-stage transition, as often observed in correlated transition-metal oxides [11,12,27,59,60]. In this limit, the strong lowenergy excitations of acoustic phonons cause the IMT to bypass the CDW-like pseudo-gap state. We therefore make a prediction that the avalanche physics which controls the CDW threshold field E T is manifested in Mott insulators in a switching field E ∆ (T ) that increases with T . We compare our results with the single-stage switching fields measured in TiS 3 [56,57] (blue circles). The agreement between theory and experiment is quite reasonable. Altogether, the E ∆ (T ) behavior in the weakand strong-coupling limits demonstrates a crossover between the quantum and thermal switching mechanisms as temperature is varied.
To conclude, we have established a quantum avalanche mechanism as a generic scenario for the nonequilibrium phase transition in correlated insulators. The proposed mechanism not only resolves the long-standing discrepancy between the observed and predicted switching fields in these materials, but also sheds new light on the origins of the crossover between the quantum and thermal scenarios in resistive transitions. The existence of the pseudo-gapped metallic states may be directly verified experimentally by using transient electrical pulses of varying duration to control the amount of hot-phonon generation. While we have discussed this mechanism with electron-phonon coupling, the avalanche may also arise via coupling to other bosonic excitations, such as the Goldstone modes associated with the order parameter responsible for the opening of the charge gap. Such a scenario may provide a more direct path for the destruction of the correlated gap. Here we have presented a minimal framework for the quantum avalanche that is fundamentally different from the Landau-Zener mechanism. The study of the avalanche with spatial inhomogeneity, and its interplay with disorder, is left for future research.

Methods
Calculation of Self-Energy and Effective Temperature: The many-body calculations in this work are based on the nonequilibrium Green's function technique approximated by the DMFT scheme [37,38,40,44], in which we limit the self-energy to be diagonal in site and orbital indices. Full details on the nonequilibrium DMFT are given in the Supplementary Information. The electron and phonon self-energies, expressed respectively as Π ≶ ep,α (ω) = −2ig 2 ep dω 2π G ≶ αα (r, ω + ω )G ≷ αα (r, ω ), are iterated to convergence. (The factor 2 in the phonon self-energy accounts for spin degeneracy.) The selfenergies without the vertex correction, i.e., Migdal approximation [61], are reasonable since we are in the weakto moderate el-ph coupling limit, as can be seen with less than 10% level shift of the bandedge in Fig. 2(b). Once we achieve full convergence, we compute the distribution functions as f el (ω) = − 1 2 α ImG < αα (r = 0, ω) α ImG R αα (r = 0, ω) , n ph (ω) = 1 2 ImD < (ω) ImD R (ω) .
(12) We define the effective temperature for electrons and phonons in terms of the first moment of the distribution, which correctly reduces to the bath temperature in the equilibrium limit, as with the step-function Θ(x). As shown in Fig. 2b, electronic distribution functions may deviate strongly from the Fermi-Dirac form, in which cases T el provides an approximate measure of nonequilibrium energy excitation.