Device Physics, TCAD, and Mathematical Modeling\n\nEvery transistor is governed by the same physics — the drift and diffusion of charge carriers through a doped crystal under electrostatic control — but no single equation is solved in practice. Device engineering is a ladder of approximations: the atomistic quantum picture is exact but unaffordable, the compact SPICE model is instant but only a calibrated fit, and the real work of technology computer-aided design (TCAD) is choosing the coarsest level that still captures the effect you care about. The map below is the spine of the whole field; everything that follows fills in one rung at a time.\n\n``svg\n\n``\n\n## 1. Physical Foundation\n\n### 1.1 Band Theory and Electronic Structure\n\n- Energy bands arise from the periodic potential of the crystal lattice — the conduction band holds empty states available for transport, the valence band holds filled states whose vacancies act as holes, and the bandgap $E_g$ separates them (Si: ~1.12 eV at 300 K).\n- Effective mass approximation — electrons and holes move as quasi-particles with a modified mass, electron $m_n^$ and hole $m_p^$, that folds the lattice potential into a single scalar.\n- Carrier statistics follow the Fermi–Dirac distribution:\n\n$$f(E) = \frac{1}{1 + \exp\left(\frac{E - E_F}{k_B T}\right)}$$\n\nIn non-degenerate semiconductors the carrier concentrations reduce to Boltzmann form:\n\n$$n = N_C \exp\left(-\frac{E_C - E_F}{k_B T}\right)$$\n\n$$p = N_V \exp\left(-\frac{E_F - E_V}{k_B T}\right)$$\n\nWhere:\n\n- $N_C$, $N_V$ = effective density of states in the conduction / valence bands\n- $E_C$, $E_V$ = conduction / valence band edges\n- $E_F$ = Fermi level\n\n### 1.2 Carrier Transport Mechanisms\n\n| Mechanism | Driving Force | Current Density |\n|-----------|---------------|-----------------|\n| Drift | Electric field $\mathbf{E}$ | $\mathbf{J} = qn\mu\mathbf{E}$ |\n| Diffusion | Concentration gradient | $\mathbf{J} = qD\nabla n$ |\n| Thermionic emission | Thermal energy over a barrier | Exponential in $\phi_B / k_B T$ |\n| Tunneling | Quantum penetration | Exponential in barrier width |\n\nThe Einstein relation ties mobility and diffusivity together, so a single measurement fixes both:\n\n$$D = \frac{k_B T}{q}\, \mu$$\n\n### 1.3 Generation and Recombination\n\nAt thermal equilibrium the mass-action law $np = n_i^2$ holds. Away from equilibrium, three mechanisms restore it: Shockley–Read–Hall (SRH) trap-assisted recombination, Auger recombination (a three-particle process that dominates at high injection), and radiative recombination (photon emission, important in direct-bandgap materials such as GaAs and InP).\n\n## 2. The Mathematical Hierarchy\n\n### 2.1 Quantum Mechanical Level (most fundamental)\n\nThe time-independent Schrödinger equation sets the states available to a confined carrier:\n\n$$\left[-\frac{\hbar^2}{2m^}\nabla^2 + V(\mathbf{r})\right]\psi = E\psi$$\n\nFor open systems — tunnel FETs, ultra-scaled MOSFETs with $L_g < 10$ nm, resonant tunneling diodes — the Non-Equilibrium Green's Function (NEGF) formalism handles contacts and coherence:\n\n$$G^R = [EI - H - \Sigma]^{-1}$$\n\nHere $H$ is the device Hamiltonian and the self-energy $\Sigma$ encodes coupling to the contacts. This is the most physically complete and the most expensive rung on the ladder.\n\n### 2.2 Boltzmann Transport Level\n\nThe Boltzmann Transport Equation (BTE) evolves the full carrier distribution in phase space and captures hot-carrier effects, velocity overshoot, and ballistic transport that the continuum models miss:\n\n$$\frac{\partial f}{\partial t} + \mathbf{v}\cdot\nabla_{\mathbf{r}} f + \frac{\mathbf{F}}{\hbar}\cdot\nabla_{\mathbf{k}} f = \left(\frac{\partial f}{\partial t}\right)_{\text{coll}}$$\n\nSolution methods: stochastic Monte Carlo particle tracking, spherical-harmonics expansion (SHE), and moment methods — the last of which is exactly what produces the drift-diffusion and hydrodynamic models below.\n\n### 2.3 Hydrodynamic / Energy-Balance Level\n\nTaking moments of the BTE with carrier energy as a variable yields an energy-balance equation whose signature feature is that the carrier temperature is allowed to decouple from the lattice, $T_n \neq T_L$:\n\n$$\frac{\partial (nw)}{\partial t} + \nabla\cdot\mathbf{S} = \mathbf{J}\cdot\mathbf{E} - \frac{n(w - w_0)}{\tau_w}$$\n\nWhere $w$ is the carrier energy density, $\mathbf{S}$ the energy flux, and $\tau_w$ the energy-relaxation time.\n\n### 2.4 Drift-Diffusion Level (the workhorse)\n\nThe overwhelming majority of production TCAD runs solve three coupled PDEs. Poisson's equation sets the electrostatics:\n\n$$\nabla\cdot(\varepsilon\nabla\psi) = -\rho = -q\,(p - n + N_D^+ - N_A^-)$$\n\nThe continuity equations conserve each carrier species:\n\n$$\frac{\partial n}{\partial t} = \frac{1}{q}\nabla\cdot\mathbf{J}_n + G_n - R_n$$\n\n$$\frac{\partial p}{\partial t} = -\frac{1}{q}\nabla\cdot\mathbf{J}_p + G_p - R_p$$\n\nAnd the current-density equations close the system, either in drift-plus-diffusion form:\n\n$$\mathbf{J}_n = q\mu_n n\,\mathbf{E} + qD_n\nabla n$$\n\n$$\mathbf{J}_p = q\mu_p p\,\mathbf{E} - qD_p\nabla p$$\n\nor, more compactly, as a gradient of the quasi-Fermi level $\mathbf{J}_n = q\mu_n n\,\nabla E_{F,n}$. The system is coupled, nonlinear, and elliptic-parabolic, and because carrier concentrations vary exponentially with potential it spans more than ten orders of magnitude across a junction — which is what makes the discretization below non-trivial.\n\n## 3. Numerical Methods\n\n### 3.1 Spatial Discretization\n\n- Finite Difference (FDM) — simple, but limited to structured rectangular grids.\n- Finite Element (FEM) — handles complex geometry through basis-function expansion and a weak variational form.\n- Finite Volume (FVM) — integrates over control volumes to guarantee local conservation, which is the natural fit for the semiconductor equations.\n\n### 3.2 Scharfetter–Gummel Discretization\n\nThe single most important trick for numerical stability: it interpolates carrier density exponentially between nodes so the current stays smooth despite huge potential swings.\n\n$$J_{n,i+\frac{1}{2}} = \frac{qD_n}{h}\left[n_i B\left(\frac{\psi_i - \psi_{i+1}}{V_T}\right) - n_{i+1} B\left(\frac{\psi_{i+1} - \psi_i}{V_T}\right)\right]$$\n\nwhere the Bernoulli function is $B(x) = x / (e^x - 1)$. It reduces to central differencing for small $\Delta\psi$ and to upwinding for large $\Delta\psi$, suppressing the spurious oscillations that a naive scheme produces. The thermal voltage $V_T = k_B T / q \approx 26$ mV at 300 K sets the scale.\n\n### 3.3 Nonlinear and Linear Solvers\n\nGummel iteration decouples the system — solve Poisson, then electron continuity, then hole continuity, and repeat to convergence. It is robust and cheap per step but converges slowly under strong coupling or high injection. Newton–Raphson solves the fully coupled linearized system $\mathbf{J}\cdot\delta\mathbf{x} = -\mathbf{F}(\mathbf{x})$ with quadratic convergence near the solution, at the cost of assembling a Jacobian and solving a larger system. In practice a hybrid strategy starts with Gummel to get close, then switches to Newton for fast final convergence. The resulting sparse, ill-conditioned Jacobians are solved with direct factorizations (PARDISO, UMFPACK) or preconditioned Krylov methods (GMRES, BiCGSTAB), with multigrid reserved for the Poisson-like blocks.\n\n## 4. Physical Models\n\n### 4.1 Mobility\n\nIndependent scattering mechanisms combine through Matthiessen's rule, $1/\mu = 1/\mu_\text{lattice} + 1/\mu_\text{impurity} + 1/\mu_\text{surface} + \cdots$. Lattice (phonon) scattering falls with temperature as $\mu_L = \mu_0 (T/300)^{-\alpha}$ ($\alpha \approx 2.4$ for Si electrons), while ionized-impurity scattering follows the Brooks–Herring model. At high field the velocity saturates via the Caughey–Thomas form:\n\n$$\mu(E) = \frac{\mu_0}{\left[1 + \left(\frac{\mu_0 E}{v_\text{sat}}\right)^\beta\right]^{1/\beta}}$$\n\nwith $v_\text{sat} \approx 10^7$ cm/s for silicon.\n\n### 4.2 Recombination\n\nShockley–Read–Hall (trap-assisted), Auger (high-density), and radiative (direct-gap) recombination each get an explicit rate:\n\n$$R_\text{SRH} = \frac{np - n_i^2}{\tau_p(n + n_1) + \tau_n(p + p_1)}$$\n\n$$R_\text{Auger} = (C_n n + C_p p)(np - n_i^2)$$\n\n$$R_\text{rad} = B(np - n_i^2)$$\n\n### 4.3 Tunneling and Quantum Corrections\n\nBand-to-band tunneling — the mechanism behind tunnel FETs and Zener breakdown — scales as $G_\text{BTBT} = A\,E^2 \exp(-B/E)$. For inversion-layer quantization in scaled MOSFETs, FinFETs, and nanowires, the density-gradient method adds a quantum potential $V_Q = -\frac{\hbar^2}{6m^}\frac{\nabla^2\sqrt{n}}{\sqrt{n}}$, while stronger confinement calls for a self-consistent 1D Schrödinger–Poisson loop that solves for subbands and iterates the quantum charge into Poisson. At high doping, bandgap narrowing $\Delta E_g = A\,N^{1/3} + B\ln(N/N_\text{ref})$ raises $n_i^2$ and feeds back into recombination.\n\n## 5. Process TCAD\n\nThe same numerical machinery models how the device is built, not just how it operates. Ion implantation is captured either by Monte Carlo trajectory tracking or by analytic Gaussian / Pearson-IV profiles. Diffusion obeys Fick's laws, $\partial C/\partial t = \nabla\cdot(D\nabla C)$, with a concentration-dependent $D$ that accounts for charged point defects. Oxidation follows the Deal–Grove relation $x_\text{ox}^2 + A\,x_\text{ox} = B(t + \tau)$, linear for thin oxides and parabolic for thick. Etch and deposition surfaces evolve by the level-set equation $\partial\phi/\partial t + v_n|\nabla\phi| = 0$, where the zero contour of $\phi$ is the moving surface.\n\n## 6. Multiphysics and Reliability\n\nReal devices are never purely electrical. Electrothermal coupling feeds Joule and recombination heating $H = \mathbf{J}\cdot\mathbf{E} + (R - G)(E_g + 3k_BT)$ into a lattice heat equation. Strain engineering shifts mobility as $\mu_\text{strained} = \mu_0(1 + \Pi\cdot\sigma)$ — the basis of strained-Si and SiGe channels. Statistical variability from random dopant fluctuations, line-edge roughness, and metal-gate granularity is swept by Monte Carlo over device instances to produce threshold-voltage distributions. And reliability models — bias-temperature instability (BTI) and hot-carrier injection (HCI) — track interface-defect generation over the device lifetime, while thermal, shot, and 1/f noise set the analog floor.\n\n## 7. Computational Architecture\n\n### 7.1 Model Hierarchy — Cost vs. Accuracy\n\n| Level | Physics captured | Governing math | Cost | Accuracy |\n|-------|------------------|----------------|------|----------|\n| NEGF | Quantum coherence | $G = [EI - H - \Sigma]^{-1}$ | Highest | Highest |\n| Monte Carlo | Full distribution function | Stochastic BTE | High | High |\n| Hydrodynamic | Carrier temperature | Hyperbolic-parabolic PDEs | Medium | Good |\n| Drift-Diffusion | Continuum transport | Elliptic-parabolic PDEs | Low | Moderate |\n| Compact | Empirical fit | Algebraic | Lowest | Calibrated |\n\n### 7.2 The TCAD ↔ Compact-Model Flow\n\nTCAD does not replace circuit simulation — it feeds it. Physics-based TCAD is calibrated against silicon measurements, then distilled into a compact model (BSIM, PSP) whose algebraic I–V equations are what SPICE actually evaluates a billion times per chip. Silicon data validates the TCAD; the compact model enables the circuit. That two-way loop — physical rigor upstream, computational speed downstream — is the reason the hierarchy at the top of this page exists at all.\n\n## 8. Reference Values\n\n| Symbol | Name | Value |\n|--------|------|-------|\n| $q$ | Elementary charge | $1.602 \times 10^{-19}$ C |\n| $k_B$ | Boltzmann constant | $1.381 \times 10^{-23}$ J/K |\n| $\hbar$ | Reduced Planck | $1.055 \times 10^{-34}$ J·s |\n| $\varepsilon_0$ | Vacuum permittivity | $8.854 \times 10^{-12}$ F/m |\n| $V_T$ | Thermal voltage (300 K) | 25.9 mV |\n\n| Silicon property (300 K) | Value |\n|--------------------------|-------|\n| Bandgap $E_g$ | 1.12 eV |\n| Intrinsic carrier density $n_i$ | $1.0 \times 10^{10}$ cm⁻³ |\n| Electron mobility $\mu_n$ | 1450 cm²/V·s |\n| Hole mobility $\mu_p$ | 500 cm²/V·s |\n| Electron saturation velocity | $1.0 \times 10^7$ cm/s |\n| Relative permittivity $\varepsilon_r$ | 11.7 |\n\nRead device physics through a quantitative lens rather than a purely qualitative one: the transistor is not a schematic symbol but a boundary-value problem, and every design decision — channel material, doping profile, gate stack, thermal budget — is ultimately a choice about which term in these equations you are willing to pay to solve exactly and which you can afford to approximate.\n