Context
PR 3771 adds the ability to tally microscopic cross sections in void materials, but only works for neutrons. Il tested it for photon dose calculations (silicon vs biological dose in void) and confirmed photons don't work. Two problems need fixing:
- Void path is neutron-only: The atom_density blocks in the tracklength/collision estimators only call
update_neutron_xs() — no photon equivalent exists.
- Latent nuclide-vs-element index bug:
score_general_ce_nonanalog calls p.photon_xs(i_nuclide) where i_nuclide is a nuclide index, but photon_xs_ is indexed by element index. This is wrong when nuclide index != element index (common with multiple isotopes).
Changes
1. Add global nuclide-to-element index mapping
include/openmc/nuclide.h — declare in namespace data:
extern vector<int> nuclide_element;
src/nuclide.cpp — define and clear:
- Add
vector<int> data::nuclide_element;
- In
nuclides_clear() (line 1218): add data::nuclide_element.clear();
src/simulation.cpp — build mapping after init_nuclide_index() loop (after line 116):
if (settings::photon_transport) {
data::nuclide_element.resize(data::nuclides.size(), C_NONE);
for (size_t i = 0; i < data::nuclides.size(); ++i) {
auto element_name = to_element(data::nuclides[i]->name_);
auto it = data::element_map.find(element_name);
if (it != data::element_map.end()) {
data::nuclide_element[i] = it->second;
}
}
}
(Already includes openmc/nuclide.h, openmc/photon.h; needs openmc/string_utils.h)
2. Add helper function for photon XS cache update
src/tallies/tally_scoring.cpp — add in anonymous namespace near top of file:
//! Ensure photon XS cache is populated for the element corresponding to a
//! nuclide. Returns element index, or C_NONE if unavailable.
int update_photon_xs_for_nuclide(Particle& p, int i_nuclide)
{
int i_element = data::nuclide_element[i_nuclide];
if (i_element != C_NONE) {
auto& micro = p.photon_xs(i_element);
if (p.E() != micro.last_E) {
data::elements[i_element]->calculate_xs(p);
}
}
return i_element;
}
3. Fix p.photon_xs(i_nuclide) index bug in score_general_ce_nonanalog
At the top of the scoring loop (after line 583), compute element index once:
int i_element = (i_nuclide >= 0 && p.type().is_photon())
? data::nuclide_element[i_nuclide] : C_NONE;
Replace 4 occurrences of p.photon_xs(i_nuclide) with p.photon_xs(i_element):
- Line 603 (SCORE_TOTAL)
- Line 627 (SCORE_SCATTER)
- Line 647 (SCORE_ABSORPTION)
- Line 1050 (COHERENT/INCOHERENT/PHOTOELECTRIC/PAIR_PROD)
4. Add photon void path in score_tracklength_tally_general (~line 2434)
Current code skips the entire block when p.material() == MATERIAL_VOID. Add an else branch for void that handles both neutrons and photons:
if (i_nuclide >= 0) {
if (p.material() != MATERIAL_VOID) {
// ... existing material path (unchanged) ...
// In the j == C_NONE sub-path, also add photon branch alongside neutron
} else {
// MATERIAL_VOID
if (!tally.multiply_density()) {
if (p.type().is_neutron()) {
// existing neutron log-union + update_neutron_xs logic
} else if (p.type().is_photon()) {
if (update_photon_xs_for_nuclide(p, i_nuclide) != C_NONE)
atom_density = 1.0;
}
}
}
}
Also add the photon branch in the existing j == C_NONE sub-path (nuclide not in material, non-void) for consistency.
5. Add photon void path in score_collision_tally (~line 2566)
Same pattern as step 4. Also wrap existing material access in a MATERIAL_VOID guard (currently missing — would crash if reached in void).
6. Add tests
tests/unit_tests/test_tally_multiply_density.py — add test function:
- Photon source in void sphere, photon_transport enabled
- Tally with nuclide filter (e.g., 'Fe56'),
multiply_density=False
- Photon scores: 'total', 'coherent', 'incoherent', 'photoelectric', 'pair-production'
- Assert all mean values > 0 in void cell
Files to modify
| File |
Change |
include/openmc/nuclide.h |
Declare data::nuclide_element |
src/nuclide.cpp |
Define + clear in nuclides_clear() |
src/simulation.cpp |
Build mapping (+ add string_utils.h include) |
src/tallies/tally_scoring.cpp |
Helper function, fix 4x photon_xs index bug, add photon void paths in 2 estimator functions |
tests/unit_tests/test_tally_multiply_density.py |
Add photon void test |
Verification
- Build OpenMC with changes
- Run existing test:
pytest tests/unit_tests/test_tally_multiply_density.py — must still pass
- Run new photon void test — must pass (non-zero photon XS scores in void)
- Run full tally test suite:
pytest tests/unit_tests/ -k tally — no regressions
Context
PR 3771 adds the ability to tally microscopic cross sections in void materials, but only works for neutrons. Il tested it for photon dose calculations (silicon vs biological dose in void) and confirmed photons don't work. Two problems need fixing:
update_neutron_xs()— no photon equivalent exists.score_general_ce_nonanalogcallsp.photon_xs(i_nuclide)wherei_nuclideis a nuclide index, butphoton_xs_is indexed by element index. This is wrong when nuclide index != element index (common with multiple isotopes).Changes
1. Add global nuclide-to-element index mapping
include/openmc/nuclide.h— declare innamespace data:src/nuclide.cpp— define and clear:vector<int> data::nuclide_element;nuclides_clear()(line 1218): adddata::nuclide_element.clear();src/simulation.cpp— build mapping afterinit_nuclide_index()loop (after line 116):(Already includes
openmc/nuclide.h,openmc/photon.h; needsopenmc/string_utils.h)2. Add helper function for photon XS cache update
src/tallies/tally_scoring.cpp— add in anonymous namespace near top of file:3. Fix
p.photon_xs(i_nuclide)index bug inscore_general_ce_nonanalogAt the top of the scoring loop (after line 583), compute element index once:
Replace 4 occurrences of
p.photon_xs(i_nuclide)withp.photon_xs(i_element):4. Add photon void path in
score_tracklength_tally_general(~line 2434)Current code skips the entire block when
p.material() == MATERIAL_VOID. Add anelsebranch for void that handles both neutrons and photons:Also add the photon branch in the existing
j == C_NONEsub-path (nuclide not in material, non-void) for consistency.5. Add photon void path in
score_collision_tally(~line 2566)Same pattern as step 4. Also wrap existing material access in a
MATERIAL_VOIDguard (currently missing — would crash if reached in void).6. Add tests
tests/unit_tests/test_tally_multiply_density.py— add test function:multiply_density=FalseFiles to modify
include/openmc/nuclide.hdata::nuclide_elementsrc/nuclide.cppnuclides_clear()src/simulation.cppstring_utils.hinclude)src/tallies/tally_scoring.cppphoton_xsindex bug, add photon void paths in 2 estimator functionstests/unit_tests/test_tally_multiply_density.pyVerification
pytest tests/unit_tests/test_tally_multiply_density.py— must still passpytest tests/unit_tests/ -k tally— no regressions