Variance reduction

MAGIC weight windows and alias-table source sampling in Nucleide.

Nucleide provides two variance-reduction utilities built on MCNP mesh tallies: MAGIC weight-window generation, weight-window emission for OpenMC and Serpent, and mesh-based source sampling with alias tables.

MAGIC weight windows

MAGIC converts a forward mesh tally into weight-window lower bounds. The algorithm is simple and deterministic:

For each energy group $g$, find the maximum flux over all volume elements:

Then the weight-window lower bound for volume element $v$ and group $g$ is

where $\varepsilon_{v,g}$ is the relative error of the tally and $\tau$ is the error tolerance. The factor of $2$ normalizes the peak-flux voxel to a lower bound of $0.5$; MAGIC is intended to be iterated so the weight windows refine over successive forward runs (Cooper & Larsen 2001)^cooper-2001.

Total vs. per-group mode

Parameters

Worked example

Consider a tally with four volume elements and total fluxes $\Phi = [1.2,\, 3.3,\, 1.6,\, 1.7]$ and relative errors all below the default tolerance of $0.5$.

1. Find the group maximum: $\Phi^{\max} = 3.3$.

2. Compute the scale factor: $1 / (2 \times 3.3) = 1 / 6.6$.

3. The lower bounds are $\text{ww} = [1.2/6.6,\, 3.3/6.6,\, 1.6/6.6,\, 1.7/6.6] = [0.182,\, 0.5,\, 0.242,\, 0.258]$.

If the tolerance were tightened to $0.15$, the element with error $0.19$ would be nulled to null_value instead.

Emitting weight windows for OpenMC and Serpent

MAGIC produces lower bounds. Two functions reformat the same MagicOutput for other transport codes — pure formatting over the meshtal mesh bounds, with no new math:

• emit_openmc_weight_windows(tally, output, ...) returns a settings.xml fragment: a <mesh> element (id, dimension, lower_left, upper_right) plus a <weight_windows> element that references it (mesh, particle_type, energy_bounds, lower_ww_bounds, upper_ww_bounds, survival_ratio, max_split, weight_cutoff). The element spelling matches what OpenMC’s own writer emits and its transport code parses (settings specification §3.29 and §3.66). Energies are converted from MeV to eV, and the window values are written with the x index varying fastest inside each energy group — the same ordering OpenMC uses internally. Paste the fragment inside the existing <settings> root.
• emit_serpent_wwin(tally, output, name, file) returns the text of a weight-window file plus the input card that references it. Serpent reads third-party window meshes in the MCNP WWINP text format through wwin <name> wf "<file>" 2 (Serpent user guide §2.2.8.3), so that is the spelling the emitter writes; the file re-parses cleanly with the workspace WWINP reader. Serpent’s own first file format is not publicly documented, so it is not produced — requesting it raises a clear error.

MAGIC computes lower bounds only. OpenMC requires upper bounds as well, so the emitter derives them as five times the lower bound (the same ratio OpenMC applies when importing MCNP files) and says so in its notes; Serpent derives its upper bounds itself at load time. Inputs that cannot be represented — negative or non-finite bounds, energy or mesh bounds that are not increasing, or tuning values outside the ranges the code accepts — raise a clear error instead of writing a partial file. Cells nulled by MAGIC (zero lower bounds) are written as zeros and called out in the notes.

Alias-table source sampling

The MeshSourceSampler uses a Walker/Vose alias table to draw particle birth voxels in $O(1)$ time per sample after an $O(n)$ construction step.

Alias-table construction

Given a discrete probability mass function $p_i$ over $n$ bins:

1. Normalize so that $\sum_i p_i = 1$.
2. Scale each probability by $n$, giving $n \bar{p}_i$.
3. Partition bins into “small” ($n\bar{p}_i < 1$) and “large” ($n\bar{p}_i \geq 1$).
4. Pair a small bin $a$ with a large bin $g$:
   • Set $\text{prob}[a] = n\bar{p}_a$.
   • Set $\text{alias}[a] = g$.
   • Reduce the large bin’s surplus: $n\bar{p}_g \leftarrow n\bar{p}_g + n\bar{p}_a - 1$.
5. Drain any remaining bins with $\text{prob} = 1$.

Sampling

To draw one index:

1. Generate $r_1, r_2 \in [0, 1)$.

2. Pick column $i = \lfloor n \, r_1 \rfloor$.

3. Return $i$ if $r_2 < \text{prob}[i]$; otherwise return $\text{alias}[i]$.

Mesh source sampling modes

Assumptions and limitations

• MAGIC requires positive flux somewhere: if every value in an energy group is $\leq 0$, the normalization divides by zero and Error::ZeroMaxFlux is returned.
• MAGIC rejects non-finite flux or error values: NaN or infinite entries in the input arrays produce Error::NonFiniteTally rather than poisoned lower bounds.
• Mesh source sampling rejects negative values: negative total_result or user_pdf entries produce Error::NegativeTally; zero values remain valid zero-probability bins.
• Error comparison is strict: a relative error exactly equal to the tolerance keeps its scaled value; only errors greater than $\tau$ are nulled.
• Alias tables are exact for the given PDF: the construction reproduces the input probabilities to machine precision.
• Voxel resolution: birth coordinates are returned as integer voxel indices; users must sample a uniform point inside the voxel themselves if sub-voxel positions are needed.
• Emission carries lower bounds only: upper bounds are derived (OpenMC) or code defaults (Serpent), and any structure the target format cannot spell is a clear error rather than a silent approximation.

Related work beyond Nucleide

Adjacent capabilities found in PyNE, OpenMC, and the literature that Nucleide does not implement:

• MAGIC iteration: the published method iterates — each forward run refines the windows for the next — and derives upper bounds from a fixed ratio to the lower bounds; Nucleide performs a single lower-bound pass.
• FW-CADIS and adjoint-based global variance reduction: hybrid deterministic/Monte Carlo methods that build weight windows from an adjoint importance solve; see Wagner et al. (2014)^wagner-2014.
• Built-in generators in OpenMC: WeightWindows and WeightWindowGenerator (MAGIC-like and FW-CADIS modes, updated on the fly during transport), including an internal random-ray adjoint solver.
• Source biasing: OpenMC supports energy/space/angle biasing of independent sources; PyNE’s mesh sampler additionally supports sub-voxel sampling modes. Nucleide’s MeshSourceSampler is voxel-level and energy-integrated.