5.2 A New Approach to Solve Local Oxidation

Using AMIGOS a new model has been developed based on the model of E. Rank [Ran90]. The key idea is the description of the local oxidation as a three component thermodynamic process involving silicon, silicon dioxide and oxidant molecules. This results in a reactive layer of finite width in contrast to the sharp interface between silicon and silicon dioxide in the conventional formulation. The numerical approximation takes advantage of this description in a finite element approach which models silicon, silicon dioxide and the reactive layer together, thus avoiding the necessity to track the interface with element edges. The smooth transition zone is a means to regularize the mathematical free boundary problem and can be selected in a way, that numerically seen, the same results appear as in case of sharp interface formulations. To distinguish between different materials a method similar to the level set method was chosen to keep the transition zone as small as possible (usually the transition zone is limited to a single element).

Figure 5.2: Domain and boundary settings

For the definition of the model consider Fig. 5.2 as computation domain $\Omega$ which consists of a pure silicon dioxide range $\Omega_{1}^{}$, an interface range $\Omega_{2}^{}$ with a mixture of silicon and silicon dioxide, a pure silicon range $\Omega_{3}^{}$ and a nitride mask $\Omega_{4}^{}$ that is defined on a separate mesh and is connected to $\Omega_{1}^{}$ via boundary $\Gamma_{4}^{}$ to transmit mechanical displacements. For the nitride mask an elastic model is used to calculate its stress-strain contribution. To describe the different phases of oxygen within the domain $\Omega_{1}^{}$  $\cup$  $\Omega_{2}^{}$  $\cup$  $\Omega_{3}^{}$ a generation/recombination rate of oxygen

$$\left(\vphantom{1-\eta(x,t)}\right. \eta \left.\vphantom{1-\eta(x,t)}\right)$$ (5.19)

is defined, where

$$\eta \left(\vphantom{\frac{C_{SiO_2}\left(x,t\right)}{C_{Si_0}}}\right. {\frac{C_{SiO_2}\left(x,t\right)}{C_{Si_0}}} \left.\vphantom{\frac{C_{SiO_2}\left(x,t\right)}{C_{Si_0}}}\right)$$ (5.20)

is a function of a normalized silicon dioxide concentration related to the C_Si0 concentration of silicon in pure crystal. $\eta$ varies between one (pure silicon dioxide) and zero (pure silicon). The function is calculated using

$${\frac{{\rm atan}(ac \cdot x)}{\pi}} \rightarrow \infty$$ (5.21)

which describes a jump in the material in case of a huge parameter ac (Fig. 5.3). In contrast to previous approaches the transition zone can now be calculated as an immediate jump from zero to one within one element. The advantage of the selected level set function is that a correct derivative of the system matrix can still be calculated that is essential for good convergence behavior of the discretized system.

Figure 5.3: Parameter dependent level set function

The generation of silicon dioxide itself is handled by a generation rate:

$${\frac{\partial C_{SiO_2}}{\partial t}}$$ (5.22)

The free oxidant diffusion in $\Omega_{1}^{}$  $\cup$  $\Omega_{2}^{}$  $\cup$  $\Omega_{3}^{}$ is described by

$${\frac{\partial C_{O}}{\partial t}} \left[\vphantom{D\left(\eta(x,t)\right)\cdot \left( grad\left( C_{O}\right) -\frac{C_O}{m}\: grad\left( m\right) \right) }\right. \left(\vphantom{\eta(x,t)}\right. \eta \left.\vphantom{\eta(x,t)}\right) \left(\vphantom{ grad\left( C_{O}\right) -\frac{C_O}{m}\: grad\left( m\right) }\right. \left(\vphantom{ C_{O}}\right. \left.\vphantom{ C_{O}}\right) {\frac{C_O}{m}} \left(\vphantom{ m}\right. \left.\vphantom{ m}\right) \left.\vphantom{ grad\left( C_{O}\right) -\frac{C_O}{m}\: grad\left( m\right) }\right) \left.\vphantom{D\left(\eta(x,t)\right)\cdot \left( grad\left( C_{O}\right) -\frac{C_O}{m}\: grad\left( m\right) \right) }\right]$$ (5.23)

with the boundary condition

$${\frac{\partial C_{O}}{\partial n}} \Gamma_{1}^{} {\frac{\partial C_{O}}{\partial n}} \Gamma_{2}^{} \Gamma_{4}^{}$$ (5.24)

The mobility of oxygen is strongly influenced by the amount of the generated silicon dioxide since compounded oxygen atoms are assumed to be immobile. The amount of generated silicon dioxide itself depends on the local concentration of already generated oxide as well as of free oxygen which reacts with silicon immediately. This effect leads to an enhancement of oxidation at the very beginning, since nearly all oxide atoms are reacting with silicon. Later on, the reaction rate slows down due to the decreasing amount of reaction partners. Now, an increasing number of oxygen molecules must diffuse through the generated oxide layer towards the interface, before the molecules can finally be incorporated into the silicon and cause the interface movement. Furthermore, the so built-in oxygen involves a local stress which is imprinted into the system. In dependence of further treatment of the mechanical behavior a volume dilatation can be calculated.

The preferred approximation of mechanical behavior of silicon and silicon dioxide can be reached by a visco-elastic model. Therefore, it suggests itself to use a Maxwell body (Fig. 5.1) to solve the mechanical equations. The advantage of the Maxwell body is, that it can handle both extreme fields, elastic as well as viscous material and of course a mixture of both - a visco-elastic behavior. Again, the level set function $\eta$ (5.20) is used to distinguish between different mechanical models:

• elastic behavior in the silicon range $\Omega_{1}^{}$
• visco-elastic behavior in the interface range $\Omega_{2}^{}$
• viscous behavior in the pure silicon dioxide layer $\Omega_{3}^{}$

The used strain relation is based on Hook's law (5.15) that is expressed in a way that the dilatational components of stress, which involve a volumetric expansion, and the deviatoric part, which only accounts for shape modification, are decoupled. To account for the nonlinear material behavior the modulus of rigidity G and the compressibility $\chi$ are modeled as functions of $\eta$.

For the volumetric expansion we solve the equilibrium condition

$$\left(\vphantom{\int\limits_V{{\cal{L}^T}\cdot{\cal{D}}\left(\eta(x,t)\right)\cdot{\cal{L}}}\cdot\:dV}\right. \int\limits_{V}^{} \cal {L} \cal {D} \left(\vphantom{\eta(x,t)}\right. \eta \left.\vphantom{\eta(x,t)}\right) \cal {L} \left.\vphantom{\int\limits_V{{\cal{L}^T}\cdot{\cal{D}}\left(\eta(x,t)\right)\cdot{\cal{L}}}\cdot\:dV}\right) \int\limits_{V}^{} \cal {L} \cal {D} \left(\vphantom{\eta(x,t)}\right. \eta \left.\vphantom{\eta(x,t)}\right) \varepsilon_{o}^{}$$ (5.25)

where u, $\varepsilon_{o}^{}$, $\cal {L}$ and $\cal {D}$ represent the displacement vector, the strain caused by silicon dioxide generation, the mechanical operator defined as {$\varepsilon$} = $\cal {L}$ ^. {u} and the elasticity matrix, respectively. The right-hand side of (5.25) can be interpreted as an energy term caused by the chemical reaction between silicon and silicon dioxide. For the dilatation effect within the oxide we assume a hydrostatic pressure term

$$\chi \varepsilon_{xx}^{} \varepsilon_{yy}^{} \varepsilon_{zz}^{} \Longrightarrow \Delta$$ (5.26)

The stress histories have to be calculated in order to get a correct stress-strain distribution within the different materials. Assuming the model suggested in [Pen91] stress within elastic material is calculated by

$$\sigma \Delta \sum_{i=1}^{n} \sigma_{i}^{} \Delta \sigma_{i}^{} \Delta \chi \varepsilon \Delta$$ (5.27)

and within viscous material by

$$\sigma \Delta \sum_{i=1}^{n} \sigma_{i}^{} \Delta {\frac{(n-i)\cdot\Delta T}{\tau}} \sigma_{i}^{} \Delta \varepsilon \Delta$$ (5.28)

where $\Delta$T is interpreted as the time-step length, $\tau$ as the relaxation rate, $\chi$ as the compressibility and G_eff = G ^. ${\frac{\tau}{\Delta T}}$ ^. $\left(\vphantom{1-e^{\frac{\Delta T}{\tau}}}\right.$1 - e^{${\frac{\Delta T}{\tau}}$}$\left.\vphantom{1-e^{\frac{\Delta T}{\tau}}}\right)$ as the effective modulus of rigidity.