From b9446ade6c5e51ae94726b96ff9ea80281c9e723 Mon Sep 17 00:00:00 2001 From: GuySten Date: Fri, 6 Mar 2026 19:26:07 +0200 Subject: [PATCH 1/2] hard fix --- src/geometry.cpp | 21 +++- src/particle.cpp | 19 ++- tests/unit_tests/test_energy_balance.py | 152 ++++++++++++++++++++++++ 3 files changed, 181 insertions(+), 11 deletions(-) create mode 100644 tests/unit_tests/test_energy_balance.py diff --git a/src/geometry.cpp b/src/geometry.cpp index ddb61385f18..e8e8a065112 100644 --- a/src/geometry.cpp +++ b/src/geometry.cpp @@ -173,12 +173,21 @@ bool find_cell_inner( } // Set the material, temperature and density multiplier. - p.material_last() = p.material(); - p.material() = c.material(p.cell_instance()); - p.sqrtkT_last() = p.sqrtkT(); - p.sqrtkT() = c.sqrtkT(p.cell_instance()); - p.density_mult_last() = p.density_mult(); - p.density_mult() = c.density_mult(p.cell_instance()); + if (p.sqrtkT() != -1) { + // if current material is set, set last material + p.material_last() = p.material(); + p.material() = c.material(p.cell_instance()); + p.sqrtkT_last() = p.sqrtkT(); + p.sqrtkT() = c.sqrtkT(p.cell_instance()); + p.density_mult_last() = p.density_mult(); + p.density_mult() = c.density_mult(p.cell_instance()); + } else { + // set both current data and last data to current data + p.material_last() = p.material() = c.material(p.cell_instance()); + p.sqrtkT_last() = p.sqrtkT() = c.sqrtkT(p.cell_instance()); + p.density_mult_last() = p.density_mult() = + c.density_mult(p.cell_instance()); + } return true; diff --git a/src/particle.cpp b/src/particle.cpp index 0a063559824..679b32fc62a 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -198,6 +198,8 @@ void Particle::event_calculate_xs() cell_last(j) = coord(j).cell(); } n_coord_last() = n_coord(); + + material_last() = material(); } // Write particle track. @@ -211,7 +213,7 @@ void Particle::event_calculate_xs() if (material() != MATERIAL_VOID) { if (settings::run_CE) { if (material() != material_last() || sqrtkT() != sqrtkT_last() || - density_mult() != density_mult_last()) { + density_mult() != density_mult_last() || true) { // If the material is the same as the last material and the // temperature hasn't changed, we don't need to lookup cross // sections again. @@ -397,10 +399,6 @@ void Particle::event_collide() // Save coordinates for tallying purposes r_last_current() = r(); - // Set last material to none since cross sections will need to be - // re-evaluated - material_last() = C_NONE; - // Set all directions to base level -- right now, after a collision, only // the base level directions are changed for (int j = 0; j < n_coord() - 1; ++j) { @@ -419,6 +417,14 @@ void Particle::event_collide() if (!model::active_tallies.empty()) score_collision_derivative(*this); + // Saving previous cell data + for (int j = 0; j < n_coord(); ++j) { + cell_last(j) = coord(j).cell(); + } + n_coord_last() = n_coord(); + + material_last() = material(); + #ifdef OPENMC_DAGMC_ENABLED history().reset(); #endif @@ -472,6 +478,8 @@ void Particle::event_revive_from_secondary() cell_last(j) = coord(j).cell(); } n_coord_last() = n_coord(); + + material_last() = material(); } pht_secondary_particles(); } @@ -696,6 +704,7 @@ void Particle::cross_reflective_bc(const Surface& surf, Direction new_u) // Reassign particle's cell and surface coord(0).cell() = cell_last(0); + material() = material_last(); surface() = -surface(); // If a reflective surface is coincident with a lattice or universe diff --git a/tests/unit_tests/test_energy_balance.py b/tests/unit_tests/test_energy_balance.py new file mode 100644 index 00000000000..d1c07b417b4 --- /dev/null +++ b/tests/unit_tests/test_energy_balance.py @@ -0,0 +1,152 @@ +import openmc +import pytest +import numpy as np + + +@pytest.fixture +def two_cell_model(): + """Simple two-cell slab model with a fixed source. + Cell1 occupies x in [-10, 0], cell2 x in [0, 10]. + """ + openmc.reset_auto_ids() + model = openmc.Model() + + m1 = openmc.Material() + m1.add_element('Li', 1.0) + m1.set_density('g/cm3', 1.0) + + m2 = openmc.Material() + m2.add_element('Al', 1.0) + m2.set_density('g/cm3', 1.0) + + xmin = openmc.XPlane(-10.0, boundary_type="reflective") + xmid = openmc.XPlane(0.0) + xmax = openmc.XPlane(10.0, boundary_type="reflective") + ymin = openmc.YPlane(-10.0, boundary_type="reflective") + ymax = openmc.YPlane(10.0, boundary_type="reflective") + zmin = openmc.ZPlane(-10.0, boundary_type="reflective") + zmax = openmc.ZPlane(10.0, boundary_type="reflective") + cell1 = openmc.Cell(fill=m1, region=+xmin & -xmid & +ymin & -ymax & +zmin & -zmax) + cell2 = openmc.Cell(fill=m2, region=+xmid & -xmax & +ymin & -ymax & +zmin & -zmax) + model.geometry = openmc.Geometry([cell1, cell2]) + + src = openmc.IndependentSource() + src.space = openmc.stats.Point((-5.0, 0.0, 0.0)) + src.particle = 'photon' + src.energy = openmc.stats.Discrete([1e6],[1.0]) + + model.settings.run_mode = 'fixed source' + model.settings.batches = 1 + model.settings.particles = 100 + model.settings.source = src + + return model, xmid, cell1, cell2, m1, m2 + +@pytest.fixture +def two_sphere_model(): + openmc.reset_auto_ids() + model = openmc.Model() + + mat = openmc.Material() + mat.add_nuclide('H2', 1.0) + mat.set_density('g/cm3', 1.0) + + surf1 = openmc.Sphere(r=0.01) + surf2 = openmc.Sphere(r=1000, boundary_type='vacuum') + cell1 = openmc.Cell(region=-surf1) + cell2 = openmc.Cell(fill=mat, region=+surf1 & -surf2) + + model.geometry = openmc.Geometry([cell1, cell2]) + + src = openmc.IndependentSource() + src.energy = openmc.stats.Discrete([1e6],[1.0]) + + model.settings.run_mode = 'fixed source' + model.settings.batches = 1 + model.settings.particles = 1000 + model.settings.source = src + + return model, cell1, cell2, mat + +def test_energy_balance_simple(two_cell_model, run_in_tmpdir): + model, xmid, cell1, cell2, *_ = two_cell_model + + tally1 = openmc.Tally() + tally1.scores = ['heating'] + tally1.filters = [openmc.CellFilter([cell1, cell2])] + + tally2 = openmc.Tally() + tally2.scores = ['current'] + tally2.filters = [openmc.EnergyFunctionFilter([0.0,20e6],[0.0,20e6]), openmc.SurfaceFilter([xmid])] + + + model.tallies = [tally1, tally2] + model.run(apply_tally_results=True) + + assert tally1.mean.sum() == pytest.approx(1e6) + + assert tally2.mean[0] == pytest.approx(tally1.mean[1]) + +def test_current_conservation(two_cell_model, run_in_tmpdir): + model, xmid, cell1, cell2, m1, m2 = two_cell_model + + energyfunc_filter = openmc.EnergyFunctionFilter([0.0,20e6],[0.0,20e6]) + + tally1 = openmc.Tally() + tally1.scores = ['current'] + tally1.filters = [energyfunc_filter, openmc.SurfaceFilter([xmid])] + + tally2 = openmc.Tally() + tally2.scores = ['current'] + tally2.filters = [energyfunc_filter, + openmc.CellFromFilter([cell1]), + openmc.CellFilter([cell2])] + + tally3 = openmc.Tally() + tally3.scores = ['current'] + tally3.filters = [energyfunc_filter, + openmc.CellFromFilter([cell2]), + openmc.CellFilter([cell1])] + + tally4 = openmc.Tally() + tally4.scores = ['current'] + tally4.filters = [energyfunc_filter, + openmc.MaterialFromFilter([m1]), + openmc.MaterialFilter([m2])] + + tally5 = openmc.Tally() + tally5.scores = ['current'] + tally5.filters = [energyfunc_filter, + openmc.MaterialFromFilter([m2]), + openmc.MaterialFilter([m1])] + + + model.tallies = [tally1, tally2, tally3, tally4, tally5] + model.run(apply_tally_results=True) + + assert pytest.approx(tally1.mean.sum()) == tally2.mean[0] + tally3.mean[0] + + assert tally2.mean[0] == pytest.approx(tally4.mean[0]) + + assert tally3.mean[0] == pytest.approx(tally5.mean[0]) + +def test_cellfrom_heating(run_in_tmpdir, two_sphere_model): + + model, cell1, cell2, mat = two_sphere_model + + tally1 = openmc.Tally() + tally1.scores = ['heating'] + tally1.filters = [openmc.CellFromFilter([cell1, cell2]), openmc.CellFilter([cell2])] + + tally2 = openmc.Tally() + tally2.scores = ['heating'] + tally2.filters = [openmc.MaterialFromFilter([mat]), openmc.MaterialFilter([mat])] + + model.tallies = [tally1, tally2] + + model.run(apply_tally_results=True) + + assert np.all(tally1.mean > 0) + + assert tally1.mean[1] == pytest.approx(tally2.mean[0]) + From 47e1d1c3ba39d27dc5cc1f6cbdf6a72a7b7cf498 Mon Sep 17 00:00:00 2001 From: GuySten Date: Fri, 6 Mar 2026 20:16:21 +0200 Subject: [PATCH 2/2] update test --- .../filter_cellfrom/results_true.dat | 28 +++++++++---------- 1 file changed, 14 insertions(+), 14 deletions(-) diff --git a/tests/regression_tests/filter_cellfrom/results_true.dat b/tests/regression_tests/filter_cellfrom/results_true.dat index 408a4965b1d..1f9d0446104 100644 --- a/tests/regression_tests/filter_cellfrom/results_true.dat +++ b/tests/regression_tests/filter_cellfrom/results_true.dat @@ -1,29 +1,29 @@ k-combined: 9.035025E-02 2.654309E-03 tally 1: -5.994069E+00 -3.594398E+00 +6.586649E+00 +4.340013E+00 tally 2: 4.559707E-04 2.080739E-08 tally 3: -5.994525E+00 -3.594945E+00 +6.587105E+00 +4.340613E+00 tally 4: 0.000000E+00 0.000000E+00 tally 5: -1.473892E+00 -2.187509E-01 +8.813120E-01 +7.829296E-02 tally 6: 5.223875E-05 2.992861E-10 tally 7: -1.473945E+00 -2.187663E-01 +8.813643E-01 +7.830214E-02 tally 8: -1.885798E+01 -3.558423E+01 +5.908648E+00 +3.492172E+00 tally 9: 7.467961E+00 5.580255E+00 @@ -34,8 +34,8 @@ tally 11: 7.468470E+00 5.581014E+00 tally 12: -1.885798E+01 -3.558423E+01 +5.908648E+00 +3.492172E+00 tally 13: 0.000000E+00 0.000000E+00 @@ -46,8 +46,8 @@ tally 15: 2.739543E-04 7.600983E-09 tally 16: -7.881296E+01 -6.221087E+02 +9.176229E+01 +8.428689E+02 tally 17: 1.051397E+02 1.106292E+03