From 86702a5e03aaf91c51b669ee85ed698e3acf4713 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Wed, 17 Nov 2021 10:57:41 -0800 Subject: [PATCH] Draft: Particle fuel modeling Provides a yaml-interface for - modeling particle fuel with arbitrary layers - multiple particle fuel types in the reactor - assigning particle fuel as children to some parent component We have been decently successful with these changes internally in that downstream plugins can see `Component.particleFuel` and perform actions based on their content. What follows is an overview of the interface, implementation, and a discussion of where to go next. Related to https://github.com/terrapower/armi/issues/228 and would support modeling the MHTGR-350 benchmark https://github.com/terrapower/armi/issues/224 This patch is submitted to kick-start a discussion on better ways to add this feature into ARMI, leveraging the domain knowledge of the ARMI developers and the "particle fuel aware" plugins internally developed at USNC Tech. Input interface --------------- ```yaml particle fuel: demo: kernel: material: UO2 id: 0 od: 0.6 Tinput: 900 Thot: 900 flags: DEPLETABLE buffer: material: SiC id: 0.6 od: 0.61 Tinput: 900 Thot: 900 ``` ```yaml matrix: shape: Circle material: Graphite Tinput: 1200 Thot: 1200 id: 0.0 od: 2.2 latticeIDs: [F] flags: DEPLETABLE particleFuelSpec: demo particleFuelPackingFraction: 0.4 ``` With this interface it's possible to define several specifications in the model and assign them to different cylindrical components. Implementation -------------- The particle fuel is stored as a child of the parent component, such that `.children` is used to dynamically find the particle fuel spec. We can't store the specification as an attribute because we have to support potentially dynamic addition and removal of children to this component. Something like `self.particleFuel = spec` that makes spec a child would also have to understand what happens if we remove the spec from the parent, e.g., `self.remove(self.particleFuel)`. What is then the outcome of `self.particleFuel` unless we always check the children of the matrix? By making the particle fuel spec a `Composite` and placing it in the `.children` of the parent, the spec is able to be written to and read from the database. Adds `Component.setParticleMultiplicity`. The method is called during block construction when the block's height is available to the matrix component (particle's parent). The multiplicity is determined from the matrix volume and target packing fraction. Note: the particle mult is, by design, for a single component. Unresolved issues ----------------- - Volume of the parent matrix is not reduced by the volume occupied by the particles. - No support for homogenizing regions that contain particle fuel - Various homogenization properties don't account for particle fuel (e.g., `Core.getHM*`) - Particle fuel is not included in some text-reporting, leading to statements like ``` [info] Nuclide categorization for cross section temperature assignments: ------------------ ------------------------------------------------------ Nuclide Category Nuclides ------------------ ------------------------------------------------------ Fuel ``` and ``` [warn] The system has no heavy metal and therefore is not a nuclear reactor. Please make sure that this is intended and not a input error. ``` - No provided routines for packing the particles into their parent. This could be facilitated with a dedicated packing plugin and an unstructured 3-D `SpatialGrid` class. It's burdensome to expect the user to define the exact location of _every_ particle every time. But, if some `Plugin` performs the packing and creates this spatial grid, the multiplicity is tackled, and you potentially avoid adding or removing particles as their parent expands or contracts. - Unsure if material modifications make their way down to the materials in the particle fuel spec. The implementation suggests it as the `matMods` argument is passed into the particle fuel YAML object constructor. But we have yet to stress test that Next steps ---------- This patch is submitted because we continue to find places where this approach does not play well with the rest of the ARMI stack. While Blocks that contain particle fuel are correctly able to compute their heavy metal mass by iterating over their children, which in turn finds the particle fuel. However, higher-level actions like `Core.getHM*` do not go down to the sub-block level, instead asking to homogenize Blocks. The homogenization methods are not yet aware of the particle fuel because they rely on `Component.getNuclides` which reports the nuclides for it's `Material`, and does not include the children. This approach is sensible because if I'm writing a neutronic input file and I can exactly model the matrix and it's particle fuel, I would expect `matrix.getNuclides` to return the nuclides for _just the matrix_. Then, being informed of the particle fuel, I can write those materials and geometry uniquely. However, codes that cannot handle particle fuel and/or work with homogenized regions (e.g., nodal diffusion) would need this homogenized data. Allowing `Component.getNuclides` to return the nuclides on the child particle fuel would support this case, but not the previous case. My speculation is that the optimal strategy lies somewhere in making the matrix object not a `Component` but a `Block` and having the ARMI Composite tree accept blocks that potentially contain blocks. I think this is valid using the API but the user interface would need some work. Having the matrix be a block would provide a better interface for homogenization and exact representation of the particle fuel. I think... Other related changes --------------------- Added `Sphere.getBoundingCircleOuterDiameter` so that the particle fuel composites can be properly added to the database. They need to be "sortable" other wise `Database3._createLayout` breaks trying to sort components. This enables the particle fuel spec to be added to and read from the HDF data file The material constructor is a more public function as it is needed both in creating `Components` but also in creating the particle fuel spec. Added `MATRIX` flag Signed-off-by: Andrew Johnson --- armi/bookkeeping/report/reportInterface.py | 3 + armi/bookkeeping/report/reportingUtils.py | 85 +++++ armi/reactor/blueprints/__init__.py | 4 + armi/reactor/blueprints/blockBlueprint.py | 8 + armi/reactor/blueprints/componentBlueprint.py | 314 +++++++++++++++--- .../blueprints/tests/test_blockBlueprints.py | 92 +++++ .../blueprints/tests/test_blueprints.py | 94 +++++- .../tests/test_componentBlueprint.py | 218 ++++++++++++ armi/reactor/components/basicShapes.py | 4 + armi/reactor/components/component.py | 56 ++++ .../reactor/components/componentParameters.py | 11 + armi/reactor/components/volumetricShapes.py | 17 + armi/reactor/flags.py | 1 + armi/reactor/particleFuel.py | 90 +++++ armi/reactor/tests/test_components.py | 7 + 15 files changed, 963 insertions(+), 41 deletions(-) create mode 100644 armi/reactor/particleFuel.py diff --git a/armi/bookkeeping/report/reportInterface.py b/armi/bookkeeping/report/reportInterface.py index d9f81546f..1c27ba801 100644 --- a/armi/bookkeeping/report/reportInterface.py +++ b/armi/bookkeeping/report/reportInterface.py @@ -56,6 +56,7 @@ def interactBOL(self): runLog.important("Beginning of BOL Reports") reportingUtils.makeCoreAndAssemblyMaps(self.r, self.cs) reportingUtils.writeAssemblyMassSummary(self.r) + reportingUtils.makeParticleFuelDesignReport(self.r) if self.cs["summarizeAssemDesign"]: reportingUtils.summarizePinDesign(self.r.core) @@ -115,6 +116,8 @@ def generateDesignReport(self, generateFullCoreMap, showBlockAxMesh): ) reportingUtils.makeBlockDesignReport(self.r) + reportingUtils.makeParticleFuelDesignReport(self.r) + def interactEOL(self): """Adds the data to the report, and generates it""" b = self.o.r.core.getFirstBlock(Flags.FUEL) diff --git a/armi/bookkeeping/report/reportingUtils.py b/armi/bookkeeping/report/reportingUtils.py index 0289ef430..b8630e2b8 100644 --- a/armi/bookkeeping/report/reportingUtils.py +++ b/armi/bookkeeping/report/reportingUtils.py @@ -25,6 +25,7 @@ import time import tabulate from copy import copy +from math import nan import numpy @@ -56,6 +57,8 @@ Operator_MasterMachine = "Master Machine:" Operator_Date = "Date and Time:" Operator_CaseDescription = "Case Description:" +# Convert a value in centimeters to micrometers +CM_TO_MICRO_METER = 10000 def writeWelcomeHeaders(o, cs): @@ -1002,3 +1005,85 @@ def makeCoreAndAssemblyMaps(r, cs, generateFullCoreMap=False, showBlockAxMesh=Tr COMPONENT_INFO = "Component Information" + + +def makeParticleFuelDesignReport(r): + """Add reports for particle fuel specification and usage + + Parameters + ---------- + r : Reactor + Reactor containing particle fuel designs. If none are found, no + work is done. + + """ + if not r.blueprints.particleFuelDesigns: + return + + descriptions = {} + + for name, design in r.blueprints.particleFuelDesigns.items(): + rows = _particleFuelSpecToTable(design) + descriptions[name] = { + "spec": rows, + "usage": [], + } + + # Find all usages of each specification + for blockDesign in r.blueprints.blockDesigns.values(): + for componentDesign in blockDesign.values(): + thisParticleFuel = componentDesign.particleFuelSpec + if thisParticleFuel is None: + continue + descriptions[thisParticleFuel]["usage"].append( + ( + f"{blockDesign.name} {componentDesign.name}", + f"{componentDesign.particleFuelPackingFraction}", + ) + ) + + for name, subdata in descriptions.items(): + grp = report.data.Table(f"{name} Specification") + report.setData( + "Layer", + "Material, Outer diameter (μm), Thickness (μm)", + group=grp, + reports=report.DESIGN, + ) + for layerName, layerData in subdata["spec"]: + report.setData(layerName, layerData, grp, report.DESIGN) + + grp = report.data.Table(f"{name} Usage") + if subdata["usage"]: + report.setData("Block and Component name", "Packing fraction") + for desc, pf in subdata["usage"]: + report.setData(desc, pf, grp, report.DESIGN) + else: + msg = f"Particle fuel specification {name} not used" + report.setData("WARNING", msg, grp, report.DESIGN) + runLog.warning(msg) + + return + + +def _particleFuelSpecToTable(spec) -> list: + rows = [] + names = [] + prevOD = None + for layer in sorted(spec.values(), key=lambda ring: ring.od): + od = layer.od + if prevOD is not None: + thickness = round(CM_TO_MICRO_METER * 0.5 * (od - prevOD)) + else: + thickness = "-" + names.append(layer.name) + rows.append( + "{}, {}, {}".format( + layer.material, + round(od * CM_TO_MICRO_METER, 6), + thickness, + ) + ) + prevOD = od + + return list(zip(names, rows)) diff --git a/armi/reactor/blueprints/__init__.py b/armi/reactor/blueprints/__init__.py index 367526cef..11bc69514 100644 --- a/armi/reactor/blueprints/__init__.py +++ b/armi/reactor/blueprints/__init__.py @@ -104,6 +104,7 @@ from armi.reactor.blueprints.blockBlueprint import BlockKeyedList from armi.reactor.blueprints import isotopicOptions from armi.reactor.blueprints.gridBlueprint import Grids, Triplet +from armi.reactor.blueprints.componentBlueprint import ParticleFuelKeyedList context.BLUEPRINTS_IMPORTED = True context.BLUEPRINTS_IMPORT_CONTEXT = "".join(traceback.format_stack()) @@ -185,6 +186,9 @@ class Blueprints(yamlize.Object, metaclass=_BlueprintsPluginCollector): ) systemDesigns = yamlize.Attribute(key="systems", type=Systems, default=None) gridDesigns = yamlize.Attribute(key="grids", type=Grids, default=None) + particleFuelDesigns = yamlize.Attribute( + key="particle fuel", type=ParticleFuelKeyedList, default=None + ) # These are used to set up new attributes that come from plugins. Defining its # initial state here to make pylint happy diff --git a/armi/reactor/blueprints/blockBlueprint.py b/armi/reactor/blueprints/blockBlueprint.py index 2ac8081af..0593aeb4e 100644 --- a/armi/reactor/blueprints/blockBlueprint.py +++ b/armi/reactor/blueprints/blockBlueprint.py @@ -168,6 +168,14 @@ def construct( b.autoCreateSpatialGrids() except (ValueError, NotImplementedError) as e: runLog.warning(str(e), single=True) + # check if particle fuel exists and set particle mult + # Note: these changes occur here during block construction instead of during + # component contruction because component parent parameters (i.e., height) + # are needed for volume calculations + for component in b.getChildren(): + if component.particleFuel: + component.setParticleMultiplicity() + return b def _getGridDesign(self, blueprint): diff --git a/armi/reactor/blueprints/componentBlueprint.py b/armi/reactor/blueprints/componentBlueprint.py index 32cfffb63..56c4405a5 100644 --- a/armi/reactor/blueprints/componentBlueprint.py +++ b/armi/reactor/blueprints/componentBlueprint.py @@ -17,6 +17,8 @@ Special logic is required for handling component links. """ +from math import isclose, inf + import six import yamlize @@ -24,8 +26,22 @@ from armi import materials from armi.reactor import components from armi.reactor.flags import Flags +from armi.reactor.particleFuel import ParticleFuel from armi.utils import densityTools -from armi.nucDirectory import nuclideBases +from armi.nucDirectory import nuclideBases, nucDir + + +class ComponentParticleFuel(yamlize.Object): + specifier = yamlize.Attribute(type=str) + packingFraction = yamlize.Attribute(type=float) + + @packingFraction.validator + def packingFraction(self, packingFraction): + if not 0 < packingFraction < 1: + raise ValueError( + "Packing fraction must be between 0 and 1, exclusive. Got " + "{}".format(packingFraction) + ) class ComponentDimension(yamlize.Object): @@ -144,6 +160,8 @@ def shape(self, shape): # pylint: disable=no-self-use; reason=yamlize requireme orientation = yamlize.Attribute(type=str, default=None) mergeWith = yamlize.Attribute(type=str, default=None) area = yamlize.Attribute(type=float, default=None) + particleFuelSpec = yamlize.Attribute(type=str, default=None) + particleFuelPackingFraction = yamlize.Attribute(type=float, default=None) def construct(self, blueprint, matMods): """Construct a component""" @@ -167,6 +185,7 @@ def construct(self, blueprint, matMods): if component.hasFlags(Flags.DEPLETABLE): # depletable components, whether auto-derived or explicitly flagged need expanded nucs _insertDepletableNuclideKeys(component, blueprint) + return component def _conformKwargs(self, blueprint, matMods): @@ -189,6 +208,16 @@ def _conformKwargs(self, blueprint, matMods): # Don't pass these to the component constructor. These are used to # override the flags derived from the type, if present. continue + elif attr.name == "particleFuelSpec": + design = blueprint.particleFuelDesigns[val] + value = design.construct(blueprint, matMods) + elif attr.name == "particleFuelPackingFraction": + value = attr.get_value(self) + if value <= 0 or value >= 1: + raise ValueError( + f"Packing fraction {value} not allowed: must be between " + "0 and 1 exclusive", + ) else: value = attr.get_value(self) @@ -202,48 +231,77 @@ def _conformKwargs(self, blueprint, matMods): return kwargs def _constructMaterial(self, blueprint, matMods): - nucsInProblem = blueprint.allNuclidesInProblem - # make material with defaults - mat = materials.resolveMaterialClassByName(self.material)() - - if self.isotopics is not None: - # Apply custom isotopics before processing input mods so - # the input mods have the final word - blueprint.customIsotopics.apply(mat, self.isotopics) - - # add mass fraction custom isotopics info, since some material modifications need to see them - # e.g. in the base Material.applyInputParams - matMods.update( - { - "customIsotopics": { - k: v.massFracs for k, v in blueprint.customIsotopics.items() - } - } - ) - if len(matMods) > 1: - # don't apply if only customIsotopics is in there - try: - # update material with updated input params from blueprints file. - mat.applyInputParams(**matMods) - except TypeError: - # This component does not accept material modification inputs of the names passed in - # Keep going since the modification could work for another component - pass - - expandElementals(mat, blueprint) - - missing = set(mat.p.massFrac.keys()).difference(nucsInProblem) + return constructMaterial(self.material, blueprint, matMods, self.isotopics) + + +def constructMaterial(name, blueprint, matMods, isotopics): + """Create a material from blueprint, applying material modifications as necessary + + Parameters + ---------- + name : str + Name of this material. ARMI must know how to resolve the material + class given the string name + blueprint : Blueprints + Various detailed information, such as nuclides to model + matMods : dict + Material modifications to be applied to this material + isotopics : dict or none + Isotopics to apply to this material, if given. Must conform to + ``custom isotopics`` blueprints specification. + + Returns + ------- + armi.materials.Material + + Raises + ------ + ValueError + Nuclides found in the material aren't modelled in the ``nuclide flags`` + section of the blueprints problem - if missing: - raise ValueError( - "The nuclides {} are present in material {} by compositions, but are not " - "specified in the `nuclide flags` section of the input file. " - "They need to be added, or custom isotopics need to be applied.".format( - missing, mat - ) + """ + # make material with defaults + mat = materials.resolveMaterialClassByName(name)() + + if isotopics is not None: + # Apply custom isotopics before processing input mods so + # the input mods have the final word + blueprint.customIsotopics.apply(mat, isotopics) + + # add mass fraction custom isotopics info, since some material modifications need to see them + # e.g. in the base Material.applyInputParams + matMods.update( + { + "customIsotopics": { + k: v.massFracs for k, v in blueprint.customIsotopics.items() + } + } + ) + if len(matMods) > 1: + # don't apply if only customIsotopics is in there + try: + # update material with updated input params from blueprints file. + mat.applyInputParams(**matMods) + except TypeError: + # This component does not accept material modification inputs of the names passed in + # Keep going since the modification could work for another component + pass + + expandElementals(mat, blueprint) + + missing = set(mat.p.massFrac.keys()).difference(blueprint.allNuclidesInProblem) + + if missing: + raise ValueError( + "The nuclides {} are present in material {} by compositions, but are not " + "specified in the `nuclide flags` section of the input file. " + "They need to be added, or custom isotopics need to be applied.".format( + missing, mat ) + ) - return mat + return mat def expandElementals(mat, blueprint): @@ -310,3 +368,179 @@ def _insertDepletableNuclideKeys(c, blueprint): dimName, yamlize.Attribute(name=dimName, type=ComponentDimension, default=None), ) + + +class _ParticleFuelLayer(yamlize.Object): + """Component-like specification for a single layer in particle fuel""" + + name = yamlize.Attribute(key="name", type=str) + material = yamlize.Attribute(type=str) + innerDiam = yamlize.Attribute(key="id", type=float, default=0) + od = yamlize.Attribute(key="od", type=float) + Tinput = yamlize.Attribute(type=float) + Thot = yamlize.Attribute(type=float) + # Need this to pick up flags like depletable + flags = yamlize.Attribute(type=str, default=None) + + def construct(self, blueprint, matMods): + """Construct a sphere and assign to a parent""" + # Very similar to ComponentBlueprint.construct, maybe share? + runLog.debug(f"Constructing particle fuel layer {self.name}") + kwargs = self._conformKwargs(blueprint, matMods) + component = components.factory("sphere", [], kwargs) + + # the component __init__ calls setType(), which gives us our initial guess at + # what the flags should be. + if self.flags is not None: + # override the flags from __init__ with the ones from the blueprint + component.p.flags = Flags.fromString(self.flags) + else: + # potentially add the DEPLETABLE flag. Don't do this if we set flags + # explicitly. WARNING: If you add flags explicitly, it will + # turn off depletion so be sure to add depletable to your list of flags + # if you expect depletion + if any(nuc in blueprint.activeNuclides for nuc in component.getNuclides()): + component.p.flags |= Flags.DEPLETABLE + + if component.hasFlags(Flags.DEPLETABLE): + # depletable components, whether auto-derived or explicitly flagged + # need expanded nucs + _insertDepletableNuclideKeys(component, blueprint) + + if any(nucDir.isHeavyMetal(nucName) for nucName in component.getNuclides()): + component.p.flags |= Flags.FUEL + + return component + + def _conformKwargs(self, blueprint, matMods) -> dict: + """Return dictionary of arguments to help with component construction""" + kwargs = {} + for attr in self.attributes: + key = attr.name + if key == "innerDiam": + key = "id" + elif key == "flags": + continue + val = attr.get_value(self) + if key == "material": + value = constructMaterial(val, blueprint, matMods, None) + else: + value = attr.get_value(self) + kwargs[key] = value + return kwargs + + +class ParticleFuelSpec(yamlize.KeyedList): + """Specification for a single particle fuel type""" + + item_type = _ParticleFuelLayer + key_attr = _ParticleFuelLayer.name + name = yamlize.Attribute(type=str) + + def construct(self, blueprint, matMods): + """Produce a particle fuel instance that can be attached to a component""" + bounds = set() + + layers = sorted(self.values(), key=lambda spec: spec.od) + spec = ParticleFuel(self.name) + + prevInner = -inf + for layer in layers: + comp = layer.construct(blueprint, matMods) + innerDim = comp.getDimension("id") + if innerDim < 0: + raise ValueError( + f"{layer.name} inner diameter must be non-negative, got " + f"{innerDim}" + ) + if innerDim <= prevInner: + inners = [ring.innerDiam for ring in layers] + raise ValueError( + f"{self.name} has inconsistent inner diameters: not " + f"increasing {inners}" + ) + od = comp.getDimension("od") + if innerDim >= od: + raise ValueError( + f"{layer.name} outer diameter must be greater than inner diameter. " + f"Got {od} {innerDim}" + ) + + bounds.update((innerDim, od)) + spec.add(comp) + prevInner = innerDim + + nLayers = len(self) + if len(bounds) != nLayers + 1: + names = [m.name for m in spec] + raise ValueError( + f"{self.name} does not have consistent boundaries. Bounds: " + f"{sorted(bounds)}, compositions: {names}" + ) + + if not isclose(min(bounds), 0): + raise ValueError( + f"Particle fuel {self.name} does not start at radius of zero" + ) + + return spec + + +class ParticleFuelKeyedList(yamlize.KeyedList): + """Keyed list to enable the ``particleFuel`` section + + Structure + --------- + + ..code:: yaml + + particle fuel: + : + : # not strictly in increasing order + material: + id: + od: + Tinput: + Thot: + flags: + + Example + ------- + + ..code:: yaml + + particle fuel: + TRISO: + kernel: + material: UO2 + id: 0.6 + od: 0.61 + Tinput: 900 + Thot: 900 + flags: DEPLETABLE + buffer: + material: SiC + id: 0.6 + od: 0.62 + Tinput: 900 + Thot: 900 + ... + + This enables components to have a ``particle fuel`` section in the + blueprints file, indicating a material like a matrix is filled with particle + fuel + + ..code:: yaml + + : + : + name: matrix + material: Graphite + # more options + particleFuelSpec: TRISO + particleFuelPackingFraction: 0.4 # 40% packing + + """ + + item_type = ParticleFuelSpec + key_attr = ParticleFuelSpec.name diff --git a/armi/reactor/blueprints/tests/test_blockBlueprints.py b/armi/reactor/blueprints/tests/test_blockBlueprints.py index 767e5d6fb..84592ef40 100644 --- a/armi/reactor/blueprints/tests/test_blockBlueprints.py +++ b/armi/reactor/blueprints/tests/test_blockBlueprints.py @@ -14,6 +14,7 @@ """Tests for block blueprints.""" import unittest import io +import math from armi.reactor import blueprints from armi import settings @@ -318,6 +319,97 @@ def test_explicitFlags(self): self.assertTrue(a2.hasFlags(Flags.FUEL | Flags.TEST, exact=True)) +FULL_BP_PARTICLES = """ +blocks: + block: &block + flags: fuel test + duct: + shape: Hexagon + material: Graphite + Tinput: 600 + Thot: 600 + ip: 20 + op: 21 + matrix: + flags: matrix + shape: Circle + material: Graphite + Tinput: 600.0 + Thot: 600.0 + id: 0.0 + od: {matrix_od} + particleFuelSpec: dual + particleFuelPackingFraction: {pf} +assemblies: + assembly: &assembly_a + specifier: IC + blocks: [*block] + height: [{height}] + axial mesh points: [1] + xs types: [A] +particle fuel: + dual: + kernel: + material: UraniumOxide + Thot: 900 + Tinput: 900 + id: 0.0 + od: {kernel_od} + flags: DEPLETABLE + shell: + material: SiC + Tinput: 800 + Thot: 800 + id: {kernel_od} + od: {shell_od} +nuclide flags: + U235: &nuc_flags {{burn: false, xs: true}} + U238: *nuc_flags + C: *nuc_flags + SI: *nuc_flags + O: *nuc_flags +""" + + +class TestBlockWithParticles(unittest.TestCase): + + PF = 0.4 + HEIGHT = 10.0 + MATRIX_OD = 1.0 + KERNEL_OD = 0.1 + SHELL_OD = 0.2 + spec = FULL_BP_PARTICLES.format( + pf=PF, + height=HEIGHT, + matrix_od=MATRIX_OD, + kernel_od=KERNEL_OD, + shell_od=SHELL_OD, + ) + + def setUp(self): + self.cs = settings.Settings() + + bp = blueprints.Blueprints.load(self.spec) + assem = bp.constructAssem(self.cs, "assembly") + fuelBlock = assem.getFirstBlock(Flags.FUEL) + self.matrix = fuelBlock.getComponent(Flags.MATRIX) + + def test_particle_mult(self): + + self.assertIsNotNone(self.matrix.particleFuel) + + maxtrixRadius = self.MATRIX_OD / 2 + particleRadius = self.SHELL_OD / 2 + + nominalMatrixVolume = math.pi * maxtrixRadius ** 2 * self.HEIGHT + singleParticleVolume = 4 * math.pi * particleRadius ** 3 / 3 + + expectedMult = round(nominalMatrixVolume * self.PF / singleParticleVolume) + + for component in self.matrix.particleFuel.layers: + self.assertEqual(expectedMult, component.p.mult) + + if __name__ == "__main__": # import sys;sys.argv = ['', 'Test.testName'] unittest.main() diff --git a/armi/reactor/blueprints/tests/test_blueprints.py b/armi/reactor/blueprints/tests/test_blueprints.py index 54b395cb4..984d9a6c9 100644 --- a/armi/reactor/blueprints/tests/test_blueprints.py +++ b/armi/reactor/blueprints/tests/test_blueprints.py @@ -28,7 +28,11 @@ from armi.utils import directoryChangers from armi.utils import textProcessors from armi.reactor.blueprints.isotopicOptions import NuclideFlags, CustomIsotopics -from armi.reactor.blueprints.componentBlueprint import ComponentBlueprint +from armi.reactor.blueprints.componentBlueprint import ( + ComponentBlueprint, + ParticleFuelKeyedList, + ComponentParticleFuel, +) class TestBlueprints(unittest.TestCase): @@ -469,6 +473,94 @@ def test_withoutBlocks(self): self.assertEqual(1, len(a)) self.assertEqual(1, len(a[0])) + def test_loadComponentWithoutParticleFuel(self): + structure = "{material: UO2, Tinput: 700, Thot: 700, name: fuel, shape: circle}" + bp = ComponentBlueprint.load(structure) + self.assertIsNone(bp.particleFuelSpec) + + def test_loadComponentWithParticleFuel(self): + structure = ( + "{material: UO2, Tinput: 700, Thot: 700, name: fuel, " + "shape: circle, particleFuelSpec: triso, " + "particleFuelPackingFraction: 0.4}" + ) + bp = ComponentBlueprint.load(structure) + self.assertEqual(bp.particleFuelSpec, "triso") + self.assertEqual(bp.particleFuelPackingFraction, 0.4) + + def test_loadManyParticleFuelSpecs(self): + spec = """ +particle a: + kernel: + material: UO2 + id: 0.6 + od: 0.61 + Tinput: 900 + Thot: 900 + flags: DEPLETABLE + buffer: + material: SiC + id: 0.61 + od: 0.62 + Tinput: 800 + Thot: 800 +particle b: + ball: + material: Graphite + id: 0.7 + od: 0.8 + Tinput: 600.5 + Thot: 600.5 +""" + particles = ParticleFuelKeyedList.load(spec) + self.assertIn("particle a", particles) + self.assertIn("particle b", particles) + self.assertEqual(len(particles), 2) + + speca = particles["particle a"] + self.assertIn("kernel", speca) + self.assertIn("buffer", speca) + self.assertEqual(len(speca), 2) + + kernel = speca["kernel"] + kExpected = ( + ("material", "UO2"), + ("Thot", 900), + ("Tinput", 900), + ("od", 0.61), + # id stored as innerDiam on intermediate object + ("innerDiam", 0.6), + ("flags", "DEPLETABLE"), + ) + for key, value in kExpected: + self.assertEqual(getattr(kernel, key), value, msg=key) + + buffer = speca["buffer"] + bExpected = ( + ("material", "SiC"), + ("innerDiam", 0.61), + ("od", 0.62), + ("Thot", 800), + ("Tinput", 800), + ) + for key, value in bExpected: + self.assertEqual(getattr(buffer, key), value, msg=key) + + specb = particles["particle b"] + self.assertIn("ball", specb) + self.assertEqual(len(specb), 1) + + ball = specb["ball"] + expected = ( + ("material", "Graphite"), + ("Thot", 600.5), + ("Tinput", 600.5), + ("innerDiam", 0.7), + ("od", 0.8), + ) + for key, value in expected: + self.assertEqual(getattr(ball, key), value, msg=key) + if __name__ == "__main__": # import sys;sys.argv = ['', 'TestBlueprints.test_nuclides']] diff --git a/armi/reactor/blueprints/tests/test_componentBlueprint.py b/armi/reactor/blueprints/tests/test_componentBlueprint.py index 5ba7b8f56..b09282a57 100644 --- a/armi/reactor/blueprints/tests/test_componentBlueprint.py +++ b/armi/reactor/blueprints/tests/test_componentBlueprint.py @@ -20,8 +20,11 @@ import unittest from armi import settings +from armi.reactor.components import Component +from armi.reactor.particleFuel import ParticleFuel from armi.reactor import blueprints from armi.reactor.flags import Flags +from armi.nucDirectory import nucDir class TestComponentBlueprint(unittest.TestCase): @@ -339,6 +342,221 @@ def test_componentInitializationThoriumNoBurnCustomIsotopics(self): self.assertIn(nuc, a[0][0].getNuclides()) +class TestCompsWithParticleFuel(unittest.TestCase): + componentString = """ +blocks: + block: &block + duct: + shape: Hexagon + material: Graphite + Tinput: 600 + Thot: 600 + ip: 20 + op: 21 + component: + flags: MATRIX + shape: Circle + material: SiC + Tinput: 600.0 + Thot: 600.0 + id: 0.0 + od: 0.8660 + particleFuelSpec: dual + particleFuelPackingFraction: {pf} +assemblies: + assembly: &assembly_a + specifier: IC + blocks: [*block] + height: [1.0] + axial mesh points: [1] + xs types: [A] +particle fuel: + dual: + kernel: + material: UraniumOxide + Thot: 900 + Tinput: 900 + id: {innerID} + od: {innerOD} + flags: DEPLETABLE + shell: + material: SiC + Tinput: 800 + Thot: 800 + id: {outerID} + od: {outerOD} +nuclide flags: + U235: {{burn: false, xs: true}} + U238: {{burn: false, xs: true}} + C: {{burn: false, xs: true}} + SI: {{burn: false, xs: true}} + O: {{burn: false, xs: true}} +""" + DEF_PF = 0.4 + DEF_INNER_ID = 0.0 + DEF_INNER_OD = 0.6 + DEF_OUTER_ID = 0.6 + DEF_OUTER_OD = 0.7 + + def render(self, **kwargs): + kwargs.setdefault("pf", self.DEF_PF) + kwargs.setdefault("innerID", self.DEF_INNER_ID) + kwargs.setdefault("innerOD", self.DEF_INNER_OD) + kwargs.setdefault("outerID", self.DEF_OUTER_ID) + kwargs.setdefault("outerOD", self.DEF_OUTER_OD) + return self.componentString.format(**kwargs) + + def getCompWithParticleFuel(self, **kwargs) -> Component: + spec = self.render(**kwargs) + bp = blueprints.Blueprints.load(spec) + cs = settings.Settings() + assem = bp.constructAssem(cs, "assembly") + return assem[0][1] + + def test_valid(self): + comp = self.getCompWithParticleFuel() + packedSpec = comp.particleFuel + self.assertIsNotNone(packedSpec) + self.assertEqual(packedSpec.name, "dual") + self.assertEqual(comp.p.get("packingFractionBOL"), self.DEF_PF) + spheres = packedSpec.layers + self.assertEqual(len(spheres), 2) + + inner, outer = spheres + self.assertEqual(inner.getDimension("id"), self.DEF_INNER_ID) + self.assertEqual(inner.getDimension("od"), self.DEF_INNER_OD) + self.assertEqual(inner.material.name, "Uranium Oxide") + self.assertTrue(inner.hasFlags(Flags.DEPLETABLE)) + + self.assertEqual(outer.getDimension("id"), self.DEF_OUTER_ID) + self.assertEqual(outer.getDimension("od"), self.DEF_OUTER_OD) + self.assertEqual(outer.material.name, "Silicon Carbide") + + def test_invalidPF(self): + # yamlize intercepts the exception raised during validation + # and instead raises an yamlize error. We don't care for the + # specific type of exception, just that one is raised indicating + # invalid packing fraction + with self.assertRaisesRegex( + Exception, "Packing fraction.*must be between 0 and 1" + ): + self.getCompWithParticleFuel(pf=-1) + + def test_invalidLayerDims(self): + spec = self.render(outerOD=0.95 * self.DEF_OUTER_ID) + bp = blueprints.Blueprints.load(spec) + cs = settings.Settings() + with self.assertRaisesRegex( + ValueError, + "^shell outer diameter must be greater than inner", + ): + bp.constructAssem(cs, "assembly") + + spec = self.render(innerID=-1) + bp = blueprints.Blueprints.load(spec) + with self.assertRaisesRegex( + ValueError, "^kernel inner diameter must be non-negative" + ): + bp.constructAssem(cs, "assembly") + + def test_specWithGaps(self): + # Define a spec where the outer bound of one layer + # is not the inner bound of the next layer + # zero = inner_id < inner_od < outer_id < outer_od + # where we should have inner_od == outer_id + outerID = 1.1 * self.DEF_INNER_OD + spec = self.render(outerID=outerID, outerOD=1.1 * outerID) + + bp = blueprints.Blueprints.load(spec) + cs = settings.Settings() + + with self.assertRaisesRegex(ValueError, "consistent boundaries"): + bp.constructAssem(cs, "assembly") + + def test_specWithOverlaps(self): + # Define a spec where the two layers overlap + # zero = inner_id < outer_id < inner_od < outer_od + outerID = 0.5 * (self.DEF_INNER_ID + self.DEF_INNER_OD) + spec = self.render(outerID=outerID) + + bp = blueprints.Blueprints.load(spec) + cs = settings.Settings() + + with self.assertRaisesRegex(ValueError, "consistent boundaries"): + bp.constructAssem(cs, "assembly") + + def test_noZeroID(self): + # Define a spec that does not start at zero + innerID = 0.5 * (self.DEF_INNER_ID + self.DEF_INNER_OD) + spec = self.render(innerID=innerID) + + bp = blueprints.Blueprints.load(spec) + cs = settings.Settings() + + with self.assertRaisesRegex( + ValueError, ".*dual does not start at radius of zero" + ): + bp.constructAssem(cs, "assembly") + + def test_weirdOrdering(self): + # Define a specification where both rings start at r=0 but have + # different outer diameters + spec = self.render( + outerID=self.DEF_INNER_ID, + ) + + bp = blueprints.Blueprints.load(spec) + cs = settings.Settings() + + with self.assertRaisesRegex(ValueError, "inconsistent inner diameters"): + bp.constructAssem(cs, "assembly") + + def test_fuelFlag(self): + # Test if FUEL flag was added to particle kernel + fuelBlock = self.getCompWithParticleFuel().parent + + hasHeavyMetal = False + for child in fuelBlock: + spec = child.particleFuel + if spec is not None: + for component in spec.layers: + if any( + nucDir.isHeavyMetal(nucName) + for nucName in component.getNuclides() + ): + self.assertIn(Flags.FUEL, component.p.flags) + hasHeavyMetal = True + else: + self.assertNotIn(Flags.FUEL, component.p.flags) + self.assertTrue(hasHeavyMetal) + + def test_particleParent(self): + # check the parent/child relationship for a component and particle fuel + comp = self.getCompWithParticleFuel() + spec = comp.particleFuel + + self.assertIs(spec.parent, comp) + + for layer in spec.layers: + self.assertIs(layer.parent, spec) + self.assertIn(layer, spec) + + def test_multipleParticleFuelChildren(self): + comp = self.getCompWithParticleFuel() + + comp.add(ParticleFuel("extra")) + + with self.assertRaisesRegex(ValueError, "^Multiple"): + comp.particleFuel + + def test_noParticleFuel(self): + comp = self.getCompWithParticleFuel() + comp.remove(comp.particleFuel) + + with self.assertRaisesRegex(ValueError, "^No"): + comp.particleFuel + + if __name__ == "__main__": # import sys;sys.argv = ['', 'TestComponentBlueprint.test_componentInitializationAmericiumCustomIsotopics'] unittest.main() diff --git a/armi/reactor/components/basicShapes.py b/armi/reactor/components/basicShapes.py index d1d3d2073..0b88f1124 100644 --- a/armi/reactor/components/basicShapes.py +++ b/armi/reactor/components/basicShapes.py @@ -47,6 +47,8 @@ def __init__( isotopics=None, mergeWith=None, components=None, + particleFuelSpec=None, + particleFuelPackingFraction: float = 0, ): ShapedComponent.__init__( self, @@ -57,6 +59,8 @@ def __init__( isotopics=isotopics, mergeWith=mergeWith, components=components, + particleFuelSpec=particleFuelSpec, + particleFuelPackingFraction=particleFuelPackingFraction, ) self._linkAndStoreDimensions( components, od=od, id=id, mult=mult, modArea=modArea diff --git a/armi/reactor/components/component.py b/armi/reactor/components/component.py index 3bae8ee96..285f2e781 100644 --- a/armi/reactor/components/component.py +++ b/armi/reactor/components/component.py @@ -142,6 +142,8 @@ class ComponentType(composites.CompositeModelType): "name", "components", "area", + "particleFuelSpec", + "particleFuelPackingFraction", ) def __new__(cls, name, bases, attrs): @@ -181,6 +183,8 @@ class Component(composites.Composite, metaclass=ComponentType): Temperature in C to which dimensions were thermally-expanded upon input. material : str or material.Material The material object that makes up this component and give it its thermo-mechanical properties. + particleFuel : armi.reactor.component.particleFuel.ParticleFuel or None + Particle fuel specification """ DIMENSION_NAMES = tuple() # will be assigned by ComponentType @@ -216,6 +220,8 @@ def __init__( isotopics="", mergeWith="", components=None, + particleFuelSpec=None, + particleFuelPackingFraction: float = 0, ): if components and name in components: raise ValueError( @@ -233,6 +239,38 @@ def __init__( self.setType(name) self.p.mergeWith = mergeWith self.p.customIsotopicsName = isotopics + if particleFuelSpec is not None: + self.p.packingFractionBOL = particleFuelPackingFraction + self.add(particleFuelSpec) + + @property + def particleFuel(self): + if self.p.packingFractionBOL == 0: + return None + # Import here to avoid circular logic of having + # Component -> ParticleFuel -> VolumetricShapes -> Component + from armi.reactor.particleFuel import ParticleFuel + + # We can't store the particle fuel as an attribute because + # we have to support potentially dynamic addition and removal + # of children to this component. Something like + # self.particleFuel = spec + # that makes spec a child would also have to understand what + # happens if we remove the spec from the parent, e.g., + # self.remove(self.particleFuel) + # self.particleFuel + specs = self.getChildren( + deep=False, + includeMaterials=False, + generationNum=1, + predicate=lambda child: isinstance(child, ParticleFuel), + ) + if len(specs) == 1: + return specs[0] + raise ValueError( + f"{'Multiple' if specs else 'No'} particle fuel specifications " + f"found on {self} despite non-zero packing fraction" + ) @property def temperatureInC(self): @@ -1170,6 +1208,24 @@ def getPitchData(self): "Please implement if this component type can be a pitch defining component." ) + def setParticleMultiplicity(self): + """Estimate and set the multiplicity for the child ParticleFuel components.""" + + targetVolumeFraction = self.p.packingFractionBOL + + # Note: the matrix volume currently includes particle volume + matrixVolume = self.getVolume() + + # calculate particle volume + singleParticleVolume = 0.0 + for layer in self.particleFuel.layers: + singleParticleVolume += layer.getVolume() + + mult = round(targetVolumeFraction * matrixVolume / singleParticleVolume) + + for component in self.particleFuel.layers: + component.setDimension("mult", mult) + class ShapedComponent(Component): """A component with well-defined dimensions.""" diff --git a/armi/reactor/components/componentParameters.py b/armi/reactor/components/componentParameters.py index de7198fd1..011125a52 100644 --- a/armi/reactor/components/componentParameters.py +++ b/armi/reactor/components/componentParameters.py @@ -121,6 +121,17 @@ def getComponentParameterDefinitions(): description="Pin number of this component in some mesh. Starts at 1.", default=None, ) + + pb.defParam( + "packingFractionBOL", + default=0.0, + units=None, + description=( + "Packing fraction between zero and one for particle fuel on " + "this component. A value of zero indicates no packed particles" + ), + ) + return pDefs diff --git a/armi/reactor/components/volumetricShapes.py b/armi/reactor/components/volumetricShapes.py index 3006f8fa3..4fbdd03ca 100644 --- a/armi/reactor/components/volumetricShapes.py +++ b/armi/reactor/components/volumetricShapes.py @@ -67,9 +67,26 @@ def getComponentVolume(self): od = self.getDimension("od") iD = self.getDimension("id") mult = self.getDimension("mult") + mult = mult if mult is not None else 1 vol = mult * 4.0 / 3.0 * math.pi * ((od / 2.0) ** 3 - (iD / 2.0) ** 3) return vol + def getBoundingCircleOuterDiameter(self, Tc=None, cold=False): + """Return outer diameter of the sphere given thermal conditions + + Parameters + ---------- + Tc : float, optional + Temperature in C + cold : bool, optional + + Returns + ------- + float + Outer diameter + """ + return self.getDimension("od", Tc=Tc, cold=cold) + class Cube(ShapedComponent): """ diff --git a/armi/reactor/flags.py b/armi/reactor/flags.py index 73de9e27a..95ddedae1 100644 --- a/armi/reactor/flags.py +++ b/armi/reactor/flags.py @@ -276,6 +276,7 @@ class Flags(Flag): STRUCTURE = auto() DEPLETABLE = auto() + MATRIX = auto() @classmethod def fromStringIgnoreErrors(cls, typeSpec): diff --git a/armi/reactor/particleFuel.py b/armi/reactor/particleFuel.py new file mode 100644 index 000000000..6f4b052df --- /dev/null +++ b/armi/reactor/particleFuel.py @@ -0,0 +1,90 @@ +from typing import Tuple, Any, Union + +from armi.reactor.composites import Composite +from armi.reactor.components.volumetricShapes import Sphere + + +class ParticleFuel(Composite): + """Composite structure representing concentric spheres of particle fuel + + Each layer of the particle fuel is expected to be added via :meth:`add` + rather than ``__init__`` because we need this to be constructable via + blueprint file and the heirarchy written to / read from the database file. + + Parameters + ---------- + name : str + Name of this specification, typically something that is easily + recognizable by an engineering team, e.g., ``"AGR TRISO"`` + + Attributes + ---------- + layers : tuple of :class:`~armi.reactor.components.Sphere` + Each layer of the specification in order of increasing outer + diameter + + """ + + def __init__(self, name: str): + super().__init__(name) + self._layers = None + + def __repr__(self) -> str: + return f"<{type(self).__qualname__}({self.name}, {self.layers}>" + + def _checkChildAndClearCacheBeforeAdding(self, obj: Union[Any, Sphere]): + if isinstance(obj, Sphere): + self._layers = None + return + raise TypeError(f"Cannot add non-spherical layer {obj} to {self}") + + def add(self, obj: Union[Any, Sphere]): + """Add a component layer to this spec""" + self._checkChildAndClearCacheBeforeAdding(obj) + return super().add(obj) + + def remove(self, obj: Union[Any, Sphere]): + """Remove a layer from this spec""" + self._checkChildAndClearCacheBeforeAdding(obj) + return super().remove(obj) + + def insert(self, index: int, obj: Union[Any, Sphere]): + """Insert a layer in this spec + + .. note:: + + Ordering of layers is not based on ordering + of the children as ARMI defines them, but the + outer diameter of each individual layer + + """ + self._checkChildAndClearCacheBeforeAdding(obj) + return super().insert(index, obj) + + def append(self, obj: Union[Any, Sphere]): + """Append a layer to the children of this spec + + .. note:: + + Ordering of layers is not based on ordering + of the children as ARMI defines them, but the + outer diameter of each individual layer + + """ + self._checkChildAndClearCacheBeforeAdding(obj) + return super().append(obj) + + @property + def layers(self) -> Tuple[Sphere]: + if self._layers is None: + self._layers = tuple(sorted(self)) + return self._layers + + def __lt__(self, other: Union["ParticleFuel", Any]) -> bool: + """Is the outer diameter of the this spec less than that of another spec""" + if isinstance(other, type(self)): + mine = self.layers + theirs = other.layers + if mine is not None and theirs is not None: + return mine[-1].getDimension("od") < theirs[-1].getDimension("od") + return super().__lt__(other) diff --git a/armi/reactor/tests/test_components.py b/armi/reactor/tests/test_components.py index 62142d503..eb9306262 100644 --- a/armi/reactor/tests/test_components.py +++ b/armi/reactor/tests/test_components.py @@ -103,6 +103,8 @@ def test_componentInitializationAndDuplication(self): thisAttrs["name"] = "banana{}".format(i) if "modArea" in thisAttrs: thisAttrs["modArea"] = None + if "particleFuelSpec" in thisAttrs: + thisAttrs["particleFuelSpec"] = None component = components.factory(name, [], thisAttrs) duped = copy.deepcopy(component) for key, val in component.p.items(): @@ -841,6 +843,11 @@ def test_getVolume(self): def test_thermallyExpands(self): self.assertFalse(self.component.THERMAL_EXPANSION_DIMS) + def test_getBoundingCircleOuterDiameter(self): + od = self.component.getDimension("od") + actual = self.component.getBoundingCircleOuterDiameter() + self.assertEqual(actual, od) + class TestTorus(TestShapedComponent): componentCls = Torus