Drift–diffusion transport is the continuum model that converts electric fields, carrier-density gradients, generation, recombination, and contact injection into semiconductor current and charge evolution. Its useful object is not one current formula but a coupled system: Poisson electrostatics determines the field, electron and hole continuity equations conserve particles, constitutive flux laws connect current to electrochemical-potential gradients, and material plus boundary models close the problem. It is the workhorse of device simulation when carriers remain near local equilibrium and transport lengths are long enough for mobility and diffusion descriptions to be meaningful.

<svg viewBox="0 0 760 470" xmlns="http://www.w3.org/2000/svg" font-family="Arial,sans-serif"><rect x="0" y="0" width="760" height="470" rx="18" fill="#0d1117"/><text x="380" y="34" fill="#e6edf3" font-size="21" font-weight="700" text-anchor="middle">Drift–diffusion is a closed transport system</text><text x="380" y="60" fill="#8b949e" font-size="12" text-anchor="middle">Field, charge, flux, and continuity must agree simultaneously</text><rect x="55" y="115" width="175" height="205" rx="14" fill="#161b22" stroke="#58a6ff" stroke-width="2"/><text x="142" y="150" fill="#79c0ff" font-size="15" font-weight="700" text-anchor="middle">Poisson</text><text x="142" y="200" fill="#e6edf3" font-size="13" text-anchor="middle">∇·(ε∇φ) = −ρ</text><text x="142" y="260" fill="#c9d1d9" font-size="11" text-anchor="middle">potential and field</text><path d="M230 215H285" stroke="#3fb950" stroke-width="4"/><polygon points="285,215 272,207 272,223" fill="#3fb950"/><rect x="285" y="90" width="190" height="255" rx="14" fill="#161b22" stroke="#3fb950" stroke-width="2"/><text x="380" y="130" fill="#7ee787" font-size="15" font-weight="700" text-anchor="middle">Flux laws</text><text x="380" y="185" fill="#e6edf3" font-size="12" text-anchor="middle">Jₙ = qμₙnE + qDₙ∇n</text><text x="380" y="225" fill="#e6edf3" font-size="12" text-anchor="middle">Jₚ = qμₚpE − qDₚ∇p</text><text x="380" y="290" fill="#c9d1d9" font-size="11" text-anchor="middle">mobility + diffusion + statistics</text><path d="M475 215H530" stroke="#d29922" stroke-width="4"/><polygon points="530,215 517,207 517,223" fill="#d29922"/><rect x="530" y="115" width="175" height="205" rx="14" fill="#161b22" stroke="#d29922" stroke-width="2"/><text x="617" y="150" fill="#e3b341" font-size="15" font-weight="700" text-anchor="middle">Continuity</text><text x="617" y="205" fill="#e6edf3" font-size="12" text-anchor="middle">∂n/∂t + flux div.</text><text x="617" y="245" fill="#e6edf3" font-size="12" text-anchor="middle">= G − R</text><text x="617" y="285" fill="#c9d1d9" font-size="11" text-anchor="middle">carrier conservation</text><path d="M617 320V385H142V320" fill="none" stroke="#a371f7" stroke-width="4"/><polygon points="142,320 134,333 150,333" fill="#a371f7"/><text x="380" y="415" fill="#d2a8ff" font-size="12" text-anchor="middle">n, p set charge; φ sets field; all equations close together</text></svg>

The model separates conservation laws from constitutive assumptions. Electron and hole continuity equations are particle balances and remain valid far beyond simple drift–diffusion. The current relations are closures derived by reducing the Boltzmann transport equation under assumptions about scattering, local distributions, and moments. Poisson's equation is the electrostatic closure. Keeping these roles distinct makes model extensions intelligible: changing mobility alters the constitutive law, adding traps alters source terms and charge, and adopting hydrodynamic transport adds higher moments without replacing carrier conservation.

Electric drift and concentration diffusion are two representations of electrochemical driving. For an isothermal, nondegenerate semiconductor under a conventional sign choice, $\mathbf J_n=q\mu_nn\mathbf E+qD_n\nabla n$ and $\mathbf J_p=q\mu_pp\mathbf E-qD_p\nabla p$, with $\mathbf E=-\nabla\phi$. Electron conventional current points opposite electron particle motion, which explains signs that otherwise look asymmetric. The safest practice is to derive each current from charge times particle flux and verify equilibrium cancellation rather than memorize isolated signs across software packages.

The Einstein relation links mobility and diffusion only within a statistical regime. In the Maxwell–Boltzmann, isothermal limit, $D_n/\mu_n=D_p/\mu_p=k_BT/q=V_T$. It guarantees that diffusion opposing an equilibrium density gradient cancels electric drift. Degenerate carriers require a generalized Einstein relation involving derivatives of density with respect to chemical potential. Hot carriers, nonlocal transport, magnetic fields, and anisotropic bands can require tensor or energy-dependent coefficients. Assigning mobility and diffusion independently can violate detailed balance and create spurious equilibrium current.

Quasi-Fermi levels provide the most stable physical interpretation of current. Under suitable conventions, electron and hole currents are proportional to carrier density, mobility, and gradients of their respective quasi-Fermi energies or electrochemical potentials. At thermal equilibrium the two quasi-Fermi levels collapse to one constant Fermi level, so current vanishes even though electric-field and concentration-gradient terms may each be large. Under bias their splitting measures nonequilibrium. MIT device-physics notes emphasize this gradient form because it unifies drift and diffusion and makes contact boundary conditions clearer.

Carrier continuity turns local imbalance into storage or flux divergence. A consistent convention gives $\partial_t n-(1/q)\nabla\cdot\mathbf J_n=G_n-R_n$ and $\partial_t p+(1/q)\nabla\cdot\mathbf J_p=G_p-R_p$. Integrating over a control volume relates carrier-number change to terminal flux plus net generation. In steady state the time derivative vanishes, but current need not be spatially constant if generation, recombination, or exchange between carrier populations occurs. Total conventional current can remain conserved when electron and hole components trade through pair recombination.

Poisson coupling makes transport nonlinear even when every isolated equation looks familiar. Charge density $\rho=q(p-n+N_D^+-N_A^-)+\rho_{fixed}+\rho_{trap}$ bends the bands, the field changes drift, potential changes carrier statistics, and carrier changes feed back into charge. Mobility and recombination may also depend on field, temperature, and density. A low residual for one subequation does not establish a self-consistent solution. Potential, both carrier equations, trap occupation, contact currents, and any thermal equation must satisfy one declared convergence standard.

A sign-and-unit ledger prevents the most expensive class of implementation errors. Record whether $q$ means positive elementary charge, whether currents are conventional or particle fluxes, whether quasi-Fermi variables use volts or electron-volts, and whether recombination is positive for carrier loss. Current density has amperes per square meter, number flux has inverse square-meter seconds, $G$ and $R$ have inverse cubic-meter seconds, and doping is a number density unless multiplied by $q$. A solver can converge smoothly with a missing charge factor, so dimensional tests are part of verification.

<svg viewBox="0 0 760 470" xmlns="http://www.w3.org/2000/svg" font-family="Arial,sans-serif"><rect x="0" y="0" width="760" height="470" rx="18" fill="#0d1117"/><text x="380" y="34" fill="#e6edf3" font-size="21" font-weight="700" text-anchor="middle">Equilibrium is exact cancellation, not absence of driving terms</text><text x="380" y="60" fill="#8b949e" font-size="12" text-anchor="middle">A constant electrochemical potential makes net current zero</text><rect x="55" y="110" width="650" height="245" rx="16" fill="#161b22" stroke="#30363d" stroke-width="2"/><path d="M95 295C210 275 280 150 380 155C490 160 555 280 665 290" fill="none" stroke="#58a6ff" stroke-width="4"/><text x="135" y="275" fill="#79c0ff" font-size="12">band edge E꜀(x)</text><line x1="95" y1="220" x2="665" y2="220" stroke="#d2a8ff" stroke-width="4"/><text x="580" y="207" fill="#d2a8ff" font-size="12">constant Eꜰ</text><path d="M120 380H315" stroke="#f85149" stroke-width="8"/><polygon points="315,380 296,368 296,392" fill="#f85149"/><text x="215" y="415" fill="#ff7b72" font-size="12" text-anchor="middle">electric drift</text><path d="M640 380H445" stroke="#3fb950" stroke-width="8"/><polygon points="445,380 464,368 464,392" fill="#3fb950"/><text x="545" y="415" fill="#7ee787" font-size="12" text-anchor="middle">concentration diffusion</text><text x="380" y="450" fill="#e6edf3" font-size="13" text-anchor="middle">Jₙ,drift + Jₙ,diffusion = 0 when ∇Eꜰₙ = 0</text></svg>

Thermal equilibrium is a stringent zero-current benchmark. With no illumination, imposed current, or time-varying drive, solve Poisson and carrier statistics so both quasi-Fermi levels are spatial constants. Drift and diffusion currents should cancel to discretization tolerance at every face. A nonzero equilibrium current commonly reveals inconsistent Einstein relations, band-edge interpolation, contact statistics, or flux discretization. This test is stronger than observing equal terminal currents because local errors may cancel globally.

Local equilibrium is the central closure assumption. Drift–diffusion treats each carrier population as sufficiently relaxed that a density, temperature, and quasi-Fermi level characterize the relevant distribution locally. Momentum relaxation is assumed fast compared with spatial and temporal variation of these macroscopic fields. The model can remain useful far from global equilibrium while failing when a local equilibrium distribution is not established. Ballistic channels, sharp energy filtering, velocity overshoot, and strongly nonthermal injection expose that limit.

Mobility is a model of momentum loss rather than a universal material constant. Low-field mobility depends on phonons, ionized impurities, neutral defects, alloy disorder, interfaces, carrier density, and temperature. Device models often combine scattering mechanisms through Matthiessen-like rules, but simple inverse-rate addition is approximate when mechanisms interact. Calibration must specify crystal orientation, stress, doping, temperature, and extraction method. A mobility fitted to one transistor geometry can absorb contact resistance or quantum confinement and fail when transferred elsewhere.

High-field transport requires velocity saturation or a higher-order model. The low-field relation $v_d=\mu E$ cannot grow without bound. Empirical field-dependent mobility or velocity-saturation laws limit drift speed and reproduce long-channel current trends. Nonlocal velocity overshoot depends on carrier energy history and cannot be captured by a purely local field law. Hydrodynamic or energy-transport models add carrier temperature or energy flux, while Monte Carlo and Boltzmann solvers resolve distribution dynamics at greater cost.

Surface mobility must distinguish normal confinement from lateral driving field. In MOS inversion layers, vertical effective field presses carriers toward an interface and changes roughness and phonon scattering, while lateral field drives channel current. Compact mobility formulas combine doping, effective field, temperature, and velocity saturation. Using total field magnitude can conflate these roles around corners. Density-gradient or quantum corrections shift the carrier centroid, which changes the effective field and therefore mobility as well as electrostatics.

Degenerate statistics change both density and transport response. Maxwell–Boltzmann approximations are accurate when quasi-Fermi levels lie several $k_BT$ from relevant band edges. Heavy doping, strong accumulation, and low temperature require Fermi–Dirac integrals. COMSOL's current semiconductor documentation explicitly distinguishes these regimes. Degeneracy modifies the generalized Einstein factor, incomplete ionization, screening, and thermoelectric response. Switching statistics without recalibrating mobility and bandgap narrowing can double-count or omit density effects.

Bandgap narrowing and incomplete ionization alter more than equilibrium charge. Heavy doping shifts effective band edges and intrinsic concentration, influencing pn-junction built-in potential, recombination, and quasi-Fermi relations. Donors and acceptors may not be fully ionized, especially at cryogenic temperature, and their occupation can depend on local potential. Abruptly changing these models with a threshold creates derivative discontinuities that hurt Newton convergence. Parameterizations must match the selected statistics and material composition.

Heterostructure transport needs thermodynamic and interface consistency. Electron affinity, bandgap, density of states, permittivity, mobility, and recombination parameters may jump across a material interface. Thermionic emission, tunneling, or interface resistance can replace simple continuity of carrier quasi-Fermi level. The normal total current must still balance with interface storage or recombination. Naively smoothing band offsets changes barrier transmission; imposing both carrier density and flux can overconstrain the interface. A heterojunction model must state which quantities are continuous and why.

<svg viewBox="0 0 760 470" xmlns="http://www.w3.org/2000/svg" font-family="Arial,sans-serif"><rect x="0" y="0" width="760" height="470" rx="18" fill="#0d1117"/><text x="380" y="34" fill="#e6edf3" font-size="21" font-weight="700" text-anchor="middle">Generation and recombination reshape carrier-current paths</text><text x="380" y="60" fill="#8b949e" font-size="12" text-anchor="middle">Continuity converts local pair exchange into current divergence</text><rect x="70" y="105" width="620" height="270" rx="16" fill="#161b22" stroke="#30363d" stroke-width="2"/><circle cx="220" cy="225" r="24" fill="#58a6ff"/><circle cx="540" cy="225" r="24" fill="#f85149"/><text x="220" y="231" fill="#0d1117" font-size="20" font-weight="700" text-anchor="middle">−</text><text x="540" y="231" fill="#0d1117" font-size="20" font-weight="700" text-anchor="middle">+</text><path d="M315 125V195" stroke="#d29922" stroke-width="5"/><polygon points="315,195 306,179 324,179" fill="#d29922"/><path d="M445 195V125" stroke="#d29922" stroke-width="5"/><polygon points="445,125 436,141 454,141" fill="#d29922"/><text x="380" y="115" fill="#e3b341" font-size="13" text-anchor="middle">optical / thermal generation G</text><path d="M245 225H515" stroke="#a371f7" stroke-width="4" stroke-dasharray="8 5"/><text x="380" y="212" fill="#d2a8ff" font-size="12" text-anchor="middle">recombination R</text><path d="M220 255V330" stroke="#58a6ff" stroke-width="4"/><path d="M540 255V330" stroke="#f85149" stroke-width="4"/><text x="220" y="352" fill="#79c0ff" font-size="12" text-anchor="middle">∇·Jₙ = −q(G−R)</text><text x="540" y="352" fill="#ff7b72" font-size="12" text-anchor="middle">∇·Jₚ = +q(G−R)</text><text x="380" y="425" fill="#c9d1d9" font-size="12" text-anchor="middle">Electron and hole currents may vary while total terminal current remains conserved.</text></svg>

Shockley–Read–Hall recombination represents trap-assisted exchange. A common rate is $R_{SRH}=(np-n_i^2)/[\tau_p(n+n_1)+\tau_n(p+p_1)]$, with trap energy embedded in $n_1$ and $p_1$. Lifetimes are effective parameters tied to defect type, density, capture cross sections, and temperature, not immutable bulk constants. Interface traps require surface or distributed boundary treatment. The sign should reverse under net generation conditions, and equilibrium $np=n_i^2$ should make the net rate vanish.

Radiative and Auger processes dominate in different density regimes. Band-to-band radiative recombination often scales with $B(np-n_i^2)$ and produces photons, central to LEDs and direct-gap solar cells. Auger recombination grows roughly cubically with carrier density through coefficients multiplying $n$ or $p$ times $np-n_i^2$, becoming important in heavy injection. Coefficients depend on material, temperature, degeneracy, and band structure. Adding all published rates without checking overlapping calibration can double-count measured lifetime behavior.

Generation must be spatially, spectrally, and dimensionally consistent. Optical generation derives from absorbed photon flux, not optical power density alone. Reflection, interference, polarization, complex refractive index, and wavelength-dependent absorption determine where pairs appear. Impact ionization depends strongly on field and carrier energy, while thermal generation follows detailed balance with recombination. Mapping an optical solution onto an electrical mesh must conserve generated pairs. In two-dimensional simulations, current and generation outputs require a declared out-of-plane depth.

Avalanche multiplication pushes local drift–diffusion toward its validity edge. Local-field ionization coefficients assume carrier energy responds instantaneously to field; dead-space and nonlocal history matter in short high-field regions. Generated pairs feed Poisson and current, creating strong positive feedback and possible breakdown branches. Continuation, current control, and external-circuit coupling may be needed to follow a stable operating path. A stationary solver failure is not by itself a physical breakdown criterion, and convergence achieved by excessive damping is not proof of correct avalanche physics.

Trap dynamics introduce memory and additional state variables. Occupation evolves through capture and emission rather than always following steady equilibrium. Trapped charge shifts threshold, modulates recombination, and can produce hysteresis, random telegraph signals, bias-temperature instability, and persistent photoconductivity. A stationary trap model assumes observation time is long relative to its kinetics. Transient simulation needs consistent initial occupation and charge conservation when a carrier enters or leaves a trap. Broad trap distributions can span many decades of time.

Heat couples back through nearly every transport coefficient. Joule heating, recombination heat, Thomson/Peltier terms, and optical absorption raise lattice temperature. Temperature changes mobility, intrinsic density, bandgap, ionization, diffusion, and recombination. Solving an electrical model at fixed ambient temperature can overpredict current or miss thermal runaway. The heat source must avoid double-counting electrochemical work, and thermal boundary resistance may dominate temperature. Coupled electrothermal convergence should include both terminal power balance and spatial heat-flux balance.

Thermoelectric transport needs gradients beyond the elementary Einstein picture. Temperature gradients drive Seebeck currents and modify carrier diffusion through density-of-states and band-edge temperature dependence. A quasi-Fermi-gradient formulation can organize these terms, but transport coefficients must satisfy compatible thermodynamics. Ignoring thermodiffusion while allowing strong self-heating may violate equilibrium in a nonuniform temperature field. Energy-transport models become preferable when carrier and lattice temperatures differ.

Magnetic fields turn scalar mobility into a tensor response. The Lorentz force produces Hall current and magnetoresistance, rotating carrier flux relative to electric and electrochemical-potential gradients. Electron and hole Hall factors need not equal one because scattering is energy dependent. The transport tensor must preserve nonnegative entropy production in its symmetric part, while its antisymmetric Hall part changes direction without dissipation. A scalar field-dependent mobility cannot reproduce these effects.

<svg viewBox="0 0 760 470" xmlns="http://www.w3.org/2000/svg" font-family="Arial,sans-serif"><rect x="0" y="0" width="760" height="470" rx="18" fill="#0d1117"/><text x="380" y="34" fill="#e6edf3" font-size="21" font-weight="700" text-anchor="middle">Contacts determine injection, extraction, and reference potential</text><text x="380" y="60" fill="#8b949e" font-size="12" text-anchor="middle">A boundary voltage alone does not specify carrier exchange</text><rect x="45" y="115" width="170" height="240" rx="12" fill="#161b22" stroke="#d29922" stroke-width="3"/><text x="130" y="150" fill="#e3b341" font-size="14" font-weight="700" text-anchor="middle">Ohmic contact</text><text x="130" y="205" fill="#e6edf3" font-size="11" text-anchor="middle">potential + carrier equilibrium</text><text x="130" y="250" fill="#c9d1d9" font-size="11" text-anchor="middle">low injection barrier</text><rect x="295" y="115" width="170" height="240" rx="12" fill="#161b22" stroke="#58a6ff" stroke-width="3"/><text x="380" y="150" fill="#79c0ff" font-size="14" font-weight="700" text-anchor="middle">Schottky contact</text><path d="M325 280L380 175L435 280" fill="none" stroke="#58a6ff" stroke-width="4"/><text x="380" y="315" fill="#c9d1d9" font-size="11" text-anchor="middle">thermionic barrier current</text><rect x="545" y="115" width="170" height="240" rx="12" fill="#161b22" stroke="#a371f7" stroke-width="3"/><text x="630" y="150" fill="#d2a8ff" font-size="14" font-weight="700" text-anchor="middle">Selective contact</text><text x="630" y="205" fill="#e6edf3" font-size="11" text-anchor="middle">one carrier extracted</text><text x="630" y="250" fill="#c9d1d9" font-size="11" text-anchor="middle">surface velocity / barrier</text><text x="380" y="410" fill="#e6edf3" font-size="13" text-anchor="middle">Terminal current must equal the boundary flux integrated with one sign convention.</text><text x="380" y="442" fill="#8b949e" font-size="11" text-anchor="middle">Contact resistance and external circuits may belong outside the semiconductor domain.</text></svg>

Ohmic contacts impose more physics than a fixed voltage. An ideal ohmic contact sets electrostatic reference and carrier populations consistent with local doping, temperature, band structure, and applied electrochemical potential, allowing majority carriers to enter without a limiting barrier. Heavy contact doping may justify this approximation but can require degeneracy and bandgap narrowing. Imposing equilibrium minority density under strong injection can artificially absorb carriers. Real contact resistivity, current crowding, and metal spreading resistance may need explicit boundary or circuit elements.

Schottky contacts require a barrier-current relation. Metal work function, semiconductor electron affinity, interface states, image-force lowering, tunneling, and interfacial layers determine injection. Thermionic-emission boundary flux depends exponentially on barrier and quasi-Fermi separation; thermionic-field emission matters for heavily doped barriers. Pinning can decouple the barrier from ideal work-function difference. Prescribing both carrier density and thermionic current overconstrains the boundary unless the formulation reconciles them.

Surface recombination is a boundary flux, not a volume lifetime. Electron and hole exchange at an interface can be expressed through surface recombination velocities and trap occupancy. Converting it to a volumetric rate by dividing by an arbitrary mesh-cell width makes the physics mesh dependent. Passivation changes capture kinetics and fixed charge simultaneously, so a fitted velocity may not transfer across bias or injection. Global continuity should count the surface flux with the same sign as bulk recombination.

Insulating boundaries block normal carrier flux but may still carry electrostatic charge. Setting $\mathbf J_n\cdot\mathbf n=\mathbf J_p\cdot\mathbf n=0$ prevents carrier crossing. Poisson can simultaneously impose normal displacement from fixed surface charge or a dielectric interface. Conflating electrical insulation with zero electric field removes surface-charge physics. Symmetry boundaries use the same zero-normal-flux form only when geometry, material, sources, and solution are truly mirror symmetric.

Periodic boundaries require compatible potential drop and carrier driving. A periodic unit cell under zero macroscopic bias identifies both values and fluxes across paired faces. A driven periodic conductor may use an affine potential or quasi-Fermi offset, not simply identical potential. Net generation and recombination must be compatible with periodic carrier balance. Periodic electrostatics also needs charge neutrality or a compensating background and a gauge. Otherwise the apparent steady state violates the integrated equations.

Initial conditions matter even when only the final DC point is desired. Transient carrier and trap states choose the basin approached by a nonlinear system with hysteresis or multiple steady branches. Equilibrium initialization is effective at zero bias; continuation from a neighboring solved bias is usually better for a sweep. Arbitrary tiny carrier densities can create enormous logarithms and unphysical space charge. A DC solution obtained by pseudo-time stepping should be checked independently for steady residual and path dependence.

Terminal current includes displacement current in transient operation. Conduction current from electrons and holes need not be spatially constant during charge storage. Maxwell displacement current $\partial\mathbf D/\partial t$ completes total current continuity in the electroquasistatic regime. Omitting it distorts capacitance, switching current, and high-frequency admittance. Integrating charge change and terminal total current provides a strong transient conservation check. At frequencies where wave propagation matters, full Maxwell coupling replaces the quasistatic approximation.

<svg viewBox="0 0 760 470" xmlns="http://www.w3.org/2000/svg" font-family="Arial,sans-serif"><rect x="0" y="0" width="760" height="470" rx="18" fill="#0d1117"/><text x="380" y="34" fill="#e6edf3" font-size="21" font-weight="700" text-anchor="middle">Scharfetter–Gummel flux resolves drift and diffusion together</text><text x="380" y="60" fill="#8b949e" font-size="12" text-anchor="middle">Exponential fitting preserves equilibrium across a cell</text><circle cx="155" cy="230" r="24" fill="#58a6ff"/><circle cx="605" cy="230" r="24" fill="#58a6ff"/><text x="155" y="236" fill="#0d1117" font-size="13" font-weight="700" text-anchor="middle">i</text><text x="605" y="236" fill="#0d1117" font-size="13" font-weight="700" text-anchor="middle">i+1</text><path d="M179 230H581" stroke="#30363d" stroke-width="8"/><path d="M180 295C285 295 360 140 580 155" fill="none" stroke="#d29922" stroke-width="4"/><text x="380" y="125" fill="#e3b341" font-size="12" text-anchor="middle">potential and carrier variation inside cell</text><path d="M250 350H510" stroke="#3fb950" stroke-width="7"/><polygon points="510,350 492,338 492,362" fill="#3fb950"/><text x="380" y="385" fill="#7ee787" font-size="13" text-anchor="middle">one conservative face flux Jᵢ₊½</text><rect x="245" y="410" width="270" height="38" rx="8" fill="#161b22" stroke="#a371f7"/><text x="380" y="435" fill="#d2a8ff" font-size="12" text-anchor="middle">Bernoulli function handles large Δψ stably</text></svg>

Naive centered differencing can fail when drift dominates diffusion. A cell Peclet number compares electrostatic potential drop with thermal voltage. At large values, centered carrier-density gradients can create negative concentrations, oscillations, and nonphysical current. Upwinding stabilizes drift but adds artificial diffusion and may destroy exact thermal equilibrium. The discretization should treat field and density coupling as one flux rather than two unrelated approximations, particularly across depletion regions and high barriers.

Scharfetter–Gummel discretization exponentially fits the cell problem. Assuming approximately constant field and coefficients across an edge, solve the one-dimensional drift–diffusion relation analytically to obtain a Bernoulli-function flux. It remains conservative, preserves the discrete equilibrium relation, and handles large potential drops more robustly than centered differences. Stable evaluation near zero uses a series or expm1-type implementation to avoid cancellation. Strong coefficient variation, multidimensional anisotropy, degeneracy, and abrupt heterojunctions require generalized fluxes rather than blind reuse of the elementary formula.

Log-density variables enforce positivity but change nonlinear conditioning. Solving for $\ln n$ and $\ln p$ prevents negative densities and spans many decades common in depleted devices. Quasi-Fermi variables similarly align unknowns with electrochemical driving and equilibrium. Density variables may be simpler for finite volumes and charge conservation. Each formulation has different Jacobian scaling and boundary transformations. Switching variables is not merely cosmetic; convergence, interpolation, and stopping norms must be interpreted in physical density and current afterward.

Gummel iteration exploits the system's physical block structure. Solve Poisson with fixed carriers, then electron continuity, then hole continuity, updating models and repeating. Damping or nonlinear Poisson variants improve robustness. Gummel is inexpensive per step and often tolerant of poor initial guesses, but can converge slowly under strong coupling, high injection, avalanche, or self-heating. Convergence must evaluate the original coupled residual, not only the change between damped iterates, because heavy damping can make updates small while equations remain unsatisfied.

Newton's method trades a coupled Jacobian for rapid local convergence. Assemble derivatives of Poisson, continuity, recombination, mobility, statistics, and boundary fluxes with respect to all unknowns. A correct Jacobian yields near-quadratic convergence close to a nonsingular solution. Line searches, trust regions, voltage continuation, and positivity-aware variables globalize the method. An approximate Jacobian that omits strong field or recombination derivatives may behave worse than Gummel. Automatic differentiation helps consistency but does not repair nondifferentiable empirical models.

Continuation is a physical route through a difficult nonlinear landscape. Ramp contact voltage, illumination, doping, interface charge, avalanche strength, or quantum correction from an easier solved state. Adaptive step size grows after easy convergence and shrinks near sharp response. Current-controlled continuation can pass voltage turning points that defeat a simple voltage sweep. The followed branch depends on circuit and stability; numerical continuation can trace mathematically unstable states that an experiment never occupies. Record direction and step history when hysteresis exists.

Scaling must accommodate densities spanning many orders of magnitude. Normalize potential by thermal voltage, length by a device or Debye scale, density by a representative doping, and current by a compatible flux. Row and variable scaling prevent Poisson residual units from overwhelming continuity residuals in a combined norm. Absolute tolerances protect near-zero currents; relative tolerances control large signals. Scaling a residual for linear algebra is distinct from defining physical convergence. Always translate the final tolerances back into volts, charge, particle balance, and terminal current.

Linear solver structure changes across nonlinear formulations. A decoupled Poisson block may be symmetric positive definite after anchoring, but the full Newton Jacobian is generally nonsymmetric and indefinite. GMRES or direct sparse factorization is common; block preconditioners approximate Poisson and electron/hole Schur complements. Algebraic multigrid effective for Poisson may struggle with advective continuity blocks unless tailored. Reordering and scaling affect fill and robustness. Solver choice must follow the assembled matrix, not the elliptic label attached to one subequation.

Time integration must resolve both storage and stiff reaction. Backward Euler is robust and dissipative; higher-order backward differentiation or implicit Runge–Kutta improves accuracy for smooth transients. Explicit stepping is restricted by diffusion, drift, dielectric relaxation, and reaction scales. Adaptive methods need error estimates in variables that reflect terminal observables, not only dominant majority density. Discontinuous voltage steps create mathematical high-frequency content; a physically finite ramp often yields a more meaningful and numerically tractable response.

<svg viewBox="0 0 760 470" xmlns="http://www.w3.org/2000/svg" font-family="Arial,sans-serif"><rect x="0" y="0" width="760" height="470" rx="18" fill="#0d1117"/><text x="380" y="34" fill="#e6edf3" font-size="21" font-weight="700" text-anchor="middle">Model validity is set by competing length and time scales</text><text x="380" y="60" fill="#8b949e" font-size="12" text-anchor="middle">Drift–diffusion lives between equilibrium and ballistic transport</text><line x1="90" y1="235" x2="670" y2="235" stroke="#30363d" stroke-width="8"/><circle cx="160" cy="235" r="20" fill="#3fb950"/><circle cx="380" cy="235" r="24" fill="#58a6ff"/><circle cx="600" cy="235" r="20" fill="#f85149"/><text x="160" y="190" fill="#7ee787" font-size="14" font-weight="700" text-anchor="middle">Local equilibrium</text><text x="160" y="285" fill="#c9d1d9" font-size="11" text-anchor="middle">L ≫ mean free path</text><text x="380" y="180" fill="#79c0ff" font-size="15" font-weight="700" text-anchor="middle">Drift–diffusion</text><text x="380" y="285" fill="#c9d1d9" font-size="11" text-anchor="middle">mobility + diffusion + continuity</text><text x="600" y="190" fill="#ff7b72" font-size="14" font-weight="700" text-anchor="middle">Nonlocal / ballistic</text><text x="600" y="285" fill="#c9d1d9" font-size="11" text-anchor="middle">L ≲ relaxation length</text><path d="M185 350H575" stroke="#d29922" stroke-width="4" stroke-dasharray="9 6"/><text x="380" y="385" fill="#e3b341" font-size="12" text-anchor="middle">field, density, energy, Debye, mean-free-path, and transit scales</text><text x="380" y="430" fill="#e6edf3" font-size="12" text-anchor="middle">A good fit does not extend a model beyond the assumptions used to close it.</text></svg>

Validity is governed by scale separation rather than device generation labels. Compare mean free path and energy-relaxation length with channel length, barrier width, and field-variation scale; compare momentum, energy, recombination, dielectric, transit, and drive times. A nominally nanoscale device may contain diffusive reservoirs and a ballistic constriction, requiring hybrid treatment. Conversely, a large device can develop a sharp high-field region beyond local closure. Mesh refinement cannot cure a continuum-model validity failure.

Quantum confinement can be corrected approximately without becoming quantum transport. Density-gradient and effective-potential models shift carrier density away from interfaces and raise confinement energy while retaining drift–diffusion current. Their calibration depends on effective mass, orientation, boundary conditions, and dimensionality. Self-consistent Poisson–Schrödinger supplies subband charge more directly but still needs a transport occupation model. Neither approach captures coherent tunneling, interference, or contact mode injection in the NEGF sense.

Tunneling must enter as a transfer mechanism consistent with continuity. Band-to-band, trap-assisted, Fowler–Nordheim, and direct tunneling models create generation terms or boundary/interface fluxes. Their exponential sensitivity to field, barrier shape, effective mass, and band alignment makes mesh and electrostatic accuracy decisive. Depositing pair generation at the wrong spatial location can violate energy or current balance. Combining a nonlocal tunneling path with local impact ionization requires careful avoidance of double counting.

Hydrodynamic transport adds carrier energy when local mobility is insufficient. Energy-balance equations evolve carrier temperature or mean energy, and flux laws include energy gradients and temperature-dependent relaxation. They can reproduce velocity overshoot and hot-carrier effects more efficiently than a full Boltzmann solver. Closure coefficients still come from kinetic assumptions or calibration, boundary conditions for energy are difficult, and numerical stiffness increases. A more elaborate model is not automatically more predictive without verified energy-relaxation data.

Boltzmann, Monte Carlo, and NEGF define distinct escalation paths. Deterministic Boltzmann solvers resolve distribution functions in phase space; ensemble Monte Carlo samples semiclassical trajectories and scattering; nonequilibrium Green functions treat quantum-coherent states and contact injection. Each adds information that drift–diffusion integrates out, at substantial computational and calibration cost. Cross-model comparison should hold band structure, geometry, contacts, and scattering assumptions as consistent as possible. Disagreement then diagnoses closure limits rather than arbitrary parameter differences.

Compact models are reductions of transport, not replacements for physical validation. MOSFET, diode, solar-cell, and LED compact equations encode selected drift–diffusion behavior into terminal relations for circuits. Parameters can be extracted from measurement or numerical simulation. A compact model may conserve charge and reproduce I–V while hiding internal field, self-heating, or breakdown mechanisms. Use detailed transport to establish parameter dependence and validity range, then verify the reduced model across bias, geometry, temperature, and frequency.

Device examples emphasize different portions of the same system. A long-channel MOSFET emphasizes field-dependent channel charge and mobility; a pn diode emphasizes minority diffusion and depletion electrostatics; a bipolar transistor emphasizes injection and recombination; a solar cell emphasizes optical generation and selective extraction; an LED emphasizes radiative recombination and current crowding; a power device emphasizes high field, heating, and avalanche. A generic solver needs model switches, but each switch must be tied to evidence and not enabled merely because it exists.

<svg viewBox="0 0 760 470" xmlns="http://www.w3.org/2000/svg" font-family="Arial,sans-serif"><rect x="0" y="0" width="760" height="470" rx="18" fill="#0d1117"/><text x="380" y="34" fill="#e6edf3" font-size="21" font-weight="700" text-anchor="middle">Verification closes three independent evidence loops</text><text x="380" y="60" fill="#8b949e" font-size="12" text-anchor="middle">Equation satisfaction, numerical convergence, and experimental meaning differ</text><circle cx="380" cy="145" r="72" fill="#161b22" stroke="#3fb950" stroke-width="3"/><text x="380" y="137" fill="#7ee787" font-size="14" font-weight="700" text-anchor="middle">Conservation</text><text x="380" y="162" fill="#c9d1d9" font-size="11" text-anchor="middle">charge and terminal flux</text><circle cx="190" cy="315" r="72" fill="#161b22" stroke="#58a6ff" stroke-width="3"/><text x="190" y="307" fill="#79c0ff" font-size="14" font-weight="700" text-anchor="middle">Numerics</text><text x="190" y="332" fill="#c9d1d9" font-size="11" text-anchor="middle">mesh, timestep, residual</text><circle cx="570" cy="315" r="72" fill="#161b22" stroke="#d29922" stroke-width="3"/><text x="570" y="307" fill="#e3b341" font-size="14" font-weight="700" text-anchor="middle">Validation</text><text x="570" y="332" fill="#c9d1d9" font-size="11" text-anchor="middle">I–V, C–V, optics, heat</text><path d="M332 198L240 265M428 198L520 265M262 315H498" stroke="#a371f7" stroke-width="3" stroke-dasharray="7 5"/><text x="380" y="430" fill="#e6edf3" font-size="12" text-anchor="middle">Agreement in one loop cannot substitute for missing evidence in another.</text></svg>

Global conservation is the first nonnegotiable verification target. Integrate each continuity equation over the domain and compare stored-carrier change, contact particle flux, bulk generation–recombination, and surface exchange. Add electron, hole, and displacement currents with consistent terminal orientation. In steady two-terminal dark operation, total current should agree at both contacts within declared tolerance. Exact global balance does not prove local accuracy, but imbalance immediately exposes source signs, boundary flux, nonlinear convergence, or postprocessing errors.

Manufactured solutions verify code paths that analytical devices do not cover. Choose smooth potential and positive carrier fields, derive Poisson sources, continuity sources, and boundary data from the implemented equations, then recover them on a mesh sequence. Exercise variable mobility, recombination derivatives, heterointerfaces, each contact type, and transient storage separately. Measure potential, density, quasi-Fermi, current, and conservation error. Expected convergence rates should appear before solver tolerance or roundoff dominates.

Equilibrium, resistor, and low-injection diode limits form a compact benchmark ladder. Equilibrium tests exact drift–diffusion cancellation. A uniformly doped bar under small bias tests Ohm's law $J=q(\mu_nn+\mu_pp)E$. A long neutral region with injected minority carriers tests exponential diffusion length $L=\sqrt{D\tau}$. An ideal long diode tests the Shockley exponential only within its assumptions. Moving through this ladder isolates electrostatics, flux, continuity, recombination, and contact defects before attempting a full transistor.

Mesh studies must resolve Debye layers, depletion edges, optical absorption, and transport gradients. Refine geometry and source projection consistently while tightening algebraic tolerances. Compare terminal current, stored charge, recombination integral, peak field, and a local current profile on at least three credible meshes. Pointwise field at an ideal sharp corner may not converge, so round the physical geometry or use an integrated output. Changing a numerical interface width with the mesh changes the physical model and invalidates an order estimate.

Bias-step convergence is separate from spatial convergence. A coarse voltage sweep can skip snapback, hysteresis, threshold structure, or sharp recombination changes even if every point is fully converged. Repeat with smaller continuation steps and both sweep directions. For transient ramps, refine time step and input waveform together. Interpolating a sparse I–V curve may conceal negative differential resistance or convergence branch changes. Store state hashes and predecessor bias so a result's continuation history is reproducible.

Validation requires outputs filtered through the experiment. Compare terminal I–V with series resistance and instrument compliance, C–V with frequency and trap response, luminescence with optical extraction, temperature with sensor placement, and transient current with circuit parasitics. Internal carrier density is rarely measured directly. Calibrating mobility, lifetime, and contact resistance to the same curve used for validation is parameter fitting, not independent prediction. Reserve geometries, temperatures, biases, or observables for validation.

Sensitivity and identifiability should precede aggressive calibration. Mobility, lifetime, contact resistance, interface charge, doping, dimensions, and temperature can compensate each other in terminal curves. Local derivatives or adjoints show which outputs respond to which parameters, while profile likelihood or Bayesian analysis reveals correlated uncertainty. A parameter with little sensitivity cannot be reliably extracted. Spatially resolved or frequency-dependent measurements can break degeneracies that DC I–V cannot.

Uncertainty propagation distinguishes numerical precision from predictive confidence. Mesh and solver errors may be below one percent while uncertain mobility, trap density, geometry, and contact barrier produce orders-of-magnitude current variation. Sample physically correlated parameters and preserve constraints such as positive lifetimes. Report distributions or intervals for decision outputs, not only a best-fit contour. Model-form uncertainty—local transport versus nonlocal, or one recombination law versus another—requires comparison across plausible closures rather than parameter sampling alone.

A diagnostic workflow should identify the failed layer before changing parameters. Separate model validity, boundary closure, discretization, nonlinear solution, linear algebra, and measurement mapping. Negative density points to variable or flux treatment; equilibrium current points to statistics or discretization inconsistency; unequal steady terminal currents point to continuity or convergence; mesh-dependent surface recombination points to a volume conversion error; bias-path dependence can be physical hysteresis or branch-selection failure. Each symptom demands a targeted invariant, not arbitrary damping.

Declare device geometry, materials, temperature, doping, traps, and reference energies
  -> Write Poisson, electron continuity, hole continuity, current, generation, and recombination signs
  -> Check local-equilibrium, diffusive-length, field, degeneracy, and quasistatic validity scales
  -> Assign contact injection, insulating/symmetry, interface, optical, thermal, and circuit boundaries
  -> Choose density, log-density, or quasi-Fermi variables and conservative spatial fluxes
  -> Scale variables and assemble Gummel blocks or the coupled Newton residual and Jacobian
  -> Continue from equilibrium through bias, illumination, temperature, or model strength
  -> Require coupled residual, positivity, global carrier balance, and terminal-current agreement
  -> Run manufactured, equilibrium, resistor, diode, mesh, timestep, and bias-step benchmarks
  -> Compare declared observables through circuit, optical, thermal, and instrument forward models
  -> Archive equations, parameters, boundary map, mesh, tolerances, branch history, and hashes

Drift–diffusion troubleshooting becomes systematic when every numerical symptom is mapped back to a conserved quantity or closure assumption. A tiny update with a large original residual means damping has hidden nonconvergence. A tiny algebraic residual with wrong equilibrium current means the discrete flux or statistics are inconsistent. A stable mesh sequence with wrong experiment points toward material, contact, heating, measurement, or validity errors. A measured fit that changes wildly under parameter perturbation indicates poor identifiability rather than a uniquely characterized device.

The one-dimensional steady continuity equation provides a transparent sign test. Integrating it across a slab states that the difference between outgoing and incoming carrier currents equals the integrated net recombination or generation with a charge-dependent sign. If $G=R$, each carrier current is constant. When recombination transfers one electron and one hole out of their mobile populations, electron and hole current components change oppositely along the device while their conventional sum remains constant. Plotting cumulative source integrals beside face currents localizes imbalance to a cell or boundary.

The minority-carrier diffusion equation is a controlled reduction of the full system. In a quasi-neutral region with negligible electric field perturbation, low injection, constant $D$ and lifetime, the excess minority density satisfies a second-order equation with diffusion length $L=\sqrt{D\tau}$. Its exponential solution explains diode injection profiles and collection probability. Near depletion fields, high injection, spatially varying lifetime, degeneracy, or significant majority perturbation, the reduction fails and the coupled equations must be restored.

The long-channel charge-sheet approximation is another reduction with a clear domain. It integrates inversion charge normal to the MOS interface and transports that sheet laterally using a gradual-channel field. This yields intuitive MOSFET current formulas and compact-model structure. It loses accuracy near source/drain junctions, short-channel barriers, two-dimensional fringing, velocity overshoot, and strong self-heating. Comparing it with two-dimensional drift–diffusion separates geometric field effects from mobility assumptions.

Solar-cell drift–diffusion couples optics, electrostatics, and selective contacts. Spectral absorption creates $G(\mathbf x,\lambda)$; minority diffusion and depletion drift collect carriers; bulk and surface recombination set loss; contacts extract one carrier preferentially. The current–voltage curve yields short-circuit current, open-circuit voltage, fill factor, and efficiency only after optical input power and area are defined. Sesame and nextnano documentation use coupled Poisson and continuity equations for this class of calculation, making solar cells valuable end-to-end benchmarks.

LED simulation reverses much of the photovoltaic causal chain. Contact injection creates electron and hole populations, transport brings them into an active region, radiative and nonradiative rates set internal quantum efficiency, and optical extraction converts emitted photons to measured power. Current crowding, polarization charge, heterobarriers, self-heating, and Auger loss create strong spatial coupling. Matching total light output while misplacing recombination can give the wrong thermal and reliability prediction, so spatial emission evidence matters.

Power-device simulation stresses high-field closure and electrothermal feedback. Drift regions trade breakdown voltage against on-resistance; junction curvature concentrates field; avalanche generates carriers; conductivity modulation changes charge; heating lowers mobility and can raise leakage. External circuit impedance selects whether breakdown settles, snaps back, or runs away. A voltage-driven stationary solve without circuit or thermal coupling may follow an irrelevant branch. Verification should include blocking-state charge balance, on-state current balance, breakdown mesh sensitivity, and total electrical-to-thermal power.

Cryogenic drift–diffusion needs more than changing $T$ in thermal voltage. Dopant freeze-out, incomplete ionization, band tails, degenerate statistics, field-assisted ionization, trap kinetics, and mobility all change. Very long relaxation and recombination times challenge steady assumptions, while tiny intrinsic densities challenge floating-point scale. Local equilibrium may fail in short channels even when low lattice temperature suggests small thermal velocity. Calibration must use cryogenic-specific data rather than extrapolating room-temperature formulas.

Mixed-dimensional devices require careful source and current normalization. A two-dimensional cross-section may report amperes per meter of assumed depth; an axisymmetric model integrates around $2\pi r$; a sheet material carries density per area and current per length. Coupling a two-dimensional channel to three-dimensional contacts or optical generation requires conservative dimensional transfer. An unexplained width multiplier can make an otherwise correct I–V curve numerically arbitrary.

Reproducible transport studies preserve the whole closure stack. Record band parameters, density of states, statistics, mobility components and combination rule, all recombination/generation models, ionization and narrowing, contact relations, interface conditions, quantum or thermal corrections, circuit elements, meshes, variable formulation, flux scheme, scaling, nonlinear and linear tolerances, continuation path, and output integration. Parameter names alone are insufficient because software versions can change defaults and formulations.

The final interpretation should distinguish density, particle flux, conventional current, electrostatic field, quasi-Fermi driving, and measured terminal response. Density can be enormous where mobility is low, current can be constant while its electron and hole shares change, and a steep electric field can coexist with zero equilibrium current. Quasi-Fermi gradients identify dissipative driving more directly than band bending alone. Terminal data combine the internal solution with contact, displacement, circuit, optical, and thermal mappings.

Read drift–diffusion transport through a conservation-closure-and-validity lens rather than a drift-term-plus-diffusion-term lens.