Allowing users to set powerDensity instead of power (terrapower#1395)

armi/bookkeeping/db/tests/test_databaseInterface.py

@@ -129,6 +129,7 @@ def setUp(self):
         self.td = directoryChangers.TemporaryDirectoryChanger()
         self.td.__enter__()
         cs = settings.Settings(os.path.join(TEST_ROOT, "armiRun.yaml"))
+        cs = cs.modified(newSettings={"power": 0.0, "powerDensity": 9e4})
         self.o, cs = getSimpleDBOperator(cs)
         self.r = self.o.r
 
@@ -139,6 +140,9 @@ def test_metaData_endSuccessfully(self):
         def goodMethod(cycle, node):
             pass
 
+        # the power should start at zero
+        self.assertEqual(self.r.core.p.power, 0)
+
         self.o.interfaces.append(MockInterface(self.o.r, self.o.cs, goodMethod))
         with self.o:
             self.o.operate()
@@ -159,6 +163,9 @@ def goodMethod(cycle, node):
             self.assertIn("settings", h5["inputs"])
             self.assertIn("blueprints", h5["inputs"])
 
+        # after operating, the power will be greater than zero
+        self.assertGreater(self.r.core.p.power, 1e9)
+
     def test_metaDataEndFail(self):
         def failMethod(cycle, node):
             if cycle == 0 and node == 1:

armi/conftest.py

@@ -21,7 +21,7 @@
 
 Tests must be invoked via pytest for this to have any affect, for example::
 
-    $ pytest -n6 framework/armi
+    $ pytest -n 6 armi
 
 """
 import os

armi/operators/operator.py

@@ -355,10 +355,13 @@ def _cycleLoop(self, cycle, startingCycle):
         if halt:
             return False
 
+        # read total core power from settings (power or powerDensity)
+        basicPower = self.cs["power"] or (
+            self.cs["powerDensity"] * self.r.core.getHMMass()
+        )
+
         for timeNode in range(startingNode, int(self.burnSteps[cycle])):
-            self.r.core.p.power = (
-                self.powerFractions[cycle][timeNode] * self.cs["power"]
-            )
+            self.r.core.p.power = self.powerFractions[cycle][timeNode] * basicPower
             self.r.p.capacityFactor = (
                 self.r.p.availabilityFactor * self.powerFractions[cycle][timeNode]
             )
@@ -373,7 +376,7 @@ def _cycleLoop(self, cycle, startingCycle):
             else:
                 powFrac = self.powerFractions[cycle][timeNode - 1]
 
-            self.r.core.p.power = powFrac * self.cs["power"]
+            self.r.core.p.power = powFrac * basicPower
             self._timeNodeLoop(cycle, timeNode)
 
         self.interactAllEOC(self.r.p.cycle)

armi/operators/settingsValidation.py

@@ -20,25 +20,25 @@
 dialogues in the GUI. They say things like: "Your ___ setting has the value ___, which
 is impossible. Would you like to switch to ___?"
 """
-import re
+import itertools
 import os
+import re
 import shutil
-import itertools
 
 from armi import context
 from armi import getPluginManagerOrFail
 from armi import runLog
-from armi.utils import pathTools
-from armi.utils.mathematics import expandRepeatedFloats
+from armi.physics import neutronics
 from armi.reactor import geometry
 from armi.reactor import systemLayoutInput
-from armi.physics import neutronics
-from armi.utils import directoryChangers
 from armi.settings.settingsIO import (
     prompt,
     RunLogPromptCancel,
     RunLogPromptUnresolvable,
 )
+from armi.utils import directoryChangers
+from armi.utils import pathTools
+from armi.utils.mathematics import expandRepeatedFloats
 
 
 class Query:
@@ -468,8 +468,15 @@ def _willBeCopiedFrom(fName):
         )
 
         self.addQuery(
-            lambda: not self.cs["power"],
-            "No power level set. You must always start by importing a base settings file.",
+            lambda: not self.cs["power"] and not self.cs["powerDensity"],
+            "No power or powerDensity set. You must always start by importing a base settings file.",
+            "",
+            self.NO_ACTION,
+        )
+
+        self.addQuery(
+            lambda: self.cs["power"] > 0 and self.cs["powerDensity"] > 0,
+            "The power and powerDensity are both set, please note the power will be used as the truth.",
             "",
             self.NO_ACTION,
         )

armi/operators/snapshots.py

@@ -73,6 +73,7 @@ def _mainOperate(self):
             # need to update reactor power after the database load
             # this is normally handled in operator._cycleLoop
             self.r.core.p.power = self.cs["power"]
+            self.r.core.p.powerDensity = self.cs["powerDensity"]
 
             halt = self.interactAllBOC(self.r.p.cycle)
             if halt:

armi/physics/neutronics/globalFlux/globalFluxInterface.py

@@ -127,6 +127,7 @@ def _checkEnergyBalance(self):
             )
             / units.WATTS_PER_MW
         )
+        self.r.core.setPowerIfNecessary()
         specifiedPower = (
             self.r.core.p.power / units.WATTS_PER_MW / self.r.core.powerMultiplier
         )

armi/reactor/reactorParameters.py

@@ -366,6 +366,14 @@ def defineCoreParameters():
             "nuclear power generated in the core.",
         )
 
+        pb.defParam(
+            "powerDensity",
+            units=f"{units.WATTS}/{units.GRAMS}",
+            description="BOL Power density of the reactor core, in units of Watts per"
+            "grams of Heavy Metal Mass. After the BOL, the power parameter will be set, "
+            "and this will entirely overridden by that.",
+        )
+
         pb.defParam(
             "powerDecay",
             units=units.WATTS,

armi/reactor/reactors.py

@@ -395,8 +395,8 @@ def refAssem(self):
         The reference assembly is defined as the center-most assembly with a FUEL flag,
         if any are present, or the center-most of any assembly otherwise.
 
-        Warnings
-        ========
+        Warning
+        -------
         The convenience of this property should be weighed against it's somewhat
         arbitrary nature for any particular client. The center-most fueled assembly is
         not particularly representative of the state of the core as a whole.
@@ -443,6 +443,20 @@ def summarizeReactorStats(self):
             )
         )
 
+    def setPowerFromDensity(self):
+        """Set the power from the powerDensity."""
+        self.p.power = self.p.powerDensity * self.getHMMass()
+
+    def setPowerIfNecessary(self):
+        """Set the core power, from the power density.
+
+        If the power density is set, but the power isn't, we set the calculate the
+        total heavy metal mass of the reactor, and set the total power. Which will
+        then be the real source of truth again.
+        """
+        if self.p.power == 0 and self.p.powerDensity > 0:
+            self.setPowerFromDensity()
+
     def setBlockMassParams(self):
         """Set the parameters kgHM and kgFis for each block and calculate Pu fraction."""
         for b in self.getBlocks():
@@ -733,7 +747,9 @@ def getNumRings(self, indexBased=False):
         """
         Returns the number of rings in this reactor. Based on location so indexing will start at 1.
 
-        WARNING: If you loop through range(maxRing) then ring+1 is the one you want!!
+        Warning
+        -------
+        If you loop through range(maxRing) then ring+1 is the one you want!
 
         Parameters
         ----------
@@ -742,7 +758,6 @@ def getNumRings(self, indexBased=False):
 
         When circular ring shuffling is activated, this changes interpretation.
         Developers plan on making this another method for the secondary interpretation.
-
         """
         if self.circularRingList and not indexBased:
             return max(self.circularRingList)
@@ -1299,8 +1314,8 @@ def getFirstAssembly(self, typeSpec=None, exact=False):
         """
         Gets the first assembly in the reactor.
 
-        Warnings
-        --------
+        Warning
+        -------
         This function should be used with great care. There are **very** few
         circumstances in which one wants the "first" of a given sort of assembly,
         `whichever that may happen to be`. Precisely which assembly is returned is

armi/reactor/tests/test_reactors.py

@@ -1098,7 +1098,7 @@ def test_buildManualZonesEmpty(self):
         self.assertEqual(len(list(self.r.core.zones)), 0)
 
     def test_getNuclideCategories(self):
-        # test that nuclides are categorized correctly
+        """Test that nuclides are categorized correctly."""
         self.r.core.getNuclideCategories()
         self.assertIn("coolant", self.r.core._nuclideCategories)
         self.assertIn("structure", self.r.core._nuclideCategories)
@@ -1107,6 +1107,27 @@ def test_getNuclideCategories(self):
         self.assertIn("FE56", self.r.core._nuclideCategories["structure"])
         self.assertIn("U235", self.r.core._nuclideCategories["fuel"])
 
+    def test_setPowerIfNecessary(self):
+        self.assertAlmostEqual(self.r.core.p.power, 0)
+        self.assertAlmostEqual(self.r.core.p.powerDensity, 0)
+
+        # to start, this method shouldn't do anything
+        self.r.core.setPowerIfNecessary()
+        self.assertAlmostEqual(self.r.core.p.power, 0)
+
+        # take the powerDensity when needed
+        self.r.core.p.power = 0
+        self.r.core.p.powerDensity = 1e9
+        mass = self.r.core.getHMMass()
+        self.r.core.setPowerIfNecessary()
+        self.assertAlmostEqual(self.r.core.p.power, 1e9 * mass)
+
+        # don't take the powerDensity when not needed
+        self.r.core.p.power = 3e9
+        self.r.core.p.powerDensity = 2e9
+        self.r.core.setPowerIfNecessary()
+        self.assertAlmostEqual(self.r.core.p.power, 3e9)
+
 
 class CartesianReactorTests(ReactorTests):
     def setUp(self):

armi/settings/fwSettings/globalSettings.py

@@ -91,6 +91,7 @@
 CONF_PHYSICS_FILES = "savePhysicsFiles"
 CONF_PLOTS = "plots"
 CONF_POWER = "power"
+CONF_POWER_DENSITY = "powerDensity"
 CONF_POWER_FRACTIONS = "powerFractions"
 CONF_PROFILE = "profile"
 CONF_REALLY_SMALL_RUN = "reallySmallRun"
@@ -627,6 +628,14 @@ def defineSettings() -> List[setting.Setting]:
             "setting the powerFractions setting.",
             schema=vol.All(vol.Coerce(float), vol.Range(min=0)),
         ),
+        setting.Setting(
+            CONF_POWER_DENSITY,
+            default=0.0,
+            label="Reactor Thermal Power Density (W/HMM)",
+            description="Thermal power of the Reactor, per gram of Heavy metal "
+            "mass. Ignore this setting if the `power` setting is non-zero.",
+            schema=vol.All(vol.Coerce(float), vol.Range(min=0)),
+        ),
         setting.Setting(
             CONF_REMOVE_PER_CYCLE, default=3, label="Move per Cycle", description="None"
         ),

doc/release/0.2.rst

@@ -12,6 +12,7 @@ What's new in ARMI
 #. Broad cleanup of ``Parameters``: filled in all empty units and descriptions, removed unused params. (`PR#1345 <https://github.com/terrapower/armi/pull/1345>`_)
 #. Removed redundant ``Material.name`` variable. (`PR#1335 <https://github.com/terrapower/armi/pull/1335>`_)
 #. Moved the ``Reactor`` assembly number from the global scope to a ``Parameter``. (`PR#1383 <https://github.com/terrapower/armi/pull/1383>`_)
+#. Added ``powerDensity`` as a high-level alternative to ``power`` to configure a Reactor. (`PR#1395 <https://github.com/terrapower/armi/pull/1395>`_)
 #. Added SHA1 hashes of XS control files to the welcome text. (`PR#1334 <https://github.com/terrapower/armi/pull/1334>`_)
 #. Add python 3.11 to ARMI's CI testing GH actions! (`PR#1341 <https://github.com/terrapower/armi/pull/1341>`_)
 #. Put back ``avgFuelTemp`` block parameter. (`PR#1362 <https://github.com/terrapower/armi/pull/1362>`_)