Getting Started

Every change re-runs the simulation automatically (with a short debounce) in a background worker thread, so the interface stays responsive even for long, high-resolution anneals.

Interface Guide

Profile Library

Copy and paste these directly into the Concentration input field to test different scenarios.

Material Systems

Four thermoelectric hosts are modeled. Each has its own transport-property fit (see Physics Model for the equations and Default Parameters for every constant) and its own allowed temperature window, enforced by the Thermal Gradient card and tied to that material's real stability limit rather than one generic range.

A. PbSe (Lead Selenide)

A rock-salt IV–VI narrow-gap semiconductor and, alongside PbTe, one of the leading mid-to-high-temperature thermoelectrics: strong band convergence and intrinsically low lattice thermal conductivity give optimized bulk/nanostructured samples zT ≈ 1.3–1.8 in the 800–950 K range. Dopants here are treated as simple substitutional/interstitial point defects (Na, Ag, Cu, Cl).

Melting point	1338–1351 K (1065–1078°C)
Tool's allowed T range	300–1300 K — ~45–50 K margin below melting
Band gap Eg (model)	0.28 eV
Available dopants	Na, Ag, Cu, Cl, Custom

B. SnSe (Tin Selenide)

A layered orthorhombic IV–VI semiconductor, famous since Zhao et al. (Nature, 2014) for ultralow single-crystal lattice thermal conductivity along the b-axis, giving a record zT ≈ 2.6–2.8 near 800–950 K. It undergoes a structural Pnma→Cmcm phase transition around 750 K that further improves performance in the high-T phase; this tool adds a small extra thermal-conductivity term above 750 K to reflect that (see Physics Model, part D).

Melting point	1134 K
Tool's allowed T range	300–1100 K — keeps the 750 K transition and the near-melting zT peak accessible
Band gap Eg (model)	0.86 eV
Available dopants	Na, Ag, K, Custom

C. MgAgSb — α-phase only

A low-to-mid-temperature p-type thermoelectric valued for a high power factor and low lattice thermal conductivity near room temperature (reported peak zT ≈ 1.2–1.4 near 475 K), attractive for waste-heat and wearable applications below ~200°C. MgAgSb exists in three polymorphs; only the low-temperature α-phase has this favorable transport behavior. It transforms to a β-phase (and eventually γ) above roughly 560 K, with different, unmodeled transport — so this tool restricts simulation to the α-phase window only, unlike the other three materials where the ceiling is set by melting.

α→β transition	≈560–573 K (true melting point is much higher, but irrelevant — the modeled phase doesn't survive that far)
Tool's allowed T range	300–560 K — the entire range stays in the α-phase
Band gap Eg (model)	0.18 eV
Available dopants	Ag, Zn, Ni, Custom

D. Bi₂(Se,Te)₃

The benchmark commercial thermoelectric for near-room-temperature Peltier cooling and power generation — essentially every commercial TE module uses a Bi₂Te₃-based alloy. Alloying with Bi₂Se₃ (n-type) or Sb₂Te₃ (p-type) reduces lattice thermal conductivity below that of pure Bi₂Te₃; this tool's parameters target that optimized-alloy regime (k_lat ≈ 0.6 W/m/K at 300 K, vs. ≈1.2 W/m/K reported for unalloyed single-crystal Bi₂Te₃). Its narrow gap causes significant bipolar carrier excitation even near 300 K — the main reason its useful range is so much lower than the other three materials here, and why the field keeps searching for wider-gap alternatives for higher-T service. This tool's simplified single-parabolic-band model reaches zT ≈ 0.9–1.1 near room temperature (optimized, nanostructured real samples reach higher, up to ≈1.2–1.4, via mechanisms — band convergence, nanostructuring — this simple model doesn't capture).

Melting point	858 K (585°C)
Tool's allowed T range	300–800 K
Band gap Eg (model)	0.14 eV — narrowest of the four, dominant limitation
Available dopants	Cu, Ag, I, Custom

Cu's D₀/E_a (1.7×10⁻&sup6; m²/s, 0.25 eV) come from a 2-point Arrhenius fit to two independent literature diffusion measurements — single-crystal fast-axis diffusion ≈1×10⁻¹⁰ m²/s at 300 K, and an effective sintered-pellet value ≈4×10⁻&sup9; m²/s at 483 K. It's the only dopant preset in this tool anchored to real measurements rather than estimated by analogy; see Default Parameters for the full table and its caveats.

Physics Model

A. Conservation Law

Dopant concentration $C(x,t)$ obeys a 1D continuity equation — the flux $J$ defined below is the only thing that moves mass around, and the zero-flux boundaries (part E, Numerical Method) mean nothing enters or leaves the leg:

B. The Soret Effect

Dopants migrate under the temperature gradient via thermodiffusion (thermomigration). The net atomic flux combines ordinary Fickian diffusion with a thermal-gradient (Soret) term:

The direction of net migration is set by the sign of the Heat of Transport ($Q^*$):

The diffusivity itself follows an Arrhenius law with prefactor $D_0$ and activation energy $E_a$ (both per-dopant, set in the sidebar or a preset):

Note the model assumes complete ionization: the simulated dopant concentration $C$ is used directly as the free-carrier concentration in every transport equation below, with no separate ionization-fraction or trap-state treatment.

C. Mobility, Conductivity & Seebeck Coefficient

Carrier mobility follows a power law in temperature, with material-specific prefactor $\mu_0$ and exponent $m$ (values in Default Parameters):

The Seebeck coefficient starts from a Pisarenko-like baseline at 300 K (material constants $A$, $B$, $C_{ref}$), then gets a linear high-T boost before bipolar suppression is applied:

$s$ is the per-material "Seebeck slope" listed in Default Parameters. For SnSe specifically, $S_T$ is additionally capped at 550 µV/K before bipolar suppression, to prevent runaway growth at low doping/high T that the simple linear boost would otherwise produce — the other three materials didn't need this in testing but could in principle at extreme enough inputs.

D. Thermal Conductivity & Figure of Merit (zT)

Total thermal conductivity sums a material-specific lattice term, an electronic term from the Wiedemann–Franz law (Lorenz number $L_0=2.4\times10^{-8}$ W·Ω/K²), a bipolar term, and — for SnSe above 750 K only — a phase-transition term:

which combines with $S_T$ (after bipolar suppression, part below) and $\sigma$ into the figure of merit:

where $S$ is the Seebeck coefficient, $\sigma$ the electrical conductivity, and $\kappa_{tot}$ the total thermal conductivity (lattice + electronic + bipolar [+ phase, SnSe only]). The Avg ZT badge on the ZT tab averages this value across the leg at the final simulated (or truncated) time.

E. High-Temperature Accuracy (Bipolar Correction)

To prevent unrealistic $zT$ predictions at high temperatures (where standard single-band models fail), this simulation uses a bipolar correction calibrated to typical TE data. The activation factor is:

• Effective Band Gap: Each system uses its own $E_g$ estimate (Default Parameters) to calculate minority-carrier activation — the smaller $E_g$ is, the lower the temperature at which bipolar effects start eating into $S$ and $zT$. This is why Bi₂Te₃ ($E_g=0.14$ eV) rolls over so much earlier than PbSe ($E_g=0.28$ eV).
• Seebeck Suppression: As $T$ rises, $\xi$ grows and the $\frac{1}{1+5\xi}$ factor pulls $S$ down, reproducing the Seebeck "rollover" seen in real experiments.
• Thermal Conductivity: The $\kappa_{bi}$ term above grows with the same $\xi$, further lowering $zT$ toward realistic limits at high T.

Default Parameters

Every constant behind the Material Systems and Physics Model sections above, for direct inspection or reproduction outside this tool. Diffusion parameters (D₀, E_a, Q*) are order-of-magnitude estimates consistent with known diffusion trends (fast interstitial dopants like Cu/Ag vs. slower substitutional ones like Cl/I) rather than independently verified literature values for each exact dopant–host pair, with the single exception noted below.

Dopant Presets

Only Bi₂Te₃:Cu is anchored to real measurements (a 2-point Arrhenius fit across 300–483 K — see Material Systems, part D). Every other row, across all four systems, is an order-of-magnitude estimate. Q* values were not independently found in literature for any of these exact dopant–host pairs; their magnitudes (0.10–0.30 eV) match the right order of magnitude for reported thermomigration heat-of-transport values but should be treated as illustrative, not authoritative.

Transport-Model Constants

$\mu(T)=\mu_0(T/300)^{-m}$; $\kappa_{lat}(T) = \kappa_{lat}(300\text{K})\times(300/T)$ for all four systems (i.e. $\propto 1/T$). SnSe additionally caps $S_T$ at 550 µV/K and adds $\kappa_{phase}=0.1\times\frac{T-750}{50}$ for $T>750$ K. Universal constants shared by all four systems: Lorenz number $L_0=2.4\times10^{-8}$ W·Ω/K², bipolar prefactor 20 (i.e. $\kappa_{bi}=20\,(T/300)\,\xi$), Boltzmann constant $k_B = 8.617333\times10^{-5}$ eV/K.

Numerical Method

The tool solves the 1D continuity equation with an explicit finite-volume scheme, offloaded to a background Web Worker so the UI never blocks while a run is in progress. Space is divided into $N$ nodes $x_i = i\Delta x$ spanning the full leg length exactly ($x_0=0$, $x_{N-1}=L$).

Discrete Update

Interior nodes update from the flux imbalance across a full cell of width $\Delta x$:

The two boundary nodes ($i=0$ and $i=N-1$) sit exactly on the domain edge, so they only own a half cell — node 0 owns $[0,\Delta x/2]$, not $[-\Delta x/2, \Delta x/2]$, since there's no leg on the other side. With the zero-flux condition $J(0)=J(L)=0$ imposed at the true edges, that gives:

The factor of 2 (rather than the interior nodes' implicit factor of 1) is what makes the scheme conserve the true integrated mass $\int_0^L C\,dx$ exactly, to floating-point precision, rather than a quantity that implicitly treats the leg as one grid spacing longer than entered. This was verified by direct calculation: without the factor of 2, the boundary cells behave as if they were full-width, so the total simulated mass drifts by $\mathcal{O}(1/N)$ — about 1–2% at the default 70-point grid, rising to about 10% at the coarsest 10-point floor. With it, ∫C dx is conserved exactly regardless of grid resolution or how strong the Soret term is.

Upwind Advection

Each face flux $J_{i+1/2}$ combines a centered Fickian term with an upwinded Soret (advective) term — the concentration used for the advective part comes from whichever side the local drift is flowing from, not a simple average of both sides:

$C_{upwind}=C_i$ if $v_{i+1/2}\ge0$ (drift toward $+x$), else $C_{i+1}$. A centered average for the Soret term is only oscillation-free while the local cell Péclet number $\big(\frac{Q^*}{k_BT^2}\frac{\partial T}{\partial x}\big)\Delta x$ stays below ≈2; grid coarsening (below) enlarges $\Delta x$ exactly when that number is already large from a big $Q^*$ or steep gradient, so a centered scheme can and did produce negative concentrations in that regime. Upwinding trades a bit of numerical smearing for guaranteed monotonicity. A $\max(0,\cdot)$ floor is also applied to every update as a belt-and-suspenders backstop.

Adaptive Grid & Stability

• Adaptive grid: The spatial grid starts at $N=70$ points and is coarsened in steps of 5 (down to a floor of $N=10$) if the step budget (4,000,000 steps) would otherwise be exceeded — this keeps very fast dopants tractable without changing the physics. If the budget still can't cover the full requested anneal duration even at the coarsest grid, the run is truncated and a banner (plus a note in the CSV export) says so along with how much simulated time you actually got — it is never left unstated.
• CFL stability: The time step obeys both a diffusive and an advective limit, and the smaller of the two is used: $\Delta t = \min(\Delta t_{diff}, \Delta t_{adv})$, with
  $\Delta t_{diff} = 0.4\,\dfrac{\Delta x^2}{D_{max}} \qquad \Delta t_{adv} = 0.25\,\dfrac{\Delta x}{v_{max}}$
  The marginal (exactly-stable) limits are $\Delta x^2/(2D)$ and $\Delta x/(2v)$ respectively — both figures above already include the boundary nodes' tighter half-cell requirement, so the 0.4 and 0.25 factors are true safety margins (not just margins for the interior nodes), not marginal values themselves.
• Snapshot sampling: Roughly 40 evenly spaced snapshots are stored across the run for the 3D Timeline view and CSV export; the 2D Final view compares only the first (t=0) and last snapshot.
• Boundary conditions: Both ends of the leg are closed (zero-flux, $J=0$) — dopant mass is neither absorbed nor emitted at the contacts, and (per the half-cell treatment above) the true integrated mass is conserved exactly, not approximately.

Limits & Validation Rules

Every hard bound the tool enforces, and why. Bounds on T_Hot, T_Cold, and Leg Length are enforced by clamping the input field itself (you'll see the number in the box change); the gradient cap is enforced by refusing to run and explaining why, never by silently changing a number you typed — the same is true of the Smart Input concentration formula (invalid syntax blocks the run with an inline error rather than falling back to something else).

Modeling Assumptions

Units & Symbols

FAQ & Troubleshooting