whaeck
commented
1 day ago I wanted to reply to this but other things needed my attention. I looked through this and it looks solid to me. Interestingly enough, this issue was brought to my attention by another person around the same time you put up the PR. There's probably more of these things in the source code. It might have made sense over 30 years ago, but not today.
I'll run a few more tests over the next few weeks and will most likely merge this into a feature branch on this repository, prior to merging it into the develop branch - which is the staging area for new releases.
I found in my PhD studies that the Chi^2 sampling scheme in PURR is pretty severely biased. It uses a pre-tabulated set of twenty chi^2 samples for one through four degrees of freedom. Obviously, choosing these at random leads to a mean of the sample mean rather than the true mean of 1, 2, 3, and 4.
This PR is partly just to make this known.
I tried generating a whole nuclear data library using this patch and couldn't pick up on any appreciable differences between the original version of the code and this one. This appears to be the case because despite the samples being biased, the widths sampled according to Chi^2 R.V.'s always appear as ratios, so the ratio terms themselves are not biased. Moreover, using the ptables as a stochastic dilution factor for pointwise grid is also normalizing out potential biases. Perhaps this would explain some of the historical disagreement between NJOY and PREPRO as Red Cullen has pointed out before, suggesting that the samples should be normalized and treated as a multiplying factor. Of course, there are other reasons for this such as decusping the interpolation scheme. I can't find the paper at the moment but will find it if anyone requests.
If you plot the CDF of the cross section distribution, you can observe for many nuclides that PURR is missing the tails by quite a bit. The small set of samples essentially fails to hit rare events where the widths are very large or very small, and this governs the tails of the cross-section distribution.
Here is a plot of U235's XS multiplier term at its first URR energy point. You can see the tails are cut off by the default implementation.
Here is the U238 total XS multiplier at the first URR energy point with and without this fix. Again it is a small effect but visible in the tails.
maybe there would be some discernible difference for deep shielding problems of neutrons in the URR range, but aside from that, I'm not sure this makes much of a difference.
I'm just raising this PR for posterity to see, in case someone else ever gets worried about the strange Chi^2 sampling method in PURR.