Merge pull request #1536 from paulromano/model-run-fix

Fix for Model.run (cleaned up version of #1498)

openmc/executor.py

@@ -1,6 +1,6 @@
 from collections.abc import Iterable
-import subprocess
 from numbers import Integral
+import subprocess
 
 import openmc
 

openmc/model/model.py

@@ -1,5 +1,6 @@
 from collections.abc import Iterable
 from pathlib import Path
+import time
 
 import openmc
 from openmc.checkvalue import check_type, check_value
@@ -172,7 +173,7 @@ class Model:
             will be created.
 
         """
-        # Create directory if
+        # Create directory if required
         d = Path(directory)
         if not d.is_dir():
             d.mkdir(parents=True)
@@ -197,27 +198,39 @@ class Model:
             self.plots.export_to_xml(d)
 
     def run(self, **kwargs):
-        """Creates the XML files, runs OpenMC, and returns k-effective
+        """Creates the XML files, runs OpenMC, and returns the path to the last
+        statepoint file generated.
 
         Parameters
         ----------
         **kwargs
-            All keyword arguments are passed to :func:`openmc.run`
+            Keyword arguments passed to :func:`openmc.run`
 
         Returns
         -------
-        uncertainties.UFloat
-            Combined estimator of k-effective from the statepoint
+        Path
+            Path to the last statepoint written by this run
+            (None if no statepoint was written)
 
         """
+
         self.export_to_xml()
 
+        # Setting tstart here ensures we don't pick up any pre-existing statepoint
+        # files in the output directory
+        tstart = time.time()
+        last_statepoint = None
+
         openmc.run(**kwargs)
 
-        n = self.settings.batches
-        if self.settings.statepoint is not None:
-            if 'batches' in self.settings.statepoint:
-                n = self.settings.statepoint['batches'][-1]
-
-        with openmc.StatePoint('statepoint.{}.h5'.format(n)) as sp:
-            return sp.k_combined
+        # Get output directory and return the last statepoint written by this run
+        if self.settings.output and 'path' in self.settings.output:
+            output_dir = Path(self.settings.output['path'])
+        else:
+            output_dir = Path.cwd()
+        for sp in output_dir.glob('statepoint.*.h5'):
+            mtime = sp.stat().st_mtime
+            if mtime >= tstart:  # >= allows for poor clock resolution
+                tstart = mtime
+                last_statepoint = sp
+        return last_statepoint

openmc/search.py

@@ -51,7 +51,9 @@ def _search_keff(guess, target, model_builder, model_args, print_iterations,
     model = model_builder(guess, **model_args)
 
     # Run the model and obtain keff
-    keff = model.run(output=print_output)
+    sp_filepath = model.run(output=print_output)
+    with openmc.StatePoint(sp_filepath) as sp:
+        keff = sp.k_combined
 
     # Record the history
     guesses.append(guess)

openmc/settings.py

@@ -813,10 +813,9 @@ class Settings:
     def _create_keff_trigger_subelement(self, root):
         if self._keff_trigger is not None:
             element = ET.SubElement(root, "keff_trigger")
-
-            for key in self._keff_trigger:
+            for key, value in sorted(self._keff_trigger.items()):
                 subelement = ET.SubElement(element, key)
-                subelement.text = str(self._keff_trigger[key]).lower()
+                subelement.text = str(value).lower()
 
     def _create_energy_mode_subelement(self, root):
         if self._energy_mode is not None:
@@ -839,8 +838,7 @@ class Settings:
     def _create_output_subelement(self, root):
         if self._output is not None:
             element = ET.SubElement(root, "output")
-
-            for key, value in self._output.items():
+            for key, value in sorted(self._output.items()):
                 subelement = ET.SubElement(element, key)
                 if key in ('summary', 'tallies'):
                     subelement.text = str(value).lower()

openmc/statepoint.py

@@ -21,7 +21,7 @@ class StatePoint:
 
     Parameters
     ----------
-    filename : str
+    filepath : str or Path
         Path to file to load
     autolink : bool, optional
         Whether to automatically link in metadata from a summary.h5 file and
@@ -115,7 +115,8 @@ class StatePoint:
 
     """
 
-    def __init__(self, filename, autolink=True):
+    def __init__(self, filepath, autolink=True):
+        filename = str(filepath)  # in case it's a Path
         self._f = h5py.File(filename, 'r')
         self._meshes = {}
         self._filters = {}

tests/regression_tests/trigger_statepoint_restart/inputs_true.dat

@@ -0,0 +1,35 @@
+<?xml version='1.0' encoding='utf-8'?>
+<geometry>
+  <cell id="1" material="1" region="-1" universe="1" />
+  <surface boundary="vacuum" coeffs="0.0 0.0 0.0 10.0" id="1" type="sphere" />
+</geometry>
+<?xml version='1.0' encoding='utf-8'?>
+<materials>
+  <material depletable="true" id="1">
+    <density units="g/cm3" value="4.5" />
+    <nuclide ao="1.0" name="U235" />
+  </material>
+</materials>
+<?xml version='1.0' encoding='utf-8'?>
+<settings>
+  <run_mode>eigenvalue</run_mode>
+  <particles>200</particles>
+  <batches>10</batches>
+  <inactive>5</inactive>
+  <keff_trigger>
+    <threshold>0.004</threshold>
+    <type>std_dev</type>
+  </keff_trigger>
+  <trigger>
+    <active>true</active>
+    <max_batches>1000</max_batches>
+    <batch_interval>1</batch_interval>
+  </trigger>
+  <verbosity>1</verbosity>
+</settings>
+<?xml version='1.0' encoding='utf-8'?>
+<tallies>
+  <tally id="1">
+    <scores>flux</scores>
+  </tally>
+</tallies>

tests/regression_tests/trigger_statepoint_restart/results_true.dat

@@ -0,0 +1,5 @@
+k-combined:
+2.916922E-01 3.293799E-03
+tally 1:
+6.184423E+01
+4.789617E+02

tests/regression_tests/trigger_statepoint_restart/test.py

@@ -0,0 +1,132 @@
+import glob
+import os
+
+import openmc
+import pytest
+
+from tests.testing_harness import PyAPITestHarness
+from tests.regression_tests import config
+
+
+@pytest.fixture
+def model():
+
+    # Materials
+    mat = openmc.Material()
+    mat.set_density('g/cm3', 4.5)
+    mat.add_nuclide('U235', 1.0)
+    materials = openmc.Materials([mat])
+
+    # Geometry
+    sph = openmc.Sphere(r=10.0, boundary_type='vacuum')
+    cell = openmc.Cell(fill=mat, region=-sph)
+    geometry = openmc.Geometry([cell])
+
+    # Settings
+    settings = openmc.Settings()
+    settings.run_mode = 'eigenvalue'
+    settings.batches = 10
+    settings.inactive = 5
+    settings.particles = 200
+    # Choose a sufficiently low threshold to trigger after more than 10 batches.
+    # 0.004 seems to take 13 batches.
+    settings.keff_trigger = {'type': 'std_dev', 'threshold': 0.004}
+    settings.trigger_max_batches = 1000
+    settings.trigger_batch_interval = 1
+    settings.trigger_active = True
+    settings.verbosity = 1 # to test that this works even with no output
+ 
+    # Tallies
+    t = openmc.Tally()
+    t.scores = ['flux']
+    tallies = openmc.Tallies([t])
+ 
+    # Put it all together
+    model = openmc.model.Model(materials=materials,
+                               geometry=geometry,
+                               settings=settings,
+                               tallies=tallies)
+    return model
+
+
+class TriggerStatepointRestartTestHarness(PyAPITestHarness):
+    def __init__(self, statepoint, model=None):
+        super().__init__(statepoint, model)
+        self._restart_sp = None
+        self._final_sp = None
+        # store the statepoint filename pattern separately to sp_name so we can reuse it
+        self._sp_pattern = self._sp_name
+
+    def _test_output_created(self):
+        """Make sure statepoint files have been created."""
+        spfiles = sorted(glob.glob(self._sp_pattern))
+        assert len(spfiles) == 2, \
+            'Two statepoint files should have been created'
+        if not self._final_sp:
+            # First non-restart run
+            self._restart_sp = spfiles[0]
+            self._final_sp = spfiles[1]
+        else:
+            # Second restart run
+            assert spfiles[1] == self._final_sp, \
+                'Final statepoint names were different'
+        # Use the final_sp as the sp_name for the 'standard' results tests
+        self._sp_name = self._final_sp
+
+    def execute_test(self):
+        """
+        Perform initial and restart runs using the model.run method,
+        Check all inputs and outputs which should be the same as those
+        generated using the normal PyAPITestHarness update methods.
+        """
+        try:
+            args = {'openmc_exec': config['exe'], 'event_based': config['event']}
+            if config['mpi']:
+                args['mpi_args'] = [config['mpiexec'], '-n', config['mpi_np']]
+
+            # First non-restart run
+            spfile = self._model.run(**args)
+            sp_batchno_1 = 0
+            print('Last sp file: %s' % spfile)
+            assert spfile
+            with openmc.StatePoint(spfile) as sp:
+                 sp_batchno_1 = sp.current_batch
+                 k_combined_1 = sp.k_combined
+            assert sp_batchno_1 > 10
+            print('Last batch no = %d' % sp_batchno_1)
+            self._write_inputs(self._get_inputs())
+            self._compare_inputs()
+            self._test_output_created()
+            self._write_results(self._get_results())
+            self._compare_results()
+
+            # Second restart run
+            restart_spfile = glob.glob(os.path.join(os.getcwd(), self._restart_sp))
+            assert len(restart_spfile) == 1
+            args['restart_file'] = restart_spfile[0]
+            spfile = self._model.run(**args)
+            sp_batchno_2 = 0
+            assert spfile
+            with openmc.StatePoint(spfile) as sp:
+                 sp_batchno_2 = sp.current_batch
+                 k_combined_2 = sp.k_combined
+            assert sp_batchno_2 > 10
+            assert sp_batchno_1 == sp_batchno_2, \
+                'Different final batch number after restart'
+            # need str() here as uncertainties.ufloat instances are always different
+            assert str(k_combined_1) == str(k_combined_2), \
+                'Different final k_combined after restart'
+            self._write_inputs(self._get_inputs())
+            self._compare_inputs()
+            self._test_output_created()
+            self._write_results(self._get_results())
+            self._compare_results()
+        finally:
+            self._cleanup()
+
+
+def test_trigger_statepoint_restart(model):
+    # Assuming we converge within 1000 batches, the statepoint filename
+    # should include the batch number padded by at least one '0'.
+    harness = TriggerStatepointRestartTestHarness('statepoint.0*.h5', model)
+    harness.main()

tests/unit_tests/test_filters.py

@@ -173,10 +173,10 @@ def test_first_moment(run_in_tmpdir, box_model):
     for t in box_model.tallies:
         t.estimator = 'analog'
 
-    box_model.run()
+    sp_name = box_model.run()
 
     # Check that first moment matches the score from the plain tally
-    with openmc.StatePoint('statepoint.10.h5') as sp:
+    with openmc.StatePoint(sp_name) as sp:
         # Get scores from tally without expansion filters
         flux, scatter = sp.tallies[plain_tally.id].mean.ravel()