diff --git a/Physlib/Electromagnetism/ThreeDimension/MaxwellEquations.lean b/Physlib/Electromagnetism/ThreeDimension/MaxwellEquations.lean
index f2b4c12b9c..de0927556b 100644
--- a/Physlib/Electromagnetism/ThreeDimension/MaxwellEquations.lean
+++ b/Physlib/Electromagnetism/ThreeDimension/MaxwellEquations.lean
@@ -40,6 +40,9 @@ boundary conditions, and constitutive laws for material media are outside the cu
 module.
 
 -/
+
+@[expose] public section
+
 namespace Electromagnetism
 namespace ThreeDimension
 
diff --git a/PhyslibAlpha.lean b/PhyslibAlpha.lean
index 2ab231a9a7..97d7d7fd9e 100644
--- a/PhyslibAlpha.lean
+++ b/PhyslibAlpha.lean
@@ -37,6 +37,8 @@ public import PhyslibAlpha.CondensedMatter.TightBindingChain.Current
 public import PhyslibAlpha.CondensedMatter.TightBindingChain.CurrentEigenstates
 public import PhyslibAlpha.CondensedMatter.TightBindingChain.OpenBoundary
 public import PhyslibAlpha.CondensedMatter.TightBindingChain.Uncertainty
+public import PhyslibAlpha.Electromagnetism.BoxChargeConservation
+public import PhyslibAlpha.Electromagnetism.Distributional.WireJunction
 public import PhyslibAlpha.Mathematics.Analysis.Normed.HolderDual
 public import PhyslibAlpha.Mathematics.Analysis.RealBounds
 public import PhyslibAlpha.Mathematics.Convex.Choquet.BoundaryRepresentation
diff --git a/PhyslibAlpha/Electromagnetism/BoxChargeConservation.lean b/PhyslibAlpha/Electromagnetism/BoxChargeConservation.lean
new file mode 100644
index 0000000000..346721f153
--- /dev/null
+++ b/PhyslibAlpha/Electromagnetism/BoxChargeConservation.lean
@@ -0,0 +1,279 @@
+/-
+Copyright (c) 2026 Joseph Tooby-Smith. All rights reserved.
+Released under Apache 2.0 license as described in the file LICENSE.
+Authors: Joseph Tooby-Smith
+-/
+module
+
+public import Physlib.Electromagnetism.ThreeDimension.MaxwellEquations
+public import Mathlib.MeasureTheory.Integral.DivergenceTheorem
+
+/-!
+
+# Charge conservation on a box from Maxwell's equations
+
+Let `V` be an electromagnetic potential and `J₄` a Lorentz current density such that `V` is an
+extremum of the free-space action with source `J₄`, i.e. Maxwell's equations hold. Taking the
+divergence of Ampère's law, the divergence of the curl of the magnetic field vanishes, and the
+divergence of the displacement current is, by Gauss's law, the time derivative of the charge
+density. This gives the continuity equation `∂ₜ ρ + ∇ ⬝ J = 0`.
+
+Integrating the continuity equation over a coordinate box `[a, b]` and applying the divergence
+theorem, the net current leaving the box equals minus the rate of change of the charge it
+contains. In particular, if no charge accumulates inside the box, the currents leaving through
+its six faces sum to zero.
+
+This is motivated by Kirchhoff's current law, but is not that law: Kirchhoff's current law is a
+statement about the branch currents meeting at the nodes of a lumped-element circuit, whereas
+the results here are the integral form of charge conservation on a single box.
+
+## Main results
+
+- `Space.integral_div_box` : the divergence theorem on a coordinate box in `Space`.
+- `Electromagnetism.ThreeDimension.continuityEquation` : the continuity equation.
+- `Electromagnetism.LorentzCurrentDensity.boxFaceCurrent`,
+  `Electromagnetism.LorentzCurrentDensity.boxOutwardCurrent` : the current through a face of a
+  coordinate box, and the net current leaving the box.
+- `Electromagnetism.ThreeDimension.boxOutwardCurrent_eq` : the integral form of charge
+  conservation.
+- `Electromagnetism.ThreeDimension.boxOutwardCurrent_eq_zero_of_steady` : the net current
+  leaving a box in which no charge accumulates is zero.
+
+## Contents
+
+- A. Vector calculus on `Space`
+  - A.1. Time derivatives and the divergence
+  - A.2. Smoothness of the curl
+  - A.3. The divergence theorem on a coordinate box
+- B. The continuity equation
+  - B.1. The continuity equation from Gauss's and Ampère's laws
+  - B.2. The continuity equation for an electromagnetic potential
+- C. Charge conservation on a box
+  - C.1. Currents through a box
+  - C.2. The integral form of charge conservation
+  - C.3. Steady charge in a box
+
+-/
+
+@[expose] public section
+
+open Space Time MeasureTheory Set
+
+namespace Space
+
+/-!
+
+## A. Vector calculus on `Space`
+
+### A.1. Time derivatives and the divergence
+
+-/
+
+lemma div_differentiable_time {d} (f : Time → Space d → EuclideanSpace ℝ (Fin d))
+    (hf : ContDiff ℝ 2 ↿f) (x : Space d) :
+    Differentiable ℝ (fun t => (∇ ⬝ f t) x) := by
+  have hfi (i : Fin d) : ContDiff ℝ 2 ↿(fun t x => f t x i) :=
+    (ContinuousLinearMap.contDiff (𝕜 := ℝ) (EuclideanSpace.proj i)).comp hf
+  simp only [div]
+  exact Differentiable.fun_sum fun i _ => space_deriv_differentiable_time (hfi i) x
+
+/-- The divergence and the time derivative commute. -/
+lemma time_deriv_div_commute {d} (f : Time → Space d → EuclideanSpace ℝ (Fin d))
+    (hf : ContDiff ℝ 2 ↿f) (t : Time) (x : Space d) :
+    ∂ₜ (fun t => (∇ ⬝ f t) x) t = (∇ ⬝ fun x => ∂ₜ (fun t => f t x) t) x := by
+  have hfi (i : Fin d) : ContDiff ℝ 2 ↿(fun t x => f t x i) :=
+    (ContinuousLinearMap.contDiff (𝕜 := ℝ) (EuclideanSpace.proj i)).comp hf
+  simp only [div]
+  rw [Time.deriv_eq, fderiv_fun_sum fun i _ =>
+    (space_deriv_differentiable_time (hfi i) x).differentiableAt]
+  simp only [FunLike.coe_sum, Finset.sum_apply, ← Time.deriv_eq]
+  refine Finset.sum_congr rfl fun i _ => ?_
+  rw [time_deriv_comm_space_deriv (hfi i)]
+  congr
+  funext y
+  rw [Time.deriv_euclid]
+  exact fun t => ((hf.differentiable (by simp)).comp (f := fun t => (t, y)) (by fun_prop)) t
+
+/-!
+
+### A.2. Smoothness of the curl
+
+-/
+
+lemma curl_contDiff {n : WithTop ℕ∞} (f : Space → EuclideanSpace ℝ (Fin 3))
+    (hf : ContDiff ℝ (n + 1) f) : ContDiff ℝ n (∇ ⨯ f) := by
+  rw [contDiff_euclidean]
+  intro i
+  have hfi (j : Fin 3) : ContDiff ℝ (n + 1) (fun x => f x j) :=
+    (ContinuousLinearMap.contDiff (𝕜 := ℝ) (EuclideanSpace.proj j)).comp hf
+  have hd (j k : Fin 3) : ContDiff ℝ n (fun x => ∂[j] (fun x => f x k) x) :=
+    (contDiff_apply ℝ ℝ j).comp (deriv_contDiff (hfi k))
+  exact (hd _ _).sub (hd _ _)
+
+/-!
+
+### A.3. The divergence theorem on a coordinate box
+
+The box `[a, b]` is given in coordinates, `a b : Fin (n + 1) → ℝ`. Its faces normal to the
+`i`-th axis are parametrized, as in Mathlib's divergence theorem, by the remaining coordinates
+`y ∈ [a ∘ i.succAbove, b ∘ i.succAbove]` through `i.insertNth c y` with `c = a i` or `c = b i`.
+
+The divergence theorem below integrates over the box in coordinates, `x ∈ [a, b]`, rather than
+over the image of the box in `Space`.
+
+-/
+
+lemma div_mk_eq_sum_fderiv {n : ℕ} (f : Space n → EuclideanSpace ℝ (Fin n))
+    (hf : Differentiable ℝ f) (x : Fin n → ℝ) :
+    (∇ ⬝ f) ⟨x⟩ = ∑ i, fderiv ℝ (fun y : Fin n → ℝ => f ⟨y⟩ i) x (Pi.single i 1) := by
+  refine Finset.sum_congr rfl fun i _ => ?_
+  have hfi : Differentiable ℝ (fun y : Space n => f y i) :=
+    (EuclideanSpace.proj (𝕜 := ℝ) i).differentiable.comp hf
+  rw [Space.deriv_eq_fderiv_basis, fderiv_fun_comp x (hfi _) (mk_differentiable x), fderiv_mk]
+  simp only [ContinuousLinearMap.coe_comp, Function.comp_apply]
+  congr 1
+  ext j
+  simp [equivPi, basis_apply, Pi.single_apply, eq_comm]
+
+/-- The divergence theorem on a coordinate box in `Space`. -/
+lemma integral_div_box {n : ℕ} (f : Space (n + 1) → EuclideanSpace ℝ (Fin (n + 1)))
+    (hf : ContDiff ℝ 1 f) (a b : Fin (n + 1) → ℝ) (hle : a ≤ b) :
+    ∫ x in Icc a b, (∇ ⬝ f) ⟨x⟩ =
+      ∑ i, ((∫ y in Icc (a ∘ i.succAbove) (b ∘ i.succAbove), f ⟨i.insertNth (b i) y⟩ i) -
+        ∫ y in Icc (a ∘ i.succAbove) (b ∘ i.succAbove), f ⟨i.insertNth (a i) y⟩ i) := by
+  have hfi (i : Fin (n + 1)) : ContDiff ℝ 1 (fun y : Fin (n + 1) → ℝ => f ⟨y⟩ i) :=
+    (EuclideanSpace.proj (𝕜 := ℝ) i).contDiff.comp (hf.comp mk_contDiff)
+  simp only [div_mk_eq_sum_fderiv f (hf.differentiable (by simp))]
+  refine integral_divergence_of_hasFDerivAt_off_countable' a b hle (fun i y => f ⟨y⟩ i)
+    (fun i y => fderiv ℝ (fun y : Fin (n + 1) → ℝ => f ⟨y⟩ i) y) ∅ countable_empty
+    (fun i => (hfi i).continuous.continuousOn)
+    (fun y _ i => ((hfi i).differentiable (by simp) y).hasFDerivAt) ?_
+  refine ContinuousOn.integrableOn_compact isCompact_Icc (Continuous.continuousOn ?_)
+  exact continuous_finsetSum _ fun i _ =>
+    ((hfi i).continuous_fderiv (by simp)).clm_apply continuous_const
+
+end Space
+
+namespace Electromagnetism
+namespace ThreeDimension
+
+open ElectromagneticPotential ContDiff
+
+/-!
+
+## B. The continuity equation
+
+### B.1. The continuity equation from Gauss's and Ampère's laws
+
+-/
+
+/-- Fields `E`, `B`, `ρ` and `J` obeying Gauss's law for the electric field and Ampère's law
+obey the continuity equation. -/
+lemma continuity_of_gauss_ampere (𝓕 : FreeSpace)
+    {E B J : Time → Space → EuclideanSpace ℝ (Fin 3)} {ρ : Time → Space → ℝ}
+    (hE : ContDiff ℝ 2 ↿E) (hB : ∀ t, ContDiff ℝ 2 (B t))
+    (hGauss : ∀ t x, (∇ ⬝ E t) x = ρ t x / 𝓕.ε₀)
+    (hAmpere : ∀ t x, (∇ ⨯ B t) x = 𝓕.μ₀ • J t x + 𝓕.μ₀ • 𝓕.ε₀ • ∂ₜ (fun t => E t x) t)
+    (t : Time) (x : Space) :
+    ∂ₜ (fun t => ρ t x) t + (∇ ⬝ J t) x = 0 := by
+  -- Solve Ampère's law for `J` and Gauss's law for `ρ`.
+  have hJ : J t = 𝓕.μ₀⁻¹ • (∇ ⨯ B t) - 𝓕.ε₀ • fun x => ∂ₜ (fun t => E t x) t := by
+    funext y
+    rw [Pi.sub_apply, Pi.smul_apply, Pi.smul_apply, hAmpere, smul_add, smul_smul, smul_smul,
+      inv_mul_cancel₀ 𝓕.μ₀_ne_zero, one_smul, one_smul, add_sub_cancel_right]
+  have hρ : (fun t => ρ t x) = fun t => 𝓕.ε₀ * (∇ ⬝ E t) x := by
+    funext s
+    rw [hGauss, mul_div_cancel₀ _ 𝓕.ε₀_ne_zero]
+  have hcurl : Differentiable ℝ (∇ ⨯ B t) :=
+    (Space.curl_contDiff (n := 1) _ (hB t)).differentiable (by simp)
+  have hdE : Differentiable ℝ fun x => ∂ₜ (fun t => E t x) t :=
+    time_deriv_differentiable_space hE t
+  -- The divergence of the curl vanishes; the displacement current carries `-∂ₜ ρ`.
+  rw [hJ, sub_eq_add_neg, ← neg_smul, div_add _ _ (hcurl.const_smul _) (hdE.const_smul _),
+    div_smul _ _ hcurl, div_smul _ _ hdE, div_of_curl_eq_zero _ (hB t), hρ, Time.deriv_eq,
+    fderiv_const_mul ((div_differentiable_time E hE x) t)]
+  simp [← Time.deriv_eq, time_deriv_div_commute E hE]
+
+/-!
+
+### B.2. The continuity equation for an electromagnetic potential
+
+-/
+
+variable {𝓕 : FreeSpace} (V : ElectromagneticPotential 3) (J₄ : LorentzCurrentDensity 3)
+
+lemma magneticField_contDiff_space (hV : ContDiff ℝ ∞ V) (t : Time) :
+    ContDiff ℝ 2 (V.magneticField 𝓕.c t) := by
+  have hV3 : ContDiff ℝ (2 + 1) V := hV.of_le (WithTop.coe_le_coe.mpr le_top)
+  simp only [magneticField_eq_3D]
+  exact Space.curl_contDiff (n := 2) _ (by fun_prop)
+
+/-- The continuity equation: Maxwell's equations imply local conservation of charge. -/
+theorem continuityEquation (t : Time) (x : Space)
+    (h : IsExtrema 𝓕 V J₄) (hV : ContDiff ℝ ∞ V) (hJ : ContDiff ℝ ∞ J₄) :
+    ∂ₜ (fun t => J₄.chargeDensity 𝓕.c t x) t + (∇ ⬝ J₄.currentDensity 𝓕.c t) x = 0 :=
+  continuity_of_gauss_ampere 𝓕 (electricField_contDiff (hV.of_le (WithTop.coe_le_coe.mpr le_top)))
+    (magneticField_contDiff_space V hV) (fun t x => gaussLawElectric V J₄ t x h hV hJ)
+    (fun t x => ampereLaw V J₄ t x h hV hJ) t x
+
+/-!
+
+## C. Charge conservation on a box
+
+### C.1. Currents through a box
+
+-/
+
+/-- The electric current at time `t` through the face `x i = s` of the coordinate box `[a, b]`,
+counted positively in the direction of increasing `x i`. -/
+noncomputable def _root_.Electromagnetism.LorentzCurrentDensity.boxFaceCurrent
+    (c : SpeedOfLight) (J : LorentzCurrentDensity 3) (t : Time) (a b : Fin 3 → ℝ) (i : Fin 3)
+    (s : ℝ) : ℝ :=
+  ∫ y in Icc (a ∘ i.succAbove) (b ∘ i.succAbove), J.currentDensity c t ⟨i.insertNth s y⟩ i
+
+/-- The net electric current at time `t` leaving the coordinate box `[a, b]` through its six
+faces. -/
+noncomputable def _root_.Electromagnetism.LorentzCurrentDensity.boxOutwardCurrent
+    (c : SpeedOfLight) (J : LorentzCurrentDensity 3) (t : Time) (a b : Fin 3 → ℝ) : ℝ :=
+  ∑ i, (J.boxFaceCurrent c t a b i (b i) - J.boxFaceCurrent c t a b i (a i))
+
+/-!
+
+### C.2. The integral form of charge conservation
+
+-/
+
+/-- The net current leaving the coordinate box `[a, b]` equals minus the rate of change of the
+charge inside it. -/
+lemma boxOutwardCurrent_eq (t : Time) (a b : Fin 3 → ℝ) (hle : a ≤ b)
+    (h : IsExtrema 𝓕 V J₄) (hV : ContDiff ℝ ∞ V) (hJ : ContDiff ℝ ∞ J₄) :
+    J₄.boxOutwardCurrent 𝓕.c t a b =
+      - ∫ x in Icc a b, ∂ₜ (fun t => J₄.chargeDensity 𝓕.c t ⟨x⟩) t := by
+  have hJt : ContDiff ℝ 1 (J₄.currentDensity 𝓕.c t) :=
+    (LorentzCurrentDensity.currentDensity_ContDiff
+      (hJ.of_le (WithTop.coe_le_coe.mpr le_top))).comp (f := fun x => (t, x)) (by fun_prop)
+  rw [LorentzCurrentDensity.boxOutwardCurrent]
+  simp only [LorentzCurrentDensity.boxFaceCurrent]
+  rw [← integral_div_box _ hJt a b hle, ← integral_neg]
+  congr 1
+  funext x
+  exact eq_neg_of_add_eq_zero_right (continuityEquation V J₄ t ⟨x⟩ h hV hJ)
+
+/-!
+
+### C.3. Steady charge in a box
+
+-/
+
+/-- If no charge accumulates inside the coordinate box `[a, b]`, the net current leaving through
+its six faces is zero. This is motivated by Kirchhoff's current law, with the box enclosing a
+node of a circuit. -/
+lemma boxOutwardCurrent_eq_zero_of_steady (t : Time) (a b : Fin 3 → ℝ) (hle : a ≤ b)
+    (h : IsExtrema 𝓕 V J₄) (hV : ContDiff ℝ ∞ V) (hJ : ContDiff ℝ ∞ J₄)
+    (hsteady : ∀ x ∈ Icc a b, ∂ₜ (fun t => J₄.chargeDensity 𝓕.c t ⟨x⟩) t = 0) :
+    J₄.boxOutwardCurrent 𝓕.c t a b = 0 := by
+  rw [boxOutwardCurrent_eq V J₄ t a b hle h hV hJ, setIntegral_eq_zero_of_forall_eq_zero hsteady,
+    neg_zero]
+
+end ThreeDimension
+end Electromagnetism
diff --git a/PhyslibAlpha/Electromagnetism/Distributional/WireJunction.lean b/PhyslibAlpha/Electromagnetism/Distributional/WireJunction.lean
new file mode 100644
index 0000000000..96ed13f01b
--- /dev/null
+++ b/PhyslibAlpha/Electromagnetism/Distributional/WireJunction.lean
@@ -0,0 +1,270 @@
+/-
+Copyright (c) 2026 Joseph Tooby-Smith. All rights reserved.
+Released under Apache 2.0 license as described in the file LICENSE.
+Authors: Joseph Tooby-Smith
+-/
+module
+
+public import Physlib.Electromagnetism.Distributional.Dynamics.IsExtrema
+public import Physlib.SpaceAndTime.TimeAndSpace.ConstantTimeDist
+
+/-!
+
+# Junctions of thin wires
+
+Let `A` be a distributional electromagnetic potential and `J` a distributional Lorentz current
+density such that Maxwell's equations hold, `IsExtrema 𝓕 A J`. Taking the divergence of
+Ampère's law and using Gauss's law gives the continuity equation `∂ₜ ρ + ∇ ⬝ J = 0` as an
+equation of distributions.
+
+Distributions allow currents to be carried by thin wires. A junction is modelled by finitely
+many straight, semi-infinite wires leaving the origin in the directions `u k`, carrying steady
+currents `I k`. The divergence of such a current density is `(∑ k, I k) δ₀`: the wires end at
+the origin, where the currents must come from somewhere. If no charge accumulates, the
+continuity equation forces `∑ k, I k = 0`.
+
+This is motivated by Kirchhoff's current law, of which it is a special case: a single node with
+steady currents, rather than the branch currents at every node of a lumped-element circuit.
+
+## Main results
+
+- `Electromagnetism.DistElectromagneticPotential.continuityEquation` : the continuity equation
+  for distributions.
+- `Electromagnetism.DistLorentzCurrentDensity.IsWireJunction` : the current density of a
+  junction of thin wires.
+- `Electromagnetism.DistLorentzCurrentDensity.IsWireJunction.distSpaceDiv_currentDensity` : the
+  divergence of the current density of a junction.
+- `Electromagnetism.DistLorentzCurrentDensity.IsWireJunction.sum_currents_eq_zero` : the
+  currents leaving a junction of thin wires sum to zero.
+
+## Contents
+
+- A. Analysis on `Space`
+  - A.1. Time derivatives and the divergence of distributions
+  - A.2. Schwartz functions along a ray
+  - A.3. Integrating test functions over time
+- B. The continuity equation
+- C. Junctions of thin wires
+  - C.1. The current density of a junction
+  - C.2. The currents at a junction sum to zero
+
+-/
+
+@[expose] public section
+
+namespace Space
+
+open MeasureTheory Set Filter SchwartzMap
+
+/-!
+
+## A. Analysis on `Space`
+
+### A.1. Time derivatives and the divergence of distributions
+
+-/
+
+/-- The divergence and the time derivative of distributions commute. -/
+lemma distTimeDeriv_distSpaceDiv {d} (f : (Time × Space d) →d[ℝ] EuclideanSpace ℝ (Fin d)) :
+    distTimeDeriv (distSpaceDiv f) = distSpaceDiv (distTimeDeriv f) := by
+  ext ε
+  rw [distTimeDeriv_apply', distSpaceDiv_apply_eq_sum_distSpaceDeriv,
+    distSpaceDiv_apply_eq_sum_distSpaceDeriv]
+  simp [apply_fderiv_eq_distTimeDeriv, ← distTimeDeriv_commute_distSpaceDeriv]
+
+/-!
+
+### A.2. Schwartz functions along a ray
+
+The ray from the origin in the direction `w` is parametrized by `s ↦ s • w`.
+
+-/
+
+/-- The restriction of a Schwartz function to a line through the origin is a Schwartz
+function. -/
+lemma exists_schwartz_eq_comp_ray {d} (χ : 𝓢(Space d, ℝ)) {w : Space d} (hw : w ≠ 0) :
+    ∃ g : 𝓢(ℝ, ℝ), ∀ s : ℝ, g s = χ (s • w) := by
+  let γ : ℝ →L[ℝ] Space d := (ContinuousLinearMap.id ℝ ℝ).smulRight w
+  have hγ : AntilipschitzWith ‖w‖₊⁻¹ γ := γ.antilipschitz_of_bound fun s => by
+    simp only [γ, ContinuousLinearMap.smulRight_apply, ContinuousLinearMap.id_apply, norm_smul]
+    rw [NNReal.coe_inv, coe_nnnorm, mul_comm ‖s‖, ← mul_assoc,
+      inv_mul_cancel₀ (norm_ne_zero_iff.mpr hw), one_mul]
+  exact ⟨compCLMOfAntilipschitz ℝ γ.hasTemperateGrowth hγ χ, fun s => rfl⟩
+
+/-- The fundamental theorem of calculus along the ray in direction `w`. -/
+lemma integral_Ioi_fderiv_ray {d} (ψ : 𝓢(Space d, ℝ)) {w : Space d} (hw : w ≠ 0) :
+    ∫ s in Ioi (0 : ℝ), fderiv ℝ ψ (s • w) w = - ψ 0 := by
+  obtain ⟨g, hg⟩ := exists_schwartz_eq_comp_ray ψ hw
+  obtain ⟨g', hg'⟩ := exists_schwartz_eq_comp_ray
+    (SchwartzMap.evalCLM ℝ (Space d) ℝ w (fderivCLM ℝ (Space d) ℝ ψ)) hw
+  have h := integral_Ioi_of_hasDerivAt_of_tendsto (f := g) (f' := g') (m := 0) (a := 0)
+    g.continuous.continuousWithinAt (fun s _ => ?_) g'.integrable.integrableOn
+    (g.tendsto_cocompact.mono_left atTop_le_cocompact)
+  · simpa [hg, hg'] using h
+  · have hd := (ψ.differentiableAt (x := s • w)).hasFDerivAt.comp_hasDerivAt s
+      ((hasDerivAt_id s).smul_const w)
+    rw [show (⇑g) = fun s => ψ (s • w) from funext hg, hg']
+    rw [one_smul] at hd
+    exact hd
+
+lemma integral_Ioi_fderiv_ray_basis {d} (ψ : 𝓢(Space d, ℝ)) {u : EuclideanSpace ℝ (Fin d)}
+    (hu : u ≠ 0) :
+    ∑ i, (∫ s in Ioi (0 : ℝ), fderiv ℝ ψ (s • basis.repr.symm u) (basis i)) * u i = - ψ 0 := by
+  have hw : basis.repr.symm u ≠ 0 := by simpa using hu
+  have hint (i : Fin d) : Integrable fun s : ℝ =>
+      fderiv ℝ ψ (s • basis.repr.symm u) (basis i) * u i := by
+    obtain ⟨g, hg⟩ := exists_schwartz_eq_comp_ray
+      (SchwartzMap.evalCLM ℝ (Space d) ℝ (basis i) (fderivCLM ℝ (Space d) ℝ ψ)) hw
+    exact (g.integrable.mul_const (u i)).congr (Eventually.of_forall fun s => by simp [hg])
+  simp only [← integral_mul_const]
+  rw [← integral_finsetSum _ fun i _ => (hint i).integrableOn, ← integral_Ioi_fderiv_ray ψ hw]
+  congr 1
+  funext s
+  rw [show fderiv ℝ ψ (s • basis.repr.symm u) (basis.repr.symm u) =
+      fderiv ℝ ψ (s • basis.repr.symm u) (∑ i, u i • basis i) by rw [basis.sum_repr_symm]]
+  simp [mul_comm]
+
+/-!
+
+### A.3. Integrating test functions over time
+
+-/
+
+/-- Integrating over time commutes with differentiating in space. -/
+lemma timeIntegralSchwartz_fderiv_space {d} (η : 𝓢(Time × Space d, ℝ)) (i : Fin d)
+    (x : Space d) :
+    timeIntegralSchwartz (SchwartzMap.evalCLM ℝ (Time × Space d) ℝ (0, basis i)
+      (fderivCLM ℝ (Time × Space d) ℝ η)) x = fderiv ℝ (timeIntegralSchwartz η) x (basis i) := by
+  have h := congrArg (fun f => f η)
+    (constantTime_distSpaceDeriv i (Physlib.Distribution.diracDelta ℝ x))
+  simpa [distSpaceDeriv_apply', constantTime_apply, distDeriv_apply,
+    Physlib.Distribution.fderivD_apply] using h
+
+/-- Some test function has a non-zero time integral at the origin. -/
+lemma exists_integral_time_ne_zero {d} :
+    ∃ η : 𝓢(Time × Space d, ℝ), ∫ t, η (t, 0) ≠ 0 := by
+  let b : ContDiffBump (0 : Time × Space d) := ⟨1, 2, one_pos, one_lt_two⟩
+  refine ⟨b.hasCompactSupport.toSchwartzMap b.contDiff, ne_of_gt ?_⟩
+  have hsupp : HasCompactSupport fun t : Time => b (t, 0) :=
+    HasCompactSupport.intro (isCompact_closedBall (0 : Time) 2) fun t ht =>
+      b.zero_of_le_dist (by
+        rw [Metric.mem_closedBall, not_le, dist_zero_right] at ht
+        simpa [b, Prod.dist_eq] using ht.le)
+  refine Continuous.integral_pos_of_hasCompactSupport_nonneg_nonzero (x := 0)
+    (by fun_prop) hsupp (fun t => b.nonneg) (b.pos_of_mem_ball (by simp [b])).ne'
+
+end Space
+
+namespace Electromagnetism
+open Space SchwartzMap
+
+namespace DistElectromagneticPotential
+
+/-!
+
+## B. The continuity equation
+
+-/
+
+/-- The continuity equation for distributions: Maxwell's equations imply local conservation of
+charge. -/
+theorem continuityEquation {d} {𝓕 : FreeSpace} (A : DistElectromagneticPotential d)
+    (J : DistLorentzCurrentDensity d) (h : IsExtrema 𝓕 A J) (ε : 𝓢(Time × Space d, ℝ)) :
+    distTimeDeriv (J.chargeDensity 𝓕.c) ε + distSpaceDiv (J.currentDensity 𝓕.c) ε = 0 := by
+  obtain ⟨hG, hA⟩ := (isExtrema_iff_vectorPotential A J).mp h
+  set E := A.electricField 𝓕.c
+  set V := A.vectorPotential 𝓕.c
+  -- Ampère's law, differentiated in the `i`-th direction.
+  have hAi (i : Fin d) : 𝓕.μ₀ * 𝓕.ε₀ * distSpaceDeriv i (distTimeDeriv E) ε i +
+      (∑ x, distSpaceDeriv i (distSpaceDeriv x (distSpaceDeriv x V)) ε i
+        - ∑ x, distSpaceDeriv i (distSpaceDeriv x (distSpaceDeriv i V)) ε x) +
+      𝓕.μ₀ * distSpaceDeriv i (J.currentDensity 𝓕.c) ε i = 0 := by
+    have := hA ((SchwartzMap.evalCLM ℝ (Time × Space d) ℝ (0, basis i))
+      ((fderivCLM ℝ (Time × Space d) ℝ) ε)) i
+    simp only [apply_fderiv_eq_distSpaceDeriv] at this
+    simp only [PiLp.neg_apply, neg_sub_neg, Finset.sum_neg_distrib, Finset.sum_sub_distrib] at this
+    linear_combination -this
+  -- The curl-curl terms cancel once summed over `i`.
+  have hcurl : ∑ i, (∑ x, distSpaceDeriv i (distSpaceDeriv x (distSpaceDeriv x V)) ε i
+      - ∑ x, distSpaceDeriv i (distSpaceDeriv x (distSpaceDeriv i V)) ε x) = 0 := by
+    rw [Finset.sum_sub_distrib, sub_eq_zero, Finset.sum_comm]
+    refine Finset.sum_congr rfl fun i _ => Finset.sum_congr rfl fun x _ => ?_
+    rw [distSpaceDeriv_commute]
+  -- Gauss's law gives the time derivative of the charge density.
+  have hρ : J.chargeDensity 𝓕.c = 𝓕.ε₀ • distSpaceDiv E := by
+    ext η
+    rw [_root_.smul_apply, hG, smul_eq_mul, ← mul_assoc,
+      mul_one_div_cancel 𝓕.ε₀_ne_zero, one_mul]
+  have hsum := Finset.sum_eq_zero (s := Finset.univ) fun i (_ : i ∈ Finset.univ) => hAi i
+  rw [hρ, map_smul, _root_.smul_apply, distTimeDeriv_distSpaceDiv, smul_eq_mul,
+    distSpaceDiv_apply_eq_sum_distSpaceDeriv, distSpaceDiv_apply_eq_sum_distSpaceDeriv]
+  rw [Finset.sum_add_distrib, Finset.sum_add_distrib, ← Finset.mul_sum, ← Finset.mul_sum,
+    hcurl] at hsum
+  apply mul_left_cancel₀ 𝓕.μ₀_ne_zero
+  linear_combination hsum
+
+end DistElectromagneticPotential
+
+namespace DistLorentzCurrentDensity
+
+open MeasureTheory Set
+
+/-!
+
+## C. Junctions of thin wires
+
+### C.1. The current density of a junction
+
+The `k`-th wire is the ray `s ↦ s • u k`, `s > 0`. A current `I k` along it has current density
+`I k` times the unit tangent times the arc-length measure on the ray; in the parameter `s` this is
+`I k • u k` times the Lebesgue measure, whatever the length of `u k`. The currents are steady, so
+the pairing with a test function `η` only sees its time integral.
+
+-/
+
+/-- The current density of `J` is that of steady currents `I k` flowing away from the origin
+along thin, straight, semi-infinite wires in the directions `u k`. -/
+def IsWireJunction {d} (c : SpeedOfLight) (J : DistLorentzCurrentDensity d) {ι : Type}
+    [Fintype ι] (u : ι → EuclideanSpace ℝ (Fin d)) (I : ι → ℝ) : Prop :=
+  ∀ η, J.currentDensity c η =
+    ∑ k, (I k * ∫ s in Ioi (0 : ℝ), timeIntegralSchwartz η (basis.repr.symm (s • u k))) • u k
+
+/-- The wires of a junction meet at the origin, where they act as a source of strength the
+total current leaving along them. -/
+lemma IsWireJunction.distSpaceDiv_currentDensity {d} {c : SpeedOfLight}
+    {J : DistLorentzCurrentDensity d} {ι : Type} [Fintype ι] {u : ι → EuclideanSpace ℝ (Fin d)}
+    {I : ι → ℝ} (hJ : J.IsWireJunction c u I) (hu : ∀ k, u k ≠ 0) (η : 𝓢(Time × Space d, ℝ)) :
+    distSpaceDiv (J.currentDensity c) η = (∑ k, I k) * ∫ t, η (t, 0) := by
+  have hJ' : ∀ η, J.currentDensity c η = _ := hJ
+  rw [distSpaceDiv_apply_eq_sum_distSpaceDeriv]
+  simp only [distSpaceDeriv_apply', PiLp.neg_apply, hJ']
+  simp only [timeIntegralSchwartz_fderiv_space, map_smul, WithLp.ofLp_sum, WithLp.ofLp_smul,
+    Finset.sum_apply, Pi.smul_apply, smul_eq_mul]
+  rw [Finset.sum_neg_distrib, Finset.sum_comm, Finset.sum_mul, ← Finset.sum_neg_distrib]
+  refine Finset.sum_congr rfl fun k _ => ?_
+  simp only [mul_assoc, ← Finset.mul_sum, integral_Ioi_fderiv_ray_basis _ (hu k),
+    timeIntegralSchwartz_apply]
+  ring
+
+/-!
+
+### C.2. The currents at a junction sum to zero
+
+-/
+
+/-- If Maxwell's equations hold and no charge accumulates, the steady currents leaving a junction
+of thin wires sum to zero. This is motivated by Kirchhoff's current law, of which it is the
+special case of a single node. -/
+lemma IsWireJunction.sum_currents_eq_zero {d} {𝓕 : FreeSpace} {J : DistLorentzCurrentDensity d}
+    {ι : Type} [Fintype ι] {u : ι → EuclideanSpace ℝ (Fin d)} {I : ι → ℝ}
+    (hJ : J.IsWireJunction 𝓕.c u I) (hu : ∀ k, u k ≠ 0) (A : DistElectromagneticPotential d)
+    (h : DistElectromagneticPotential.IsExtrema 𝓕 A J)
+    (hρ : distTimeDeriv (J.chargeDensity 𝓕.c) = 0) :
+    ∑ k, I k = 0 := by
+  obtain ⟨η, hη⟩ := exists_integral_time_ne_zero (d := d)
+  have hcont := DistElectromagneticPotential.continuityEquation A J h η
+  rw [hρ, _root_.zero_apply, zero_add, hJ.distSpaceDiv_currentDensity hu] at hcont
+  exact (mul_eq_zero.mp hcont).resolve_right hη
+
+end DistLorentzCurrentDensity
+end Electromagnetism