From e120de61e4797875562ebb4aff854065169838fb Mon Sep 17 00:00:00 2001 From: Claude Date: Fri, 2 Oct 2026 06:11:23 +0000 Subject: [PATCH 1/3] Use direction of travel when looking up weight window mesh bin When a geometry surface coincides with a weight window mesh plane, the surface checkpoint looks up the weight window with the particle sitting exactly on the plane. Structured mesh index lookups resolve an exact boundary position to the lower-coordinate element, so particles crossing in the +x/+y/+z direction were assigned the window of the element they were leaving. This produced asymmetric splitting/rouletting and a direction-dependent loss of FOM (reported on Discourse for a gamma shielding problem using ADVANTG weight windows whose mesh planes match the geometry). Nudge the lookup position by TINY_BIT along the direction of travel, as is already done for tracklength mesh tallies, so the particle is assigned to the element it is entering. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_014RN7JroEkrbkLKVHbMDap9 --- src/weight_windows.cpp | 9 ++- .../test_ww_surface_checkpoint.py | 81 +++++++++++++++++++ 2 files changed, 88 insertions(+), 2 deletions(-) create mode 100644 tests/unit_tests/weightwindows/test_ww_surface_checkpoint.py diff --git a/src/weight_windows.cpp b/src/weight_windows.cpp index 114531bf754..47f52af6ca6 100644 --- a/src/weight_windows.cpp +++ b/src/weight_windows.cpp @@ -8,6 +8,7 @@ #include "openmc/tensor.h" +#include "openmc/constants.h" #include "openmc/error.h" #include "openmc/file_utils.h" #include "openmc/hdf5_interface.h" @@ -289,9 +290,13 @@ std::pair WeightWindows::get_weight_window( if (E < energy_bounds_.front() || E > energy_bounds_.back()) return {false, {}}; - // Get mesh index for particle's position + // Get mesh index for particle's position. The position is nudged along the + // direction of travel so that a particle sitting exactly on a mesh boundary + // (e.g., at a surface checkpoint where a geometry surface coincides with a + // mesh plane) is assigned to the element it is entering rather than always + // to the element on the lower-coordinate side. const auto& mesh = this->mesh(); - int mesh_bin = mesh->get_bin(p.r()); + int mesh_bin = mesh->get_bin(p.r() + TINY_BIT * p.u()); // particle is outside the weight window mesh if (mesh_bin < 0) diff --git a/tests/unit_tests/weightwindows/test_ww_surface_checkpoint.py b/tests/unit_tests/weightwindows/test_ww_surface_checkpoint.py new file mode 100644 index 00000000000..49ea09f1af4 --- /dev/null +++ b/tests/unit_tests/weightwindows/test_ww_surface_checkpoint.py @@ -0,0 +1,81 @@ +import openmc +import pytest + + +@pytest.mark.parametrize("direction", [1.0, -1.0]) +def test_ww_surface_on_mesh_boundary(run_in_tmpdir, direction): + """Weight windows applied at a surface that coincides with a weight window + mesh boundary must use the mesh element the particle is entering, + regardless of the direction of travel.""" + + # Void-like one-group material so the only weight window checks happen at + # surface crossings + groups = openmc.mgxs.EnergyGroups([0.0, 20.0e6]) + xsdata = openmc.XSdata('void', groups) + xsdata.order = 0 + xsdata.set_total([0.0]) + xsdata.set_absorption([0.0]) + xsdata.set_scatter_matrix([[[0.0]]]) + mg_library = openmc.MGXSLibrary(groups) + mg_library.add_xsdata(xsdata) + mg_library.export_to_hdf5('mgxs.h5') + + mat = openmc.Material() + mat.add_macroscopic('void') + materials = openmc.Materials([mat]) + materials.cross_sections = 'mgxs.h5' + + # Two cells separated by the plane x = 0, which is also a mesh boundary + x_min = openmc.XPlane(-1.0, boundary_type='vacuum') + x_mid = openmc.XPlane(0.0) + x_max = openmc.XPlane(1.0, boundary_type='vacuum') + yz = openmc.model.RectangularPrism(2.0, 2.0, axis='x', + boundary_type='vacuum') + left = openmc.Cell(fill=mat, region=+x_min & -x_mid & -yz) + right = openmc.Cell(fill=mat, region=+x_mid & -x_max & -yz) + geometry = openmc.Geometry([left, right]) + + # Particle is born in the window of its birth element and must be split + # five ways when it enters the element on the other side of x = 0 + birth_cell, target_cell = (left, right) if direction > 0 else (right, left) + mesh = openmc.RectilinearMesh() + mesh.x_grid = [-1.0, 0.0, 1.0] + mesh.y_grid = [-1.0, 1.0] + mesh.z_grid = [-1.0, 1.0] + birth_bounds = (0.5, 1.5) + target_bounds = (0.1, 0.2) + if direction > 0: + lower = [birth_bounds[0], target_bounds[0]] + upper = [birth_bounds[1], target_bounds[1]] + else: + lower = [target_bounds[0], birth_bounds[0]] + upper = [target_bounds[1], birth_bounds[1]] + ww = openmc.WeightWindows(mesh, lower_ww_bounds=lower, + upper_ww_bounds=upper) + + settings = openmc.Settings() + settings.energy_mode = 'multi-group' + settings.run_mode = 'fixed source' + settings.particles = 10 + settings.batches = 1 + settings.source = openmc.IndependentSource( + space=openmc.stats.Point((-0.5*direction, 0.0, 0.0)), + angle=openmc.stats.Monodirectional((direction, 0.0, 0.0)), + ) + settings.weight_windows = ww + settings.weight_window_checkpoints = {'collision': False, 'surface': True} + + tally = openmc.Tally() + tally.filters = [ + openmc.CellFilter([target_cell]), + openmc.WeightFilter([0.0, 0.5, 2.0]) + ] + tally.scores = ['flux'] + + model = openmc.Model(geometry, materials, settings, openmc.Tallies([tally])) + model.run(apply_tally_results=True) + + # All flux in the target cell should be carried by split particles with + # weight 0.2; total flux is conserved (1 cm path length per source particle) + flux = tally.mean.squeeze() + assert flux == pytest.approx([1.0, 0.0]) From 57fefbc340d6adb294b9dbb2d848ed08ea41431e Mon Sep 17 00:00:00 2001 From: Claude Date: Fri, 2 Oct 2026 20:03:44 +0000 Subject: [PATCH 2/3] Explain in the surface checkpoint test why the +x case failed Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_014RN7JroEkrbkLKVHbMDap9 --- .../weightwindows/test_ww_surface_checkpoint.py | 12 +++++++++++- 1 file changed, 11 insertions(+), 1 deletion(-) diff --git a/tests/unit_tests/weightwindows/test_ww_surface_checkpoint.py b/tests/unit_tests/weightwindows/test_ww_surface_checkpoint.py index 49ea09f1af4..911de6ade2a 100644 --- a/tests/unit_tests/weightwindows/test_ww_surface_checkpoint.py +++ b/tests/unit_tests/weightwindows/test_ww_surface_checkpoint.py @@ -6,7 +6,17 @@ def test_ww_surface_on_mesh_boundary(run_in_tmpdir, direction): """Weight windows applied at a surface that coincides with a weight window mesh boundary must use the mesh element the particle is entering, - regardless of the direction of travel.""" + regardless of the direction of travel. + + At the surface checkpoint the particle sits exactly on the plane x = 0. + Structured mesh index lookups resolve a position exactly on an element + boundary to the element on the lower-coordinate side, so when the window + was looked up at the particle's position alone, a particle moving in +x + was given the window of the element it was leaving. It was then not split + and the +x case failed, while the -x case passed only because there the + lower-coordinate element happens to be the one being entered. The lookup + now nudges the position along the direction of travel. + """ # Void-like one-group material so the only weight window checks happen at # surface crossings From 58c4fca5d859876f39e23222571b817376b699b3 Mon Sep 17 00:00:00 2001 From: Claude Date: Fri, 2 Oct 2026 20:27:20 +0000 Subject: [PATCH 3/3] Update reference results affected by weight window lookup fix Both tests apply weight windows at surface checkpoints on geometry surfaces that coincide with weight window mesh planes. In the survival biasing test, the reflective planes at x, y, z = 0 lie on the lower edge of the mesh, so particles reflecting off them were previously found outside the mesh and no window was applied; in the bootstrap test, a cavity face lies on an interior mesh plane. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_014RN7JroEkrbkLKVHbMDap9 --- .../results_true.dat | 16 ++++++++-------- .../survival_biasing/local/results_true.dat | 2 +- .../survival_biasing/shared/results_true.dat | 2 +- 3 files changed, 10 insertions(+), 10 deletions(-) diff --git a/tests/regression_tests/random_ray_auto_convert_bootstrap/results_true.dat b/tests/regression_tests/random_ray_auto_convert_bootstrap/results_true.dat index c32200aeff1..92dd0f70bce 100644 --- a/tests/regression_tests/random_ray_auto_convert_bootstrap/results_true.dat +++ b/tests/regression_tests/random_ray_auto_convert_bootstrap/results_true.dat @@ -1,9 +1,9 @@ tally 1: -8.094106E+01 -6.974487E+02 -3.181906E-04 -1.455701E-08 -9.840422E+02 -1.110377E+05 -8.515372E+01 -8.078508E+02 +8.087496E+01 +7.339588E+02 +2.446860E-04 +8.405655E-09 +9.761237E+02 +1.018073E+05 +9.815520E+01 +1.065938E+03 diff --git a/tests/regression_tests/weightwindows/survival_biasing/local/results_true.dat b/tests/regression_tests/weightwindows/survival_biasing/local/results_true.dat index 8a4d3368cc4..b8f32dc2f9f 100644 --- a/tests/regression_tests/weightwindows/survival_biasing/local/results_true.dat +++ b/tests/regression_tests/weightwindows/survival_biasing/local/results_true.dat @@ -1 +1 @@ -0d7b17d4e364bda2be21feb2e5b2aae4be69fc8ed634c4f45062e8efc61b7bd24dcf448837692fd603afe3756501fe8821d134f2d6289179618f85178ddb8793 \ No newline at end of file +7a9d895b430f0c7696fc36f73dc68179d484664a29bd168df3b28a2b0549552d6216e4e9b32de4449d8134a42a84883d7e7e7364f0b9e4cc75174e95ecaf6b1c \ No newline at end of file diff --git a/tests/regression_tests/weightwindows/survival_biasing/shared/results_true.dat b/tests/regression_tests/weightwindows/survival_biasing/shared/results_true.dat index 2e88e2f0e81..f99e3c14d80 100644 --- a/tests/regression_tests/weightwindows/survival_biasing/shared/results_true.dat +++ b/tests/regression_tests/weightwindows/survival_biasing/shared/results_true.dat @@ -1 +1 @@ -2f1c8efbf6ab0b0c7319fadd7ee8f27b545c1ebab94ebe15cf1133b54227487542966b6d798970c8b37bc19e7041e687687d5285803241bb83f4c8dc40a4d276 \ No newline at end of file +5ee91c9d9285829528fb2e53a1c5bf069318259e0d4e8768852309a921bc6dd77852de251bd6e91245f68bc33b0fdc2321cfc3b315432d822603976f74079396 \ No newline at end of file