diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index 4feb38547..b02fc763a 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -2,9 +2,9 @@ name: CI testing on: pull_request: - branches: ["*dev*", alpha-test*] + branches: ["*dev*", "*-test*"] push: - branches: ["*dev*", main, alpha-test*] + branches: ["*dev*", "main", "*-test*"] jobs: build: diff --git a/doc/source/_templates/mytallyclass.rst b/doc/source/_templates/mytallyclass.rst new file mode 100644 index 000000000..836466f1f --- /dev/null +++ b/doc/source/_templates/mytallyclass.rst @@ -0,0 +1,10 @@ +{{ objname }} +{{ underline }} + +.. currentmodule:: {{ module }} + +.. autoclass:: {{ objname }} + :members: + :inherited-members: + :undoc-members: + :show-inheritance: diff --git a/doc/source/api/modules.rst b/doc/source/api/modules.rst index b662a6121..2c8e8ffdc 100644 --- a/doc/source/api/modules.rst +++ b/doc/source/api/modules.rst @@ -36,6 +36,7 @@ Collections montepy.CommentCollection montepy.Materials montepy.Surfaces + montepy.Tallies montepy.Transforms montepy.Universes @@ -185,7 +186,78 @@ Materials montepy.Nucleus montepy.Nuclide montepy.ThermalScatteringLaw - + + + +Tallies +^^^^^^^ + +Tally Objects +~~~~~~~~~~~~~ + +.. autosummary:: + :toctree: generated + :nosignatures: + :template: mytallyclass.rst + + montepy.Tally + montepy.SurfaceTally + montepy.CellTally + montepy.DetectorTally + montepy.SurfaceCurrentTally + montepy.SurfaceFluxTally + montepy.CellFluxTally + montepy.EnergyDepositionTally + montepy.FissionEnergyDepositionTally + montepy.EnergyDetectorPulseTally + montepy.TallyType + montepy.Score + +Scoring Groups and Filters +~~~~~~~~~~~~~~~~~~~~~~~~~~ + +.. note:: + + You will rarely create these directly, rather get them from + :attr:`montepy.Tally.groups` and :attr:`montepy.Tally.filters`. + +.. autosummary:: + :toctree: generated + :nosignatures: + :template: mytallyclass.rst + + montepy.data_inputs.tally.TallyGroup + montepy.data_inputs.tally.FlatGroup + montepy.data_inputs.tally.PathGroup + montepy.data_inputs.tally.LatticeIndex + montepy.data_inputs.tally.Filter + montepy.data_inputs.tally.ParticleFilter + montepy.data_inputs.tally.SpatialFilter + +Tally Multipliers +~~~~~~~~~~~~~~~~~ + +.. note:: + + You will rarely create the ``Multiplier*`` classes directly, rather get + them from :attr:`montepy.TallyMultiplier.bins`. + +.. autosummary:: + :toctree: generated + :nosignatures: + :template: mytallyclass.rst + + montepy.TallyMultiplier + montepy.Reaction + montepy.data_inputs.tally_multiplier.ReactionExpression + montepy.data_inputs.tally_multiplier_type.ReactionOperator + montepy.MultiplierSet + montepy.SpecialMultiplierSet + montepy.SpecialMultiplier + montepy.AttenuatorLayer + montepy.AttenuatorSet + montepy.data_inputs.tally_multiplier.MultiplierScore + montepy.data_inputs.tally_multiplier.MultiplierBin Cell Modifiers diff --git a/doc/source/changelog.rst b/doc/source/changelog.rst index 18c7200d4..84cd0d31b 100644 --- a/doc/source/changelog.rst +++ b/doc/source/changelog.rst @@ -10,6 +10,14 @@ MontePy Changelog ============= +#Next Version# +-------------- + +**Features Added** + +* Added a :class:`~montepy.Tally`/:class:`~montepy.TallyMultiplier` object model for reading, creating, and editing MCNP tally (``F``) and tally multiplier (``FM``) inputs (:issue:`11`). + + 1.6.0b1 -------------- diff --git a/doc/source/developing.rst b/doc/source/developing.rst index 703659fee..350776e72 100644 --- a/doc/source/developing.rst +++ b/doc/source/developing.rst @@ -8,8 +8,8 @@ Developer's Reference ===================== MontePy can be thought of as having two layers: the syntax, and the semantic layers. -The syntax layers handle the boring syntax things: like multi-line cards, and comments, etc. -The semantic layer takes this information and makes sense of it, like what the material number in a cell card is. +The syntax layers handle the boring syntax things: like multi-line inputs, and comments, etc. +The semantic layer takes this information and makes sense of it, like what the material number in a cell input is. .. note:: diff --git a/doc/source/devguide/cell_modifier.rst b/doc/source/devguide/cell_modifier.rst index 8092da983..072888eb4 100644 --- a/doc/source/devguide/cell_modifier.rst +++ b/doc/source/devguide/cell_modifier.rst @@ -2,7 +2,7 @@ :description lang=en: How to implement CellModifierInput subclasses in MontePy for data inputs that modify cells. -Data Cards that Modify Cells: :class:`~montepy.data_inputs.cell_modifier.CellModifierInput` +Data Inputs that Modify Cells: :class:`~montepy.data_inputs.cell_modifier.CellModifierInput` ============================================================================================ This is a subclass of :class:`~montepy.data_inputs.data_input.DataInputAbstract` that is meant to handle data inputs that specify information about, and modify cells. diff --git a/doc/source/devguide/collections_dev.rst b/doc/source/devguide/collections_dev.rst index 682ead3aa..660225a32 100644 --- a/doc/source/devguide/collections_dev.rst +++ b/doc/source/devguide/collections_dev.rst @@ -87,7 +87,7 @@ During init the inputs' "name word" (e.g., ``M3``, ``kcode``, ``f7:n``) is valid Conceptually these names can contain up to four sections. This information is stored in an instance of :class:`~montepy.input_parser.syntax_node.ClassifierNode`. -#. A ``prefix_modifier`` this modifies the whole card with a special character such as ``*tr5`` +#. A ``prefix_modifier`` this modifies the whole input with a special character such as ``*tr5`` #. A ``Prefix``, which is a series of letters that identifies the type such as ``m`` #. A ``number``, which numbers it. These must be an unsigned integer. #. A particle classifier such as ``:n,p``. @@ -117,7 +117,7 @@ How to Add an Object to :class:`~montepy.MCNP_Problem` The :class:`~montepy.MCNP_Problem` automatically consumes problem level data inputs, and adds them to itself. -Cards this would be appropriate for would be things like ``mode`` and ``kcode``. +Inputs this would be appropriate for would be things like ``mode`` and ``kcode``. To do this it uses the dictionary ``inputs_to_property`` in the ``__load_data_inputs_to_object`` method. To add a problem level data Object you need to diff --git a/doc/source/devguide/mcnp_object.rst b/doc/source/devguide/mcnp_object.rst index 036596741..793f5dc82 100644 --- a/doc/source/devguide/mcnp_object.rst +++ b/doc/source/devguide/mcnp_object.rst @@ -7,7 +7,7 @@ Input: :class:`~montepy.mcnp_object.MCNP_Object` ================================================= -All classes that represent a single input card *must* subclass this. +All classes that represent a single input *must* subclass this. For example: some children are: :class:`~montepy.Cell`, :class:`~montepy.Surface`. How to ``__init__`` diff --git a/doc/source/devguide/pointers_generators.rst b/doc/source/devguide/pointers_generators.rst index 8a1fde3c2..8ff4acb3b 100644 --- a/doc/source/devguide/pointers_generators.rst +++ b/doc/source/devguide/pointers_generators.rst @@ -42,11 +42,11 @@ if the surfaces did know, this would be bidirectional. So how do we decide which direction to point? In general we should default to MCNP. -So a cell borrows a surface because a cell card in MCNP +So a cell borrows a surface because a cell input in MCNP references surface numbers, and not vice versa. The exception to this is the case of inputs that modify another object. -For example the ``MT`` card modifies its parent ``M`` card. +For example the ``MT`` input modifies its parent ``M`` input. In general the parent object should own its children modifiers. This is an area of new development, and this may change. diff --git a/doc/source/guide/cells.rst b/doc/source/guide/cells.rst index 288ad4ffa..cebbeaae3 100644 --- a/doc/source/guide/cells.rst +++ b/doc/source/guide/cells.rst @@ -66,10 +66,10 @@ Their importances will all be set to 0. Setting How Cell Data Gets Displayed in the Input File ------------------------------------------------------ -Much of the cell data can show up in the cell block or the data block, like the importance card. +Much of the cell data can show up in the cell block or the data block, like the importance input. These are referred to MontePy as "cell modifiers". You can change how these cell modifiers are printed with :attr:`~montepy.MCNP_Problem.print_in_data_block`. -This acts like a dictionary where the key is the MCNP card name. +This acts like a dictionary where the key is the MCNP input name. So to make cell importance data show up in the cell block just run: ``problem.print_in_data_block["imp"] = False``. diff --git a/doc/source/guide/reading_writing.rst b/doc/source/guide/reading_writing.rst index 718a1115a..e65272a3d 100644 --- a/doc/source/guide/reading_writing.rst +++ b/doc/source/guide/reading_writing.rst @@ -181,7 +181,7 @@ Information Kept #. The optional message block at the beginning of the problem (it's a niche feature; check out section :manual63:`4.4.1` of the user manual) #. The problem title #. ``C`` style comments (e.g., ``C this is a banana``) -#. (Almost) all MCNP inputs (cards). Only the read input is discarded. +#. (Almost) all MCNP inputs. Only the read input is discarded. #. Dollar sign comments (e.g., ``1 0 $ this is a banana``) #. Other user formatting and spaces. If extra spaces between values are given the space will be expanded or shortened to try to keep the position of the next value in the same spot as the length of the first value changes. @@ -195,11 +195,11 @@ Information Kept Information Lost ^^^^^^^^^^^^^^^^ -#. Read cards. These are handled properly, but when written out these cards themselves will disappear. - When MontePy encounters a read card it notes the file in the card, and then discard the card. +#. Read inputs. These are handled properly, but when written out these inputs themselves will disappear. + When MontePy encounters a read input it notes the file in the input, and then discards the input. It will then read these extra files and append their contents to the appropriate block. - So If you were to write out a problem that used the read card in the surface block the surface - cards in that file from the read card will appear at the end of the new surface block in the newly written file. + So If you were to write out a problem that used the read input in the surface block the surface + inputs in that file from the read input will appear at the end of the new surface block in the newly written file. .. note:: diff --git a/doc/source/guide/tallies.rst b/doc/source/guide/tallies.rst new file mode 100644 index 000000000..2069d8fca --- /dev/null +++ b/doc/source/guide/tallies.rst @@ -0,0 +1,611 @@ +.. meta:: + :description lang=en: + Working with MCNP tallies in MontePy: the Tally object model, building tallies, cloning, tally multipliers, and building reaction expressions with Python operators. + +Tallies +======= + +.. testsetup:: * + + import montepy + problem = montepy.read_input("tests/inputs/test.imcnp") + +This guide covers how to inspect and build tallies: what a tally +scores, which cells or surfaces it covers, and how its tally +multiplier modifies it. Every tally (``F``) input and tally multiplier +(``FM``) input in a problem is available as a real object through +``problem.tallies``, addressable by number like any other MontePy +collection. + +The Tally Object Hierarchy +--------------------------- + +Every tally in a problem is stored in ``problem.tallies``, a +:class:`~montepy.Tallies` collection, and can be accessed by its number like +any other collection in MontePy. + +.. testcode:: + + tally = problem.tallies[4] + print(tally) + +.. testoutput:: + + CellFluxTally: 4 + +MontePy picks the class for you based on the tally's type digit (e.g., the ``4`` in +``F4``). +:class:`~montepy.Tally` is the base class, and it has one subclass for every tally +type: + +.. list-table:: + :header-rows: 1 + + * - Quantity Tallied + - Type Digit + - MontePy Class + - Shorthand Alias + * - Surface current + - ``F1`` + - :class:`~montepy.SurfaceCurrentTally` + - ``montepy.F1Tally`` + * - Average surface flux + - ``F2`` + - :class:`~montepy.SurfaceFluxTally` + - ``montepy.F2Tally`` + * - Cell flux + - ``F4`` + - :class:`~montepy.CellFluxTally` + - ``montepy.F4Tally`` + * - Point/ring detector flux + - ``F5`` + - :class:`~montepy.DetectorTally` + - ``montepy.F5Tally`` + * - Energy deposition + - ``F6`` + - :class:`~montepy.EnergyDepositionTally` + - ``montepy.F6Tally`` + * - Fission energy deposition + - ``F7`` + - :class:`~montepy.FissionEnergyDepositionTally` + - ``montepy.F7Tally`` + * - Pulse height (energy deposition in a detector) + - ``F8`` + - :class:`~montepy.EnergyDetectorPulseTally` + - ``montepy.F8Tally`` + +The Shorthand Alias is just another name for the same class, e.g. +``montepy.F4Tally`` is :class:`~montepy.CellFluxTally`. + +Underneath these, there are two intermediate classes worth knowing about: +:class:`~montepy.SurfaceTally` for tallies that score on surfaces +(i.e., ``F1`` and ``F2``), and :class:`~montepy.CellTally` for tallies that score in +cells (i.e., ``F4``, ``F6``, ``F7``, and ``F8``). +You can check :attr:`~montepy.Tally.tally_type` to get the +:class:`~montepy.TallyType` for any tally. + +.. doctest:: + + >>> tally.tally_type + + +Geometry Filters +------------------ + +The cells or surfaces a tally scores over are its geometry filter, exposed as +:attr:`~montepy.Tally.groups`, a list of +:class:`~montepy.data_inputs.tally.TallyGroup` objects. + +This section covers flat cell/surface lists, where every entry is a +:class:`~montepy.data_inputs.tally.FlatGroup`. +A ``FlatGroup`` knows the numbers it covers, and whether they form a single averaged +bin or separate bins. +Repeated structures and lattices use the ``<``/``[...]``/``U=`` syntax instead, which +produces :class:`~montepy.data_inputs.tally.PathGroup` entries; see +`Universe and Lattice Paths`_ below. + +In the simple, "flat" case, a tally just lists cells or surfaces one after another, +and MCNP creates a separate bin for each one: + +.. testcode:: + + tally = problem.tallies[4] + for group in tally.groups: + cells = [problem.cells[n] for n in group.old_numbers] + print(*cells, group.is_grouped) + +.. testoutput:: + + Cell: 1 False + Cell: 2 False + Cell: 3 False + +``is_grouped`` is ``False`` for every group here, since ``F4:n 1 2 3`` has no +parentheses: each of cells 1, 2, and 3 gets its own separate bin. +``old_numbers`` is looked up by hand above to show how it relates to the raw numbers +on the card, but once a group is linked to a problem you don't need to do that +lookup yourself: :attr:`~montepy.data_inputs.tally.FlatGroup.cells_or_surfaces` gives +you the resolved objects directly. + +.. doctest:: + + >>> [c.number for c in tally.groups[0].cells_or_surfaces] + [1] + +Wrapping cells or surfaces in parentheses instead unions them into a single bin, +averaged for normalized tally types like ``F2``/``F4``/``F6``/``F7``, or summed for +``F1``/``F8``, rather than reported separately. +A tally can mix flat entries and multiple parenthesized groups on the same input, and +each group becomes its own entry in ``groups``: + +.. testcode:: + + grouped = montepy.CellFluxTally("f14:n (1 2) (3)") + for group in grouped.groups: + cells = [problem.cells[n] for n in group.old_numbers] + print(*cells, group.is_grouped) + +.. testoutput:: + + Cell: 1 Cell: 2 True + Cell: 3 True + +Notice that the second group here still has ``is_grouped`` set to ``True``, even +though it only covers a single cell: what matters is whether the parentheses were +there, not how many cells are inside them. Without the parentheses, ``f14:n 1 2 3`` +would instead produce three separate, ungrouped bins, exactly like the cell flux +tally example above. +To build flat and grouped bins like these from scratch instead of reading them from +an existing tally, see `Building Tallies from Scratch`_ below. + +For cell tallies, :attr:`~montepy.CellTally.cells` gives you the +flattened, deduplicated set of every cell the tally touches, across all of its groups. +Surface tallies have the equivalent +:attr:`~montepy.SurfaceTally.surfaces`. + +.. doctest:: + + >>> print(tally.cells) + Cells: [1, 2, 3] + +Use ``cells``/``surfaces`` when you just want to know what's involved. +Use ``groups`` when the bin structure itself matters. + +Checking whether a specific cell is scored by a tally works the way you'd expect: + +.. doctest:: + + >>> problem.cells[1] in tally + True + >>> problem.cells[99] in tally + False + +A tally can also end with a trailing ``T``, for "total": an extra bin that's the +union of every other bin on the input, rather than a scoring region of its own. +Because it isn't really its own region, MontePy doesn't represent it as another +entry in ``groups``. +Instead it's a separate flag, :attr:`~montepy.Tally.include_total`: + +.. testcode:: + + totaled = montepy.CellFluxTally("f24:n (1 2) (3) T") + for group in totaled.groups: + cells = [problem.cells[n] for n in group.old_numbers] + print(*cells, group.is_grouped) + print(totaled.include_total) + +.. testoutput:: + + Cell: 1 Cell: 2 True + Cell: 3 True + True + +Notice that ``groups`` only has the two real bins; the ``T`` never shows up as a +third entry there, no matter how many bins came before it. +Unlike ``scores`` and ``filters``, ``include_total`` is settable, so you can turn +total-bin reporting on (or off) for any tally, including one you build from scratch +with :func:`~montepy.CellTally.add_cell`/:func:`~montepy.CellTally.add_group` (see +`Building Tallies from Scratch`_ below): + +.. doctest:: + + >>> totaled.include_total = False + >>> totaled.include_total + False + >>> "T" in totaled.mcnp_str() + False + +Building Tallies from Scratch +------------------------------- + +It's also possible to build a tally from scratch. +Instantiate the concrete subclass for the quantity you want to score (e.g. +:class:`~montepy.CellFluxTally` for flux), give it a number, and add cells or +surfaces to it. + +.. testcode:: + + new_tally = montepy.CellFluxTally(number=104) + new_tally.add_cell(problem.cells[1]) + new_tally.add_cell(problem.cells[2]) + problem.tallies.append(new_tally) + +:func:`~montepy.CellTally.add_cell` adds a cell as its own separate +bin. +If you want a group of cells averaged into a single bin instead, use +:func:`~montepy.CellTally.add_group`: + +.. testcode:: + + new_tally.add_group([problem.cells[1], problem.cells[3]]) + +.. doctest:: + + >>> for group in new_tally.groups: + ... print(group.old_numbers, group.is_grouped) + [1] False + [2] False + [1, 3] True + +:class:`~montepy.SurfaceTally` has the matching +:func:`~montepy.SurfaceTally.add_surface` and +:func:`~montepy.SurfaceTally.add_group`. + +Each of these has a matching removal method: +:func:`~montepy.CellTally.remove_cell`/:func:`~montepy.SurfaceTally.remove_surface` +removes the single-item bin that ``add_cell``/``add_surface`` would have created, and +:func:`~montepy.CellTally.remove_group` removes any group outright, whether it came +from ``add_cell``, ``add_group``, or ``add_path_group``. +A cell or surface only drops out of :attr:`~montepy.CellTally.cells`/ +:attr:`~montepy.SurfaceTally.surfaces` once no remaining group references it. + +.. doctest:: + + >>> new_tally.remove_cell(problem.cells[2]) + >>> for group in new_tally.groups: + ... print(group.old_numbers, group.is_grouped) + [1] False + [1, 3] True + >>> problem.cells[2] in new_tally.cells + False + +Scores and Filters +-------------------- + +.. note:: + + "Score" and "filter" aren't MCNP terms. + MontePy borrows this vocabulary from + `OpenMC `_, since it's more general than + anything MCNP itself uses, and describes the same underlying concepts. + See `OpenMC's tallies user guide + `_ for more on this + terminology. + +Sometimes you don't need the raw group structure, you just want to know what physical +quantity a tally is measuring, and for which particles. +:attr:`~montepy.Tally.scores` and +:attr:`~montepy.Tally.filters` give you a shallow, read-only summary +of this. + +.. doctest:: + + >>> tally.scores + [] + +Every tally type has a default score: +:class:`~montepy.Score` is an enum with entries like ``FLUX``, +``CURRENT``, and ``ENERGY_DEPOSITION``. +This default is overridden if the tally has a tally multiplier attached; see +`Tally Multipliers`_ below. +``filters`` returns a list of :class:`~montepy.data_inputs.tally.ParticleFilter` and +:class:`~montepy.data_inputs.tally.SpatialFilter` objects, whichever apply: + +.. testcode:: + + f2 = problem.tallies[2] + for filt in f2.filters: + print(filt) + +.. testoutput:: + + ParticleFilter(':p') + SpatialFilter([FlatGroup([1005], grouped=False)]) + +.. note:: + + These aren't a new way to write tallies to the input file, they're just a + convenient way to read what's already there. + ``scores`` and ``filters`` aren't settable. + +Cloning Tallies +----------------- + +Like most MontePy objects, tallies support +:func:`~montepy.Tally.clone`, which copies a tally and gives it a +new, unused number of the same tally type. + +.. testcode:: + + clone = tally.clone() + +.. doctest:: + + >>> clone.number + 14 + +Tallies also have something the other objects don't: :func:`~montepy.Tally.clone_as`. +This copies a tally's scoring geometry into a *different* tally type. +It's handy when you already have flux tallied over a set of cells and you also want +the heating in those same cells, without retyping the cell list: + +.. testcode:: + + heating = tally.clone_as(montepy.EnergyDepositionTally) + +.. doctest:: + + >>> print(heating) + EnergyDepositionTally: 16 + >>> print(heating.cells) + Cells: [1, 2, 3] + >>> heating.scores + [] + +``clone_as`` only allows conversions within the same family: surface tallies +(i.e., ``F1`` and ``F2``) convert to other surface tallies, and cell tallies +(i.e., ``F4``, ``F6``, ``F7``, and ``F8``) convert to other cell tallies. +``F5`` point/ring detectors are their own family, since they don't have cells or +surfaces to carry over. +Trying to cross families raises a ``ValueError``. + +Tally Multipliers +------------------- + +A tally multiplier input multiplies a tally's flux or current by a cross section, turning a plain +flux tally into a reaction rate, a heating rate, or similar. +MontePy represents this as a :class:`~montepy.TallyMultiplier`, +linked to its tally through :attr:`~montepy.Tally.multiplier`. + +.. testcode:: + + fm = montepy.TallyMultiplier("fm4 (1.0 26 16 103)") + problem.tallies.append(fm) + +.. doctest:: + + >>> tally.multiplier is fm + True + +An ``FMn`` input is linked to its tally purely by number, the same way an ``MTn`` +thermal scattering input gets linked to material ``n``. +You can append the ``TallyMultiplier`` and its ``Tally`` to the problem in either +order, and MontePy will connect them once both are present. +Because that link is the whole point of an ``FMn`` card, MontePy keeps it +consistent for you: deleting a tally that has a linked multiplier removes the +now-orphaned ``FM`` card from the problem too (with a warning), and renumbering a +tally renumbers its linked multiplier to match. + +Since a ``TallyMultiplier``'s number always tracks its parent tally's rather than +being independently assignable, it has its own +:func:`~montepy.TallyMultiplier.clone`, separate from +:func:`~montepy.Tally.clone`: pass the tally to attach the clone to (it takes that +tally's number), or omit it for a detached, unregistered copy. + +.. testcode:: + + fm_clone = fm.clone(problem.tallies[14]) + +.. doctest:: + + >>> fm_clone.number + 14 + >>> problem.tallies[14].multiplier is fm_clone + True + +The bulk of a tally multiplier input is its :attr:`~montepy.TallyMultiplier.bins`, +a list of :class:`~montepy.data_inputs.tally_multiplier.MultiplierBin`. +Each bin holds one or more +:class:`~montepy.MultiplierSet` or +:class:`~montepy.SpecialMultiplierSet` terms, and +optionally an :class:`~montepy.AttenuatorSet`. + +.. testcode:: + + term = fm.bins[0].terms[0] + +.. doctest:: + + >>> term.constant + 1.0 + >>> term.material + 26 + +Once a tally has a multiplier, its :attr:`~montepy.Tally.scores` +switches from the generic default to a list of +:class:`~montepy.data_inputs.tally_multiplier.MultiplierScore`, one for every output +bin the multiplier defines. This gives you the full recipe (constant, material, +reaction, and any attenuator) for every column of the tally's output. + +.. doctest:: + + >>> tally.scores + [MultiplierScore(constant=1.0, material=26, reaction=ReactionExpression(Reaction(16), ReactionOperator.MULTIPLY, Reaction(103)), kind=None, attenuator=None)] + +Building Reaction Expressions +------------------------------ + +A reaction list on a tally multiplier input, like ``16 103`` in ``FM4 (1.0 26 16 +103)``, is really a small expression: multiply reaction 16 by reaction 103. +Rather than making you build this out of strings, MontePy lets you write it as an +actual Python expression, using ``*`` for multiply, ``+`` for add, and ``-`` for +subtract. + +.. note:: + + Raw MT numbers like ``16`` are hard to read and easy to mistype. + :class:`~montepy.Reaction` has named constants for the common ones, e.g. + ``Reaction.N_2N`` for MT 16; use those instead of ``Reaction(16)`` whenever a name + is available. + Raw numbers are still there for the reactions that don't have one. + +.. doctest:: + + >>> from montepy import Reaction + >>> expr = Reaction.N_2N * Reaction.N_P + >>> expr == term.reactions[0] + True + >>> fm.mcnp_str() + 'fm4 (1.0 26 16 103)' + +That's exactly the reaction expression already attached to ``tally`` above, built by +hand, matching the actual MCNP text of ``fm``. +Python's own operator precedence already does the right thing for a longer list too, +exactly like the MCNP manual says a reaction list should work, so you don't need to +write unnecessary parentheses to get the correct grouping. + +.. doctest:: + + >>> bigger_expr = Reaction.N_2N * Reaction.N_P + Reaction.N_D + >>> bigger_expr.left == expr + True + >>> bigger_expr.operator + + +You can go one step further and build a whole +:class:`~montepy.MultiplierSet` with ``&``, joining a +material number to a reaction expression: + +.. doctest:: + + >>> built = 26 & Reaction.N_2N * Reaction.N_P + >>> built == term + True + >>> fm.mcnp_str() + 'fm4 (1.0 26 16 103)' + +You can scale the constant afterward with ``*``: + +.. doctest:: + + >>> scaled = 1.5 * built + >>> scaled.constant + 1.5 + +:class:`~montepy.AttenuatorSet` supports the same kind of +chaining with ``&``, for building up multiple attenuating layers: + +.. doctest:: + + >>> from montepy import AttenuatorSet, AttenuatorLayer + >>> attenuator = AttenuatorSet(1.0, [AttenuatorLayer(26, 0.5)]) + >>> attenuator = attenuator & AttenuatorLayer(27, 0.3, is_atom_density=False) + >>> len(attenuator.layers) + 2 + >>> fm2 = montepy.TallyMultiplier("fm104:n (1.0 -1 26 0.5 27 -0.3)") + >>> attenuator == fm2.bins[0].attenuator + True + >>> fm2.mcnp_str() + 'fm104:n (1.0 -1 26 0.5 27 -0.3)' + +These operators build and compare :class:`~montepy.MultiplierSet`/ +:class:`~montepy.data_inputs.tally_multiplier.ReactionExpression` values, but a +:class:`~montepy.TallyMultiplier` needs one more step to actually add them to a card: +wrap each term (and optional attenuator) in a +:class:`~montepy.data_inputs.tally_multiplier.MultiplierBin`, then use +:func:`~montepy.TallyMultiplier.add_bin` (and :func:`~montepy.TallyMultiplier.remove_bin` +to take one back out). ``bins`` itself stays a read-only view, the same way ``groups`` +does for a :class:`~montepy.Tally`. + +.. testcode:: + + from montepy.data_inputs.tally_multiplier import MultiplierBin + + fm3 = montepy.TallyMultiplier(number=4) + fm3.add_bin(MultiplierBin([built])) + +.. doctest:: + + >>> fm3.bins[0].terms[0] == built + True + >>> fm3.mcnp_str() + 'FM4 1.0 26 16 103 ' + +Universe and Lattice Paths +----------------------------- + +This section is for the less common case: tallying a specific cell inside a specific +universe or lattice, using the ``<`` path syntax. +Most tallies don't need this. + +A tally group written with ``<`` becomes a +:class:`~montepy.data_inputs.tally.PathGroup` instead of a ``FlatGroup``. +Each step in the chain is still a ``FlatGroup``, accessible through +:attr:`~montepy.data_inputs.tally.PathGroup.levels`, ordered from innermost to +outermost. + +In the simplest case, a chain is just cell numbers separated by ``<``, read +innermost-to-outermost, e.g. "cell 2, inside cell 5": + +.. testcode:: + + simple_path = montepy.CellFluxTally("f104:n (2 < 5)") + for level in simple_path.groups[0].levels: + print(level.old_numbers) + +.. testoutput:: + + [2] + [5] + +A level can narrow further with a lattice-element index in brackets +(``[i j k]``), and a universe designator (``u=n``) can stand in for an entire level, +meaning "any cell in universe n": + +.. testcode:: + + path_tally = montepy.CellFluxTally("f114:n (u=1 < 2[0 0 0] < 5)") + for level in path_tally.groups[0].levels: + print(level.old_numbers, level.universe_spec) + +.. testoutput:: + + [] 1 + [2] None + [5] None + +The first level has no cell number at all, just a universe designator +(``u=1``), meaning "any cell in universe 1". +The second level narrows that down to lattice element ``[0 0 0]`` of cell 2, and the +third level says that whole path has to live inside cell 5. +That lattice element is available directly through +:attr:`~montepy.data_inputs.tally.FlatGroup.lattice_indices`, parallel to +``old_numbers``, as a list of :class:`~montepy.data_inputs.tally.LatticeIndex`: + +.. doctest:: + + >>> path_tally.groups[0].levels[1].lattice_indices[0].dimensions + [0, 0, 0] + +You can also build a path group from scratch with +:func:`~montepy.CellTally.add_path_group` and +:func:`~montepy.data_inputs.tally.PathGroup.inside`: + +.. testcode:: + + scratch_tally = montepy.CellFluxTally(number=304) + pg = scratch_tally.add_path_group(problem.cells[1]) + pg.inside(problem.cells[2]) + +.. doctest:: + + >>> for level in pg.levels: + ... print(level.old_numbers) + [1] + [2] + +References +---------- + +* :manual63:`5.9.1` +* :manual63:`5.9.7` diff --git a/doc/source/guide/universes.rst b/doc/source/guide/universes.rst index 51a586fe9..d9eff46de 100644 --- a/doc/source/guide/universes.rst +++ b/doc/source/guide/universes.rst @@ -97,7 +97,7 @@ You can also easy apply a transform to the filling universe with: MCNP supports some rather complicated cell filling systems. Mainly the ability to fill a cell with different universes for every lattice site, - and to create an "anonymous transform" in the fill card. + and to create an "anonymous transform" in the fill input. MontePy can understand and manipulate fills with these features in the input. However, generating these from scratch may be cumbersome. diff --git a/doc/source/starting.rst b/doc/source/starting.rst index 9b6d43f31..90b6530cc 100644 --- a/doc/source/starting.rst +++ b/doc/source/starting.rst @@ -13,7 +13,7 @@ Getting Started with MontePy MontePy is a Python API for reading, editing, and writing MCNP input files. The library provides a semantic interface for working with input files ("MCNP problems"). It does not run MCNP, nor does it parse MCNP output files. -It understands that the second entry on a cell card is the material number, +It understands that the second entry on a cell input is the material number, and will link the cell with its material object. .. note:: @@ -140,3 +140,4 @@ You can install this with ``pip install montepy[demo-present]``. guide/cells guide/materials guide/universes + guide/tallies diff --git a/montepy/__init__.py b/montepy/__init__.py index 33fa1e188..c957bec4a 100644 --- a/montepy/__init__.py +++ b/montepy/__init__.py @@ -26,6 +26,35 @@ from montepy.data_inputs import Mode from montepy.data_inputs.thermal_scattering import ThermalScatteringLaw from montepy.data_inputs.data_parser import parse_data +from montepy.data_inputs.tally import ( + Tally, + SurfaceTally, + CellTally, + DetectorTally, + SurfaceCurrentTally, + SurfaceFluxTally, + CellFluxTally, + EnergyDepositionTally, + FissionEnergyDepositionTally, + EnergyDetectorPulseTally, + F1Tally, + F2Tally, + F4Tally, + F5Tally, + F6Tally, + F7Tally, + F8Tally, +) +from montepy.data_inputs.tally_multiplier import ( + TallyMultiplier, + Reaction, + MultiplierSet, + AttenuatorLayer, + AttenuatorSet, + SpecialMultiplierSet, +) +from montepy.data_inputs.tally_type import TallyType, Score +from montepy.data_inputs.tally_multiplier_type import SpecialMultiplier # geometry from montepy.geometry_operators import Operator @@ -50,6 +79,7 @@ from montepy.universes import Universes from montepy.surface_collection import Surfaces from montepy.transforms import Transforms +from montepy.tallies import Tallies import montepy.exceptions import montepy.errors # deprecated diff --git a/montepy/cell.py b/montepy/cell.py index ed9b666ea..8a0647fce 100644 --- a/montepy/cell.py +++ b/montepy/cell.py @@ -711,6 +711,21 @@ def surfaces(self): ) return self._surfaces + @property + def tallies(self): + """Generator of tallies in the parent problem that score this cell. + + Yields + ------ + Tally + + .. versionadded:: 1.6.0b2 + """ + if self._problem: + for t in self._problem.tallies: + if self in t: + yield t + @property @needs_full_ast def parameters(self) -> dict[str, str]: diff --git a/montepy/data_inputs/__init__.py b/montepy/data_inputs/__init__.py index 60536b330..92fc386ce 100644 --- a/montepy/data_inputs/__init__.py +++ b/montepy/data_inputs/__init__.py @@ -5,3 +5,7 @@ from .material import Material from .thermal_scattering import ThermalScatteringLaw from .data_parser import parse_data +from .material import Material +from .tally import Tally +from .tally_multiplier import TallyMultiplier +from .thermal_scattering import ThermalScatteringLaw diff --git a/montepy/data_inputs/data_parser.py b/montepy/data_inputs/data_parser.py index 96629cf9c..e7ee8cf3f 100644 --- a/montepy/data_inputs/data_parser.py +++ b/montepy/data_inputs/data_parser.py @@ -11,6 +11,8 @@ lattice_input, material, mode, + tally, + tally_multiplier, thermal_scattering, universe_input, volume, @@ -23,6 +25,8 @@ lattice_input.LatticeInput, material.Material, mode.Mode, + tally.Tally, + tally_multiplier.TallyMultiplier, thermal_scattering.ThermalScatteringLaw, transform.Transform, volume.Volume, @@ -67,6 +71,8 @@ def parse_data( return data_input.ForbiddenDataInput(input) DataClass = PREFIX_MATCHES.get(prefix) if DataClass is not None: + if DataClass is tally.Tally: + return tally.Tally.from_input(input, jit_parse=jit_parse) if issubclass(DataClass, montepy.data_inputs.cell_modifier.CellModifierInput): return DataClass(input, problem=problem, jit_parse=jit_parse) return DataClass(input, jit_parse=jit_parse) diff --git a/montepy/data_inputs/tally.py b/montepy/data_inputs/tally.py new file mode 100644 index 000000000..2f1f556c5 --- /dev/null +++ b/montepy/data_inputs/tally.py @@ -0,0 +1,1435 @@ +# Copyright 2024, Battelle Energy Alliance, LLC All Rights Reserved. +from __future__ import annotations +import copy +from typing import Generator + +import montepy +from montepy.cells import Cells +from montepy.surface_collection import Surfaces +from montepy.data_inputs.data_input import DataInputAbstract +from montepy.data_inputs import tally_multiplier +from montepy.data_inputs.tally_type import Score, TallyType +from montepy.exceptions import MalformedInputError, NumberConflictError +from montepy.input_parser.tally_parser import TallyParser +from montepy.input_parser import syntax_node +from montepy.numbered_mcnp_object import Numbered_MCNP_Object +import montepy.types as ty +from montepy.utilities import * +from montepy.mcnp_object import InitInput + +_TALLY_TYPE_MODULUS = 10 + + +def _make_value_node(value_type, default, padding=" ", never_pad=False): + """Build a fresh :class:`~montepy.input_parser.syntax_node.ValueNode`. + + Mirrors :func:`montepy.mcnp_object.MCNP_Object._generate_default_node`, + duplicated here since ``FlatGroup``/``PathGroup`` aren't ``MCNP_Object`` + subclasses and so don't inherit that helper. + """ + padding_node = syntax_node.PaddingNode(padding) if padding else None + if default is None: + return syntax_node.ValueNode(default, value_type, padding_node, never_pad) + return syntax_node.ValueNode(str(default), value_type, padding_node, never_pad) + + +class LatticeIndex: + """A lattice element index ``[i j k]`` in a tally path specification. + + .. versionadded:: 1.6.0b2 + + Parameters + ---------- + dimensions : list + List of :class:`int` (single element index) or + ``tuple[int, int]`` (range ``i1:i2``). + """ + + __slots__ = ("_dimensions",) + + @args_checked + def __init__(self, dimensions: list[ty.Integral | tuple[ty.Integral, ty.Integral]]): + self._dimensions = list(dimensions) + + @property + def dimensions(self): + """The list of indices or (start, end) ranges.""" + return list(self._dimensions) + + def __repr__(self): + return f"LatticeIndex({self._dimensions})" + + +class TallyGroup: + """Abstract base for a tally scoring group. + + .. versionadded:: 1.6.0b2 + """ + + def __contains__(self, item) -> bool: + raise NotImplementedError + + +class Filter: + """Abstract analog of an OpenMC-style tally filter. + + .. versionadded:: 1.6.0b2 + """ + + +class ParticleFilter(Filter): + """Filters a tally to the particle types in its classifier (e.g. ``:n,p``). + + A thin wrapper around the underlying + :class:`~montepy.input_parser.syntax_node.ParticleNode`, mirroring + :class:`~montepy.Mode`'s design: particle *membership* is what matters, not + order. Two filters with the same particles compare equal no matter what order + they were written in (``:n,p`` == ``:p,n``); order is only ever meaningful when + the underlying node formats itself back to MCNP text. + + .. versionadded:: 1.6.0b2 + + Parameters + ---------- + particles : montepy.input_parser.syntax_node.ParticleNode, list[montepy.Particle], set[montepy.Particle] + The parsed node backing this filter's particles, or a plain collection of + particles to build one from. + """ + + __slots__ = ("_node",) + + @args_checked + def __init__( + self, + particles: ( + syntax_node.ParticleNode | list[montepy.Particle] | set[montepy.Particle] + ), + ): + if isinstance(particles, syntax_node.ParticleNode): + self._node = particles + else: + token = ",".join(p.value for p in particles) + self._node = syntax_node.ParticleNode("particle_filter", token) + + @property + def particles(self) -> set[montepy.Particle]: + """The particles this tally is restricted to.""" + return set(self._node.particles) + + def __eq__(self, other): + if not isinstance(other, ParticleFilter): + return NotImplemented + return frozenset(self.particles) == frozenset(other.particles) + + def __repr__(self): + return f"ParticleFilter({self._node.format()!r})" + + +class SpatialFilter(Filter): + """Filters a tally to its scoring bins (cells, surfaces, or paths). + + A thin wrapper around a :class:`~montepy.Tally`'s :attr:`~montepy.Tally.groups`. + + .. versionadded:: 1.6.0b2 + + Parameters + ---------- + groups : list[TallyGroup] + The scoring bins. + """ + + __slots__ = ("_groups",) + + @args_checked + def __init__(self, groups: list[TallyGroup]): + self._groups = list(groups) + + @property + def groups(self): + """The scoring bins (:class:`~montepy.data_inputs.tally.TallyGroup` objects) this filter covers.""" + return list(self._groups) + + def __eq__(self, other): + return isinstance(other, SpatialFilter) and self._groups == other._groups + + def __repr__(self): + return f"SpatialFilter({self._groups})" + + +class FlatGroup(TallyGroup): + """A flat list of cells/surfaces to score over. + + Used as both a top-level bin and as an individual level in a + :class:`PathGroup` chain. + + .. versionadded:: 1.6.0b2 + + Parameters + ---------- + numbers : list[int] + Cell or surface numbers. + lattice_indices : list[LatticeIndex | None], optional + Lattice indices parallel to ``numbers``. + is_grouped : bool + ``True`` = parenthesized union (one averaged bin); + ``False`` = separate bins. + universe_spec : int, optional + Universe number if ``U=N`` syntax was used. + """ + + __slots__ = ( + "_old_numbers", + "_lattice_indices", + "_is_grouped", + "_universe_spec", + "_cells_or_surfaces", + "_node", + "_number_nodes", + ) + + @args_checked + def __init__( + self, + numbers: list[ty.Integral], + lattice_indices: list[LatticeIndex | None] | None = None, + *, + is_grouped: bool, + universe_spec: ty.Integral | None = None, + node: syntax_node.SyntaxNodeBase | None = None, + number_nodes: list | None = None, + ): + self._old_numbers = list(numbers) + self._lattice_indices = lattice_indices or [None] * len(self._old_numbers) + self._is_grouped = is_grouped + self._universe_spec = universe_spec + self._cells_or_surfaces = [] + self._node = node + self._number_nodes = list(number_nodes) if number_nodes else [] + + @property + def old_numbers(self): + """The original cell/surface numbers as read.""" + return list(self._old_numbers) + + @property + def is_grouped(self): + """``True`` if entries form a union bin (parenthesized).""" + return self._is_grouped + + @property + def universe_spec(self): + """Universe number if ``U=N`` syntax was used, else ``None``.""" + return self._universe_spec + + @property + def lattice_indices(self): + """The :class:`LatticeIndex` for each number, parallel to :attr:`old_numbers`. + + .. versionadded:: 1.6.0b3 + """ + return list(self._lattice_indices) + + @property + def cells_or_surfaces(self): + """The resolved :class:`~montepy.Cell`/:class:`~montepy.Surface` objects + this group covers, once linked to a problem. Empty if not yet linked. + + .. versionadded:: 1.6.0b3 + """ + return list(self._cells_or_surfaces) + + def _current_numbers(self): + """The numbers to write out: live cell/surface numbers if linked, + else the numbers as originally read.""" + if self._cells_or_surfaces: + return [obj.number for obj in self._cells_or_surfaces] + return list(self._old_numbers) + + @property + def node(self): + """The syntax node for this group. + + If this was generated by parsing, that node (patched to reflect + current numbers) is returned. If this was created from scratch, a + new node is generated the first time this is accessed. + + .. versionadded:: 1.6.0b3 + """ + self._ensure_has_node() + self._update_node() + return self._node + + def _ensure_has_node(self): + if self._node is not None: + return + body = syntax_node.ListNode("flat group body") + if self._universe_spec is not None: + body.append(_make_value_node(str, "u", padding=None)) + body.append(syntax_node.PaddingNode("=")) + spec_node = _make_value_node(int, self._universe_spec, padding=None) + body.append(spec_node) + numbers = self._current_numbers() + self._number_nodes = [] + for i, num in enumerate(numbers): + num_node = _make_value_node(int, num, padding=" ") + self._number_nodes.append(num_node) + body.append(num_node) + idx = self._lattice_indices[i] if i < len(self._lattice_indices) else None + if idx is not None: + body.append(_build_lattice_node(idx)) + if self._is_grouped: + wrapped = syntax_node.ListNode("tally group") + wrapped.append(syntax_node.PaddingNode("(")) + for n in body.nodes: + wrapped.append(n) + wrapped.append(syntax_node.PaddingNode(")")) + self._node = wrapped + else: + self._node = body + + def _update_node(self): + numbers = self._current_numbers() + for num_node, num in zip(self._number_nodes, numbers): + if num_node.value != num: + num_node.value = num + + def __contains__(self, item) -> bool: + if self._cells_or_surfaces: + return item in self._cells_or_surfaces + for num in self._old_numbers: + if hasattr(item, "old_number") and item.old_number == num: + return True + if hasattr(item, "number") and item.number == num: + return True + return False + + def __repr__(self): + return f"FlatGroup({self._old_numbers}, grouped={self._is_grouped})" + + +class PathGroup(TallyGroup): + """A universe-path group for repeated-structures tallies. + + .. versionadded:: 1.6.0b2 + + Parameters + ---------- + levels : list[FlatGroup] + Levels in the ``<`` chain, innermost (scored) first. + """ + + __slots__ = ("_levels", "_node") + + @args_checked + def __init__( + self, + levels: list[FlatGroup], + node: syntax_node.SyntaxNodeBase | None = None, + ): + self._levels = list(levels) + self._node = node + + @property + def levels(self): + """FlatGroup levels, innermost (scored) first.""" + return list(self._levels) + + @property + def node(self): + """The syntax node for this path group. + + .. versionadded:: 1.6.0b3 + """ + self._ensure_has_node() + self._update_node() + return self._node + + def _ensure_has_node(self): + if self._node is not None: + for level in self._levels: + level._ensure_has_node() + return + body = syntax_node.ListNode("flat group body") + for i, level in enumerate(self._levels): + level._ensure_has_node() + level._update_node() + if isinstance(level.node, syntax_node.ListNode): + for n in level.node.nodes: + body.append(n) + else: + body.append(level.node) + if i != len(self._levels) - 1: + body.append(_make_value_node(str, "<", padding=" ")) + wrapped = syntax_node.ListNode("tally group") + wrapped.append(syntax_node.PaddingNode("(")) + for n in body.nodes: + wrapped.append(n) + wrapped.append(syntax_node.PaddingNode(")")) + self._node = wrapped + + def _update_node(self): + for level in self._levels: + level._update_node() + + @args_checked + def inside( + self, + *cells_or_surfaces: montepy.Cell | montepy.Surface, + lattice: list[ty.Integral] | None = None, + ) -> PathGroup: + """Append an outer level and return self for chaining. + + Parameters + ---------- + cells_or_surfaces : Cell | Surface + Objects at this level. + lattice : list[int], optional + Lattice index dimensions for the first element. + + Returns + ------- + PathGroup + ``self``, for method chaining. + """ + numbers = [obj.number for obj in cells_or_surfaces] + indices = [None] * len(numbers) + if lattice is not None and numbers: + indices[0] = LatticeIndex(lattice) + is_grouped = len(cells_or_surfaces) > 1 + level = FlatGroup(numbers, indices, is_grouped=is_grouped) + level._cells_or_surfaces = list(cells_or_surfaces) + self._levels.append(level) + return self + + def __contains__(self, item) -> bool: + if not self._levels: + return False + return item in self._levels[0] + + def __repr__(self): + return f"PathGroup(levels={len(self._levels)})" + + +def _build_lattice_node(index: LatticeIndex) -> syntax_node.ListNode: + """Build a fresh ``ListNode("lattice phrase")`` from a :class:`LatticeIndex`. + + Inverse of :func:`_parse_lattice_phrase`. Freshly-built dimensions are + always space-separated (the common ``[i j k]`` single-element-index + form); the data model doesn't distinguish that from a comma-separated + list once parsed, so there's no lossless "original style" to preserve + here anyway. + """ + node = syntax_node.ListNode("lattice phrase") + node.append(syntax_node.PaddingNode("[")) + dims = index.dimensions + for i, dim in enumerate(dims): + if isinstance(dim, tuple): + item = syntax_node.ListNode("lattice range") + item.append(_make_value_node(int, dim[0], padding=None)) + item.append(syntax_node.PaddingNode(":")) + item.append(_make_value_node(int, dim[1], padding=None)) + node.append(item) + else: + node.append(_make_value_node(int, dim, padding=None)) + if i != len(dims) - 1: + node.append(syntax_node.PaddingNode(" ")) + node.append(syntax_node.PaddingNode("]")) + return node + + +def _parse_lattice_phrase(lattice_node) -> LatticeIndex: + """Parse a ``ListNode("lattice phrase")`` into a :class:`LatticeIndex`.""" + dimensions = [] + for n in lattice_node.nodes: + if isinstance(n, syntax_node.PaddingNode): + continue + if isinstance(n, syntax_node.ListNode) and n.name == "lattice range": + vals = [m.value for m in n.nodes if isinstance(m, syntax_node.ValueNode)] + if len(vals) >= 2: + dimensions.append((int(vals[0]), int(vals[1]))) + elif isinstance(n, syntax_node.ValueNode) and isinstance(n.value, (int, float)): + dimensions.append(int(n.value)) + return LatticeIndex(dimensions) + + +def _extract_numbers_with_lattice(nodes): + """Pair each numeric ValueNode with its immediately following lattice phrase. + + Returns ``(numbers, lattice_indices, number_nodes)`` where + ``lattice_indices[i]`` is a :class:`LatticeIndex` or ``None``, and + ``number_nodes[i]`` is the real, original :class:`~montepy.input_parser.syntax_node.ValueNode` + parsed for that number (kept so it can be patched in place later instead + of rebuilt, preserving the original formatting/whitespace). + + If an MCNP shortcut (e.g. ``3i``) is present, ``number_nodes`` is left + empty for the whole group: a shortcut's "virtual" expanded value nodes + aren't meant to be formatted/patched individually (only the + :class:`~montepy.input_parser.syntax_node.ShortcutNode` itself knows how + to compress back to e.g. ``3i``), so live-renumbering isn't supported for + a group that used one -- it keeps its original compressed text as-is. + """ + numbers = [] + lattice_indices = [] + number_nodes = [] + has_shortcut = False + i = 0 + while i < len(nodes): + n = nodes[i] + if isinstance(n, syntax_node.ShortcutNode): + has_shortcut = True + for inner in n.nodes: + if isinstance(inner.value, (int, float)) and not isinstance( + inner.value, bool + ): + numbers.append(int(inner.value)) + lattice_indices.append(None) + elif ( + isinstance(n, syntax_node.ValueNode) + and isinstance(n.value, (int, float)) + and not isinstance(n.value, bool) + ): + numbers.append(int(n.value)) + number_nodes.append(n) + if ( + i + 1 < len(nodes) + and isinstance(nodes[i + 1], syntax_node.ListNode) + and nodes[i + 1].name == "lattice phrase" + ): + lattice_indices.append(_parse_lattice_phrase(nodes[i + 1])) + i += 2 + continue + else: + lattice_indices.append(None) + i += 1 + if has_shortcut: + number_nodes = [] + return numbers, lattice_indices, number_nodes + + +def _extract_universe_spec_from_nodes(nodes): + """Extract universe number from nodes that may contain a universe_phrase ListNode.""" + for n in nodes: + if isinstance(n, syntax_node.ListNode) and n.name == "universe phrase": + for m in n.nodes: + if isinstance(m, syntax_node.ValueNode) and isinstance( + m.value, (int, float) + ): + return int(m.value) + elif isinstance(n, syntax_node.ListNode) and n.name == "tally group": + inner = list(n.nodes)[1:-1] + result = _extract_universe_spec_from_nodes(inner) + if result is not None: + return result + return None + + +def _parse_body_segment(nodes, *, is_grouped, node=None) -> FlatGroup: + numbers, lattice_indices, number_nodes = _extract_numbers_with_lattice(nodes) + universe_spec = _extract_universe_spec_from_nodes(nodes) + if node is None and nodes: + wrapper = syntax_node.ListNode("flat group body") + for n in nodes: + wrapper.append(n) + node = wrapper + return FlatGroup( + numbers, + lattice_indices, + is_grouped=is_grouped, + universe_spec=universe_spec, + node=node, + number_nodes=number_nodes, + ) + + +def _parse_segment_as_level(seg) -> FlatGroup: + """Parse a path segment (between ``<`` separators) into a :class:`FlatGroup` level.""" + non_pad = [n for n in seg if not isinstance(n, syntax_node.PaddingNode)] + if ( + len(non_pad) == 1 + and isinstance(non_pad[0], syntax_node.ListNode) + and non_pad[0].name == "tally group" + ): + inner_body = list(non_pad[0].nodes)[1:-1] + return _parse_body_segment(inner_body, is_grouped=True, node=non_pad[0]) + return _parse_body_segment(seg, is_grouped=False) + + +def _parse_tally_group_node(group_node) -> TallyGroup: + nodes = list(group_node.nodes) + for i, n in enumerate(nodes): + if ( + isinstance(n, syntax_node.ValueNode) + and n.padding is None + and i + 1 < len(nodes) + and not isinstance(nodes[i + 1], syntax_node.PaddingNode) + ): + n.never_pad = True + body = nodes[1:-1] if len(nodes) >= 2 else nodes + + path_sep_indices = [ + i + for i, n in enumerate(body) + if ( + isinstance(n, syntax_node.ValueNode) + and isinstance(n.value, str) + and n.value.strip() == "<" + ) + ] + + if not path_sep_indices: + return _parse_body_segment(body, is_grouped=True, node=group_node) + + segments = [] + start = 0 + for sep_i in path_sep_indices: + segments.append(body[start:sep_i]) + start = sep_i + 1 + segments.append(body[start:]) + + return PathGroup( + [_parse_segment_as_level(seg) for seg in segments], node=group_node + ) + + +def _parse_tally_numbers(tally_numbers_node) -> list[TallyGroup]: + groups = [] + for node in tally_numbers_node: + if isinstance(node, syntax_node.ValueNode): + v = node.value + if v is None: + continue + if isinstance(v, (int, float)) and not isinstance(v, bool): + wrapper = syntax_node.ListNode("flat group body") + wrapper.append(node) + groups.append( + FlatGroup( + [int(v)], is_grouped=False, node=wrapper, number_nodes=[node] + ) + ) + elif isinstance(node, syntax_node.ListNode) and node.name == "tally group": + groups.append(_parse_tally_group_node(node)) + return groups + + +def _link_multiplier_to_tally(self, fm): + fm._link_to_parent(self) + if self._problem is not None: + self._problem.tallies.append(fm) + + +class Tally(DataInputAbstract, Numbered_MCNP_Object): + """Base class for MCNP F-card tallies (F1, F2, F4, F5, F6, F7, F8). + + Use :meth:`from_input` as a factory to create the appropriate subclass + when reading from a file. + + .. versionadded:: 1.6.0b2 + """ + + _POINTER_ATTRS = set() + _DEFAULT_SCORES = () + _KEYS_TO_PRESERVE = {"_multiplier"} + + @staticmethod + def _parser(): + return TallyParser() + + def _init_blank(self): + super()._init_blank() + self._old_number = self._generate_default_node(int, -1) + self._groups = [] + self._include_total = False + self._multiplier = None + + def _jit_light_init(self, input): + super()._jit_light_init(input) + self._old_number = self._input_number + + def _parse_tree(self): + super()._parse_tree() + num = self._input_number + self._old_number = copy.deepcopy(num) + self._number = num + self._parse_tally_body() + + def _generate_default_tree(self, **kwargs): + ret = {} + ret["start_pad"] = syntax_node.PaddingNode() + ret["classifier"] = syntax_node.ClassifierNode() + ret["classifier"].prefix = syntax_node.ValueNode( + self._class_prefix().upper(), str, padding=None, never_pad=True + ) + # A non-negative placeholder: ValueNode._reverse_engineer_formatting + # reserves a leading sign column for any token starting with "-", + # which would otherwise permanently corrupt this node's formatting + # once a real (positive) tally number is assigned to it. + ret["classifier"].number = self._generate_default_node(int, 1) + ret["keyword"] = syntax_node.ValueNode(None, str, padding=None) + tally_numbers = syntax_node.ListNode("tally numbers") + end_node = syntax_node.ValueNode(None, str) + ret["data"] = syntax_node.SyntaxNode( + "tally list", {"tally": tally_numbers, "end": end_node} + ) + ret["parameters"] = syntax_node.ParametersNode() + self._tree = syntax_node.SyntaxNode("blank data tree", ret) + + @args_checked + def __init__( + self, + input: InitInput = None, + number: ty.PositiveInt = None, + *, + jit_parse: bool = True, + ): + Numbered_MCNP_Object.__init__(self, input, number, jit_parse=jit_parse) + + @staticmethod + def _class_prefix() -> str: + return "f" + + @staticmethod + def _has_number() -> bool: + return True + + @staticmethod + def _has_classifier() -> int: + return 1 + + @staticmethod + def _parent_collections(): + return () + + def _parse_tally_body(self): + if self._input is None: + return + num = self._input_number.value + try: + TallyType(num % _TALLY_TYPE_MODULUS) + except ValueError as e: + raise MalformedInputError(self._input, f"Invalid tally type digit: {e}") + tally_list = self._tree["data"] + end_node = tally_list["end"] + self._include_total = ( + end_node.value is not None and str(end_node.value).upper() == "T" + ) + self._groups = _parse_tally_numbers(tally_list["tally"]) + + def _number_validator(self, number): + tally_type = getattr(type(self), "_TALLY_TYPE", None) + if tally_type is not None and number % _TALLY_TYPE_MODULUS != tally_type.value: + raise ValueError( + f"Cannot change tally type via number setter; " + f"expected last digit {tally_type.value}, " + f"got {number % _TALLY_TYPE_MODULUS}." + ) + super()._number_validator(number) + if self._multiplier is not None: + self._multiplier.number = number + + @make_prop_val_node("_old_number") + def old_number(self): + """The tally number as read from the input file.""" + pass + + @property + @needs_full_ast + def tally_type(self) -> TallyType | None: + """The MCNP tally type (e.g. ``TallyType.CELL_FLUX`` for F4).""" + return getattr(type(self), "_TALLY_TYPE", None) + + @property + @needs_full_ast + def groups(self) -> list[TallyGroup]: + """The list of :class:`~montepy.data_inputs.tally.TallyGroup` objects defining what is scored.""" + return list(self._groups) + + @property + @needs_full_ast + def include_total(self) -> bool: + """``True`` if a total bin (T) is appended.""" + return self._include_total + + @include_total.setter + @args_checked + @needs_full_cst + def include_total(self, value: bool): + self._include_total = value + + @make_prop_pointer( + "_multiplier", + tally_multiplier.TallyMultiplier, + validator=_link_multiplier_to_tally, + ) + def multiplier(self) -> tally_multiplier.TallyMultiplier: + """The ``FM`` tally-multiplier card linked to this tally, if any. + + Returns + ------- + TallyMultiplier + """ + pass + + @property + @needs_full_ast + def scores(self) -> list[Score] | list[tally_multiplier.MultiplierScore]: + """The physical quantities this tally scores, e.g. ``[Score.FLUX]`` for F4. + + This is just the quantity implied by the tally type digit, unless an + ``FM`` tally-multiplier card is linked (see :attr:`multiplier`), in + which case this returns one :class:`~montepy.data_inputs.tally_multiplier.MultiplierScore` + per output bin the multiplier defines instead. + """ + if self.multiplier is not None: + return [score for bin_ in self.multiplier.bins for score in bin_.scores] + return list(self._DEFAULT_SCORES) + + @property + @needs_full_ast + def filters(self) -> list[Filter]: + """A shallow analog of OpenMC's tally filters. + + Defaults to a :class:`~montepy.data_inputs.tally.ParticleFilter` (from + :attr:`particle_classifiers`) and a + :class:`~montepy.data_inputs.tally.SpatialFilter` (from :attr:`groups`), + whichever are present. + """ + filters = [] + if self.particle_classifiers: + filters.append(ParticleFilter(self._classifier.particles)) + if self._groups: + filters.append(SpatialFilter(self._groups)) + return filters + + @needs_full_ast + def __contains__(self, item) -> bool: + for group in self._groups: + if item in group: + return True + return False + + @staticmethod + def _dispatch_class(input, num: int) -> type[Tally]: + """The :class:`Tally` subclass for a tally number.""" + try: + tally_type = TallyType(num % _TALLY_TYPE_MODULUS) + except ValueError as e: + raise MalformedInputError( + input, f"Tally type digit {num % _TALLY_TYPE_MODULUS} is not valid." + ) from e + # _TALLY_TYPE_MAP's keys are exactly TallyType's members (both + # defined by hand in lockstep in this module), so this can never miss. + return _TALLY_TYPE_MAP[tally_type] + + @classmethod + def from_input(cls, input, *, jit_parse: bool = True) -> Tally: + """Factory: create the appropriate :class:`Tally` subclass from an input. + + Parameters + ---------- + input : Input | str + The raw MCNP input object. + jit_parse : bool + Whether to defer full parsing. + + Returns + ------- + Tally + An instance of the correct subclass for the tally type digit. + """ + try: + bare_tree = Tally._peek_light_parse(input) + number_node = bare_tree.nodes["classifier"].number + if number_node is None: + raise ValueError("Tally classifier has no number.") + subclass = Tally._dispatch_class(input, number_node.value) + except Exception: + # The JIT light parser isn't fully robust and can fail on valid + # syntax. Fall back to building a real Tally: its own + # JIT-with-fallback-to-full-parse handling in _parse_input will + # reliably determine the number instead of guessing, and gives + # proper file/line context on error. + base = Tally(input, jit_parse=True) + subclass = Tally._dispatch_class(input, base._number.value) + + return subclass(input, jit_parse=jit_parse) + + def link_to_problem(self, problem, *, deepcopy=False): + super().link_to_problem(problem) + + def _update_values(self): + tally_numbers_node = self._tree["data"]["tally"] + tally_numbers_node.nodes.clear() + for group in self._groups: + node = group.node + if isinstance(node, syntax_node.ListNode) and node.name != "tally group": + for n in node.nodes: + tally_numbers_node.nodes.append(n) + else: + tally_numbers_node.nodes.append(node) + end_node = self._tree["data"]["end"] + end_node.value = "T" if self._include_total else None + + @staticmethod + def _align_to_type(tally_type: TallyType, start: int) -> int: + """The smallest number ``>= start`` whose last digit matches ``tally_type``.""" + aligned = start - (start % 10) + tally_type.value + if aligned < start: + aligned += 10 + return aligned + + def _next_number_for_type( + self, tally_type: TallyType, starting_number, step + ) -> int: + """Finds the next free tally number matching ``tally_type``'s digit. + + Note + ---- + This probes with :meth:`~montepy.numbered_object_collection.NumberedObjectCollection.check_number` + rather than delegating to :meth:`~montepy.numbered_object_collection.NumberedObjectCollection.request_number`, + because that method tracks a single collection-wide + ``_last_assigned_number`` ratchet that isn't digit-aware: a prior + request for one tally-type digit (e.g. ``clone()`` landing on 124) + pushes that ratchet past 124, so a later request for a *different* + digit (e.g. ``clone_as`` targeting type 6) would start its search + from >134 instead of the correctly-aligned 6, and drift to a number + that still doesn't end in 6. + """ + collection = self._problem.tallies if self._problem else None + if collection is not None: + start = ( + starting_number + if starting_number is not None + else collection.starting_number + ) + step = step if step is not None else collection.step + else: + start = starting_number if starting_number is not None else 1 + step = step if step is not None else 1 + candidate = self._align_to_type(tally_type, start) + while True: + if collection is not None: + try: + collection.check_number(candidate) + return candidate + except NumberConflictError: + pass + elif candidate != self.number: + return candidate + candidate += step * 10 + + @staticmethod + def _tally_category(cls: type[Tally]) -> type[Tally] | None: + """Which of {SurfaceTally, CellTally, DetectorTally} ``cls`` belongs to.""" + for category in (SurfaceTally, CellTally, DetectorTally): + if issubclass(cls, category): + return category + return None + + @args_checked + @needs_full_cst + def clone( + self, + starting_number: ty.PositiveInt = None, + step: ty.PositiveInt = None, + ) -> Tally: + """Clone this tally with a new number. + + Note that the clone does **not** carry over a linked ``FM`` + multiplier (see :attr:`multiplier`) -- a multiplier is a companion + card tied to this exact tally number, not something that + meaningfully transfers to a renumbered copy. + + See :meth:`~montepy.numbered_mcnp_object.Numbered_MCNP_Object.clone`. + """ + ret = copy.deepcopy(self) + ret._multiplier = None + new_number = self._next_number_for_type(self.tally_type, starting_number, step) + if self._problem: + ret.link_to_problem(self._problem) + ret.number = new_number + self._problem.tallies.append(ret) + else: + ret.number = new_number + return ret + + @args_checked + @needs_full_cst + def clone_as( + self, + new_type: TallyType | type[Tally], + starting_number: ty.PositiveInt = None, + step: ty.PositiveInt = None, + ) -> Tally: + """Clone this tally as a different tally type, keeping the same scoring geometry. + + For example, this can turn an F4 cell-flux tally into an F6 + energy-deposition tally scoring the same cells: + + .. code-block:: python + + from montepy.data_inputs.tally import F6Tally + from montepy.data_inputs.tally_type import TallyType + + heating = flux_tally.clone_as(F6Tally) + # or, equivalently: + heating = flux_tally.clone_as(TallyType.ENERGY_DEPOSITION) + + Only conversions within the same tally category are allowed: + F1/F2 (surface-based) convert freely among each other, as do + F4/F6/F7/F8 (cell-based); F5 (point/ring detector) has no + cell/surface geometry to carry over and can't be converted to or + from. + + Parameters + ---------- + new_type : TallyType, type[Tally] + The target tally type, either as a :class:`TallyType` member or + as a :class:`Tally` subclass (e.g. ``montepy.F6Tally``). + starting_number : int + The starting number to request for the new object's number. + step : int + The step size to use to find a new valid number. + + Returns + ------- + Tally + A new tally of the requested type, with the same scoring groups. + + Note + ---- + The clone does **not** carry over a linked ``FM`` multiplier (see + :attr:`multiplier`) -- a multiplier is a companion card tied to this + exact tally number, not something that meaningfully transfers to a + retyped/renumbered copy. + """ + if isinstance(new_type, TallyType): + # _TALLY_TYPE_MAP's keys are exactly TallyType's members (both + # defined by hand in lockstep in this module), so this can never + # miss. + target_cls = _TALLY_TYPE_MAP[new_type] + elif isinstance(new_type, type) and issubclass(new_type, Tally): + target_cls = new_type + else: + raise TypeError( + f"new_type must be a TallyType or a Tally subclass, got {new_type!r}." + ) + + if self._tally_category(type(self)) != self._tally_category(target_cls): + raise ValueError( + f"Cannot clone a {type(self).__name__} (tally type " + f"{self.tally_type}) as a {target_cls.__name__} (tally type " + f"{target_cls._TALLY_TYPE}); incompatible tally categories." + ) + + if target_cls is type(self): + return self.clone(starting_number, step) + + ret = copy.deepcopy(self) + ret.__class__ = target_cls + ret._multiplier = None + new_number = self._next_number_for_type( + target_cls._TALLY_TYPE, starting_number, step + ) + if self._problem: + ret.link_to_problem(self._problem) + ret.number = new_number + self._problem.tallies.append(ret) + else: + ret.number = new_number + return ret + + def __str__(self): + try: + return f"{type(self).__name__}: {self.number}" + except Exception: + return f"{type(self).__name__}: (unparsed)" + + def __repr__(self): + try: + ttype = getattr(type(self), "_TALLY_TYPE", None) + ngroups = len(getattr(self, "_groups", [])) + return f"TALLY: {self.number}, type: {ttype}, groups: {ngroups}" + except Exception: + return "TALLY: (unparsed)" + + +class SurfaceTally(Tally): + """Intermediate class for tallies that score on surfaces (F1, F2). + + .. versionadded:: 1.6.0b2 + """ + + def _init_blank(self): + super()._init_blank() + self._surfaces = Surfaces() + + @property + @needs_full_ast + def surfaces(self) -> Surfaces: + """The surfaces this tally scores over.""" + return self._surfaces + + @args_checked + @needs_full_cst + def add_surface(self, surface: montepy.Surface) -> None: + """Add a single surface as a separate scoring bin. + + Parameters + ---------- + surface : Surface + The surface to add. + """ + group = FlatGroup([surface.number], is_grouped=False) + group._cells_or_surfaces = [surface] + self._groups.append(group) + if surface not in self._surfaces: + self._surfaces.append(surface) + + @args_checked + @needs_full_cst + def add_group(self, surfaces: list[montepy.Surface] | set[montepy.Surface]) -> None: + """Add surfaces as a single union (averaged) bin. + + Parameters + ---------- + surfaces : list[Surface], set[Surface] + The surfaces to group. + """ + surfaces = list(surfaces) + numbers = [s.number for s in surfaces] + group = FlatGroup(numbers, is_grouped=True) + group._cells_or_surfaces = list(surfaces) + self._groups.append(group) + for s in surfaces: + if s not in self._surfaces: + self._surfaces.append(s) + + @args_checked + @needs_full_cst + def add_path_group(self, *surfaces: montepy.Surface) -> PathGroup: + """Add a universe-path group rooted at the given surfaces. + + Returns the :class:`~montepy.data_inputs.tally.PathGroup` for chaining via :meth:`~montepy.data_inputs.tally.PathGroup.inside`. + + Parameters + ---------- + surfaces : Surface + The innermost-level surfaces. + + Returns + ------- + PathGroup + The new path group (already appended). + """ + numbers = [s.number for s in surfaces] + is_grouped = len(surfaces) > 1 + first_level = FlatGroup(numbers, is_grouped=is_grouped) + first_level._cells_or_surfaces = list(surfaces) + pg = PathGroup([first_level]) + self._groups.append(pg) + return pg + + @args_checked + @needs_full_cst + def remove_surface(self, surface: montepy.Surface) -> None: + """Remove the single-surface scoring bin added via :meth:`add_surface`. + + Parameters + ---------- + surface : Surface + The surface to remove. + """ + for group in self._groups: + if ( + isinstance(group, FlatGroup) + and not group.is_grouped + and group.cells_or_surfaces == [surface] + ): + self.remove_group(group) + return + raise ValueError( + f"No single-surface scoring group found for surface {surface.number}." + ) + + @args_checked + @needs_full_cst + def remove_group(self, group: TallyGroup) -> None: + """Remove a scoring group previously added via ``add_surface``/``add_group``/``add_path_group``. + + A surface is only dropped from :attr:`surfaces` if no other + remaining group still references it. + + Parameters + ---------- + group : TallyGroup + The group to remove. + """ + self._groups.remove(group) + for surface in list(self._surfaces): + if surface not in self: + self._surfaces.remove(surface) + + def link_to_problem(self, problem, *, deepcopy=False): + super().link_to_problem(problem) + if problem is not None and not hasattr(self, "_not_parsed"): + # Rebuild from scratch: a deepcopy (e.g. from clone()) carries + # stale copied Surface objects that must be discarded, not + # merged with the live ones from `problem`. + self._surfaces = Surfaces() + for group in self._groups: + self._link_group_surfaces(group, problem) + + def _link_group_surfaces(self, group, problem): + if isinstance(group, FlatGroup): + group._cells_or_surfaces = [] + for num in group._old_numbers: + try: + # Use _surfaces directly to avoid triggering __relink_objs via the property. + surface = problem._surfaces[num] + except KeyError: + continue + group._cells_or_surfaces.append(surface) + try: + self._surfaces[surface.number] + except KeyError: + self._surfaces.append(surface) + elif isinstance(group, PathGroup): + for level in group._levels: + self._link_group_surfaces(level, problem) + + +class CellTally(Tally): + """Intermediate class for tallies that score in cells (F4, F6, F7, F8). + + .. versionadded:: 1.6.0b2 + """ + + def _init_blank(self): + super()._init_blank() + self._cells = Cells() + + @property + @needs_full_ast + def cells(self) -> Cells: + """The cells this tally scores in.""" + return self._cells + + @args_checked + @needs_full_cst + def add_cell(self, cell: montepy.Cell) -> None: + """Add a single cell as a separate scoring bin. + + Parameters + ---------- + cell : Cell + The cell to add. + """ + group = FlatGroup([cell.number], is_grouped=False) + group._cells_or_surfaces = [cell] + self._groups.append(group) + if cell not in self._cells: + self._cells.append(cell) + + @args_checked + @needs_full_cst + def add_group(self, cells: list[montepy.Cell] | set[montepy.Cell]) -> None: + """Add cells as a single union (averaged) bin. + + Parameters + ---------- + cells : list[Cell], set[Cell] + The cells to group. + """ + cells = list(cells) + numbers = [c.number for c in cells] + group = FlatGroup(numbers, is_grouped=True) + group._cells_or_surfaces = list(cells) + self._groups.append(group) + for c in cells: + if c not in self._cells: + self._cells.append(c) + + @args_checked + @needs_full_cst + def add_path_group(self, *cells: montepy.Cell) -> PathGroup: + """Add a universe-path group rooted at the given cells. + + Returns the :class:`~montepy.data_inputs.tally.PathGroup` for chaining via :meth:`~montepy.data_inputs.tally.PathGroup.inside`. + + Parameters + ---------- + cells : Cell + The innermost-level cells. + + Returns + ------- + PathGroup + The new path group (already appended). + """ + numbers = [c.number for c in cells] + is_grouped = len(cells) > 1 + first_level = FlatGroup(numbers, is_grouped=is_grouped) + first_level._cells_or_surfaces = list(cells) + pg = PathGroup([first_level]) + self._groups.append(pg) + return pg + + @args_checked + @needs_full_cst + def remove_cell(self, cell: montepy.Cell) -> None: + """Remove the single-cell scoring bin added via :meth:`add_cell`. + + Parameters + ---------- + cell : Cell + The cell to remove. + """ + for group in self._groups: + if ( + isinstance(group, FlatGroup) + and not group.is_grouped + and group.cells_or_surfaces == [cell] + ): + self.remove_group(group) + return + raise ValueError(f"No single-cell scoring group found for cell {cell.number}.") + + @args_checked + @needs_full_cst + def remove_group(self, group: TallyGroup) -> None: + """Remove a scoring group previously added via ``add_cell``/``add_group``/``add_path_group``. + + A cell is only dropped from :attr:`cells` if no other remaining + group still references it. + + Parameters + ---------- + group : TallyGroup + The group to remove. + """ + self._groups.remove(group) + for cell in list(self._cells): + if cell not in self: + self._cells.remove(cell) + + def link_to_problem(self, problem, *, deepcopy=False): + super().link_to_problem(problem) + if problem is not None and not hasattr(self, "_not_parsed"): + # Rebuild from scratch: a deepcopy (e.g. from clone()) carries + # stale copied Cell objects that must be discarded, not merged + # with the live ones from `problem`. + self._cells = Cells() + for group in self._groups: + self._link_group_cells(group, problem) + + def _link_group_cells(self, group, problem): + if isinstance(group, FlatGroup): + group._cells_or_surfaces = [] + for num in group._old_numbers: + try: + # Use _cells directly to avoid triggering __relink_objs via the property. + cell = problem._cells[num] + except KeyError: + continue + group._cells_or_surfaces.append(cell) + try: + self._cells[cell.number] + except KeyError: + self._cells.append(cell) + elif isinstance(group, PathGroup): + if group._levels: + self._link_group_cells(group._levels[0], problem) + + +class DetectorTally(Tally): + """F5: point/ring detector tally. + + .. versionadded:: 1.6.0b2 + """ + + _TALLY_TYPE = TallyType.DETECTOR + _DEFAULT_SCORES = (Score.FLUX,) + + +# ── Concrete subclasses ──────────────────────────────────────────────────────── + + +class SurfaceCurrentTally(SurfaceTally): + """F1: surface current tally. + + .. versionadded:: 1.6.0b2 + """ + + _TALLY_TYPE = TallyType.CURRENT + _DEFAULT_SCORES = (Score.CURRENT,) + + +class SurfaceFluxTally(SurfaceTally): + """F2: average surface flux tally. + + .. versionadded:: 1.6.0b2 + """ + + _TALLY_TYPE = TallyType.SURFACE_FLUX + _DEFAULT_SCORES = (Score.FLUX,) + + +class CellFluxTally(CellTally): + """F4: cell flux tally. + + .. versionadded:: 1.6.0b2 + """ + + _TALLY_TYPE = TallyType.CELL_FLUX + _DEFAULT_SCORES = (Score.FLUX,) + + +class EnergyDepositionTally(CellTally): + """F6: energy deposition tally. + + .. versionadded:: 1.6.0b2 + """ + + _TALLY_TYPE = TallyType.ENERGY_DEPOSITION + _DEFAULT_SCORES = (Score.ENERGY_DEPOSITION,) + + +class FissionEnergyDepositionTally(CellTally): + """F7: fission energy deposition tally. + + .. versionadded:: 1.6.0b2 + """ + + _TALLY_TYPE = TallyType.FISSION_ENERGY_DEPOSITION + _DEFAULT_SCORES = (Score.FISSION_ENERGY_DEPOSITION,) + + +class EnergyDetectorPulseTally(CellTally): + """F8: energy-detector pulse height tally. + + .. versionadded:: 1.6.0b2 + """ + + _TALLY_TYPE = TallyType.ENERGY_DETECTOR_PULSE + _DEFAULT_SCORES = (Score.PULSE_HEIGHT,) + + +_TALLY_TYPE_MAP: dict[TallyType, type[Tally]] = { + TallyType.CURRENT: SurfaceCurrentTally, + TallyType.SURFACE_FLUX: SurfaceFluxTally, + TallyType.CELL_FLUX: CellFluxTally, + TallyType.DETECTOR: DetectorTally, + TallyType.ENERGY_DEPOSITION: EnergyDepositionTally, + TallyType.FISSION_ENERGY_DEPOSITION: FissionEnergyDepositionTally, + TallyType.ENERGY_DETECTOR_PULSE: EnergyDetectorPulseTally, +} + +# ── Convenience aliases ──────────────────────────────────────────────────────── + +F1Tally = SurfaceCurrentTally +F2Tally = SurfaceFluxTally +F4Tally = CellFluxTally +F5Tally = DetectorTally +F6Tally = EnergyDepositionTally +F7Tally = FissionEnergyDepositionTally +F8Tally = EnergyDetectorPulseTally diff --git a/montepy/data_inputs/tally_multiplier.py b/montepy/data_inputs/tally_multiplier.py new file mode 100644 index 000000000..69364f9a9 --- /dev/null +++ b/montepy/data_inputs/tally_multiplier.py @@ -0,0 +1,1723 @@ +# Copyright 2024, Battelle Energy Alliance, LLC All Rights Reserved. +from __future__ import annotations +import copy +import warnings +from typing import Union + +import montepy +from montepy.data_inputs.data_input import DataInputAbstract +from montepy.data_inputs.tally_multiplier_type import ( + ReactionOperator, + SpecialMultiplier, +) +from montepy.exceptions import MalformedInputWarning +from montepy.input_parser.tally_parser import TallyParser +from montepy.input_parser import syntax_node +from montepy.numbered_mcnp_object import Numbered_MCNP_Object +import montepy.types as ty +from montepy.utilities import * +from montepy.mcnp_object import InitInput + +_SPECIAL_KIND_MAP = { + -1: SpecialMultiplier.INVERSE_WEIGHT, + -2: SpecialMultiplier.INVERSE_VELOCITY, + -3: SpecialMultiplier.FIRST_INTERACTION_XS, +} +_SPECIAL_KIND_MAP_INVERSE = {v: k for k, v in _SPECIAL_KIND_MAP.items()} + + +def _make_value_node(value_type, default, padding=" ", never_pad=False): + padding_node = syntax_node.PaddingNode(padding) if padding else None + if default is None: + return syntax_node.ValueNode(default, value_type, padding_node, never_pad) + return syntax_node.ValueNode(str(default), value_type, padding_node, never_pad) + + +def _coerce(value) -> ReactionExpression: + """Wrap a bare ``int`` reaction number in a :class:`Reaction`, or pass through.""" + if isinstance(value, ReactionExpression): + return value + if isinstance(value, ty.Integral): + return Reaction(int(value)) + raise TypeError(f"Cannot combine a reaction expression with {value!r}.") + + +class ReactionExpression: + """A binary expression tree for one FM reaction list. + + Modeled on :class:`montepy.surfaces.half_space.HalfSpace` + (``left``/``operator``/``right``, ``&``/``|``/``~`` there → ``*``/``+``/``-`` + here). Python's own operator precedence (``*`` binds tighter than + ``+``/``-``) gives MCNP's "multiply first" reaction-list rule for free: + ``Reaction(16) * Reaction(103) + Reaction(104)`` builds ``(16*103) + 104`` + with no custom precedence-climbing code. + + .. versionadded:: 1.6.0b2 + """ + + def __init__( + self, + left: ReactionExpression, + operator: ReactionOperator, + right: ReactionExpression, + ): + self._left = left + self._operator = operator + self._right = right + self._node = None + + @property + def node(self) -> syntax_node.ListNode: + """The syntax node for this reaction expression. + + .. versionadded:: 1.6.0b3 + """ + self._ensure_has_node() + return self._node + + def _ensure_has_node(self): + if self._node is not None: + return + node = syntax_node.ListNode("reaction expr") + for n in self.left.node.nodes: + node.append(n) + if self._operator != ReactionOperator.MULTIPLY: + symbol = ":" if self._operator == ReactionOperator.ADD else "#" + node.append(_make_value_node(str, symbol, padding=" ")) + for n in self.right.node.nodes: + node.append(n) + self._node = node + + @make_prop_pointer("_left") + def left(self): + """The left side of this expression.""" + pass + + @make_prop_pointer("_operator") + def operator(self): + """The :class:`~montepy.data_inputs.tally_multiplier_type.ReactionOperator` joining the two sides.""" + pass + + @make_prop_pointer("_right") + def right(self): + """The right side of this expression.""" + pass + + def __mul__(self, other) -> ReactionExpression: + return ReactionExpression(self, ReactionOperator.MULTIPLY, _coerce(other)) + + def __add__(self, other) -> ReactionExpression: + return ReactionExpression(self, ReactionOperator.ADD, _coerce(other)) + + def __sub__(self, other) -> ReactionExpression: + return ReactionExpression(self, ReactionOperator.SUBTRACT, _coerce(other)) + + def __rmul__(self, other) -> ReactionExpression: + # Not a plain alias to __mul__: that would put ``self`` on the left, + # reversing the operand order implied by e.g. ``16 * Reaction(103)``. + return ReactionExpression(_coerce(other), ReactionOperator.MULTIPLY, self) + + def __radd__(self, other) -> ReactionExpression: + return ReactionExpression(_coerce(other), ReactionOperator.ADD, self) + + def __rsub__(self, other) -> ReactionExpression: + return ReactionExpression(_coerce(other), ReactionOperator.SUBTRACT, self) + + def __rand__( + self, material: Union[ty.Integral, "montepy.Material"] + ) -> MultiplierSet: + """``material_or_number & reaction_expr`` -> a one-term :class:`MultiplierSet`. + + Defined here so both leaves and composite trees support it via + inheritance: ``mat1 & Reaction.CAPTURE`` and + ``26 & (Reaction.CAPTURE - Reaction.INELASTIC_SCATTER)`` both work, + building a single-reaction :class:`MultiplierSet` with + ``constant=1.0``. Scale it with ``*`` afterwards (see + ``MultiplierSet.__rmul__``) or drop it straight into a + :class:`~montepy.data_inputs.tally_multiplier.MultiplierBin`'s + ``terms`` list. + """ + return MultiplierSet(1.0, material, [self]) + + def __eq__(self, other): + if not isinstance(other, ReactionExpression): + return NotImplemented + return ( + self._left == other._left + and self._operator == other._operator + and self._right == other._right + ) + + def __repr__(self): + return f"ReactionExpression({self._left!r}, {self._operator}, {self._right!r})" + + +class Reaction(ReactionExpression): + """A leaf reaction number: a single ENDF (MT) or special (R) reaction. + + Does **not** call ``ReactionExpression.__init__`` — mirrors + ``UnitHalfSpace``, which holds independent leaf state rather than being a + degenerate composite node pointing at itself. Inherits + ``__mul__``/``__add__``/``__sub__``/``__rand__`` from + :class:`~montepy.data_inputs.tally_multiplier.ReactionExpression` unchanged; + only ``__eq__``/``__repr__`` need leaf-specific overrides. + + Common reaction numbers are available as ready-to-use class attributes, + e.g. ``Reaction.CAPTURE`` for MT 102, (n,gamma) radiative capture, so you + don't need to remember raw MT numbers. These are not exhaustive or + closed — any other MT/reaction number still works via ``Reaction(n)`` + directly; the class attributes are just a convenience. The full list is + intentionally large (500+ constants, covering every officially-assigned + ENDF-6 MT), so use your editor's autocomplete or search rather than + scanning the alphabetical attribute list below. Roughly, the categories + are: + + * **Everyday physics** — ``TOTAL``, ``ELASTIC``, ``INELASTIC_SCATTER``, + ``N_2N``, ``N_3N``, ``FISSION``, ``CAPTURE``, ``N_P``, ``N_D``, ``N_T``, + ``N_HE3``, ``N_ALPHA`` — the handful of reactions most tallies actually + use. + * **MCNP's own aliases** (negative numbers: ``TOTAL_MCNP``, + ``ABSORPTION``, ``ELASTIC_MCNP``, ``HEATING``, ``PHOTON_PRODUCTION``, + ``FISSION_MCNP``) — computed directly from transport data, not a + single ENDF MT channel. + * **NJOY-derived quantities** — ``RADIATION_DAMAGE`` and its + ``RADIATION_DAMAGE_*`` siblings (HEATR displacement-damage energies, + split the same way ENDF splits total/elastic/inelastic/capture), plus + ``AVERAGE_LETHARGY``, ``INVERSE_VELOCITY``, ``WEIGHTING_FLUX``, + ``PHOTON_HEATING``, ``KINEMATIC_KERMA``, + ``FISSION_STEADY_STATE_SPECTRUM``, and ``FISSION_DELAYED_SPECTRUM`` + (from GROUPR/HEATR/DTFR) — not official ENDF-6 reaction numbers, but + real values NJOY computes and MCNP libraries carry. + * **Everything else in ENDF-6 Appendix B** — redundant summary/total + cross sections (``NONELASTIC``, ``TOTAL_ABSORPTION``, ...); partial + fission chances (``FISSION_FIRST_CHANCE`` .. ``FISSION_FOURTH_CHANCE``); + exclusive multi-particle-emission channels (``N_2N_D``, ``N_N_ALPHA``, + ... roughly MT 11-200); total-particle-production sums + (``TOTAL_NEUTRON_PRODUCTION``, ...); fission nu-bar/yield/decay data + (``NU_TOTAL``, ``NU_DELAYED``, ...); photo-/electro-atomic data for + incident photons/electrons (``TOTAL_ATOMIC_INTERACTION``, + ``PHOTON_COHERENT_SCATTERING``, atomic subshells ``SUBSHELL_K`` .. + ``SUBSHELL_Q3``, ...); and discrete-level exit channels (``_L00`` .. + ``_L48``, plus a trailing ``_CONTINUUM``) for ``INELASTIC_SCATTER``, + ``N_P``, ``N_D``, ``N_T``, ``N_HE3``, ``N_ALPHA``, and ``N_2N``. + + Non-standard, library-specific MT numbers that aren't covered by any + named constant — e.g. IRDF-II's non-standard capture reaction, MT 11102 + — still work via plain ``Reaction(n)``: MontePy doesn't validate MT + numbers against any of these lists. + + .. versionadded:: 1.6.0b2 + """ + + def __init__(self, number: int): + self._number = number + self._left = None + self._operator = None + self._right = None + self._node = None + + @property + def number(self) -> int: + """The raw ENDF (MT) or special (R) reaction number.""" + return self._number + + def _ensure_has_node(self): + if self._node is None: + node = syntax_node.ListNode("reaction expr") + node.append(_make_value_node(int, self._number, padding=" ")) + self._node = node + + def __eq__(self, other): + if not isinstance(other, Reaction): + return NotImplemented + return self._number == other._number + + def __repr__(self): + return f"Reaction({self._number})" + + +# Common reaction numbers as ready-to-use Reaction instances, attached here +# rather than in the class body above, since `Reaction` isn't bound as a +# name until the class statement finishes executing. + +# ENDF MT numbers (positive; direct ENDF cross-section channel) +Reaction.TOTAL = Reaction(1) +Reaction.ELASTIC = Reaction(2) +Reaction.INELASTIC_SCATTER = Reaction(4) +Reaction.N_2N = Reaction(16) +Reaction.N_3N = Reaction(17) +Reaction.FISSION = Reaction(18) +Reaction.CAPTURE = Reaction(102) +Reaction.N_P = Reaction(103) +Reaction.N_D = Reaction(104) +Reaction.N_T = Reaction(105) +Reaction.N_HE3 = Reaction(106) +Reaction.N_ALPHA = Reaction(107) + +# NJOY HEATR radiation-damage-energy family. +Reaction.RADIATION_DAMAGE = Reaction(444) +Reaction.RADIATION_DAMAGE_ELASTIC = Reaction(445) +Reaction.RADIATION_DAMAGE_INELASTIC = Reaction(446) +Reaction.RADIATION_DAMAGE_DISAPPEARANCE = Reaction(447) + +# Other NJOY-module-specific "MT" identifiers that aren't official ENDF-6 +# reaction numbers (transcribed from the NJOY2016 manual, LA-UR-17-20093): +# GROUPR's special mtd values for slowing-down-moment group constants, +# HEATR's mtk values for heating diagnostics, and DTFR's special edit MTs +# for its DTF-format multigroup edits. +Reaction.AVERAGE_LETHARGY = Reaction(258) +Reaction.INVERSE_VELOCITY = Reaction(259) +Reaction.WEIGHTING_FLUX = Reaction(300) +Reaction.PHOTON_HEATING = Reaction(442) +Reaction.KINEMATIC_KERMA = Reaction(443) +Reaction.FISSION_STEADY_STATE_SPECTRUM = Reaction(470) +Reaction.FISSION_DELAYED_SPECTRUM = Reaction(471) + +# MCNP's own special reaction-number aliases (negative; computed directly +# from transport data, not a single ENDF MT channel). +Reaction.TOTAL_MCNP = Reaction(-1) +Reaction.ABSORPTION = Reaction(-2) +Reaction.ELASTIC_MCNP = Reaction(-3) +Reaction.HEATING = Reaction(-4) +Reaction.PHOTON_PRODUCTION = Reaction(-5) +Reaction.FISSION_MCNP = Reaction(-6) + +# The remaining constants below cover every other officially-assigned +# ENDF-6 MT number from Appendix B of the ENDF-6 Formats Manual, transcribed +# directly from that appendix. Skipped: MT numbers Appendix B marks +# "(Unassigned)" or "Not allowed in Version 6"; MT 6-9, 26, 31, 39, 40, 46-49, +# 120, 465-466 (old Version-5-only assignments); MT 301-450 (a formulaic +# "MT=300+reaction" energy-release/KERMA transform over the other reaction +# MTs, not itself a set of individually-assigned reactions); MT 451 (File 1 +# heading/title metadata, not a cross section); and MT 851-870 ("Lumped +# reaction covariances", a covariance grouping rather than a reaction). + +# Redundant summary / total-type cross sections. +Reaction.NONELASTIC = Reaction(3) +Reaction.ANYTHING = Reaction(5) +Reaction.TOTAL_CONTINUUM = Reaction(10) +Reaction.TOTAL_ABSORPTION = Reaction(27) + +# Partial (chance) fission cross sections; sum to Reaction.FISSION. +Reaction.FISSION_FIRST_CHANCE = Reaction(19) +Reaction.FISSION_SECOND_CHANCE = Reaction(20) +Reaction.FISSION_THIRD_CHANCE = Reaction(21) +Reaction.FISSION_FOURTH_CHANCE = Reaction(38) + +# Exclusive multi-particle-emission channels (single/few discrete exit +# channels, as opposed to the summed families like N_P above). +Reaction.N_2N_D = Reaction(11) +Reaction.N_N_ALPHA = Reaction(22) +Reaction.N_N_3ALPHA = Reaction(23) +Reaction.N_2N_ALPHA = Reaction(24) +Reaction.N_3N_ALPHA = Reaction(25) +Reaction.N_N_P = Reaction(28) +Reaction.N_N_2ALPHA = Reaction(29) +Reaction.N_2N_2ALPHA = Reaction(30) +Reaction.N_N_D = Reaction(32) +Reaction.N_N_T = Reaction(33) +Reaction.N_N_HE3 = Reaction(34) +Reaction.N_N_D_2ALPHA = Reaction(35) +Reaction.N_N_T_2ALPHA = Reaction(36) +Reaction.N_4N = Reaction(37) +Reaction.N_2N_P = Reaction(41) +Reaction.N_3N_P = Reaction(42) +Reaction.N_N_2P = Reaction(44) +Reaction.N_N_P_ALPHA = Reaction(45) +Reaction.N_2ALPHA = Reaction(108) +Reaction.N_3ALPHA = Reaction(109) +Reaction.N_2P = Reaction(111) +Reaction.N_P_ALPHA = Reaction(112) +Reaction.N_T_2ALPHA = Reaction(113) +Reaction.N_D_2ALPHA = Reaction(114) +Reaction.N_P_D = Reaction(115) +Reaction.N_P_T = Reaction(116) +Reaction.N_D_ALPHA = Reaction(117) + +# Neutron disappearance (capture-like absorption, excludes fission). +Reaction.NEUTRON_DISAPPEARANCE = Reaction(101) + +# Resonance-parameter data (File 2); incident neutrons only. +Reaction.RESONANCE_PARAMETERS = Reaction(151) + +# High-energy multi-particle-emission open channels, allocated to cover +# all reaction channels (within +/-10 mb) up to 60 MeV incident energy. +Reaction.N_5N = Reaction(152) +Reaction.N_6N = Reaction(153) +Reaction.N_2N_T = Reaction(154) +Reaction.N_T_ALPHA = Reaction(155) +Reaction.N_4N_P = Reaction(156) +Reaction.N_3N_D = Reaction(157) +Reaction.N_N_D_ALPHA = Reaction(158) +Reaction.N_2N_P_ALPHA = Reaction(159) +Reaction.N_7N = Reaction(160) +Reaction.N_8N = Reaction(161) +Reaction.N_5N_P = Reaction(162) +Reaction.N_6N_P = Reaction(163) +Reaction.N_7N_P = Reaction(164) +Reaction.N_4N_ALPHA = Reaction(165) +Reaction.N_5N_ALPHA = Reaction(166) +Reaction.N_6N_ALPHA = Reaction(167) +Reaction.N_7N_ALPHA = Reaction(168) +Reaction.N_4N_D = Reaction(169) +Reaction.N_5N_D = Reaction(170) +Reaction.N_6N_D = Reaction(171) +Reaction.N_3N_T = Reaction(172) +Reaction.N_4N_T = Reaction(173) +Reaction.N_5N_T = Reaction(174) +Reaction.N_6N_T = Reaction(175) +Reaction.N_2N_HE3 = Reaction(176) +Reaction.N_3N_HE3 = Reaction(177) +Reaction.N_4N_HE3 = Reaction(178) +Reaction.N_3N_2P = Reaction(179) +Reaction.N_3N_2ALPHA = Reaction(180) +Reaction.N_3N_P_ALPHA = Reaction(181) +Reaction.N_D_T = Reaction(182) +Reaction.N_N_P_D = Reaction(183) +Reaction.N_N_P_T = Reaction(184) +Reaction.N_N_D_T = Reaction(185) +Reaction.N_N_P_HE3 = Reaction(186) +Reaction.N_N_D_HE3 = Reaction(187) +Reaction.N_N_T_HE3 = Reaction(188) +Reaction.N_N_T_ALPHA = Reaction(189) +Reaction.N_2N_2P = Reaction(190) +Reaction.N_P_HE3 = Reaction(191) +Reaction.N_D_HE3 = Reaction(192) +Reaction.N_HE3_ALPHA = Reaction(193) +Reaction.N_4N_2P = Reaction(194) +Reaction.N_4N_2ALPHA = Reaction(195) +Reaction.N_4N_P_ALPHA = Reaction(196) +Reaction.N_3P = Reaction(197) +Reaction.N_N_3P = Reaction(198) +Reaction.N_3N_2P_ALPHA = Reaction(199) +Reaction.N_5N_2P = Reaction(200) + +# Redundant total-particle-production cross sections (derived files). +Reaction.TOTAL_NEUTRON_PRODUCTION = Reaction(201) +Reaction.TOTAL_GAMMA_PRODUCTION = Reaction(202) +Reaction.TOTAL_PROTON_PRODUCTION = Reaction(203) +Reaction.TOTAL_DEUTERON_PRODUCTION = Reaction(204) +Reaction.TOTAL_TRITON_PRODUCTION = Reaction(205) +Reaction.TOTAL_HE3_PRODUCTION = Reaction(206) +Reaction.TOTAL_ALPHA_PRODUCTION = Reaction(207) +Reaction.TOTAL_PI_PLUS_PRODUCTION = Reaction(208) +Reaction.TOTAL_PI_ZERO_PRODUCTION = Reaction(209) +Reaction.TOTAL_PI_MINUS_PRODUCTION = Reaction(210) +Reaction.TOTAL_MU_PLUS_PRODUCTION = Reaction(211) +Reaction.TOTAL_MU_MINUS_PRODUCTION = Reaction(212) +Reaction.TOTAL_KAON_PLUS_PRODUCTION = Reaction(213) +Reaction.TOTAL_KAON_ZERO_LONG_PRODUCTION = Reaction(214) +Reaction.TOTAL_KAON_ZERO_SHORT_PRODUCTION = Reaction(215) +Reaction.TOTAL_KAON_MINUS_PRODUCTION = Reaction(216) +Reaction.TOTAL_ANTIPROTON_PRODUCTION = Reaction(217) +Reaction.TOTAL_ANTINEUTRON_PRODUCTION = Reaction(218) + +# Elastic-scattering slowing-down moments (derived files only). +Reaction.AVERAGE_COSINE_ELASTIC = Reaction(251) +Reaction.AVERAGE_LOG_ENERGY_DECREMENT_ELASTIC = Reaction(252) +Reaction.AVERAGE_ENERGY_DECREMENT_RATIO_ELASTIC = Reaction(253) + +# Fission nu-bar, yield, and decay data. +Reaction.NU_TOTAL = Reaction(452) +Reaction.FISSION_YIELD_INDEPENDENT = Reaction(454) +Reaction.NU_DELAYED = Reaction(455) +Reaction.NU_PROMPT = Reaction(456) +Reaction.RADIOACTIVE_DECAY_DATA = Reaction(457) +Reaction.FISSION_ENERGY_RELEASE = Reaction(458) +Reaction.FISSION_YIELD_CUMULATIVE = Reaction(459) +Reaction.DELAYED_FISSION_PHOTONS = Reaction(460) + +# Photo-/electro-atomic interaction data (incident photons/electrons +# only, never incident neutrons; kept plainly-named for that reason). +Reaction.TOTAL_CHARGED_PARTICLE_STOPPING_POWER = Reaction(500) +Reaction.TOTAL_ATOMIC_INTERACTION = Reaction(501) +Reaction.PHOTON_COHERENT_SCATTERING = Reaction(502) +Reaction.PHOTON_INCOHERENT_SCATTERING = Reaction(504) +Reaction.IMAGINARY_SCATTERING_FACTOR = Reaction(505) +Reaction.REAL_SCATTERING_FACTOR = Reaction(506) +Reaction.PAIR_PRODUCTION_ELECTRON_FIELD = Reaction(515) +Reaction.PAIR_PRODUCTION_TOTAL = Reaction(516) +Reaction.PAIR_PRODUCTION_NUCLEAR_FIELD = Reaction(517) +Reaction.IONIZATION_TOTAL = Reaction(522) +Reaction.PHOTOEXCITATION = Reaction(523) +Reaction.LARGE_ANGLE_SCATTERING = Reaction(525) +Reaction.TOTAL_ELECTRO_ATOMIC_SCATTERING = Reaction(526) +Reaction.ELECTRO_ATOMIC_BREMSSTRAHLUNG = Reaction(527) +Reaction.ELECTRO_ATOMIC_EXCITATION = Reaction(528) +Reaction.ATOMIC_RELAXATION_DATA = Reaction(533) + +# Atomic-subshell photoelectric/electro-atomic cross sections. +Reaction.SUBSHELL_K = Reaction(534) +Reaction.SUBSHELL_L1 = Reaction(535) +Reaction.SUBSHELL_L2 = Reaction(536) +Reaction.SUBSHELL_L3 = Reaction(537) +Reaction.SUBSHELL_M1 = Reaction(538) +Reaction.SUBSHELL_M2 = Reaction(539) +Reaction.SUBSHELL_M3 = Reaction(540) +Reaction.SUBSHELL_M4 = Reaction(541) +Reaction.SUBSHELL_M5 = Reaction(542) +Reaction.SUBSHELL_N1 = Reaction(543) +Reaction.SUBSHELL_N2 = Reaction(544) +Reaction.SUBSHELL_N3 = Reaction(545) +Reaction.SUBSHELL_N4 = Reaction(546) +Reaction.SUBSHELL_N5 = Reaction(547) +Reaction.SUBSHELL_N6 = Reaction(548) +Reaction.SUBSHELL_N7 = Reaction(549) +Reaction.SUBSHELL_O1 = Reaction(550) +Reaction.SUBSHELL_O2 = Reaction(551) +Reaction.SUBSHELL_O3 = Reaction(552) +Reaction.SUBSHELL_O4 = Reaction(553) +Reaction.SUBSHELL_O5 = Reaction(554) +Reaction.SUBSHELL_O6 = Reaction(555) +Reaction.SUBSHELL_O7 = Reaction(556) +Reaction.SUBSHELL_O8 = Reaction(557) +Reaction.SUBSHELL_O9 = Reaction(558) +Reaction.SUBSHELL_P1 = Reaction(559) +Reaction.SUBSHELL_P2 = Reaction(560) +Reaction.SUBSHELL_P3 = Reaction(561) +Reaction.SUBSHELL_P4 = Reaction(562) +Reaction.SUBSHELL_P5 = Reaction(563) +Reaction.SUBSHELL_P6 = Reaction(564) +Reaction.SUBSHELL_P7 = Reaction(565) +Reaction.SUBSHELL_P8 = Reaction(566) +Reaction.SUBSHELL_P9 = Reaction(567) +Reaction.SUBSHELL_P10 = Reaction(568) +Reaction.SUBSHELL_P11 = Reaction(569) +Reaction.SUBSHELL_Q1 = Reaction(570) +Reaction.SUBSHELL_Q2 = Reaction(571) +Reaction.SUBSHELL_Q3 = Reaction(572) + +# (n,n') to individual discrete excited levels of the residual nucleus, +# and the continuum remainder; sum to Reaction.INELASTIC_SCATTER (MT 4). +# No L00/MT 50 (ground state) for incident neutrons -- that's elastic +# scattering, Reaction.ELASTIC. +Reaction.INELASTIC_SCATTER_L01 = Reaction(51) +Reaction.INELASTIC_SCATTER_L02 = Reaction(52) +Reaction.INELASTIC_SCATTER_L03 = Reaction(53) +Reaction.INELASTIC_SCATTER_L04 = Reaction(54) +Reaction.INELASTIC_SCATTER_L05 = Reaction(55) +Reaction.INELASTIC_SCATTER_L06 = Reaction(56) +Reaction.INELASTIC_SCATTER_L07 = Reaction(57) +Reaction.INELASTIC_SCATTER_L08 = Reaction(58) +Reaction.INELASTIC_SCATTER_L09 = Reaction(59) +Reaction.INELASTIC_SCATTER_L10 = Reaction(60) +Reaction.INELASTIC_SCATTER_L11 = Reaction(61) +Reaction.INELASTIC_SCATTER_L12 = Reaction(62) +Reaction.INELASTIC_SCATTER_L13 = Reaction(63) +Reaction.INELASTIC_SCATTER_L14 = Reaction(64) +Reaction.INELASTIC_SCATTER_L15 = Reaction(65) +Reaction.INELASTIC_SCATTER_L16 = Reaction(66) +Reaction.INELASTIC_SCATTER_L17 = Reaction(67) +Reaction.INELASTIC_SCATTER_L18 = Reaction(68) +Reaction.INELASTIC_SCATTER_L19 = Reaction(69) +Reaction.INELASTIC_SCATTER_L20 = Reaction(70) +Reaction.INELASTIC_SCATTER_L21 = Reaction(71) +Reaction.INELASTIC_SCATTER_L22 = Reaction(72) +Reaction.INELASTIC_SCATTER_L23 = Reaction(73) +Reaction.INELASTIC_SCATTER_L24 = Reaction(74) +Reaction.INELASTIC_SCATTER_L25 = Reaction(75) +Reaction.INELASTIC_SCATTER_L26 = Reaction(76) +Reaction.INELASTIC_SCATTER_L27 = Reaction(77) +Reaction.INELASTIC_SCATTER_L28 = Reaction(78) +Reaction.INELASTIC_SCATTER_L29 = Reaction(79) +Reaction.INELASTIC_SCATTER_L30 = Reaction(80) +Reaction.INELASTIC_SCATTER_L31 = Reaction(81) +Reaction.INELASTIC_SCATTER_L32 = Reaction(82) +Reaction.INELASTIC_SCATTER_L33 = Reaction(83) +Reaction.INELASTIC_SCATTER_L34 = Reaction(84) +Reaction.INELASTIC_SCATTER_L35 = Reaction(85) +Reaction.INELASTIC_SCATTER_L36 = Reaction(86) +Reaction.INELASTIC_SCATTER_L37 = Reaction(87) +Reaction.INELASTIC_SCATTER_L38 = Reaction(88) +Reaction.INELASTIC_SCATTER_L39 = Reaction(89) +Reaction.INELASTIC_SCATTER_L40 = Reaction(90) +Reaction.INELASTIC_SCATTER_CONTINUUM = Reaction(91) + +# (n,p) to individual discrete levels (L00 = ground state) and the +# continuum remainder; sum to Reaction.N_P (MT 103). +Reaction.N_P_L00 = Reaction(600) +Reaction.N_P_L01 = Reaction(601) +Reaction.N_P_L02 = Reaction(602) +Reaction.N_P_L03 = Reaction(603) +Reaction.N_P_L04 = Reaction(604) +Reaction.N_P_L05 = Reaction(605) +Reaction.N_P_L06 = Reaction(606) +Reaction.N_P_L07 = Reaction(607) +Reaction.N_P_L08 = Reaction(608) +Reaction.N_P_L09 = Reaction(609) +Reaction.N_P_L10 = Reaction(610) +Reaction.N_P_L11 = Reaction(611) +Reaction.N_P_L12 = Reaction(612) +Reaction.N_P_L13 = Reaction(613) +Reaction.N_P_L14 = Reaction(614) +Reaction.N_P_L15 = Reaction(615) +Reaction.N_P_L16 = Reaction(616) +Reaction.N_P_L17 = Reaction(617) +Reaction.N_P_L18 = Reaction(618) +Reaction.N_P_L19 = Reaction(619) +Reaction.N_P_L20 = Reaction(620) +Reaction.N_P_L21 = Reaction(621) +Reaction.N_P_L22 = Reaction(622) +Reaction.N_P_L23 = Reaction(623) +Reaction.N_P_L24 = Reaction(624) +Reaction.N_P_L25 = Reaction(625) +Reaction.N_P_L26 = Reaction(626) +Reaction.N_P_L27 = Reaction(627) +Reaction.N_P_L28 = Reaction(628) +Reaction.N_P_L29 = Reaction(629) +Reaction.N_P_L30 = Reaction(630) +Reaction.N_P_L31 = Reaction(631) +Reaction.N_P_L32 = Reaction(632) +Reaction.N_P_L33 = Reaction(633) +Reaction.N_P_L34 = Reaction(634) +Reaction.N_P_L35 = Reaction(635) +Reaction.N_P_L36 = Reaction(636) +Reaction.N_P_L37 = Reaction(637) +Reaction.N_P_L38 = Reaction(638) +Reaction.N_P_L39 = Reaction(639) +Reaction.N_P_L40 = Reaction(640) +Reaction.N_P_L41 = Reaction(641) +Reaction.N_P_L42 = Reaction(642) +Reaction.N_P_L43 = Reaction(643) +Reaction.N_P_L44 = Reaction(644) +Reaction.N_P_L45 = Reaction(645) +Reaction.N_P_L46 = Reaction(646) +Reaction.N_P_L47 = Reaction(647) +Reaction.N_P_L48 = Reaction(648) +Reaction.N_P_CONTINUUM = Reaction(649) + +# (n,d) to individual discrete levels (L00 = ground state) and the +# continuum remainder; sum to Reaction.N_D (MT 104). +Reaction.N_D_L00 = Reaction(650) +Reaction.N_D_L01 = Reaction(651) +Reaction.N_D_L02 = Reaction(652) +Reaction.N_D_L03 = Reaction(653) +Reaction.N_D_L04 = Reaction(654) +Reaction.N_D_L05 = Reaction(655) +Reaction.N_D_L06 = Reaction(656) +Reaction.N_D_L07 = Reaction(657) +Reaction.N_D_L08 = Reaction(658) +Reaction.N_D_L09 = Reaction(659) +Reaction.N_D_L10 = Reaction(660) +Reaction.N_D_L11 = Reaction(661) +Reaction.N_D_L12 = Reaction(662) +Reaction.N_D_L13 = Reaction(663) +Reaction.N_D_L14 = Reaction(664) +Reaction.N_D_L15 = Reaction(665) +Reaction.N_D_L16 = Reaction(666) +Reaction.N_D_L17 = Reaction(667) +Reaction.N_D_L18 = Reaction(668) +Reaction.N_D_L19 = Reaction(669) +Reaction.N_D_L20 = Reaction(670) +Reaction.N_D_L21 = Reaction(671) +Reaction.N_D_L22 = Reaction(672) +Reaction.N_D_L23 = Reaction(673) +Reaction.N_D_L24 = Reaction(674) +Reaction.N_D_L25 = Reaction(675) +Reaction.N_D_L26 = Reaction(676) +Reaction.N_D_L27 = Reaction(677) +Reaction.N_D_L28 = Reaction(678) +Reaction.N_D_L29 = Reaction(679) +Reaction.N_D_L30 = Reaction(680) +Reaction.N_D_L31 = Reaction(681) +Reaction.N_D_L32 = Reaction(682) +Reaction.N_D_L33 = Reaction(683) +Reaction.N_D_L34 = Reaction(684) +Reaction.N_D_L35 = Reaction(685) +Reaction.N_D_L36 = Reaction(686) +Reaction.N_D_L37 = Reaction(687) +Reaction.N_D_L38 = Reaction(688) +Reaction.N_D_L39 = Reaction(689) +Reaction.N_D_L40 = Reaction(690) +Reaction.N_D_L41 = Reaction(691) +Reaction.N_D_L42 = Reaction(692) +Reaction.N_D_L43 = Reaction(693) +Reaction.N_D_L44 = Reaction(694) +Reaction.N_D_L45 = Reaction(695) +Reaction.N_D_L46 = Reaction(696) +Reaction.N_D_L47 = Reaction(697) +Reaction.N_D_L48 = Reaction(698) +Reaction.N_D_CONTINUUM = Reaction(699) + +# (n,t) to individual discrete levels (L00 = ground state) and the +# continuum remainder; sum to Reaction.N_T (MT 105). +Reaction.N_T_L00 = Reaction(700) +Reaction.N_T_L01 = Reaction(701) +Reaction.N_T_L02 = Reaction(702) +Reaction.N_T_L03 = Reaction(703) +Reaction.N_T_L04 = Reaction(704) +Reaction.N_T_L05 = Reaction(705) +Reaction.N_T_L06 = Reaction(706) +Reaction.N_T_L07 = Reaction(707) +Reaction.N_T_L08 = Reaction(708) +Reaction.N_T_L09 = Reaction(709) +Reaction.N_T_L10 = Reaction(710) +Reaction.N_T_L11 = Reaction(711) +Reaction.N_T_L12 = Reaction(712) +Reaction.N_T_L13 = Reaction(713) +Reaction.N_T_L14 = Reaction(714) +Reaction.N_T_L15 = Reaction(715) +Reaction.N_T_L16 = Reaction(716) +Reaction.N_T_L17 = Reaction(717) +Reaction.N_T_L18 = Reaction(718) +Reaction.N_T_L19 = Reaction(719) +Reaction.N_T_L20 = Reaction(720) +Reaction.N_T_L21 = Reaction(721) +Reaction.N_T_L22 = Reaction(722) +Reaction.N_T_L23 = Reaction(723) +Reaction.N_T_L24 = Reaction(724) +Reaction.N_T_L25 = Reaction(725) +Reaction.N_T_L26 = Reaction(726) +Reaction.N_T_L27 = Reaction(727) +Reaction.N_T_L28 = Reaction(728) +Reaction.N_T_L29 = Reaction(729) +Reaction.N_T_L30 = Reaction(730) +Reaction.N_T_L31 = Reaction(731) +Reaction.N_T_L32 = Reaction(732) +Reaction.N_T_L33 = Reaction(733) +Reaction.N_T_L34 = Reaction(734) +Reaction.N_T_L35 = Reaction(735) +Reaction.N_T_L36 = Reaction(736) +Reaction.N_T_L37 = Reaction(737) +Reaction.N_T_L38 = Reaction(738) +Reaction.N_T_L39 = Reaction(739) +Reaction.N_T_L40 = Reaction(740) +Reaction.N_T_L41 = Reaction(741) +Reaction.N_T_L42 = Reaction(742) +Reaction.N_T_L43 = Reaction(743) +Reaction.N_T_L44 = Reaction(744) +Reaction.N_T_L45 = Reaction(745) +Reaction.N_T_L46 = Reaction(746) +Reaction.N_T_L47 = Reaction(747) +Reaction.N_T_L48 = Reaction(748) +Reaction.N_T_CONTINUUM = Reaction(749) + +# (n,He3) to individual discrete levels (L00 = ground state) and the +# continuum remainder; sum to Reaction.N_HE3 (MT 106). +Reaction.N_HE3_L00 = Reaction(750) +Reaction.N_HE3_L01 = Reaction(751) +Reaction.N_HE3_L02 = Reaction(752) +Reaction.N_HE3_L03 = Reaction(753) +Reaction.N_HE3_L04 = Reaction(754) +Reaction.N_HE3_L05 = Reaction(755) +Reaction.N_HE3_L06 = Reaction(756) +Reaction.N_HE3_L07 = Reaction(757) +Reaction.N_HE3_L08 = Reaction(758) +Reaction.N_HE3_L09 = Reaction(759) +Reaction.N_HE3_L10 = Reaction(760) +Reaction.N_HE3_L11 = Reaction(761) +Reaction.N_HE3_L12 = Reaction(762) +Reaction.N_HE3_L13 = Reaction(763) +Reaction.N_HE3_L14 = Reaction(764) +Reaction.N_HE3_L15 = Reaction(765) +Reaction.N_HE3_L16 = Reaction(766) +Reaction.N_HE3_L17 = Reaction(767) +Reaction.N_HE3_L18 = Reaction(768) +Reaction.N_HE3_L19 = Reaction(769) +Reaction.N_HE3_L20 = Reaction(770) +Reaction.N_HE3_L21 = Reaction(771) +Reaction.N_HE3_L22 = Reaction(772) +Reaction.N_HE3_L23 = Reaction(773) +Reaction.N_HE3_L24 = Reaction(774) +Reaction.N_HE3_L25 = Reaction(775) +Reaction.N_HE3_L26 = Reaction(776) +Reaction.N_HE3_L27 = Reaction(777) +Reaction.N_HE3_L28 = Reaction(778) +Reaction.N_HE3_L29 = Reaction(779) +Reaction.N_HE3_L30 = Reaction(780) +Reaction.N_HE3_L31 = Reaction(781) +Reaction.N_HE3_L32 = Reaction(782) +Reaction.N_HE3_L33 = Reaction(783) +Reaction.N_HE3_L34 = Reaction(784) +Reaction.N_HE3_L35 = Reaction(785) +Reaction.N_HE3_L36 = Reaction(786) +Reaction.N_HE3_L37 = Reaction(787) +Reaction.N_HE3_L38 = Reaction(788) +Reaction.N_HE3_L39 = Reaction(789) +Reaction.N_HE3_L40 = Reaction(790) +Reaction.N_HE3_L41 = Reaction(791) +Reaction.N_HE3_L42 = Reaction(792) +Reaction.N_HE3_L43 = Reaction(793) +Reaction.N_HE3_L44 = Reaction(794) +Reaction.N_HE3_L45 = Reaction(795) +Reaction.N_HE3_L46 = Reaction(796) +Reaction.N_HE3_L47 = Reaction(797) +Reaction.N_HE3_L48 = Reaction(798) +Reaction.N_HE3_CONTINUUM = Reaction(799) + +# (n,alpha) to individual discrete levels (L00 = ground state) and the +# continuum remainder; sum to Reaction.N_ALPHA (MT 107). +Reaction.N_ALPHA_L00 = Reaction(800) +Reaction.N_ALPHA_L01 = Reaction(801) +Reaction.N_ALPHA_L02 = Reaction(802) +Reaction.N_ALPHA_L03 = Reaction(803) +Reaction.N_ALPHA_L04 = Reaction(804) +Reaction.N_ALPHA_L05 = Reaction(805) +Reaction.N_ALPHA_L06 = Reaction(806) +Reaction.N_ALPHA_L07 = Reaction(807) +Reaction.N_ALPHA_L08 = Reaction(808) +Reaction.N_ALPHA_L09 = Reaction(809) +Reaction.N_ALPHA_L10 = Reaction(810) +Reaction.N_ALPHA_L11 = Reaction(811) +Reaction.N_ALPHA_L12 = Reaction(812) +Reaction.N_ALPHA_L13 = Reaction(813) +Reaction.N_ALPHA_L14 = Reaction(814) +Reaction.N_ALPHA_L15 = Reaction(815) +Reaction.N_ALPHA_L16 = Reaction(816) +Reaction.N_ALPHA_L17 = Reaction(817) +Reaction.N_ALPHA_L18 = Reaction(818) +Reaction.N_ALPHA_L19 = Reaction(819) +Reaction.N_ALPHA_L20 = Reaction(820) +Reaction.N_ALPHA_L21 = Reaction(821) +Reaction.N_ALPHA_L22 = Reaction(822) +Reaction.N_ALPHA_L23 = Reaction(823) +Reaction.N_ALPHA_L24 = Reaction(824) +Reaction.N_ALPHA_L25 = Reaction(825) +Reaction.N_ALPHA_L26 = Reaction(826) +Reaction.N_ALPHA_L27 = Reaction(827) +Reaction.N_ALPHA_L28 = Reaction(828) +Reaction.N_ALPHA_L29 = Reaction(829) +Reaction.N_ALPHA_L30 = Reaction(830) +Reaction.N_ALPHA_L31 = Reaction(831) +Reaction.N_ALPHA_L32 = Reaction(832) +Reaction.N_ALPHA_L33 = Reaction(833) +Reaction.N_ALPHA_L34 = Reaction(834) +Reaction.N_ALPHA_L35 = Reaction(835) +Reaction.N_ALPHA_L36 = Reaction(836) +Reaction.N_ALPHA_L37 = Reaction(837) +Reaction.N_ALPHA_L38 = Reaction(838) +Reaction.N_ALPHA_L39 = Reaction(839) +Reaction.N_ALPHA_L40 = Reaction(840) +Reaction.N_ALPHA_L41 = Reaction(841) +Reaction.N_ALPHA_L42 = Reaction(842) +Reaction.N_ALPHA_L43 = Reaction(843) +Reaction.N_ALPHA_L44 = Reaction(844) +Reaction.N_ALPHA_L45 = Reaction(845) +Reaction.N_ALPHA_L46 = Reaction(846) +Reaction.N_ALPHA_L47 = Reaction(847) +Reaction.N_ALPHA_L48 = Reaction(848) +Reaction.N_ALPHA_CONTINUUM = Reaction(849) + +# (n,2n) to individual discrete levels (L00 = ground state) and the +# continuum remainder; sum to Reaction.N_2N (MT 16). +Reaction.N_2N_L00 = Reaction(875) +Reaction.N_2N_L01 = Reaction(876) +Reaction.N_2N_L02 = Reaction(877) +Reaction.N_2N_L03 = Reaction(878) +Reaction.N_2N_L04 = Reaction(879) +Reaction.N_2N_L05 = Reaction(880) +Reaction.N_2N_L06 = Reaction(881) +Reaction.N_2N_L07 = Reaction(882) +Reaction.N_2N_L08 = Reaction(883) +Reaction.N_2N_L09 = Reaction(884) +Reaction.N_2N_L10 = Reaction(885) +Reaction.N_2N_L11 = Reaction(886) +Reaction.N_2N_L12 = Reaction(887) +Reaction.N_2N_L13 = Reaction(888) +Reaction.N_2N_L14 = Reaction(889) +Reaction.N_2N_L15 = Reaction(890) +Reaction.N_2N_CONTINUUM = Reaction(891) + + +class AttenuatorLayer: + """One layer of an FM attenuator set: ``m px``. + + .. versionadded:: 1.6.0b2 + + Parameters + ---------- + material : int + Material number identified on an ``Mm`` card. + areal_density : float + Density times thickness of the attenuating layer. Always stored + positive (like :attr:`montepy.Material.is_atom_fraction`'s + convention) — the atom-vs-mass distinction is carried separately in + ``is_atom_density``, not via sign. + is_atom_density : bool + ``True`` if ``areal_density`` is an atom density (atoms/barn-cm), + ``False`` if it's a mass density (g/cm3). Corresponds to the sign + of the raw MCNP ``px`` value (positive = atom, negative = mass). + """ + + __slots__ = ("_material", "_areal_density", "_is_atom_density") + + @args_checked + def __init__( + self, + material: ty.Integral, + areal_density: ty.Real, + is_atom_density: bool = True, + ): + self._material = material + self._areal_density = abs(areal_density) + self._is_atom_density = is_atom_density + + @property + def material(self) -> int: + """The material number for this layer, from an ``Mm`` card.""" + return self._material + + @property + def areal_density(self) -> float: + """Density times thickness of this layer (always positive).""" + return self._areal_density + + @property + def is_atom_density(self) -> bool: + """``True`` if :attr:`areal_density` is an atom density, ``False`` if mass density.""" + return self._is_atom_density + + def __eq__(self, other): + if not isinstance(other, AttenuatorLayer): + return NotImplemented + return ( + self._material == other._material + and self._areal_density == other._areal_density + and self._is_atom_density == other._is_atom_density + ) + + def __repr__(self): + return f"AttenuatorLayer({self._material}, {self._areal_density}, is_atom_density={self._is_atom_density})" + + +class AttenuatorSet: + """An FM attenuator set: ``c -1 m1 px1 m2 px2 ...``. + + Models the thin-shield line-of-sight attenuation factor + ``exp(-sum(sigma_i * px_i))``. Layers chain via ``&`` (mirrors + :class:`~montepy.HalfSpace`'s own ``&``; layers stack multiplicatively + in the exponent, like an intersection of independent attenuating + conditions): + + .. code-block:: python + + attenuator = AttenuatorSet(1.0, [AttenuatorLayer(3, 0.05)]) + attenuator = attenuator & AttenuatorLayer(4, 0.1, is_atom_density=False) + + .. versionadded:: 1.6.0b2 + """ + + __slots__ = ("_constant", "_layers", "_node") + + @args_checked + def __init__(self, constant: ty.Real, layers: list[AttenuatorLayer]): + self._constant = constant + self._layers = list(layers) + self._node = None + + @property + def constant(self) -> float: + """The scalar constant ``c`` for this attenuator set.""" + return self._constant + + @property + def layers(self) -> list[AttenuatorLayer]: + """The attenuating layers, in order.""" + return list(self._layers) + + @property + def node(self) -> syntax_node.ListNode: + """The syntax node for this attenuator set. + + .. versionadded:: 1.6.0b3 + """ + self._ensure_has_node() + return self._node + + def _ensure_has_node(self): + if self._node is not None: + return + node = syntax_node.ListNode("attenuator set") + node.append(_make_value_node(float, self._constant)) + node.append(_make_value_node(int, -1)) + for layer in self._layers: + node.append(_make_value_node(int, layer.material)) + density = ( + layer.areal_density if layer.is_atom_density else -layer.areal_density + ) + node.append(_make_value_node(float, density)) + self._node = node + + @args_checked + def __and__(self, other: Union[AttenuatorLayer, "AttenuatorSet"]) -> AttenuatorSet: + """Return a new :class:`AttenuatorSet` with ``other``'s layer(s) appended.""" + if isinstance(other, AttenuatorLayer): + new_layers = [other] + else: + new_layers = other.layers + return AttenuatorSet(self._constant, self._layers + new_layers) + + def __eq__(self, other): + if not isinstance(other, AttenuatorSet): + return NotImplemented + return self._constant == other._constant and self._layers == other._layers + + def __repr__(self): + return f"AttenuatorSet({self._constant}, {self._layers!r})" + + +class MultiplierSet: + """An FM multiplier set: ``c m (reaction list 1) (reaction list 2) ...``. + + .. versionadded:: 1.6.0b2 + + Parameters + ---------- + constant : float + The scalar constant ``c``. If negative (type-4 tallies only), MCNP + replaces ``|c|`` with ``|c|`` times the tallying cell's atom density. + material : int, optional + Material number from an ``Mm`` card. ``None``/``0`` means "the + material of the current cell." + reactions : list[montepy.data_inputs.tally_multiplier.ReactionExpression] + One entry per output bin this set creates (MCNP creates one bin per + reaction list, per FM spec footnote 4). + """ + + __slots__ = ("_constant", "_material", "_reactions", "_node", "_material_node") + + @args_checked + def __init__( + self, + constant: ty.Real, + material: Union[ty.Integral, "montepy.Material", None], + reactions: list[ReactionExpression], + ): + self._constant = constant + self._material = material + self._reactions = list(reactions) + self._node = None + self._material_node = None + + @property + def constant(self) -> float: + """The scalar constant ``c`` for this multiplier set.""" + return self._constant + + @property + def material(self) -> int | None: + """The material number, or ``None`` for "current cell's material". + + Resolved live from the linked :class:`~montepy.Material` if this set + was built from one (e.g. ``mat & Reaction.CAPTURE``), so a later + renumber of that material is reflected here too. + """ + if isinstance(self._material, montepy.Material): + return self._material.number + return self._material + + @property + def reactions(self) -> list[ReactionExpression]: + """One :class:`~montepy.data_inputs.tally_multiplier.ReactionExpression` per output bin this set creates.""" + return list(self._reactions) + + @property + def node(self) -> syntax_node.ListNode: + """The syntax node for this multiplier set. + + .. versionadded:: 1.6.0b3 + """ + self._ensure_has_node() + if ( + self._material_node is not None + and self._material_node.value != self.material + ): + self._material_node.value = self.material + return self._node + + def _ensure_has_node(self): + if self._node is not None: + return + node = syntax_node.ListNode("multiplier set") + node.append(_make_value_node(float, self._constant)) + if self.material is not None: + self._material_node = _make_value_node(int, self.material) + node.append(self._material_node) + if len(self._reactions) == 1: + node.append(self._reactions[0].node) + else: + for reaction in self._reactions: + group = syntax_node.ListNode("tally group") + group.append(syntax_node.PaddingNode("(")) + for n in reaction.node.nodes: + group.append(n) + group.append(syntax_node.PaddingNode(")")) + node.append(group) + self._node = node + + def __rmul__(self, constant: ty.Real) -> MultiplierSet: + """``1.5 * (mat1 & Reaction.CAPTURE)`` sets the constant. + + Completes the DSL alongside ``ReactionExpression.__rand__``: + ``material & reaction`` builds a :class:`MultiplierSet` with + ``constant=1.0``, and this lets you scale it afterwards, mirroring + ``ReactionExpression.__rmul__``'s int-first convenience + (``16 * Reaction(103)``). + """ + return MultiplierSet(constant, self._material, self._reactions) + + def __eq__(self, other): + if not isinstance(other, MultiplierSet): + return NotImplemented + return ( + self._constant == other._constant + and self.material == other.material + and self._reactions == other._reactions + ) + + def __repr__(self): + return f"MultiplierSet({self._constant}, {self.material}, {self._reactions!r})" + + +class SpecialMultiplierSet: + """An FM special multiplier set: ``c k``. + + .. versionadded:: 1.6.0b2 + + Parameters + ---------- + constant : float + The scalar constant ``c``. + kind : SpecialMultiplier + Which special multiplier option (``k``) this is. + """ + + __slots__ = ("_constant", "_kind", "_node") + + @args_checked + def __init__(self, constant: ty.Real, kind: SpecialMultiplier): + self._constant = constant + self._kind = kind + self._node = None + + @property + def constant(self) -> float: + """The scalar constant ``c`` for this special multiplier set.""" + return self._constant + + @property + def kind(self) -> SpecialMultiplier: + """Which special multiplier option this is.""" + return self._kind + + @property + def node(self) -> syntax_node.ListNode: + """The syntax node for this special multiplier set. + + .. versionadded:: 1.6.0b3 + """ + self._ensure_has_node() + return self._node + + def _ensure_has_node(self): + if self._node is not None: + return + node = syntax_node.ListNode("special multiplier set") + node.append(_make_value_node(float, self._constant)) + node.append(_make_value_node(int, _SPECIAL_KIND_MAP_INVERSE[self._kind])) + self._node = node + + def __eq__(self, other): + if not isinstance(other, SpecialMultiplierSet): + return NotImplemented + return self._constant == other._constant and self._kind == other._kind + + def __repr__(self): + return f"SpecialMultiplierSet({self._constant}, {self._kind})" + + +class MultiplierScore: + """One physical output bin's full recipe: exactly one FM-derived tally score. + + Exactly one of ``reaction``/``kind`` is non-``None``: ``reaction`` for a + bin coming from a :class:`MultiplierSet`, ``kind`` for one coming from a + :class:`SpecialMultiplierSet`. + + .. versionadded:: 1.6.0b2 + """ + + __slots__ = ("_constant", "_material", "_reaction", "_kind", "_attenuator") + + def __init__( + self, + constant: ty.Real, + material: ty.Integral | None, + reaction: ReactionExpression | None, + kind: SpecialMultiplier | None, + attenuator: AttenuatorSet | None, + ): + self._constant = constant + self._material = material + self._reaction = reaction + self._kind = kind + self._attenuator = attenuator + + @property + def constant(self) -> float: + """The scalar constant ``c`` for this score.""" + return self._constant + + @property + def material(self) -> int | None: + """The material number, or ``None``.""" + return self._material + + @property + def reaction(self) -> ReactionExpression | None: + """The reaction expression for this score, if from a :class:`MultiplierSet`.""" + return self._reaction + + @property + def kind(self) -> SpecialMultiplier | None: + """The special multiplier kind for this score, if from a :class:`SpecialMultiplierSet`.""" + return self._kind + + @property + def attenuator(self) -> AttenuatorSet | None: + """The attenuator applied to this score, if any, inherited from the parent :class:`MultiplierBin`.""" + return self._attenuator + + def __eq__(self, other): + if not isinstance(other, MultiplierScore): + return NotImplemented + return ( + self._constant == other._constant + and self._material == other._material + and self._reaction == other._reaction + and self._kind == other._kind + and self._attenuator == other._attenuator + ) + + def __repr__(self): + return ( + f"MultiplierScore(constant={self._constant}, material={self._material}, " + f"reaction={self._reaction!r}, kind={self._kind}, attenuator={self._attenuator!r})" + ) + + +class MultiplierBin: + """One top-level ``(bin set k)`` group of an FM card. + + .. versionadded:: 1.6.0b2 + + Parameters + ---------- + terms : list[MultiplierSet | SpecialMultiplierSet] + The multiplier/special-multiplier sets in this bin set. + attenuator : AttenuatorSet, optional + The attenuator set for this bin set, if any -- applies to every bin + the ``terms`` produce. + """ + + __slots__ = ("_terms", "_attenuator", "_node") + + def __init__( + self, + terms: list[MultiplierSet | SpecialMultiplierSet], + attenuator: AttenuatorSet | None = None, + ): + self._terms = list(terms) + self._attenuator = attenuator + self._node = None + + @property + def terms(self) -> list[MultiplierSet | SpecialMultiplierSet]: + """The multiplier/special-multiplier sets in this bin set.""" + return list(self._terms) + + @property + def attenuator(self) -> AttenuatorSet | None: + """The attenuator set for this bin set, if any.""" + return self._attenuator + + @property + def node(self) -> syntax_node.ListNode: + """The syntax node for this bin's own content (not including the + outer parens a sibling bin or multi-term structure may require -- + that's decided by :class:`TallyMultiplier`, which knows about + sibling bins). + + .. versionadded:: 1.6.0b3 + """ + self._ensure_has_node() + items = list(self._terms) + if self._attenuator is not None: + items.append(self._attenuator) + for item in items: + if item._node is not None: + item.node + return self._node + + def _own_item_count(self): + return len(self._terms) + (1 if self._attenuator is not None else 0) + + def _ensure_has_node(self): + if self._node is not None: + return + items = list(self._terms) + if self._attenuator is not None: + items.append(self._attenuator) + node = syntax_node.ListNode("bin body") + wrap_each = len(items) > 1 + for item in items: + item_node = item.node + if wrap_each: + group = syntax_node.ListNode("tally group") + group.append(syntax_node.PaddingNode("(")) + for n in item_node.nodes: + group.append(n) + group.append(syntax_node.PaddingNode(")")) + node.append(group) + else: + for n in item_node.nodes: + node.append(n) + self._node = node + + @property + def scores(self) -> list[MultiplierScore]: + """Flatten this bin set's terms into one :class:`MultiplierScore` per actual output bin.""" + if not self._terms and self._attenuator is not None: + return [ + MultiplierScore( + self._attenuator.constant, None, None, None, self._attenuator + ) + ] + result = [] + for term in self._terms: + if isinstance(term, MultiplierSet): + if not term.reactions: + result.append( + MultiplierScore( + term.constant, term.material, None, None, self._attenuator + ) + ) + for reaction in term.reactions: + result.append( + MultiplierScore( + term.constant, + term.material, + reaction, + None, + self._attenuator, + ) + ) + elif isinstance(term, SpecialMultiplierSet): + result.append( + MultiplierScore( + term.constant, None, None, term.kind, self._attenuator + ) + ) + return result + + @classmethod + def from_items(cls, items: list) -> MultiplierBin: + """Parse a bin set's flat CST items into a :class:`MultiplierBin`.""" + original = list(items) + items = _non_padding(items) + nested_groups = [n for n in items if _is_group(n)] + if nested_groups and len(nested_groups) == len(items): + # Multiple sibling multiplier/attenuator sets, each individually + # parenthesized (FM parenthesization rule 2). + raw_terms = [_parse_term(_group_body(g)) for g in nested_groups] + else: + # A single term, inline (rule 1) -- items are that term's content. + raw_terms = [_parse_term(items)] + terms = [] + attenuator = None + for term in raw_terms: + if isinstance(term, AttenuatorSet): + attenuator = term + else: + terms.append(term) + result = cls(terms, attenuator) + node = syntax_node.ListNode("bin body") + for n in original: + node.append(n) + result._node = node + return result + + def __eq__(self, other): + if not isinstance(other, MultiplierBin): + return NotImplemented + return self._terms == other._terms and self._attenuator == other._attenuator + + def __repr__(self): + return f"MultiplierBin(terms={self._terms!r}, attenuator={self._attenuator!r})" + + +def _is_group(node) -> bool: + return isinstance(node, syntax_node.ListNode) and node.name == "tally group" + + +def _non_padding(nodes): + return [n for n in nodes if not isinstance(n, syntax_node.PaddingNode)] + + +def _group_body(group_node): + """Strip the leading/trailing paren ``ValueNode``s off a "tally group" ``ListNode``.""" + nodes = list(group_node.nodes) + return _non_padding(nodes[1:-1]) if len(nodes) >= 2 else _non_padding(nodes) + + +def _numeric_value(node): + return node.value + + +def _parse_term(items: list) -> MultiplierSet | SpecialMultiplierSet | AttenuatorSet: + """Parse one multiplier/special-multiplier/attenuator set's flat content.""" + items = _non_padding(items) + constant = _numeric_value(items[0]) + rest = items[1:] + if not rest: + return MultiplierSet(constant, None, []) + + first = rest[0] + first_val = ( + _numeric_value(first) if isinstance(first, syntax_node.ValueNode) else None + ) + + if first_val in _SPECIAL_KIND_MAP and len(rest) == 1: + return SpecialMultiplierSet(constant, _SPECIAL_KIND_MAP[first_val]) + + if first_val == -1 and len(rest) > 1: + layer_tokens = rest[1:] + layers = [] + for i in range(0, len(layer_tokens), 2): + mat = int(_numeric_value(layer_tokens[i])) + px = _numeric_value(layer_tokens[i + 1]) + layers.append(AttenuatorLayer(mat, abs(px), px >= 0)) + return AttenuatorSet(constant, layers) + + # A multiplier set: material, then reaction list(s). + material = int(first_val) + reaction_items = rest[1:] + nested = [n for n in reaction_items if _is_group(n)] + if nested: + reactions = [_parse_reaction_expr(_group_body(g)) for g in nested] + elif reaction_items: + reactions = [_parse_reaction_expr(reaction_items)] + else: + reactions = [] + return MultiplierSet(constant, material, reactions) + + +def _parse_reaction_expr(items: list) -> ReactionExpression: + """Resolve a flat reaction-list token stream into a :class:`ReactionExpression`. + + Applies "multiply first" precedence in Python (fold space-separated + numbers with ``*``, then fold ``:``/``#``-separated groups with + ``+``/``-``) rather than encoding it in the grammar -- mirrors how + ``tally.py``'s own ``_parse_tally_group_node`` does its real + interpretation as a second pass over a loosely structured CST. + """ + original = list(items) + items = _non_padding(items) + groups: list[tuple[ReactionOperator | None, list[Reaction]]] = [] + current_op = None + current_nums = [] + for item in items: + val = item.value + if isinstance(val, str) and val.strip() in (":", "#"): + groups.append((current_op, current_nums)) + current_op = ( + ReactionOperator.ADD + if val.strip() == ":" + else ReactionOperator.SUBTRACT + ) + current_nums = [] + else: + current_nums.append(Reaction(int(val))) + groups.append((current_op, current_nums)) + + def fold_multiply(nums): + expr = nums[0] + for n in nums[1:]: + expr = expr * n + return expr + + expr = fold_multiply(groups[0][1]) + for op, nums in groups[1:]: + term = fold_multiply(nums) + expr = expr + term if op == ReactionOperator.ADD else expr - term + node = syntax_node.ListNode("reaction expr") + for n in original: + node.append(n) + expr._node = node + return expr + + +def _mark_never_pad(nodes): + for i, n in enumerate(nodes): + if ( + isinstance(n, syntax_node.ValueNode) + and n.padding is None + and i + 1 < len(nodes) + and not isinstance(nodes[i + 1], syntax_node.PaddingNode) + ): + n.never_pad = True + if isinstance(n, syntax_node.ListNode): + _mark_never_pad(n.nodes) + + +def _parse_multiplier_bins(tally_numbers_node) -> list[MultiplierBin]: + """Parse a full FM card's ``tally numbers`` CST node into its bin sets.""" + raw = list(tally_numbers_node) + _mark_never_pad(raw) + items = _non_padding(raw) + top_groups = [n for n in items if _is_group(n)] + if not top_groups: + # FM parenthesization rule 3: the whole card is one bin set with one + # term, and no parens are used anywhere. + return [MultiplierBin.from_items(items)] + return [MultiplierBin.from_items(_group_body(g)) for g in top_groups] + + +class TallyMultiplier(DataInputAbstract, Numbered_MCNP_Object): + """An ``FMn`` tally multiplier card. + + Multiplies tally ``n``'s flux/current by a cross-section-derived + response function. Must be paired with a + :class:`~montepy.data_inputs.tally.Tally` of the same number -- see + :attr:`parent_tally`. + + .. versionadded:: 1.6.0b2 + """ + + _KEYS_TO_PRESERVE = {"_parent_tally"} + + @staticmethod + def _parser(): + return TallyParser() + + def _init_blank(self): + super()._init_blank() + self._old_number = self._generate_default_node(int, -1) + self._bins = [] + self._parsed_bins = [] + self._include_total = False + self._cumulative = False + self._parent_tally = None + + def _jit_light_init(self, input): + super()._jit_light_init(input) + self._old_number = self._input_number + + def _parse_tree(self): + super()._parse_tree() + num = self._input_number + self._old_number = copy.deepcopy(num) + self._number = num + self._parse_multiplier_body() + + def _generate_default_tree(self, **kwargs): + ret = {} + ret["start_pad"] = syntax_node.PaddingNode() + ret["classifier"] = syntax_node.ClassifierNode() + ret["classifier"].prefix = syntax_node.ValueNode( + self._class_prefix().upper(), str, padding=None, never_pad=True + ) + # A non-negative placeholder: ValueNode._reverse_engineer_formatting + # reserves a leading sign column for any token starting with "-", + # which would otherwise permanently corrupt this node's formatting + # once a real (positive) tally number is assigned to it. + ret["classifier"].number = self._generate_default_node(int, 1) + ret["keyword"] = syntax_node.ValueNode(None, str, padding=None) + tally_numbers = syntax_node.ListNode("tally numbers") + end_node = syntax_node.ValueNode(None, str) + ret["data"] = syntax_node.SyntaxNode( + "tally list", {"tally": tally_numbers, "end": end_node} + ) + ret["parameters"] = syntax_node.ParametersNode() + self._tree = syntax_node.SyntaxNode("blank data tree", ret) + + @args_checked + def __init__( + self, + input: InitInput = None, + number: ty.PositiveInt = None, + *, + jit_parse: bool = True, + ): + Numbered_MCNP_Object.__init__(self, input, number, jit_parse=jit_parse) + + @staticmethod + def _class_prefix() -> str: + return "fm" + + @staticmethod + def _has_number() -> bool: + return True + + @staticmethod + def _has_classifier() -> int: + return 1 + + @staticmethod + def _parent_collections(): + return () + + def _parse_multiplier_body(self): + if self._input is None: + return + tally_list = self._tree["data"] + end_node = tally_list["end"] + end_val = str(end_node.value).upper() if end_node.value is not None else None + self._include_total = end_val == "T" + self._cumulative = end_val == "C" + self._bins = _parse_multiplier_bins(tally_list["tally"]) + self._parsed_bins = list(self._bins) + + @make_prop_val_node("_old_number") + def old_number(self): + """The FM number as read from the input file.""" + pass + + @property + @needs_full_ast + def bins(self) -> list[MultiplierBin]: + """The bin sets (top-level parenthesized groups) of this FM card.""" + return list(self._bins) + + @args_checked + @needs_full_cst + def add_bin(self, bin_: MultiplierBin) -> None: + """Add a bin set to this FM card. + + Parameters + ---------- + bin_ : MultiplierBin + The bin set to add. + """ + self._bins.append(bin_) + + @args_checked + @needs_full_cst + def remove_bin(self, bin_: MultiplierBin) -> None: + """Remove a bin set from this FM card. + + Parameters + ---------- + bin_ : MultiplierBin + The bin set to remove. + """ + self._bins.remove(bin_) + + @property + @needs_full_ast + def include_total(self) -> bool: + """``True`` if a total bin (``T``) is appended.""" + return self._include_total + + @include_total.setter + @args_checked + @needs_full_cst + def include_total(self, value: bool): + self._include_total = value + if value: + self._cumulative = False + + @property + @needs_full_ast + def cumulative(self) -> bool: + """``True`` if the bins are cumulative (``C``), with the last being the total.""" + return self._cumulative + + @cumulative.setter + @args_checked + @needs_full_cst + def cumulative(self, value: bool): + self._cumulative = value + if value: + self._include_total = False + + @property + def parent_tally(self): + """The :class:`~montepy.data_inputs.tally.Tally` this multiplier is linked to.""" + return self._parent_tally + + def _link_to_parent(self, tally: "montepy.data_inputs.tally.Tally"): + if tally.multiplier is not None: + warnings.warn( + f"Multiple FM inputs were specified for tally: {self.old_number}.", + MalformedInputWarning, + ) + self._parent_tally = tally + + def link_to_problem(self, problem, *, deepcopy=False): + super().link_to_problem(problem) + + @args_checked + @needs_full_cst + def clone( + self, tally: "montepy.data_inputs.tally.Tally" = None + ) -> "TallyMultiplier": + """Create an independent copy of this ``FM`` card. + + Unlike the generic :meth:`~montepy.numbered_mcnp_object.Numbered_MCNP_Object.clone`, + a ``TallyMultiplier`` has no independent number or collection of its + own -- its number always matches its parent tally's -- so this is a + bespoke override. + + Parameters + ---------- + tally : Tally + The tally to link the clone to. Its number is copied onto the + clone, and the clone is registered as that tally's + :attr:`~montepy.data_inputs.tally.Tally.multiplier`. If omitted, + a detached, unregistered clone is returned instead. + + Returns + ------- + TallyMultiplier + The cloned ``FM`` card. + """ + ret = copy.deepcopy(self) + ret._parent_tally = None + if tally is not None: + ret.number = tally.number + ret._old_number.value = tally.number + tally.multiplier = ret + return ret + + def _update_values(self): + if self._bins != self._parsed_bins: + tally_numbers_node = self._tree["data"]["tally"] + tally_numbers_node.nodes.clear() + wrap_each = len(self._bins) > 1 + for bin_ in self._bins: + bin_node = bin_.node + if wrap_each or bin_._own_item_count() > 1: + group = syntax_node.ListNode("tally group") + group.append(syntax_node.PaddingNode("(")) + for n in bin_node.nodes: + group.append(n) + group.append(syntax_node.PaddingNode(")")) + tally_numbers_node.nodes.append(group) + else: + for n in bin_node.nodes: + tally_numbers_node.nodes.append(n) + self._parsed_bins = list(self._bins) + else: + for bin_ in self._bins: + if bin_._node is not None: + bin_.node + end_node = self._tree["data"]["end"] + if self._include_total: + end_node.value = "T" + elif self._cumulative: + end_node.value = "C" + else: + end_node.value = None + + def __str__(self): + try: + return f"{type(self).__name__}: {self.number}" + except Exception: + return f"{type(self).__name__}: (unparsed)" + + def __repr__(self): + try: + nbins = len(getattr(self, "_bins", [])) + return f"TALLY MULTIPLIER: {self.number}, bins: {nbins}" + except Exception: + return "TALLY MULTIPLIER: (unparsed)" diff --git a/montepy/data_inputs/tally_multiplier_type.py b/montepy/data_inputs/tally_multiplier_type.py new file mode 100644 index 000000000..89eb20a30 --- /dev/null +++ b/montepy/data_inputs/tally_multiplier_type.py @@ -0,0 +1,38 @@ +# Copyright 2024, Battelle Energy Alliance, LLC All Rights Reserved. + +from enum import unique, Enum + + +@unique +class ReactionOperator(Enum): + """The combinator between reaction numbers in an FM reaction list. + + See MCNP manual section 5.9.7, footnote 4: a space means multiply, a + colon means add, and a pound sign means subtract, with multiply binding + tighter than add/subtract. + + .. versionadded:: 1.6.0b2 + """ + + MULTIPLY = " " + ADD = ":" + SUBTRACT = "#" + + +@unique +class SpecialMultiplier(Enum): + """The ``c k`` special-multiplier flags (FM spec footnote 2). + + A closed, fixed 3-value set with no arithmetic use case, unlike + reaction numbers, so unlike :class:`~montepy.data_inputs.tally_multiplier.Reaction` + this is a real :class:`~enum.Enum`. + + .. versionadded:: 1.6.0b2 + """ + + INVERSE_WEIGHT = -1 + """The tally is multiplied by 1/weight; the tally is the number of tracks (or collisions for F5).""" + INVERSE_VELOCITY = -2 + """The tally is multiplied by 1/velocity; the tally is the neutron population or removal lifetime.""" + FIRST_INTERACTION_XS = -3 + """The tally is multiplied by the microscopic cross section of the first interaction.""" diff --git a/montepy/data_inputs/tally_type.py b/montepy/data_inputs/tally_type.py new file mode 100644 index 000000000..6dc064688 --- /dev/null +++ b/montepy/data_inputs/tally_type.py @@ -0,0 +1,40 @@ +# Copyright 2024, Battelle Energy Alliance, LLC All Rights Reserved. + +from enum import unique, Enum + + +@unique +class TallyType(Enum): + """The MCNP tally type, i.e. the last digit of an F-card's number. + + .. versionadded:: 1.6.0b2 + """ + + CURRENT = 1 + SURFACE_FLUX = 2 + CELL_FLUX = 4 + DETECTOR = 5 + ENERGY_DEPOSITION = 6 + FISSION_ENERGY_DEPOSITION = 7 + ENERGY_DETECTOR_PULSE = 8 + + +@unique +class Score(Enum): + """The physical quantity a :class:`~montepy.Tally` scores. + + A shallow analog of OpenMC's tally scores: for MontePy this is just the + quantity implied by the tally type digit (e.g. F4 always scores + :class:`Score.FLUX`). If an FM tally-multiplier card is linked to the + tally, :attr:`~montepy.Tally.scores` returns a list of + :class:`~montepy.data_inputs.tally_multiplier.MultiplierScore` instead of + this enum -- see :attr:`~montepy.Tally.multiplier`. + + .. versionadded:: 1.6.0b2 + """ + + CURRENT = 1 + FLUX = 2 + ENERGY_DEPOSITION = 6 + FISSION_ENERGY_DEPOSITION = 7 + PULSE_HEIGHT = 8 diff --git a/montepy/input_parser/tally_parser.py b/montepy/input_parser/tally_parser.py index a8d363ced..442bb73d3 100644 --- a/montepy/input_parser/tally_parser.py +++ b/montepy/input_parser/tally_parser.py @@ -1,5 +1,6 @@ # Copyright 2024, Battelle Energy Alliance, LLC All Rights Reserved. from montepy.input_parser.data_parser import DataParser +from montepy.input_parser.tokens import TallyLexer from montepy.input_parser import syntax_node @@ -13,6 +14,7 @@ class TallyParser(DataParser): """ debugfile = None + _lexer_class = TallyLexer @_("introduction tally_specification") def tally(self, p): @@ -33,7 +35,7 @@ def tally_specification(self, p): "tally list", {"tally": p.tally_numbers, "end": text} ) - @_("PARTICLE", "PARTICLE padding") + @_("PARTICLE", "PARTICLE padding", "TEXT", "TEXT padding") def end_phrase(self, p): """A non-zero number with or without padding. @@ -41,38 +43,147 @@ def end_phrase(self, p): ------- ValueNode a float ValueNode + + Notes + ----- + ``T`` (total) happens to lex as ``PARTICLE`` (the triton letter), but + ``C`` (cumulative, FM cards only) isn't a particle letter and lexes + as ``TEXT`` instead -- both alternatives are needed here. """ return self._flush_phrase(p, str) + # tally_numbers uses fresh rules (no number_sequence) so that the inherited + # `number_sequence → "(" number_sequence ")"` production is unreachable from + # the `tally` start symbol and cannot create a shift/reduce conflict with + # lparen_phrase inside tally_group. @_( - "tally_numbers tally_numbers", - "number_sequence", + "tally_flat_item", "tally_group", - "tally_numbers padding", + "tally_numbers tally_flat_item", + "tally_numbers tally_group", ) def tally_numbers(self, p): if hasattr(p, "tally_numbers"): ret = p.tally_numbers - ret.nodes["right"] += p.padding - return ret - if hasattr(p, "tally_numbers1"): - return syntax_node.SyntaxNode("tally tree", {"left": p[0], "right": p[1]}) + item = p.tally_flat_item if hasattr(p, "tally_flat_item") else p.tally_group + else: + ret = syntax_node.ListNode("tally numbers") + item = p[0] + # type() is, not isinstance(): ShortcutNode is a ListNode subclass + # and must stay unflattened, same as number_sequence in parser_base.py. + if type(item) is syntax_node.ListNode and item.name != "tally group": + for node in item.nodes: + ret.append(node) else: - left = syntax_node.PaddingNode(None) - right = syntax_node.PaddingNode(None) - return syntax_node.SyntaxNode( - "tally set", {"left": left, "tally": p[0], "right": right} - ) + ret.append(item) + return ret + + @_("number_phrase", "null_phrase", "shortcut_phrase") + def tally_flat_item(self, p): + return p[0] + + @_("lparen_phrase tally_group_body rparen_phrase") + def tally_group(self, p): + ret = syntax_node.ListNode("tally group") + ret.append(p.lparen_phrase) + for node in p.tally_group_body.nodes: + ret.append(node) + ret.append(p.rparen_phrase) + return ret + + @_("tally_group_item", "tally_group_body tally_group_item") + def tally_group_body(self, p): + if hasattr(p, "tally_group_body"): + ret = p.tally_group_body + else: + ret = syntax_node.ListNode("tally group body") + ret.append(p.tally_group_item) + return ret @_( - '"(" number_sequence ")"', - '"(" padding number_sequence ")"', + "number_phrase", + "null_phrase", + "shortcut_phrase", + "path_sep", + "lattice_phrase", + "universe_phrase", + "tally_group", + "reaction_operator", ) - def tally_group(self, p): - left = syntax_node.PaddingNode(p[0]) + def tally_group_item(self, p): + return p[0] + + @_("PARTICLE_SPECIAL", "PARTICLE_SPECIAL padding") + def path_sep(self, p): + return self._flush_phrase(p, str) + + @_('":"', '":" padding', "COMPLEMENT", "COMPLEMENT padding") + def reaction_operator(self, p): + """An FM tally-multiplier reaction-list operator: ``:`` (add) or ``#`` (subtract). + + Returns + ------- + ValueNode + a str ValueNode holding the raw operator symbol. + """ + return self._flush_phrase(p, str) + + @_('"[" lattice_body "]"', '"[" lattice_body "]" padding') + def lattice_phrase(self, p): + ret = syntax_node.ListNode("lattice phrase") + ret.append(syntax_node.PaddingNode(p[0])) + for node in p.lattice_body.nodes: + ret.append(node) if hasattr(p, "padding"): - left.append(p.padding) - right = syntax_node.PaddingNode(p[-1]) - return syntax_node.SyntaxNode( - "tally set", {"left": left, "tally": p.number_sequence, "right": right} - ) + ret.append(syntax_node.PaddingNode(p[2])) + ret.append(p.padding) + else: + ret.append(syntax_node.PaddingNode(p[2])) + return ret + + @_( + "lattice_item", + "lattice_body lattice_item", + 'lattice_body "," lattice_item', + 'lattice_body "," padding lattice_item', + ) + def lattice_body(self, p): + if hasattr(p, "lattice_body"): + ret = p.lattice_body + else: + ret = syntax_node.ListNode("lattice body") + if hasattr(p, "padding"): + # lattice_body "," padding lattice_item + ret.append(syntax_node.PaddingNode(p[1])) + ret.append(p.padding) + elif len(p) > 1 and isinstance(p[1], str) and p[1] == ",": + # lattice_body "," lattice_item + ret.append(syntax_node.PaddingNode(p[1])) + ret.append(p.lattice_item) + return ret + + @_( + "number_phrase", + "null_phrase", + 'number_phrase ":" number_phrase', + 'null_phrase ":" number_phrase', + 'number_phrase ":" null_phrase', + 'null_phrase ":" null_phrase', + ) + def lattice_item(self, p): + if len(p) > 1: + ret = syntax_node.ListNode("lattice range") + ret.append(p[0]) + ret.append(syntax_node.PaddingNode(p[1])) + ret.append(p[2]) + return ret + return p[0] + + @_("KEYWORD equals_sign number_phrase", "PARTICLE equals_sign number_phrase") + def universe_phrase(self, p): + token_val = p.KEYWORD if hasattr(p, "KEYWORD") else p.PARTICLE + ret = syntax_node.ListNode("universe phrase") + ret.append(syntax_node.ValueNode(token_val, str)) + ret.append(p.equals_sign) + ret.append(p.number_phrase) + return ret diff --git a/montepy/input_parser/tokens.py b/montepy/input_parser/tokens.py index d749d2e77..051f63e05 100644 --- a/montepy/input_parser/tokens.py +++ b/montepy/input_parser/tokens.py @@ -462,6 +462,25 @@ def TEXT(self, t): return t +class TallyLexer(DataLexer): + """A lexer for tally inputs. + + Adds ``[`` and ``]`` as literals so lattice index syntax like + ``[0 0 0]`` tokenizes correctly instead of being consumed by FILE_PATH. + Tally inputs never legitimately contain a ``READ``-style file path, and the + tally grammar never references the ``file_atom``/``file_name`` productions + that are the only consumers of a FILE_PATH token, so FILE_PATH is disabled + outright here (rather than merely narrowed) with a regex that can never + match, ``(?!)``. This is stricter than necessary for lattice brackets alone, + but safer: it guarantees no tally input can ever be silently misread as a + file path. + """ + + tokens = DataLexer.tokens + literals = DataLexer.literals | {"[", "]"} + FILE_PATH = r"(?!)" + + class SurfaceLexer(MCNP_Lexer): """A lexer for Surface inputs. diff --git a/montepy/mcnp_problem.py b/montepy/mcnp_problem.py index 7a3ffd06a..c69d1526e 100644 --- a/montepy/mcnp_problem.py +++ b/montepy/mcnp_problem.py @@ -7,7 +7,7 @@ import os import warnings -from montepy.data_inputs import mode, transform +from montepy.data_inputs import mode, tally as tally_mod, tally_multiplier, transform from montepy._cell_data_control import CellDataPrintController from montepy.utilities import * from montepy.cell import Cell @@ -17,6 +17,7 @@ from montepy.materials import Material, Materials from montepy.surfaces import surface, surface_builder from montepy.surface_collection import Surfaces +from montepy.tallies import Tallies import montepy.types as ty # weird way to avoid circular imports @@ -101,6 +102,7 @@ class MCNP_Problem: surface.Surface: Surfaces, Material: Materials, transform.Transform: Transforms, + tally_mod.Tally: Tallies, Universe: Universes, } @@ -419,6 +421,20 @@ def universes(self): """ return self._universes + @property + def tallies(self): + """A collection of the Tally objects in this problem. + + Returns + ------- + Tallies + a collection of the tally objects, ordered by the order + they appeared in the input file. + + .. versionadded:: 1.6.0b2 + """ + return self._tallies + @property def transforms(self): """The collection of transform objects in this problem. @@ -518,6 +534,10 @@ def parse_input( self._materials.append(obj, insert_in_data=False) elif isinstance(obj, transform.Transform): self._transforms.append(obj, insert_in_data=False) + elif isinstance( + obj, (tally_mod.Tally, tally_multiplier.TallyMultiplier) + ): + self._tallies.append(obj, insert_in_data=False) elif isinstance( obj, montepy.data_inputs.cell_modifier.CellModifierInput ): @@ -847,8 +867,12 @@ def parse( self.data_inputs.append(obj) if isinstance(obj, Material): self._materials.append(obj, insert_in_data=False) - if isinstance(obj, transform.Transform): + elif isinstance(obj, transform.Transform): self._transforms.append(obj, insert_in_data=False) + elif isinstance( + obj, (tally_mod.Tally, tally_multiplier.TallyMultiplier) + ): + self._tallies.append(obj, insert_in_data=False) return obj def full_parse(self): diff --git a/montepy/tallies.py b/montepy/tallies.py new file mode 100644 index 000000000..ea5e14bb7 --- /dev/null +++ b/montepy/tallies.py @@ -0,0 +1,94 @@ +# Copyright 2024, Battelle Energy Alliance, LLC All Rights Reserved. +import warnings + +import montepy +from montepy.exceptions import MalformedInputError, MalformedInputWarning +from montepy.numbered_object_collection import NumberedDataObjectCollection +from montepy.utilities import * + + +class Tallies(NumberedDataObjectCollection): + """ + A container of multiple :class:`~montepy.data_inputs.tally.Tally` instances. + + :param objects: the list of tallies to start with if needed + :type objects: list + + .. versionadded:: 1.6.0b2 + """ + + def __init__(self, objects=None, problem=None): + super().__init__(montepy.data_inputs.tally.Tally, objects, problem) + self._fm_queue = {} + self._multipliers = [] + + @property + def multipliers(self): + """The :class:`~montepy.data_inputs.tally_multiplier.TallyMultiplier` instances + held by this collection. + + Unlike the :class:`~montepy.data_inputs.tally.Tally` instances in this + collection, these are not stored in a + :class:`~montepy.numbered_object_collection.NumberedObjectCollection`, + as a ``TallyMultiplier``'s number is not an independent identity: it + always matches its parent tally's number, and two ``FM`` cards can + transiently share a number before one is linked. + + Returns + ------- + list + the tally multipliers ("FM" cards) in this problem. + """ + return list(self._multipliers) + + @args_checked + def append( + self, + obj: "montepy.data_inputs.tally.Tally | montepy.data_inputs.tally_multiplier.TallyMultiplier", + insert_in_data: bool = True, + **kwargs, + ): + if isinstance(obj, montepy.data_inputs.tally.Tally): + if obj.number in self._fm_queue: + fm = self._fm_queue.pop(obj.number) + fm._link_to_parent(obj) + obj._multiplier = fm + super().append(obj, insert_in_data=insert_in_data, **kwargs) + elif isinstance(obj, montepy.data_inputs.tally_multiplier.TallyMultiplier): + if obj not in self._multipliers: + self._multipliers.append(obj) + try: + tally = self[obj._old_number.value] + obj._link_to_parent(tally) + tally._multiplier = obj + except KeyError: + self._fm_queue[obj._old_number.value] = obj + if ( + insert_in_data + and self._problem is not None + and obj not in self._problem.data_inputs + ): + self._problem.data_inputs.append(obj) + + def _delete_hook(self, obj, **kwargs): + fm = getattr(obj, "_multiplier", None) + if fm is not None: + if fm in self._multipliers: + self._multipliers.remove(fm) + if self._problem is not None and fm in self._problem.data_inputs: + self._problem.data_inputs.remove(fm) + fm._parent_tally = None + warnings.warn( + f"Tally multiplier (FM) card for tally {obj.number} was removed " + "because its parent tally was deleted.", + MalformedInputWarning, + ) + super()._delete_hook(obj, **kwargs) + + def finalize_init(self, jit_parse: bool = False): + # Raise error for unflushed connection + for num, fm in self._fm_queue.items(): + raise MalformedInputError( + fm._input, + f'Tally multiplier "FM" input has no parent tally with number: {num}', + ) diff --git a/tests/inputs/test_tally.imcnp b/tests/inputs/test_tally.imcnp new file mode 100644 index 000000000..0803196bd --- /dev/null +++ b/tests/inputs/test_tally.imcnp @@ -0,0 +1,191 @@ +MCNP Test Model for Tallies +C cells +c +1 1 20 + -1000 $ dollar comment + imp:n,p=1 trcl=5 + u=1 +2 2 8 + -1005 + imp:n=1 + imp:p=0.5 + lat 1 + fill= 0:1 0:1 0:0 1 0 1 (5) +3 3 -1 + 1000 1005 -1010 + imp:n,p=1 +99 0 + 1010 + imp:n,p=0 + fill=1 +5 0 + #99 + imp:n,p=3 fill=1 (1 0.0 0.0) +c foo end comment + +C surfaces +1000 SO 1 +1005 RCC 0 1.5 -0.5 0 0 1 0.25 +1010 SO 3 + +C data +C materials +C UO2 5 atpt enriched +m1 92235.80c 5 & +92238.80c 95 +C Iron +m2 26054.80c 5.85 + 26056.80c 91.75 + 26057.80c 2.12 + 26058.80c 0.28 +C water +m3 1001.80c 2 + 8016.80c 1 +MT3 lwtr.23t +TR5 0.0 0.0 1.0 +C execution +ksrc 0 0 0 +kcode 100000 1.000 50 1050 +phys:p j 1 2j 1 +mode n p +vol NO 2J 1 1.5 J +C tallies +C F1: surface current, multiple particles +fc1 Surface current tally comment +f1:n,p 1000 +C F2: average surface flux, photons only +fc2 Average surface flux +f2:p 1005 +C F4: cell flux, simple list +fc4 Cell flux simple +f4:n 1 2 3 +C with energy bins +e4 0.625e-6 +C F4: cell flux with two groups +fc14 Cell flux groups +f14:n (1 2) (3) +C F4: cell flux with two groups and total +fc24 Cell flux groups with total +f24:n (1 2) (3) T +C F4: cell flux with interpolated group and total +fc34 Cell flux interpolated group +f34:n (1 3i 5) T +C F4: cell flux with multiple groups +fc44 Cell flux multiple groups +f44:n (1 3i 5) (7 8 9) +C F4: cell flux with multiple groups and total +fc54 Cell flux multiple groups with total +f54:n (1 3i 5) (7 8 9) T +C F4: cell flux with nested universe path +fc64 Cell 1 in universe 1 +f64:n (1<1) +C F4: deeper nested universe path +fc74 Cell nested two universes deep +f74:n (1<1< 2) +C F4: lattice element tally +fc84 Lattice element tally +f84:n (1[0 0 0]<2) +C F4: lattice element in nested universe +fc94 Lattice element nested universe tally +f94:n (1[0 0 0]<2<3) +C F6: energy deposition, photons +fc6 Energy deposition photon +f6:p 1 3 +C F7: fission energy deposition, neutrons +fc7 Fission energy deposition +f7:n 1 +C F8: pulse height tally, photons +fc8 Pulse height tally +f8:p 3 +C tally multiplier for f4 +fm4 (1.0) +C tally multiplier with reaction specification +fm14 (1.0 1 -6) +C tally multiplier with comment and multiple bins +fm44 (1.0 1 444) $ reaction rate +C tally multiplier with multiple reaction lists (separate bins) in one multiplier set +fm54:n (1.0 26 (16) (103)) +C tally multiplier with multiple multiplier sets (each parenthesized, whole set wrapped) +fm64:n ((1.0 26 16) (2.0 27 102)) +C tally multiplier reaction list: multiply (space) within one bin +fm74:n (1.0 26 16 103) +C tally multiplier reaction list: add (colon) operator +fm84:n (1.0 26 16:103) +C tally multiplier reaction list: subtract (pound) operator +fm94:n (1.0 26 16#103) +C tally multiplier reaction list precedence: (16*103) + 104 +fm104:n (1.0 26 16 103 : 104) +C tally for fm204 +f204:n 1 2 3 +C tally multiplier reaction list precedence: 16 + (103*104) +fm204:n (1.0 26 16 : 103 104) +C tally for fm214 +f214:n 1 2 3 +C tally multiplier reaction list precedence: (16*103) + (104*105) +fm214:n (1.0 26 16 103 : 104 105) +C tally for fm224 +f224:n 1 2 3 +C tally multiplier reaction list precedence: 16 - (103*104) +fm224:n (1.0 26 16 # 103 104) +C tally for fm234 +f234:n 1 2 3 +C tally multiplier reaction list precedence: (16*103) - 104 +fm234:n (1.0 26 16 103 # 104) +C tally for fm244 +f244:n 1 2 3 +C tally multiplier reaction list precedence: 16 + 103 - (104*105) +fm244:n (1.0 26 16 : 103 # 104 105) +C tally for fm254 +f254:n 1 2 3 +C tally multiplier reaction list precedence: (16*103) + 104 - (105*106) +fm254:n (1.0 26 16 103 : 104 # 105 106) +C tally for fm264 +f264:n 1 2 3 +C tally multiplier reaction list precedence: 16 - 103 + (104*105) +fm264:n (1.0 26 16 # 103 : 104 105) +C tally multiplier single-layer attenuator +fm114:n (1.0 -1 26 0.5) +C tally for fm124 +f124:n 1 2 3 +C tally multiplier multi-layer attenuator +fm124:n (1.0 -1 26 0.5 27 -0.3) +C tally for fm134 +f134:n 1 2 3 +C tally multiplier combining multiplier sets with an attenuator set +fm134:n ((1.0 26 16) (2.0 27 102) (3.0 -1 28 0.1)) +C tally for fm144 +f144:n 1 2 3 +C tally multiplier negative constant (type-4 atom density normalization) +fm144:n (-1.0 26 103) +C tally multiplier special option: 1/weight +fm154:n 1 -1 +C tally for fm164 +f164:n 1 2 3 +C tally multiplier special option: 1/velocity +fm164:n (1 -2) +C tally multiplier special option: microscopic xs of first interaction +fm174:n (1 -3) +C tally multiplier with total flag +fm184:n (1.0 26 16) (2.0 27 102) T +C tally multiplier with cumulative flag +fm194:n (1.0 26 16) (2.0 27 102) C +C tally segment, single surface +fs1 -1000 +C tally segment with total flag +fs2 -1005 t +C tally segment, multiple surfaces with total and complement +fs44 -1000 -1005 t c +C tally segment, positive surface orientation +fs54 +1000 c +F154:n,p (1 < (2[0 0 0] 2[0 1 0]) < 5) +F1464:n (1 < 2[0:1 0:1 0:0] < 5) +c --- FORM 7: Range of elements, collapsed into one averaged bin -------------- +F174:n (1 < (2[0:1 0:1 0:0]) < 5) +c --- FORM 8: Non-contiguous individual elements, separate bins --------------- +F184:n (1 < 2[0 0 0, 0 1 0] < 5) +c --- FORM 9: No brackets -- all elements, all merged into one bin ------------ +F194:n (1 < 2 < 5) +F104:n ((u=1) < 2[0 0 0] < 5) +F114:n (u=1 < 2[0 0 0] < 5) +C tally print specifier +fq4 f e diff --git a/tests/test_tally.py b/tests/test_tally.py index bfbf52441..1550bf925 100644 --- a/tests/test_tally.py +++ b/tests/test_tally.py @@ -1,10 +1,27 @@ # Copyright 2024-2025, Battelle Energy Alliance, LLC All Rights Reserved. +import io + import pytest import montepy from montepy.data_inputs.data_parser import parse_data +from montepy.data_inputs.tally import ( + EnergyDepositionTally, + F1Tally, + F4Tally, + F6Tally, + FlatGroup, + LatticeIndex, + ParticleFilter, + PathGroup, + SpatialFilter, + TallyGroup, +) +from montepy.data_inputs.tally_type import Score, TallyType from montepy.input_parser.block_type import BlockType from montepy.input_parser.mcnp_input import Input +from montepy.input_parser.tally_parser import TallyParser +from montepy.input_parser.tokens import TallyLexer class TestTallyParser: @@ -68,8 +85,515 @@ def test_de_parsing_jail(_, line): data.data +class TestTallyPathSyntax: + """Tests for complex MCNP tally path syntax (universe paths, lattice elements).""" + + _parser = TallyParser() + _lexer = TallyLexer() + + @pytest.mark.parametrize( + "line", + [ + "f64:n (1<1)", + "f74:n (1<1< 2)", + "f84:n (1[0 0 0]<2)", + "f94:n (1[0 0 0]<2<3)", + "F154:n,p (1 < (2[0 0 0] 2[0 1 0]) < 5)", + "F1464:n (1 < 2[0:1 0:1 0:0] < 5)", + "F174:n (1 < (2[0:1 0:1 0:0]) < 5)", + "F184:n (1 < 2[0 0 0, 0 1 0] < 5)", + "F184:n (1 < 2[0 0 0,0 1 0] < 5)", # comma with no trailing padding + "F194:n (1 < 2 < 5)", + "F104:n ((u=1) < 2[0 0 0] < 5)", + "F114:n (u=1 < 2[0 0 0] < 5)", + ], + ) + def test_tally_path_parsing(self, line): + result = self._parser.parse(self._lexer.tokenize(line)) + assert result is not None, f"TallyParser failed to parse: {line}" + + class TestFmesh: # this is hacky; just makes sure it doesn't crash @pytest.mark.parametrize("line", ["fmesh14:n vec=0 0 0", "fmesh14:n vec=0, 0, 0"]) def test_fmesh_parse(_, line): parse_data(line) + + +@pytest.fixture +def tally_problem(): + return montepy.read_input("tests/inputs/test_tally.imcnp") + + +class TestTallyObject: + """Tests for the Tally object model: clone, clone_as, scores, filters.""" + + def test_clone_same_type(self, tally_problem): + f4 = tally_problem.tallies[4] + clone = f4.clone() + assert clone.number != f4.number + assert clone.number % 10 == 4 + assert clone in tally_problem.tallies + assert list(clone.cells.numbers) == list(f4.cells.numbers) + + def test_clone_as_class(self, tally_problem): + f4 = tally_problem.tallies[4] + new = f4.clone_as(F6Tally) + assert isinstance(new, EnergyDepositionTally) + assert new.number % 10 == 6 + assert new in tally_problem.tallies + assert list(new.cells.numbers) == list(f4.cells.numbers) + assert new.scores == [Score.ENERGY_DEPOSITION] + + def test_clone_as_enum(self, tally_problem): + f4 = tally_problem.tallies[4] + new = f4.clone_as(TallyType.ENERGY_DEPOSITION) + assert isinstance(new, EnergyDepositionTally) + assert new.number % 10 == 6 + assert new in tally_problem.tallies + + def test_clone_as_incompatible_category(self, tally_problem): + f1 = tally_problem.tallies[1] + with pytest.raises(ValueError): + f1.clone_as(F4Tally) + + def test_clone_as_bad_type(self, tally_problem): + f4 = tally_problem.tallies[4] + with pytest.raises(TypeError): + f4.clone_as("f6") + + @pytest.mark.parametrize( + "number, expected", + [ + (1, [Score.CURRENT]), + (2, [Score.FLUX]), + # not 4: it has a linked fm4 card, tested in test_tally_multiplier.py + (34, [Score.FLUX]), + (6, [Score.ENERGY_DEPOSITION]), + (7, [Score.FISSION_ENERGY_DEPOSITION]), + (8, [Score.PULSE_HEIGHT]), + ], + ) + def test_scores_default(self, tally_problem, number, expected): + assert tally_problem.tallies[number].scores == expected + + def test_filters_default(self, tally_problem): + f1 = tally_problem.tallies[1] # f1:n,p 1000 + filters = f1.filters + assert len(filters) == 2 + particle_filter, spatial_filter = filters + assert isinstance(particle_filter, ParticleFilter) + assert isinstance(spatial_filter, SpatialFilter) + assert set(particle_filter.particles) == set(f1.particle_classifiers) + assert spatial_filter.groups == f1.groups + + def test_add_cell_before_full_parse_preserves_existing_groups(self, tally_problem): + # Regression test: add_cell/add_group/add_path_group must trigger a + # full parse *before* mutating _groups, or the mutation is silently + # lost the next time a @needs_full_ast getter forces a full parse. + f4 = tally_problem.tallies[4] # f4:n 1 2 3 + assert not f4.fully_parsed + new_cell = tally_problem.cells[1].clone() + f4.add_cell(new_cell) + assert f4.fully_parsed + numbers = list(f4.cells.numbers) + assert {1, 2, 3}.issubset(set(numbers)) + assert new_cell.number in numbers + + def test_from_input_invalid_tally_type_digit(self): + with pytest.raises(montepy.exceptions.MalformedInputError): + parse_data(Input(["f3:n 1 2 3"], BlockType.DATA)) + + def test_direct_tally_construction_invalid_digit(self): + # Tally() bypasses from_input's dispatch-time digit check, so + # _parse_tally_body must catch it too, once fully parsed. + from montepy.data_inputs.tally import Tally + + t = Tally(Input(["f3:n 1 2 3"], BlockType.DATA), jit_parse=True) + with pytest.raises(montepy.exceptions.MalformedInputError): + t.groups + + def test_jump_in_tally_numbers_is_skipped(self): + t = F4Tally(Input(["f4:n 1 J 3"], BlockType.DATA)) + assert [g.old_numbers[0] for g in t.groups] == [1, 3] + + def test_universe_spec_nested_in_bare_group(self): + # A universe designator wrapped in its own parens, inside a group + # with no `<` path separator at all: _extract_universe_spec_from_nodes + # must recurse into the nested "tally group" node to find it. + t = F4Tally(Input(["f4:n (5 (u=2) 7)"], BlockType.DATA)) + assert t.groups[0].universe_spec == 2 + + def test_from_input_bare_prefix_no_number(self): + # A tally card with no number at all: the light JIT parser succeeds + # (it doesn't validate structure), so from_input's own "has no + # number" check and except-Exception fallback both have to run, + # and the fallback's real full construction is what actually + # surfaces the parsing error. + with pytest.raises(montepy.exceptions.ParsingError): + parse_data(Input(["f"], BlockType.DATA)) + + def test_clone_as_non_tally_type_raises_type_error(self, tally_problem): + f34 = tally_problem.tallies[34] + with pytest.raises(TypeError): + f34.clone_as(str) + + def test_clone_as_custom_subclass_outside_category_raises(self, tally_problem): + from montepy.data_inputs.tally import Tally + + class CustomTally(Tally): + _TALLY_TYPE = TallyType.CELL_FLUX + + f34 = tally_problem.tallies[34] + with pytest.raises(ValueError): + f34.clone_as(CustomTally) + + def test_str_and_repr_on_uninitialized_tally(self): + # __new__ bypasses __init__ entirely, leaving no _number/_groups -- + # exactly the partially-constructed state __str__/__repr__'s except + # branches exist to report gracefully. + from montepy.data_inputs.tally import Tally + + t = Tally.__new__(Tally) + assert str(t) == "Tally: (unparsed)" + assert repr(t) == "TALLY: (unparsed)" + + def test_all_fixture_tallies_fully_parse(self, tally_problem): + # Force a full parse of every tally in the fixture, including the + # path/lattice/universe forms (104, 114, 154, 1464, 174, 184, 194) + # that no other test touches directly. + for t in tally_problem.tallies: + groups = t.groups + assert groups is not None + assert len(groups) > 0 + _ = t.include_total + + def test_include_total(self, tally_problem): + assert tally_problem.tallies[24].include_total is True + assert tally_problem.tallies[34].include_total is True + assert tally_problem.tallies[1].include_total is False + + def test_contains_linked_flat_group(self, tally_problem): + f34 = tally_problem.tallies[34] + assert tally_problem.cells[1] in f34 + assert tally_problem.cells[99] not in f34 + + def test_contains_linked_path_group(self, tally_problem): + f64 = tally_problem.tallies[64] + assert tally_problem.cells[1] in f64 + + def test_contains_forces_full_parse(self, tally_problem): + # __contains__ is @needs_full_ast: membership must be checked + # against real data, not silently read as False while still JIT. + f34 = tally_problem.tallies[34] + assert not f34.fully_parsed + assert tally_problem.cells[1] in f34 + assert f34.fully_parsed + + def test_contains_unlinked_flat_group_by_number_only(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA)) + + class FakeCellByNumber: + number = 2 + + assert FakeCellByNumber() in t + + def test_contains_unlinked_flat_group_by_old_number(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA)) + + class FakeCellByOldNumber: + old_number = 2 + + assert FakeCellByOldNumber() in t + + def test_number_validator_rejects_wrong_digit(self, tally_problem): + f34 = tally_problem.tallies[34] + with pytest.raises(ValueError): + f34.number = 16 + + def test_clone_unlinked_tally(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA)) + clone = t.clone() + assert clone.number != t.number + assert clone.number % 10 == 4 + + def test_clone_as_unlinked_tally(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA)) + new = t.clone_as(F6Tally) + assert isinstance(new, EnergyDepositionTally) + assert new.number % 10 == 6 + + def test_clone_as_same_type_delegates_to_clone(self, tally_problem): + f34 = tally_problem.tallies[34] + new = f34.clone_as(type(f34)) + assert new.number % 10 == 4 + assert new in tally_problem.tallies + + def test_clone_starting_number_needs_alignment(self, tally_problem): + # starting_number's digit (9) doesn't match the tally type's digit + # (4), forcing _align_to_type's "aligned < start" +10 branch. + f34 = tally_problem.tallies[34] + new = f34.clone(starting_number=9) + assert new.number % 10 == 4 + assert new.number >= 14 + + def test_blank_tally_construction(self): + t = F4Tally() + assert t.groups == [] + + +class TestTallyBuilders: + """Tests for the from-scratch tally-building API: add_surface, add_group, + add_path_group, and PathGroup.inside chaining.""" + + def test_surface_tally_add_surface_and_group(self, tally_problem): + f1 = tally_problem.tallies[1] # f1:n,p 1000 -- only surface 1000 so far + s = tally_problem.surfaces[1005] + f1.add_surface(s) + assert s in f1.surfaces + assert f1.groups[-1].old_numbers == [s.number] + + s_already_present = tally_problem.surfaces[1000] + s_new = tally_problem.surfaces[1010] + f1.add_group([s_already_present, s_new]) + assert f1.groups[-1].is_grouped + assert set(f1.groups[-1].old_numbers) == { + s_already_present.number, + s_new.number, + } + assert s_new in f1.surfaces + + def test_surface_tally_add_path_group_and_inside_chaining(self, tally_problem): + f1 = tally_problem.tallies[1] + s1, s2 = tally_problem.surfaces[1000], tally_problem.surfaces[1005] + pg = f1.add_path_group(s1) + assert isinstance(pg, PathGroup) + pg.inside(s2, lattice=[0, 0, 0]) + assert pg in f1.groups + assert len(pg.levels) == 2 + + # Re-linking must recurse into the newly added PathGroup too. + f1.link_to_problem(tally_problem) + assert s1 in f1.surfaces and s2 in f1.surfaces + + def test_link_group_surfaces_skips_missing_surface(self, tally_problem): + f1 = tally_problem.tallies[1] + _ = f1.surfaces # force full parse + f1._groups.append(FlatGroup([99999], is_grouped=False)) + f1.link_to_problem(tally_problem) + assert 99999 not in list(f1.surfaces.numbers) + + def test_cell_tally_add_group(self, tally_problem): + f34 = tally_problem.tallies[34] # cells 1,2,3,5 -- 99 is new + c_already_present = tally_problem.cells[1] + c_new = tally_problem.cells[99] + f34.add_group([c_already_present, c_new]) + assert f34.groups[-1].is_grouped + assert set(f34.groups[-1].old_numbers) == { + c_already_present.number, + c_new.number, + } + assert c_new in f34.cells + + def test_cell_tally_add_path_group_and_inside_chaining(self, tally_problem): + f34 = tally_problem.tallies[34] + c1, c2 = tally_problem.cells[1], tally_problem.cells[2] + pg = f34.add_path_group(c1) + assert isinstance(pg, PathGroup) + pg.inside(c2) + assert pg in f34.groups + assert len(pg.levels) == 2 + + +def verify_export(tally): + """Format ``tally`` to MCNP text, re-parse it standalone, and confirm + the result is equivalent. Mirrors the ``verify_export`` convention in + ``tests/test_surfaces.py``/``tests/test_cell_problem.py``.""" + output = tally.format_for_mcnp_input((6, 3, 0)) + joined = "\n".join(output) + assert joined == tally.mcnp_str((6, 3, 0)) + new_tally = type(tally)(joined) + assert new_tally.number == tally.number + assert new_tally.tally_type == tally.tally_type + assert new_tally.include_total == tally.include_total + assert len(new_tally.groups) == len(tally.groups) + for old_group, new_group in zip(tally.groups, new_tally.groups): + assert isinstance(new_group, type(old_group)) + if isinstance(old_group, FlatGroup): + assert old_group.is_grouped == new_group.is_grouped + return new_tally + + +def verify_prob_export(problem, tally): + """Write the whole ``problem`` out and re-read it, returning the + equivalent tally from the new problem. The only test shape that can + catch bugs in problem/collection-level registration (e.g. an FM card + silently missing from ``data_inputs``), since a per-object + :func:`verify_export` check never sees the problem at all.""" + with io.StringIO() as fh: + problem.write_problem(fh) + fh.seek(0) + new_problem = montepy.read_input(fh) + return new_problem.tallies[tally.number] + + +class TestGroupRoundTrip: + """Mutating a Tally's groups through the public API must be reflected in + mcnp_str(), not just in the in-memory Python state. These lock in the + Tally._update_values / FlatGroup/PathGroup node-generation layer.""" + + def test_unmodified_tally_round_trips_exactly(self, tally_problem): + for number in (4, 14, 24, 34, 44, 54, 64, 74, 84, 94): + tally = tally_problem.tallies[number] + before = tally.mcnp_str() + tally.full_parse() + assert tally.mcnp_str() == before + verify_export(tally) + + def test_add_cell_reflected_in_mcnp_str(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA), jit_parse=False) + cell = montepy.Cell() + cell.number = 99 + t.add_cell(cell) + assert "99" in t.mcnp_str() + assert "1" in t.mcnp_str() and "2" in t.mcnp_str() and "3" in t.mcnp_str() + verify_export(t) + + def test_blank_tally_add_cell_writes_valid_card(self): + t = F4Tally() + t.number = 4 + cell = montepy.Cell() + cell.number = 1 + t.add_cell(cell) + text = t.mcnp_str() + assert "1" in text + # the unfixed bug wrote out only 'F 4 ', with no cell number at all + assert text.strip() != "F 4" + + def test_add_group_reflected_in_mcnp_str(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA), jit_parse=False) + c1, c2 = montepy.Cell(), montepy.Cell() + c1.number = 10 + c2.number = 11 + t.add_group([c1, c2]) + text = t.mcnp_str() + assert "10" in text and "11" in text + assert "(" in text and ")" in text + verify_export(t) + + def test_renumbered_cell_reflected_in_mcnp_str(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA), jit_parse=False) + cell = montepy.Cell() + cell.number = 99 + t.add_cell(cell) + cell.number = 199 + text = t.mcnp_str() + assert "199" in text + assert "99" not in text.replace("199", "") + verify_export(t) + + def test_include_total_settable(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA), jit_parse=False) + assert not t.include_total + t.include_total = True + assert t.include_total + assert t.mcnp_str().strip().endswith("T") + t.include_total = False + assert not t.include_total + assert not t.mcnp_str().strip().endswith("T") + + def test_remove_cell_reflected_in_mcnp_str(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA), jit_parse=False) + cell = montepy.Cell() + cell.number = 99 + t.add_cell(cell) + assert "99" in t.mcnp_str() + t.remove_cell(cell) + assert "99" not in t.mcnp_str() + assert cell not in t.cells + + def test_remove_cell_keeps_cell_if_referenced_elsewhere(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA), jit_parse=False) + cell = montepy.Cell() + cell.number = 99 + other = montepy.Cell() + other.number = 98 + t.add_cell(cell) + t.add_group([cell, other]) + t.remove_cell(cell) + assert cell in t.cells + assert "99" in t.mcnp_str() + + def test_remove_cell_raises_if_not_a_single_cell_group(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA), jit_parse=False) + cell = montepy.Cell() + cell.number = 99 + with pytest.raises(ValueError): + t.remove_cell(cell) + + def test_remove_group_reflected_in_mcnp_str(self): + t = F4Tally(Input(["f4:n 1 2 3"], BlockType.DATA), jit_parse=False) + c1, c2 = montepy.Cell(), montepy.Cell() + c1.number = 10 + c2.number = 11 + t.add_group([c1, c2]) + group = t.groups[-1] + t.remove_group(group) + text = t.mcnp_str() + assert "10" not in text and "11" not in text + assert c1 not in t.cells and c2 not in t.cells + + def test_remove_surface_reflected_in_mcnp_str(self, tally_problem): + f1 = tally_problem.tallies[1] + s = tally_problem.surfaces[1005] + f1.add_surface(s) + assert "1005" in f1.mcnp_str() + f1.remove_surface(s) + assert "1005" not in f1.mcnp_str() + assert s not in f1.surfaces + + +class TestReprAndEquality: + """Smoke tests for __repr__ and __eq__-against-wrong-type on the small + standalone value objects in tally.py. Coverage only counts a line as hit + if it executes, and pytest only reprs on assertion *failure*, so these + branches need explicit exercising.""" + + def test_lattice_index(self): + li = LatticeIndex([1, (2, 3)]) + assert li.dimensions == [1, (2, 3)] + assert "LatticeIndex" in repr(li) + + def test_tally_group_contains_not_implemented(self): + with pytest.raises(NotImplementedError): + 1 in TallyGroup() + + def test_particle_filter_eq_wrong_type_and_repr(self): + obj = ParticleFilter([montepy.Particle.NEUTRON]) + assert obj != "not a filter" + assert "ParticleFilter" in repr(obj) + + def test_spatial_filter_eq_wrong_type_and_repr(self): + obj = SpatialFilter([FlatGroup([1], is_grouped=False)]) + assert obj != "not a filter" + assert "SpatialFilter" in repr(obj) + + def test_flat_group_properties_and_repr(self): + fg = FlatGroup([1, 2], is_grouped=True, universe_spec=3) + assert fg.old_numbers == [1, 2] + assert fg.is_grouped is True + assert fg.universe_spec == 3 + assert "FlatGroup" in repr(fg) + + def test_path_group_levels_and_repr(self): + pg = PathGroup([FlatGroup([1], is_grouped=False)]) + assert len(pg.levels) == 1 + assert "PathGroup" in repr(pg) + + def test_path_group_empty_levels_contains(self): + assert 1 not in PathGroup([]) + + def test_tally_parent_collections(self): + from montepy.data_inputs.tally import Tally + + assert Tally._parent_collections() == () diff --git a/tests/test_tally_multiplier.py b/tests/test_tally_multiplier.py new file mode 100644 index 000000000..c7be10bff --- /dev/null +++ b/tests/test_tally_multiplier.py @@ -0,0 +1,691 @@ +# Copyright 2024-2025, Battelle Energy Alliance, LLC All Rights Reserved. +import io +import warnings + +import pytest + +import montepy +from montepy.data_inputs.data_parser import parse_data +from montepy.data_inputs.tally_multiplier import ( + AttenuatorLayer, + AttenuatorSet, + MultiplierBin, + MultiplierScore, + MultiplierSet, + Reaction, + SpecialMultiplierSet, + TallyMultiplier, +) +from montepy.data_inputs.tally_multiplier_type import ( + ReactionOperator, + SpecialMultiplier, +) +from montepy.data_inputs.tally_type import Score +from montepy.input_parser.block_type import BlockType +from montepy.input_parser.mcnp_input import Input + + +@pytest.fixture +def tally_problem(): + return montepy.read_input("tests/inputs/test_tally.imcnp") + + +def verify_export(fm): + """Format ``fm`` to MCNP text, re-parse it standalone, and confirm the + result is equivalent. Mirrors the ``verify_export`` convention in + ``tests/test_surfaces.py``/``tests/test_cell_problem.py``.""" + output = fm.format_for_mcnp_input((6, 3, 0)) + joined = "\n".join(output) + assert joined == fm.mcnp_str((6, 3, 0)) + new_fm = type(fm)(joined) + assert new_fm.number == fm.number + assert new_fm.include_total == fm.include_total + assert new_fm.cumulative == fm.cumulative + assert len(new_fm.bins) == len(fm.bins) + return new_fm + + +def verify_prob_export(problem, fm): + """Write the whole ``problem`` out and re-read it, returning the + equivalent multiplier from the new problem. The only test shape that + can catch bugs in problem/collection-level registration, since a + per-object :func:`verify_export` check never sees the problem at all.""" + with io.StringIO() as fh: + problem.write_problem(fh) + fh.seek(0) + new_problem = montepy.read_input(fh) + return new_problem.tallies[fm.number].multiplier + + +# Every "fm" line currently in tests/inputs/test_tally.imcnp, kept in sync +# with that fixture so the grammar round-trip test below actually exercises +# what's on disk. +FM_FIXTURE_LINES = [ + "fm4 (1.0)", + "fm14 (1.0 1 -6)", + "fm44 (1.0 1 444)", + "fm54:n (1.0 26 (16) (103))", + "fm64:n ((1.0 26 16) (2.0 27 102))", + "fm74:n (1.0 26 16 103)", + "fm84:n (1.0 26 16:103)", + "fm94:n (1.0 26 16#103)", + "fm104:n (1.0 26 16 103 : 104)", + "fm204:n (1.0 26 16 : 103 104)", + "fm214:n (1.0 26 16 103 : 104 105)", + "fm224:n (1.0 26 16 # 103 104)", + "fm234:n (1.0 26 16 103 # 104)", + "fm244:n (1.0 26 16 : 103 # 104 105)", + "fm254:n (1.0 26 16 103 : 104 # 105 106)", + "fm264:n (1.0 26 16 # 103 : 104 105)", + "fm114:n (1.0 -1 26 0.5)", + "fm124:n (1.0 -1 26 0.5 27 -0.3)", + "fm134:n ((1.0 26 16) (2.0 27 102) (3.0 -1 28 0.1))", + "fm144:n (-1.0 26 103)", + "fm154:n 1 -1", + "fm164:n (1 -2)", + "fm174:n (1 -3)", + "fm184:n (1.0 26 16) (2.0 27 102) T", + "fm194:n (1.0 26 16) (2.0 27 102) C", +] + + +class TestGrammar: + """Locks in the tally_parser.py grammar fix: reaction-list operators + (space/colon/COMPLEMENT) and the C (cumulative) end flag must parse + outside of lattice-index brackets. + """ + + @pytest.mark.parametrize("line", FM_FIXTURE_LINES) + def test_fm_line_full_parses(self, line): + data = parse_data(Input([line], BlockType.DATA)) + data.full_parse() # must not raise ParsingError + assert isinstance(data, TallyMultiplier) + + def test_mixed_particle_classifier(self): + # fm4 has no particle designator; fm54:n does. DataInputAbstract's + # __enforce_name rejects a classifier mismatch, so both must be + # legal -- regression test for _has_classifier() == 1, not 0. + for line in ["fm4 (1.0)", "fm54:n (1.0 26 (16) (103))"]: + data = parse_data(Input([line], BlockType.DATA)) + data.full_parse() + + +class TestReactionExpressionPrecedence: + @pytest.mark.parametrize( + "number, expected", + [ + (104, Reaction(16) * Reaction(103) + Reaction(104)), + (204, Reaction(16) + Reaction(103) * Reaction(104)), + (214, Reaction(16) * Reaction(103) + Reaction(104) * Reaction(105)), + (224, Reaction(16) - Reaction(103) * Reaction(104)), + (234, Reaction(16) * Reaction(103) - Reaction(104)), + (244, Reaction(16) + Reaction(103) - Reaction(104) * Reaction(105)), + ( + 254, + Reaction(16) * Reaction(103) + + Reaction(104) + - Reaction(105) * Reaction(106), + ), + (264, Reaction(16) - Reaction(103) + Reaction(104) * Reaction(105)), + ], + ) + def test_precedence_matches_hand_built_expression( + self, tally_problem, number, expected + ): + fm = tally_problem.tallies[number].multiplier + [reaction] = fm.bins[0].terms[0].reactions + assert reaction == expected + + +class TestOperatorOverloading: + def test_multiply_binds_tighter_than_add(self): + expr = Reaction(16) * Reaction(103) + Reaction(104) + assert expr.left == Reaction(16) * Reaction(103) + assert expr.operator == ReactionOperator.ADD + assert expr.right == Reaction(104) + + def test_int_first_forms(self): + assert 16 * Reaction(103) == Reaction(16) * Reaction(103) + assert 16 + Reaction(103) == Reaction(16) + Reaction(103) + + def test_rsub_swaps_operand_order(self): + # 16 - Reaction(103) must be Reaction(16) - Reaction(103), not the reverse. + expr = 16 - Reaction(103) + assert expr.left == Reaction(16) + assert expr.right == Reaction(103) + + def test_rand_builds_multiplier_set(self): + built = 26 & Reaction.CAPTURE + assert built == MultiplierSet(1.0, 26, [Reaction.CAPTURE]) + + def test_rand_with_material_object(self): + mat = montepy.Material() + mat.number = 26 + built = mat & Reaction.CAPTURE + assert built == MultiplierSet(1.0, 26, [Reaction.CAPTURE]) + + def test_rmul_scales_constant(self): + built = 1.5 * (26 & Reaction.CAPTURE) + assert built == MultiplierSet(1.5, 26, [Reaction.CAPTURE]) + + def test_named_reaction_constants_compose(self): + expr = Reaction.TOTAL - Reaction.CAPTURE - Reaction.INELASTIC_SCATTER + assert expr.left == Reaction.TOTAL - Reaction.CAPTURE + assert expr.right == Reaction.INELASTIC_SCATTER + assert expr.operator == ReactionOperator.SUBTRACT + + def test_named_reaction_constants_are_unique(self): + constants = { + name: obj + for name, obj in vars(Reaction).items() + if isinstance(obj, Reaction) + } + # sanity: this should have picked up more than just a handful, + # confirming the Appendix B transcription actually landed. + assert len(constants) > 400 + numbers_seen = {} + for name, reaction in constants.items(): + if reaction.number in numbers_seen: + pytest.fail( + f"Reaction.{name} (MT {reaction.number}) collides with " + f"Reaction.{numbers_seen[reaction.number]}" + ) + numbers_seen[reaction.number] = name + + +class TestAttenuator: + def test_and_chains_layers(self): + att = AttenuatorSet(1.0, [AttenuatorLayer(3, 0.05)]) + att = att & AttenuatorLayer(4, 0.1, is_atom_density=False) + assert att.layers == [ + AttenuatorLayer(3, 0.05), + AttenuatorLayer(4, 0.1, is_atom_density=False), + ] + + def test_and_chains_attenuator_sets(self): + left = AttenuatorSet(1.0, [AttenuatorLayer(3, 0.05)]) + right = AttenuatorSet(1.0, [AttenuatorLayer(4, 0.1)]) + combined = left & right + assert combined.layers == [AttenuatorLayer(3, 0.05), AttenuatorLayer(4, 0.1)] + + def test_single_layer_from_fixture(self, tally_problem): + fm = tally_problem.tallies[114].multiplier + assert fm.bins[0].attenuator == AttenuatorSet(1.0, [AttenuatorLayer(26, 0.5)]) + + def test_multi_layer_from_fixture_normalizes_sign(self, tally_problem): + fm = tally_problem.tallies[124].multiplier + attenuator = fm.bins[0].attenuator + assert attenuator.layers[0] == AttenuatorLayer(26, 0.5, is_atom_density=True) + assert attenuator.layers[1] == AttenuatorLayer(27, 0.3, is_atom_density=False) + # sign is normalized away -- areal_density is always positive + assert attenuator.layers[1].areal_density == 0.3 + + +class TestCompanionCardLinking: + def test_multiplier_is_linked(self, tally_problem): + assert tally_problem.tallies[4].multiplier is not None + assert ( + tally_problem.tallies[4].multiplier.parent_tally is tally_problem.tallies[4] + ) + + def test_multiplier_is_none_without_fm(self, tally_problem): + assert tally_problem.tallies[1].multiplier is None + assert tally_problem.tallies[2].multiplier is None + assert tally_problem.tallies[6].multiplier is None + + def test_orphaned_fm_raises(self): + problem = montepy.MCNP_Problem(None) + fm = TallyMultiplier(Input(["fm999 (1.0)"], BlockType.DATA), jit_parse=False) + problem.tallies.append(fm) + with pytest.raises(montepy.exceptions.MalformedInputError): + problem.tallies.finalize_init() + + def test_appended_fm_registers_in_data_inputs(self): + problem = montepy.MCNP_Problem(None) + tally = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + problem.tallies.append(tally) + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + problem.tallies.append(fm) + assert fm in problem.data_inputs + assert fm in problem.tallies.multipliers + assert tally.multiplier is fm + assert fm.parent_tally is tally + + def test_multiplier_setter_registers_in_data_inputs(self): + problem = montepy.MCNP_Problem(None) + tally = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + problem.tallies.append(tally) + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + tally.multiplier = fm + assert fm in problem.data_inputs + assert fm.parent_tally is tally + + def test_appended_fm_survives_full_problem_export(self): + problem = montepy.MCNP_Problem(None) + problem.title = "test problem" + tally = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + problem.tallies.append(tally) + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + problem.tallies.append(fm) + with io.StringIO() as fh: + problem.write_problem(fh) + fh.seek(0) + new_problem = montepy.read_input(fh) + new_fm = new_problem.tallies[4].multiplier + assert new_fm is not None + assert new_fm.bins[0].terms[0].material == 26 + + +class TestScoresIntegration: + def test_scores_reflects_multiplier_as_first_touch(self, tally_problem): + # This must be the *first* @needs_full_ast touch on this object -- + # regression test for the _KEYS_TO_PRESERVE bug where the FM link + # vanished the moment full_parse() ran. + f4 = tally_problem.tallies[4] + assert not f4.fully_parsed + scores = f4.scores + assert scores == [MultiplierScore(1.0, None, None, None, None)] + + def test_scores_falls_back_to_default_without_multiplier(self, tally_problem): + assert tally_problem.tallies[2].scores == [Score.FLUX] + + def test_scores_for_multi_bin_multiplier(self, tally_problem): + f64 = tally_problem.tallies[64] + scores = f64.scores + assert scores == [ + MultiplierScore(1.0, 26, Reaction(16), None, None), + MultiplierScore(2.0, 27, Reaction(102), None, None), + ] + + def test_scores_for_special_multiplier(self, tally_problem): + f154 = tally_problem.tallies[154] + assert f154.scores == [ + MultiplierScore(1.0, None, None, SpecialMultiplier.INVERSE_WEIGHT, None) + ] + + def test_scores_carry_attenuator(self, tally_problem): + f134 = tally_problem.tallies[134] + attenuator = f134.multiplier.bins[0].attenuator + assert f134.scores == [ + MultiplierScore(1.0, 26, Reaction(16), None, attenuator), + MultiplierScore(2.0, 27, Reaction(102), None, attenuator), + ] + + def test_clone_as_does_not_carry_multiplier(self, tally_problem): + from montepy.data_inputs.tally import F6Tally + + f4 = tally_problem.tallies[4] + assert f4.multiplier is not None + new = f4.clone_as(F6Tally) + assert new.multiplier is None + assert new.scores == [Score.ENERGY_DEPOSITION] + + def test_clone_does_not_carry_multiplier(self, tally_problem): + f4 = tally_problem.tallies[4] + clone = f4.clone() + assert clone.multiplier is None + + +class TestFlags: + def test_include_total_and_cumulative_flags(self, tally_problem): + assert tally_problem.tallies[184].multiplier.include_total is True + assert tally_problem.tallies[184].multiplier.cumulative is False + assert tally_problem.tallies[194].multiplier.include_total is False + assert tally_problem.tallies[194].multiplier.cumulative is True + + def test_include_total_settable(self): + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + assert not fm.include_total + fm.include_total = True + assert fm.include_total + assert not fm.cumulative + assert fm.mcnp_str().strip().endswith("T") + + def test_cumulative_settable(self): + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + fm.cumulative = True + assert fm.cumulative + assert not fm.include_total + assert fm.mcnp_str().strip().endswith("C") + + def test_include_total_and_cumulative_are_mutually_exclusive(self): + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + fm.include_total = True + fm.cumulative = True + assert fm.cumulative + assert not fm.include_total + fm.include_total = True + assert fm.include_total + assert not fm.cumulative + + def test_attenuator_only_bin_scores(self, tally_problem): + fm = tally_problem.tallies[114].multiplier + attenuator = fm.bins[0].attenuator + assert tally_problem.tallies[114].scores == [ + MultiplierScore(1.0, None, None, None, attenuator) + ] + + +class TestParseTermEdgeCases: + def test_material_with_no_reaction_list(self): + fm = TallyMultiplier(Input(["fm999 (1.0 26)"], BlockType.DATA), jit_parse=False) + assert fm.bins[0].terms[0] == MultiplierSet(1.0, 26, []) + + def test_coerce_type_error(self): + with pytest.raises(TypeError): + Reaction(16) + "not a reaction" + + +class TestDuplicateFmCards: + def test_duplicate_fm_cards_warn(self): + problem = montepy.MCNP_Problem(None) + tally = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + problem.tallies.append(tally) + fm1 = TallyMultiplier(Input(["fm4 (1.0)"], BlockType.DATA), jit_parse=False) + fm2 = TallyMultiplier(Input(["fm4 (2.0)"], BlockType.DATA), jit_parse=False) + problem.tallies.append(fm1) + with pytest.warns(montepy.exceptions.MalformedInputWarning): + problem.tallies.append(fm2) + + def test_deleting_tally_cascades_to_multiplier(self): + problem = montepy.MCNP_Problem(None) + tally = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + problem.tallies.append(tally) + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + problem.tallies.append(fm) + with pytest.warns(montepy.exceptions.MalformedInputWarning): + del problem.tallies[4] + assert fm not in problem.data_inputs + assert fm not in problem.tallies.multipliers + assert fm.parent_tally is None + + def test_deleting_tally_without_multiplier_does_not_warn(self): + problem = montepy.MCNP_Problem(None) + tally = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + problem.tallies.append(tally) + with warnings.catch_warnings(): + warnings.simplefilter("error") + del problem.tallies[4] + assert tally not in problem.data_inputs + + def test_renumbering_tally_syncs_multiplier_number(self): + problem = montepy.MCNP_Problem(None) + tally = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + problem.tallies.append(tally) + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + problem.tallies.append(fm) + tally.number = 14 + assert fm.number == 14 + + def test_deleting_tally_removes_fm_from_full_problem_export(self): + # Regression test for the collection/problem-level cascade: a + # per-object mcnp_str() check can't see whether the FM actually made + # it out of the *problem's* data_inputs on a full write. + problem = montepy.MCNP_Problem(None) + problem.title = "test problem" + tally = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + problem.tallies.append(tally) + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + problem.tallies.append(fm) + with pytest.warns(montepy.exceptions.MalformedInputWarning): + del problem.tallies[4] + with io.StringIO() as fh: + problem.write_problem(fh) + fh.seek(0) + written = fh.read() + fh.seek(0) + new_problem = montepy.read_input(fh) + assert "fm4" not in written.lower() + assert 4 not in new_problem.tallies.numbers + assert new_problem.tallies.multipliers == [] + + def test_renumbering_tally_syncs_multiplier_across_full_problem_export(self): + problem = montepy.MCNP_Problem(None) + problem.title = "test problem" + tally = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + problem.tallies.append(tally) + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + problem.tallies.append(fm) + tally.number = 14 + # 4 -> 14 widens the field; MCNP's historically column-based format + # warns on that regardless of this plan's changes (see test_integration.py). + with pytest.warns(montepy.exceptions.LineExpansionWarning): + new_fm = verify_prob_export(problem, fm) + assert new_fm is not None + assert new_fm.number == 14 + assert new_fm.parent_tally.number == 14 + + +class TestClone: + def test_clone_to_new_tally_registers_and_links(self): + problem = montepy.MCNP_Problem(None) + tally4 = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + tally14 = parse_data(Input(["f14:n 4 5 6"], BlockType.DATA)) + problem.tallies.append(tally4) + problem.tallies.append(tally14) + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + problem.tallies.append(fm) + + new_fm = fm.number + clone = fm.clone(tally14) + + assert clone is not fm + assert clone.number == 14 + assert clone.parent_tally is tally14 + assert tally14.multiplier is clone + assert clone in problem.data_inputs + assert clone in problem.tallies.multipliers + # the original is untouched + assert fm.number == new_fm + assert fm.parent_tally.number == 4 + + def test_clone_without_tally_is_detached(self): + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + clone = fm.clone() + assert clone is not fm + assert clone.parent_tally is None + assert clone.bins[0].terms[0].material == 26 + + def test_clone_to_tally_with_existing_multiplier_warns(self): + problem = montepy.MCNP_Problem(None) + tally4 = parse_data(Input(["f4:n 1 2 3"], BlockType.DATA)) + tally14 = parse_data(Input(["f14:n 4 5 6"], BlockType.DATA)) + problem.tallies.append(tally4) + problem.tallies.append(tally14) + fm4 = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + fm14 = TallyMultiplier( + Input(["fm14 (2.0 27 102)"], BlockType.DATA), jit_parse=False + ) + problem.tallies.append(fm4) + problem.tallies.append(fm14) + with pytest.warns(montepy.exceptions.MalformedInputWarning): + fm4.clone(tally14) + + +class TestBlankConstruction: + def test_blank_tally_multiplier_construction(self): + fm = TallyMultiplier() + assert fm.bins == [] + assert str(fm) + assert repr(fm) + + def test_str_and_repr_on_uninitialized_tally_multiplier(self): + # __new__ bypasses __init__ entirely -- exactly the partially + # constructed state __str__/__repr__'s except branches exist to + # report gracefully. + fm = TallyMultiplier.__new__(TallyMultiplier) + assert str(fm) == "TallyMultiplier: (unparsed)" + assert repr(fm) == "TALLY MULTIPLIER: (unparsed)" + + +class TestReprAndEquality: + """Smoke tests for __repr__, __eq__-against-wrong-type, and otherwise + unexercised property getters across the tally_multiplier.py value + objects. Coverage only counts a line as hit if it executes, and pytest + only reprs on assertion *failure*, so these branches need explicit + exercising; likewise the existing tests only ever compare MultiplierScore + via == (which reads the private attributes directly), never through its + public properties.""" + + def test_reaction_expression_eq_wrong_type_and_repr(self): + expr = Reaction(16) * Reaction(103) + assert expr != "not an expression" + assert "ReactionExpression" in repr(expr) + + def test_reaction_number_and_eq_wrong_type_and_repr(self): + r = Reaction(16) + assert r.number == 16 + assert r != "not a reaction" + assert repr(r) == "Reaction(16)" + + def test_attenuator_layer_material_and_eq_wrong_type_and_repr(self): + layer = AttenuatorLayer(26, 0.5) + assert layer.material == 26 + assert layer.is_atom_density is True + assert layer != "not a layer" + assert "AttenuatorLayer" in repr(layer) + + def test_attenuator_set_constant_and_eq_wrong_type_and_repr(self): + att = AttenuatorSet(1.0, [AttenuatorLayer(26, 0.5)]) + assert att.constant == 1.0 + assert att != "not an attenuator set" + assert "AttenuatorSet" in repr(att) + + def test_multiplier_set_eq_wrong_type_and_repr(self): + ms = MultiplierSet(1.0, 26, [Reaction(16)]) + assert ms != "not a multiplier set" + assert "MultiplierSet" in repr(ms) + + def test_special_multiplier_set_eq_wrong_type_and_repr(self): + sms = SpecialMultiplierSet(1.0, SpecialMultiplier.INVERSE_WEIGHT) + assert sms == SpecialMultiplierSet(1.0, SpecialMultiplier.INVERSE_WEIGHT) + assert sms != "not a special multiplier set" + assert "SpecialMultiplierSet" in repr(sms) + + def test_multiplier_score_properties_eq_wrong_type_and_repr(self): + attenuator = AttenuatorSet(1.0, [AttenuatorLayer(26, 0.5)]) + score = MultiplierScore( + 1.0, 26, Reaction(16), SpecialMultiplier.INVERSE_WEIGHT, attenuator + ) + assert score.constant == 1.0 + assert score.material == 26 + assert score.reaction == Reaction(16) + assert score.kind == SpecialMultiplier.INVERSE_WEIGHT + assert score.attenuator == attenuator + assert score != "not a score" + assert "MultiplierScore" in repr(score) + + def test_multiplier_bin_eq_wrong_type_and_repr(self): + mb = MultiplierBin([MultiplierSet(1.0, 26, [Reaction(16)])]) + assert mb == MultiplierBin([MultiplierSet(1.0, 26, [Reaction(16)])]) + assert mb != "not a bin" + assert "MultiplierBin" in repr(mb) + + def test_parent_collections(self): + assert TallyMultiplier._parent_collections() == () + + +class TestBinRoundTrip: + """Mutating a TallyMultiplier's bins through the public API must be + reflected in mcnp_str(), not just in the in-memory Python state.""" + + @pytest.mark.parametrize("line", FM_FIXTURE_LINES) + def test_unmodified_fm_round_trips_exactly(self, line): + fm = TallyMultiplier(Input([line], BlockType.DATA), jit_parse=False) + assert fm.mcnp_str() == line + fm.full_parse() + assert fm.mcnp_str() == line + verify_export(fm) + + def test_add_bin_reflected_in_mcnp_str(self): + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + fm.add_bin(MultiplierBin([MultiplierSet(2.0, 27, [Reaction(102)])])) + text = fm.mcnp_str() + assert "27" in text and "102" in text and "2.0" in text + # two bins now -- MCNP requires each parenthesized separately. + assert text.count("(") == 2 and text.count(")") == 2 + verify_export(fm) + + def test_add_single_bin_from_scratch_has_no_extra_parens(self): + fm = TallyMultiplier() + fm.number = 4 + fm.add_bin(MultiplierBin([MultiplierSet(1.0, 26, [Reaction(16)])])) + text = fm.mcnp_str() + assert "(" not in text and ")" not in text + verify_export(fm) + + def test_add_attenuator_only_bin(self): + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + att_bin = MultiplierBin( + [], attenuator=AttenuatorSet(1.0, [AttenuatorLayer(28, 0.2)]) + ) + fm.add_bin(att_bin) + text = fm.mcnp_str() + assert "28" in text and "0.2" in text + new_fm = verify_export(fm) + assert new_fm.bins[-1].attenuator == att_bin.attenuator + assert new_fm.bins[-1].terms == [] + + def test_blank_fm_add_bin_writes_valid_card(self): + fm = TallyMultiplier() + fm.number = 4 + fm.add_bin(MultiplierBin([MultiplierSet(1.0, 26, [Reaction(16)])])) + text = fm.mcnp_str() + assert "26" in text and "16" in text + verify_export(fm) + + def test_remove_bin(self): + fm = TallyMultiplier( + Input(["fm4 (1.0 26 16)"], BlockType.DATA), jit_parse=False + ) + bin_ = MultiplierBin([MultiplierSet(2.0, 27, [Reaction(102)])]) + fm.add_bin(bin_) + fm.remove_bin(bin_) + assert bin_ not in fm.bins + text = fm.mcnp_str() + assert "27" not in text + verify_export(fm) + + def test_multiplier_set_material_object_resolves_live(self): + mat = montepy.Material() + mat.number = 26 + built = mat & Reaction.N_2N + fm = TallyMultiplier() + fm.number = 4 + fm.add_bin(MultiplierBin([built])) + assert "26" in fm.mcnp_str() + mat.number = 99 + text = fm.mcnp_str() + assert "99" in text + assert "26" not in text + assert built.material == 99