OpenMC - within Cardinal - to model neutronics for Greisheimer and Kooreman Analytical Benchmark (PHYSOR 2022)

Greisheimer and Kooreman present a multiphysics Analytical Benchmark that I would like to use OpenMC for the neutronics modeling, within Cardinal. The paper uses deterministic (S_{2}) transport, but it should be possible to use Monte Carlo for this purpose. Since the paper assumes 1D transport with vacuum conditions on the endpoints, OpenMC can handle this by using a vacuum boundary condition on each x boundary and reflective boundary conditions on the y and z boundaries. Additionally, the mesh will be N x 1 x 1, so there is only refinement in the x direction, which will allow for computing \phi(x).

The part that I am looking for feedback on is how to handle the cross sections. This paper assumes a one-group structure for the cross sections, but it also couples some temperature and density feedbacks: namely Doppler Broadening and thermal expansion (in the x direction). Equation (10) in the paper shows that the total macroscopic cross section can be written as \Sigma_{t}(x) = \frac{\Sigma_{t,0}T_{0}}{T(x)}.

Another assumption of the paper is that the temperature and neutron flux have the same shape, or T(x)=f\phi(x). This allows the total macroscopic cross to be written also as \Sigma_{t}(x) = \frac{\Sigma_{t,0}T_{0}}{f\phi(x)}. The paper also explains how to determine \Sigma_{t,0} from other chosen parameters, so the question I have is this: how can I generate the one-group cross sections as a function of position (or temperature and/or flux)? It seems like this may require some iteration and some initial guess for the distribution. The paper gives a procedure for the selecting the parameters, and even gives some canonical ones as a way to compare. I am going to use the canonical parameters in this benchmark, so this will determine \Sigma_{t,0} and a few other parameters that are derived from initial choices. For more on the paper and details of the benchmark, see my GitHub repository for this project.

For cross section variation as a function of position, it seems like the simplest thing to do is to just have a different material / MG cross section specified for each region in the x direction. For temperature variation, I don’t think there’s any way you could represent the cross section as a smooth function of T. Instead, you would have to discretize T into a suitably fine temperature mesh and then provide a different cross section for each temperature. It looks like our XSdata class can handle multiple temperatures, but I’ve never tried it myself. Maybe @agnelson can chime in with some advice.