Skip to content

Fix MaterialFromFilter dropping scores after the first collision in a cell - #4090

Open
GuySten wants to merge 2 commits into
openmc-dev:developfrom
GuySten:material-filter-fix
Open

Fix MaterialFromFilter dropping scores after the first collision in a cell#4090
GuySten wants to merge 2 commits into
openmc-dev:developfrom
GuySten:material-filter-fix

Conversation

@GuySten

@GuySten GuySten commented Aug 31, 2026

Copy link
Copy Markdown
Contributor

Fix MaterialFromFilter dropping scores after the first collision in a cell

This changes volumetric MaterialFromFilter results. Tallies using it were
previously missing scores; they will now be larger. Cell-to-cell partial
currents (a flux score with this filter) are unaffected. See below for why.

Problem

material_last was doing two unrelated jobs:

  1. the tally attribute read by MaterialFromFilter, and
  2. the cache key for the cross section lookup in event_calculate_xs().

Because a particle's energy changes at every collision, event_collide() forced
a cache miss by setting material_last() = C_NONE. From the second collision in
a cell onward the particle therefore carried material_last == C_NONE, which
matches no filter bin, so those scores silently vanished from volumetric
MaterialFromFilter tallies. Cell-to-cell partial currents (a flux score with
one of these filters) are scored at the crossing itself, where the previous
material was still correct, and were not affected.

cell_last is not reset on collision, so the equivalent CellFromFilter tally
counted those scores, and the two decompositions of the same score disagreed.

The effect is largest where it is least visible: a thick, scattering-dominated
cell loses nearly everything after the first collision.

MWE

import openmc
 
SHELL_THICKNESS = 10.0  # cm of iron; raise it to make the shortfall worse
 
water = openmc.Material()
water.set_density('g/cm3', 1.0)
water.add_nuclide('H1', 2.0)
water.add_nuclide('O16', 1.0)
 
iron = openmc.Material()
iron.set_density('g/cm3', 7.87)
iron.add_nuclide('Fe56', 1.0)
 
inner_surf = openmc.Sphere(r=5.0)
outer_surf = openmc.Sphere(r=5.0 + SHELL_THICKNESS, boundary_type='vacuum')
inner = openmc.Cell(fill=water, region=-inner_surf)
outer = openmc.Cell(fill=iron, region=+inner_surf & -outer_surf)
 
model = openmc.Model()
model.geometry = openmc.Geometry([inner, outer])
 
model.settings.run_mode = 'fixed source'
model.settings.batches = 10
model.settings.particles = 2000
model.settings.source = openmc.IndependentSource(
    space=openmc.stats.Point((0.0, 0.0, 0.0)),
    energy=openmc.stats.delta_function(2.0e6),
)
 
in_outer = openmc.CellFilter([outer])
 
total = openmc.Tally(name='total')
total.filters = [in_outer]
total.scores = ['total']
 
by_cell = openmc.Tally(name='by cell')
by_cell.filters = [in_outer, openmc.CellFromFilter([inner, outer])]
by_cell.scores = ['total']
 
by_material = openmc.Tally(name='by material')
by_material.filters = [in_outer, openmc.MaterialFromFilter([water, iron])]
by_material.scores = ['total']
 
model.tallies = openmc.Tallies([total, by_cell, by_material])
 
with openmc.StatePoint(model.run(output=False)) as sp:
    t = sp.get_tally(name='total').mean.sum()
    c = sp.get_tally(name='by cell').mean.sum()
    m = sp.get_tally(name='by material').mean.sum()
 
print(f"\n  reaction rate in outer cell        {t:12.6f}")
print(f"  summed over CellFromFilter        {c:12.6f}"
      f"   ({100 * c / t:6.2f}% of total)")
print(f"  summed over MaterialFromFilter    {m:12.6f}"
      f"   ({100 * m / t:6.2f}% of total)")
 
print("\n  Both decompositions should equal the total exactly.")
if abs(m - t) > 1e-9 * t:
    print(f"  MaterialFromFilter is short by {100 * (1 - m / t):.1f}%.")
else:
    print("  Both match.")

Changes

  • material_last is renamed to material_xs_cache, which is what it actually
    is. It is written in exactly the same places as before, so cross section
    caching behaviour is unchanged.
  • MaterialFromFilter now derives the material from cell_last at the lowest
    coordinate level instead of tracking it separately. The two "from" filters
    read one source of truth and cannot drift apart again.
  • Adds cell_instance_last, maintained on the same cadence as cell_last. A
    cell with distributed materials resolves a different material per instance, so
    the cell index alone does not determine which material a particle came from.
  • Documents the semantics of cell_last, cell_instance_last and
    material_xs_cache on their accessors, including that the first two are
    deliberately not updated on collision and that the last must never be used as
    a tally attribute.

This deliberately leaves event_calculate_xs() alone. Rerouting the cache key
was the alternative; deriving the material instead means there is no cross
section behaviour to re-verify.

Testing

tests/unit_tests/test_filter_material_from.py:

  • A water sphere inside an iron shell, one material per cell so the two
    decompositions are the same partition. Asserts both sum to the undecomposed
    total and agree bin for bin. The material sum falls well short of the total
    before this change.
  • A 2x2 lattice of one repeated cell carrying four distributed materials,
    asserting all four "from" bins are populated. Resolving the fill without the
    instance would put the whole total in the first bin.
  • A cell-to-cell partial current (flux score), asserting the material and cell
    decompositions still agree. cell_last and cell_instance_last are saved at
    the top of event_cross_surface(), before the crossing, so the derived filter
    sees the same thing the tracked attribute used to.

Both are per-history identities within a single run, so they hold to round-off
rather than statistically.

tests/regression_tests/surface_tally uses MaterialFromFilter with a flux
score, which makes it a partial current rather than a volumetric tally, so its
stored results should be unchanged. A unit test covers that path so the change
is shown not to regress it.

Notes

material_last was stale at birth as well: a fresh particle inherited whatever
the reused Particle object last held, since find_cell() writes the
pre-existing material before overwriting it. Deriving from cell_last, which is
initialised to the particle's own cells at birth, removes that too.

This is the material_last portion of #3849, split out as suggested in review
there. It does not change cell_last semantics.

Checklist

  • I have performed a self-review of my own code
  • I have run clang-format (version 18) on any C++ source files (if applicable)
  • I have followed the style guidelines for Python source files (if applicable)
  • I have made corresponding changes to the documentation (if applicable)
  • I have added tests that prove my fix is effective or that my feature works (if applicable)

@GuySten
GuySten requested a review from paulromano August 31, 2026 17:54
@GuySten
GuySten marked this pull request as ready for review August 31, 2026 18:45
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant