depletion-thermochemistry: Redox control transfer rates

[church89]

Description

This PR adds a redox control (initially thought for molten salt reactors) functionality to the TransferRates class, as a further term added to the Bateman equation to maintain the redox of the fuel salt constant during a depletion simulation, by continuously adding or removing a user-defined buffer nuclides vector.
In uranium fluoride molten salts , for example, assuming all uranium is UF4 (therefore assuming oxidation state of Uranium +4), 4 atoms of fluorine are free every time a fission event occurs, which will then bind with fission products with lower oxidation states than uranium (with an average oxidation state of approx. +3), giving rise to an excess of free fluorine anions ($F^{-}$), a cationic charge deficiency and therefore a reduction of the fuel salt (decrease of redox potential).

Assuming the oxidation states of the nuclides present in the salt are known, the average oxidation state $\bar{OS}$ can be computed multiplying the number of atoms $N_i$ of each nuclide present in the salt, by their oxidation states $OS_i$:

$$\begin{equation} \bar{OS} = \sum_i N_i OS_i \end{equation}$$

A generic buffer nuclide $N_j$ can be used to compensate the exact change of fluorine atoms and a further term added to the Bateman equation as such:

$$\frac{dN_j(t)}{dt} = \cdots - \frac{1}{OS_j}\sum_i N_i a_{ij} \cdot OS_i$$

where the first terms of the equation have been omitted.
$OS$ is the oxidation state vector and $a_{ij}$ the corresponding term of the depletion matrix.

Example

# Define fuel as UF4
uf4 = openmc.Material(name='UF4')
uf4.set_density('g/cm3', 4.5)
uf4.add_element('U', 1., enrichment=2.4)
uf4.add_element('F', 4.)

...

# my guess of oxidation states 
ox = {
   'H': +1, 'He': 0, 'Li': +1, 'Be': +2, 'B': +3,  'C': 0,   'N': +3, 'O': -2,
   'F': -1, 'Ne': 0, 'Na': +1, 'Mg': +2, 'Al': +3, 'Si': +4, 'P': -3, 'S': -2,
   'Cl':-1, 'Ar': 0, 'K':  +1, 'Ca': +2, 'Sc': +3, 'Ti': +4, 'V': +5, 'Cr':+3,
   'Mn':+2, 'Fe':+3, 'Co': +2, 'Ni': +2, 'Cu': +2, 'Zn': +2, 'Ga':+3, 'Ge':+2,
   'As':+3, 'Se':-2, 'Br': -1, 'Kr':  0, 'Rb': +1, 'Sr': +2, 'Y': +3, 'Zr':+4,
   'Nb':+6, 'Mo':+4, 'Tc': +4, 'Ru': +2, 'Rh': +3, 'Pd': +2, 'Ag':+1, 'Cd':+2,
   'In':+3, 'Sn':+2, 'Sb': +3, 'Te': +2, 'I':  -1, 'Xe':  0, 'Cs':+1, 'Ba':+2,
   'La':+3, 'Ce':+3, 'Pr': +3, 'Nd': +3, 'Pm': +3, 'Sm': +3, 'Eu':+3, 'Gd':+3,
   'Tb':+3, 'Dy':+3, 'Ho': +3, 'Er': +3, 'Tm': +3, 'Yb': +3, 'Lu':+3, 'Hf':+4,
   'Ta':+5, 'W': +4, 'Re': +4, 'Os': +4, 'Ir': +3, 'Pt': +2, 'Au':+3, 'Hg':+1,
   'Tl':+1, 'Pb':+2, 'Bi': +3, 'Po': +2, 'At': -1, 'Rn':  0, 'Fr':+1, 'Ra':+2,
   'Ac':+2, 'Th':+4, 'Pa': +4, 'U':  +4, 'Np': +4, 'Pu': +3, 'Am':+3, 'Cm':+3,
   'Bk':+3, 'Cf':+3, 'Es': +3, 'Fm': +3, 'Md': +3, 'No': +2, 'Lr':+3
   }
# Add transfer rates for main absorbers to maintain an almost constant burnup 
integrator.add_transfer_rate('uf4', ['Xe', 'Kr'], 0.1)
# Use U238 as buffer for redox control
integrator.add_redox(uf4, {'U238':1}, ox)
integrator.integrate()

Results

In case of no redox control, there is an increase of fluorine atoms that can be compensated by adding buffer U238 to the fuel salt:

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)

[gridley]

Why should this be added to OpenMC? IMO it is a little too specific and unrelated to particle transport. There are a variety of ways to solve this particular problem, potentially using the Gibbs energy minimization approach, for example.

What I hope we could add is a way to modify the conditions of the depletion run so that very specific code like this can be kept in an external package. For example it would be nice to make changes through the C API over the course of depletion. As a matter of fact this is something I've been looking into lately, for example doing branch calculations at each depletion point.

Considering how many reactor concepts are out there, IMO we should make the code general so that specifics like this can be handled on the user side.

[church89]

Hi @gridley and thanks for your comment.

I agree this can be solved in many different ways, but this implementation relies on recent TransferRates implementation to add a new term to the Bateman equation directly, i.e. modifying the nuclides vector continuously, without acting step-wise at each depletion point.

I've also pulled another PR #2693 to perform batchwise actions during a depletion calculation, if is somehow related to your branch implementation, we should probably combine the efforts.

[gridley]

After looking at the code, I now see that this is a single up-front step rather than something done continuously over the depletion run.

So, it even further solidifies my opinion that this should go into an external package for two main reasons.

1. The code mainly just modifies the TransferRates object, and it doesn't require being integrated into the OpenMC ecosystem directly.
2. The assumptions in this depletion method are not exact. For example uranium in fluouride salts exists mainly in a IV oxidation state but can be III (albeit in more limited quantities) just as easily depending on the composition.

I think we should just wait for someone else's opinion on the matter here though. Perhaps my thinking is totally out of line with how the OpenMC project should be steered. After all, why not add all the bells and whistles?

[church89]

@gridley thanks for getting back to this.

In regards to your concerns:

1. What I've added here is a further term to the Bateman equation that balances the redox of the material by adjusting the composition of a user-defined buffer nuclide. Therefore it is taken into account when solving the depletion matrix and thus I consider it as a continuous action, am I missing something?

2. The method I use here to calculate the redox potential is through oxidation states, which they are user defined. If you want to set your initial uranium valence state to 4 or for example 3.99 (in presence of some UF3) is totally up to you. As you pointed out there are other methods to do it, but through valence states is particularly handy here.

[gridley]

@church89 Well, I have been thinking about this problem for a while. I think it sounds like a useful feature for a small subpopulation of openmc users.

Are there any published papers showing that this works? Or showing the impacts of tracking the redox? Back in my day studying this stuff, I found that adding buffer has a pretty negligible impact on neutronics for the most part, although I understand why it'd be nice to be extra sure and check this.

One thing that bothers me is that there doesn't seem to be a straightforward way to track the total amount of buffer added over time, at the moment. Am I correct in thinking the current code does not do that? It would be essential for checking the mass balance is reasonable.

Lastly, I'm a little confused about one thing. The buffer addition rate is controlled by the loss or gain rate of nuclides in a material, which in general might represent values that change over time. So, this matrix has to get re-formed at every single depletion substep, correct? It isn't just using the 1rst step's reaction rates, right?

[church89]

Hi @gridley thanks again for taking a look at this:

1. Not that I know of, I derived this myself, but it's a simple oxidation state balance equation. Haven't really tested many use cases, but I can tell that for certain type of liquid fueled reactors where fuel reprocessing is a thing this might be an important feature to have.
2. Correct, at the moment it is not accounted for, but wouldn't be hard to add this bit.
3. That's right, the matrix gets re-formed at every depletion step by simply recomputing gain and loss with oxidation states.

[gridley]

Makes sense to me @church89! This is a pretty neat feature... I wonder if there's any experimental data for MSRE, as to how much buffer was added to the salt over time. Wasn't it just uranium metal that was added in the sampler bucket?

Also, it appears one pretty essential piece of machinery is missing. Don't you think we should have the addition rate and total amount of buffer added over time as an output? Currently this approach is going to continuously add mass over time, but appears to not account for the volume increase that comes along with that. It would be nice to check the mass balance.

I assume that mass density is frozen. This will make it look like other nuclides are vanishing if you add in material and account for it in the isotopic breakdown but do not adjust density to preserve total mass.

[paulromano]

@church89 Can you get this branch up-to-date (looks like there are conflicts currently)? Apologies this one has lingered!

[paulromano]

@church89 Here's an initial review for this PR. I don't see any major blocking issues here so just a few minor things to take care of. It looks like your branch has gotten out-of-sync again with the develop branch, so please update that when you get a chance as well.

Once you get the branch updated based on the comments, feel free to hit the "Allow maintainers to edit" checkbox and I can take care of any lingering items at that point and get this merged. Thanks for your contribution again and sorry for the long delay on this PR!

[paulromano]

I've stared at this a long time and it looks to me like both terms here are positive, which would mean there's really no point of treating the off-diagonal terms separately. Based on the equation, it looks like the gains should be negative?

[church89]

You're right, I guess I wanted to separate the two terms (loss and gain) for better comprehensions, but it's much more cleaner and simple to combine in one simple notation.
Also regarding the sign it's a bit misleading, I'll update the equation and the text to better explain.

[paulromano]

Sounds good, thanks @church89! Just ping me once this PR is ready for another review

[church89]

Hey @paulromano, should be good now. Unfortunately I don't find the "Allow maintainers to edit" button as the PR is created from a fork living under an organization (I think it's a current Github problem, inconveniently not easy to overcome)