04 — Divertor Heat Flux & Blanket Neutronics

This tutorial demonstrates two nuclear engineering modules:

1. Divertor Thermal Simulation — heat-flux spreading on the divertor target
2. Blanket Neutronics — tritium breeding ratio (TBR) calculation

Both are critical for DEMO-class reactor feasibility.

License: © 1998–2026 Miroslav Šotek. GNU AGPL v3.

---

import numpy as np
import matplotlib.pyplot as plt

Part A: Divertor Thermal Simulation

The divertor receives the exhaust heat from the plasma scrape-off layer. ITER's divertor must handle ~10 MW/m² steady-state, with transients up to 20 MW/m².

from scpn_fusion.core.divertor_thermal_sim import DivertorThermalSim

# ITER-like divertor parameters
sim = DivertorThermalSim(
    target_length_m=1.2,
    target_width_m=0.3,
    n_points=100,
    coolant_temp_K=423.0,  # 150°C water
)

# Sweep over power levels
powers = np.linspace(5.0, 25.0, 20)  # MW/m²
peak_temps = []
for q in powers:
    result = sim.compute_steady_state(heat_flux_mw_m2=q)
    peak_temps.append(result['peak_temp_K'])

fig, ax = plt.subplots(figsize=(8, 5))
ax.plot(powers, np.array(peak_temps) - 273.15, 'b-o', markersize=4)
ax.axhline(y=1800, color='r', linestyle='--', label='Tungsten recrystallization (1800°C)')
ax.axhline(y=3422, color='darkred', linestyle='--', label='Tungsten melting (3422°C)')
ax.set_xlabel('Heat Flux [MW/m²]')
ax.set_ylabel('Peak Surface Temperature [°C]')
ax.set_title('Divertor Target Temperature vs Heat Flux')
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

print(f"At 10 MW/m²: {peak_temps[5]-273.15:.0f}°C")
print(f"At 20 MW/m²: {peak_temps[15]-273.15:.0f}°C")

Part B: Blanket Neutronics (Tritium Breeding)

A D-T fusion reactor must breed its own tritium. The Tritium Breeding Ratio (TBR) must exceed 1.0 for fuel self-sufficiency. ITER targets TBR ≈ 1.10.

from scpn_fusion.nuclear.blanket_neutronics import BlanketNeutronics

neutronics = BlanketNeutronics()

# Compare blanket concepts
concepts = {
    'HCLL (Li-Pb)':       {'li6_enrichment': 0.30, 'blanket_thickness_m': 0.8},
    'HCPB (Li-ceramic)':  {'li6_enrichment': 0.60, 'blanket_thickness_m': 0.7},
    'WCLL (Water-cooled)': {'li6_enrichment': 0.40, 'blanket_thickness_m': 0.75},
}

results = {}
for name, params in concepts.items():
    tbr = neutronics.compute_tbr(**params)
    results[name] = tbr
    print(f"{name:25s}: TBR = {tbr:.3f} {'✓' if tbr > 1.0 else '✗'}")

# Plot
fig, ax = plt.subplots(figsize=(8, 4))
bars = ax.barh(list(results.keys()), list(results.values()), color=['#2196F3', '#4CAF50', '#FF9800'])
ax.axvline(x=1.0, color='r', linestyle='--', linewidth=2, label='TBR = 1.0 (self-sufficiency)')
ax.set_xlabel('Tritium Breeding Ratio')
ax.set_title('Blanket Concept TBR Comparison')
ax.legend()
plt.tight_layout()
plt.show()

Benchmark

%%timeit -n1 -r3
for _name, _params in concepts.items():
    neutronics.compute_tbr(**_params)

Part C: Li-6 Enrichment Sensitivity

TBR depends strongly on Li-6 enrichment. Natural lithium is only 7.5% Li-6.

enrichments = np.linspace(0.075, 0.90, 50)
tbr_values = [neutronics.compute_tbr(li6_enrichment=e, blanket_thickness_m=0.8) for e in enrichments]

fig, ax = plt.subplots(figsize=(8, 5))
ax.plot(enrichments * 100, tbr_values, 'g-', linewidth=2)
ax.axhline(y=1.0, color='r', linestyle='--', label='TBR = 1.0')
ax.axhline(y=1.1, color='orange', linestyle='--', alpha=0.7, label='TBR = 1.1 (target)')
ax.axvline(x=7.5, color='gray', linestyle=':', label='Natural Li-6 (7.5%)')
ax.set_xlabel('Li-6 Enrichment [%]')
ax.set_ylabel('Tritium Breeding Ratio')
ax.set_title('TBR vs Li-6 Enrichment (0.8m blanket)')
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

Performance Benchmarks

Timing the key computations in this notebook:

1. Single divertor steady-state solve (compute_steady_state)
2. Heat flux sweep (20-point parameter scan)
3. Single TBR calculation (compute_tbr)
4. Li-6 enrichment sweep (50-point parameter scan)

import timeit

# 1. Single divertor steady-state solve
def bench_divertor_single():
    sim.compute_steady_state(heat_flux_mw_m2=10.0)

t_div = timeit.repeat(bench_divertor_single, number=100, repeat=5)
print(f"compute_steady_state single call (100 calls):")
print(f"  Mean: {np.mean(t_div)*1000:.1f} ms +/- {np.std(t_div)*1000:.1f} ms")
print(f"  Per call: {np.mean(t_div)/100*1000:.2f} ms")

# 2. Heat flux sweep (20-point parameter scan)
def bench_divertor_sweep():
    for q in np.linspace(5.0, 25.0, 20):
        sim.compute_steady_state(heat_flux_mw_m2=q)

t_sweep = timeit.repeat(bench_divertor_sweep, number=10, repeat=3)
print(f"\nDivertor 20-point heat flux sweep (10 runs):")
print(f"  Mean: {np.mean(t_sweep)*1000:.1f} ms +/- {np.std(t_sweep)*1000:.1f} ms")
print(f"  Per run: {np.mean(t_sweep)/10*1000:.2f} ms")

# 3. Single TBR calculation
def bench_tbr_single():
    neutronics.compute_tbr(li6_enrichment=0.30, blanket_thickness_m=0.8)

t_tbr = timeit.repeat(bench_tbr_single, number=1000, repeat=5)
print(f"\ncompute_tbr single call (1000 calls):")
print(f"  Mean: {np.mean(t_tbr)*1000:.1f} ms +/- {np.std(t_tbr)*1000:.1f} ms")
print(f"  Per call: {np.mean(t_tbr)/1000*1e6:.2f} us")

# 4. Li-6 enrichment sweep (50-point parameter scan)
def bench_tbr_sweep():
    for e in np.linspace(0.075, 0.90, 50):
        neutronics.compute_tbr(li6_enrichment=e, blanket_thickness_m=0.8)

t_tbr_sweep = timeit.repeat(bench_tbr_sweep, number=10, repeat=3)
print(f"\nTBR 50-point enrichment sweep (10 runs):")
print(f"  Mean: {np.mean(t_tbr_sweep)*1000:.1f} ms +/- {np.std(t_tbr_sweep)*1000:.1f} ms")
print(f"  Per run: {np.mean(t_tbr_sweep)/10*1000:.2f} ms")

Summary

• The divertor model computes 1D/2D temperature profiles under steady-state and transient heat loads
• The neutronics model estimates TBR from Li-6 enrichment, blanket geometry, and neutron multiplication
• Both modules feed into the compact reactor optimizer (see notebook 01) for integrated design
• With Rust acceleration, these calculations can be embedded in real-time digital twin scenarios