Nuclear Reactor

AI authored to showcase ironpad capabilities.

Introduction

This notebook simulates a 2D cross-section of a nuclear reactor using one-group neutron diffusion theory with six groups of delayed-neutron precursors, live in your browser: a Rust cell compiled to WebAssembly steps a 300x300 flux grid 450 times per second while you drive the control rods.

Physics

One-group neutron diffusion equation (Stacey, Nuclear Reactor Physics, §5.1):

\frac{\partial \phi}{\partial t} = v \left[ D \nabla^2 \phi - \Sigma_a \phi + (1 - \beta)\, \nu \Sigma_f \phi + \sum_{i=1}^{6} \lambda_i C_i \right] + S

where:

  • \phi(\mathbf{r}, t): scalar neutron flux (neutrons / cm² / s)
  • v: average neutron speed (cm / s)
  • D: diffusion coefficient (cm); how far neutrons scatter between collisions
  • \nabla^2 \phi: spatial Laplacian of the flux (leakage term)
  • \Sigma_a: macroscopic absorption cross-section (1 / cm)
  • \nu \Sigma_f: fission neutron production rate (1 / cm); \nu is neutrons per fission
  • \beta = \sum_i \beta_i: total delayed-neutron fraction (≈ 0.0065 for {}^{235}\text{U})
  • S: external neutron source strength (neutrons / cm³ / s)

Read the bracket as a ledger: diffusion spreads neutrons, absorption removes them, prompt fission replaces (1 - \beta) of each generation immediately, and the \lambda_i C_i sum trickles in from precursor decay. That last term looks like a rounding error. It is the entire reason this machine is operable.

Delayed-neutron precursor kinetics (Keepin, Physics of Nuclear Kinetics, Ch. 4):

\frac{\partial C_i}{\partial t} = \beta_i \, \nu \Sigma_f \, \phi - \lambda_i \, C_i \qquad (i = 1 \ldots 6)

where:

  • C_i(\mathbf{r}, t): concentration of delayed-neutron precursor group i
  • \beta_i: delayed-neutron fraction for group i
  • \lambda_i: decay constant for group i (1 / s); how quickly precursors emit their delayed neutrons

The effective multiplication factor is computed from the material balance:

k_\text{eff} = \frac{\nu \Sigma_f \phi}{\Sigma_a \phi + \text{leakage}}

When k_\text{eff} = 1 the reactor is critical (steady-state); k < 1 is subcritical (flux decays); k > 1 is supercritical (flux grows).

Before you drive the reactor below, the next sections explain why it can be driven at all. They collapse the field equation to point kinetics, quantify the delayed-neutron timescale gap, and price reactivity the way operators do, in dollars and cents.

From a field to a point

You do not need the full PDE to understand reactor control. Integrate the flux over the whole core and what remains is point kinetics, the workhorse model of reactor dynamics:

\frac{dn}{dt} = \frac{\rho - \beta}{\Lambda}\, n + \sum_{i=1}^{6} \lambda_i c_i \qquad\qquad \frac{dc_i}{dt} = \frac{\beta_i}{\Lambda}\, n - \lambda_i c_i

Narrated plainly:

  • n(t): total neutron population; \rho = (k - 1)/k: reactivity, distance from perfect balance.
  • \Lambda: mean generation time, birth to next fission. Roughly 25 microseconds in a water-moderated core.
  • First term: the prompt ledger. Every \Lambda seconds the population multiplies by (\rho - \beta) worth of imbalance.
  • Second term: delayed neutrons arriving from precursor decay; the c_i equation is shelf bookkeeping (fission stocks the shelf, decay clears it).

That (\rho - \beta) grouping is the most load-bearing notation in reactor physics. While \rho < \beta, prompt neutrons alone form a dying chain; the reactor grows only as fast as the precursor shelf restocks, in seconds. The next cell shows the signature two-phase response: a fast prompt jump, then a slow, steerable climb.

Prompt Jump

The 0.65% that makes reactors possible

Delayed neutrons are the reason a fission reactor is a machine humans can operate. The numbers make the case:

  • A prompt neutron finishes its entire life in about 25 microseconds.
  • A precursor sits on the shelf (half-lives 0.23 s to 55.9 s, weighted average about 13 seconds) before releasing its neutron.
  • Fold that 0.65% into the average and the effective generation time stretches from 25 μs to roughly 0.085 s: a factor of about 3,400.

For a mental model: prompt-only at k = 1.001 has an e-folding period of \Lambda / \delta k = 25 milliseconds; power multiplies by e^{40} \approx 10^{17} every second. No mechanical system moves rods on that timescale. With delayed neutrons pacing it, the same k gives a period near 70 seconds: steerable with an electric motor. The prompt-only case is not worth running live, so it stays in a fence:

// Prompt-only at k = 1.001: dn/dt = (0.001 / 25e-6) * n = 40 n.
// n(t) = e^(40 t): about 5e8 by t = 0.5 s, f64 overflow by t = 17.7 s.
let mut n = 1.0_f64;
loop { n *= 1.00004; } // 1 microsecond steps, hopeless within one frame

The next cell computes these numbers straight from Keepin's six-group data. They are directional rather than precise (\Lambda varies by core design), but the ratio is the point.

Timescale Gap

Reactivity in dollars and cents

Reactivity has natural units, and they are money. Define one dollar = \beta: exactly the reactivity at which prompt neutrons alone sustain the chain. A hundredth of a dollar is a cent; small \rho is also quoted in pcm (10^{-5}). The ladder:

  • \rho < 0: subcritical; the chain dies.
  • 0 < \rho < \beta (under one dollar): delayed supercritical. Precursors set the pace; the operator is in charge.
  • \rho \geq \beta (one dollar and up): prompt critical. The 3,400x cushion evaporates and power rises on the microsecond ledger. This is the cliff; SL-1 and Chernobyl both went over it.

Operators quote rod moves in cents because the dollar scale pins the cliff at a fixed place regardless of fuel; \beta is about 0.0065 for U-235 but only 0.0021 for Pu-239, one reason plutonium-heavy cores handle more nervously. The next cell converts k_\text{eff}, the number on the live gauge below, into pcm, dollars, and cents.

Reactivity Ladder

Control-rod worth

A rod's worth is the reactivity change it causes, and it is nowhere near linear in insertion depth. First-order perturbation theory says a small absorber steals reactivity proportional to the square of the flux at its location (flux weights both the neutrons available and the importance of losing them). With a roughly cosine axial flux:

  • Differential worth (reactivity per unit of travel) follows \sin^2(\pi x): nearly zero at the ends, peaked at mid-core.
  • Integral worth follows the classic S-curve W(x) = x - \sin(2\pi x)/(2\pi).

The numbers are lopsided: the first 10% of insertion buys about 0.6% of total worth; the middle 20% buys about 39%. That is why fine power control happens with rods near mid-core. Watch for this below: the k_\text{eff} gauge barely reacts to the first slider increments and swings hardest through mid-range.

Rod Worth S-Curve

What this simulation leaves out

A model should state its own error bars. Four simplifications here matter in a real core:

  • One energy group. Real neutrons are born near 2 MeV and thermalize to 0.025 eV; cross-sections swing by orders of magnitude along the way. Production codes use 2 to hundreds of groups.
  • No xenon. Xenon-135 is a fission-product poison with a capture cross-section near 2.6 million barns, building and decaying on a 6-to-9-hour timescale; it is why real reactors fall into a "xenon pit" after a power drop. Our clock runs for seconds; xenon never shows up.
  • No thermal feedback. Rising fuel temperature Doppler-broadens U-238 resonances and promptly pushes reactivity down; moderator density feedback does the same. Those negative coefficients make water reactors forgiving. Here the feedback loop is you and the sliders.
  • Assorted pedagogy. A 2D slice of a 3D core; neutron speed scaled from 2.2 × 10⁵ cm/s to 50 cm/s so your eyes can follow the flux; explicit Euler; rod worths tuned for feel, not calibrated.

None of this changes the story: \beta is still 0.0065, precursors still take seconds, and the response to your rod moves is qualitatively the real thing.

What You're Seeing

ElementDescription
Plasma coreFuel region; colour intensity shows neutron flux
Teal ringWater/graphite reflector surrounding the core
Silver barsControl rod blades (7 total); they absorb neutrons
Left gaugek_\text{eff} (white line = criticality at k = 1)
Right gaugeTotal neutron flux relative to peak
Bottom plotsk_\text{eff} and total flux time histories

Controls

Use the sliders to adjust:

  • Control Rod Depth: how far the rods are inserted (1.0 = fully in)
  • Rod Absorption (\Sigma_a): neutron capture strength of the rod material
  • Fission XS (\nu\Sigma_f): fission production rate of the fuel
  • Fuel Absorption (\Sigma_a): parasitic neutron capture in fuel
  • Reflector Diffusion (D): how easily neutrons scatter through the reflector
  • Reflector Absorption (\Sigma_a): neutron capture rate in the reflector
  • Neutron Source: toggle the external startup source on/off

Try it: Set k_\text{eff} < 1 by inserting the rods fully, then turn off the neutron source and watch the flux decay exponentially to zero. Then pull the rods out slowly and watch subcritical multiplication stretch the settling time as k creeps toward 1.

Reactor Controls
Reactor Simulation
Instruments
Dashboard