(image) (image) [Previous] [Next]

Modeling of Defect Related Reliability Phenomena
in SiC Power-MOSFETs

4.6 Compact Physics Framework (Comphy)

The development of the Compact Physics Framework (Comphy) [203, 211] aimed at reducing the computational effort needed to extract the large number of parameters as required for detailed BTI models, i.e. the four-state NMP or the GSHR model which are partly implemented in commercial TCAD frameworks [232, 233, 234]. However, it has been demonstrated that the mean BTI degradation observed in large area MOSFETs can be explained physically with a two-state NMP model, which requires only a reduced defect parameter set. Comphy and its implementation of the two-state model have been initially applied to reproduce the degradation for a large number of process splits with gate stacks employing SiO2 and HfO2 including different gate contact materials [235]. As a large ensemble of defects can be handled efficiently, it makes the framework particularly suitable for calculating BTI in SiC-MOSFETs, allowing for the extraction of physical defect parameters by using established material parameters. In the following, the main physical models for providing the electrostatic quantities required to calculate the charge transfer kinetics and transient ∆Vth will be outlined. Also, the extension of the framework for efficiently calculating gate leakage currents employing the models presented in the previous section will be presented.

4.6.1 Electrostatic Quantities

In order to calculate the effective trap level of a defect relative to the channel carrier reservoirs, the relation of the applied gate bias to the surface potential at the channel / oxide interface needs to be known. Therefore, by assuming a uniform doping concentration and charge neutrality exists deep in the bulk semiconductor (far away from the interface), the approximation of the surface charge \( Q_\mathrm {s} \) as a function of the surface potential \( \varphi _\mathrm {s} \) can be used, which reads [31]

(4.57) \{begin}{multline} Q_\mathrm {s}\left (\varphi _\mathrm {s}\right ) = \pm \frac {\sqrt {2}k_\mathrm {B}T}{q_0 L_\mathrm {D,e}} \\ \sqrt { \frac {p_0}{n_0} \left (
e^{-q_0\varphi _\mathrm {s} / \left ( k_\mathrm {B} T \right ) } + q_0 \varphi _\mathrm {s} / \left (k_\mathrm {B} T \right ) - 1 \right ) + \left ( e^{q_0\varphi _\mathrm {s} / \left (
k_\mathrm {B} T \right ) } - q_0 \varphi _\mathrm {s} / \left ( k_\mathrm {B} T \right ) - 1 \right )}. \label {equ:surf_charge} \{end}{multline}

Thereby, \( L_\mathrm {D,e} = \sqrt {k_\mathrm {B} T \varepsilon _\mathrm {s} / \left ( n_0 q_0^2 \right )} \) denotes the electron Debye length, \( \varepsilon _\mathrm {s} \) the semiconductor permittivity and \( n_0,p_0 \) the carrier concentrations at thermal equilibrium, which are derived from the doping concentrations \( N_\mathrm {A}, N_\mathrm {D} \) and the band gap \( E_\mathrm {G} \) [203]. In the case of a charge free insulator, the potential drop across the gate stack is given by

(4.58) \{begin}{align} V_\mathrm {G} = \frac {Q_{s}\left (\varphi _\mathrm {s}\right ) W L}{C_\mathrm {ox}} + \varphi _\mathrm {s} + \frac {\Delta E_\mathrm {W}} {q_0} \label
{equ:vg_pot_drop_ideal} \{end}{align}

with the work-function difference \( \Delta E_\mathrm {W} \) between the gate and channel materials. As the inverse relation \( \varphi _\mathrm {s} \left ( V_\mathrm {G} \right ) \) is required to calculate the channel electrostatic quantities like Fermi-Level \( E_\mathrm {F} \) and carrier densities \( n,p \) with employing the Joyce-Dixon approximation [236], (4.58) is solved numerically for \( \varphi _\mathrm {s} \) employing an iterative scheme. Figure 4.9 shows a comparison of the approximations used in Comphy to a full numerical solution of the Poisson equation across a poly-Si/SiO2/SiC stack with a Finite-Volume method and Fermi-Dirac statistics.


Figure 4.9: Comparison of the numeric solution of the Poisson equation for \( \varphi _\mathrm {s}\left (V_\mathrm {G}\right ) \) of a n-MOSCAP with a charge free oxide shows good agreement with the efficient approximations used in Comphy (left). The full numeric solution of the Poisson equation employing Fermi-Dirac statistics is shown in a band diagram for strong accumulation (right).

4.6.2 Threshold Voltage Shift

With the surface potential derived in the previous section, the Fermi-Levels in the channel and gate are known and can be used to calculate the transition rates (4.29) and (4.30) for each input tuple \( \left ( V_\mathrm {G}, t, T \right ) \). These allow the computation of the transient defect occupancy \( p(t) \) using (4.9). By using a simple charge sheet approximation [237, 69], the perturbation of the potential due to \( N \) defects reads

(4.59) \{begin}{align} \Delta \varphi = - q_0 C_\mathrm {ox} \sum _i^{N} p_i \bigg ( 1 - \frac {x_\mathrm {T,i}}{t_\mathrm {ox}} \bigg ) \label {equ:CSA} \{end}{align}

with \( x_\mathrm {T} \) being the distance of the defect from the channel/oxide interface.

Often a solution of the occupation for AC signals is sought in order to calculate ∆Vth within circuit simulations. Therefore, an efficient numerical solution for the occupations given by a two-state Master equation for arbitrary shaped periodic two-level bias signals has been proposed by Giering et. al [238]. For the special case of a periodic digital gate AC signal with period \( t_\mathrm {p} = 1/f = t^\mathrm {H}+t^\mathrm {L} \) at biases \( V_\mathrm {G}^\mathrm {H} \) at high level and \( V_\mathrm {G}^\mathrm {L} \) at low level an analytical expression, evaluated after \( n \) periods, can be found by [203]

(4.60) \{begin}{align} p_2\left (t_0 + n t_\mathrm {p} \right ) = p_2\left (t_0\right ) P_1\left ( t_0 + t_\mathrm {p}, t_0 \right )^n + \frac {1-P_1\left ( t_0 + t_\mathrm {p}, t_0
\right )^n}{1-P_1\left ( t_0 + t_\mathrm {p}, t_0 \right )} P_2\left ( t_0 + t_\mathrm {p}, t_0 \right ) \label {equ:AC} \{end}{align}


(4.61–4.62) \{begin}{align} P_1\left ( t_0 + t_\mathrm {p}, t_0 \right ) &= \mathrm {exp} \bigg ( -\left (k_\mathrm {12,H} + k_\mathrm {21,H} \right ) t^\mathrm {H} - \left
(k_\mathrm {12,L} + k_\mathrm {21,L} \right ) t^\mathrm {L} \bigg ),\\ P_2\left ( t_0 + t_\mathrm {p}, t_0 \right ) &= - \frac {k_\mathrm {12,H} }{\left (k_\mathrm {12,H} + k_\mathrm
{21,H}\right )} \bigg ( \mathrm {e}^{- \left (k_\mathrm {12,H} + k_\mathrm {21,H} \right ) t^\mathrm {H}} - 1 \bigg ) \notag \\ &\phantom {=} - \frac {k_\mathrm {12,L} }{\left (k_\mathrm
{12,L} + k_\mathrm {21,L}\right )} \bigg ( \mathrm {e}^{- \left (k_\mathrm {12,L} + k_\mathrm {21,L} \right ) t^\mathrm {L}} - 1 \bigg ) \mathrm {e}^{- \left (k_\mathrm {12,H} + k_\mathrm
{21,H} \right ) t^\mathrm {H}}. \label {equ:AC_analytic} \{end}{align}

The availability of an analytic expression to calculate the defect occupation for periodic AC signals allows for efficient extrapolation of the device degradation ∆Vth and therefore \( \Delta R_\mathrm {on} \) under operating conditions as will be shown in Section 5.1.6.

4.6.3 Gate Leakage Current Computation


Figure 4.10: To account for the amorphous nature of the most common gate insulators, the defects in Comphy are sampled in three dimensional slabs, with boundaries on the gate contact and channel interface. The volume \( V = a^2 t_\mathrm {diel} \) of the slabs is calculated to contain a given mean number of defects \( \overline {N} \) for a density \( N_\mathrm {T} \), with the actual number of traps \( N \) being Poisson distributed within the \( M \) slabs. (taken from [CSJ2])

Irrespective of the fact that the Poisson equation is only solved in one dimension within the Comphy framework, the derivation of the correct spatial tunneling distance \( d_\mathrm {T} \) requires a three dimensional defect distribution. Therefore, based on the input quantities defect density NT and average defect number \( \overline {N} \) that should be sampled, a volume with \( V_i = A_i t_\mathrm {diel} \) is defined that fulfills \( V = \overline {N} / N_\mathrm {T} \). The dielectric thickness is fixed by the gate stack input quantities, and \( A_i = a^2 \) a square area. Within this volume, the defects are sampled using a Monte-Carlo scheme according to their uniform spatial and normally distributed energetic parameters. To avoid low defect densities, in particular those which lead to \( d_\mathrm {T} \) below 1 nm, defects that result in smaller \( d_\mathrm {T} \) to any neighboring defect are removed from the drawn sample and re-drawn, i.e. rejection sampling is implemented. The actual number of defects \( N \) within a slab is further sampled based on a Poisson distribution, as such is observed in small area devices [239, 152], with probability of \( N \) defects within the slab given by

(4.63) \{begin}{align} P(N) = \frac {\mathrm {e}^{-\overline {N}}\overline {N}^N}{N!}. \{end}{align}

For the energetic defect parameters (ET ,ER ) samples are drawn based on a normal distribution with ( \( \sigma _{E_\text {T}} \), \( \sigma _{E_\mathrm {R}} \)). The expectation values and variations of the current density and threshold voltage shifts can then be calculated from repeating the simulation on \( M \) slabs, as schematically shown for the simulation space in Figure 4.10.

The computation of the currents in steady state is based on solving the coupled system of equations (4.10) for zero derivative, which reads

(4.64) \{begin}{align} k_{\mathrm {in},i} \big (\vec {f}\big ) \left ( 1 - f_i \right ) = k_{\mathrm {out},i} \big (\vec {f}\big ) f_i \label {equ:coupled_master_equation_ss}

with the total in and out rates \( k_{\mathrm {in},i} \) and \( k_{\mathrm {out},i} \) at each defect \( i \), which consist of reservoir and defect interaction rates. Therefore, for each step \( k \) with the input tuple ( \( t_k \), \( V_\mathrm {G,k} \), \( T_k \)) all rates as described in Section 4.4 are computed and (4.64) is solved with a Newton scheme. For the transient case, the equation is discretized with an implicit Euler scheme, resulting in the occupations as

(4.65) \{begin}{align} f_{i,k} = \left [ \left ( 1-f_{i,k} \right ) k_{\mathrm {in},i,k-1}\big (\vec {f}_{k-1}\big ) - f_{i,k} k_{\mathrm {out},i,k-1}\big (\vec {f}_{k-1}\big )
\right ] \left (t_k - t_{k-1}\right ) + f_{i,k-1}. \label {equ:impl_eul} \{end}{align}

The discrete equation can then again be solved with a Newton scheme and the resulting occupations together with the individual rate contributions. With the occupations, the TAT current (4.56) can be evaluated ant the total gate current sums up from the contributions of the TAT and band-to-band current \( I_\mathrm {G,TE} \) (Tsu-Esaki) (4.48) to a total gate current

(4.66) \{begin}{align} I_\mathrm {G,tot} = I_\mathrm {G,TAT} + I_\mathrm {G,TE}. \label {equ:gate_leak_tot} \{end}{align}

With the full modeling framework at hand to compute both TAT currents via single and multiple steps, the significance of a charge transition between two contacts via multiple defects in terms of its contribution to a total leakage current will be evaluated in the following section.