Fixing the affect of custom isotopics on material density

armi/reactor/blueprints/blockBlueprint.py

@@ -28,6 +28,7 @@
 from armi.reactor.composites import Composite
 from armi.reactor.converters import blockConverters
 from armi.reactor.flags import Flags
+from armi.settings.fwSettings.globalSettings import CONF_INPUT_HEIGHTS_HOT
 
 
 def _configureGeomOptions():
@@ -144,7 +145,11 @@ def construct(
             filteredMaterialInput, byComponentMatModKeys = self._filterMaterialInput(
                 materialInput, componentDesign
             )
-            c = componentDesign.construct(blueprint, filteredMaterialInput)
+            c = componentDesign.construct(
+                blueprint,
+                filteredMaterialInput,
+                cs[CONF_INPUT_HEIGHTS_HOT],
+            )
             components[c.name] = c
 
             # check that the mat mods for this component are valid options

armi/reactor/blueprints/componentBlueprint.py

@@ -186,7 +186,7 @@ def shape(self, shape):
     mergeWith = yamlize.Attribute(type=str, default=None)
     area = yamlize.Attribute(type=float, default=None)
 
-    def construct(self, blueprint, matMods):
+    def construct(self, blueprint, matMods, inputHeightsConsideredHot):
         """Construct a component or group.
 
         .. impl:: User-defined on material alterations are applied here.
@@ -205,6 +205,17 @@ def construct(self, blueprint, matMods):
             The ``applyInputParams()`` method of that material class is then called,
             passing in the associated material modifications data, which the material
             class can then use to modify the isotopics as necessary.
+
+        Parameters
+        ----------
+        blueprint : Blueprints
+            Blueprints object containing various detailed information, such as nuclides to model
+
+        matMods : dict
+            Material modifications to apply to the component.
+
+        inputHeightsConsideredHot : bool
+            See the case setting of the same name.
         """
         runLog.debug("Constructing component {}".format(self.name))
         kwargs = self._conformKwargs(blueprint, matMods)
@@ -214,7 +225,9 @@ class can then use to modify the isotopics as necessary.
             constructedObject = composites.Composite(self.name)
             for groupedComponent in group:
                 componentDesign = blueprint.componentDesigns[groupedComponent.name]
-                component = componentDesign.construct(blueprint, matMods=dict())
+                component = componentDesign.construct(
+                    blueprint, {}, inputHeightsConsideredHot
+                )
                 # override free component multiplicity if it's set based on the group definition
                 component.setDimension("mult", groupedComponent.mult)
                 _setComponentFlags(component, self.flags, blueprint)
@@ -229,29 +242,33 @@ class can then use to modify the isotopics as necessary.
                 constructedObject.material.getTD()
             )
 
-        # set the custom density for non-custom material components after construction
-        self.setCustomDensity(constructedObject, blueprint, matMods)
+        self._setComponentCustomDensity(
+            constructedObject,
+            blueprint,
+            matMods,
+            inputHeightsConsideredHot,
+        )
 
         return constructedObject
 
-    def setCustomDensity(self, constructedComponent, blueprint, matMods):
+    def _setComponentCustomDensity(
+        self, comp, blueprint, matMods, inputHeightsConsideredHot
+    ):
         """Apply a custom density to a material with custom isotopics but not a 'custom material'."""
         if self.isotopics is None:
             # No custom isotopics specified
             return
 
-        density = blueprint.customIsotopics[self.isotopics].density
-        if density is None:
+        densityFromCustomIsotopic = blueprint.customIsotopics[self.isotopics].density
+        if densityFromCustomIsotopic is None:
             # Nothing to do
             return
 
-        if density <= 0:
+        if densityFromCustomIsotopic <= 0:
             runLog.error(
                 "A zero or negative density was specified in a custom isotopics input. "
                 "This is not permitted, if a 0 density material is needed, use 'Void'. "
-                "The component is {} and the isotopics entry is {}.".format(
-                    constructedComponent, self.isotopics
-                )
+                f"The component is {comp} and the isotopics entry is {self.isotopics}."
             )
             raise ValueError(
                 "A zero or negative density was specified in the custom isotopics for a component"
@@ -261,30 +278,40 @@ def setCustomDensity(self, constructedComponent, blueprint, matMods):
         if not isinstance(mat, materials.Custom):
             # check for some problem cases
             if "TD_frac" in matMods.keys():
-                runLog.warning(
-                    "Both TD_frac and a custom density (custom isotopics) has been specified for "
-                    "material {}. The custom density will override the density calculated using "
-                    "TD_frac.".format(self.material)
+                runLog.error(
+                    f"Both TD_frac and a custom isotopic with density {blueprint.customIsotopics[self.isotopics]} "
+                    f"has been specified for material {self.material}. This is an overspecification."
                 )
             if not mat.density(Tc=self.Tinput) > 0:
                 runLog.error(
-                    "A custom density has been assigned to material '{}', which has no baseline "
+                    f"A custom density has been assigned to material '{self.material}', which has no baseline "
                     "density. Only materials with a starting density may be assigned a density. "
-                    "This comes up e.g. if isotopics are assigned to 'Void'.".format(
-                        self.material
-                    )
+                    "This comes up e.g. if isotopics are assigned to 'Void'."
                 )
                 raise ValueError(
                     "Cannot apply custom densities to materials without density."
                 )
 
-            densityRatio = density / constructedComponent.density()
-            constructedComponent.changeNDensByFactor(densityRatio)
+            # Apply a density scaling to account for the temperature change between Tinput
+            # Thot. There may be a better place in the initialization to determine
+            # if the block height will be interpreted as hot dimensions, which would
+            # allow us to not have to pass the case settings down this far
+            dLL = comp.material.linearExpansionFactor(
+                Tc=comp.temperatureInC, T0=comp.inputTemperatureInC
+            )
+            if inputHeightsConsideredHot:
+                f = 1.0 / (1 + dLL) ** 2
+            else:
+                f = 1.0 / (1 + dLL) ** 3
+
+            scaledDensity = comp.density() / f
+            densityRatio = densityFromCustomIsotopic / scaledDensity
+            comp.changeNDensByFactor(densityRatio)
 
             runLog.important(
                 "A custom material density was specified in the custom isotopics for non-custom "
-                "material {}. The component density has been altered to "
-                "{}.".format(mat, constructedComponent.density()),
+                f"material {mat}. The component density has been altered to "
+                f"{comp.density() * scaledDensity} at temperature {comp.temperatureInC} C",
                 single=True,
             )
 

armi/reactor/blueprints/isotopicOptions.py

@@ -394,6 +394,7 @@ def apply(self, material):
         material.massFrac = dict(self.massFracs)
         if self.density is not None:
             if not isinstance(material, materials.Custom):
+                # TODO I think this warning should go away?
                 runLog.important(
                     "A custom density or number densities has been specified for non-custom "
                     "material {}. The material object's density will not be updated to prevent unintentional "

armi/reactor/blueprints/tests/test_blueprints.py

@@ -32,6 +32,7 @@
 from armi.utils import directoryChangers
 from armi.utils import textProcessors
 from armi.reactor.blueprints.gridBlueprint import saveToStream
+from armi.settings.fwSettings.globalSettings import CONF_INPUT_HEIGHTS_HOT
 
 
 class TestBlueprints(unittest.TestCase):
@@ -612,7 +613,7 @@ def test_topLevelComponentInput(self):
         # which is required during construction of a component
         design._resolveNuclides(cs)
         componentDesign = design.componentDesigns["freefuel"]
-        topComponent = componentDesign.construct(design, matMods=dict())
+        topComponent = componentDesign.construct(design, {}, cs[CONF_INPUT_HEIGHTS_HOT])
         self.assertEqual(topComponent.getDimension("od", cold=True), 4.0)
         self.assertGreater(topComponent.getVolume(), 0.0)
         self.assertGreater(topComponent.getMass("U235"), 0.0)

armi/reactor/blueprints/tests/test_customIsotopics.py

@@ -309,7 +309,12 @@ class TestCustomIsotopics(unittest.TestCase):
     @classmethod
     def setUpClass(cls):
         cs = settings.Settings()
-        cs = cs.modified(newSettings={CONF_XS_KERNEL: "MC2v2"})
+        cs = cs.modified(
+            newSettings={
+                CONF_XS_KERNEL: "MC2v2",
+                "inputHeightsConsideredHot": False,
+            }
+        )
 
         cls.bp = blueprints.Blueprints.load(cls.yamlString)
         cls.a = cls.bp.constructAssem(cs, name="fuel a")
@@ -356,14 +361,19 @@ def test_densitiesAppliedToNonCustomMaterials(self):
         # A block with custom density set via number density
         fuel8 = self.a[8].getComponent(Flags.FUEL)
 
+        dLL = fuel2.material.linearExpansionFactor(Tc=600, T0=25)
+        # the exponent here is 3 because inputHeightsConsideredHot = False.
+        # if inputHeightsConsideredHot were True, then we would use a factor of 2 instead
+        f = 1 / ((1 + dLL) ** 3)
+
         # Check that the density is set correctly on the custom density block,
         # and that it is not the same as the original
-        self.assertAlmostEqual(19.1, fuel2.density())
+        self.assertAlmostEqual(19.1 * f, fuel2.density())
         self.assertNotAlmostEqual(fuel0.density(), fuel2.density(), places=2)
         # Check that the custom density block has the correct material
         self.assertEqual("UZr", fuel2.material.name)
         # Check that the block with only number densities set has a new density
-        self.assertAlmostEqual(19.1, fuel8.density())
+        self.assertAlmostEqual(19.1 * f, fuel8.density())
         # original material density should not be changed after setting a custom density component,
         # so a new block without custom isotopics and density should have the same density as the original
         self.assertAlmostEqual(fuel6.density(), fuel0.density())