openmc-dev / openmc

OpenMC Monte Carlo Code
https://docs.openmc.org
Other
748 stars 482 forks source link

Rare bug for a randomly generated nuclear data file #1699

Open makeclean opened 3 years ago

makeclean commented 3 years ago

We've (through Ander Gray) come across a bug which results in a failure to sample an inelastic reaction. In this instance, its due to the fact that at the the elastic scattering cross section is 0 at the energy of the particle, the logic flow then breaks down. The original initiating issue is that with implicit capture on, the in physics.cpp

if (settings::survival_biasing) {
    // Determine weight absorbed in survival biasing
    p.wgt_absorb_ = p.wgt_ * p.neutron_xs_[i_nuclide].absorption /
          p.neutron_xs_[i_nuclide].total;

In this instance neutron_xs_[i_nuclide].absorption == neutron_xs_[i_nuclide].total and thus the p.wgt_absorb_ is the entire particle weight and p.wgt_ is 0. The check if the particle is dead happens only in the non-implicit capture branch. So my question is, is this truely a bug? I feel like a nuclide physically can't have an elastic scattering cross section of 0 can it? Which points to a nuclear data processesing error upstream.

paulromano commented 3 years ago

Hmm, yeah I'm no physicist but I seem to recall from my studies of neutron interaction theory that elastic scattering is always possible. There have been problems (some fairly recent) in evaluations of negative elastic scattering cross sections, which processing codes typically zero out or set to some very low value. In light of what you're seeing, obviously a very small positive value would cause less problems, but even so I feel like we should probably add some kind of check here. Can you guys share / point to the evaluation that is leading to this issue?

makeclean commented 3 years ago

problem_child.tar.gz

Here is the full problem, the particular nuclide causing the problem is ni60

paulromano commented 3 years ago

Yikes, I see the 2 keV wide hole in that cross section. Is the original ENDF file using SLBW resonance formalism? In any event, I think we should ought to fix this "bug". One possible solution that should prevent problems whether or not implicit capture is on would be to check for a non-zero scattering cross section before calling scatter, something like

if (p.neutron_xs_[i_nuclide].total > p.neutron_xs_[i_nuclide].absorption) {
  scatter(p, i_nuclide);
}

where currently there is no conditional around the call to scatter.

AnderGray commented 3 years ago

Hi both, the perturbed files produced by sandy if helps randomEndf.zip