Use of custom isotopics with thermal expansion turned on does not correctly preserve mass

[keckler]

#1745 added the ability to combine a custom isotopic density vector with the material properties of a material class. I reviewed that PR and thought that everything looked good, but now I've been trying to use it on a real case and I've found a problem with mass conservation.

It appears that the component properly thermally expands. But during the thermal expansion, the component's number density is not adjusted away from the value specified in custom isotopics definition. This means that the density and the volume of the component become out of sync with each other, and as the component expands it gains mass.

I've attached a simple example input that you can use to reproduce the problem. It has only one single real component which uses a custom isotopic as well as a material class. Follow these steps to see the problem:

1. Note that Tinput == Thot in the component
2. Run the case
3. Open iPython
4. import armi; armi.configure()
5. from armi.bookkeeping.db import databaseFactory
6. db = databaseFactory("anl-afci-177.h5", 'r')
7. with db: r = db.load(0,0)
8. r.core[0][0].getMass()
9. Note the mass as 26.567467
10. Edit the blueprints with Thot: 250
11. Rerun the case
12. db = databaseFactory("anl-afci-177.h5", 'r')
13. with db: r = db.load(0,0)
14. r.core[0][0].getMass()
15. Note the mass as 26.7636707

These masses should be exactly equal, but they are not. I have tried the same thing without using a custom isotopic in the component definition to confirm that this is not a problem with the axial expansion routine itself, and the masses were the same in that case.

@bsculac @onufer

example.zip

[bsculac]

I replicated the problem, the output from the test case shows that the component with custom isotopics is initialized to the same mass, specifically the info message for the application of custom isotopics is the same:
[impt] A custom material density was specified in the custom isotopics for non-custom material <Material: MOX>. The component density has been altered to 26.567467529113692.

However, the BOL mass report indicates the difference:

250 C:

25 C:

So the error is occurring somewhere between component construction and the BOL mass report.

Outputs attached.

25.txt
250.txt

[jakehader]

I am not stepping through the debugger, but the mass density of ~26 g/cc seems odd and is set here -

armi/armi/reactor/blueprints/componentBlueprint.py

Line 278 in 41ce854

densityRatio = density / constructedComponent.density()

[Impt] A custom density or number densities has been specified for non-custom material <Material: MOX>. The material object's density will not be updated to prevent unintentional density changes across the model. Only custom materials may have a density specified.

I feel like this message is saying something cryptic about the behavior where the density once set on the component is not being scaled with temperature changes, so when volume changes due to thermal expansion the density is not being reduced to conserve mass. Is the material class assigned when Custom Isotopics are assigned set to Custom?

armi/armi/materials/custom.py

Line 27 in 41ce854

class Custom(Material):

If so, the Custom material object inherits from the base Material and the linear expansion is 0.0:

armi/armi/materials/material.py

Line 213 in 41ce854

def linearExpansionPercent(self, Tk: float = None, Tc: float = None) -> float:

On my phone so I didn't run the example, but just adding some bread crumbs to look at. I'd need to get a single component isolated in ipython, check the material object, and then see how the linear expansion is being applied. Separately, the initial density that is being set / derived seems wrong (?). Not sure what the custom density is being set to since I cannot open the zip file on my phone.

[keckler]

Don't focus on the value of ~26, it is just some contrived example that I built with placeholder values.

[jakehader]

Don't focus on the value of ~26, it is just some contrived example that I built with placeholder values.

I understand now, thanks! Update - the material is MOX and not custom. Still debugging by I think I know where the issue is.

[bsculac]

I feel like this message is saying something cryptic about the behavior where the density once set on the component is not being scaled with temperature changes, so when volume changes due to thermal expansion the density is not being reduced to conserve mass. Is the material class assigned when Custom Isotopics are assigned set to Custom?

The intent of that message is to indicate to the user that the base material density is not being changed, only the component that is using custom isotopics. So a second component using the same material but not the same custom isotopics will have a different density.

[jakehader]

@keckler / @bsculac - This is not necessarily a good fix, but can you test this out - #1822?

[bsculac]

@jakehader Even with #1822 the issue would persist, I think. The example blueprints have inputHeightsConsideredHot: False so something is being missed when the cold dimensions are expanded.

[bsculac]

https://github.com/terrapower/armi/blob/main/armi/reactor/converters/axialExpansionChanger.py#L112 This line in the axial expansion changer is most likely what is increasing the density of the components. I talked with @albeanth , my understanding is that the axial expansion changer assumes the heights provided in the blueprints are for hot heights, so the cold density of the component is increased such that the expanded/hot density will match the input density. So, #1822 could fix this but I'm still confused as to why the axial expansion changer will always increase the density. The function linked above even states that it should only be called when the inputHeightsConsideredHot setting is used, though it doesn't say whether that setting should be true or false. https://github.com/terrapower/armi/blob/main/armi/reactor/converters/axialExpansionChanger.py#L268-L283

[albeanth]

my understanding is that the axial expansion changer assumes the heights provided in the blueprints are for hot heights, so the cold density of the component is increased such that the expanded/hot density will match the input density.

The axial expansion changer doesn't assume anything about the blueprints height. That assumption of hot heights in the blueprints is in the building of the blueprints (i.e., armi/reactor/blueprints) and is a legacy assumption that hasn't been changed. So to correct for cold heights in the blueprints, the expandColdDimsToHot method gets called. This corrects the relationship between number density and volume (really the height).

I'm still confused as to why the axial expansion changer will always increase the density

It only increases the number density to account for the cold heights. I.e., when inputHeightsConsideredHot is false.

should only be called when the inputHeightsConsideredHot setting is used

This isn't the best description, I agree. Rather, it should say "when the setting inputHeightsConsideredHot is False" (since False is not the default).

The confusion stems from 1) a lack of documentation 😞 , and 2) hacking together solutions to fit the legacy way of having hot heights in the blueprints. It'd probably be clearer to just pull the plug and refactor the blueprints system to have inputted heights in the blueprints match Tinput. That's just not something we've prioritized.

[bsculac]

Thanks for the clarification. Since blueprints assume hot densities, I'm more on board with the fix in #1822 until we do something like what is discussed at #767 (comment).

I tested the fix in #1822 and the block masses come out the same for the 25C and 250C cases. The BOL mass report shows a discrepancy still, but that may be from having frozen sodium (the 25C case has more mass). Chris is going to be testing on a non-toy problem using #1822.

[jakehader]

Hi all - thanks for the discussion. I didn't intend #1822 to be the final solution but it peaked my curiosity. @bsculac - I know you originally made the PR that this ticket is about so you're welcome to take over the branch and we can think about a better way to resolve this without passing in a case setting down to the component construct method.

I got to a similar conclusion as you and @albeanth did but didn't get to write down all the nuances so I'm glad you spoke about it and wrote it down. I think that assuming hot dimensions for blocks is a modeling choice and we have historically had this because there was no axial expansion consideration for ARMI until it was implemented. Now that it is, there may be benchmark edge cases where this modeling choice is still convenient. It's hard to say we should remove it but doing so would definitely simplify the logic for expansion and make documenting the behavior or intended behavior much more straight forward.

Let me know if you'd like to drive that straw man PR / branch to closure or please provide me some feedback on how we can better account for this density scaling when Tinput != Thot. I could use some help crafting some better unit testing for this too. I think we should add testing around mass fractions and other areas of customIsotopics since it is useful but not stress tested all the time.

[albeanth]

@jakehader see #1824. Let's chat a little more about simplifying this beast 💪