Repository: terrapower/armi

Issue:
Issue when custom isotopics are combined with material classes
When #1745 was implemented, some problems were quickly found (see #1818 ). #1822 was put through shortly after to fix those problems. However it appears that #1822 did not fix the whole problem, in particular when custom isotopics are combined with material classes that are children of `Fluid`. In these cases, the density of the `Fluid`-based components is not initialized correctly to the `Thot` temperature. Instead, the hot initialized fluid component stays at the density that was entered on the `custom isotopic`, when that density should actually correspond to the `Tinput` temperature. We need to make sure that the density is initialized properly to `Thot`.

This is only an issue with fluid components. The capability still appears to be working correctly with solid components.

I will put up a PR shortly, which will include a new unit test that demonstrates the problem.

Discussion:

Solution patch:
diff --git a/armi/reactor/blueprints/componentBlueprint.py b/armi/reactor/blueprints/componentBlueprint.py
index c871d8ca..54676cbe 100644
--- a/armi/reactor/blueprints/componentBlueprint.py
+++ b/armi/reactor/blueprints/componentBlueprint.py
@@ -290,20 +290,29 @@ class ComponentBlueprint(yamlize.Object):
                     "Cannot apply custom densities to materials without density."
                 )
 
-            # 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
+            # Apply a density scaling to account for the temperature change between
+            # Tinput and Thot
+            if isinstance(mat, materials.Fluid):
+                densityRatio = densityFromCustomIsotopic / mat.density(
+                    Tc=comp.inputTemperatureInC
+                )
             else:
-                f = 1.0 / (1 + dLL) ** 3
+                # for solids we need to consider if the input heights are hot or
+                # cold, in order to get the density correct.
+                # 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 = mat.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
 
-            scaledDensity = comp.density() / f
-            densityRatio = densityFromCustomIsotopic / scaledDensity
             comp.changeNDensByFactor(densityRatio)
 
             runLog.important(
diff --git a/doc/release/0.5.rst b/doc/release/0.5.rst
index 33051845..11a81876 100644
--- a/doc/release/0.5.rst
+++ b/doc/release/0.5.rst
@@ -21,6 +21,7 @@ Bug Fixes
 ---------
 #. Fixing check for jagged arrays during ``_writeParams``. (`PR#2051 <https://github.com/terrapower/armi/pull/2051>`_)
 #. Fixing BP-section ignoring tool in ``PassiveDBLoadPlugin``. (`PR#2055 <https://github.com/terrapower/armi/pull/2055>`_)
+#. Fixing number densities when custom isotopics are combined with Fluid components. (`PR#2071 <https://github.com/terrapower/armi/pull/2071>`_)
 #. Fixing scaling of volume-integrated parameters on edge assemblies. (`PR#2060 <https://github.com/terrapower/armi/pull/2060>`_)
 #. Fixing strictness of ``HexGrid`` rough equality check. (`PR#2058 <https://github.com/terrapower/armi/pull/2058>`_)
 #. TBD


Test patch:
diff --git a/armi/reactor/blueprints/tests/test_customIsotopics.py b/armi/reactor/blueprints/tests/test_customIsotopics.py
index 27ce95ed..c79a2343 100644
--- a/armi/reactor/blueprints/tests/test_customIsotopics.py
+++ b/armi/reactor/blueprints/tests/test_customIsotopics.py
@@ -19,6 +19,7 @@ from logging import DEBUG
 import yamlize
 
 from armi import runLog, settings
+from armi.materials import Fluid, Sodium
 from armi.physics.neutronics.settings import (
     CONF_MCNP_LIB_BASE,
     CONF_NEUTRONICS_KERNEL,
@@ -42,6 +43,7 @@ nuclide flags:
     AL: {burn: false, xs: true}
     FE: {burn: false, xs: true}
     C: {burn: false, xs: true}
+    NA: {burn: false, xs: true}
     DUMP2: {burn: true, xs: true}
     DUMP1: {burn: true, xs: true}
     LFP35: {burn: true, xs: true}
@@ -118,6 +120,11 @@ custom isotopics:
         C: 0.3
         density: 7.0
 
+    sodium custom isotopics:
+        input format: mass fractions
+        NA: 1
+        density: 666
+
 """
 
     yamlGoodBlocks = r"""
@@ -141,6 +148,25 @@ blocks:
             mult: 1.0
             od: 10.0
 
+        sodium1:
+            shape: Circle
+            material: Sodium
+            Tinput: 25
+            Thot: 600
+            id: 0
+            mult: 1
+            od: 1
+
+        sodium2:
+            shape: Circle
+            material: Sodium
+            isotopics: sodium custom isotopics
+            Tinput: 25
+            Thot: 600
+            id: 0
+            mult: 1
+            od: 1
+
     uranium fuel from isotopic mass fractions : &block_1
         fuel:
             <<: *basic_fuel
@@ -381,6 +407,35 @@ assemblies:
         self.assertEqual(fuel6.material.name, fuel0.material.name)
         self.assertEqual("UZr", fuel0.material.name)
 
+    def test_densitiesAppliedToNonCustomMaterialsFluid(self):
+        """
+        Ensure that a density can be set in custom isotopics for components using library materials,
+        specifically in the case of a fluid component. In this case, inputHeightsConsideredHot
+        does not matter, and the material has a zero dLL value.
+        """
+        # The template block
+        sodium1 = self.a[0].getComponentByName("sodium1")
+        sodium2 = self.a[0].getComponentByName("sodium2")
+
+        self.assertEqual(sodium1.material.name, "Sodium")
+        self.assertEqual(sodium2.material.name, "Sodium")
+        self.assertTrue(isinstance(sodium1.material, Fluid))
+        self.assertTrue(isinstance(sodium2.material, Fluid))
+        self.assertEqual(sodium1.p.customIsotopicsName, "")
+        self.assertEqual(sodium2.p.customIsotopicsName, "sodium custom isotopics")
+
+        # show that, even though the two components have the same material class
+        # and the same temperatures, their densities are different
+        self.assertNotEqual(sodium1.density(), sodium2.density())
+
+        # show that sodium1 has a density from the material class, while sodium2
+        # has a density from the blueprint and adjusted from Tinput -> Thot
+        s = Sodium()
+        self.assertAlmostEqual(sodium1.density(), s.density(Tc=600))
+        self.assertAlmostEqual(
+            sodium2.density(), s.density(Tc=600) * (666 / s.density(Tc=25))
+        )
+
     def test_customDensityLogsAndErrors(self):
         """Test that the right warning messages and errors are emitted when applying custom densities."""
         # Check for warnings when specifying both TD_frac and custom isotopics
diff --git a/armi/reactor/blueprints/tests/test_reactorBlueprints.py b/armi/reactor/blueprints/tests/test_reactorBlueprints.py
index 2e107517..0fa76239 100644
--- a/armi/reactor/blueprints/tests/test_reactorBlueprints.py
+++ b/armi/reactor/blueprints/tests/test_reactorBlueprints.py
@@ -143,7 +143,12 @@ class TestReactorBlueprints(unittest.TestCase):
 
     def test_materialDataSummary(self):
         """Test that the material data summary for the core is valid as a printout to the stdout."""
-        expectedMaterialData = [("Custom", "ARMI"), ("HT9", "ARMI"), ("UZr", "ARMI")]
+        expectedMaterialData = [
+            ("Custom", "ARMI"),
+            ("HT9", "ARMI"),
+            ("Sodium", "ARMI"),
+            ("UZr", "ARMI"),
+        ]
         core, _sfp, _evst = self._setupReactor()
         materialData = reactorBlueprint.summarizeMaterialData(core)
         for actual, expected in zip(materialData, expectedMaterialData):