From d55c77a65fc11aa78c327b879e224b054e649caf Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 09:37:07 -0500
Subject: [PATCH 01/19] Pass depletion process time to deplete.Results.save

Method openmc.deplete.Results.save now accepts a new argument,
proc_time, to represent the time each process spent depleting
all materials. This value is summed across all processes to
get the total cpu time depleting. Time in transport simulations
should not be included.

The value written to depletion_results.h5 is
proc_time / (n_proc * n_burn_materials)
to get a look at the average time required to deplete a material.
A new dataset, "depletion time" is created in this file, and has
one fewer rows than other data set, like eigenvalues. This corresponds
to there being one fewer round of depletion simulations than transport.

The proc_time is also read in from the depletion file in the
class method Results.from_hdf5(), with one notable exception.
If data from the last step is requested, then the proc_time
attribute will be a numpy array of nans, reflecting the lack
of depletion data.
---
 openmc/deplete/__init__.py |  2 ++
 openmc/deplete/results.py  | 37 +++++++++++++++++++++++++++++++++----
 2 files changed, 35 insertions(+), 4 deletions(-)

diff --git a/openmc/deplete/__init__.py b/openmc/deplete/__init__.py
index e49c4e69cf..32588588ef 100644
--- a/openmc/deplete/__init__.py
+++ b/openmc/deplete/__init__.py
@@ -10,9 +10,11 @@ try:
     from mpi4py import MPI
     comm = MPI.COMM_WORLD
     have_mpi = True
+    mpi_sum = MPI.SUM
 except ImportError:
     comm = DummyCommunicator()
     have_mpi = False
+    mpi_sum = lambda x: x
 
 from .nuclide import *
 from .chain import *
diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index 6170a44643..e37ccb5407 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -10,7 +10,7 @@ from warnings import warn
 import numpy as np
 import h5py
 
-from . import comm, have_mpi
+from . import comm, have_mpi, mpi_sum
 from .reaction_rates import ReactionRates
 
 _VERSION_RESULTS = (1, 0)
@@ -55,6 +55,7 @@ class Results(object):
         self.power = None
         self.rates = None
         self.volume = None
+        self.proc_time = None
 
         self.mat_to_ind = None
         self.nuc_to_ind = None
@@ -245,6 +246,10 @@ class Results(object):
         handle.create_dataset("power", (1, n_stages), maxshape=(None, n_stages),
                               dtype='float64')
 
+        handle.create_dataset(
+            "depletion time", (1, 1), maxshape=(None, 1),
+            dtype="float64")
+
     def _to_hdf5(self, handle, index):
         """Converts results object into an hdf5 object.
 
@@ -268,6 +273,7 @@ class Results(object):
         eigenvalues_dset = handle["/eigenvalues"]
         time_dset = handle["/time"]
         power_dset = handle["/power"]
+        proc_time_dset = handle["/depletion time"]
 
         # Get number of results stored
         number_shape = list(number_dset.shape)
@@ -296,6 +302,13 @@ class Results(object):
             power_shape[0] = new_shape
             power_dset.resize(power_shape)
 
+            # keep depletion time matrix with one less index
+            # corresponding to one fewer calculation
+            if self.proc_time:
+                proc_shape = list(proc_time_dset.shape)
+                proc_shape[0] = new_shape - 1
+                proc_time_dset.resize(proc_shape)
+
         # If nothing to write, just return
         if len(self.mat_to_ind) == 0:
             return
@@ -314,6 +327,10 @@ class Results(object):
         if comm.rank == 0:
             time_dset[index, :] = self.time
             power_dset[index, :] = self.power
+            if self.proc_time is not None:
+                proc_time_dset[index - 1] = (
+                    self.proc_time / (comm.size * self.n_hdf5_mats)
+                )
 
     @classmethod
     def from_hdf5(cls, handle, step):
@@ -324,8 +341,7 @@ class Results(object):
         handle : h5py.File or h5py.Group
             An HDF5 file or group type to load from.
         step : int
-            What step is this?
-
+            Index for depletion step
         """
         results = cls()
 
@@ -334,12 +350,20 @@ class Results(object):
         eigenvalues_dset = handle["/eigenvalues"]
         time_dset = handle["/time"]
         power_dset = handle["/power"]
+        proc_time_dset = handle["/depletion time"]
 
         results.data = number_dset[step, :, :, :]
         results.k = eigenvalues_dset[step, :]
         results.time = time_dset[step, :]
         results.power = power_dset[step, :]
 
+        # depletion time data has one fewer index, corresponding
+        # to one fewer calculation stage
+        if step < proc_time_dset.shape[0]:
+            results.proc_time = proc_time_dset[step, :]
+        else:
+            results.proc_time = np.array([np.nan])
+
         # Reconstruct dictionaries
         results.volume = OrderedDict()
         results.mat_to_ind = OrderedDict()
@@ -375,7 +399,7 @@ class Results(object):
         return results
 
     @staticmethod
-    def save(op, x, op_results, t, power, step_ind):
+    def save(op, x, op_results, t, power, step_ind, proc_time):
         """Creates and writes depletion results to disk
 
         Parameters
@@ -392,6 +416,8 @@ class Results(object):
             Power during time step
         step_ind : int
             Step index.
+        proc_time: float or None
+            time spent depleting materials.
 
         """
         # Get indexing terms
@@ -422,6 +448,9 @@ class Results(object):
         results.rates = [r.rates for r in op_results]
         results.time = t
         results.power = power
+        results.proc_time = proc_time
+        if step_ind and results.proc_time is not None:
+            results.proc_time = comm.reduce(proc_time, op=mpi_sum)
 
         results.export_to_hdf5("depletion_results.h5", step_ind)
 

From a31e6ced65d8a04c3ffe271ad5a9cf1c0572d82d Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 09:41:29 -0500
Subject: [PATCH 02/19] Add ResultsList.get_depletion_time method

Returns a numpy array with the time spent depleting materials
between transport steps. This array has shape
(N-1, 1), where N is the number of transport steps
---
 openmc/deplete/results_list.py | 15 +++++++++++++++
 1 file changed, 15 insertions(+)

diff --git a/openmc/deplete/results_list.py b/openmc/deplete/results_list.py
index 9fb6eec86c..3442a8ce1e 100644
--- a/openmc/deplete/results_list.py
+++ b/openmc/deplete/results_list.py
@@ -103,3 +103,18 @@ class ResultsList(list):
             eigenvalue[i] = result.k[0]
 
         return time, eigenvalue
+
+    def get_depletion_time(self):
+        """Return an array of time spent depleting materials
+
+        ..note::
+
+            Will have one fewer row than number of other methods,
+            like :meth:`get_eigenvalues`, because no depletion
+            is performed at the final transport stage
+
+        """
+        times = np.empty((len(self) - 1, 1))
+        for ix, res in enumerate(self[:-1]):
+            times[ix] = res.proc_time[0]
+        return times

From 2ecca7bff696f31aeb18f57c657f20751e120376 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 09:43:17 -0500
Subject: [PATCH 03/19] Add wrapped timed_deplete method for providing process
 time

Simple wrapper function that passes all arguments to
openmc.deplete.integrator.deplete and returns a two-item
tuple. The first item is the process time spent depleting
all materials (can be on one of many MPI processes).
The second item is what is returned directly from the
deplete function.
---
 openmc/deplete/integrator/cram.py | 19 +++++++++++++++++++
 1 file changed, 19 insertions(+)

diff --git a/openmc/deplete/integrator/cram.py b/openmc/deplete/integrator/cram.py
index 8a8940df8b..f5da6d696b 100644
--- a/openmc/deplete/integrator/cram.py
+++ b/openmc/deplete/integrator/cram.py
@@ -53,6 +53,25 @@ def deplete(chain, x, rates, dt, print_out=True, matrix_func=None):
     return x_result
 
 
+def timed_deplete(*args, **kwargs):
+    """Wrapper over :func:`deplete` that also returns process time
+
+    All arguments and keyword arguments are passed onto
+    :func:`deplete` directly.
+
+    Returns
+    -------
+    proc_time: float
+        Process time required to return from deplete
+    results: list of numpy arrays
+        Output from :func:`deplete` call
+    """
+
+    start = time.time()
+    results = deplete(*args, **kwargs)
+    return time.time() - start, results
+
+
 def _cram_wrapper(chain, n0, rates, dt, matrix_func=None):
     """Wraps depletion matrix creation / CRAM solve for multiprocess execution
 

From e291725b6b54bb670f9859f619b5e7782f61e17e Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 09:47:36 -0500
Subject: [PATCH 04/19] Modify deplete.integrator.predictor to store, write
 process time

Use new timed_deplete wrapper func to obtain process time.
Pass process time to modifed openmc.deplete.Results.save
to be reduced and written to depletion file.

test_deplete_predictor is modified to ensure that the
shape of the resulting "depletion time" vector is correct
---
 openmc/deplete/integrator/predictor.py     | 11 +++++++----
 tests/unit_tests/test_deplete_predictor.py |  4 ++++
 2 files changed, 11 insertions(+), 4 deletions(-)

diff --git a/openmc/deplete/integrator/predictor.py b/openmc/deplete/integrator/predictor.py
index 6872ae6ef6..b30484c695 100644
--- a/openmc/deplete/integrator/predictor.py
+++ b/openmc/deplete/integrator/predictor.py
@@ -3,7 +3,7 @@
 import copy
 from collections.abc import Iterable
 
-from .cram import deplete
+from .cram import timed_deplete
 from ..results import Results
 
 
@@ -53,6 +53,8 @@ def predictor(operator, timesteps, power=None, power_density=None,
     if not isinstance(power, Iterable):
         power = [power]*len(timesteps)
 
+    proc_time = None
+
     # Generate initial conditions
     with operator as vec:
         # Initialize time and starting index
@@ -74,7 +76,7 @@ def predictor(operator, timesteps, power=None, power_density=None,
                 op_results = [operator(x[0], p)]
 
                 # Create results, write to disk
-                Results.save(operator, x, op_results, [t, t + dt], p, i_res + i)
+                Results.save(operator, x, op_results, [t, t + dt], p, i_res + i, proc_time)
             else:
                 # Get initial concentration
                 x = [operator.prev_res[-1].data[0]]
@@ -89,7 +91,8 @@ def predictor(operator, timesteps, power=None, power_density=None,
                 op_results[0].rates *= ratio_power[0]
 
             # Deplete for full timestep
-            x_end = deplete(chain, x[0], op_results[0].rates, dt, print_out)
+            proc_time, x_end = timed_deplete(
+                chain, x[0], op_results[0].rates, dt, print_out)
 
             # Advance time, update vector
             t += dt
@@ -100,4 +103,4 @@ def predictor(operator, timesteps, power=None, power_density=None,
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), proc_time)
diff --git a/tests/unit_tests/test_deplete_predictor.py b/tests/unit_tests/test_deplete_predictor.py
index b39cc7dcac..0af320eb96 100644
--- a/tests/unit_tests/test_deplete_predictor.py
+++ b/tests/unit_tests/test_deplete_predictor.py
@@ -35,3 +35,7 @@ def test_predictor(run_in_tmpdir):
 
     assert y1[2] == approx(s2[0])
     assert y2[2] == approx(s2[1])
+
+    # Test structure of depletion time dataset
+
+    assert res.get_depletion_time().shape == (len(dt), 1)

From cb3e89bf95fc54ec61f8eaeb914b31e2c8dcf915 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 10:20:50 -0500
Subject: [PATCH 05/19] Always reshape depletion process time dataset

The conditional to reshape/extend the "depletion time"
dataset in the exporting of a deplete.Results object
has been removed. This caused the last exporting to
hang when running with MPI
---
 openmc/deplete/results.py | 7 +++----
 1 file changed, 3 insertions(+), 4 deletions(-)

diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index e37ccb5407..4b4b91baea 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -304,10 +304,9 @@ class Results(object):
 
             # keep depletion time matrix with one less index
             # corresponding to one fewer calculation
-            if self.proc_time:
-                proc_shape = list(proc_time_dset.shape)
-                proc_shape[0] = new_shape - 1
-                proc_time_dset.resize(proc_shape)
+            proc_shape = list(proc_time_dset.shape)
+            proc_shape[0] = new_shape - 1
+            proc_time_dset.resize(proc_shape)
 
         # If nothing to write, just return
         if len(self.mat_to_ind) == 0:

From da23e4646d73fa4a5ab0a71bea40fe4f398ad5da Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 10:25:10 -0500
Subject: [PATCH 06/19] Reduce in deplete.Results.save only if proc_time is not
 None

---
 openmc/deplete/results.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index 4b4b91baea..0d6e839154 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -448,7 +448,7 @@ class Results(object):
         results.time = t
         results.power = power
         results.proc_time = proc_time
-        if step_ind and results.proc_time is not None:
+        if results.proc_time is not None:
             results.proc_time = comm.reduce(proc_time, op=mpi_sum)
 
         results.export_to_hdf5("depletion_results.h5", step_ind)

From 3a963b4183efbf52193c02393ea9ddd6585cafc2 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 11:49:17 -0500
Subject: [PATCH 07/19] Modify depletion time matrix size; Remove index - 1 in
 Results.save

The index - 1 was causing values to be written in a correctly sized
matrix, but in an offset pattern, e.g. first value written
to last position, second value to first position, etc.
The ResultsList.get_depletion_time method will already
trim this array in the retrieval process.
---
 openmc/deplete/results.py | 6 ++----
 1 file changed, 2 insertions(+), 4 deletions(-)

diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index 0d6e839154..f7d7550c35 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -302,10 +302,8 @@ class Results(object):
             power_shape[0] = new_shape
             power_dset.resize(power_shape)
 
-            # keep depletion time matrix with one less index
-            # corresponding to one fewer calculation
             proc_shape = list(proc_time_dset.shape)
-            proc_shape[0] = new_shape - 1
+            proc_shape[0] = new_shape
             proc_time_dset.resize(proc_shape)
 
         # If nothing to write, just return
@@ -327,7 +325,7 @@ class Results(object):
             time_dset[index, :] = self.time
             power_dset[index, :] = self.power
             if self.proc_time is not None:
-                proc_time_dset[index - 1] = (
+                proc_time_dset[index] = (
                     self.proc_time / (comm.size * self.n_hdf5_mats)
                 )
 

From f1157e90d43d9016058565627eeb8cb2756fe010 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 12:08:33 -0500
Subject: [PATCH 08/19] Store and test storage of process time for cecm, celi,
 epc_rk4

---
 openmc/deplete/integrator/cecm.py        | 16 +++++++++----
 openmc/deplete/integrator/celi.py        | 12 +++++-----
 openmc/deplete/integrator/cf4.py         | 30 ++++++++++++++----------
 openmc/deplete/integrator/epc_rk4.py     | 26 ++++++++++++--------
 tests/unit_tests/test_deplete_cecm.py    |  6 +++++
 tests/unit_tests/test_deplete_celi.py    |  6 +++++
 tests/unit_tests/test_deplete_cf4.py     |  6 +++++
 tests/unit_tests/test_deplete_epc_rk4.py |  6 +++++
 8 files changed, 74 insertions(+), 34 deletions(-)

diff --git a/openmc/deplete/integrator/cecm.py b/openmc/deplete/integrator/cecm.py
index b1bcec448d..e463bca575 100644
--- a/openmc/deplete/integrator/cecm.py
+++ b/openmc/deplete/integrator/cecm.py
@@ -3,7 +3,7 @@
 import copy
 from collections.abc import Iterable
 
-from .cram import deplete
+from .cram import timed_deplete
 from ..results import Results
 
 
@@ -69,6 +69,8 @@ def cecm(operator, timesteps, power=None, power_density=None, print_out=True):
 
         chain = operator.chain
 
+        proc_time = None
+
         for i, (dt, p) in enumerate(zip(timesteps, power)):
             # Get beginning-of-timestep concentrations and reaction rates
             # Avoid doing first transport run if already done in previous
@@ -94,7 +96,8 @@ def cecm(operator, timesteps, power=None, power_density=None, print_out=True):
                 op_results[0].rates *= ratio_power[0]
 
             # Deplete for first half of timestep
-            x_middle = deplete(chain, x[0], op_results[0].rates, dt/2, print_out)
+            proc_time, x_middle = timed_deplete(
+                chain, x[0], op_results[0].rates, dt/2, print_out)
 
             # Get middle-of-timestep reaction rates
             x.append(x_middle)
@@ -102,10 +105,13 @@ def cecm(operator, timesteps, power=None, power_density=None, print_out=True):
 
             # Deplete for full timestep using beginning-of-step materials
             # and middle-of-timestep reaction rates
-            x_end = deplete(chain, x[0], op_results[1].rates, dt, print_out)
+            pt_end, x_end = timed_deplete(
+                chain, x[0], op_results[1].rates, dt, print_out)
 
             # Create results, write to disk
-            Results.save(operator, x, op_results, [t, t + dt], p, i_res + i)
+            Results.save(
+                operator, x, op_results, [t, t + dt], p, i_res + i,
+                proc_time + pt_end)
 
             # Advance time, update vector
             t += dt
@@ -116,4 +122,4 @@ def cecm(operator, timesteps, power=None, power_density=None, print_out=True):
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
diff --git a/openmc/deplete/integrator/celi.py b/openmc/deplete/integrator/celi.py
index e7b8014990..b1fe4d5e89 100644
--- a/openmc/deplete/integrator/celi.py
+++ b/openmc/deplete/integrator/celi.py
@@ -3,7 +3,7 @@
 import copy
 from collections.abc import Iterable
 
-from .cram import deplete
+from .cram import timed_deplete
 from ..results import Results
 
 
@@ -88,7 +88,7 @@ def celi(operator, timesteps, power=None, power_density=None,
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
 
 
 def celi_inner(operator, vec, p, i, i_res, t, dt, print_out):
@@ -149,19 +149,19 @@ def celi_inner(operator, vec, p, i, i_res, t, dt, print_out):
         op_results[0].rates *= ratio_power[0]
 
     # Deplete to end
-    x_new = deplete(chain, x[0], op_results[0].rates, dt, print_out)
+    proc_time, x_new = timed_deplete(chain, x[0], op_results[0].rates, dt, print_out)
     x.append(x_new)
     op_results.append(operator(x[1], p))
 
     # Deplete with two matrix exponentials
     rates = list(zip(op_results[0].rates, op_results[1].rates))
-    x_end = deplete(chain, x[0], rates, dt, print_out,
+    time_1, x_end = timed_deplete(chain, x[0], rates, dt, print_out,
                     matrix_func=_celi_f1)
-    x_end = deplete(chain, x_end, rates, dt, print_out,
+    time_2, x_end = timed_deplete(chain, x_end, rates, dt, print_out,
                     matrix_func=_celi_f2)
 
     # Create results, write to disk
-    Results.save(operator, x, op_results, [t, t + dt], p, i_res + i)
+    Results.save(operator, x, op_results, [t, t + dt], p, i_res + i, proc_time + time_1 + time_2)
 
     # return updated time and vectors
     return x_end, t + dt, op_results[0]
diff --git a/openmc/deplete/integrator/cf4.py b/openmc/deplete/integrator/cf4.py
index fecb917517..37d7023a84 100644
--- a/openmc/deplete/integrator/cf4.py
+++ b/openmc/deplete/integrator/cf4.py
@@ -3,7 +3,7 @@
 import copy
 from collections.abc import Iterable
 
-from .cram import deplete
+from .cram import timed_deplete
 from ..results import Results
 
 
@@ -119,34 +119,38 @@ def cf4(operator, timesteps, power=None, power_density=None, print_out=True):
                 op_results[0].rates *= ratio_power[0]
 
             # Step 1: deplete with matrix 1/2*A(y0)
-            x_new = deplete(chain, x[0], op_results[0].rates, dt, print_out,
-                            matrix_func=_cf4_f1)
+            time_1, x_new = timed_deplete(
+                chain, x[0], op_results[0].rates, dt, print_out,
+                matrix_func=_cf4_f1)
             x.append(x_new)
             op_results.append(operator(x_new, p))
 
             # Step 2: deplete with matrix 1/2*A(y1)
-            x_new = deplete(chain, x[0], op_results[1].rates, dt, print_out,
-                            matrix_func=_cf4_f1)
+            time_2, x_new = timed_deplete(
+                chain, x[0], op_results[1].rates, dt, print_out,
+                matrix_func=_cf4_f1)
             x.append(x_new)
             op_results.append(operator(x_new, p))
 
             # Step 3: deplete with matrix -1/2*A(y0)+A(y2)
             rates = list(zip(op_results[0].rates, op_results[2].rates))
-            x_new = deplete(chain, x[1], rates, dt, print_out,
-                            matrix_func=_cf4_f2)
+            time_3, x_new = timed_deplete(
+                chain, x[1], rates, dt, print_out, matrix_func=_cf4_f2)
             x.append(x_new)
             op_results.append(operator(x_new, p))
 
             # Step 4: deplete with two matrix exponentials
             rates = list(zip(op_results[0].rates, op_results[1].rates,
                              op_results[2].rates, op_results[3].rates))
-            x_end = deplete(chain, x[0], rates, dt, print_out,
-                            matrix_func=_cf4_f3)
-            x_end = deplete(chain, x_end, rates, dt, print_out,
-                            matrix_func=_cf4_f4)
+            time_4, x_end = timed_deplete(
+                chain, x[0], rates, dt, print_out, matrix_func=_cf4_f3)
+            time_5, x_end = timed_deplete(
+                chain, x_end, rates, dt, print_out, matrix_func=_cf4_f4)
 
             # Create results, write to disk
-            Results.save(operator, x, op_results, [t, t + dt], p, i_res + i)
+            Results.save(
+                operator, x, op_results, [t, t + dt], p, i_res + i,
+                sum((time_1, time_2, time_3, time_4, time_5)))
 
             # Advance time, update vector
             t += dt
@@ -157,4 +161,4 @@ def cf4(operator, timesteps, power=None, power_density=None, print_out=True):
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
diff --git a/openmc/deplete/integrator/epc_rk4.py b/openmc/deplete/integrator/epc_rk4.py
index 58c4d08f7e..c57c71204a 100644
--- a/openmc/deplete/integrator/epc_rk4.py
+++ b/openmc/deplete/integrator/epc_rk4.py
@@ -3,7 +3,7 @@
 import copy
 from collections.abc import Iterable
 
-from .cram import deplete
+from .cram import timed_deplete
 from ..results import Results
 
 
@@ -105,30 +105,35 @@ def epc_rk4(operator, timesteps, power=None, power_density=None, print_out=True)
                 op_results[0].rates *= ratio_power[0]
 
             # Step 1: deplete with matrix 1/2*A(y0)
-            x_new = deplete(chain, x[0], op_results[0].rates, dt, print_out,
-                            matrix_func=_rk4_f1)
+            time_1, x_new = timed_deplete(
+                chain, x[0], op_results[0].rates, dt, print_out,
+                matrix_func=_rk4_f1)
             x.append(x_new)
             op_results.append(operator(x[1], p))
 
             # Step 2: deplete with matrix 1/2*A(y1)
-            x_new = deplete(chain, x[0], op_results[1].rates, dt, print_out,
-                            matrix_func=_rk4_f1)
+            time_2, x_new = timed_deplete(
+                chain, x[0], op_results[1].rates, dt, print_out,
+                matrix_func=_rk4_f1)
             x.append(x_new)
             op_results.append(operator(x[2], p))
 
             # Step 3: deplete with matrix A(y2)
-            x_new = deplete(chain, x[0], op_results[2].rates, dt, print_out)
+            time_3, x_new = timed_deplete(
+                chain, x[0], op_results[2].rates, dt, print_out)
             x.append(x_new)
             op_results.append(operator(x[3], p))
 
             # Step 4: deplete with matrix 1/6*A(y0)+1/3*A(y1)+1/3*A(y2)+1/6*A(y3)
             rates = list(zip(op_results[0].rates, op_results[1].rates,
                              op_results[2].rates, op_results[3].rates))
-            x_end = deplete(chain, x[0], rates, dt, print_out,
-                            matrix_func=_rk4_f4)
+            time_4, x_end = timed_deplete(
+                chain, x[0], rates, dt, print_out, matrix_func=_rk4_f4)
 
             # Create results, write to disk
-            Results.save(operator, x, op_results, [t, t + dt], p, i_res + i)
+            Results.save(
+                operator, x, op_results, [t, t + dt], p, i_res + i,
+                sum((time_1, time_2, time_3, time_4)))
 
             # Advance time, update vector
             t += dt
@@ -139,4 +144,5 @@ def epc_rk4(operator, timesteps, power=None, power_density=None, print_out=True)
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
+        Results.save(
+            operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
diff --git a/tests/unit_tests/test_deplete_cecm.py b/tests/unit_tests/test_deplete_cecm.py
index 466a2eec74..d143658af3 100644
--- a/tests/unit_tests/test_deplete_cecm.py
+++ b/tests/unit_tests/test_deplete_cecm.py
@@ -35,3 +35,9 @@ def test_cecm(run_in_tmpdir):
 
     assert y1[2] == approx(s2[0])
     assert y2[2] == approx(s2[1])
+
+    # Test structure of depletion time dataset
+
+    dep_time = res.get_depletion_time()
+    assert dep_time.shape == (len(dt), 1)
+    assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_celi.py b/tests/unit_tests/test_deplete_celi.py
index 0b23b1a438..bb479f6cc2 100644
--- a/tests/unit_tests/test_deplete_celi.py
+++ b/tests/unit_tests/test_deplete_celi.py
@@ -35,3 +35,9 @@ def test_celi(run_in_tmpdir):
 
     assert y1[2] == approx(s2[0])
     assert y2[2] == approx(s2[1])
+
+    # Test structure of depletion time dataset
+
+    dep_time = res.get_depletion_time()
+    assert dep_time.shape == (len(dt), 1)
+    assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_cf4.py b/tests/unit_tests/test_deplete_cf4.py
index 784d228beb..5b7d08c229 100644
--- a/tests/unit_tests/test_deplete_cf4.py
+++ b/tests/unit_tests/test_deplete_cf4.py
@@ -35,3 +35,9 @@ def test_cf4(run_in_tmpdir):
 
     assert y1[2] == approx(s2[0])
     assert y2[2] == approx(s2[1])
+
+    # Test structure of depletion time dataset
+
+    dep_time = res.get_depletion_time()
+    assert dep_time.shape == (len(dt), 1)
+    assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_epc_rk4.py b/tests/unit_tests/test_deplete_epc_rk4.py
index dd15bf8e1d..0623f9e6e8 100644
--- a/tests/unit_tests/test_deplete_epc_rk4.py
+++ b/tests/unit_tests/test_deplete_epc_rk4.py
@@ -35,3 +35,9 @@ def test_epc_rk4(run_in_tmpdir):
 
     assert y1[2] == approx(s2[0])
     assert y2[2] == approx(s2[1])
+
+    # Test structure of depletion time dataset
+
+    dep_time = res.get_depletion_time()
+    assert dep_time.shape == (len(dt), 1)
+    assert all(dep_time > 0)

From 6583834bd4a8460d708be6492025ca65500e9323 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 13:06:54 -0500
Subject: [PATCH 09/19] Store depletion time for leqi, si_celi, si_leqi
 routines

---
 openmc/deplete/integrator/leqi.py        | 25 +++++++++++----------
 openmc/deplete/integrator/si_celi.py     | 19 +++++++++-------
 openmc/deplete/integrator/si_leqi.py     | 28 ++++++++++++++----------
 tests/unit_tests/test_deplete_leqi.py    |  6 +++++
 tests/unit_tests/test_deplete_si_celi.py | 12 ++++++++++
 tests/unit_tests/test_deplete_si_leqi.py |  6 +++++
 6 files changed, 66 insertions(+), 30 deletions(-)

diff --git a/openmc/deplete/integrator/leqi.py b/openmc/deplete/integrator/leqi.py
index 66b0c253ad..9e1ed56a4e 100644
--- a/openmc/deplete/integrator/leqi.py
+++ b/openmc/deplete/integrator/leqi.py
@@ -5,7 +5,7 @@ from collections.abc import Iterable
 from itertools import repeat
 
 from .celi import celi_inner
-from .cram import deplete
+from .cram import timed_deplete
 from ..results import Results
 
 
@@ -137,22 +137,24 @@ def leqi(operator, timesteps, power=None, power_density=None, print_out=True):
 
             inputs = list(zip(op_res_last.rates, op_results[0].rates,
                               repeat(dt_l), repeat(dt)))
-            x_new = deplete(chain, x[0], inputs, dt, print_out,
-                            matrix_func=_leqi_f1)
-            x_new = deplete(chain, x_new, inputs, dt, print_out,
-                            matrix_func=_leqi_f2)
+            time_1, x_new = timed_deplete(
+                chain, x[0], inputs, dt, print_out, matrix_func=_leqi_f1)
+            time_2, x_new = timed_deplete(
+                chain, x_new, inputs, dt, print_out, matrix_func=_leqi_f2)
             x.append(x_new)
             op_results.append(operator(x[1], p))
 
             inputs = list(zip(op_res_last.rates, op_results[0].rates,
                               op_results[1].rates, repeat(dt_l), repeat(dt)))
-            x_new = deplete(chain, x[0], inputs, dt, print_out,
-                            matrix_func=_leqi_f3)
-            x_new = deplete(chain, x_new, inputs, dt, print_out,
-                            matrix_func=_leqi_f4)
+            time_3, x_new = timed_deplete(
+                chain, x[0], inputs, dt, print_out, matrix_func=_leqi_f3)
+            time_4, x_new = timed_deplete(
+                chain, x_new, inputs, dt, print_out, matrix_func=_leqi_f4)
 
             # Create results, write to disk
-            Results.save(operator, x, op_results, [t, t+dt], p, i_res+i)
+            Results.save(
+                operator, x, op_results, [t, t+dt], p, i_res+i,
+                sum((time_1, time_2, time_3, time_4)))
 
             # update results
             op_res_last = copy.deepcopy(op_results[0])
@@ -164,4 +166,5 @@ def leqi(operator, timesteps, power=None, power_density=None, print_out=True):
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
+        Results.save(
+            operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
diff --git a/openmc/deplete/integrator/si_celi.py b/openmc/deplete/integrator/si_celi.py
index 4f34bace6d..dcddf5ea39 100644
--- a/openmc/deplete/integrator/si_celi.py
+++ b/openmc/deplete/integrator/si_celi.py
@@ -3,7 +3,7 @@
 import copy
 from collections.abc import Iterable
 
-from .cram import deplete
+from .cram import timed_deplete
 from ..results import Results
 from ..abc import OperatorResult
 from .celi import _celi_f1, _celi_f2
@@ -94,7 +94,8 @@ def si_celi(operator, timesteps, power=None, power_density=None,
                                              i, i_res, t, dt, print_out, m)
 
         # Create results for last point, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
+        Results.save(
+            operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
 
 
 def si_celi_inner(operator, x, op_results, p, i, i_res, t, dt, print_out, m=10):
@@ -136,7 +137,8 @@ def si_celi_inner(operator, x, op_results, p, i, i_res, t, dt, print_out, m=10):
     chain = operator.chain
 
     # Deplete to end
-    x_new = deplete(chain, x[0], op_results[0].rates, dt, print_out)
+    proc_time, x_new = timed_deplete(
+        chain, x[0], op_results[0].rates, dt, print_out)
     x.append(x_new)
 
     for j in range(m + 1):
@@ -150,14 +152,15 @@ def si_celi_inner(operator, x, op_results, p, i, i_res, t, dt, print_out, m=10):
             op_res_bar = OperatorResult(k, rates)
 
         rates = list(zip(op_results[0].rates, op_res_bar.rates))
-        x_new = deplete(chain, x[0], rates, dt, print_out,
-                        matrix_func=_celi_f1)
-        x_new = deplete(chain, x_new, rates, dt, print_out,
-                        matrix_func=_celi_f2)
+        time_1, x_new = timed_deplete(
+            chain, x[0], rates, dt, print_out, matrix_func=_celi_f1)
+        time_2, x_new = timed_deplete(
+            chain, x_new, rates, dt, print_out, matrix_func=_celi_f2)
+        proc_time += time_1 + time_2
 
     # Create results, write to disk
     op_results.append(op_res_bar)
-    Results.save(operator, x, op_results, [t, t+dt], p, i_res+i)
+    Results.save(operator, x, op_results, [t, t+dt], p, i_res+i, proc_time)
 
     # return updated time and vectors
     return [x_new], t + dt, [op_res_bar]
diff --git a/openmc/deplete/integrator/si_leqi.py b/openmc/deplete/integrator/si_leqi.py
index 05002880a1..b4b586feed 100644
--- a/openmc/deplete/integrator/si_leqi.py
+++ b/openmc/deplete/integrator/si_leqi.py
@@ -6,7 +6,7 @@ from itertools import repeat
 
 from .si_celi import si_celi_inner
 from .leqi import _leqi_f1, _leqi_f2, _leqi_f3, _leqi_f4
-from .cram import deplete
+from .cram import timed_deplete
 from ..results import Results
 from ..abc import OperatorResult
 
@@ -112,12 +112,14 @@ def si_leqi(operator, timesteps, power=None, power_density=None,
             # Perform remaining LE/QI
             inputs = list(zip(op_res_last.rates, op_results[0].rates,
                               repeat(dt_l), repeat(dt)))
-            x_new = deplete(chain, x[0], inputs, dt, print_out,
-                            matrix_func=_leqi_f1)
-            x_new = deplete(chain, x_new, inputs, dt, print_out,
-                            matrix_func=_leqi_f2)
+            proc_time, x_new = timed_deplete(
+                chain, x[0], inputs, dt, print_out, matrix_func=_leqi_f1)
+            time_1, x_new = timed_deplete(
+                chain, x_new, inputs, dt, print_out, matrix_func=_leqi_f2)
             x.append(x_new)
 
+            proc_time += time_1
+
             # Loop on inner
             for j in range(m + 1):
                 op_res = operator(x_new, p)
@@ -131,14 +133,17 @@ def si_leqi(operator, timesteps, power=None, power_density=None,
 
                 inputs = list(zip(op_res_last.rates, op_results[0].rates,
                                   op_res_bar.rates, repeat(dt_l), repeat(dt)))
-                x_new = deplete(chain, x[0], inputs, dt, print_out,
-                                matrix_func=_leqi_f3)
-                x_new = deplete(chain, x_new, inputs, dt, print_out,
-                                matrix_func=_leqi_f4)
+                time_1, x_new = timed_deplete(
+                    chain, x[0], inputs, dt, print_out, matrix_func=_leqi_f3)
+                time_2, x_new = timed_deplete(
+                    chain, x_new, inputs, dt, print_out, matrix_func=_leqi_f4)
+
+                proc_time += time_1 + time_2
 
             # Create results, write to disk
             op_results.append(op_res_bar)
-            Results.save(operator, x, op_results, [t, t+dt], p, i_res+i)
+            Results.save(
+                operator, x, op_results, [t, t+dt], p, i_res+i, proc_time)
 
             # update results
             x = [x_new]
@@ -148,4 +153,5 @@ def si_leqi(operator, timesteps, power=None, power_density=None,
             dt_l = dt
 
         # Create results for last point, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res+len(timesteps))
+        Results.save(
+                operator, x, op_results, [t, t], p, i_res+len(timesteps), None)
diff --git a/tests/unit_tests/test_deplete_leqi.py b/tests/unit_tests/test_deplete_leqi.py
index a1a32e092c..97a841f6ce 100644
--- a/tests/unit_tests/test_deplete_leqi.py
+++ b/tests/unit_tests/test_deplete_leqi.py
@@ -35,3 +35,9 @@ def test_leqi(run_in_tmpdir):
 
     assert y1[2] == approx(s2[0])
     assert y2[2] == approx(s2[1])
+
+    # Test structure of depletion time dataset
+
+    dep_time = res.get_depletion_time()
+    assert dep_time.shape == (len(dt), 1)
+    assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_si_celi.py b/tests/unit_tests/test_deplete_si_celi.py
index 8e84bcba8a..c064172976 100644
--- a/tests/unit_tests/test_deplete_si_celi.py
+++ b/tests/unit_tests/test_deplete_si_celi.py
@@ -35,3 +35,15 @@ def test_si_celi(run_in_tmpdir):
 
     assert y1[2] == approx(s2[0])
     assert y2[2] == approx(s2[1])
+
+    # Test structure of depletion time dataset
+
+    dep_time = res.get_depletion_time()
+    assert dep_time.shape == (len(dt), 1)
+    assert all(dep_time > 0)
+
+    # Test structure of depletion time dataset
+
+    dep_time = res.get_depletion_time()
+    assert dep_time.shape == (len(dt), 1)
+    assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_si_leqi.py b/tests/unit_tests/test_deplete_si_leqi.py
index 4814d8ab86..1003a73161 100644
--- a/tests/unit_tests/test_deplete_si_leqi.py
+++ b/tests/unit_tests/test_deplete_si_leqi.py
@@ -35,3 +35,9 @@ def test_si_leqi(run_in_tmpdir):
 
     assert y1[2] == approx(s2[0])
     assert y2[2] == approx(s2[1])
+
+    # Test structure of depletion time dataset
+
+    dep_time = res.get_depletion_time()
+    assert dep_time.shape == (len(dt), 1)
+    assert all(dep_time > 0)

From df07ecbcfafc6daf2f6712c3b055319954bf0e2b Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 13:52:35 -0500
Subject: [PATCH 10/19] Pass proc_time of None for test_deplete_integrator.py

---
 tests/unit_tests/test_deplete_integrator.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/tests/unit_tests/test_deplete_integrator.py b/tests/unit_tests/test_deplete_integrator.py
index 8cade45b14..c90d7122d1 100644
--- a/tests/unit_tests/test_deplete_integrator.py
+++ b/tests/unit_tests/test_deplete_integrator.py
@@ -75,8 +75,8 @@ def test_results_save(run_in_tmpdir):
 
     op_result1 = [OperatorResult(k, rates) for k, rates in zip(eigvl1, rate1)]
     op_result2 = [OperatorResult(k, rates) for k, rates in zip(eigvl2, rate2)]
-    Results.save(op, x1, op_result1, t1, 0, 0)
-    Results.save(op, x2, op_result2, t2, 0, 1)
+    Results.save(op, x1, op_result1, t1, 0, 0, None)
+    Results.save(op, x2, op_result2, t2, 0, 1, None)
 
     # Load the files
     res = ResultsList("depletion_results.h5")

From bd0a9a18e4a4e92f04a6bfb9c6a845f470d35ec5 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 14:09:34 -0500
Subject: [PATCH 11/19] Modify Results.from_hdf5 to support for files w/o
 depletion time

Needed because ./tests/regression_tests/deplete/test_reference.h5
does not contain this dataset. Now, a check is performed
to see if "depletion file" is present in the root of the hdf5 file.
Otherwise, the proc_time is set to the array of nans
---
 openmc/deplete/results.py | 12 ++++++------
 1 file changed, 6 insertions(+), 6 deletions(-)

diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index f7d7550c35..8211a59b14 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -347,18 +347,18 @@ class Results(object):
         eigenvalues_dset = handle["/eigenvalues"]
         time_dset = handle["/time"]
         power_dset = handle["/power"]
-        proc_time_dset = handle["/depletion time"]
 
         results.data = number_dset[step, :, :, :]
         results.k = eigenvalues_dset[step, :]
         results.time = time_dset[step, :]
         results.power = power_dset[step, :]
 
-        # depletion time data has one fewer index, corresponding
-        # to one fewer calculation stage
-        if step < proc_time_dset.shape[0]:
-            results.proc_time = proc_time_dset[step, :]
-        else:
+        if "depletion time" in handle:
+            proc_time_dset = handle["/depletion time"]
+            if step < proc_time_dset.shape[0]:
+                results.proc_time = proc_time_dset[step, :]
+
+        if results.proc_time is None:
             results.proc_time = np.array([np.nan])
 
         # Reconstruct dictionaries

From 71d31d6ea67dff85911cb7a8f82b31befbd9f33e Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 21 Jun 2019 14:14:35 -0500
Subject: [PATCH 12/19] Update test_deplete_predictor to check depletion time

Due to where the results are saved for the predictor routine,
at the top of the for loop vs. the end, the positions
are incremented by one in the depletion time. This causes
a value of zero to be returned as the first item in the array
from get_depletion_time.

This save was moved in commit 87aed7a0b095b3902156e2a385bc46db2ba7f5ac
of PR #1014 to make the restart feature easier.
---
 openmc/deplete/integrator/predictor.py     | 2 +-
 tests/unit_tests/test_deplete_predictor.py | 4 +++-
 2 files changed, 4 insertions(+), 2 deletions(-)

diff --git a/openmc/deplete/integrator/predictor.py b/openmc/deplete/integrator/predictor.py
index b30484c695..54a6258dc4 100644
--- a/openmc/deplete/integrator/predictor.py
+++ b/openmc/deplete/integrator/predictor.py
@@ -103,4 +103,4 @@ def predictor(operator, timesteps, power=None, power_density=None,
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), proc_time)
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
diff --git a/tests/unit_tests/test_deplete_predictor.py b/tests/unit_tests/test_deplete_predictor.py
index 0af320eb96..80b060cc9e 100644
--- a/tests/unit_tests/test_deplete_predictor.py
+++ b/tests/unit_tests/test_deplete_predictor.py
@@ -38,4 +38,6 @@ def test_predictor(run_in_tmpdir):
 
     # Test structure of depletion time dataset
 
-    assert res.get_depletion_time().shape == (len(dt), 1)
+    dep_time = res.get_depletion_time()
+    assert dep_time.shape == (len(dt), 1)
+    assert all(dep_time > 0)

From 875ba94825b766ee56afa7e62bbf90ea024d7aca Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Wed, 26 Jun 2019 17:47:02 -0500
Subject: [PATCH 13/19] Set proc_time argument in Results.save to default to
 None

Adjust final calls to integrator routines accordingly.

Minor documentation changes
---
 openmc/deplete/integrator/cecm.py      | 4 +---
 openmc/deplete/integrator/celi.py      | 2 +-
 openmc/deplete/integrator/cf4.py       | 4 ++--
 openmc/deplete/integrator/epc_rk4.py   | 2 +-
 openmc/deplete/integrator/leqi.py      | 2 +-
 openmc/deplete/integrator/predictor.py | 2 +-
 openmc/deplete/integrator/si_celi.py   | 2 +-
 openmc/deplete/integrator/si_leqi.py   | 2 +-
 openmc/deplete/results.py              | 8 +++++---
 9 files changed, 14 insertions(+), 14 deletions(-)

diff --git a/openmc/deplete/integrator/cecm.py b/openmc/deplete/integrator/cecm.py
index e463bca575..758302b027 100644
--- a/openmc/deplete/integrator/cecm.py
+++ b/openmc/deplete/integrator/cecm.py
@@ -69,8 +69,6 @@ def cecm(operator, timesteps, power=None, power_density=None, print_out=True):
 
         chain = operator.chain
 
-        proc_time = None
-
         for i, (dt, p) in enumerate(zip(timesteps, power)):
             # Get beginning-of-timestep concentrations and reaction rates
             # Avoid doing first transport run if already done in previous
@@ -122,4 +120,4 @@ def cecm(operator, timesteps, power=None, power_density=None, print_out=True):
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
diff --git a/openmc/deplete/integrator/celi.py b/openmc/deplete/integrator/celi.py
index b1fe4d5e89..2e80c8b739 100644
--- a/openmc/deplete/integrator/celi.py
+++ b/openmc/deplete/integrator/celi.py
@@ -88,7 +88,7 @@ def celi(operator, timesteps, power=None, power_density=None,
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
 
 
 def celi_inner(operator, vec, p, i, i_res, t, dt, print_out):
diff --git a/openmc/deplete/integrator/cf4.py b/openmc/deplete/integrator/cf4.py
index 37d7023a84..8b7f3a05a0 100644
--- a/openmc/deplete/integrator/cf4.py
+++ b/openmc/deplete/integrator/cf4.py
@@ -150,7 +150,7 @@ def cf4(operator, timesteps, power=None, power_density=None, print_out=True):
             # Create results, write to disk
             Results.save(
                 operator, x, op_results, [t, t + dt], p, i_res + i,
-                sum((time_1, time_2, time_3, time_4, time_5)))
+                time_1 + time_2 + time_3 + time_4 + time_5)
 
             # Advance time, update vector
             t += dt
@@ -161,4 +161,4 @@ def cf4(operator, timesteps, power=None, power_density=None, print_out=True):
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
diff --git a/openmc/deplete/integrator/epc_rk4.py b/openmc/deplete/integrator/epc_rk4.py
index c57c71204a..1f0fbe2a84 100644
--- a/openmc/deplete/integrator/epc_rk4.py
+++ b/openmc/deplete/integrator/epc_rk4.py
@@ -145,4 +145,4 @@ def epc_rk4(operator, timesteps, power=None, power_density=None, print_out=True)
 
         # Create results, write to disk
         Results.save(
-            operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
+            operator, x, op_results, [t, t], p, i_res + len(timesteps))
diff --git a/openmc/deplete/integrator/leqi.py b/openmc/deplete/integrator/leqi.py
index 9e1ed56a4e..465139abf4 100644
--- a/openmc/deplete/integrator/leqi.py
+++ b/openmc/deplete/integrator/leqi.py
@@ -167,4 +167,4 @@ def leqi(operator, timesteps, power=None, power_density=None, print_out=True):
 
         # Create results, write to disk
         Results.save(
-            operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
+            operator, x, op_results, [t, t], p, i_res + len(timesteps))
diff --git a/openmc/deplete/integrator/predictor.py b/openmc/deplete/integrator/predictor.py
index 54a6258dc4..a32ef12b9d 100644
--- a/openmc/deplete/integrator/predictor.py
+++ b/openmc/deplete/integrator/predictor.py
@@ -103,4 +103,4 @@ def predictor(operator, timesteps, power=None, power_density=None,
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
diff --git a/openmc/deplete/integrator/si_celi.py b/openmc/deplete/integrator/si_celi.py
index dcddf5ea39..26087c1f2c 100644
--- a/openmc/deplete/integrator/si_celi.py
+++ b/openmc/deplete/integrator/si_celi.py
@@ -95,7 +95,7 @@ def si_celi(operator, timesteps, power=None, power_density=None,
 
         # Create results for last point, write to disk
         Results.save(
-            operator, x, op_results, [t, t], p, i_res + len(timesteps), None)
+            operator, x, op_results, [t, t], p, i_res + len(timesteps))
 
 
 def si_celi_inner(operator, x, op_results, p, i, i_res, t, dt, print_out, m=10):
diff --git a/openmc/deplete/integrator/si_leqi.py b/openmc/deplete/integrator/si_leqi.py
index b4b586feed..7d01782bae 100644
--- a/openmc/deplete/integrator/si_leqi.py
+++ b/openmc/deplete/integrator/si_leqi.py
@@ -154,4 +154,4 @@ def si_leqi(operator, timesteps, power=None, power_density=None,
 
         # Create results for last point, write to disk
         Results.save(
-                operator, x, op_results, [t, t], p, i_res+len(timesteps), None)
+                operator, x, op_results, [t, t], p, i_res+len(timesteps))
diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index 8211a59b14..30a2dcf8ac 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -47,6 +47,8 @@ class Results(object):
         Number of stages in simulation.
     data : numpy.ndarray
         Atom quantity, stored by stage, mat, then by nuclide.
+    proc_time: int
+        Amount of process time spent depleting materials
 
     """
     def __init__(self):
@@ -247,7 +249,7 @@ class Results(object):
                               dtype='float64')
 
         handle.create_dataset(
-            "depletion time", (1, 1), maxshape=(None, 1),
+            "depletion time", (1,), maxshape=(None,),
             dtype="float64")
 
     def _to_hdf5(self, handle, index):
@@ -396,7 +398,7 @@ class Results(object):
         return results
 
     @staticmethod
-    def save(op, x, op_results, t, power, step_ind, proc_time):
+    def save(op, x, op_results, t, power, step_ind, proc_time=None):
         """Creates and writes depletion results to disk
 
         Parameters
@@ -413,7 +415,7 @@ class Results(object):
             Power during time step
         step_ind : int
             Step index.
-        proc_time: float or None
+        proc_time : float or None
             time spent depleting materials.
 
         """

From 8deda73647971637051f2d31ea29cef4a1809355 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Wed, 26 Jun 2019 17:50:40 -0500
Subject: [PATCH 14/19] Update Results and ResultsList doc on proc_time

Process time that is stored in the dataset is the
average time spent depleting a material across
all materials and processes.
---
 openmc/deplete/results.py      | 7 +++++--
 openmc/deplete/results_list.py | 2 +-
 2 files changed, 6 insertions(+), 3 deletions(-)

diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index 30a2dcf8ac..8c47706b75 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -48,7 +48,8 @@ class Results(object):
     data : numpy.ndarray
         Atom quantity, stored by stage, mat, then by nuclide.
     proc_time: int
-        Amount of process time spent depleting materials
+        Average time spent depleting a material across all
+        materials and processes
 
     """
     def __init__(self):
@@ -416,7 +417,9 @@ class Results(object):
         step_ind : int
             Step index.
         proc_time : float or None
-            time spent depleting materials.
+            Total process time spent depleting materials. This may
+            be process-dependent and will be reduced across MPI
+            processes.
 
         """
         # Get indexing terms
diff --git a/openmc/deplete/results_list.py b/openmc/deplete/results_list.py
index 3442a8ce1e..58dcf0c59d 100644
--- a/openmc/deplete/results_list.py
+++ b/openmc/deplete/results_list.py
@@ -105,7 +105,7 @@ class ResultsList(list):
         return time, eigenvalue
 
     def get_depletion_time(self):
-        """Return an array of time spent depleting materials
+        """Return an array of the average time to deplete a material
 
         ..note::
 

From 544cbd0c39dc1fb9d4dcaa6edd21a39c4114e225 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Wed, 26 Jun 2019 17:52:27 -0500
Subject: [PATCH 15/19] Use unittest.Mock for MPI when collecting proc_time

---
 openmc/deplete/__init__.py | 5 +++--
 openmc/deplete/results.py  | 4 ++--
 2 files changed, 5 insertions(+), 4 deletions(-)

diff --git a/openmc/deplete/__init__.py b/openmc/deplete/__init__.py
index 32588588ef..08338c3f74 100644
--- a/openmc/deplete/__init__.py
+++ b/openmc/deplete/__init__.py
@@ -5,16 +5,17 @@ openmc.deplete
 A depletion front-end tool.
 """
 
+from unittest import Mock
+
 from .dummy_comm import DummyCommunicator
 try:
     from mpi4py import MPI
     comm = MPI.COMM_WORLD
     have_mpi = True
-    mpi_sum = MPI.SUM
 except ImportError:
     comm = DummyCommunicator()
     have_mpi = False
-    mpi_sum = lambda x: x
+    MPI = Mock()
 
 from .nuclide import *
 from .chain import *
diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index 8c47706b75..bdee89b6ca 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -10,7 +10,7 @@ from warnings import warn
 import numpy as np
 import h5py
 
-from . import comm, have_mpi, mpi_sum
+from . import comm, have_mpi, MPI
 from .reaction_rates import ReactionRates
 
 _VERSION_RESULTS = (1, 0)
@@ -452,7 +452,7 @@ class Results(object):
         results.power = power
         results.proc_time = proc_time
         if results.proc_time is not None:
-            results.proc_time = comm.reduce(proc_time, op=mpi_sum)
+            results.proc_time = comm.reduce(proc_time, op=MPI.SUM)
 
         results.export_to_hdf5("depletion_results.h5", step_ind)
 

From 2efc7e8b71606b59708868bbc0e1267856fd6297 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Thu, 27 Jun 2019 10:29:56 -0500
Subject: [PATCH 16/19] Fix issue with vector proc_time in reading
 depletion_summary

Previous changes made proc_time a vector, not an (N, 1) array.
The reading of this data was not updated accordingly, until now
---
 openmc/deplete/__init__.py     | 2 +-
 openmc/deplete/results.py      | 2 +-
 openmc/deplete/results_list.py | 2 +-
 3 files changed, 3 insertions(+), 3 deletions(-)

diff --git a/openmc/deplete/__init__.py b/openmc/deplete/__init__.py
index 0c3a85d071..b5eb2dec2e 100644
--- a/openmc/deplete/__init__.py
+++ b/openmc/deplete/__init__.py
@@ -8,7 +8,7 @@ from sys import exit
 
 from h5py import get_config
 
-from unittest import Mock
+from unittest.mock import Mock
 
 from .dummy_comm import DummyCommunicator
 
diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py
index bdee89b6ca..c61664d42f 100644
--- a/openmc/deplete/results.py
+++ b/openmc/deplete/results.py
@@ -359,7 +359,7 @@ class Results(object):
         if "depletion time" in handle:
             proc_time_dset = handle["/depletion time"]
             if step < proc_time_dset.shape[0]:
-                results.proc_time = proc_time_dset[step, :]
+                results.proc_time = proc_time_dset[step]
 
         if results.proc_time is None:
             results.proc_time = np.array([np.nan])
diff --git a/openmc/deplete/results_list.py b/openmc/deplete/results_list.py
index 58dcf0c59d..c5dd1a8926 100644
--- a/openmc/deplete/results_list.py
+++ b/openmc/deplete/results_list.py
@@ -116,5 +116,5 @@ class ResultsList(list):
         """
         times = np.empty((len(self) - 1, 1))
         for ix, res in enumerate(self[:-1]):
-            times[ix] = res.proc_time[0]
+            times[ix] = res.proc_time
         return times

From e139310894ae2871917636e6ed5354b64fc2378d Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Tue, 2 Jul 2019 14:22:36 -0500
Subject: [PATCH 17/19] Add predictor-specific logic in get_depletion_time

The predictor saves BOS number densities, eigenvalues, and
other results before the depletion event, unlike all
other schemes, which have at least one depletion event between
the BOS transport solution and the depletion solve.
This causes the first value in the resulting depletion process
time vector to be zero, and the last value to be non-zero,
opposite to other schemes.
---
 openmc/deplete/integrator/predictor.py |  2 +-
 openmc/deplete/results_list.py         | 10 +++++++++-
 2 files changed, 10 insertions(+), 2 deletions(-)

diff --git a/openmc/deplete/integrator/predictor.py b/openmc/deplete/integrator/predictor.py
index a32ef12b9d..b30484c695 100644
--- a/openmc/deplete/integrator/predictor.py
+++ b/openmc/deplete/integrator/predictor.py
@@ -103,4 +103,4 @@ def predictor(operator, timesteps, power=None, power_density=None,
         op_results = [operator(x[0], power[-1])]
 
         # Create results, write to disk
-        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps))
+        Results.save(operator, x, op_results, [t, t], p, i_res + len(timesteps), proc_time)
diff --git a/openmc/deplete/results_list.py b/openmc/deplete/results_list.py
index c5dd1a8926..39cce893d1 100644
--- a/openmc/deplete/results_list.py
+++ b/openmc/deplete/results_list.py
@@ -115,6 +115,14 @@ class ResultsList(list):
 
         """
         times = np.empty((len(self) - 1, 1))
-        for ix, res in enumerate(self[:-1]):
+        # Need special logic because the predictor
+        # writes EOS values for step i as BOS values
+        # for step i+1
+        # The first proc_time may be zero
+        if self[0].proc_time > 0.0:
+            items = self[:-1]
+        else:
+            items = self[1:]
+        for ix, res in enumerate(items):
             times[ix] = res.proc_time
         return times

From 6732f5908d8e105e7f65338729c9b45231892146 Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 5 Jul 2019 08:02:32 -0500
Subject: [PATCH 18/19] Add depletion time to io formats

---
 docs/source/io_formats/depletion_results.rst | 3 +++
 1 file changed, 3 insertions(+)

diff --git a/docs/source/io_formats/depletion_results.rst b/docs/source/io_formats/depletion_results.rst
index 4f92384a0a..3c782b1d90 100644
--- a/docs/source/io_formats/depletion_results.rst
+++ b/docs/source/io_formats/depletion_results.rst
@@ -25,6 +25,9 @@ The current version of the depletion results file format is 1.0.
              nuclides, number of reactions).
            - **time** (*double[][2]*) -- Time in [s] at beginning/end of each
              step.
+           - **depletion time** (*double[]*) -- Average process time in [s] 
+             spent depleting a material across all burnable materials and,
+             if applicable, MPI processes.
 
 **/materials/<id>/**
 

From 6d1b6fc23074f2a7b189d7e3a2afd94ab1d48e4a Mon Sep 17 00:00:00 2001
From: Andrew Johnson <drewej@protonmail.com>
Date: Fri, 5 Jul 2019 13:45:26 -0500
Subject: [PATCH 19/19] Return vector from ResultsList.get_depletion_time

---
 openmc/deplete/integrator/epc_rk4.py        | 2 +-
 openmc/deplete/integrator/leqi.py           | 2 +-
 openmc/deplete/results_list.py              | 9 ++++++++-
 tests/unit_tests/test_deplete_cecm.py       | 2 +-
 tests/unit_tests/test_deplete_celi.py       | 2 +-
 tests/unit_tests/test_deplete_cf4.py        | 2 +-
 tests/unit_tests/test_deplete_epc_rk4.py    | 2 +-
 tests/unit_tests/test_deplete_integrator.py | 4 ++--
 tests/unit_tests/test_deplete_leqi.py       | 2 +-
 tests/unit_tests/test_deplete_predictor.py  | 2 +-
 tests/unit_tests/test_deplete_si_celi.py    | 8 +-------
 tests/unit_tests/test_deplete_si_leqi.py    | 2 +-
 12 files changed, 20 insertions(+), 19 deletions(-)

diff --git a/openmc/deplete/integrator/epc_rk4.py b/openmc/deplete/integrator/epc_rk4.py
index 1f0fbe2a84..373c323e42 100644
--- a/openmc/deplete/integrator/epc_rk4.py
+++ b/openmc/deplete/integrator/epc_rk4.py
@@ -133,7 +133,7 @@ def epc_rk4(operator, timesteps, power=None, power_density=None, print_out=True)
             # Create results, write to disk
             Results.save(
                 operator, x, op_results, [t, t + dt], p, i_res + i,
-                sum((time_1, time_2, time_3, time_4)))
+                time_1 + time_2 + time_3 + time_4)
 
             # Advance time, update vector
             t += dt
diff --git a/openmc/deplete/integrator/leqi.py b/openmc/deplete/integrator/leqi.py
index 465139abf4..19b1025e5c 100644
--- a/openmc/deplete/integrator/leqi.py
+++ b/openmc/deplete/integrator/leqi.py
@@ -154,7 +154,7 @@ def leqi(operator, timesteps, power=None, power_density=None, print_out=True):
             # Create results, write to disk
             Results.save(
                 operator, x, op_results, [t, t+dt], p, i_res+i,
-                sum((time_1, time_2, time_3, time_4)))
+                time_1 + time_2 + time_3 + time_4)
 
             # update results
             op_res_last = copy.deepcopy(op_results[0])
diff --git a/openmc/deplete/results_list.py b/openmc/deplete/results_list.py
index 9f256eb3db..79203a6318 100644
--- a/openmc/deplete/results_list.py
+++ b/openmc/deplete/results_list.py
@@ -115,8 +115,15 @@ class ResultsList(list):
             like :meth:`get_eigenvalues`, because no depletion
             is performed at the final transport stage
 
+        Returns
+        -------
+
+        times : :class:`numpy.ndarray`
+            Vector of average time to deplete a single material
+            across all processes and materials.
+
         """
-        times = np.empty((len(self) - 1, 1))
+        times = np.empty(len(self) - 1)
         # Need special logic because the predictor
         # writes EOS values for step i as BOS values
         # for step i+1
diff --git a/tests/unit_tests/test_deplete_cecm.py b/tests/unit_tests/test_deplete_cecm.py
index d143658af3..dd8b769fe6 100644
--- a/tests/unit_tests/test_deplete_cecm.py
+++ b/tests/unit_tests/test_deplete_cecm.py
@@ -39,5 +39,5 @@ def test_cecm(run_in_tmpdir):
     # Test structure of depletion time dataset
 
     dep_time = res.get_depletion_time()
-    assert dep_time.shape == (len(dt), 1)
+    assert dep_time.shape == (len(dt), )
     assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_celi.py b/tests/unit_tests/test_deplete_celi.py
index bb479f6cc2..47d3319fdb 100644
--- a/tests/unit_tests/test_deplete_celi.py
+++ b/tests/unit_tests/test_deplete_celi.py
@@ -39,5 +39,5 @@ def test_celi(run_in_tmpdir):
     # Test structure of depletion time dataset
 
     dep_time = res.get_depletion_time()
-    assert dep_time.shape == (len(dt), 1)
+    assert dep_time.shape == (len(dt), )
     assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_cf4.py b/tests/unit_tests/test_deplete_cf4.py
index 5b7d08c229..0c6199c3e5 100644
--- a/tests/unit_tests/test_deplete_cf4.py
+++ b/tests/unit_tests/test_deplete_cf4.py
@@ -39,5 +39,5 @@ def test_cf4(run_in_tmpdir):
     # Test structure of depletion time dataset
 
     dep_time = res.get_depletion_time()
-    assert dep_time.shape == (len(dt), 1)
+    assert dep_time.shape == (len(dt), )
     assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_epc_rk4.py b/tests/unit_tests/test_deplete_epc_rk4.py
index 0623f9e6e8..ea687bfd5d 100644
--- a/tests/unit_tests/test_deplete_epc_rk4.py
+++ b/tests/unit_tests/test_deplete_epc_rk4.py
@@ -39,5 +39,5 @@ def test_epc_rk4(run_in_tmpdir):
     # Test structure of depletion time dataset
 
     dep_time = res.get_depletion_time()
-    assert dep_time.shape == (len(dt), 1)
+    assert dep_time.shape == (len(dt), )
     assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_integrator.py b/tests/unit_tests/test_deplete_integrator.py
index 451b703a30..40090264ae 100644
--- a/tests/unit_tests/test_deplete_integrator.py
+++ b/tests/unit_tests/test_deplete_integrator.py
@@ -76,8 +76,8 @@ def test_results_save(run_in_tmpdir):
 
     op_result1 = [OperatorResult(k, rates) for k, rates in zip(eigvl1, rate1)]
     op_result2 = [OperatorResult(k, rates) for k, rates in zip(eigvl2, rate2)]
-    Results.save(op, x1, op_result1, t1, 0, 0, None)
-    Results.save(op, x2, op_result2, t2, 0, 1, None)
+    Results.save(op, x1, op_result1, t1, 0, 0)
+    Results.save(op, x2, op_result2, t2, 0, 1)
 
     # Load the files
     res = ResultsList("depletion_results.h5")
diff --git a/tests/unit_tests/test_deplete_leqi.py b/tests/unit_tests/test_deplete_leqi.py
index 97a841f6ce..a26329cf4a 100644
--- a/tests/unit_tests/test_deplete_leqi.py
+++ b/tests/unit_tests/test_deplete_leqi.py
@@ -39,5 +39,5 @@ def test_leqi(run_in_tmpdir):
     # Test structure of depletion time dataset
 
     dep_time = res.get_depletion_time()
-    assert dep_time.shape == (len(dt), 1)
+    assert dep_time.shape == (len(dt), )
     assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_predictor.py b/tests/unit_tests/test_deplete_predictor.py
index 80b060cc9e..6af6f8bca4 100644
--- a/tests/unit_tests/test_deplete_predictor.py
+++ b/tests/unit_tests/test_deplete_predictor.py
@@ -39,5 +39,5 @@ def test_predictor(run_in_tmpdir):
     # Test structure of depletion time dataset
 
     dep_time = res.get_depletion_time()
-    assert dep_time.shape == (len(dt), 1)
+    assert dep_time.shape == (len(dt), )
     assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_si_celi.py b/tests/unit_tests/test_deplete_si_celi.py
index c064172976..71679439c5 100644
--- a/tests/unit_tests/test_deplete_si_celi.py
+++ b/tests/unit_tests/test_deplete_si_celi.py
@@ -39,11 +39,5 @@ def test_si_celi(run_in_tmpdir):
     # Test structure of depletion time dataset
 
     dep_time = res.get_depletion_time()
-    assert dep_time.shape == (len(dt), 1)
-    assert all(dep_time > 0)
-
-    # Test structure of depletion time dataset
-
-    dep_time = res.get_depletion_time()
-    assert dep_time.shape == (len(dt), 1)
+    assert dep_time.shape == (len(dt), )
     assert all(dep_time > 0)
diff --git a/tests/unit_tests/test_deplete_si_leqi.py b/tests/unit_tests/test_deplete_si_leqi.py
index 1003a73161..1396065c0f 100644
--- a/tests/unit_tests/test_deplete_si_leqi.py
+++ b/tests/unit_tests/test_deplete_si_leqi.py
@@ -39,5 +39,5 @@ def test_si_leqi(run_in_tmpdir):
     # Test structure of depletion time dataset
 
     dep_time = res.get_depletion_time()
-    assert dep_time.shape == (len(dt), 1)
+    assert dep_time.shape == (len(dt), )
     assert all(dep_time > 0)