Home Knowledge Base Semiconductor Manufacturing Process: Multi-Physics Coupling & Mathematical Modeling

Semiconductor Manufacturing Process: Multi-Physics Coupling & Mathematical Modeling

Keywords: multi physics coupling, multiphysics modeling, coupled simulation, process simulation, transport phenomena, heat transfer plasma coupling, electromagnetic plasma


Semiconductor Manufacturing Process: Multi-Physics Coupling & Mathematical Modeling

1. Overview: Why Multi-Physics Coupling Matters

Semiconductor fabrication involves hundreds of process steps where multiple physical phenomena occur simultaneously and interact nonlinearly. At the 3nm node and below, these couplings become critical—small perturbations propagate across physics domains, affecting yield, uniformity, and device performance.

2. Key Processes and Their Coupled Physics

2.1 Plasma Etching (RIE, ICP, CCP)

Coupled domains:

Coupling chain:

RF Power → EM Fields → Electron Heating → Plasma Density → Sheath Voltage
    ↓                                            ↓
Ion Energy Distribution ← ─────────────────────────┘
    ↓
Surface Bombardment + Radical Flux → Etch Rate & Profile
    ↓
Feature Geometry Evolution → Local Field Modification (feedback)

2.2 Chemical Vapor Deposition (CVD/ALD)

Coupled domains:

2.3 Thermal Processing (RTP, Annealing)

Coupled domains:

2.4 EUV Lithography

Coupled domains:

3. Mathematical Framework: Governing Equations

3.1 Electromagnetics (Plasma Systems)

For RF-driven plasma, the time-harmonic Maxwell's equations:

$$abla \times \left(\mu_r^{-1} abla \times \mathbf{E}\right) - k_0^2 \epsilon_r \mathbf{E} = -j\omega\mu_0 \mathbf{J}_{ext}$$

The plasma permittivity encodes the coupling to electron density:

$$\epsilon_r = 1 - \frac{\omega_{pe}^2}{\omega(\omega + j u_m)}$$

Where the plasma frequency is:

$$\omega_{pe} = \sqrt{\frac{n_e e^2}{m_e \epsilon_0}}$$

Key parameters:

u_m$ — electron-neutral collision frequency

Note: This creates a strong nonlinear coupling: the EM field depends on plasma density, which in turn depends on power absorption from the EM field.

3.2 Plasma Transport (Drift-Diffusion Approximation)

Electron continuity equation:

$$\frac{\partial n_e}{\partial t} + abla \cdot \boldsymbol{\Gamma}_e = S_e$$

Electron flux:

$$\boldsymbol{\Gamma}_e = -\mu_e n_e \mathbf{E} - D_e abla n_e$$

Electron energy density equation:

$$\frac{\partial n_\epsilon}{\partial t} + abla \cdot \boldsymbol{\Gamma}_\epsilon + \mathbf{E} \cdot \boldsymbol{\Gamma}_e = S_\epsilon - \sum_j \varepsilon_j R_j$$

Where:

Ion transport (for multiple species $i$):

$$\frac{\partial n_i}{\partial t} + abla \cdot \boldsymbol{\Gamma}_i = S_i$$

3.3 Neutral Gas Flow (Navier-Stokes Equations)

Continuity equation:

$$\frac{\partial \rho}{\partial t} + abla \cdot (\rho \mathbf{u}) = 0$$

Momentum equation:

$$\rho \frac{D\mathbf{u}}{Dt} = - abla p + abla \cdot \boldsymbol{\tau} + \mathbf{F}_{body}$$

Where:

Low-pressure corrections (Knudsen effects):

At low pressures where Knudsen number $Kn = \lambda/L > 0.01$, slip boundary conditions are required:

$$u_{slip} = \frac{2-\sigma}{\sigma} \lambda \left.\frac{\partial u}{\partial n}\right|_{wall}$$

Where:

3.4 Species Transport and Chemistry

Convection-diffusion-reaction equation:

$$\frac{\partial c_k}{\partial t} + abla \cdot (c_k \mathbf{u}) = abla \cdot (D_k abla c_k) + R_k$$

Gas-phase reaction rates:

$$R_k = \sum_j u_{kj} \, k_j(T) \prod_l c_l^{a_{lj}}$$

Where:

u_{kj}$ — stoichiometric coefficient

Surface reactions (Langmuir-Hinshelwood kinetics):

$$r_s = k_s \theta_A \theta_B$$

Surface coverage:

$$\theta_i = \frac{K_i c_i}{1 + \sum_j K_j c_j}$$

3.5 Heat Transfer

Energy equation:

$$\rho c_p \frac{\partial T}{\partial t} + \rho c_p \mathbf{u} \cdot abla T = abla \cdot (k abla T) + Q$$

Heat sources in plasma systems:

$$Q = Q_{Joule} + Q_{ion} + Q_{reaction} + Q_{radiation}$$

Joule heating (time-averaged):

$$Q_{Joule} = \frac{1}{2} \text{Re}(\mathbf{J}^* \cdot \mathbf{E})$$

Where:

3.6 Solid Mechanics (Film Stress)

Equilibrium equation:

$$abla \cdot \boldsymbol{\sigma} = 0$$

Constitutive relation with thermal strain:

$$\boldsymbol{\sigma} = \mathbf{C} : (\boldsymbol{\epsilon} - \boldsymbol{\epsilon}_{th} - \boldsymbol{\epsilon}_{intrinsic})$$

Thermal strain tensor:

$$\boldsymbol{\epsilon}_{th} = \alpha(T - T_0)\mathbf{I}$$

Where:

Stoney equation (wafer curvature from film stress):

$$\sigma_f = \frac{E_s h_s^2}{6(1- u_s)h_f}\kappa$$

Where:

u_s$ — substrate Poisson's ratio

4. Feature-Scale Modeling

At the nanometer scale within etched features, continuum assumptions break down.

4.1 Profile Evolution (Level Set Method)

The etch front $\phi(\mathbf{x},t) = 0$ evolves according to:

$$\frac{\partial \phi}{\partial t} + V_n | abla \phi| = 0$$

Local etch rate depends on coupled physics:

$$V_n = \Gamma_{ion}(E,\theta) \cdot Y_{phys}(E,\theta) + \Gamma_{rad} \cdot Y_{chem}(T) + \Gamma_{ion} \cdot \Gamma_{rad} \cdot Y_{synergy}$$

Where:

4.2 Feature-Scale Transport

Within high-aspect-ratio features, Knudsen diffusion dominates:

$$D_{Kn} = \frac{d}{3}\sqrt{\frac{8k_BT}{\pi m}}$$

Where:

View factor calculations for flux at the bottom of features:

$$\Gamma_{bottom} = \Gamma_{top} \cdot \int_{\Omega} f(\theta) \cos\theta \, d\Omega$$

4.3 Ion Angular and Energy Distribution

At the sheath-feature interface:

$$f(E, \theta) = f_E(E) \cdot f_\theta(\theta)$$

Angular distribution (from sheath collisionality):

$$f_\theta(\theta) \propto \cos^n(\theta) \exp\left(-\frac{\theta^2}{2\sigma_\theta^2}\right)$$

Where:

5. Multi-Scale Coupling Strategy

┌─────────────────────────────────────────────────────────────┐
│                    REACTOR SCALE (cm–m)                     │
│    Continuum: Navier-Stokes, Maxwell, Drift-Diffusion       │
│    Methods: FEM, FVM                                        │
└─────────────────────┬───────────────────────────────────────┘
                      │ Boundary fluxes, plasma parameters
                      ▼
┌─────────────────────────────────────────────────────────────┐
│                    FEATURE SCALE (nm–μm)                    │
│    Kinetic transport: DSMC, Angular distribution            │
│    Profile evolution: Level set, Cell-based methods         │
└─────────────────────┬───────────────────────────────────────┘
                      │ Sticking coefficients, reaction rates
                      ▼
┌─────────────────────────────────────────────────────────────┐
│                    ATOMIC SCALE (Å–nm)                      │
│    DFT: Reaction barriers, surface energies                 │
│    MD: Sputtering yields, sticking probabilities            │
│    KMC: Surface evolution, roughness                        │
└─────────────────────────────────────────────────────────────┘

Scale hierarchy:

1. Reactor scale (cm–m)

2. Feature scale (nm–μm)

3. Atomic scale (Å–nm)

6. Coupled System Structure

The full system can be written abstractly as:

$$\mathbf{M}(\mathbf{u})\frac{\partial \mathbf{u}}{\partial t} = \mathbf{F}(\mathbf{u}, abla\mathbf{u}, abla^2\mathbf{u}, t)$$

State vector:

$$\mathbf{u} = \begin{bmatrix} n_e \\ n_\epsilon \\ n_{i,k} \\ c_j \\ T \\ \mathbf{E} \\ \mathbf{u}_{gas} \\ p \\ \boldsymbol{\sigma} \\ \phi_{profile} \\ \vdots \end{bmatrix}$$

Jacobian structure reveals coupling:

$$\mathbf{J} = \frac{\partial \mathbf{F}}{\partial \mathbf{u}} = \begin{pmatrix} J_{ee} & J_{e\epsilon} & J_{ei} & J_{ec} & \cdots \\ J_{\epsilon e} & J_{\epsilon\epsilon} & J_{\epsilon i} & & \\ J_{ie} & J_{i\epsilon} & J_{ii} & & \\ J_{ce} & & & J_{cc} & \\ \vdots & & & & \ddots \end{pmatrix}$$

Off-diagonal blocks represent inter-physics coupling strengths.

7. Numerical Solution Strategies

7.1 Coupling Approaches

Monolithic (fully coupled):

Partitioned (sequential):

Hybrid approach:

7.2 Spatial Discretization

Finite Element Method (FEM) — weak form for species transport:

$$\int_\Omega w \frac{\partial c}{\partial t} \, d\Omega + \int_\Omega w (\mathbf{u} \cdot abla c) \, d\Omega + \int_\Omega abla w \cdot (D abla c) \, d\Omega = \int_\Omega w R \, d\Omega$$

SUPG Stabilization for convection-dominated problems:

$$w \rightarrow w + \tau_{SUPG} \, \mathbf{u} \cdot abla w$$

Where:

7.3 Time Integration

Stiff systems require implicit methods:

Operator splitting for multi-physics:

$$\mathbf{u}^{n+1} = \mathcal{L}_1(\Delta t) \circ \mathcal{L}_2(\Delta t) \circ \mathcal{L}_3(\Delta t) \, \mathbf{u}^n$$

Where:

8. Specific Application: ICP Etch Model

Complete coupled system summary:

Physics DomainGoverning EquationsKey Coupling Variables
EM (inductive)

abla \times ( abla \times \mathbf{E}) + k^2\epsilon_p \mathbf{E} = 0$ | $n_e \rightarrow \epsilon_p$ |

Electron transport

abla \cdot \Gamma_e = S_e$ | $\mathbf{E}_{dc}, n_e, T_e$ |

Electron energy

abla \cdot \Gamma_\epsilon = Q_{EM} - Q_{loss}$ | $T_e \rightarrow$ rate coefficients |

Ion transport

abla \cdot \Gamma_i = S_i$ | $n_e, \mathbf{E}_{dc}$ |

Neutral chemistry

abla \cdot (c_k \mathbf{u} - D_k abla c_k) = R_k$ | $T_e \rightarrow k_{diss}$ |

Gas flowNavier-Stokes$T_{gas}$
Heat transfer

abla \cdot (k abla T) + Q = 0$ | $Q_{plasma}$ |

SheathChild-Langmuir / PIC$n_e, T_e, V_{dc}$
Feature transportKnudsen + angular$\Gamma_{ion}, \Gamma_{rad}$ from reactor
Profile evolutionLevel set$V_n$ from surface kinetics

9. EUV Lithography: Stochastic Multi-Physics

At EUV wavelength (13.5 nm), photon shot noise becomes significant.

9.1 Aerial Image Formation

$$I(\mathbf{r}) = \left|\mathcal{F}^{-1}\left[\tilde{M}(\mathbf{f}) \cdot H(\mathbf{f})\right]\right|^2$$

Where:

9.2 Photon Statistics

$$N \sim \text{Poisson}(\bar{N})$$
$$\sigma_N = \sqrt{\bar{N}}$$

Where:

9.3 Resist Exposure (Stochastic Dill Model)

$$\frac{\partial [PAG]}{\partial t} = -C \cdot I \cdot [PAG] + \xi(t)$$

Where:

9.4 Line Edge Roughness (LER)

$$\sigma_{LER} \propto \sqrt{\frac{1}{\text{dose}}} \cdot \frac{1}{\text{image contrast}}$$

Note: This requires Kinetic Monte Carlo or Gillespie algorithm rather than continuum PDEs.

10. Process Optimization (Inverse Problem)

10.1 Problem Formulation

Objective: Minimize profile deviation from target

$$\min_{\mathbf{p}} J = \int_\Gamma \left|\phi(\mathbf{x}; \mathbf{p}) - \phi_{target}\right|^2 \, d\Gamma$$

Subject to physics constraints:

$$\mathbf{F}(\mathbf{u}, \mathbf{p}) = 0$$

Control parameters $\mathbf{p}$:

10.2 Adjoint Method for Efficient Gradients

Gradient computation:

$$\frac{dJ}{d\mathbf{p}} = \frac{\partial J}{\partial \mathbf{p}} - \boldsymbol{\lambda}^T \frac{\partial \mathbf{F}}{\partial \mathbf{p}}$$

Adjoint equation:

$$\left(\frac{\partial \mathbf{F}}{\partial \mathbf{u}}\right)^T \boldsymbol{\lambda} = \left(\frac{\partial J}{\partial \mathbf{u}}\right)^T$$

Where:

11. Emerging Approaches

11.1 Physics-Informed Neural Networks (PINNs)

Loss function:

$$\mathcal{L} = \mathcal{L}_{data} + \lambda \mathcal{L}_{PDE}$$

Where:

11.2 Digital Twins

Key features:

11.3 Uncertainty Quantification

Methods:

12. Mathematical Structure

The semiconductor manufacturing multi-physics problem has a characteristic mathematical structure:

1. Hierarchy of scales (atomic → feature → reactor)

2. Nonlinear coupling between physics domains

3. Stiff ODEs/DAEs

4. Moving boundaries

5. Rarefied gas effects

6. Stochastic effects

Key Physical Constants

SymbolValueDescription
$e$$1.602 \times 10^{-19}$ CElementary charge
$m_e$$9.109 \times 10^{-31}$ kgElectron mass
$\epsilon_0$$8.854 \times 10^{-12}$ F/mPermittivity of free space
$\mu_0$$4\pi \times 10^{-7}$ H/mPermeability of free space
$k_B$$1.381 \times 10^{-23}$ J/KBoltzmann constant
$N_A$$6.022 \times 10^{23}$ mol$^{-1}$Avogadro's number

Common Dimensionless Numbers

NumberDefinitionPhysical Meaning
Knudsen ($Kn$)$\lambda / L$Mean free path / characteristic length
Reynolds ($Re$)$\rho u L / \mu$Inertia / viscous forces
Péclet ($Pe$)$u L / D$Convection / diffusion
Damköhler ($Da$)$k L / u$Reaction / convection rate
Biot ($Bi$)$h L / k$Surface / bulk heat transfer

Source: ChipFoundryServicesSearch this topicAsk CFSGPT

multi physics couplingmultiphysics modelingcoupled simulationprocess simulationtransport phenomenaheat transfer plasma couplingelectromagnetic plasma

Explore 500+ Semiconductor & AI Topics

From EUV lithography to CUDA optimization — search the full knowledge base or chat with our AI assistant.