Hybrid Multigroup Transport

This fixed-source problem combines native H-1 composition with one-group neutron multigroup data on the same material. The 1 eV source lies inside the multigroup interval from 0.1 to 10 eV, while the 1 MeV source uses native H-1 physics outside that interval. Multigroup scattering represents the outgoing group at its 5.05 eV midpoint. The tally uses physical energy boundaries because particle energy remains in eV during hybrid transport. The example requires MCDC_LIB to identify a native-data library containing H1-293.6K.h5.

Full Input

import numpy as np

import mcdc

simulation = mcdc.Simulation("Hybrid multigroup sphere")

# Material with native H-1 data and a low-energy multigroup treatment
# MCDC_LIB must contain H1-293.6K.h5.
material = mcdc.Material(
    nuclide_composition={"H1": 5.0e-2},
    neutron_multigroup=mcdc.NeutronMultigroupData(
        capture=np.array([0.2]),
        scatter=np.array([[0.8]]),
        speed=np.array([1.383e6]),
        energy_grid=np.array([0.1, 10.0]),
        energy_representation="midpoint",
    ),
)

# Spherical domain
boundary = mcdc.Surface.Sphere(radius=5.0, boundary_condition="vacuum")
cell = mcdc.Cell(region=-boundary, fill=material)
simulation.set_model([cell])

# The low-energy source uses multigroup physics
multigroup_source = mcdc.Source(
    position=[0.0, 0.0, 1.0],
    isotropic=True,
    energy=1.0,
    probability=0.5,
)

# The high-energy source lies outside the multigroup grid and uses native physics
native_source = mcdc.Source(
    position=[0.0, 0.0, -1.0],
    isotropic=True,
    energy=1.0e6,
    probability=0.5,
)
simulation.set_sources([multigroup_source, native_source])

# Hybrid tallies use physical energy boundaries in eV
mesh = mcdc.MeshStructured(z=np.array([-5.0, 0.0, 5.0]))
tally = mcdc.Tally(
    name="hybrid_flux",
    mesh=mesh,
    particle_type="neutron",
    energy=np.array([0.0, 0.1, 10.0, 1.0e5, 20.0e6]),
    scores=["flux", "collision"],
)
simulation.set_tallies([tally])

simulation.settings.N_particle = 1_000
simulation.settings.N_batch = 5
simulation.settings.output_name = "hybrid_multigroup"

simulation.run()

Post-processing

import h5py
import matplotlib.pyplot as plt
import numpy as np

with h5py.File("hybrid_multigroup.h5", "r") as output:
    tally = output["tallies/hybrid_flux"]
    z = tally["grid/z"][:]
    energy = tally["grid/energy"][:]
    flux = tally["flux/mean"][:]

dz = z[1:] - z[:-1]
z_midpoint = 0.5 * (z[:-1] + z[1:])
flux = np.reshape(flux, (len(energy) - 1, len(z) - 1)) / dz

figure, axis = plt.subplots()
for index in range(len(energy) - 1):
    axis.plot(
        z_midpoint,
        flux[index],
        marker="o",
        label=f"{energy[index]:.1e}–{energy[index + 1]:.1e} eV",
    )

axis.set_xlabel("z [cm]")
axis.set_ylabel("Flux")
axis.grid()
axis.legend()
figure.tight_layout()
figure.savefig("hybrid_multigroup_flux.png", dpi=150)

How to Run

From inside examples/hybrid_multigroup:

python input.py
python process-output.py

The transport calculation writes hybrid_multigroup.h5 and the post-processing script writes hybrid_multigroup_flux.png.