aluminum etch,al metal etch,aluminum metal etch modeling,al etch modeling,aluminum chlorine etch,alcl3,metal etch plasma,aluminum plasma etch,bcl3 etch
# Aluminum Metal Etch Mathematical Modeling
1. Overview
1.1 Why Aluminum Etch Modeling is Complex
Aluminum etching (typically using $\text{Cl}_2/\text{BCl}_3$ plasmas) involves multiple coupled physical and chemical phenomena:
- Plasma generation and transport → determines species fluxes to wafer
- Ion-surface interactions → physical and chemical mechanisms
- Surface reactions → Langmuir-Hinshelwood kinetics
- Feature-scale evolution → profile development inside trenches/vias
- Redeposition and passivation → sidewall chemistry
1.2 Fundamental Reaction
The basic aluminum chlorination reaction:
$$
\text{Al} + 3\text{Cl} \rightarrow \text{AlCl}_3 \uparrow
$$
Complications requiring sophisticated modeling:
- Breaking through native $\text{Al}_2\text{O}_3$ layer (15-30 Å)
- Maintaining profile anisotropy
- Controlling selectivity to mask and underlayers
- Managing Cu residues in Al-Cu alloys
2. Kinetic and Chemical Rate Modeling
2.1 General Etch Rate Formulation
A comprehensive etch rate model combines three primary mechanisms:
$$
ER = \underbrace{k_{th} \cdot \Gamma_{Cl} \cdot f(\theta)}_{\text{thermal chemical}} + \underbrace{Y_s \cdot \Gamma_{ion} \cdot \sqrt{E_{ion}}}_{\text{physical sputtering}} + \underbrace{\beta \cdot \Gamma_{ion}^a \cdot \Gamma_{Cl}^b \cdot E_{ion}^c}_{\text{ion-enhanced (synergistic)}}
$$
Parameter Definitions:
| Symbol | Description | Units |
|--------|-------------|-------|
| $\Gamma_{Cl}$ | Neutral chlorine flux | $\text{cm}^{-2}\text{s}^{-1}$ |
| $\Gamma_{ion}$ | Ion flux | $\text{cm}^{-2}\text{s}^{-1}$ |
| $E_{ion}$ | Ion energy | eV |
| $\theta$ | Surface coverage of reactive species | dimensionless |
| $Y_s$ | Physical sputtering yield | atoms/ion |
| $\beta$ | Synergy coefficient | varies |
| $a, b, c$ | Exponents (typically 0.5-1) | dimensionless |
2.2 Surface Coverage Dynamics
The reactive site balance follows Langmuir-Hinshelwood kinetics:
$$
\frac{d\theta}{dt} = k_{ads} \cdot \Gamma_{Cl} \cdot (1-\theta) - k_{des} \cdot \theta \cdot \exp\left(-\frac{E_d}{k_B T}\right) - Y_{react}(\theta, E_{ion}) \cdot \Gamma_{ion} \cdot \theta
$$
Term-by-term breakdown:
- Term 1: $k_{ads} \cdot \Gamma_{Cl} \cdot (1-\theta)$ — Adsorption rate (proportional to empty sites)
- Term 2: $k_{des} \cdot \theta \cdot \exp(-E_d/k_B T)$ — Thermal desorption (Arrhenius)
- Term 3: $Y_{react} \cdot \Gamma_{ion} \cdot \theta$ — Ion-induced reaction/removal
Steady-State Solution ($d\theta/dt = 0$):
$$
\theta_{ss} = \frac{k_{ads} \cdot \Gamma_{Cl}}{k_{ads} \cdot \Gamma_{Cl} + k_{des} \cdot e^{-E_d/k_B T} + Y_{react} \cdot \Gamma_{ion}}
$$
2.3 Temperature Dependence
All rate constants follow Arrhenius behavior:
$$
k_i(T) = A_i \cdot \exp\left(-\frac{E_{a,i}}{k_B T}\right)
$$
Typical activation energies for aluminum etching:
- Ion-enhanced reactions: $E_a \approx 0.1 - 0.3 \text{ eV}$
- Purely thermal processes: $E_a \approx 0.5 - 1.0 \text{ eV}$
- Chlorine desorption: $E_d \approx 0.3 - 0.5 \text{ eV}$
2.4 Complete Etch Rate Expression
Combining all terms with explicit dependencies:
$$
ER(T, \Gamma_{ion}, \Gamma_{Cl}, E_{ion}) = A_1 e^{-E_1/k_B T} \Gamma_{Cl} \theta + Y_0 \Gamma_{ion} \sqrt{E_{ion}} + A_2 e^{-E_2/k_B T} \Gamma_{ion}^{0.5} \Gamma_{Cl}^{0.5} E_{ion}^{0.5}
$$
3. Ion-Surface Interaction Physics
3.1 Ion Energy Distribution Function (IEDF)
For RF-biased electrodes, the IEDF is approximately bimodal:
$$
f(E) \propto \frac{1}{\sqrt{|E - E_{dc}|}} \quad \text{for } E_{dc} - E_{rf} < E < E_{dc} + E_{rf}
$$
Key parameters:
- $E_{dc} = e \cdot V_{dc}$ — DC self-bias energy
- $E_{rf} = e \cdot V_{rf}$ — RF amplitude energy
- Peak separation: $\Delta E = 2 E_{rf}$
Collisional effects:
In collisional sheaths, charge-exchange collisions broaden the distribution:
$$
f(E) \propto \exp\left(-\frac{E}{\bar{E}}\right) \cdot \left[1 + \text{erf}\left(\frac{E - E_{dc}}{\sigma_E}\right)\right]
$$
3.2 Ion Angular Distribution Function (IADF)
The angular spread is approximately Gaussian:
$$
f(\theta) = \frac{1}{\sqrt{2\pi}\sigma_\theta} \exp\left(-\frac{\theta^2}{2\sigma_\theta^2}\right)
$$
Angular spread calculation:
$$
\sigma_\theta \approx \sqrt{\frac{k_B T_i}{e V_{sheath}}} \approx \arctan\left(\sqrt{\frac{T_i}{V_{sheath}}}\right)
$$
Typical values:
- Ion temperature: $T_i \approx 0.05 - 0.5 \text{ eV}$
- Sheath voltage: $V_{sheath} \approx 50 - 500 \text{ V}$
- Angular spread: $\sigma_\theta \approx 2° - 5°$
3.3 Physical Sputtering Yield
Yamamura Formula (Angular Dependence)
$$
Y(\theta) = Y(0°) \cdot \cos^{-f}(\theta) \cdot \exp\left[b\left(1 - \frac{1}{\cos\theta}\right)\right]
$$
Parameters for aluminum:
- $f \approx 1.5 - 2.0$
- $b \approx 0.1 - 0.3$ (depends on ion/target mass ratio)
- Maximum yield typically at $\theta \approx 60° - 70°$
Sigmund Theory (Energy Dependence)
$$
Y(E) = \frac{0.042 \cdot Q \cdot \alpha(M_2/M_1) \cdot S_n(E)}{U_s}
$$
Where:
- $S_n(E)$ = nuclear stopping power (Thomas-Fermi)
- $U_s = 3.4 \text{ eV}$ (surface binding energy for Al)
- $Q$ = dimensionless factor ($\approx 1$ for metals)
- $\alpha$ = mass-dependent parameter
- $M_1, M_2$ = projectile and target masses
Nuclear Stopping Power
$$
S_n(\epsilon) = \frac{0.5 \ln(1 + 1.2288\epsilon)}{\epsilon + 0.1728\sqrt{\epsilon} + 0.008\epsilon^{0.1504}}
$$
With reduced energy:
$$
\epsilon = \frac{M_2 E}{(M_1 + M_2) Z_1 Z_2 e^2} \cdot \frac{a_{TF}}{1}
$$
3.4 Ion-Enhanced Etching Yield
The total etch yield combines mechanisms:
$$
Y_{total} = Y_{physical} + Y_{chemical} + Y_{synergistic}
$$
Synergistic enhancement factor:
$$
\eta = \frac{Y_{total}}{Y_{physical} + Y_{chemical}} > 1
$$
For Al/Cl₂ systems, $\eta$ can exceed 10 under optimal conditions.
4. Plasma Modeling (Reactor Scale)
4.1 Species Continuity Equations
For each species $i$ (electrons, ions, neutrals):
$$
\frac{\partial n_i}{\partial t} + \nabla \cdot \vec{\Gamma}_i = S_i - L_i
$$
Flux expressions:
- Drift-diffusion: $\vec{\Gamma}_i = -D_i \nabla n_i + \mu_i n_i \vec{E}$
- Full momentum: $\vec{\Gamma}_i = n_i \vec{v}_i$ with momentum equation
Source/sink terms:
$$
S_i = \sum_j k_{ij} n_j n_e \quad \text{(ionization, dissociation)}
$$
$$
L_i = \sum_j k_{ij}^{loss} n_i n_j \quad \text{(recombination, attachment)}
$$
4.2 Electron Energy Balance
$$
\frac{\partial}{\partial t}\left(\frac{3}{2} n_e k_B T_e\right) + \nabla \cdot \vec{Q}_e = P_{abs} - P_{loss}
$$
Heat flux:
$$
\vec{Q}_e = \frac{5}{2} k_B T_e \vec{\Gamma}_e - \kappa_e \nabla T_e
$$
Power absorption (ICP):
$$
P_{abs} = \frac{1}{2} \text{Re}(\sigma_p) |E|^2
$$
Collisional losses:
$$
P_{loss} = \sum_j n_e n_j k_j \varepsilon_j
$$
Where $\varepsilon_j$ is the energy loss per collision event $j$.
4.3 Plasma Conductivity
$$
\sigma_p = \frac{n_e e^2}{m_e(\nu_m + i\omega)}
$$
Skin depth:
$$
\delta = \sqrt{\frac{2}{\omega \mu_0 \text{Re}(\sigma_p)}}
$$
4.4 Electromagnetic Field Equations
Maxwell's equations (frequency domain):
$$
\nabla \times \vec{E} = -i\omega \vec{B}
$$
$$
\nabla \times \vec{B} = \mu_0 \sigma_p \vec{E} + i\omega \mu_0 \epsilon_0 \vec{E}
$$
Wave equation:
$$
\nabla^2 \vec{E} + \left(\frac{\omega^2}{c^2} - i\omega\mu_0\sigma_p\right)\vec{E} = 0
$$
4.5 Sheath Physics
Child-Langmuir Law (Collisionless Sheath)
$$
J_{ion} = \frac{4\epsilon_0}{9}\sqrt{\frac{2e}{M}} \cdot \frac{V_s^{3/2}}{s^2}
$$
Where:
- $J_{ion}$ = ion current density
- $V_s$ = sheath voltage
- $s$ = sheath thickness
- $M$ = ion mass
Bohm Criterion
Ions must enter sheath with velocity:
$$
v_{Bohm} = \sqrt{\frac{k_B T_e}{M}}
$$
Ion flux at sheath edge:
$$
\Gamma_{ion} = n_s \cdot v_{Bohm} = 0.61 \cdot n_0 \sqrt{\frac{k_B T_e}{M}}
$$
Sheath Thickness
$$
s \approx \lambda_D \cdot \left(\frac{2 e V_s}{k_B T_e}\right)^{3/4}
$$
Debye length:
$$
\lambda_D = \sqrt{\frac{\epsilon_0 k_B T_e}{n_e e^2}}
$$
5. Feature-Scale Profile Evolution
5.1 Level Set Method
The surface is represented implicitly by $\phi(\vec{r}, t) = 0$:
$$
\frac{\partial \phi}{\partial t} + V_n |\nabla \phi| = 0
$$
Normal velocity calculation:
$$
V_n(\vec{r}) = \int_0^{E_{max}} \int_0^{\theta_{max}} Y(E, \theta_{local}) \cdot f_{IEDF}(E) \cdot f_{IADF}(\theta) \cdot \Gamma_{ion}(\vec{r}) \, dE \, d\theta
$$
Plus contributions from:
- Neutral chemical etching
- Redeposition
- Surface diffusion
5.2 Hamilton-Jacobi Formulation
$$
\frac{\partial \phi}{\partial t} + H(\nabla \phi, \vec{r}, t) = 0
$$
Hamiltonian for etch:
$$
H = V_n \sqrt{\phi_x^2 + \phi_y^2 + \phi_z^2}
$$
With $V_n$ dependent on:
- Local surface normal: $\hat{n} = -\nabla\phi / |\nabla\phi|$
- Local fluxes: $\Gamma(\vec{r})$
- Local angles: $\theta = \arccos(\hat{n} \cdot \hat{z})$
5.3 Visibility and View Factors
Direct Flux
The flux reaching a point inside a feature depends on solid angle visibility:
$$
\Gamma_{direct}(\vec{r}) = \int_{\Omega_{visible}} \Gamma_0 \cdot \cos\theta \cdot \frac{d\Omega}{\pi}
$$
Reflected/Reemitted Flux
For neutrals with sticking coefficient $s$:
$$
\Gamma_{total}(\vec{r}) = \Gamma_{direct}(\vec{r}) + (1-s) \cdot \Gamma_{reflected}(\vec{r})
$$
This leads to coupled integral equations:
$$
\Gamma(\vec{r}) = \Gamma_{plasma}(\vec{r}) + (1-s) \int_{S'} K(\vec{r}, \vec{r'}) \Gamma(\vec{r'}) dS'
$$
Kernel function:
$$
K(\vec{r}, \vec{r'}) = \frac{\cos\theta \cos\theta'}{\pi |\vec{r} - \vec{r'}|^2} \cdot V(\vec{r}, \vec{r'})
$$
Where $V(\vec{r}, \vec{r'})$ is the visibility function (1 if visible, 0 otherwise).
5.4 Aspect Ratio Dependent Etching (ARDE)
Empirical model:
$$
\frac{ER(AR)}{ER_0} = \frac{1}{1 + (AR/AR_c)^n}
$$
Where:
- $AR = \text{depth}/\text{width}$ (aspect ratio)
- $AR_c$ = critical aspect ratio (process-dependent)
- $n \approx 1 - 2$
Knudsen transport model:
$$
\Gamma_{neutral}(z) = \Gamma_0 \cdot \frac{W}{W + \alpha \cdot z}
$$
Where:
- $z$ = feature depth
- $W$ = feature width
- $\alpha$ = Clausing factor (depends on geometry and sticking)
Clausing factor for cylinder:
$$
\alpha = \frac{8}{3} \cdot \frac{1 - s}{s}
$$
6. Aluminum-Specific Phenomena
6.1 Native Oxide Breakthrough
$\text{Al}_2\text{O}_3$ (15-30 Å native oxide) requires physical sputtering:
$$
ER_{oxide} \approx Y_{\text{BCl}_3^+}(E) \cdot \Gamma_{ion}
$$
Why BCl₃ is critical:
1. Heavy $\text{BCl}_3^+$ ions provide efficient momentum transfer
2. BCl₃ scavenges oxygen chemically:
$$
2\text{BCl}_3 + \text{Al}_2\text{O}_3 \rightarrow 2\text{AlCl}_3 \uparrow + \text{B}_2\text{O}_3
$$
Breakthrough time:
$$
t_{breakthrough} = \frac{d_{oxide}}{ER_{oxide}} = \frac{d_{oxide}}{Y_{BCl_3^+} \cdot \Gamma_{ion}}
$$
6.2 Sidewall Passivation Dynamics
Anisotropic profiles require passivation of sidewalls:
$$
\frac{d\tau_{pass}}{dt} = R_{dep}(\Gamma_{redeposition}, s_{stick}) - R_{removal}(\Gamma_{ion}, \theta_{sidewall})
$$
Deposition sources:
- $\text{AlCl}_x$ redeposition from etch products
- Photoresist erosion products (C, H, O, N)
- Intentional additives: $\text{N}_2 \rightarrow \text{AlN}$ formation
Why sidewalls are protected:
At grazing incidence ($\theta \approx 85° - 90°$):
- Ion flux geometric factor: $\Gamma_{sidewall} = \Gamma_0 \cdot \cos(90° - \alpha) \approx \Gamma_0 \cdot \sin\alpha$
- For $\alpha = 5°$: $\Gamma_{sidewall} \approx 0.09 \cdot \Gamma_0$
- Sputtering yield at grazing incidence approaches zero
- Net passivation accumulates → blocks lateral etching
6.3 Notching and Charging Effects
At dielectric interfaces, differential charging causes ion deflection:
Surface charge evolution:
$$
\frac{d\sigma}{dt} = J_{ion} - J_{electron}
$$
Where:
- $\sigma$ = surface charge density (C/cm²)
- $J_{ion}$ = ion current (always positive)
- $J_{electron}$ = electron current (depends on local potential)
Local electric field:
$$
\vec{E}_{charging} = -\nabla V_{charging}
$$
Laplace equation in feature:
$$
\nabla^2 V = -\frac{\rho}{\epsilon_0} \quad \text{(with } \rho = 0 \text{ in vacuum)}
$$
Modified ion trajectory:
$$
m \frac{d^2\vec{r}}{dt^2} = e\left(\vec{E}_{sheath} + \vec{E}_{charging}\right)
$$
Result: Ions deflect toward charged surfaces → notching at feature bottom.
Mitigation strategies:
- Pulsed plasmas (allow electron neutralization)
- Low-frequency bias (time for charge equilibration)
- Conductive underlayers
6.4 Copper Residue Formation (Al-Cu Alloys)
Al-Cu alloys (0.5-4% Cu) leave Cu residues because Cu chlorides are less volatile:
Volatility comparison:
| Species | Sublimation/Boiling Point |
|---------|---------------------------|
| $\text{AlCl}_3$ | 180°C (sublimes) |
| $\text{CuCl}$ | 430°C (sublimes) |
| $\text{CuCl}_2$ | 300°C (decomposes) |
Residue accumulation rate:
$$
\frac{d[\text{Cu}]_{surface}}{dt} = x_{Cu} \cdot ER_{Al} - ER_{Cu}
$$
Where:
- $x_{Cu}$ = Cu atomic fraction in alloy
- At low temperature: $ER_{Cu} \ll x_{Cu} \cdot ER_{Al}$
Solutions:
- Elevated substrate temperature ($>$150°C)
- Increased BCl₃ fraction
- Post-etch treatments
7. Numerical Methods
7.1 Level Set Discretization
Upwind Finite Differences
Using Hamilton-Jacobi ENO (Essentially Non-Oscillatory) schemes:
$$
\phi_i^{n+1} = \phi_i^n - \Delta t \cdot H(\phi_x^-, \phi_x^+, \phi_y^-, \phi_y^+)
$$
One-sided derivatives:
$$
\phi_x^- = \frac{\phi_i - \phi_{i-1}}{\Delta x}, \quad \phi_x^+ = \frac{\phi_{i+1} - \phi_i}{\Delta x}
$$
Godunov flux for $H = V_n |\nabla\phi|$:
$$
H^{Godunov} =
\begin{cases}
V_n \sqrt{\max(\phi_x^{-,+},0)^2 + \max(\phi_y^{-,+},0)^2} & \text{if } V_n > 0 \\
V_n \sqrt{\max(\phi_x^{+,-},0)^2 + \max(\phi_y^{+,-},0)^2} & \text{if } V_n < 0
\end{cases}
$$
Reinitialization
Maintain $|\nabla\phi| = 1$ using:
$$
\frac{\partial \phi}{\partial \tau} = \text{sign}(\phi_0)(1 - |\nabla\phi|)
$$
Iterate in pseudo-time $\tau$ until convergence.
7.2 Monte Carlo Feature-Scale Simulation
Algorithm:
1. INITIALIZE surface mesh
2. FOR each time step:
a. FOR i = 1 to N_particles:
- Sample particle from IEDF, IADF
- Launch from plasma boundary
- TRACE trajectory until surface hit
- APPLY reaction probability:
* Etch (remove cell) with probability P_etch
* Reflect with probability P_reflect
* Deposit with probability P_deposit
b. UPDATE surface mesh
c. CHECK for convergence
3. OUTPUT final profile
Variance reduction techniques:
- Importance sampling: Weight particles toward features of interest
- Particle splitting: Increase statistics in critical regions
- Russian roulette: Terminate low-weight particles probabilistically
7.3 Coupled Multi-Scale Modeling
| Scale | Domain | Method | Outputs |
|-------|--------|--------|---------|
| Reactor | m | Fluid/hybrid plasma | $n_e$, $T_e$, species densities |
| Sheath | mm | PIC or fluid | IEDF, IADF, fluxes |
| Feature | nm-μm | Level set / Monte Carlo | Profile evolution |
| Atomistic | Å | MD / DFT | Yields, sticking coefficients |
Coupling strategy:
$$
\text{Reactor} \xrightarrow{\Gamma_i, f(E), f(\theta)} \text{Feature} \xrightarrow{ER(\vec{r})} \text{Reactor}
$$
7.4 Plasma Solver Discretization
Finite element for Poisson's equation:
$$
\nabla \cdot (\epsilon \nabla V) = -\rho
$$
Weak form:
$$
\int_\Omega \epsilon \nabla V \cdot \nabla w \, d\Omega = \int_\Omega \rho \, w \, d\Omega
$$
Finite volume for transport:
$$
\frac{d(n_i V_j)}{dt} = -\sum_{faces} \Gamma_i \cdot \hat{n} \cdot A + S_i V_j
$$
8. Process Window and Optimization
8.1 Response Surface Modeling
Quadratic response surface:
$$
ER = \beta_0 + \sum_{i=1}^{k} \beta_i x_i + \sum_{i=1}^{k} \beta_{ii} x_i^2 + \sum_{i