Subsections

• 6.3.1 Model Implementation
• 6.3.2 Model Evaluation
• 6.3.3 Model Discussion
• 6.3.4 Modifications

6.3 Distribution Function Based Modeling

In our group a new physics-based modeling approach was developed and verified. This approach is based on TCAD device simulators, thereby providing good degradation results on the device level which can be used for reliable device lifetime prediction. The model spans from the microscopic level of interface defect generation up to the device level. The microscopic level considers SP and MP processes relying on the carrier energy distribution functions for electrons and holes. As for the device level, the drift-diffusion device simulation technique allows to extract the device parameter degradation. This gives the unique possibility to define the life-time using the design relevant attributes instead of simply dealing with trap concentrations. This is especially important, since the trap creation is highly localized and the effective device degradation depends on the position of the traps in the channel. The full device simulation step used in this approach gives a more accurate physical link between generated traps and device parameter degradation.

Like in the models by Hess and Bravaix, it is distinguished between single-particle and a multiple-particle components. In the first realizations of this model, only electrons were considered, giving good results for n-channel MOS devices with a channel length of 0.5 µm [17,245]. Later, also secondary generated holes were considered and consolidated by the model [18], as already proposed by other authors, e.g. Moens et al. [246,247]. The degradation mechanism is driven by the acceleration integrals ( $I_\ensuremath{\mathrm{SP}}$ and $I_\ensuremath{\mathrm{MP}}$ ) which are similarly defined as in the approaches by Hess and Bravaix, compare (6.4). These carrier acceleration integrals for the SP and MP process are

$$I_\ensuremath{\mathrm{SP}}= \int_{{\ensuremath{\ensuremath{\mathc... ... v(\ensuremath{\mathcal{E}}) \ensuremath{ \mathrm{d}}\ensuremath{\mathcal{E}},$$ (6.25)

$$I_\ensuremath{\mathrm{MP}}= \int_{{\ensuremath{\ensuremath{\mathc... ... v(\ensuremath{\mathcal{E}}) \ensuremath{ \mathrm{d}}\ensuremath{\mathcal{E}}.$$ (6.26)

In contrast to the other presented modes, here the carrier energy distribution function is explicitly formulated as $f(\ensuremath{\mathcal{E}}).$ Together with the density of states $g(\ensuremath{\mathcal{E}})$ and the carrier velocity $v(\ensuremath{\mathcal{E}})$ the quantity $I(\ensuremath{\mathcal{E}})$ as used in (6.4) can be explicitly calculated. A Keldysh-like reaction cross section

$$\sigma_{\ensuremath{\mathrm{SP}}/\ensuremath{\mathrm{MP}}}(\ensur... ...rm{th}}}_{,\ensuremath{\mathrm{SP}}/\ensuremath{\mathrm{MP}}} )^{p_\mathrm{it}}$$ (6.27)

with the values $p_\mathrm{it}=11$ for SP and MP [223]. All equations exist for electrons and holes with corresponding $f_{n/p}(\ensuremath{\mathcal{E}}),$ $g_{n/p}(\ensuremath{\mathcal{E}}),$ $v_{n/p}(\ensuremath{\mathcal{E}}),$ $\sigma_{n/p}(\ensuremath{\mathcal{E}}),$ and ${\ensuremath{\ensuremath{\mathcal{E}}_\mathrm{th}}}_{,n/p},$ respectively. However, to increase the readability of the equations the subscripts are omitted here.

The interface state generation rate for the SP process in this model is assumed to be described by a first-order chemical reaction with the activation rate $P_{\ensuremath{\mathrm{SP}},\mathrm{act}},$ giving

$$\ensuremath{\ensuremath{\frac{\ensuremath{ \mathrm{d}}n_0}{\ensuremath{ \mathrm{d}}t}}} = - P_{\ensuremath{\mathrm{SP}},\mathrm{act}} n_0 .$$ (6.28)

$n_0$ is the number of non-broken interface Si-H bonds. The activation rate $P_{\ensuremath{\mathrm{SP}},\mathrm{act}}$ is modeled using the attempt frequencies $\nu_{\ensuremath{\mathrm{SP}},n/p}$ and the acceleration integral as

$$P_{\ensuremath{\mathrm{SP}},\mathrm{act}} = \nu_{\ensuremath{\mat... ...thrm{SP}},n} + \nu_{\ensuremath{\mathrm{SP}},p} I_{\ensuremath{\mathrm{SP}},p}.$$ (6.29)

The concentration $n_0$ is related to the total number of interface bonds $N_0$ and the number of broken interface bonds $N_\ensuremath{\mathrm{SP}},$ i.e. the number of interface traps activated due to the SP process, as

$$N_\ensuremath{\mathrm{SP}}= N_0 - n_0 ,$$ (6.30)

Together with the boundary condition of $n_0(t=0s) = N_0,$ the solution of (6.28) using (6.30) gives the number of interface traps created due to the SP process over time as

$$N_\ensuremath{\mathrm{SP}}= N_0 \Big( 1 - \exp (-P_{\ensuremath{\mathrm{SP}},\mathrm{act}} t) \Big) .$$ (6.31)

Like in the Bravaix approach, a truncated harmonic oscillator is employed to describe the bond energetics. However, we write the rate equations for the last bonded state in a different fashion: defined as

$$\ensuremath{\ensuremath{\frac{\ensuremath{ \mathrm{d}}n_0}{\ensuremath{ \mathrm{d}}t}}} = P_d n_1 - P_u n_0$$ (6.32)

$$\ensuremath{\ensuremath{\frac{\ensuremath{ \mathrm{d}}n_i}{\ensuremath{ \mathrm{d}}t}}} = P_d (n_{i+1} - n_i) - P_u ( n_i - n_{i-1})$$ (6.33)

$$\ensuremath{\ensuremath{\frac{\ensuremath{ \mathrm{d}}n_{N_l}}{\... ...ilde{P}_{\ensuremath{\mathrm{MP}},\mathrm{pass}} N_\ensuremath{\mathrm{MP}}^2 .$$ (6.34)

The rate equation for the topmost energy level (6.34) consists of four processes: the transition to the topmost level from the level below, the transition from the topmost level down to the level below, the desorption rate $P_{\ensuremath{\mathrm{MP}},\mathrm{act}}$ for the actual bond breaking, and the passivation rate $P_{\ensuremath{\mathrm{MP}},\mathrm{pass}}.$ The rates are defined following the Arrhenius relation,

$$P_{\ensuremath{\mathrm{MP}},\mathrm{pass}} = \nu_{\ensuremath{\ma... ...{\mathcal{E}}_{\mathrm{pass}}}{\ensuremath{\mathrm{k_B}}T_{\mathrm{L}}} \right)$$ (6.35)

$$P_{\ensuremath{\mathrm{MP}},\mathrm{act}} = \nu_{\ensuremath{\mat... ...\mathcal{E}}_{\mathrm{emi}}}{\ensuremath{\mathrm{k_B}}T_{\mathrm{L}}} \right) ,$$ (6.36)

using the energy barriers as shown in Fig. 6.2 and the attempt frequencies $\nu_{\ensuremath{\mathrm{MP}},\mathrm{pass}}$ and $\nu_{\ensuremath{\mathrm{MP}},\mathrm{act}}.$ The passivation reaction depends on the existence of a dangling bond and the existence of a hydrogen atom which can passivated, making this a second-order reaction. The density of dangling bonds is equal to the density of generated traps. No initial free hydrogen is assumed, so the available number of hydrogen atoms corresponds to the number of generated traps. This leads to the quadratic dependence $N_\ensuremath{\mathrm{MP}}^2$ in (6.34). To satisfy the dimensionality, the passivation rate is written as $\widetilde{P}_{\ensuremath{\mathrm{MP}},\mathrm{pass}} = P_{\ensuremath{\mathrm{MP}},\mathrm{pass}}/N_0.$

The great disparity between the times between bond excitation and decay on one side and hydrogen bond dissociation dictates that these two processes can be treated quasi-separately [223,245]. The steady-state within the oscillator is assumed, which is established momentarily compared to the hydrogen release/absorption. By neglecting the last two terms in (6.34) the truncated harmonic oscillator is decoupled from the bond breaking process and the occupation dynamics give

$$n_i = \left( \frac{P_u}{P_d} \right) ^i n_0 .$$ (6.37)

In a second step, as the harmonic oscillator is considered in steady-state, for the bond dissociation the first two terms in (6.34) can be omitted. While solving this differential equation two assumptions were used. The first is similar to that in [223] and claims that the concentration of virgin bonds does not change drastically. The second is the boundary condition that initially all bonds are virgin. The trap concentration due to the MP process over time becomes

$$N_\ensuremath{\mathrm{MP}}= N_0 \left[ \frac{P_{\ensuremath{\math... ...g( 1-\exp (- P_{\ensuremath{\mathrm{MP}},\mathrm{act}} t) \Big) \right]^{1/2} .$$ (6.38)

For weak stresses and/or short stress times, meaning $P_{\ensuremath{\mathrm{MP}},\mathrm{act}} t \ll 1,$ a Taylor expansion gives the approximation $1-\exp (-P_{\ensuremath{\mathrm{MP}},\mathrm{act}} t) \approx P_{\ensuremath{\mathrm{MP}},\mathrm{act}} t$ and one receives the square root time dependence like in the Bravaix model, see (6.19).

The probabilities $P_u$ and $P_d$ for the excitation and decay of the Si-H bond are defined similarly to equations (6.20) and (6.21) in the Bravaix model, using the acceleration integrals, the phonon frequency $w_e$ and the distance between the oscillator levels $\ensuremath{\hbar \omega}$ :

$$P_u = \nu_{\ensuremath{\mathrm{MP}},n} I_{\ensuremath{\mathrm{MP}... ...rac{\ensuremath{\hbar \omega}}{\ensuremath{\mathrm{k_B}}T_{\mathrm{L}}} \right)$$ (6.39)

$$P_d = \nu_{\ensuremath{\mathrm{MP}},n} I_{\ensuremath{\mathrm{MP}},n} + \nu_{\ensuremath{\mathrm{MP}},p} I_{\ensuremath{\mathrm{MP}},p} + w_e .$$ (6.40)

The excitation corresponds to a phonon absorption, meaning stimulation of the bond leading to a step up on the energy levels of the oscillator, and the decay corresponds to a phonon emission, meaning a step down one level.

The total trap concentration includes the two concurrent SP- and MP-processes and is combined using probabilities to balance the two processes as

$$N_\mathrm{it} = p_\ensuremath{\mathrm{SP}}N_{SP} + p_\ensuremath{\mathrm{MP}}N_{MP} .$$ (6.41)

Due to the position dependence of the carrier energy distribution function, and therefore the acceleration integrals, the resulting concentration of generated traps also depends on the position.

The interface states created impact the device performance by trapping and de-trapping charges. This trapping process can be simulated using an interface Shockley-Read-Hall (SRH) modeling approach like shown in Section 4.2.2. Only the trapped charges influence the output characteristics of the degraded device. In this model the interface trap is characterized by its charge state evaluated for the given coordinate along the interface and particular stress/operating conditions. The charge influences the electrostatic potential and degrades the carrier mobilities due to Coulombic scattering mechanisms. A simple interface charge induced mobility reduction model has been presented in Section 4.2.1.

6.3.1 Model Implementation

The current implementation of the model spans from the microscopic trap generation level to the device operation level. The degradation is modeled using partly existing tools and partly newly implemented modules which are used to calculate the microscopic spatially distributed damage. The implementation consists of three components shown in Fig. 6.4. A Monte Carlo simulator is used to calculate the distribution functions, the degradation is calculated using the equations presented in the previous section, and the impact on the device performance is analyzed using a drift-diffusion simulator.

Figure 6.4: The flow chart of the distribution function based hot-carrier degradation model. The Monte Carlo (MC) simulator delivers the spatially varying carrier energy distribution function $f(\ensuremath {\mathcal {E}},x)$ for given stress conditions. This information is then used to calculate the density of generated traps for all time steps. For each time step the output characteristic of the device at operation condition is calculated using the macroscopic device simulation to evaluate the degradation, i.e. $\Delta \ensuremath {V_{\mathrm {T}}}(t).$

To calculate the carrier energy distribution function the Monte Carlo method is applied using the full-band device simulator MONJU [16]. To capture the hole contribution, impact-ionization needs to be considered. Simulations are performed on the device under test for a certain stress condition, i.e. drain- and gate-voltage. The output of this model is the electron and hole energy distribution functions along the $\mathrm{Si-SiO_2}$ interface. Due to the stochastic nature of the Monte Carlo method, the convergency behavior is poor and long simulation times are required. This concerns especially high-energy tails of the distribution function, where the number of carriers is low. Therefore the computational process would take disproportionally long simulation times to obtain smooth, i.e. noiseless, results and therefore commonly noisy data has to be used.

In the next step, the distribution functions $f(\ensuremath{\mathcal{E}})$ for electrons and holes are transferred to the degradation module (see Fig. 6.4) which calculates the trap density as a function of the lateral coordinate for each stress time step. The microscopic process of the generation for MP and the SP processes is described using the acceleration integrals in equations (6.25) and (6.26) and the trap generation in equations (6.31) and (6.38). As a result, the interface trap density at each position over time is available. In the current model implementation, the initially calculated distribution function is used for all time steps, because a self-consistent re-calculation using the Monte Carlo method would be extremely time-consuming. However, the change of the distribution function during the degradation has only a small impact on the results [248]. For the calculation of the charge trapped by interface states, the SRH equations have to be solved. The operation conditions in the reference measurements used for model validation ($V_G = 2$ V) suggest that the traps are charged throughout the simulation. Hence, all interface traps are assumed as fixed interface charges.

Finally, the position dependent interface charge is then introduced into our multi-purpose device simulation tool MINIMOS-NT [120]. To gain reasonable simulation times, the drift-diffusion transport model as described in Section 4.1.2 is used.

6.3.2 Model Evaluation

The presented model is evaluated by comparing simulation results employing a set of 5 V n-MOS transistors fabricated using the same architecture but with different channel lengths. The devices are part of a 0.35 µm high-voltage mixed signal process by ams. The channel lengths are 0.5, 1.2, and 2.0 µm. The 0.5 µm device is depicted in Fig. 6.5. The devices are stressed at a gate voltage of $\ensuremath{V_{\mathrm{GS}}}= 2.0$ V and a drain voltage of $\ensuremath{V_{\mathrm{DS}}}= 6.25$ V at $25$ $^{\circ}$ C. The stress was measured for a period of 10,000 s.

Figure 6.5: Geometry and net doping concentration of the 0.5 µm n-MOSFET used for simulations. The device structure was generated using process simulation. Legend omitted due to non-disclosure agreement.

The simulation setup was used as described previously and the model was calibrated to fit the drain current degradation operating with a single parameter set. For the fitting process it was considered, that both acceleration integrals have the same functional structure. They differ only in the prefactors $\nu$ (see (6.29), (6.39), and (6.40)). This leads to $I_\ensuremath{\mathrm{SP}}= I_\ensuremath{\mathrm{MP}}$ for each carrier type. Hence, the resulting degradation in Fig. 6.6 shows a good agreement between simulations and measurements.

Figure 6.6: $I_{\mathrm {D,lin}}$ degradation in devices with different channel lengths. The simulation result using only electrons, only holes, and using both components are compared to the measurement results [18].

Figure 6.7: The figures show the acceleration integral for electrons and holes in the 0.5 (a), 1.2 (c), and 2.0 µm (e) devices and the total $N_\mathrm {it}$ concentration as well as contributions to $N_\mathrm {it}$ produced only by holes also for 0.5 (b), 1.2 (d), and 2.0 µm (f) channel MOSFETs [18]. All the quantities are plotted as functions of the lateral coordinate. The source corresponds to the abscissa origin. The $N_\mathrm {it}$ concentrations are shown after 10 and 10,000 s stress times.

For a proper model evaluation, the contributions from electron and hole components are plotted separately. As can be seen, the 0.5 µm device is the only one in which the full degradation can be represented using only the electron induced degradation component. This tendency was also observed in [18] and is already mentioned in Section 6.3. This can be investigated in more detail looking at the acceleration integral and the interface state density $N_\mathrm {it}$ along the channel, both are shown in Fig. 6.7.

First, one can see that the worst degradation happens at the drain (right) end of the device, where the electrons accelerated along the channel have reached the highest energy. This part of the degradation is caused by the electron SP process. Note that the damage is partly located outside the channel. The hole contribution is caused by secondary holes generated by impact-ionization. Consequently, they are accelerated towards the source side (left), gain energy, and reach the maximum energy within the channel area. One can clearly see that the hole component of the acceleration integral is always shifted towards the source and is situated within the channel area. This explains why the hole contribution, which is characterized by a much lower portion of $N_\mathrm{it},$ still plays a relevant role for the entire $\ensuremath{I_{\mathrm{D,lin}}}$ degradation.

An important analysis using this model has been done by Starkov et al. [249] comparing measurements and simulation results for the worst-case hot-carrier degradation. In the n-MOS transistor, the substrate current $\ensuremath{I_{\mathrm{Sub}}}$ is utilized as a criterion for the worst-case condition, reflecting the impact-ionization generation [223]. For the p-MOS transistor the gate current has been used as an indicator for maximum hot-carrier degradation [224]. In the simulations, the maxima of the acceleration integrals have been used to compare the severity of the ongoing degradation process. The impact has been compared with measurement data from the 0.5 µm devices over a wide range of varying $\ensuremath {V_{\mathrm {DS}}}$ and $\ensuremath{V_{\mathrm{GS}}}.$ The results are depicted as a color map in Fig. 6.8. The figure reflects the correct tendency of the bias dependent hot-carrier stress.

Figure 6.8: Color map plotted over the gate ( $\ensuremath {V_{\mathrm {GS}}}$ ) and drain ( $\ensuremath {V_{\mathrm {DS}}}$ ) voltages. The measured substrate current of the n-MOS (a) and of the gate current of the p-MOS (b) transistor are compared to the maximum value of the acceleration integral, shown in (c) and (d), respectively. (Figures taken from [249])

6.3.3 Model Discussion

This model tries to address the whole hierarchy of physical phenomena, taking information from the carrier energy distribution function to model the microscopic degradation processes, generating interface traps, and then simulating the influence of the traps on the device behavior. The evaluations performed on various levels show good agreement in comparison to the measurement data.

Unfortunately, the complexity of this approach, which considers the distribution function for electrons and holes, results in severe limitations regarding simulation time and flexibility of the simulation setup. The Monte Carlo simulations are very time consuming and it is only reasonable to calculate the set of distribution functions once, i.e. for the virgin device. Changes of the stress conditions during the stress cycle cannot be captured straightforward, and hence a new Monte Carlo simulation would be required for each step. Another method to calculate the carrier distribution functions is the spherical harmonics expansion [124,125]. A device simulator based on this approach is currently being developed at our institute and first results have been published [250]. However, while this approach is computationally much more efficient than Monte Carlo simulations, it still requires a huge amount of computers working memory, thereby posing some limitations in simulations of real devices. Therefore, other simplified approaches to overcome these problems are worthy to be discussed. A method which is solely based on drift-diffusion is the topic of the upcoming section Section 6.4.

A weakness of the current model might be that so far only two limiting cases related to hydrogen desorption are considered, i.e. from the ground state (SP mechanism) and from the topmost energy level (MP mechanism) of the oscillator. An extension considering hydrogen desorption from all energy levels would represent the actual dissociation process more precisely. For this purpose, the rate equations must be extended and some new barriers need to be introduced for this. This approach has already been suggested by McMahon et al. [231].

Apart from the bond breakage at the interface, the model does not yet contribute to possible oxide bulk traps or charges [251]. However, up to this point the model seems to be able to represent the experimental data reasonably well and the matter whether bulk oxide traps contribute to hot-carrier degradation or not is unresolved. In fact, in the intimately related degradation mode, in bias temperature instability, just trapping/de-trapping in the oxide bulk is responsible for the recoverable component of the damage [252,253,254]. However, hot-carrier degradation demonstrates, if at all, rather weak recovery. This suggest that bulk oxide traps do not play a substantial role in hot-carrier degradation. At the same time, recent measurements have proven the existence of the threshold voltage turn-around effect that can be explained by oxide charges [248]. The next improvement steps can, therefore, be oriented towards the inclusion of oxide traps into the hot-carrier degradation model.

6.3.4 Modifications

There are some modification of our model that would make the formalism more complete. First, also the single particle process can be considered to undergo a passivation reaction, introducing the rate $P_{\ensuremath{\mathrm{SP}},\mathrm{pass}},$ extending (6.28) to

$$\ensuremath{\ensuremath{\frac{\ensuremath{ \mathrm{d}}n_0}{\ensu... ...tilde{P}_{\ensuremath{\mathrm{SP}},\mathrm{pass}} N_\ensuremath{\mathrm{SP}}^2.$$ (6.42)

The formulation uses $\widetilde{P}_{\ensuremath{\mathrm{SP}},\mathrm{pass}}$ instead of $P_{\ensuremath{\mathrm{SP}},\mathrm{pass}}$ to keep the units consistent.

The second modification suggested, concerns the fundamental splitting between the SP and the MP process introduced in (6.41). By using just a single trap density $N_\mathrm {it}$ instead of the two separate ones ( $N_\ensuremath{\mathrm{SP}}$ and $N_\ensuremath{\mathrm{MP}}$ ), it is possible to combine the rate equations (6.32) - (6.34) and (6.42) leading to

$$\ensuremath{\ensuremath{\frac{\ensuremath{ \mathrm{d}}n_0}{\ensu... ...}} n_0 + \widetilde{P}_{\ensuremath{\mathrm{SP}},\mathrm{pass}} N_\mathrm{it}^2$$ (6.43)

$$\ensuremath{\ensuremath{\frac{\ensuremath{ \mathrm{d}}n_i}{\ensuremath{ \mathrm{d}}t}}} = P_d (n_{i+1} - n_i) - P_u ( n_i - n_{i-1})$$ (6.44)

$$\ensuremath{\ensuremath{\frac{\ensuremath{ \mathrm{d}}n_{N_l}}{\... ...N_l} + \widetilde{P}_{\ensuremath{\mathrm{MP}},\mathrm{pass}} N_\mathrm{it}^2 .$$ (6.45)

The model as presented evaluates the concentration of the broken bonds $N_\mathrm {it}$ in the range of $[0 \dots N_0].$ For a convenient handling, a fraction $f = \left.N_\mathrm{it}\middle/N_0\right.$ can be used defining the relative bond breakage. For this, one can define the total number of passivated bonds (on all energy levels of the oscillator) as

$$N = \sum_{i=0}^{N_l} n_i = n_0 \sum_{i=0}^{N_l} \left( \frac{P_u}{P_d} \right) ^ i = n_0 k,$$ (6.46)

where the geometric series describing $k$ can be evaluated as

$$k = \frac{\displaystyle \left( \frac{P_u}{P_d} \right) ^{N_l+1} -1}{\displaystyle \left( \frac{P_u}{P_d} \right) - 1} .$$ (6.47)

At time $t=0$ s no interface traps have been generated and the total number of bound traps is, therefore, equal to the maximal available number of bonds $N_0.$ The relation between $N,$ $N_0,$ and $N_\mathrm {it}$ is thus

$$N = N_0 - N_\mathrm{it}$$ (6.48)

and inserting $f = \left.N_\mathrm{it}\middle/N_0\right.$ gives

$$N = \left( 1-f \right) N_0 .$$ (6.49)

Considering the oscillator in steady-state and by combining the SP and MP rates from equations (6.31) and (6.38), the resulting bond generation rate can then be written as

$$\ensuremath{\ensuremath{\frac{\ensuremath{ \mathrm{d}}f}{\ensure... ...d}{P_u} \right) ^{N_l} P_{\ensuremath{\mathrm{MP}},\mathrm{pass}} \right] f^2 .$$ (6.50)

Note that the rates were substituted as $P_{\ensuremath{\mathrm{SP}},\mathrm{pass}} = \tilde{P}_{\ensuremath{\mathrm{SP}},\mathrm{pass}} N_0$ and $P_{\ensuremath{\mathrm{MP}},\mathrm{pass}} = \tilde{P}_{\ensuremath{\mathrm{MP}},\mathrm{pass}} N_0.$