Random Ray Eigenvalue Flux Normalization Change

[jtramm]

Problem Description

Currently the random ray solver does not take much care to normalize the flux solution when in eigenvalue mode. This is because the flux is the eigenvector, so its magnitude is not particularly meaningful. It is currently left up to the user to normalize their eigenvalue fluxes to whatever level they want/need.

The lack of normalization causes a few problems:

1. Tallies may be reported at different magnitudes depending on settings, leading to confusion if results appear to be changing (despite the shape being consistent).

2. Eigenvalue fluxes are not directly comparable to Monte Carlo tally solutions generated by OpenMC, due to the normalization differences.

It is therefore highly preferable to normalize in a consistent manner

The Fix

Normalization in a consistent manner is easy. However, normalization so as to match Monte Carlo (which lists flux in terms of "per starting source particle") is trickier, as in random ray we are natively computing flux "per cm of particle travel". Thus, we will naively be off by a factor equivalent to the mean free path of a neutron in that problem, which is also not trivial to compute.

It looks like a good way is to instead normalize the random ray fluxes to the total number of fissions produced in the system, such that fluxes are tabulated "per fission produced". We can then multiple the fluxes by the eigenvalue to convert these quantities into the expected "per source particle". To test, I used a 1D slab problem (shared by Paul Ferney and @GiudGiud of INL - thanks guys!). The problem has some fuel at the far left, then some air, then a large shield wall, then some more air.

To test the new PR, I ran multigroup Monte Carlo to generate a reference solution given a set of XS data, then ran the current branch random ray solver and then the random ray solver from this PR. As can be seen, the random ray solver now matches very closely with the Monte Carlo solver. The Monte Carlo solver flux gives zero once it gets a few tens of cm into the shield wall due to a lack of weight windows in these tests, but clearly the fluxes are now normalized correctly.

Notably, the shapes of the random ray flux distributions haven't changed, but the magnitude now is consistent with MC so will make for a much nicer user experience.

These changes resulted in scaling of many tally values, so all the k-eff tests for random ray had to be updated.

Added Support for Fission Particle Disabling

I also added in support for the existing settings::create_fission_neutrons setting, which can be set to False so as to prevent both the Monte Carlo and now random ray solvers from creating secondary neutrons (or computing the fission source) when in fixed source mode.

Checklist

• I have performed a self-review of my own code
• I have run clang-format (version 15) on any C++ source files (if applicable)
• I have followed the style guidelines for Python source files (if applicable)
• I have made corresponding changes to the documentation (if applicable)
• I have added tests that prove my fix is effective or that my feature works (if applicable)

[Copilot]

Pull Request Overview

This PR implements consistent eigenvalue flux normalization for the random ray solver to address scaling issues and improve compatibility with Monte Carlo solutions. The changes normalize random ray flux to "per fission produced" and then scale by the eigenvalue to match Monte Carlo's "per starting source particle" convention.

• Implements flux normalization in eigenvalue mode to match Monte Carlo solver output magnitude
• Adds support for the settings::create_fission_neutrons flag in random ray mode
• Updates all regression test reference values due to the flux scaling changes

Reviewed Changes

Copilot reviewed 11 out of 11 changed files in this pull request and generated 1 comment.

File	Description
src/random_ray/flat_source_domain.cpp	Implements flux normalization logic and fission neutron creation flag support
src/random_ray/linear_source_domain.cpp	Adds fission neutron creation flag support for linear source domains
Multiple test result files	Updates reference values to reflect the new flux normalization

---

_{Tip: Customize your code reviews with copilot-instructions.md. Create the file or learn how to get started.}

[paulromano]

Thanks for this PR @jtramm. I have one question below about the normalization factor. Also, what happens for fixed source mode? I'm assuming in that case it's easier to normalize since we know the source strength?

[paulromano]

Does this make sense?

norm_factor = k_eff_new / total_fission_neutrons
= (k_eff_ * (fission_rate_new / fission_rate_old)) / (fission_rate_new * simulation_volume_)
= k_eff_ / (fission_rate_old * simulation_volume_)

[jtramm]

It does seem a little odd. I wonder if (effectively) scaling by some property of the previous iteration as I'm proposing creates extra noise, whereas perhaps there is a way that only uses the current iteration's properties and thus might have lower variance? It seems this would be inevitable though as k_new is inherently a function of the ratio of fluxes between two batches, so you have to incorporate previous batch info into this. I suppose k_new / F_new = k_old / F_old may be in inherent property of power iteration and eigenvalues & the magnitude of eigenvectors (F)? Any other ideas for scaling the fluxes?

[jtramm]

Fixed source mode uses a separate method for normalizing fluxes, which is more easily done, as we can pretty easily compute the total external source term present in the simulation and normalize against the user specified total strength. The random ray fixed source solver has been computing fluxes with the same magnitudes as MGMC for awhile now, it is just k-eff mode that was still an issue.