Repository: openmc-dev/openmc

Issue:
`Universe.plot` seed does not change color behavior.
## Bug Description
I noticed today that setting the `seed` parameter when calling `Universe.plot` doesn't affect the colors in a plot.

## Steps to Reproduce
Run the following in a notebook for the assembly model located at `examples/assembly` in the repo.

<img width="1067" alt="Screenshot 2023-08-13 at 2 57 47 PM" src="https://github.com/openmc-dev/openmc/assets/4563941/06afb0a2-7d19-43e2-8d3b-4c6470a15c1c">

Discussion:

Patch:
diff --git a/docs/source/io_formats/settings.rst b/docs/source/io_formats/settings.rst
index 1016fbaa3..36027ef2f 100644
--- a/docs/source/io_formats/settings.rst
+++ b/docs/source/io_formats/settings.rst
@@ -345,6 +345,15 @@ either "false" or "true".
 
   *Default*: false
 
+-----------------------
+``<plot_seed>`` Element
+-----------------------
+
+The ``<plot_seed>`` element is used to set the seed for the pseudorandom number
+generator during generation of colors in plots.
+
+  *Default*: 1
+
 ---------------------
 ``<ptables>`` Element
 ---------------------
@@ -1175,4 +1184,3 @@ mesh-based weight windows.
 
   The ``weight_windows_file`` element has no attributes and contains the path to
   a weight windows HDF5 file to load during simulation initialization.
-
diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index e781e154f..d205bf571 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -512,10 +512,11 @@ class Results(list):
         if math.isclose(time, times[ix], rel_tol=rtol, abs_tol=atol):
             return ix
 
+        closest = min(times, key=lambda t: abs(time - t))
         raise ValueError(
-            "A value of {} {} was not found given absolute and "
-            "relative tolerances {} and {}.".format(
-                time, time_units, atol, rtol)
+            f"A value of {time} {time_units} was not found given absolute and "
+            f"relative tolerances {atol} and {rtol}. Closest time is {closest} "
+            f"{time_units}."
         )
 
     def export_to_materials(
diff --git a/openmc/material.py b/openmc/material.py
index 95dae032e..4ef16ab01 100644
--- a/openmc/material.py
+++ b/openmc/material.py
@@ -642,7 +642,8 @@ class Material(IDManagerMixin):
     def add_element(self, element: str, percent: float, percent_type: str = 'ao',
                     enrichment: Optional[float] = None,
                     enrichment_target: Optional[str] = None,
-                    enrichment_type: Optional[str] = None):
+                    enrichment_type: Optional[str] = None,
+                    cross_sections: Optional[str] = None):
         """Add a natural element to the material
 
         Parameters
@@ -669,6 +670,8 @@ class Material(IDManagerMixin):
             Default is: 'ao' for two-isotope enrichment; 'wo' for U enrichment
 
             .. versionadded:: 0.12
+        cross_sections : str, optional
+            Location of cross_sections.xml file.
 
         Notes
         -----
@@ -747,7 +750,8 @@ class Material(IDManagerMixin):
                                       percent_type,
                                       enrichment,
                                       enrichment_target,
-                                      enrichment_type):
+                                      enrichment_type,
+                                      cross_sections):
             self.add_nuclide(*nuclide)
 
     def add_elements_from_formula(self, formula: str, percent_type: str = 'ao',
@@ -1457,7 +1461,7 @@ class Material(IDManagerMixin):
         if "cfg" in elem.attrib:
             cfg = elem.get("cfg")
             return Material.from_ncrystal(cfg, material_id=mat_id)
-     
+
         mat = cls(mat_id)
         mat.name = elem.get('name')
 
diff --git a/openmc/settings.py b/openmc/settings.py
index 7937cd2e7..c4be8699f 100644
--- a/openmc/settings.py
+++ b/openmc/settings.py
@@ -133,6 +133,8 @@ class Settings:
         Number of particles per generation
     photon_transport : bool
         Whether to use photon transport.
+    plot_seed : int
+       Initial seed for randomly generated plot colors.
     ptables : bool
         Determine whether probability tables are used.
     resonance_scattering : dict
@@ -270,6 +272,7 @@ class Settings:
         self._confidence_intervals = None
         self._electron_treatment = None
         self._photon_transport = None
+        self._plot_seed = None
         self._ptables = None
         self._seed = None
         self._survival_biasing = None
@@ -504,6 +507,16 @@ class Settings:
         cv.check_type('photon transport', photon_transport, bool)
         self._photon_transport = photon_transport
 
+    @property
+    def plot_seed(self):
+        return self._plot_seed
+
+    @plot_seed.setter
+    def plot_seed(self, seed):
+        cv.check_type('random plot color seed', seed, Integral)
+        cv.check_greater_than('random plot color seed', seed, 0)
+        self._plot_seed = seed
+
     @property
     def seed(self) -> int:
         return self._seed
@@ -1118,6 +1131,11 @@ class Settings:
             element = ET.SubElement(root, "photon_transport")
             element.text = str(self._photon_transport).lower()
 
+    def _create_plot_seed_subelement(self, root):
+        if self._plot_seed is not None:
+            element = ET.SubElement(root, "plot_seed")
+            element.text = str(self._plot_seed)
+
     def _create_ptables_subelement(self, root):
         if self._ptables is not None:
             element = ET.SubElement(root, "ptables")
@@ -1491,6 +1509,11 @@ class Settings:
         if text is not None:
             self.photon_transport = text in ('true', '1')
 
+    def _plot_seed_from_xml_element(self, root):
+        text = get_text(root, 'plot_seed')
+        if text is not None:
+            self.plot_seed = int(text)
+
     def _ptables_from_xml_element(self, root):
         text = get_text(root, 'ptables')
         if text is not None:
@@ -1706,6 +1729,7 @@ class Settings:
         self._create_energy_mode_subelement(element)
         self._create_max_order_subelement(element)
         self._create_photon_transport_subelement(element)
+        self._create_plot_seed_subelement(element)
         self._create_ptables_subelement(element)
         self._create_seed_subelement(element)
         self._create_survival_biasing_subelement(element)
@@ -1802,6 +1826,7 @@ class Settings:
         settings._energy_mode_from_xml_element(elem)
         settings._max_order_from_xml_element(elem)
         settings._photon_transport_from_xml_element(elem)
+        settings._plot_seed_from_xml_element(elem)
         settings._ptables_from_xml_element(elem)
         settings._seed_from_xml_element(elem)
         settings._survival_biasing_from_xml_element(elem)
diff --git a/openmc/universe.py b/openmc/universe.py
index c490ac298..d4743e263 100644
--- a/openmc/universe.py
+++ b/openmc/universe.py
@@ -370,8 +370,8 @@ class Universe(UniverseBase):
 
         Returns
         -------
-        matplotlib.image.AxesImage
-            Resulting image
+        matplotlib.axes.Axes
+            Axes containing resulting image
 
         """
         import matplotlib.image as mpimg
@@ -423,7 +423,7 @@ class Universe(UniverseBase):
             model = openmc.Model()
             model.geometry = openmc.Geometry(self)
             if seed is not None:
-                model.settings.seed = seed
+                model.settings.plot_seed = seed
 
             # Determine whether any materials contains macroscopic data and if
             # so, set energy mode accordingly
@@ -483,7 +483,7 @@ class Universe(UniverseBase):
             # add legend showing which colors represent which material
             # or cell if that was requested
             if legend:
-                if plot.colors is None:
+                if plot.colors == {}:
                     raise ValueError("Must pass 'colors' dictionary if you "
                                      "are adding a legend via legend=True.")
 
@@ -522,7 +522,8 @@ class Universe(UniverseBase):
                 axes.legend(handles=patches, **legend_kwargs)
 
             # Plot image and return the axes
-            return axes.imshow(img, extent=(x_min, x_max, y_min, y_max), **kwargs)
+            axes.imshow(img, extent=(x_min, x_max, y_min, y_max), **kwargs)
+            return axes
 
     def add_cell(self, cell):
         """Add a cell to the universe.
diff --git a/openmc/volume.py b/openmc/volume.py
index 72fe89fbc..dddd30a52 100644
--- a/openmc/volume.py
+++ b/openmc/volume.py
@@ -122,6 +122,10 @@ class VolumeCalculation:
                 raise ValueError('Could not automatically determine bounding box '
                                  'for stochastic volume calculation.')
 
+        if np.isinf(self.lower_left).any() or np.isinf(self.upper_right).any():
+            raise ValueError('Lower-left and upper-right bounding box '
+                             'coordinates must be finite.')
+
     @property
     def ids(self):
         return self._ids
diff --git a/src/settings.cpp b/src/settings.cpp
index 1c039e9ba..ace2887a6 100644
--- a/src/settings.cpp
+++ b/src/settings.cpp
@@ -22,6 +22,7 @@
 #include "openmc/mesh.h"
 #include "openmc/message_passing.h"
 #include "openmc/output.h"
+#include "openmc/plot.h"
 #include "openmc/random_lcg.h"
 #include "openmc/simulation.h"
 #include "openmc/source.h"
@@ -390,6 +391,12 @@ void read_settings_xml(pugi::xml_node root)
     }
   }
 
+  // Copy plotting random number seed if specified
+  if (check_for_node(root, "plot_seed")) {
+    auto seed = std::stoll(get_node_value(root, "plot_seed"));
+    model::plotter_seed = seed;
+  }
+
   // Copy random number seed if specified
   if (check_for_node(root, "seed")) {
     auto seed = std::stoll(get_node_value(root, "seed"));
diff --git a/src/volume_calc.cpp b/src/volume_calc.cpp
index 0b1ea3575..b58da4750 100644
--- a/src/volume_calc.cpp
+++ b/src/volume_calc.cpp
@@ -534,8 +534,17 @@ int openmc_calculate_volumes()
 
       // Display domain volumes
       for (int j = 0; j < vol_calc.domain_ids_.size(); j++) {
-        write_message(4, "{}{}: {} +/- {} cm^3", domain_type,
-          vol_calc.domain_ids_[j], results[j].volume[0], results[j].volume[1]);
+        std::string cell_name {""};
+        if (vol_calc.domain_type_ == VolumeCalculation::TallyDomain::CELL) {
+          int cell_idx = model::cell_map[vol_calc.domain_ids_[j]];
+          cell_name = model::cells[cell_idx]->name();
+          if (cell_name.size())
+            cell_name.insert(0, " "); // prepend space for formatting
+        }
+
+        write_message(4, "{}{}{}: {} +/- {} cm^3", domain_type,
+          vol_calc.domain_ids_[j], cell_name, results[j].volume[0],
+          results[j].volume[1]);
       }
 
       // Write volumes to HDF5 file


Test patch:
diff --git a/tests/unit_tests/test_settings.py b/tests/unit_tests/test_settings.py
index d75a50477..79f2038b2 100644
--- a/tests/unit_tests/test_settings.py
+++ b/tests/unit_tests/test_settings.py
@@ -23,6 +23,7 @@ def test_export_to_xml(run_in_tmpdir):
     s.surf_source_write = {'surface_ids': [2], 'max_particles': 200}
     s.confidence_intervals = True
     s.ptables = True
+    s.plot_seed = 100
     s.survival_biasing = True
     s.cutoff = {'weight': 0.25, 'weight_avg': 0.5, 'energy_neutron': 1.0e-5,
                 'energy_photon': 1000.0, 'energy_electron': 1.0e-5,
@@ -82,6 +83,7 @@ def test_export_to_xml(run_in_tmpdir):
     assert s.surf_source_write == {'surface_ids': [2], 'max_particles': 200}
     assert s.confidence_intervals
     assert s.ptables
+    assert s.plot_seed == 100
     assert s.seed == 17
     assert s.survival_biasing
     assert s.cutoff == {'weight': 0.25, 'weight_avg': 0.5,
@@ -108,7 +110,7 @@ def test_export_to_xml(run_in_tmpdir):
                                       'energy_min': 1.0, 'energy_max': 1000.0,
                                       'nuclides': ['U235', 'U238', 'Pu239']}
     assert s.create_fission_neutrons
-    assert not s.create_delayed_neutrons 
+    assert not s.create_delayed_neutrons
     assert s.log_grid_bins == 2000
     assert not s.photon_transport
     assert s.electron_treatment == 'led'
diff --git a/tests/unit_tests/test_universe.py b/tests/unit_tests/test_universe.py
index f4cf43d91..630e66df3 100644
--- a/tests/unit_tests/test_universe.py
+++ b/tests/unit_tests/test_universe.py
@@ -61,14 +61,17 @@ def test_plot(run_in_tmpdir, sphere_model):
     }
 
     for basis in ('xy', 'yz', 'xz'):
-        pincell.geometry.root_universe.plot(
+        plot = pincell.geometry.root_universe.plot(
             colors=mat_colors,
             color_by="material",
             legend=True,
             pixels=(10, 10),
             basis=basis,
-            outline=True
+            outline=True,
+            axis_units='m'
         )
+        assert plot.xaxis.get_label().get_text() == f'{basis[0]} [m]'
+        assert plot.yaxis.get_label().get_text() == f'{basis[1]} [m]'
 
     # model with no inf values in bounding box
     m = sphere_model.materials[0]
@@ -77,7 +80,7 @@ def test_plot(run_in_tmpdir, sphere_model):
     colors = {m: 'limegreen'}
 
     for basis in ('xy', 'yz', 'xz'):
-        univ.plot(
+        plot = univ.plot(
             colors=colors,
             color_by="cell",
             legend=False,
@@ -85,6 +88,17 @@ def test_plot(run_in_tmpdir, sphere_model):
             basis=basis,
             outline=False
         )
+        assert plot.xaxis.get_label().get_text() == f'{basis[0]} [cm]'
+        assert plot.yaxis.get_label().get_text() == f'{basis[1]} [cm]'
+
+    msg = "Must pass 'colors' dictionary if you are adding a legend via legend=True."
+    # This plot call should fail as legend is True but colors is None
+    with pytest.raises(ValueError, match=msg):
+        univ.plot(
+            color_by="cell",
+            legend=True,
+            pixels=100,
+        )
 
 
 def test_get_nuclides(uo2):
diff --git a/tests/unit_tests/test_volume.py b/tests/unit_tests/test_volume.py
new file mode 100644
index 000000000..f49bb3501
--- /dev/null
+++ b/tests/unit_tests/test_volume.py
@@ -0,0 +1,15 @@
+import numpy as np
+import pytest
+
+import openmc
+
+
+def test_infinity_handling():
+    surf1 = openmc.Sphere(boundary_type="vacuum")
+    cell1 = openmc.Cell(region=-surf1)
+
+    lower_left = (-2, -np.inf, -2)
+    upper_right = (np.inf, 2, 2)
+
+    with pytest.raises(ValueError, match="must be finite"):
+        openmc.VolumeCalculation([cell1], 100, lower_left, upper_right)