crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 06:05:58 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Merge pull request #2433 from paulromano/photon-from-ace-fix

Browse source 
Handle zero photon cross sections in `IncidentPhoton.from_ace`
This commit is contained in:
Ethan Peterson 2023-03-24 12:58:28 -04:00 • committed by GitHub
commit 7200d42644
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
2 changed files with 13 additions and 6 deletions
Download patch file Download diff file Expand all files Collapse all files

9
openmc/data/photon.py
View file

@ -514,13 +514,13 @@ class IncidentPhoton(EqualityMixin):
# Read each reaction
data = cls(Z)
for mt in (502, 504, 515, 522, 525):
for mt in (502, 504, 517, 522, 525):
data.reactions[mt] = PhotonReaction.from_ace(ace, mt)
# Get heating cross sections [eV-barn] from factors [eV per collision]
# by multiplying with total xs
data.reactions[525].xs.y *= sum([data.reactions[mt].xs.y for mt in
(502, 504, 515, 522)])
(502, 504, 517, 522)])
# Compton profiles
n_shell = ace.nxs[5]
@ -1000,7 +1000,7 @@ class PhotonReaction(EqualityMixin):
elif mt == 504:
# Incoherent scattering
idx = ace.jxs[1] + n
elif mt == 515:
elif mt == 517:
# Pair production
idx = ace.jxs[1] + 4*n
elif mt == 522:
@ -1021,6 +1021,9 @@ class PhotonReaction(EqualityMixin):
else:
nonzero = (xs != 0.0)
xs[nonzero] = np.exp(xs[nonzero])
# Replace zero elements to small non-zero to enable log-log
xs[~nonzero] = np.exp(-500.0)
rx.xs = Tabulated1D(energy, xs, [n], [5])
# Get form factors for incoherent/coherent scattering

10
src/photon.cpp
View file

@ -91,9 +91,13 @@ PhotonInteraction::PhotonInteraction(hid_t group)
close_group(rgroup);
// Read pair production
rgroup = open_group(group, "pair_production_electron");
read_dataset(rgroup, "xs", pair_production_electron_);
close_group(rgroup);
if (object_exists(group, "pair_production_electron")) {
rgroup = open_group(group, "pair_production_electron");
read_dataset(rgroup, "xs", pair_production_electron_);
close_group(rgroup);
} else {
pair_production_electron_ = xt::zeros_like(energy_);
}
// Read pair production
if (object_exists(group, "pair_production_nuclear")) {