Add distributed cell densities

[nuclearkevin]

Description

This PR adds a new distributed cell quantity called a density multiplier which allows for the modification of per-cell densities (cell density = cell multiplier * material density) to better account for multiphysics feedback. This is an improvement compared to the previous approach used in coupling codes (e.g. Cardinal) where the user would need to add a new material per cell to get density feedback (as feedback could only be specified on a per-material basis), which is highly labour intensive for most reactor systems. The approach taken here is similar to the method used to get per-cell temperatures: two doubles are added to the GeometryState class which tracks the current cell density multiplier and the previous cell density multipliers. These are used to correct the sampled material density during cross section calculations, distance to next collision calculations, and tally scoring. The density multipliers are initialized to unity, they are set indirectly with the following OpenMC C-API functions (necessary for Cardinal):

int openmc_cell_set_density(int32_t index, double rho, const int32_t* instance, bool set_contained = false);

int openmc_cell_get_density(int32_t index, const int32_t* instance, double* rho);

A new attribute for the Cell class has been added in the Python API called density which is treated almost identically to temperature. This is exported to the cell's xml element and used to set the cell density multiplier on the C++ side by dividing by the material density (in g/cc). Cell densities can be imported / exported as properties and loaded from xml elements for backups / restarts with multiphysics data.

Closes #2900

cc'ing @aprilnovak and @lewisgross1296

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)

[pshriwise]

Thanks @nuclearkevin! This one's been long awaited I know. Some initial thoughts here.

[nuclearkevin]

Thanks for the initial review @pshriwise! I've addressed your comments, added documentation for the new openmc.lib and C-API functions, and added unit tests for setting / getting density multipliers.

I've been looking around to see if cell temperatures have any regression tests (so I can mirror them for cell density multipliers), but I can't seem to find any. Are per-cell temperatures regression tested at all in OpenMC, or is that left up to coupling applications?

[pshriwise]

Thanks for the initial review @pshriwise! I've addressed your comments, added documentation for the new openmc.lib and C-API functions, and added unit tests for setting / getting density multipliers.

Nice! Thanks for the quick response and the additional work!

I've been looking around to see if cell temperatures have any regression tests (so I can mirror them for cell density multipliers), but I can't seem to find any. Are per-cell temperatures regression tested at all in OpenMC, or is that left up to coupling applications?

There isn't a dedicated test for this via the Python API, but the functionality in transport is tested in the cpp_driver tests. There may be others, I'd have to poke around more.

Thought: Rather than setting density multipliers on materials, what if we densities directly, but store them as multipliers internally relative to the first entry. I think this would make for a more intuitive API as the inputs/outputs will have physical meaning.

[nuclearkevin]

Thought: Rather than setting density multipliers on materials, what if we densities directly, but store them as multipliers internally relative to the first entry. I think this would make for a more intuitive API as the inputs/outputs will have physical meaning.

It would make it quite a bit nicer to have a set / get cell density function rather then a density multiplier. The reason I decided against going that route is the possibility for unit confusion. Unlike temperature where we can reasonably expect people to be OK with the units being Kelvin, there are plenty of valid density unit combinations that are used depending on the code using the C-API. I elected to force the user to handle the unit conversion on their side before setting the multiplier in OpenMC instead of handling the mess of possible units, but if we're OK with only allowing a single density unit (e.g. g/cc) I can change the API over.

Thinking about this some more, I'll change the API over so users specify a density (in g/cm^3). It's easier to understand, and this way people know what they need to convert to.

[paulromano]

Thanks @nuclearkevin. Given that this has been commonly requested by users, it would be desirable to have an interface from the Python API. I could see two ways of going about this:

1. Add a density_multiplier attribute to Cell that gets passed as-is to the functionality you already have here.
2. Add a density attribute to Cell that overrides the density of the fill material. This would result in a multiplier on the C++ side equal to the ratio of the density attribute to the density of the material fill.

I tend to like the second option because it mirrors what we do for temperature (cell temperatures override material temperatures), but open to others' thoughts here.

[nuclearkevin]

@paulromano sure thing! I'm also in favour of option number two as well since that's more intuitive for users, and it's similar to how cell densities are specified here in the C-API.

[nuclearkevin]

@paulromano I've added a pure Python API equivalent for cell densities and a series of Python regression & unit tests for the capability. I believe this PR is ready for the next round of reviews now.

It also appears as if the activation test is still rather flaky.

[pshriwise]

Looking good @nuclearkevin!

A few more comments from me here, but I think this is close from my perspective!

[pshriwise]

A few small items here. But also, some additional thoughts on design for your consideration.

Upon first looking at this PR, it was nice to see that the density multiplier is applied in a reasonably small number of places. I'd like the option/intention that multipliers can be optionally applied to material atom densities to be part of the Material class itself. The implementation has been nagging at me a bit and I'm curious about your thoughts on the following options.

1. Create a new method for obtaining atom_densities that reflects the intention to support multipliers.

double Material::atom_density(int i, double rho_multiplier=1.0) {
  return atom_density_(i) * rho_multiplier;
}

(mixing declaration and definition here a bit I know)

2. Adjust the Material::atom_density_ values to be relative to the material reference density (specific density?) by adding the following loop at the end of Material::normalize_density

  // set specific/relative atom densities
  for (int i=0; i < atom_density_.size(); i++) {
    atom_density_(i) /= density_gpcc_;
  }

(and probably update the name of that attribute accordingly)

Then the Material::atom_density method would look like

// return density at nominal value
double Material::atom_density(int i) {
  return atom_density_(i) * density_gpcc_;
}

// return atom density for a provided macroscopic density (gpcc)
double Material::atom_density(int i, double density) {
  return atom_density_(i) * density;
}

This would allow for storage of absolute densities on the Cell class, which is great for debugging (if needed) and file I/O, but also potentially complicates higher level code or results in a more complicated internal function.

I think in the end I prefer option 1 b/c it still reflects this capability on the Material object itself, separates data from implementation, and doesn't affect the higher level code as much as option 2 would, but I'm curious about your thoughts @kevinm387.

[nuclearkevin]

@pshriwise I'm on board with making density multiplier support more obvious from the Material side, and my preference would be the first option as well (as it's minimally invasive).

[pshriwise]

A couple of minor items here to tie up, but I'm pretty happy with where it's at.

[pshriwise]

LGTM! Thanks for your work on this @nuclearkevin!