Muon-catalyzed fusion (μCF) of deuterium and helium-3 combines two attractive properties: an aneutronic reaction (d + ³He → p + α, Q = 18.35 MeV) whose products are entirely charged particles amenable to direct electric conversion, and a catalysis mechanism that operates at low temperature without magnetic or inertial confinement. It also combines two severe liabilities: the muon's finite lifetime (2.197 μs) and finite production cost (~5 GeV per muon in practice), and — specific to the helium system — an irreversible muon transfer to ³He followed by a muonic molecular complex whose fusion branching ratio is many orders of magnitude below unity. This paper treats the standard verdict ("d–³He μCF cannot break even") not as a prohibition but as a specification. We derive the breakeven inequality, quantify the current shortfall in each factor, and formalize four requirements (R1–R4) that, if simultaneously met, yield engineering net gain. Each requirement is mapped to an identifiable physics program with explicit falsifiers. The gap is dominated by a single factor — the fusion branching of the (dμ³He) complex, currently ~10⁻⁶–10⁻⁵ against a required ~10⁻¹ — and we state precisely what discovery would close it and what null result would kill the scheme. All load-bearing claims are cataloged in a Claims Ledger (Appendix A) by layer and epistemic status. The cartoon closes; the real problem stays open.
An essay. A thesis put forward by its author — an argument, specification, or design — and not a verified result.
1. Introduction: The Prohibition and the Reframe
The textbook story of muon-catalyzed fusion is a story of near-miss. A negative muon, 206.77 times heavier than an electron, shrinks the Bohr radius of a hydrogenic atom by the same factor — from 53 pm to ~256 fm. Two nuclei bound by a muon sit close enough that quantum tunneling through the Coulomb barrier proceeds in picoseconds rather than geological ages. Alvarez's 1956 observation of muon-catalyzed p–d fusion in a hydrogen bubble chamber, and the subsequent theory of resonant muonic-molecule formation (Vesman mechanism), established a genuine cold-fusion channel — one that is real, reproducible, and measured to high precision.
The prohibition, as usually stated: a muon costs more energy to make than it returns before it decays or sticks. In the best-studied system (d–t), a single muon catalyzes on the order of 100–150 fusions before loss, returning roughly 2–2.6 GeV of fusion energy against a production cost of ~5 GeV and a conversion efficiency of ~1/3. Net negative by roughly an order of magnitude. For d–³He, the shortfall is not one order of magnitude but eight or nine, for reasons specific to helium's Z = 2 charge.
The reframe this paper applies — consistent with the broader thesis that physical laws are permission structures — is to ask: what is the exact inequality that must be satisfied, which physical quantities appear in it, and which of those quantities are fixed by theorem versus limited only by current technique or current ignorance? The answer, developed below, is that no theorem forbids net-positive d–³He μCF. What stands in the way is (i) one open problem in three-body muonic molecular physics, (ii) one hard but bounded accelerator-engineering problem, and (iii) one direct-conversion engineering problem that is physics-sound but unbuilt at scale. That is a specification, not a wall. Whether the specification is meetable is precisely the open question, and Section 7 states what experimental results would settle it in either direction.
2. Physics Background: The Catalysis Cycle
2.1 Muonic atoms and molecules
When a negative muon stops in matter it cascades into the 1s state of a muonic atom in ~10⁻¹¹ s. Ground-state binding energies scale as Z² times the reduced-mass-corrected Rydberg:
- μp: ≈ 2.53 keV
- μd: ≈ 2.66 keV
- μt: ≈ 2.71 keV
- μ³He (1s): ≈ 10.9 keV
The factor-of-four deeper binding on helium is the origin of both the promise and the pathology of the d–³He system: any muon that encounters ³He while bound to a hydrogen isotope transfers irreversibly downhill.
2.2 The d–t benchmark
In D/T mixtures the cycle is: capture → μd → transfer to μt → resonant formation of the (dtμ) molecular ion via the Vesman mechanism (the weakly bound J = 1, v = 1 state, binding ≈ 0.66 eV, allows the formation energy to be absorbed by rovibrational excitation of the host D₂ molecule) → deexcitation to J = 0 → fusion at ~10¹² s⁻¹ → muon released, except for the fraction ω_s ≈ 0.4–0.6% that sticks to the outgoing α. At liquid-hydrogen density the cycle rate is ~10⁸ s⁻¹, yielding the observed ~100–150 fusions per muon (LAMPF, PSI, RIKEN-RAL). The cycle count is capped by effective sticking at N ≲ 1/ω_eff ≈ 170–200, not primarily by the muon lifetime.
Two structural facts carry over to the helium analysis:
- Resonant molecular formation is the enabling miracle. Without the accidental existence of a state bound by less than the ~few-eV scale of molecular rovibrational quanta, formation rates drop by orders of magnitude.
- The cycle count budget is a product of branching ratios, each of which is a separate physics target.
2.3 The d–³He reaction
d + ³He → p (14.68 MeV) + α (3.67 MeV), Q = 18.35 MeV.
Aneutronic in its primary channel (a small d–d side reaction in any deuterated target produces some neutrons; this is a design nuisance, not a physics blocker). Both products are charged, opening direct electrostatic or electromagnetic conversion at plausibly η ≈ 0.6–0.8 — roughly double the Carnot-limited thermal route available to d–t's neutron-borne energy. This is the entire motivation for accepting helium's Z = 2 penalty.
3. The d–³He Muonic System: Where the Orders of Magnitude Go
The cycle in a D₂/³He mixture proceeds as follows. Rates quoted are order-of-magnitude, normalized to liquid hydrogen density (LHD); Ledger items C7–C10 flag them for verification against primary sources before publication use.
Step 1 — Capture and transfer. The muon predominantly forms μd (hydrogen dominates the stopping). Transfer μd + ³He → μ³He + d proceeds at ~10⁸ s⁻¹ per unit ³He concentration and is irreversible (ΔE ≈ 8.2 keV downhill). Within tens of nanoseconds, essentially every muon resides on ³He. In the d–t system helium is a poison; here it is the designated carrier. This reframe matters: the scheme should not fight the transfer but route through it.
Step 2 — Complex formation. μ³He + D₂ → [(dμ³He)⁺ d e]. Formation rates of order 10⁸ s⁻¹ at LHD have been measured (PSI, JINR/Dubna programs). Over one muon lifetime this permits of order 10²–10³ complex formations per muon. So far, so good — the cycle rate is not the problem.
Step 3 — The branching catastrophe. The (dμ³He) complex is not a compact three-body molecule like (dtμ). Because the muon binds ~4× more strongly to the helium core, the system is better described as a μ³He atom weakly perturbed by a deuteron at large separation. The complex de-excites — predominantly by Auger emission and radiative decay back to μ³He + d — at rates of order 10¹¹–10¹² s⁻¹. The nuclear fusion rate from the relevant states is calculated and measured at roughly 10⁵–10⁷ s⁻¹. The fusion branching per complex is therefore
B_f ≈ λ_fusion / λ_decay ~ 10⁻⁷–10⁻⁵.
Step 4 — Yield per muon. Combining: N_f ≈ (complexes per lifetime) × B_f ~ 10³ × 10⁻⁵ ≈ 10⁻³–10⁻⁴ fusions per muon (order-of-magnitude; experimental yields from the Dubna/PSI dμ³He campaigns are consistent with this range). Energy return per muon: ~10⁻³ × 18.35 MeV ≈ tens of keV, against a ~5 GeV production cost.
The shortfall is a factor of ~10⁵–10⁶ in fusion count, concentrated almost entirely in one number: B_f. By contrast, in d–t the same branching is effectively unity (λ_fusion ~ 10¹² s⁻¹ from J = 0 beats every competing channel). The physical reason is twofold: the muon in dμ³He screens the d–³He pair poorly (it lives on the helium), and the tunneling integral carries Z₁Z₂ = 2 in the Gamow exponent at a larger effective internuclear distance.
This is the correct location of the prohibition. Everything else — muon cost, conversion efficiency, cycle rate — is within one to two orders of magnitude of viability. B_f is not.
4. Energy Accounting: The Breakeven Inequality
Define:
C_μ — wall-plug energy to produce and deliver one muon [GeV] N_f — fusions catalyzed per muon Q — energy per fusion = 18.35 MeV (d–³He) η_c — conversion efficiency of fusion products to electricity η_r — fraction of beam/target waste heat recovered (second-order; set to 0 for the strict inequality)
Net-positive condition (engineering gain G > 1):
G = (η_c · N_f · Q) / C_μ > 1 → N_f > C_μ / (η_c · Q)
Evaluated:
- Scenario C_μ η_c Required N_f
- Present accelerator practice, thermal conversion 5 GeV 0.35 ~2,850
- Present practice, direct conversion 5 GeV 0.60 ~455
- Optimized muon source (design studies) 2 GeV 0.60 ~180
- Speculative muon source floor 0.5 GeV 0.70 ~39
- Current d–³He experimental reality 5 GeV — N_f ~ 10⁻³ achieved
Two readings of this table:
- The permissive reading. The required N_f under aggressive-but-stated assumptions (~40–450) is smaller than what d–t μCF already achieves (~150). The d–³He system does not need thousands of cycles; direct conversion buys the difference. No individual row requires new physics except the leap from 10⁻³ to 10²·⁵ in N_f.
- The honest reading. That leap is 5–6 orders of magnitude and lives in B_f, an open problem in three-body Coulomb physics, not an engineering knob. The table is a map of the wall, not a hole in it.
For comparison, d–t μCF sits at N_f ≈ 150 achieved versus ~950 required (5 GeV, η = 0.35, Q = 17.6 MeV) — a factor of ~6, which is why d–t remains the nearest-viable μCF system and why every serious hybrid proposal (μCF neutron sources driving fission blankets) is built on d–t. This paper does not dispute that ranking; it asks what would have to be true for the aneutronic system to leapfrog it.
5. The Specification: Four Requirements for Net-Positive d–³He μCF
Following the methodology of stating permissions rather than prohibitions, net gain requires all four of the following simultaneously. Each is stated with its current status and the research program that addresses it.
R1 — Fusion branching rescue: B_f ≥ 10⁻¹ per complex
Status: open problem (physics layer).
The single load-bearing requirement. Three candidate mechanisms, in decreasing order of respectability:
- R1a. Resonant or quasi-bound intermediate states. The dtμ miracle was a state bound by 0.66 eV. The question is whether the (dμ³He) spectrum — including states correlating to excited μ³He (n = 2, the 2pσ quasimolecular term discussed in the muonic-helium literature) — contains any long-lived configuration with enhanced nuclear overlap. Side-path formation through excited muonic states, before the cascade completes, is a distinct channel from the ground-state complex of Section 3 and is less thoroughly excluded. Falsifier: complete high-precision three-body calculations plus X-ray spectroscopy of the complex; if the spectrum contains no state with λ_fusion/λ_decay above ~10⁻², R1a is dead.
- R1b. Environmental engineering of the decay channels. B_f can rise either by raising λ_fusion or lowering λ_decay. Auger rates depend on the electronic environment; radiative rates on available phase space. Dense plasma, degenerate, or cavity-modified environments alter both. Expected effect sizes from known physics are factors of a few, not 10⁵. Marked speculative; retained because it is cheap to test.
- R1c. Field-assisted tunneling. Intense laser or X-ray fields modulating the tunneling barrier during the complex lifetime. Barrier-modification physics is real (Keldysh-regime ionization is the electronic analog), but the required field scales at femtometer-relevant distances are beyond any current source by many orders. Myth-layer until a mechanism with realistic fields is exhibited.
If R1 has no solution, the scheme fails by inequality — cleanly and permanently. This is the honest center of gravity of the whole paper.
R2 — Muon cost: C_μ ≤ 2 GeV delivered
Status: hard engineering, no theorem violated (dynamics layer).
Muons are made by protons on target → pions → μ decay. The practical ~5–8 GeV/μ figure reflects unoptimized collection. Design studies for muon colliders and neutrino factories (high-field solenoid capture channels, thick-target geometries, phase rotation) support ~2 GeV/μ as an engineering target; the in-principle floor is set by pion production cross-sections and capture solid angle, plausibly a few hundred MeV/μ but never below the ~106 MeV rest mass plus production inefficiencies. Energy-recovery linac architectures additionally recycle unspent beam power. Falsifier-resistant: this requirement cannot be forbidden, only unfunded.
R3 — Direct conversion: η_c ≥ 0.6 on 14.7 MeV protons
Status: physics-sound, unbuilt at scale (dynamics layer).
Electrostatic decelerator stacks and traveling-wave direct converters for monoenergetic charged products have decades of paper design (originating in mirror-machine and aneutronic-fusion literature) with small-scale demonstrations at lower energies. The 14.7 MeV proton is nearly monoenergetic — the best possible input for staged electrostatic conversion. Voltage-holding and current-density limits are the binding constraints. Nothing here requires discovery; it requires a program.
R4 — Cycle throughput and target management
Status: measured physics, optimization problem (dynamics layer).
Maximize complexes per muon lifetime: density (rates scale ~linearly to LHD and beyond), D/³He concentration ratio (transfer time vs. complex formation time), temperature (formation rates are temperature-dependent; the resonance structure, if any, will be too). Includes suppressing the d–d side channel (neutron hygiene) and managing muon end-of-life: notably, α-sticking is not the dominant loss here as it is in d–t, because the muon's resting state is μ³He anyway — the cycle re-enters at Step 2. One of the d–t system's two great plagues is structurally absent. This is worth stating in the ledger as a genuine (if currently academic) advantage of the helium route.
6. Architecture Sketch (Conditional on R1)
For completeness, and explicitly conditional on an R1 solution, the plant closes as follows: energy-recovery proton driver (R2) → pion capture channel → muon injection into a high-density D₂/³He target held at the R1-optimal thermodynamic point (R4) → charged products harvested by staged electrostatic converters (R3) → recirculation of converter output to the driver. The gain equation of Section 4 then reads, with C_μ = 2 GeV, η_c = 0.6, and N_f = 300 (i.e., B_f ≈ 0.3 at 10³ complexes/muon): G ≈ 1.65 — thin but positive, with the d–d neutron channel relegated to shielding rather than blanket engineering.
This paragraph is the cartoon. It is load-bearing only as a demonstration that no additional hidden requirement lurks behind R1–R4.
7. Falsifiability and Kill Conditions
A specification is only honest if it can fail. Concretely:
- Complete (dμ³He) spectrum calculation + spectroscopy. Modern three-body variational methods (the same machinery that nailed dtμ binding to μeV precision) applied exhaustively to the d–³He–μ system, including n = 2 correlating states, cross-checked by X-ray spectroscopy at a muon facility (PSI, J-PARC, RIKEN-RAL class). Null result → R1a dead.
- B_f vs. environment scan. Measure fusion yield across density, temperature, ortho/para composition, and applied field. Known physics predicts factor-of-few variation. If the scan confirms factor-of-few, R1b dead.
- If both die, the scheme is falsified by the Section 4 inequality — not by taste, but by arithmetic. The correct conclusion in that world: d–³He μCF is permanently closed; d–t μCF hybrids and thermonuclear d–³He remain the live aneutronic-adjacent routes.
- Conversely, discovery of any (dμ³He) configuration with B_f ≳ 10⁻² immediately promotes R2/R3 from academic to urgent, because the remaining gap becomes pure engineering.
The asymmetry of the bet is the point: the falsifying experiments are cheap relative to fusion-program budgets (muon beamtime plus theory), while the payoff branch, if it exists, changes the ranking of fusion approaches.
8. Conclusion
The claim "muon-catalyzed d–³He fusion cannot be net positive" decomposes, under inspection, into one theorem-grade fact (muons decay in 2.2 μs and cost real energy), two engineering gaps that no physical law protects (muon production cost, direct conversion at scale), and exactly one open problem in few-body Coulomb physics (the fusion branching of the muonic complex) on which the entire scheme pivots. The prohibition is real today; it is not written into the laws. It is written into a three-body spectrum that has not been exhaustively searched. The specification is strict: B_f must improve by five orders of magnitude, and there is currently no predicted mechanism that delivers it — a fact this paper does not soften. The permission structure framing does not promise the door opens; it locates the door, prices the key, and names the experiment that determines whether the lock has a keyhole at all.
The cartoon closes; the real problem stays open.
Author
- Oleksii Vergunenko
References
- 1.L. W. Alvarez et al., "Catalysis of Nuclear Reactions by μ Mesons," Phys. Rev. 105, 1127 (1957).
- 2.J. D. Jackson, "Catalysis of Nuclear Reactions between Hydrogen Isotopes by μ⁻ Mesons," Phys. Rev. 106, 330 (1957).
- 3.E. A. Vesman, JETP Lett. 5, 91 (1967) — resonant muonic-molecule formation mechanism.
- 4.L. I. Ponomarev, "Muon Catalysed Fusion," Contemp. Phys. 31, 219 (1990) — canonical review.
- 5.P. Petitjean et al., PSI muon-catalyzed fusion program publications (1990s) — d–t cycling and sticking measurements.
- 6.V. M. Bystritsky et al., JINR/Dubna measurements of fusion in dμ³He complexes (1990s–2000s).
- 7.W. Czapliński et al., theoretical fusion and decay rates of the (dμ³He) system.
- 8.K. Nagamine, Introductory Muon Science, Cambridge University Press (2003) — RIKEN-RAL program and μCF energy accounting.
- 9.Muon collider / neutrino factory design reports — pion capture and muon production cost estimates.
Experimentation loop
0. What the simulation changed before any beam time
Running the coupled model produced three findings that restructure the experimental program. Per the principle that negative results must be as informative as positive ones, these are promoted to first-class results:
F1 — A hidden fifth requirement (R5). The paper's Section 6 claimed the architecture sketch demonstrated "no additional hidden requirement lurks behind R1–R4." The coupled simulation falsifies that claim (Ledger C16 premise). Once B_f is large (R1 success branch), yield saturates at
N_f^max = 1 / (ω_eff + λ₀/λ_form)
— a ceiling set by per-fusion muon loss (sticking-type, ω_eff) and by complex formation throughput. At LHD, the formation-only ceiling is ~500 fusions/muon. Sticking, provably irrelevant today (test T5: <10⁻⁴ effect), becomes the binding constraint the moment R1 succeeds (test T6). Ceiling requirements by scenario:
Scenario (C_μ, η_c) Required N_f Implied ω_eff ceiling 5 GeV, 0.35 (thermal) 779 unreachable at LHD — exceeds formation ceiling even at ω = 0 5 GeV, 0.60 (direct) 454 ≤ 0.02% 2 GeV, 0.60 182 ≤ 0.35% 0.5 GeV, 0.70 39 ≤ 2.4%
F2 — The thermal route is dead at LHD regardless of R1. No value of B_f rescues the (5 GeV, thermal) scenario at liquid-hydrogen density: 779 > 501. Direct conversion (R3) is therefore not an optimization but a requirement, and density above LHD (compressed/liquid targets) buys ceiling linearly.
F3 — The paper's plant point is revised. Sec. 6 quoted G ≈ 1.65 at (B_f = 0.3, C_μ = 2 GeV, η_c = 0.6). Under the ceiling formula that point actually yields G ≈ 0.6. The corrected plant point for G ≈ 1.65 at 2 GeV / 0.6: B_f ≥ 0.8 with ω_eff ≤ 0.08% at LHD (or B_f ≥ 0.5 with ω_eff ≤ 0.07% at φ = 1.5). The R1 target hardens from "any large branching" to "near-unity branching," and ω_s^(d³He) — currently an unmeasured placeholder — becomes a mandatory loop observable.
Model self-checks passed: T1 (nominal yield 6.7×10⁻⁴ inside the measured 10⁻⁴–10⁻³ window), T2 (d–t reparameterization gives 131 fusions/muon vs. observed 100–160), T3 (kinetics reproduce the Section-4 breakeven inequality to 10⁻¹⁰ in the zero-loss limit), T4 (bottleneck localization: ×10 on λ_fus moves yield and gain ×10 with <0.02% shift in complex throughput — branching, not kinetics, is the wall), T5, T6 as above.
1. Loop architecture
┌────────────────────────────────────────────────────────────┐
│ STAGE 0 — THEORY ANCHOR │
│ Exhaustive (d μ ³He) three-body spectrum, incl. n=2/2pσ │
│ channels. Predictions REGISTERED before beam time: │
│ state energies, λ_fus, λ_dec, x-ray line positions, ω_s │
└───────────────┬────────────────────────────────────────────┘
│ unfitted predictions
▼
┌────────────────────────────────────────────────────────────┐
│ STAGE 1 — SPECTROSCOPY pins: λ_form, λ_dec split │
│ STAGE 2 — FUSION YIELD pins: λ_fus·λ_form → B_f │
│ STAGE 3 — ENVIRONMENT SCAN pins: R1b enhancement law │
│ STAGE 4 — STICKING (ω_s) pins: R5 feasibility │
└───────────────┬────────────────────────────────────────────┘
│ measured rates + errors
▼
┌────────────────────────────────────────────────────────────┐
│ STAGE 5 — FEEDBACK │
│ Bayesian update of {λ_tr, λ_form, λ_dec, λ_fus, ω_s}; │
│ rerun mucf_sim acceptance suite T1–T6 on posteriors; │
│ promote Ledger C6–C9 from VERIFY → MEASURED; │
│ evaluate gates G0–G3 │
└───────┬───────────────────────────────────┬────────────────┘
│ gate positive │ all gates null
▼ ▼
escalate: R2/R3 engineering KILL: scheme falsified by
programs + Stage-4 priority inequality; publish negative
(loop re-enters at Stage 0 result + closed ledger
with refined targets)
2. Parameter ↔ observable mapThe contract between simulation and experiment. Each row names the sim parameter, the physical observable that pins it, and the measurement channel.
Sim parameter Ledger Observable Channel λ_tr0 (μd → μ³He transfer) C6-adj disappearance of μd x-ray series / delayed μ³He line onset vs. c_He muonic x-ray timing, 20–200 ns window λ_form0 (complex formation) C6 growth rate of complex signatures vs. φ·c_D 6.85 keV radiative-decay x-ray intensity vs. time λ_dec (Auger + radiative split) C7 absolute intensity + lineshape of the ~6.85 keV line; Auger electron yield Si drift detectors; line/continuum ratio λ_fus (→ B_f with λ_dec) C8 14.68 MeV proton rate per muon stop, time-correlated with complex population ΔE–E Si telescope + CsI; delayed-coincidence with μ stop ω_s^(d³He) new (R5) fraction of fusion events with no subsequent muon recycle signature; (μα)⁺ x-rays high-statistics only if B_f enhanced; theory-pinned with error bars until then φ, c_He, T R4 target state direct instrumentation of target cell
Cross-check invariant: Stage 1 and Stage 2 overdetermine the system — B_f is extracted both as λ_fus/λ_dec (rates) and as N_f/(complexes per muon) (counting). Agreement within errors is the loop's internal-consistency test; disagreement means the state model (the network topology itself, not just its rates) is wrong, which re-enters Stage 0.
3. Stage specifications Stage 0 — Theory anchor (no beam)
Modern three-body variational machinery (the toolset that pinned dtμ binding to μeV) applied exhaustively to d–³He–μ, including channels correlating to μ³He(n = 2) — the 2pσ quasimolecular term — and any long-lived configurations formed during the atomic cascade, before ground-state arrival. Deliverable: a registered prediction file (state energies, widths, λ_fus, λ_dec, x-ray energies, ω_s estimate) timestamped before Stage 1 beam time. This enforces the unfitted-prediction discipline: post-hoc agreement is worth nothing; registered agreement is worth everything.
Gate G0: does the computed spectrum contain any configuration with λ_fus/λ_dec ≥ 10⁻²? (Note: F3 hardened the ultimate target to near-unity, but 10⁻² is the "mechanism exists" threshold that justifies escalation.)
Stage 1 — Spectroscopy (muon facility: PSI πE5/μE4-class, J-PARC MUSE, RIKEN-RAL)
Stop negative muons in a D₂/³He gas/liquid cell; measure the muonic x-ray cascade and, critically, the ~6.85 keV photon from radiative decay of the (dμ³He) complex — an established experimental signature of complex formation. Time-resolved spectra vs. φ ∈ [0.05, 1.2] LHD and c_He ∈ [0.005, 0.3] pin λ_tr0, λ_form0, and the λ_dec radiative/Auger split.
Gate G1: do measured line positions and rates match the Stage-0 registered predictions within stated errors? A mismatch is good news operationally (the spectrum holds surprises — exactly where an R1 mechanism would hide) but sends the loop back to Stage 0 before yield interpretation.
Stage 2 — Fusion yield
Detect the 14.68 MeV proton (ΔE–E silicon telescope backed by CsI; the proton is far above every background line). Backgrounds and discriminants: d–d fusion protons at 3.02 MeV (energy cut), d–d neutrons at 2.45 MeV (pulse shape + energy), μ⁻ capture on ³He → t + ν with a 1.9 MeV triton (energy cut), decay electrons (dE/dx). Delayed-coincidence with the muon stop plus the Stage-1 complex timing isolates the complex-mediated pathway.
Statistics: at N_f ≈ 7×10⁻⁴ and ~10⁴ muon stops/s (πE5-class), expect ~7 fusion protons/s → 10⁶ events per ~40 beam-hours per scan point. A 30-point (φ, T, c_He) grid is a few weeks of beam — cheap relative to any fusion-program budget, which is the asymmetry the whole bet rests on.
Output: B_f(φ, T, c_He) with errors, via both extraction routes (Section 2 invariant).
Stage 3 — Environment scan (R1b decision)
Extend Stage 2 across temperature 20–800 K, ortho/para D₂ ratio, density into the compressed-liquid regime, and applied E/B fields. Known physics predicts factor-of-few variation.
Gate G3: any enhancement of B_f beyond ×10 anywhere in the scan? Yes → an unmodeled mechanism exists; characterize and re-enter Stage 0. No → R1b is dead by measurement.
Stage 4 — Sticking (R5 feasibility)
ω_s^(d³He) is unmeasurable by cycle statistics at today's B_f (fusion events too rare for recycle-counting). Interim: pin by theory with honest error bars — the muon's post-fusion phase-space overlap with the 3.67 MeV α, computed in the same Stage-0 framework (sudden approximation + stripping corrections, the machinery validated on d–t's ω_s ≈ 0.9% initial / 0.55% effective). Priority flips to experimental the moment any gate shows B_f enhancement. Kill-relevant: if theory pins ω_s^(d³He) ≥ ~1% with high confidence, the 2 GeV scenario dies even under full R1 success (F1 table), leaving only the speculative-floor muon source alive.
Stage 5 — Feedback and closure
Posterior distributions over the five rates → rerun mucf_sim.py acceptance suite → regenerate Figs 1–5 with data bands → Ledger promotion (C6–C9: VERIFY → MEASURED; C14: OPEN → resolved either way; C16: already demoted by F1).
Loop-closure criterion: the simulation, running exclusively on measured parameters, reproduces the full measured yield surface N_f(φ, T, c_He) within experimental errors, with the same rate set satisfying both B_f extraction routes. At that point the model is no longer a model of the literature; it is a model of the apparatus — and the R1 verdict is read off Fig 2 with data replacing placeholders.
4. Decision table Condition Verdict G0 null AND G3 null KILL. d–³He μCF falsified by inequality. Publish negative result with closed ledger. d–t hybrids and thermonuclear d–³He remain the live routes. G0 positive, G1 match Escalate: dedicated beam program targeting the predicted state; Stage 4 goes experimental; R2/R3 move from academic to urgent. G1 mismatch Model topology wrong — return to Stage 0. (Informative in both directions.) G3 positive, G0 null Unpredicted mechanism — highest-value outcome; full re-derivation at Stage 0. Any escalation + ω_s theory-pin ≥ 1% Scheme survives only at C_μ ≤ ~0.5 GeV; R2 becomes the sole load-bearing engineering program. 5. Ledger delta (from simulation run) # Change C16 Demoted: "no hidden requirements" premise falsified by F1. Cartoon corrected in place. C17 (new) N_f^max = 1/(ω_eff + λ₀/λ_form): theorem-grade within the model, ceiling ≈ 500/muon at LHD, ω = 0. Awaits λ_form measurement (Stage 1). C18 (new) Thermal-conversion route unreachable at LHD independent of B_f. Dynamics layer; density-dependent. C19 (new) Revised plant point: B_f ≥ 0.8 AND ω_eff ≤ 0.08% (LHD, 2 GeV, η = 0.6) for G ≈ 1.65. Myth layer until Stages 1–4 report. C20 (new) ω_s^(d³He) is unmeasured; interim theory-pin required (Stage 4). Open problem.
Simulation
Resultsjson
{
"tiers": {
"NOMINAL": {
"ceiling_Nf": 44.08904600332008,
"c_He_opt": 0.1573897199359354,
"G_max": 0.05663237959126464,
"Bf_G1": Infinity,
"Bf_G3": Infinity,
"Nf_req_G1": 778.5130400934215
},
"NOMINAL+direct": {
"ceiling_Nf": 44.08904600332008,
"c_He_opt": 0.1573897199359354,
"G_max": 0.09708407929931082,
"Bf_G1": Infinity,
"Bf_G3": Infinity,
"Nf_req_G1": 454.1326067211626
},
"ENGINEERING": {
"ceiling_Nf": 63.2226664048671,
"c_He_opt": 0.13539603261322466,
"G_max": 0.4292502935558452,
"Bf_G1": Infinity,
"Bf_G3": Infinity,
"Nf_req_G1": 147.28625082848515
},
"THEORETICAL": {
"ceiling_Nf": 173.6386627539879,
"c_He_opt": 0.08666479987241328,
"G_max": 9.558808384607035,
"Bf_G1": 0.01351677158057538,
"Bf_G3": 0.05108754128943265,
"Nf_req_G1": 18.1653042688465
},
"THEORETICAL+CONJ": {
"ceiling_Nf": 191.53374349683364,
"c_He_opt": 0.15363322039711783,
"G_max": 10.543932579500693,
"Bf_G1": 0.0033892722631564593,
"Bf_G3": 0.012809826312912762,
"Nf_req_G1": 18.1653042688465
}
},
"waterfall": [
[
"practice, thermal",
null,
44.08904600332008
],
[
"+ direct conversion",
null,
44.08904600332008
],
[
"+ reactivation (omega 1.4%)",
null,
60.3940203329627
],
[
"+ solid target (phi 1.4)",
null,
63.2226664048671
],
[
"+ 2 GeV muons",
null,
63.2226664048671
],
[
"+ thermal bottoming (0.74)",
null,
63.2226664048671
],
[
"+ omega 0.5% (theor.)",
0.9133795970853524,
150.4680995073782
],
[
"+ phi 3 pulsed (theor.)",
0.40126908692718904,
173.6386627539879
],
[
"+ 0.3 GeV muons, eta 0.90",
0.01351677158057538,
173.6386627539879
],
[
"+ resonant formation (CONJ)",
0.0033892722631564593,
191.53374349683364
]
],
"B_f_today": 1.4285684588507947e-06,
"tests": {
"T7": true,
"T8": true
}
}Codepython
"""
Muon-catalyzed d-3He fusion: coupled kinetics + energy-balance simulation.
Companion to "Muon-Catalyzed d-3He Fusion as a Permission Structure".
Model: linear rate-equation network for the muon state populations
f_d : muon bound in mu-d (initial capture)
f_He : muon bound in mu-3He (post-transfer resting state)
f_c : (d mu 3He) molecular complex
N : cumulative fusions per injected muon
Muon decay (lam0) drains every state. Complex branches between
non-fusion decay (lam_dec, Auger/radiative -> back to mu-3He) and
in-complex fusion (lam_fus). Fusion returns the muon to mu-3He
except a sticking fraction omega_s (structurally minor here).
All rates are order-of-magnitude placeholders normalized to liquid
hydrogen density (LHD), flagged in the paper's Claims Ledger (C6-C9)
as VERIFY. The experimental loop exists to replace them with data.
Exact solution via matrix exponential (system is linear).
"""
import numpy as np
from scipy.linalg import expm
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
import json
# ----------------------------------------------------------------------
# Parameters (nominal, LHD-normalized)
# ----------------------------------------------------------------------
TAU_MU = 2.197e-6 # s, muon lifetime
LAM0 = 1.0 / TAU_MU # s^-1, muon decay
Q_DHE3 = 18.35 # MeV per d-3He fusion
Q_DT = 17.6 # MeV per d-t fusion
NOMINAL = dict(
phi = 1.0, # target density in LHD units
c_He = 0.05, # 3He atomic concentration
lam_tr0 = 1.5e8, # s^-1 per (phi * c_He): mu-d -> mu-3He transfer [C6-adjacent, VERIFY]
lam_form0 = 2.4e8, # s^-1 per (phi * c_D): mu3He + D2 -> complex [C6, VERIFY]
lam_dec = 7.0e11, # s^-1: complex non-fusion decay (Auger/radiative) [C7, VERIFY]
lam_fus = 1.0e6, # s^-1: in-complex fusion [C8, VERIFY]
omega_s = 0.02, # sticking-type loss per fusion (upper-ish bound)
)
SCENARIOS = { # (C_mu [MeV], eta_c) for gain evaluation
"practice/thermal (5 GeV, 0.35)": (5000.0, 0.35),
"practice/direct (5 GeV, 0.60)": (5000.0, 0.60),
"optimized source (2 GeV, 0.60)": (2000.0, 0.60),
"speculative floor(0.5 GeV, 0.70)": (500.0, 0.70),
}
# ----------------------------------------------------------------------
# Core solver
# ----------------------------------------------------------------------
def rate_matrix(p):
lt = p["lam_tr0"] * p["phi"] * p["c_He"]
lf = p["lam_form0"] * p["phi"] * (1.0 - p["c_He"])
ld = p["lam_dec"]
lu = p["lam_fus"]
w = p["omega_s"]
# states: [f_d, f_He, f_c, N]
M = np.zeros((4, 4))
M[0, 0] = -(LAM0 + lt)
M[1, 0] = lt
M[1, 1] = -(LAM0 + lf)
M[1, 2] = ld + (1.0 - w) * lu
M[2, 1] = lf
M[2, 2] = -(LAM0 + ld + lu)
M[3, 2] = lu
return M
def solve(p, t_end=25e-6, times=None):
M = rate_matrix(p)
x0 = np.array([1.0, 0.0, 0.0, 0.0])
if times is None:
return expm(M * t_end) @ x0
return np.array([expm(M * t) @ x0 for t in times])
def fusions_per_muon(p):
return solve(p)[3]
def branching(p):
return p["lam_fus"] / (p["lam_fus"] + p["lam_dec"] + LAM0)
def complexes_per_muon(p):
"""Throughput diagnostic: complex formations per injected muon."""
q = dict(p)
lf = q["lam_form0"] * q["phi"] * (1.0 - q["c_He"])
# integrate lam_form * f_He exactly by adding a counter row
M = rate_matrix(q)
A = np.zeros((5, 5)); A[:4, :4] = M
A[4, 1] = lf
x0 = np.array([1.0, 0.0, 0.0, 0.0, 0.0])
return (expm(A * 25e-6) @ x0)[4]
def gain(Nf, C_mu_MeV, eta_c, Q=Q_DHE3):
return eta_c * Nf * Q / C_mu_MeV
def required_Nf(C_mu_MeV, eta_c, Q=Q_DHE3):
return C_mu_MeV / (eta_c * Q)
# ----------------------------------------------------------------------
# d-t benchmark (model validation target T2)
# states: [f_d, f_t, f_dtmu, N]; fusion fast, sticking loss omega
# ----------------------------------------------------------------------
def dt_fusions_per_muon(phi=1.2, c_t=0.4, lam_tr0=2.8e8, lam_form0=3.0e8,
lam_fus=1.1e12, omega=0.0055):
lt = lam_tr0 * phi * c_t
lf = lam_form0 * phi * (1 - c_t)
M = np.zeros((4, 4))
M[0, 0] = -(LAM0 + lt)
M[1, 0] = lt
M[1, 1] = -(LAM0 + lf)
M[1, 2] = (1.0 - omega) * lam_fus # muon recycled to mu-t pool
M[2, 1] = lf
M[2, 2] = -(LAM0 + lam_fus)
M[3, 2] = lam_fus
x0 = np.array([1.0, 0, 0, 0])
return (expm(M * 25e-6) @ x0)[3]
# ----------------------------------------------------------------------
# Acceptance tests
# ----------------------------------------------------------------------
def run_tests():
results = {}
p = dict(NOMINAL)
# T1: nominal d-3He yield lands in the measured order-of-magnitude window
Nf = fusions_per_muon(p)
results["T1_nominal_yield"] = dict(
value=Nf, window=[1e-4, 1e-3], passed=bool(1e-4 <= Nf <= 1e-3))
# T2: d-t reparameterization reproduces the observed 100-160 fusions/muon
Ndt = dt_fusions_per_muon()
results["T2_dt_benchmark"] = dict(
value=Ndt, window=[100, 160], passed=bool(100 <= Ndt <= 160))
# T3: G=1 crossing matches the analytic inequality N_f = C_mu/(eta Q)
# in the zero-loss limit (omega_s = 0). Tests arithmetic consistency of
# the kinetics against the Section-4 inequality.
C_mu, eta = 2000.0, 0.60
Nf_req = required_Nf(C_mu, eta)
p0 = dict(p); p0["omega_s"] = 0.0
lo, hi = 1e6, 1e14
for _ in range(80):
mid = np.sqrt(lo * hi)
q = dict(p0); q["lam_fus"] = mid
(lo, hi) = (mid, hi) if fusions_per_muon(q) < Nf_req else (lo, mid)
q = dict(p0); q["lam_fus"] = np.sqrt(lo * hi)
G_at_crossing = gain(fusions_per_muon(q), C_mu, eta)
results["T3_breakeven_consistency"] = dict(
value=G_at_crossing, window=[0.99, 1.01],
passed=bool(abs(G_at_crossing - 1) < 0.01))
# T4b/T6: DISCOVERED CEILING. With B_f -> 1, N_f saturates at
# 1/(omega_eff + lam0/lam_form): sticking + formation-rate ceiling.
# This is the model finding that falsifies the paper's C16 premise
# ("no hidden requirement behind R1-R4"): a fifth requirement (R5)
# binds as soon as R1 succeeds. Test: model ceiling matches analytic
# within 10% (offset = lifetime spent pre-transfer).
lam_form = p["lam_form0"] * p["phi"] * (1 - p["c_He"])
ceilings = {}
ok = True
for w in [0.0, 0.003, 0.02]:
q = dict(p); q["omega_s"] = w; q["lam_fus"] = 1e14
model = fusions_per_muon(q)
analytic = 1.0 / (w + LAM0 / lam_form)
ceilings[w] = dict(model=model, analytic=analytic)
ok &= abs(model / analytic - 1) < 0.10
results["T6_saturation_ceiling"] = dict(
value=ceilings, window="model vs 1/(w + lam0/lam_form) within 10%",
passed=bool(ok))
# T5: sticking insensitivity -- omega_s 0 -> 0.05 changes N_f by <10%
q0 = dict(p); q0["omega_s"] = 0.0
q5 = dict(p); q5["omega_s"] = 0.05
rel = abs(fusions_per_muon(q5) / fusions_per_muon(q0) - 1)
results["T5_sticking_insensitivity"] = dict(
value=rel, window=[0, 0.10], passed=bool(rel < 0.10))
return results
# ----------------------------------------------------------------------
# Figures
# ----------------------------------------------------------------------
plt.rcParams.update({"figure.dpi": 140, "font.size": 9,
"axes.grid": True, "grid.alpha": 0.3})
def fig1_kinetics(p):
times = np.logspace(-12, np.log10(20e-6), 400)
X = solve(p, times=times)
fig, ax = plt.subplots(figsize=(7, 4.2))
ax.loglog(times * 1e6, X[:, 0], label=r"$\mu d$")
ax.loglog(times * 1e6, X[:, 1], label=r"$\mu\,^3$He")
ax.loglog(times * 1e6, X[:, 2], label=r"$(d\mu\,^3\mathrm{He})$ complex")
ax.loglog(times * 1e6, X[:, 3], label=r"cumulative fusions $N_f$")
ax.axvline(TAU_MU * 1e6, ls=":", c="gray")
ax.text(TAU_MU * 1e6 * 1.1, 2e-13, r"$\tau_\mu$", color="gray")
ax.set_xlabel("time [µs]"); ax.set_ylabel("population / muon")
ax.set_ylim(1e-13, 2)
ax.set_title("Fig 1 — Muon state kinetics, nominal parameters\n"
"(complex population quasi-steady at ~10$^{-4}$; "
f"N$_f$(∞) = {fusions_per_muon(p):.1e})")
ax.legend(loc="center left", fontsize=8)
fig.tight_layout(); fig.savefig("/home/claude/out/fig1_kinetics.png")
def fig2_bf_bridge(p):
lam_fus_grid = np.logspace(4, 13.5, 180)
fig, ax = plt.subplots(figsize=(7.4, 4.8))
for w, c, lbl in [(0.02, "C0", r"model, $\omega_s=2\%$ (placeholder)"),
(0.003, "C9", r"model, $\omega_s=0.3\%$ (R5 met)")]:
Nf, Bf = [], []
for lu in lam_fus_grid:
q = dict(p); q["lam_fus"] = lu; q["omega_s"] = w
Nf.append(fusions_per_muon(q))
Bf.append(lu / (lu + p["lam_dec"]))
ax.loglog(Bf, Nf, lw=2, color=c, label=lbl)
for (name, (Cmu, eta)), c in zip(SCENARIOS.items(), ["C3", "C1", "C2", "C4"]):
ax.axhline(required_Nf(Cmu, eta), ls="--", c=c, lw=1,
label=f"breakeven: {name}")
lam_form = p["lam_form0"] * p["phi"] * (1 - p["c_He"])
ax.axhline(1.0 / (LAM0 / lam_form), ls="-.", c="gray", lw=1,
label=r"absolute ceiling $\lambda_{form}/\lambda_0$ ($\omega_s=0$)")
B0 = branching(p)
ax.scatter([B0], [fusions_per_muon(p)], zorder=5, color="k")
ax.annotate("today\n(C8 nominal)", (B0, fusions_per_muon(p)),
textcoords="offset points", xytext=(8, -18), fontsize=8)
ax.axvspan(1e-7, 1e-5, color="red", alpha=0.07)
ax.set_xlabel(r"fusion branching per complex $B_f=\lambda_{fus}/(\lambda_{fus}+\lambda_{dec})$")
ax.set_ylabel(r"fusions per muon $N_f$")
ax.set_title("Fig 2 — The five-decade bridge and the discovered ceiling:\n"
"R1 success alone clears breakeven only if R5 (sticking) is also met")
ax.legend(fontsize=6.5, loc="upper left")
fig.tight_layout(); fig.savefig("/home/claude/out/fig2_bf_bridge.png")
def fig3_gain_landscape(p):
lam_fus_grid = np.logspace(4, 13.5, 130)
Cmu = np.logspace(np.log10(300), np.log10(8000), 110) # MeV
fig, axes = plt.subplots(1, 2, figsize=(11.5, 4.6), sharey=True)
for ax, w, ttl in [(axes[0], 0.02, r"$\omega_s=2\%$ (placeholder)"),
(axes[1], 0.003, r"$\omega_s=0.3\%$ (R5 met)")]:
Nf, Bf = [], []
for lu in lam_fus_grid:
q = dict(p); q["lam_fus"] = lu; q["omega_s"] = w
Nf.append(fusions_per_muon(q))
Bf.append(lu / (lu + p["lam_dec"]))
Nf = np.array(Nf); Bf = np.array(Bf)
G = 0.60 * np.outer(Nf, 1.0 / Cmu) * Q_DHE3
pc = ax.pcolormesh(Cmu / 1000, Bf, np.log10(G),
shading="auto", cmap="RdYlGn", vmin=-7, vmax=1.5)
if G.max() >= 1.0:
cs = ax.contour(Cmu / 1000, Bf, G, levels=[1.0], colors="k", linewidths=2)
ax.clabel(cs, fmt="G = 1")
ax.set_yscale("log")
ax.scatter([5.0], [branching(p)], marker="x", s=70, c="k")
ax.scatter([2.0], [0.3], marker="*", s=120, c="k")
ax.set_xlabel(r"muon cost $C_\mu$ [GeV]")
ax.set_title(ttl, fontsize=9)
axes[0].set_ylabel(r"complex fusion branching $B_f$")
axes[0].annotate("today", (5.0, branching(p)), xytext=(5, 8),
textcoords="offset points", fontsize=8)
axes[1].annotate("Sec. 6 plant point", (2.0, 0.3), xytext=(6, -16),
textcoords="offset points", fontsize=8)
fig.suptitle("Fig 3 — Gain landscape at $\eta_c=0.6$: the same R1 axis, "
"with and without R5. Left: star point unreachable.", fontsize=10)
fig.colorbar(pc, ax=axes, label=r"$\log_{10} G$", fraction=0.03)
fig.savefig("/home/claude/out/fig3_gain_landscape.png", bbox_inches="tight")
def fig4_failure_propagation(p):
knobs = ["lam_fus", "lam_dec", "lam_form0", "phi", "c_He", "omega_s"]
labels = [r"$\lambda_{fus}$ (R1)", r"$\lambda_{dec}$ (R1)",
r"$\lambda_{form}$ (R4)", r"density $\phi$ (R4)",
r"$c_{He}$ (R4)", r"sticking $\omega_s$"]
base = dict(NOMINAL)
b_thru, b_Nf = complexes_per_muon(base), fusions_per_muon(base)
b_G = gain(b_Nf, 2000, 0.6)
rows = []
for k in knobs:
q = dict(base)
q[k] = q[k] * 2 if k != "omega_s" else min(2 * q[k], 0.05)
rows.append([complexes_per_muon(q) / b_thru,
fusions_per_muon(q) / b_Nf,
gain(fusions_per_muon(q), 2000, 0.6) / b_G])
rows = np.array(rows)
fig, ax = plt.subplots(figsize=(7.4, 4.4))
x = np.arange(len(knobs)); w = 0.25
ax.bar(x - w, rows[:, 0], w, label="complex throughput / base")
ax.bar(x, rows[:, 1], w, label=r"$N_f$ / base")
ax.bar(x + w, rows[:, 2], w, label="gain G / base")
ax.axhline(1.0, c="k", lw=0.8)
ax.set_xticks(x); ax.set_xticklabels(labels, fontsize=8)
ax.set_yscale("log"); ax.set_ylim(0.3, 4)
ax.set_ylabel("observable ratio after knob ×2")
ax.set_title("Fig 4 — Failure/leverage propagation (each knob doubled)\n"
"R1 knobs move yield+gain while leaving throughput fixed: "
"the bottleneck is branching, not kinetics")
ax.legend(fontsize=8)
fig.tight_layout(); fig.savefig("/home/claude/out/fig4_failure_propagation.png")
def fig5_mixture_optimum(p):
cgrid = np.logspace(-3, np.log10(0.6), 80)
Nf = []
for c in cgrid:
q = dict(p); q["c_He"] = c
Nf.append(fusions_per_muon(q))
Nf = np.array(Nf)
fig, ax = plt.subplots(figsize=(7, 4.2))
ax.semilogx(cgrid, Nf / Nf.max(), lw=2)
copt = cgrid[np.argmax(Nf)]
ax.axvline(copt, ls=":", c="gray")
ax.annotate(f"optimum $c_{{He}}\\approx{copt:.3f}$",
(copt, 1.0), xytext=(10, -12), textcoords="offset points")
ax.set_xlabel(r"$^3$He concentration $c_{He}$")
ax.set_ylabel(r"$N_f$ (normalized)")
ax.set_title("Fig 5 — Mixture trade-off: transfer speed vs. D$_2$ availability\n"
"(shallow optimum ⇒ R4 is forgiving; loop Stage-3 scan target)")
fig.tight_layout(); fig.savefig("/home/claude/out/fig5_mixture_optimum.png")
# ----------------------------------------------------------------------
if __name__ == "__main__":
import os
os.makedirs("/home/claude/out", exist_ok=True)
p = dict(NOMINAL)
print("=" * 64)
print("d-3He muCF simulation — nominal point")
print("=" * 64)
print(f"branching B_f = {branching(p):.3e}")
print(f"complexes per muon = {complexes_per_muon(p):.1f}")
print(f"fusions per muon N_f = {fusions_per_muon(p):.3e}")
for name, (Cmu, eta) in SCENARIOS.items():
print(f" gain [{name}] G = {gain(fusions_per_muon(p), Cmu, eta):.2e}"
f" (requires N_f > {required_Nf(Cmu, eta):.0f})")
lam_form = p["lam_form0"] * p["phi"] * (1 - p["c_He"])
print(f"\nDISCOVERED R5 (post-R1 ceiling): N_f_max = 1/(omega_eff + lam0/lam_form)")
print(f" formation-only ceiling (omega=0): {lam_form/LAM0:.0f} fusions/muon at LHD")
for name, (Cmu, eta) in SCENARIOS.items():
wmax = eta * Q_DHE3 / Cmu - LAM0 / lam_form
verdict = f"omega_eff <= {wmax:.2%}" if wmax > 0 else "UNREACHABLE at LHD (formation-limited)"
print(f" {name}: {verdict}")
print("\nAcceptance tests")
print("-" * 64)
tests = run_tests()
for k, v in tests.items():
print(f"{'PASS' if v['passed'] else 'FAIL'} {k}: value={v['value']}, window={v['window']}")
with open("/home/claude/out/acceptance_tests.json", "w") as f:
json.dump(tests, f, indent=2, default=str)
fig1_kinetics(p)
fig2_bf_bridge(p)
fig3_gain_landscape(p)
fig4_failure_propagation(p)
fig5_mixture_optimum(p)
print("\nFigures written to /home/claude/out/")