depletion.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Pincell Depletion\n",
    "This notebook is intended to introduce the reader to the depletion interface contained in OpenMC. It is recommended that you are moderately familiar with building models using the OpenMC Python API. The earlier examples are excellent starting points, as this notebook will not focus heavily on model building.\n",
    "\n",
    "If you have a real power reactor, the fuel composition is constantly changing as fission events produce energy, remove some fissile isotopes, and produce fission products. Other reactions, like $(n, \\alpha)$ and $(n, \\gamma)$ will alter the composition as well. Furthermore, some nuclides undergo spontaneous decay with widely ranging frequencies. Depletion is the process of modeling this behavior.\n",
    "\n",
    "In this notebook, we will model a simple fuel pin in an infinite lattice using the Python API. We will then build and examine some of the necessary components for performing depletion analysis. Then, we will use the depletion interface in OpenMC to simulate the fuel pin producing power over several months. Lastly, we will wrap up with some helpful tips to improve the fidelity of depletion simulations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import math\n",
    "import openmc"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Build the Geometry\n",
    "\n",
    "Much of this section is borrowed from the \"Modeling a Pin-Cell\" example. If you find yourself not understanding some aspects of this section, feel free to refer to that example, as some details may be glossed over for brevity.\n",
    "\n",
    "First, we will create our fuel, cladding, and water materials to represent a typical PWR."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "fuel = openmc.Material(name=\"uo2\")\n",
    "fuel.add_element(\"U\", 1, percent_type=\"ao\", enrichment=4.25)\n",
    "fuel.add_element(\"O\", 2)\n",
    "fuel.set_density(\"g/cc\", 10.4)\n",
    "\n",
    "clad = openmc.Material(name=\"clad\")\n",
    "clad.add_element(\"Zr\", 1)\n",
    "clad.set_density(\"g/cc\", 6)\n",
    "\n",
    "water = openmc.Material(name=\"water\")\n",
    "water.add_element(\"O\", 1)\n",
    "water.add_element(\"H\", 2)\n",
    "water.set_density(\"g/cc\", 1.0)\n",
    "water.add_s_alpha_beta(\"c_H_in_H2O\")\n",
    "materials = openmc.Materials([fuel, clad, water])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, we are going to use the `openmc.model.pin` function to build our pin cell. The `pin` function anticipates concentric cylinders and materials to fill the inner regions. One additional material is needed than the number of cylinders to cover the domain outside the final ring. \n",
    "\n",
    "To do this, we define two radii for the outer radius of our fuel pin, and the outer radius of the cladding."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "radii = [0.42, 0.45]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using these radii, we define concentric `ZCylinder` objects. So long as the cylinders are concentric and increasing in radius, any orientation can be used. We also take advantage of the fact that the `openmc.Materials` object is a subclass of the `list` object to assign materials to the regions defined by the surfaces."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "pin_surfaces = [openmc.ZCylinder(r=r) for r in radii]\n",
    "pin_univ = openmc.model.pin(pin_surfaces, materials)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The first material, in our case `fuel`, is placed inside the first cylinder in the inner-most region. The second material, `clad`, fills the space between our cylinders, while `water` is placed outside the last ring. The `pin` function returns an `openmc.Universe` object, and has some additional features we will mention later. Finally, we need to place the fuel pin universe in a bounding cell."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "bound_box = openmc.model.RectangularPrism(1.24, 1.24, boundary_type=\"reflective\")\n",
    "root_cell = openmc.Cell(fill=pin_univ, region=-bound_box)\n",
    "geometry = openmc.Geometry([root_cell])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To ensure our geometry looks right, let's plot it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x7f8a97aecf70>"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPUAAADnCAYAAADGrxD1AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAARIUlEQVR4nO3db6xk9V3H8feH3WyNKSx/pICAQONqQo2h9bqmD6AatpEayjaptqwhsgazUeSBaUCXkDSGPuCvRZLWKMEmSLQLJaluulQstEQelMo2InExsBdiw1IKFpCkIS2S/fpgzt07d/acM2dmzjlzzu98XsmGmTu/O+fHct/8fufce2cUEZhZOk5Y9gTMrF6O2iwxjtosMY7aLDGO2iwxjtosMZuXPYEi2qrgfcuehVlHrfLDiDg976HORs374IS7lj0Js246+nG+V/SYt99miXHUZolx1GaJcdRmiXHUZolx1GaJcdRmienu96mtPXX9Sr1qeh5biKMeipJwb9i1WsPTB3d+eVvxAAffGkedopyA6wi3jFDhMQqDd+iNcNQpmIi46YBnlRd8buiOvBaOuq/GQu5axFVMhn5c5A58bo66T3oecpnxyB34Yhx117UScrcufzvwxdQStaTLgLuBTcC9EXFrwbhPAg8BvxoRB+s4drKyzuoLuTjcSw7cuPizBzxx+S0lI+arsTBwx11Ii77ut6RNwPPAR4EjwFPAroh4dmLcicABYAtw3bSotU0xyN+nri3mjf9d6wh3XvnBz1+l44ajH+e7EbGS91gdK/V2YDUiXgSQtA/YCTw7Me5zwG3ADTUcMz0Lx9ydiCdJG+ezaORrq7fjzldH1GcDL43dPwL82vgASR8Czo2IA5Ic9biFYl4PuUsRTzM98mqFOu58jV8ok3QC8Hlgd4Wxe4A9AOS++lJC5o65nyGXGY98nsAd90Z1nFN/GPjziPjN7P6NABFxS3Z/K/AC8KPsU84E3gCuKDuvTvacesGYUwm5io2BV690CHGXnVPXEfVmRhfKLgVeZnSh7Hcj4lDB+MeB6wd5oSxmjTm9VXke827Pj8WdYNiNXiiLiHclXQc8wuhbWl+KiEOSbgYORsT+RY/RezOvzsNblcsUb8/LaxXi+l2Hk1+1Jy28UjclmZV6ptXZMVc169Y8tVW70e13U3of9Uyrs2Oe1yxxp3Su7ajbNuPq7JgXtx73MFbtsqj9ckZ1CmYIejTYQddDgou/diPH/iOUjc3OtSsM7SWv1HWZKWZvtZs015a8Z6u2t99NmyFox9yeqlvyPobt7XeTKgXtrfYyVN2Sb9iOJ8Ar9bwqX912zF1QZdXu09Vxb7/r5u12L6W0Hff2u07ebvfWULbjXqlnUTFox9x9M23HO7hie6Wug4NOysZVu2BMT1dsR12Fg05SqmE76mkcdNJSDNtRl3HQg5Ba2I66iIMelJTCdtR5HPQgzRR2hznquTjoVFUJG6Y/vEyOetLUVdpBp25a2F3fhjvqcQ7aMn0O21GvcdA2oa9hO2pw0Faoctgd4qinctBDV+niWYdWa0dduko7aBtZDzvnsY5tw4cd9czvmGHW/W34sKMu5VXaNurLNny4UXvbbXMoC7sr2/BhRu2gbQGVzq+XaJhRm9Wim9vw4UXtVdpq0OVt+LCidtBWo65uw4cVtVkjurUNryVqSZdJek7SqqS9OY9/RtKzkp6R9Jik8+o47ky8SlsDurhaLxy1pE3AF4GPARcCuyRdODHs34GViPhl4CHg9kWPa9Yt3Vmt61iptwOrEfFiRLwD7AN2jg+IiG9FxNvZ3SeBc2o4bnVepa1BlS6ataiOqM8GXhq7fyT7WJFrgK/XcNwaOGirR9k2HGh1tW71Qpmkq4AV4I6Cx/dIOijpIG/VdFD/fLe1avmrdR1RvwycO3b/nOxjG0jaAdwEXBERP8l7ooi4JyJWImKFrTXMrJRXaatXV1brOqJ+Ctgm6QJJW4Argf3jAyR9EPgbRkG/VsMxq/EqbUux3NV64agj4l3gOuAR4L+AByPikKSbJV2RDbsDeC/wFUlPS9pf8HQt8SptzejCar25jieJiIeBhyc+9tmx2zvqOM5sk/IqbcsUTL5d5tpqfee+bY0eeYA/UeZV2pq17NU6zai9SltHtXFunWbUZp2wnF/TGljU3npbO5a5BU8vam+9rVPa//ZWelEX8ipt7Zq6WjdkQFGbdUxDW/C0ovbW2zqp3S14WlEX8tbblmMZW/CBRG3WUQ1swdOJ2ltv67T2tuDpRF3IW29brra34AOI2mxYHLXZstV8Xp1G1D6ftl5o57w6jagL+XzauqHN8+rEozYbHkdtlhhHbdaqgqtiNV4s63/UhRfJfD5t3VJ0Xl33xbL+R21mGzhqs8Q4arPEOGqzxDhqs8Q4arPWNfttrX5H7W9nWc+08W2tfkdtZsdx1GaJcdRmiXHUZompJWpJl0l6TtKqpL05j79H0gPZ49+RdH4dxzWz4y0ctaRNwBeBjwEXArskXTgx7BrgzYj4eeAu4LZFj2tm+epYqbcDqxHxYkS8A+wDdk6M2Qncl91+CLhUkmo4tplNqCPqs4GXxu4fyT6WOyYi3gXeAk6r4dhmNmHzsicwTtIeYA8Apy93LmZ9VcdK/TJw7tj9c7KP5Y6RtBnYCrw++UQRcU9ErETECltrmJnZANUR9VPANkkXSNoCXAnsnxizH7g6u/3bwDcjoqE38jQbtoW33xHxrqTrgEeATcCXIuKQpJuBgxGxH/hb4H5Jq8AbjMI3swbUck4dEQ8DD0987LNjt38M/E4dxzKzcv6JMrPEOGqzxDhqs8T0PuooeLkIX1u3oep31II7v7wt94EnLr+l9emYTRMBT1x+6/EfJ7hzX97X8uz6HbWZHcdRmyXGUZslxlGbJcZRmyUmiaj9bS2zdf2P2t/Wsp5o49tZkELUZraBozZLTPJR+7zahiaZqPMvlvm82rqh6Hy6CWlEXXixzKzb6r5IBqlEbWbHDCJqn1fbkCQVtc+rrYvaPJ+GlKL2ebX1TBPn05BS1FN4C25DkVzU3oJbl7S99YbUovYW3Hqiqa03pBb1FN6C2xAkGbW34NYFy9h6Q4pRT9mCe7W2ZWty6w0pRl3Kq7W1Y1mrNCQcddGroZilLs2ovQW3jmp66w2pRp3xBTNbhmVuvWHBqCWdKukbkg5n/zwlZ8xFkr4t6ZCkZyR9epFjVp+cV2vrljZWaVh8pd4LPBYR24DHsvuT3gZ+LyI+AFwG/KWkkxc8bmVera1Ny16lYfGodwL3ZbfvAz4xOSAino+Iw9nt7wOvAacveNxqvFpbR7S1SsPiUZ8REa9kt38AnFE2WNJ2YAvwwoLHnYlXa2tDF1ZpqBC1pEcl/WfOn53j4yIioPj7SJLOAu4Hfj8ijhaM2SPpoKSDvDXjv0nhgb1a23K1uUoDbJ42ICJ2FD0m6VVJZ0XEK1m0rxWMOwk4ANwUEU+WHOse4B4AbVOtuQWB0OTMeOLyW7j4azeiyYfMZtCVVRoW337vB67Obl8N/NPkAElbgK8CfxcRDy14vPmUrtbehttiyoJue5WGxaO+FfiopMPAjuw+klYk3ZuN+RRwCbBb0tPZn4sWPO5cyn7KzNtwS4Wio1/N2qY44a6anzTghl2rhQ9ecuDGmg9oqVvWKn3043w3IlbyHkv6J8ryeLW21A0r6uzcuuxbXA7bquraufSaYUUNvmhmtehq0DDEqDPehluqhhm1t+G2gC6v0jDUqKHSNtxh26SuBw1DjjpTvA33+bVt1KWfGisz7KhLt+EjXq2tiq6s0jD0qMHbcKukD9vuNY46M20b7rCHqy/b7jWOGipswx32UE0LumurNDjqdQ7bJvQxaHDUGzlsy/Q1aHDUx5v6zpkOO3V9DhocdaHyd/hw2Knq20WxPI46T4XvXzvs9FQJuuurNDjqYg57UFIJGhx1OYc9CCkFDY56OoedtNSCBkddjcNOUopBg6OubsawHXd3rf33STFoGNqridYh4Ppdh3PeGOD4gX6TgO6p+i2rrgftVxOt09iK7e14v1Rdnbse9DReqRdRadUerdiAV+0lWfsST2m7XbZSO+pFeTveaalstyd5+92kShfQRgN9Ea09VS+GQf+CnsYrdV2yFRuovGqDt+R1q7rVhvWf7+9j0N5+t6nydnw02Fvy+szyyxh9X529/W5T5avjo8Heki9u1q1234Oexit1k+ZYtcFb8qpm2WpD/1fncd5+L9NM59qjT3Dc5eaJGfp57lyksaglnQo8AJwP/DfwqYh4s2DsScCzwD9GxHVTnzuVqNfMtGqPPmEtbnDg41+ms7yIQUqr87gmo74deCMibpW0FzglIv6sYOzdwOnZ+OFFDXOs2uufONTVe9ZV+djnJbg6j2sy6ueAX4+IVySdBTweEb+YM+5XgBuAfwZWBhv1mhrihnQDn3dVhvRjXtNk1P8bESdntwW8uXZ/bMwJwDeBq4AdOOp1c8c9+uSUAl8kZBhOzGvKot487ZMlPQqcmfPQTeN3IiIk5f0f4lrg4Yg4oilfeZL2AHuA0UY9dcq+COeKW2Nf/BsDh+5HPrmWzPtif0OLuYrGt9+S/h64GDgKvBfYAvxVROwtfe4hrNSTFlq5Nz5R1yKvK+JjzzfwmJvcft8BvD52oezUiPjTkvG78fZ7utriXn/CycjH1RV82ZdSXS+7O/SY1zQZ9WnAg8DPAd9j9C2tNyStAH8YEX8wMX43jrq6sbihrsCPP0hZ8LNo6vWyx38yb+gxr/EPn6SglcC7wyGXW+hCmXWExr64Ew3cIdfDUfdRSeCjh/sR+eQvvDjkejjqvtNEDDmRj4YtN/S831hzxM1w1KmZjBwKQ1//lHqCL/tVUwfcHkc9BHmhr5kS/Cwcbjc46qErC956ya98YpYYR22WGEdtlhhHbZYYR22WGEdtlhhHbZYYR22WGEdtlpju/j619D+MXnihLj8D/LDG56ub57eYoc3vvIjIfSW/zkZdN0kHi36pvAs8v8V4fuu8/TZLjKM2S8yQor5n2ROYwvNbjOeXGcw5tdlQDGmlNhuEZKOWdKqkb0g6nP3zlJKxJ0k6IukLXZqfpIskfVvSIUnPSPp0C/O6TNJzklazN2iYfPw9kh7IHv+OpPObntOM8/uMpGezv6/HJJ3XpfmNjfukpMheI79eEZHkH+B2YG92ey9wW8nYu4F/AL7QpfkBvwBsy27/LPAKcHKDc9oEvAC8n9HbI/0HcOHEmGuBv85uXwk80OLfWZX5/Qbw09ntP+ra/LJxJwL/CjzJ6M0tap1Hsis1sBO4L7t9H/CJvEHZ2+yeAfxLO9M6Zur8IuL5iDic3f4+8BrNvnXgdmA1Il6MiHeAfdk8x43P+yHgUk1758MW5xcR34qIt7O7TwLntDS3SvPLfA64DfhxE5NIOeozIuKV7PYPGIW7QfY2u38BXN/mxDJT5zdO0nZG//d/ocE5nQ28NHb/SPax3DER8S7wFnBag3PKPXYmb37jrgG+3uiMNpo6P0kfAs6NiANNTaLXLzzY5tvsLml+a89zFnA/cHVEHK13lmmSdBWwAnxk2XNZky0inwd2N3mcXkcdETuKHpP0qqSzYv1tdl/LGfZh4GJJ15K9za6kH8WUt9ltcX5IOgk4ANwUEU/WMa8SLwPnjt0/J/tY3pgjkjYDW4HXG57X5LHX5M0PSTsY/Y/zIxHxk5bmBtPndyLwS8Dj2SJyJrBf0hURcbC2WbR1EaHtP8AdbLwQdfuU8btp90LZ1Pkx2m4/BvxJS3PaDLwIXMD6hZ4PTIz5YzZeKHuwxb+zKvP7IKNTlG1L+JqbOr+J8Y/TwIWyVv+lW/4LPi0L4jDwKKP3zobRluzenPFtRz11fsBVwP8BT4/9uajhef0W8HwWxk3Zx24Grshu/xTwFWAV+Dfg/S3/d502v0eBV8f+vvZ3aX4TYxuJ2j9RZpaYlK9+mw2SozZLjKM2S4yjNkuMozZLjKM2S4yjNkuMozZLzP8D3rsAXKYwN8gAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 258.065x264.901 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "geometry.root_universe.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lastly we construct our settings. For the sake of time, a relatively low number of particles will be used."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "settings = openmc.Settings()\n",
    "settings.particles = 1000\n",
    "settings.inactive = 10\n",
    "settings.batches = 50"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The depletion interface relies on `OpenMC` to perform the transport simulation and obtain reaction rates and other important information. Normally, we would need to export XML files before running OpenMC, but the depletion interface takes care of this for us."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "However, we must first add one bit of information: the volume of our fuel. In order to translate the reaction rates obtained by `openmc` to meaningful units for depletion, we have to normalize them to a correct power. This requires us to know, or be able to calculate, how much fuel is in our problem. Correctly setting the volumes is a critical step, and can lead to incorrect answers, as the fuel is over- or under-depleted due to poor normalization.\n",
    "\n",
    "For our problem, we can assign the \"volume\" to be the cross-sectional area of our fuel. This is identical to modeling our fuel pin inside a box with height of 1 cm."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "fuel.volume = math.pi * radii[0] ** 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setting up for depletion\n",
    "\n",
    "The OpenMC depletion interface can be accessed from the `openmc.deplete` module, and has a variety of classes that will help us."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "import openmc.deplete"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In order to run the depletion calculation we need the following information:\n",
    "\n",
    "1. Nuclide decay, fission yield, and reaction data\n",
    "2. Operational power or power density\n",
    "3. Desired depletion schedule\n",
    "4. Desired time integration scheme\n",
    "\n",
    "The first item is necessary to determine the paths by which nuclides transmute over the depletion simulation. This includes spontaneous decay, fission product yield distributions, and nuclides produced through neutron-reactions. For example,\n",
    "* Te129 decays to I129 with a half life of ~70 minutes\n",
    "* A fission event for U-235 produces fission products like Xe135 according to a distribution\n",
    "* For thermal problems, Am241 will produce metastable Am242 about 8% of the time during an $(n,\\gamma)$ reaction. The other 92% of capture reactions will produce ground state Am242\n",
    "\n",
    "These data are often distributed with other nuclear data, like incident neutron cross sections with ENDF/B-VII.\n",
    "OpenMC uses the [`openmc.deplete.Chain`](https://docs.openmc.org/en/latest/pythonapi/generated/openmc.deplete.Chain.html#openmc.deplete.Chain) to collect represent the various decay and transmutation pathways in a single object.\n",
    "While a complete `Chain` can be created using nuclear data files, users may prefer to download pre-generated XML-representations instead.\n",
    "Such files can be found at https://openmc.org/depletion-chains/ and include full and compressed chains, with capture branching ratios derived using PWR- or SFR-spectra.\n",
    "\n",
    "For this problem, we will be using a much smaller depletion chain that contains very few nuclides. In a realistic problem, over 1000 isotopes may be included in the depletion chain."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "OrderedDict([('I135', 0),\n",
       "             ('Xe135', 1),\n",
       "             ('Xe136', 2),\n",
       "             ('Cs135', 3),\n",
       "             ('Gd157', 4),\n",
       "             ('Gd156', 5),\n",
       "             ('U234', 6),\n",
       "             ('U235', 7),\n",
       "             ('U238', 8)])"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "chain = openmc.deplete.Chain.from_xml(\"./chain_simple.xml\")\n",
    "chain.nuclide_dict"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The primary entry point for depletion is the `openmc.deplete.Operator`. It relies on the `openmc.deplete.Chain` and helper classes to run `openmc`, retrieve and normalize reaction rates, and other perform other tasks. For a thorough description, please see the full API documentation."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will create our Operator using the geometry and settings from above, and our simple chain file. The materials are read in automatically using the `materials.xml` file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = openmc.Model(geometry=geometry, settings=settings)\n",
    "operator = openmc.deplete.CoupledOperator(model, \"./chain_simple.xml\")"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "We will then simulate our fuel pin operating at linear power of 174 W/cm, or 174 W given a unit height for our problem."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "power = 174"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For this problem, we will take depletion step sizes of 30 days, and instruct OpenMC to re-run a transport simulation every 30 days until we have modeled the problem over a six month cycle. The depletion interface expects the time to be given in seconds, so we will have to convert. Note that these values are not cumulative."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "time_steps = [30] * 6"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And lastly, we will use the basic predictor, or forward Euler, time integration scheme. Other, more advanced methods are provided to the user through `openmc.deplete`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "integrator = openmc.deplete.PredictorIntegrator(operator, time_steps, power, timestep_units='d')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To perform the simulation, we use the `integrate` method, and let `openmc` take care of the rest."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                                %%%%%%%%%%%%%%%\n",
      "                           %%%%%%%%%%%%%%%%%%%%%%%%\n",
      "                        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\n",
      "                      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\n",
      "                    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\n",
      "                   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\n",
      "                                    %%%%%%%%%%%%%%%%%%%%%%%%\n",
      "                                     %%%%%%%%%%%%%%%%%%%%%%%%\n",
      "                 ###############      %%%%%%%%%%%%%%%%%%%%%%%%\n",
      "                ##################     %%%%%%%%%%%%%%%%%%%%%%%\n",
      "                ###################     %%%%%%%%%%%%%%%%%%%%%%%\n",
      "                ####################     %%%%%%%%%%%%%%%%%%%%%%\n",
      "                #####################     %%%%%%%%%%%%%%%%%%%%%\n",
      "                ######################     %%%%%%%%%%%%%%%%%%%%\n",
      "                #######################     %%%%%%%%%%%%%%%%%%\n",
      "                 #######################     %%%%%%%%%%%%%%%%%\n",
      "                 ######################     %%%%%%%%%%%%%%%%%\n",
      "                  ####################     %%%%%%%%%%%%%%%%%\n",
      "                    #################     %%%%%%%%%%%%%%%%%\n",
      "                     ###############     %%%%%%%%%%%%%%%%\n",
      "                       ############     %%%%%%%%%%%%%%%\n",
      "                          ########     %%%%%%%%%%%%%%\n",
      "                                      %%%%%%%%%%%\n",
      "\n",
      "                 | The OpenMC Monte Carlo Code\n",
      "       Copyright | 2011-2022 MIT, UChicago Argonne LLC, and contributors\n",
      "         License | https://docs.openmc.org/en/latest/license.html\n",
      "         Version | 0.13.1\n",
      "        Git SHA1 | 33bc948f4b855c037975f16d16091fe4ecd12de3\n",
      "       Date/Time | 2022-10-05 19:02:05\n",
      "  OpenMP Threads | 2\n",
      "\n",
      " Reading settings XML file...\n",
      " Reading cross sections XML file...\n",
      " Reading materials XML file...\n",
      " Reading geometry XML file...\n",
      " Reading U234 from /home/pshriwise/data/xs/openmc/nndc_hdf5/U234.h5\n",
      " Reading U235 from /home/pshriwise/data/xs/openmc/nndc_hdf5/U235.h5\n",
      " Reading U238 from /home/pshriwise/data/xs/openmc/nndc_hdf5/U238.h5\n",
      " Reading O16 from /home/pshriwise/data/xs/openmc/nndc_hdf5/O16.h5\n",
      " Reading O17 from /home/pshriwise/data/xs/openmc/nndc_hdf5/O17.h5\n",
      " Reading U236 from /home/pshriwise/data/xs/openmc/nndc_hdf5/U236.h5\n",
      " Reading Zr90 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Zr90.h5\n",
      " Reading Zr91 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Zr91.h5\n",
      " Reading Zr92 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Zr92.h5\n",
      " Reading Zr94 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Zr94.h5\n",
      " Reading Zr96 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Zr96.h5\n",
      " Reading H1 from /home/pshriwise/data/xs/openmc/nndc_hdf5/H1.h5\n",
      " Reading H2 from /home/pshriwise/data/xs/openmc/nndc_hdf5/H2.h5\n",
      " Reading c_H_in_H2O from /home/pshriwise/data/xs/openmc/nndc_hdf5/c_H_in_H2O.h5\n",
      " Minimum neutron data temperature: 294 K\n",
      " Maximum neutron data temperature: 294 K\n",
      " Preparing distributed cell instances...\n",
      " Reading plot XML file...\n",
      " Writing summary.h5 file...\n",
      "[openmc.deplete] t=0.0 s, dt=2592000 s, source=174\n",
      " Reading I135 from /home/pshriwise/data/xs/openmc/nndc_hdf5/I135.h5\n",
      " Reading Xe135 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Xe135.h5\n",
      " Reading Xe136 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Xe136.h5\n",
      " Reading Cs135 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Cs135.h5\n",
      " Reading Gd157 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Gd157.h5\n",
      " Reading Gd156 from /home/pshriwise/data/xs/openmc/nndc_hdf5/Gd156.h5\n",
      " Maximum neutron transport energy: 20000000 eV for U235\n",
      " Initializing source particles...\n",
      "\n",
      " ====================>     K EIGENVALUE SIMULATION     <====================\n",
      "\n",
      "  Bat./Gen.      k            Average k\n",
      "  =========   ========   ====================\n",
      "        1/1    1.49721\n",
      "        2/1    1.49226\n",
      "        3/1    1.44838\n",
      "        4/1    1.47163\n",
      "        5/1    1.41429\n",
      "        6/1    1.42100\n",
      "        7/1    1.43656\n",
      "        8/1    1.41211\n",
      "        9/1    1.45570\n",
      "       10/1    1.33919\n",
      "       11/1    1.47822\n",
      "       12/1    1.47425    1.47624 +/- 0.00198\n",
      "       13/1    1.39919    1.45056 +/- 0.02571\n",
      "       14/1    1.35785    1.42738 +/- 0.02946\n",
      "       15/1    1.38972    1.41985 +/- 0.02403\n",
      "       16/1    1.47812    1.42956 +/- 0.02189\n",
      "       17/1    1.50905    1.44092 +/- 0.02171\n",
      "       18/1    1.48299    1.44618 +/- 0.01952\n",
      "       19/1    1.49089    1.45114 +/- 0.01792\n",
      "       20/1    1.50543    1.45657 +/- 0.01692\n",
      "       21/1    1.41183    1.45251 +/- 0.01584\n",
      "       22/1    1.49487    1.45604 +/- 0.01488\n",
      "       23/1    1.45853    1.45623 +/- 0.01369\n",
      "       24/1    1.45628    1.45623 +/- 0.01268\n",
      "       25/1    1.37824    1.45103 +/- 0.01290\n",
      "       26/1    1.40298    1.44803 +/- 0.01243\n",
      "       27/1    1.48608    1.45027 +/- 0.01189\n",
      "       28/1    1.46773    1.45124 +/- 0.01125\n",
      "       29/1    1.46127    1.45177 +/- 0.01066\n",
      "       30/1    1.38587    1.44847 +/- 0.01063\n",
      "       31/1    1.42115    1.44717 +/- 0.01020\n",
      "       32/1    1.45934    1.44772 +/- 0.00974\n",
      "       33/1    1.40481    1.44586 +/- 0.00949\n",
      "       34/1    1.47602    1.44711 +/- 0.00917\n",
      "       35/1    1.47601    1.44827 +/- 0.00887\n",
      "       36/1    1.44354    1.44809 +/- 0.00853\n",
      "       37/1    1.54954    1.45185 +/- 0.00902\n",
      "       38/1    1.46876    1.45245 +/- 0.00872\n",
      "       39/1    1.51863    1.45473 +/- 0.00872\n",
      "       40/1    1.51529    1.45675 +/- 0.00866\n",
      "       41/1    1.39841    1.45487 +/- 0.00858\n",
      "       42/1    1.58314    1.45888 +/- 0.00923\n",
      "       43/1    1.51678    1.46063 +/- 0.00911\n",
      "       44/1    1.44328    1.46012 +/- 0.00886\n",
      "       45/1    1.42208    1.45903 +/- 0.00867\n",
      "       46/1    1.50118    1.46020 +/- 0.00850\n",
      "       47/1    1.50796    1.46150 +/- 0.00837\n",
      "       48/1    1.44889    1.46116 +/- 0.00816\n",
      "       49/1    1.48715    1.46183 +/- 0.00797\n",
      "       50/1    1.49559    1.46267 +/- 0.00782\n",
      " Creating state point statepoint.50.h5...\n",
      "\n",
      " =======================>     TIMING STATISTICS     <=======================\n",
      "\n",
      " Total time for initialization     = 1.1968e+00 seconds\n",
      "   Reading cross sections          = 1.1906e+00 seconds\n",
      " Total time in simulation          = 1.3040e+01 seconds\n",
      "   Time in transport only          = 1.3031e+01 seconds\n",
      "   Time in inactive batches        = 2.2501e+00 seconds\n",
      "   Time in active batches          = 1.0790e+01 seconds\n",
      "   Time synchronizing fission bank = 4.3167e-03 seconds\n",
      "     Sampling source sites         = 3.9236e-03 seconds\n",
      "     SEND/RECV source sites        = 3.6606e-04 seconds\n",
      "   Time accumulating tallies       = 2.4592e-04 seconds\n",
      "   Time writing statepoints        = 1.4237e-03 seconds\n",
      " Total time for finalization       = 5.2430e-05 seconds\n",
      " Total time elapsed                = 1.4251e+01 seconds\n",
      " Calculation Rate (inactive)       = 4444.21 particles/second\n",
      " Calculation Rate (active)         = 3707.27 particles/second\n",
      "\n",
      " ============================>     RESULTS     <============================\n",
      "\n",
      " k-effective (Collision)     = 1.46857 +/- 0.00652\n",
      " k-effective (Track-length)  = 1.46267 +/- 0.00782\n",
      " k-effective (Absorption)    = 1.46527 +/- 0.00455\n",
      " Combined k-effective        = 1.46478 +/- 0.00392\n",
      " Leakage Fraction            = 0.00000 +/- 0.00000\n",
      "\n",
      " Creating state point openmc_simulation_n0.h5...\n",
      "[openmc.deplete] t=2592000.0 s, dt=2592000 s, source=174\n",
      " Maximum neutron transport energy: 20000000 eV for U235\n",
      " Initializing source particles...\n",
      "\n",
      " ====================>     K EIGENVALUE SIMULATION     <====================\n",
      "\n",
      "  Bat./Gen.      k            Average k\n",
      "  =========   ========   ====================\n",
      "        1/1    1.53241\n",
      "        2/1    1.40513\n",
      "        3/1    1.45199\n",
      "        4/1    1.48296\n",
      "        5/1    1.45516\n",
      "        6/1    1.49568\n",
      "        7/1    1.41651\n",
      "        8/1    1.41637\n",
      "        9/1    1.41213\n",
      "       10/1    1.43048\n",
      "       11/1    1.43207\n",
      "       12/1    1.51388    1.47297 +/- 0.04091\n",
      "       13/1    1.44591    1.46395 +/- 0.02528\n",
      "       14/1    1.50325    1.47378 +/- 0.02040\n",
      "       15/1    1.44463    1.46795 +/- 0.01684\n",
      "       16/1    1.45528    1.46583 +/- 0.01391\n",
      "       17/1    1.39473    1.45568 +/- 0.01554\n",
      "       18/1    1.41153    1.45016 +/- 0.01454\n",
      "       19/1    1.36173    1.44033 +/- 0.01616\n",
      "       20/1    1.48965    1.44527 +/- 0.01527\n",
      "       21/1    1.42710    1.44361 +/- 0.01391\n",
      "       22/1    1.41558    1.44128 +/- 0.01291\n",
      "       23/1    1.39316    1.43758 +/- 0.01244\n",
      "       24/1    1.46698    1.43968 +/- 0.01171\n",
      "       25/1    1.42137    1.43846 +/- 0.01097\n",
      "       26/1    1.44859    1.43909 +/- 0.01028\n",
      "       27/1    1.39159    1.43630 +/- 0.01005\n",
      "       28/1    1.51030    1.44041 +/- 0.01033\n",
      "       29/1    1.58514    1.44802 +/- 0.01239\n",
      "       30/1    1.33291    1.44227 +/- 0.01309\n",
      "       31/1    1.38694    1.43963 +/- 0.01272\n",
      "       32/1    1.46008    1.44056 +/- 0.01217\n",
      "       33/1    1.39234    1.43847 +/- 0.01181\n",
      "       34/1    1.40957    1.43726 +/- 0.01138\n",
      "       35/1    1.47587    1.43881 +/- 0.01102\n",
      "       36/1    1.43015    1.43847 +/- 0.01059\n",
      "       37/1    1.38802    1.43661 +/- 0.01036\n",
      "       38/1    1.33953    1.43314 +/- 0.01057\n",
      "       39/1    1.48062    1.43478 +/- 0.01033\n",
      "       40/1    1.44053    1.43497 +/- 0.00998\n",
      "       41/1    1.46435    1.43592 +/- 0.00970\n",
      "       42/1    1.49916    1.43789 +/- 0.00960\n",
      "       43/1    1.46589    1.43874 +/- 0.00934\n",
      "       44/1    1.40238    1.43767 +/- 0.00913\n",
      "       45/1    1.44733    1.43795 +/- 0.00887\n",
      "       46/1    1.35951    1.43577 +/- 0.00889\n",
      "       47/1    1.45303    1.43623 +/- 0.00866\n",
      "       48/1    1.43982    1.43633 +/- 0.00843\n",
      "       49/1    1.53025    1.43874 +/- 0.00855\n",
      "       50/1    1.45253    1.43908 +/- 0.00834\n",
      " Creating state point statepoint.50.h5...\n",
      "\n",
      " =======================>     TIMING STATISTICS     <=======================\n",
      "\n",
      " Total time for initialization     = 0.0000e+00 seconds\n",
      "   Reading cross sections          = 0.0000e+00 seconds\n",
      " Total time in simulation          = 1.2879e+01 seconds\n",
      "   Time in transport only          = 1.2870e+01 seconds\n",
      "   Time in inactive batches        = 2.1932e+00 seconds\n",
      "   Time in active batches          = 1.0686e+01 seconds\n",
      "   Time synchronizing fission bank = 4.6210e-03 seconds\n",
      "     Sampling source sites         = 4.2406e-03 seconds\n",
      "     SEND/RECV source sites        = 3.5471e-04 seconds\n",
      "   Time accumulating tallies       = 2.2210e-04 seconds\n",
      "   Time writing statepoints        = 2.5594e-03 seconds\n",
      " Total time for finalization       = 5.4901e-05 seconds\n",
      " Total time elapsed                = 1.2893e+01 seconds\n",
      " Calculation Rate (inactive)       = 4559.54 particles/second\n",
      " Calculation Rate (active)         = 3743.15 particles/second\n",
      "\n",
      " ============================>     RESULTS     <============================\n",
      "\n",
      " k-effective (Collision)     = 1.43941 +/- 0.00730\n",
      " k-effective (Track-length)  = 1.43908 +/- 0.00834\n",
      " k-effective (Absorption)    = 1.43889 +/- 0.00475\n",
      " Combined k-effective        = 1.43891 +/- 0.00466\n",
      " Leakage Fraction            = 0.00000 +/- 0.00000\n",
      "\n",
      " Creating state point openmc_simulation_n1.h5...\n",
      "[openmc.deplete] t=5184000.0 s, dt=2592000 s, source=174\n",
      " Maximum neutron transport energy: 20000000 eV for U235\n",
      " Initializing source particles...\n",
      "\n",
      " ====================>     K EIGENVALUE SIMULATION     <====================\n",
      "\n",
      "  Bat./Gen.      k            Average k\n",
      "  =========   ========   ====================\n",
      "        1/1    1.35648\n",
      "        2/1    1.41942\n",
      "        3/1    1.36608\n",
      "        4/1    1.46674\n",
      "        5/1    1.47455\n",
      "        6/1    1.40320\n",
      "        7/1    1.36282\n",
      "        8/1    1.41123\n",
      "        9/1    1.45872\n",
      "       10/1    1.43235\n",
      "       11/1    1.43551\n",
      "       12/1    1.42639    1.43095 +/- 0.00456\n",
      "       13/1    1.38410    1.41533 +/- 0.01584\n",
      "       14/1    1.40010    1.41152 +/- 0.01183\n",
      "       15/1    1.45122    1.41946 +/- 0.01212\n",
      "       16/1    1.46395    1.42688 +/- 0.01237\n",
      "       17/1    1.50374    1.43786 +/- 0.01516\n",
      "       18/1    1.41657    1.43520 +/- 0.01340\n",
      "       19/1    1.39988    1.43127 +/- 0.01245\n",
      "       20/1    1.48773    1.43692 +/- 0.01248\n",
      "       21/1    1.50428    1.44304 +/- 0.01285\n",
      "       22/1    1.40603    1.43996 +/- 0.01213\n",
      "       23/1    1.44173    1.44009 +/- 0.01115\n",
      "       24/1    1.47605    1.44266 +/- 0.01064\n",
      "       25/1    1.41864    1.44106 +/- 0.01004\n",
      "       26/1    1.42455    1.44003 +/- 0.00944\n",
      "       27/1    1.43444    1.43970 +/- 0.00888\n",
      "       28/1    1.39999    1.43749 +/- 0.00866\n",
      "       29/1    1.45318    1.43832 +/- 0.00823\n",
      "       30/1    1.46514    1.43966 +/- 0.00792\n",
      "       31/1    1.49160    1.44213 +/- 0.00793\n",
      "       32/1    1.34013    1.43750 +/- 0.00887\n",
      "       33/1    1.48829    1.43971 +/- 0.00876\n",
      "       34/1    1.47141    1.44103 +/- 0.00849\n",
      "       35/1    1.39193    1.43906 +/- 0.00838\n",
      "       36/1    1.39036    1.43719 +/- 0.00826\n",
      "       37/1    1.44249    1.43739 +/- 0.00795\n",
      "       38/1    1.48357    1.43904 +/- 0.00784\n",
      "       39/1    1.40641    1.43791 +/- 0.00765\n",
      "       40/1    1.42094    1.43734 +/- 0.00741\n",
      "       41/1    1.40773    1.43639 +/- 0.00723\n",
      "       42/1    1.45978    1.43712 +/- 0.00704\n",
      "       43/1    1.42308    1.43669 +/- 0.00683\n",
      "       44/1    1.47645    1.43786 +/- 0.00673\n",
      "       45/1    1.51314    1.44001 +/- 0.00688\n",
      "       46/1    1.48159    1.44117 +/- 0.00679\n",
      "       47/1    1.48586    1.44238 +/- 0.00671\n",
      "       48/1    1.38670    1.44091 +/- 0.00669\n",
      "       49/1    1.53726    1.44338 +/- 0.00697\n",
      "       50/1    1.44556    1.44344 +/- 0.00680\n",
      " Creating state point statepoint.50.h5...\n",
      "\n",
      " =======================>     TIMING STATISTICS     <=======================\n",
      "\n",
      " Total time for initialization     = 0.0000e+00 seconds\n",
      "   Reading cross sections          = 0.0000e+00 seconds\n",
      " Total time in simulation          = 1.2995e+01 seconds\n",
      "   Time in transport only          = 1.2986e+01 seconds\n",
      "   Time in inactive batches        = 2.2152e+00 seconds\n",
      "   Time in active batches          = 1.0780e+01 seconds\n",
      "   Time synchronizing fission bank = 4.3730e-03 seconds\n",
      "     Sampling source sites         = 3.9823e-03 seconds\n",
      "     SEND/RECV source sites        = 3.6477e-04 seconds\n",
      "   Time accumulating tallies       = 2.3376e-04 seconds\n",
      "   Time writing statepoints        = 2.6625e-03 seconds\n",
      " Total time for finalization       = 6.8176e-03 seconds\n",
      " Total time elapsed                = 1.3016e+01 seconds\n",
      " Calculation Rate (inactive)       = 4514.3 particles/second\n",
      " Calculation Rate (active)         = 3710.67 particles/second\n",
      "\n",
      " ============================>     RESULTS     <============================\n",
      "\n",
      " k-effective (Collision)     = 1.43781 +/- 0.00627\n",
      " k-effective (Track-length)  = 1.44344 +/- 0.00680\n",
      " k-effective (Absorption)    = 1.43243 +/- 0.00527\n",
      " Combined k-effective        = 1.43604 +/- 0.00493\n",
      " Leakage Fraction            = 0.00000 +/- 0.00000\n",
      "\n",
      " Creating state point openmc_simulation_n2.h5...\n",
      "[openmc.deplete] t=7776000.0 s, dt=2592000 s, source=174\n",
      " Maximum neutron transport energy: 20000000 eV for U235\n",
      " Initializing source particles...\n",
      "\n",
      " ====================>     K EIGENVALUE SIMULATION     <====================\n",
      "\n",
      "  Bat./Gen.      k            Average k\n",
      "  =========   ========   ====================\n",
      "        1/1    1.52145\n",
      "        2/1    1.45763\n",
      "        3/1    1.38715\n",
      "        4/1    1.49491\n",
      "        5/1    1.45451\n",
      "        6/1    1.41857\n",
      "        7/1    1.36661\n",
      "        8/1    1.35323\n",
      "        9/1    1.32219\n",
      "       10/1    1.45559\n",
      "       11/1    1.56651\n",
      "       12/1    1.37746    1.47198 +/- 0.09452\n",
      "       13/1    1.41360    1.45252 +/- 0.05794\n",
      "       14/1    1.46543    1.45575 +/- 0.04110\n",
      "       15/1    1.36700    1.43800 +/- 0.03645\n",
      "       16/1    1.41995    1.43499 +/- 0.02991\n",
      "       17/1    1.48070    1.44152 +/- 0.02611\n",
      "       18/1    1.46150    1.44402 +/- 0.02275\n",
      "       19/1    1.34499    1.43302 +/- 0.02288\n",
      "       20/1    1.47485    1.43720 +/- 0.02089\n",
      "       21/1    1.50424    1.44329 +/- 0.01985\n",
      "       22/1    1.41394    1.44085 +/- 0.01829\n",
      "       23/1    1.43248    1.44020 +/- 0.01683\n",
      "       24/1    1.38295    1.43611 +/- 0.01611\n",
      "       25/1    1.42331    1.43526 +/- 0.01503\n",
      "       26/1    1.44462    1.43585 +/- 0.01407\n",
      "       27/1    1.43415    1.43575 +/- 0.01321\n",
      "       28/1    1.48582    1.43853 +/- 0.01276\n",
      "       29/1    1.44884    1.43907 +/- 0.01209\n",
      "       30/1    1.34656    1.43444 +/- 0.01236\n",
      "       31/1    1.44433    1.43492 +/- 0.01177\n",
      "       32/1    1.48062    1.43699 +/- 0.01141\n",
      "       33/1    1.37226    1.43418 +/- 0.01126\n",
      "       34/1    1.31991    1.42942 +/- 0.01179\n",
      "       35/1    1.48850    1.43178 +/- 0.01155\n",
      "       36/1    1.47446    1.43342 +/- 0.01122\n",
      "       37/1    1.49478    1.43569 +/- 0.01103\n",
      "       38/1    1.35729    1.43289 +/- 0.01099\n",
      "       39/1    1.37560    1.43092 +/- 0.01079\n",
      "       40/1    1.41980    1.43055 +/- 0.01043\n",
      "       41/1    1.36832    1.42854 +/- 0.01029\n",
      "       42/1    1.41065    1.42798 +/- 0.00997\n",
      "       43/1    1.36825    1.42617 +/- 0.00984\n",
      "       44/1    1.37304    1.42461 +/- 0.00967\n",
      "       45/1    1.40813    1.42414 +/- 0.00940\n",
      "       46/1    1.39351    1.42329 +/- 0.00918\n",
      "       47/1    1.46699    1.42447 +/- 0.00900\n",
      "       48/1    1.41202    1.42414 +/- 0.00877\n",
      "       49/1    1.45476    1.42493 +/- 0.00858\n",
      "       50/1    1.38327    1.42388 +/- 0.00842\n",
      " Creating state point statepoint.50.h5...\n",
      "\n",
      " =======================>     TIMING STATISTICS     <=======================\n",
      "\n",
      " Total time for initialization     = 0.0000e+00 seconds\n",
      "   Reading cross sections          = 0.0000e+00 seconds\n",
      " Total time in simulation          = 1.3031e+01 seconds\n",
      "   Time in transport only          = 1.3022e+01 seconds\n",
      "   Time in inactive batches        = 2.2077e+00 seconds\n",
      "   Time in active batches          = 1.0824e+01 seconds\n",
      "   Time synchronizing fission bank = 4.5578e-03 seconds\n",
      "     Sampling source sites         = 4.1598e-03 seconds\n",
      "     SEND/RECV source sites        = 3.7168e-04 seconds\n",
      "   Time accumulating tallies       = 2.3355e-04 seconds\n",
      "   Time writing statepoints        = 2.7490e-03 seconds\n",
      " Total time for finalization       = 4.9850e-05 seconds\n",
      " Total time elapsed                = 1.3045e+01 seconds\n",
      " Calculation Rate (inactive)       = 4529.68 particles/second\n",
      " Calculation Rate (active)         = 3695.64 particles/second\n",
      "\n",
      " ============================>     RESULTS     <============================\n",
      "\n",
      " k-effective (Collision)     = 1.42700 +/- 0.00561\n",
      " k-effective (Track-length)  = 1.42388 +/- 0.00842\n",
      " k-effective (Absorption)    = 1.43128 +/- 0.00502\n",
      " Combined k-effective        = 1.42959 +/- 0.00429\n",
      " Leakage Fraction            = 0.00000 +/- 0.00000\n",
      "\n",
      " Creating state point openmc_simulation_n3.h5...\n",
      "[openmc.deplete] t=10368000.0 s, dt=2592000 s, source=174\n",
      " Maximum neutron transport energy: 20000000 eV for U235\n",
      " Initializing source particles...\n",
      "\n",
      " ====================>     K EIGENVALUE SIMULATION     <====================\n",
      "\n",
      "  Bat./Gen.      k            Average k\n",
      "  =========   ========   ====================\n",
      "        1/1    1.38925\n",
      "        2/1    1.42976\n",
      "        3/1    1.41226\n",
      "        4/1    1.39211\n",
      "        5/1    1.50564\n",
      "        6/1    1.47182\n",
      "        7/1    1.40604\n",
      "        8/1    1.47053\n",
      "        9/1    1.38703\n",
      "       10/1    1.40210\n",
      "       11/1    1.49465\n",
      "       12/1    1.38331    1.43898 +/- 0.05567\n",
      "       13/1    1.36733    1.41510 +/- 0.04004\n",
      "       14/1    1.43836    1.42091 +/- 0.02891\n",
      "       15/1    1.28111    1.39295 +/- 0.03582\n",
      "       16/1    1.46121    1.40433 +/- 0.03138\n",
      "       17/1    1.33978    1.39511 +/- 0.02808\n",
      "       18/1    1.41161    1.39717 +/- 0.02440\n",
      "       19/1    1.39461    1.39689 +/- 0.02152\n",
      "       20/1    1.47066    1.40426 +/- 0.02062\n",
      "       21/1    1.42437    1.40609 +/- 0.01874\n",
      "       22/1    1.46622    1.41110 +/- 0.01782\n",
      "       23/1    1.37032    1.40797 +/- 0.01669\n",
      "       24/1    1.41432    1.40842 +/- 0.01546\n",
      "       25/1    1.38921    1.40714 +/- 0.01445\n",
      "       26/1    1.40164    1.40679 +/- 0.01352\n",
      "       27/1    1.42757    1.40802 +/- 0.01276\n",
      "       28/1    1.32544    1.40343 +/- 0.01288\n",
      "       29/1    1.45468    1.40613 +/- 0.01247\n",
      "       30/1    1.49292    1.41047 +/- 0.01261\n",
      "       31/1    1.38744    1.40937 +/- 0.01204\n",
      "       32/1    1.46078    1.41171 +/- 0.01171\n",
      "       33/1    1.42757    1.41240 +/- 0.01122\n",
      "       34/1    1.44109    1.41359 +/- 0.01080\n",
      "       35/1    1.39698    1.41293 +/- 0.01038\n",
      "       36/1    1.47438    1.41529 +/- 0.01025\n",
      "       37/1    1.45001    1.41658 +/- 0.00995\n",
      "       38/1    1.42864    1.41701 +/- 0.00960\n",
      "       39/1    1.44516    1.41798 +/- 0.00931\n",
      "       40/1    1.46250    1.41946 +/- 0.00912\n",
      "       41/1    1.36719    1.41778 +/- 0.00898\n",
      "       42/1    1.32118    1.41476 +/- 0.00920\n",
      "       43/1    1.48217    1.41680 +/- 0.00915\n",
      "       44/1    1.46043    1.41808 +/- 0.00897\n",
      "       45/1    1.47425    1.41969 +/- 0.00886\n",
      "       46/1    1.46925    1.42107 +/- 0.00872\n",
      "       47/1    1.46076    1.42214 +/- 0.00854\n",
      "       48/1    1.34478    1.42010 +/- 0.00856\n",
      "       49/1    1.38493    1.41920 +/- 0.00839\n",
      "       50/1    1.34114    1.41725 +/- 0.00841\n",
      " Creating state point statepoint.50.h5...\n",
      "\n",
      " =======================>     TIMING STATISTICS     <=======================\n",
      "\n",
      " Total time for initialization     = 0.0000e+00 seconds\n",
      "   Reading cross sections          = 0.0000e+00 seconds\n",
      " Total time in simulation          = 1.3126e+01 seconds\n",
      "   Time in transport only          = 1.3116e+01 seconds\n",
      "   Time in inactive batches        = 2.2380e+00 seconds\n",
      "   Time in active batches          = 1.0888e+01 seconds\n",
      "   Time synchronizing fission bank = 4.1247e-03 seconds\n",
      "     Sampling source sites         = 3.7224e-03 seconds\n",
      "     SEND/RECV source sites        = 3.7812e-04 seconds\n",
      "   Time accumulating tallies       = 2.2789e-04 seconds\n",
      "   Time writing statepoints        = 3.0614e-03 seconds\n",
      " Total time for finalization       = 5.4060e-05 seconds\n",
      " Total time elapsed                = 1.3140e+01 seconds\n",
      " Calculation Rate (inactive)       = 4468.2 particles/second\n",
      " Calculation Rate (active)         = 3673.84 particles/second\n",
      "\n",
      " ============================>     RESULTS     <============================\n",
      "\n",
      " k-effective (Collision)     = 1.41574 +/- 0.00666\n",
      " k-effective (Track-length)  = 1.41725 +/- 0.00841\n",
      " k-effective (Absorption)    = 1.42868 +/- 0.00429\n",
      " Combined k-effective        = 1.42755 +/- 0.00464\n",
      " Leakage Fraction            = 0.00000 +/- 0.00000\n",
      "\n",
      " Creating state point openmc_simulation_n4.h5...\n",
      "[openmc.deplete] t=12960000.0 s, dt=2592000 s, source=174\n",
      " Maximum neutron transport energy: 20000000 eV for U235\n",
      " Initializing source particles...\n",
      "\n",
      " ====================>     K EIGENVALUE SIMULATION     <====================\n",
      "\n",
      "  Bat./Gen.      k            Average k\n",
      "  =========   ========   ====================\n",
      "        1/1    1.41046\n",
      "        2/1    1.38415\n",
      "        3/1    1.41886\n",
      "        4/1    1.49276\n",
      "        5/1    1.51088\n",
      "        6/1    1.47652\n",
      "        7/1    1.35478\n",
      "        8/1    1.35145\n",
      "        9/1    1.49048\n",
      "       10/1    1.44388\n",
      "       11/1    1.46152\n",
      "       12/1    1.46663    1.46408 +/- 0.00255\n",
      "       13/1    1.47189    1.46668 +/- 0.00299\n",
      "       14/1    1.39703    1.44927 +/- 0.01754\n",
      "       15/1    1.41472    1.44236 +/- 0.01524\n",
      "       16/1    1.42886    1.44011 +/- 0.01265\n",
      "       17/1    1.45836    1.44272 +/- 0.01100\n",
      "       18/1    1.44504    1.44301 +/- 0.00953\n",
      "       19/1    1.35868    1.43364 +/- 0.01259\n",
      "       20/1    1.36823    1.42710 +/- 0.01302\n",
      "       21/1    1.44284    1.42853 +/- 0.01186\n",
      "       22/1    1.52754    1.43678 +/- 0.01362\n",
      "       23/1    1.44989    1.43779 +/- 0.01257\n",
      "       24/1    1.48022    1.44082 +/- 0.01202\n",
      "       25/1    1.40460    1.43840 +/- 0.01145\n",
      "       26/1    1.36586    1.43387 +/- 0.01163\n",
      "       27/1    1.35671    1.42933 +/- 0.01183\n",
      "       28/1    1.54670    1.43585 +/- 0.01292\n",
      "       29/1    1.41323    1.43466 +/- 0.01228\n",
      "       30/1    1.41355    1.43361 +/- 0.01170\n",
      "       31/1    1.56992    1.44010 +/- 0.01288\n",
      "       32/1    1.35745    1.43634 +/- 0.01284\n",
      "       33/1    1.42131    1.43569 +/- 0.01229\n",
      "       34/1    1.45210    1.43637 +/- 0.01179\n",
      "       35/1    1.36098    1.43335 +/- 0.01170\n",
      "       36/1    1.45621    1.43423 +/- 0.01128\n",
      "       37/1    1.47355    1.43569 +/- 0.01095\n",
      "       38/1    1.44393    1.43598 +/- 0.01055\n",
      "       39/1    1.52203    1.43895 +/- 0.01061\n",
      "       40/1    1.43143    1.43870 +/- 0.01025\n",
      "       41/1    1.44187    1.43880 +/- 0.00991\n",
      "       42/1    1.27834    1.43379 +/- 0.01083\n",
      "       43/1    1.38906    1.43243 +/- 0.01058\n",
      "       44/1    1.41543    1.43193 +/- 0.01028\n",
      "       45/1    1.39747    1.43095 +/- 0.01003\n",
      "       46/1    1.38755    1.42974 +/- 0.00982\n",
      "       47/1    1.40362    1.42904 +/- 0.00958\n",
      "       48/1    1.32535    1.42631 +/- 0.00971\n",
      "       49/1    1.38622    1.42528 +/- 0.00952\n",
      "       50/1    1.40968    1.42489 +/- 0.00928\n",
      " Creating state point statepoint.50.h5...\n",
      "\n",
      " =======================>     TIMING STATISTICS     <=======================\n",
      "\n",
      " Total time for initialization     = 0.0000e+00 seconds\n",
      "   Reading cross sections          = 0.0000e+00 seconds\n",
      " Total time in simulation          = 1.3098e+01 seconds\n",
      "   Time in transport only          = 1.3089e+01 seconds\n",
      "   Time in inactive batches        = 2.2179e+00 seconds\n",
      "   Time in active batches          = 1.0881e+01 seconds\n",
      "   Time synchronizing fission bank = 4.5930e-03 seconds\n",
      "     Sampling source sites         = 4.1842e-03 seconds\n",
      "     SEND/RECV source sites        = 3.7955e-04 seconds\n",
      "   Time accumulating tallies       = 2.4737e-04 seconds\n",
      "   Time writing statepoints        = 2.7292e-03 seconds\n",
      " Total time for finalization       = 5.2371e-05 seconds\n",
      " Total time elapsed                = 1.3112e+01 seconds\n",
      " Calculation Rate (inactive)       = 4508.68 particles/second\n",
      " Calculation Rate (active)         = 3676.29 particles/second\n",
      "\n",
      " ============================>     RESULTS     <============================\n",
      "\n",
      " k-effective (Collision)     = 1.42490 +/- 0.00784\n",
      " k-effective (Track-length)  = 1.42489 +/- 0.00928\n",
      " k-effective (Absorption)    = 1.42584 +/- 0.00494\n",
      " Combined k-effective        = 1.42575 +/- 0.00483\n",
      " Leakage Fraction            = 0.00000 +/- 0.00000\n",
      "\n",
      " Creating state point openmc_simulation_n5.h5...\n",
      "[openmc.deplete] t=15552000.0 (final operator evaluation)\n",
      " Maximum neutron transport energy: 20000000 eV for U235\n",
      " Initializing source particles...\n",
      "\n",
      " ====================>     K EIGENVALUE SIMULATION     <====================\n",
      "\n",
      "  Bat./Gen.      k            Average k\n",
      "  =========   ========   ====================\n",
      "        1/1    1.46692\n",
      "        2/1    1.42016\n",
      "        3/1    1.35815\n",
      "        4/1    1.38064\n",
      "        5/1    1.44504\n",
      "        6/1    1.36182\n",
      "        7/1    1.31279\n",
      "        8/1    1.47632\n",
      "        9/1    1.48287\n",
      "       10/1    1.40492\n",
      "       11/1    1.34730\n",
      "       12/1    1.46660    1.40695 +/- 0.05965\n",
      "       13/1    1.47009    1.42800 +/- 0.04036\n",
      "       14/1    1.41037    1.42359 +/- 0.02888\n",
      "       15/1    1.35040    1.40895 +/- 0.02673\n",
      "       16/1    1.41945    1.41070 +/- 0.02190\n",
      "       17/1    1.54521    1.42992 +/- 0.02668\n",
      "       18/1    1.45638    1.43322 +/- 0.02334\n",
      "       19/1    1.38313    1.42766 +/- 0.02132\n",
      "       20/1    1.48971    1.43386 +/- 0.02006\n",
      "       21/1    1.43673    1.43412 +/- 0.01814\n",
      "       22/1    1.38255    1.42983 +/- 0.01711\n",
      "       23/1    1.47338    1.43318 +/- 0.01609\n",
      "       24/1    1.43212    1.43310 +/- 0.01490\n",
      "       25/1    1.47334    1.43578 +/- 0.01413\n",
      "       26/1    1.44139    1.43613 +/- 0.01322\n",
      "       27/1    1.47399    1.43836 +/- 0.01262\n",
      "       28/1    1.45271    1.43916 +/- 0.01192\n",
      "       29/1    1.41352    1.43781 +/- 0.01136\n",
      "       30/1    1.40013    1.43593 +/- 0.01094\n",
      "       31/1    1.44858    1.43653 +/- 0.01042\n",
      "       32/1    1.38137    1.43402 +/- 0.01025\n",
      "       33/1    1.37217    1.43133 +/- 0.01015\n",
      "       34/1    1.36865    1.42872 +/- 0.01007\n",
      "       35/1    1.35863    1.42592 +/- 0.01005\n",
      "       36/1    1.41942    1.42567 +/- 0.00966\n",
      "       37/1    1.47196    1.42738 +/- 0.00945\n",
      "       38/1    1.50806    1.43026 +/- 0.00956\n",
      "       39/1    1.38029    1.42854 +/- 0.00938\n",
      "       40/1    1.32528    1.42510 +/- 0.00969\n",
      "       41/1    1.43717    1.42549 +/- 0.00938\n",
      "       42/1    1.40786    1.42494 +/- 0.00910\n",
      "       43/1    1.45468    1.42584 +/- 0.00887\n",
      "       44/1    1.32809    1.42296 +/- 0.00907\n",
      "       45/1    1.37289    1.42153 +/- 0.00892\n",
      "       46/1    1.39638    1.42083 +/- 0.00870\n",
      "       47/1    1.45351    1.42172 +/- 0.00851\n",
      "       48/1    1.37909    1.42059 +/- 0.00836\n",
      "       49/1    1.41817    1.42053 +/- 0.00814\n",
      "       50/1    1.36996    1.41927 +/- 0.00803\n",
      " Creating state point statepoint.50.h5...\n",
      "\n",
      " =======================>     TIMING STATISTICS     <=======================\n",
      "\n",
      " Total time for initialization     = 0.0000e+00 seconds\n",
      "   Reading cross sections          = 0.0000e+00 seconds\n",
      " Total time in simulation          = 1.3234e+01 seconds\n",
      "   Time in transport only          = 1.3225e+01 seconds\n",
      "   Time in inactive batches        = 2.2456e+00 seconds\n",
      "   Time in active batches          = 1.0989e+01 seconds\n",
      "   Time synchronizing fission bank = 4.2442e-03 seconds\n",
      "     Sampling source sites         = 3.8368e-03 seconds\n",
      "     SEND/RECV source sites        = 3.8075e-04 seconds\n",
      "   Time accumulating tallies       = 2.4804e-04 seconds\n",
      "   Time writing statepoints        = 2.7128e-03 seconds\n",
      " Total time for finalization       = 5.0970e-05 seconds\n",
      " Total time elapsed                = 1.3248e+01 seconds\n",
      " Calculation Rate (inactive)       = 4453.13 particles/second\n",
      " Calculation Rate (active)         = 3640.13 particles/second\n",
      "\n",
      " ============================>     RESULTS     <============================\n",
      "\n",
      " k-effective (Collision)     = 1.42429 +/- 0.00588\n",
      " k-effective (Track-length)  = 1.41927 +/- 0.00803\n",
      " k-effective (Absorption)    = 1.42275 +/- 0.00472\n",
      " Combined k-effective        = 1.42368 +/- 0.00460\n",
      " Leakage Fraction            = 0.00000 +/- 0.00000\n",
      "\n",
      " Creating state point openmc_simulation_n6.h5...\n"
     ]
    }
   ],
   "source": [
    "integrator.integrate()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Processing the outputs\n",
    "\n",
    "The depletion simulation produces a few output files. First, the statepoint files from each individual transport simulation are written to `openmc_simulation_n<N>.h5`, where `<N>` indicates the current depletion step. Any tallies that we defined in `tallies.xml` will be included in these files across our simulations. We have 7 such files, one for each our of 6 depletion steps and the initial state."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "c5g7.h5\t\t\t openmc_simulation_n2.h5  openmc_simulation_n6.h5\n",
      "depletion_results.h5\t openmc_simulation_n3.h5  statepoint.50.h5\n",
      "openmc_simulation_n0.h5  openmc_simulation_n4.h5  summary.h5\n",
      "openmc_simulation_n1.h5  openmc_simulation_n5.h5\n"
     ]
    }
   ],
   "source": [
    "!ls *.h5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `depletion_results.h5` file contains information that is aggregated over all time steps through depletion. This includes the multiplication factor, as well as concentrations. We can process this file using the `openmc.deplete.Results` object"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "results = openmc.deplete.Results(\"./depletion_results.h5\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "time, k = results.get_keff()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "time /= (24 * 60 * 60)  # convert back to days from seconds"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1.46477526, 0.00392422],\n",
       "       [1.43890826, 0.00465543],\n",
       "       [1.43604211, 0.00493173],\n",
       "       [1.42958815, 0.00429288],\n",
       "       [1.42754967, 0.00463556],\n",
       "       [1.42575323, 0.00483288],\n",
       "       [1.42368305, 0.00459654]])"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "k"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The first column of `k` is the value of `k-combined` at each point in our simulation, while the second column contains the associated uncertainty. We can plot this using `matplotlib`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "from matplotlib import pyplot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEGCAYAAAB2EqL0AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAAkpElEQVR4nO3deXRV9b338fc38wSZw5gBE1QcEYJKFUWsvWpvq7ZWa9tlvfVWrRbt6m1r+/Q+1z63t+uWjmqtY0txKtqi0tYOapXBAQcQUHBgMmEQSQhjCAlJ+D5/nJOQxIScYM7ZJ8nntVZWTvbe5+TDTuDDbw+/Y+6OiIhIbxKCDiAiIgODCkNERCKiwhARkYioMEREJCIqDBERiUhS0AGipaCgwMvKyoKOISIyoCxbtmy7uxd2t27QFkZZWRlLly4NOoaIyIBiZtU9rYvZISkzm21mNWa2qpftpphZi5ldGv76HDNb0eGj0cwujkloERFpF8tzGHOA8w+3gZklArOAp9uWufsCd5/o7hOBGUBDx/UiIhIbMSsMd18M7Ohls5nAY0BND+svBf7u7g39mU1ERHoXN1dJmdkY4BLgrsNs9nlg7mFe4xozW2pmS2tra/s7oojIkBY3hQHcCtzs7ge7W2lmo4ATgad6egF3v9fdK929srCw25P8IiJyhOLpKqlK4BEzAygALjSzFnefH15/GfCEuzcHlE9EZEiLm8Jw93Ftj81sDvBkh7IAuAL4XoxjiYhIWMwKw8zmAtOBAjPbDNwCJAO4+929PLcMKAYWRTeliIj0JGaF4e5X9GHbq7p8XQWM6edI3br8niUAPHrt1Fh8OxGRASOeTnqLiEgcU2GIiEhEVBgiIhIRFYaIiEREhSEiIhFRYYiISERUGCIiEhEVhoiIRESFISIiEVFhiIhIRFQYIiISERWGiIhERIUhIiIRUWF046B70BFEROKOCqOL3fubWb5xFxvrGoKOIiISV1QYXaQnJ9Lqzp0L1wUdRUQkrqgwukhJSqBoWCrzlm1m806NMkRE2qgwujEqO50EM+5auD7oKCIicUOF0Y3UpAQumzKWPyzdxPu79gcdR0QkLqgwevC16RUA3LNIowwREVBh9GhMTjqXTh7L3Nc2sW1PY9BxREQCp8I4jOunV9B60Lln0Yago4iIBE6FcRjFeRl85pQxPPxKNTV7NcoQkaFNhdGLG86poLn1IL95/r2go4iIBEqF0YuygkwunjiGB5dUU1ffFHQcEZHAqDAicMOMChpbWvnNCxpliMjQpcKIQHlhFp86aTQPvFTFzn0Hgo4jIhIIFUaEvj6jgobmVma/qFGGiAxNKowIHT1iGBeeMIo5L1axu6E56DgiIjGnwuiDr8+oYG9TC797SaMMERl6VBh9MGHUcP7l+BHMfuE99jRqlCEiQ0vMCsPMZptZjZmt6mW7KWbWYmaXdlhWYmZPm9nbZvaWmZVFK+ej107l0Wun9rh+5ozx7Gls4YGXqqIVQUQkLsVyhDEHOP9wG5hZIjALeLrLqgeAn7r7BOBUoCYaASNxwphsPj6hiN+88B71TS1BxRARibmYFYa7LwZ29LLZTOAxOhSCmR0HJLn7M+HXqXf3QN/ZaOaM8exqaObBJdVBxhARiam4OYdhZmOAS4C7uqw6GthlZo+b2XIz+2l4JNLda1xjZkvNbGltbW3Usp5cnMP0Ywq57/kNNBzQKENEhoa4KQzgVuBmdz/YZXkSMA34FjAFOAq4qrsXcPd73b3S3SsLCwujGDU0ytix7wAPv7wxqt9HRCRexFNhVAKPmFkVcClwp5ldDGwGVrj7BndvAeYDk4IK2WZyaS7Txhdwz+IN7D/QGnQcEZGoi5vCcPdx7l7m7mXAPOB6d58PvAbkmFnbkGEG8FYwKTu78dzxbK9vYu6rGmWIyOAXy8tq5wJLgGPMbLOZXW1m15nZdYd7nru3Ejoc9ayZvQkYcF/0E/duSlkeU4/K5+5F62ls1ihDRAa3pFh9I3e/og/bXtXl62eAk/o7U3+48dzxXHHfy/xh6SaunFoWdBwRkaiJm0NSA9XpR+Vxalkedy1cT1OLRhkiMnipMD4iM+PGc8ezdXcj85ZtDjqOiEjUqDD6wRkV+UwqyeHOBes50NL1qmARkcFBhdEP2kYZW3bt54nlGmWIyOCkwugnZx9dyMljs7ljwTqaWzXKEJHBR4XRT9pGGZt27Gf+8i1BxxER6XcqjH4049gijh89nF8vWEeLRhkiMsioMPpR2yijqq6Bv7zxftBxRET6lQqjn503YQTHjhzGr55bR+tBDzqOiEi/UWH0s4SE0ChjQ+0+/vrm1qDjiIj0GxVGFJx//EjGF2Xxq2fXclCjDBEZJFQYUZCQYMw8dzxra+r5x+oPgo4jItIvVBhR8skTR3FUYSa3a5QhIoOECiNKEhOMmTMqeOeDvTzz9rag44iIfGQqjCj61EmjKcvP4PZn1+KuUYaIDGwqjChKSkzghnMqWP3+Hp57pyboOCIiH4kKI8ouPmUMxXnpGmWIyICnwoiy5MQEbphewcrNu1m0pjboOCIiR0yFEQOfmTSWMTnp3KZRhogMYCqMGEhJSuBr08tZvnEXL66rCzqOiMgRUWHEyOcqxzJyeBq3PbtGowwRGZBUGDGSmpTI16aX81rVTl7esCPoOCIifabCiKHLpxRTNCyV259dG3QUEZE+U2HEUFpyIteeXc6SDXW8+p5GGSIysKgwYuwLp5ZQkJXCr57TKENEBhYVRoylpyRyzVlH8fza7Syr3hl0HBGRiKkwAvDF00rJy9QoQ0QGFhVGADJTk/j3aeNY+G4tKzftCjqOiEhEVBgBuXJqGTkZyRpliMiAocIISFZqElefMY5/vl3Dqi27g44jItIrFUaAvnxGGcPSkjTKEJEBIWaFYWazzazGzFb1st0UM2sxs0s7LGs1sxXhjz9HP21sDE9L5itnjOOp1dt4e+ueoOOIiBxWLEcYc4DzD7eBmSUCs4Cnu6za7+4Twx+fjlK+QHzljHFkpSZxx3Prgo4iInJYMSsMd18M9HZ780zgMWDIvD1ddkYyV32sjL+t2sqabXuDjiMi0qO4OYdhZmOAS4C7ulmdZmZLzexlM7v4MK9xTXi7pbW1A+fNiq4+cxzpyYkaZYhIXIubwgBuBW5294PdrCt190rgC8CtZlbe3Qu4+73uXunulYWFhVGM2r9yM1O4cmoZf3njfdbV1AcdR0SkW/FUGJXAI2ZWBVwK3Nk2mnD3LeHPG4CFwCnBRIyef582jrSkRO5coFGGiMSnuCkMdx/n7mXuXgbMA6539/lmlmtmqQBmVgCcAbwVYNSoKMhK5UunlzB/xRaqtu8LOo6IyIfE8rLaucAS4Bgz22xmV5vZdWZ2XS9PnQAsNbOVwALgx+4+6AoD4KtnHUVyYgK/1ihDROJQUqy+kbtf0Ydtr+rw+CXgxGhkijdFw9L4wmklPLCkmhvPHU9xXkbQkURE2sXNISkJue7schITjDsXapQhIvFFhRFnRgxP4/NTipm3bDObdzYEHUdEpJ0KIw5dd3boquG7F60POImIyCERF4aZ/SCKOaSD0TnpfK6ymD+8tpmtu/cHHUdEBOjbCOO/zGyWmd1nZl8zs9yopRK+dnY5B925Z9GGoKOIiAB9KwwHGoGngGLgJTM7OSqphOK8DD47aSy/f3UjNXsag44jItKnwnjH3W9x93nu/n+Ai4BfRimXANefU07rQeeexRpliEjw+lIY281sctsX7r4GGDgTNg1ApfmZXDxxDA+/Uk3t3qag44jIENeXwrgReMjMHjKzm83sYeC9KOWSsBvOKedAy0F+87xGGSISrIgLw91XAhOBueFFC4CI796WI3NUYRafPnk0Dyyppq5eowwRCU6f7sNw9yZ3/6u7z3L337i7ZsmLga/PqKCxpZXfvqABnYgERzfuDQAVRcP45ImjuP+lKnY1HAg6jogMUX0qDDM773BfS/TMnDGefQdama1RhogEpK8jjFm9fC1RcszIYVxwwkh+92IVu/c3Bx1HRIagj3pIyvolhUTk6zMq2NvUwpwXq4KOIiJDUESFYWa/M7PZQImZzQ4/lhg7fnQ25x03gt++sIG9jRpliEhsRTrCmAPcD+wMf74/WoHk8G6cMZ49jS08sKQ66CgiMsREVBjuvsjdFwF7OzyG0PxSEkMnjs1mxrFF3Pf8BuqbWoKOIyJDSF/PYXS9plPXeAZg5owKdjU089DLGmWISOz09ca90w/3tcTGKSW5nHV0Ifct3kDDgZ5HGZffs4TL71kSw2QiMpjpxr0B6qZzK6jbd4Dfv7Ix6CgiMkT0uTDM7OhoBJG+mVyaxxkV+dy9aAONza1BxxGRIeBIRhhfNbPr+z2J9NmNM8azvb6Jua9qlCEi0XckhVEHXGtmt5nZV8xsUn+HksicdlQ+p43L4+5F6zXKEJGoO5LCyAfmAf8XWA+c2a+JpE9uOnc82/Y08celm4KOIiKD3JGOMC4FfgiUAy/0ayLpk6nl+VSW5nLnwvU0tWiUISLR0+fCcPcfA18F/guNMAJnZtx47ni27m7ksWVbgo4jIoPYEV1W6+6vuvvu8F3ft/d3KOmbaeMLmFicw68XrKO59WDQcURkkDqiwgif8P5d+PEn+jeS9JWZcdO549myaz9PvK5RhohEx5HeuHcQaHsnnxn9lEU+gunHFHLimGzuWLCOFo0yRCQKjrQwGoBsM0sGSiJ5Qnha9BozW9XLdlPMrMXMLu2yfLiZbTazO44w86DWdi5j444G/rTi/aDj9DtNcyISvIgLw8z+p8OXtxA64X0n8PsIX2IOcH4v3yOR0Lv4Pd3N6h8CiyP8XkPSxycUMWHUcO5YsI7Wg5pIWET6V19GGGPM7AoAd28hdC/GeHd/MpInu/tiYEcvm80EHgNqOi40s8nACLovEgkLncuo4L3t+3jyjcE3yhCRYPWlMK4FrjGzU81sCvAc8JP+CmJmY4BLgLu6LE8Afg58q7++12D2ieNGcsyIYfzquXW4a5QhIv2n18IwswfM7BvAVOAG4F5C/6hf7O5/68cstwI3u3vXM7bXA39z980RZL3GzJaa2dLa2tp+jDZwJCQYM8+tYF1NPTv26e1KRKT/JEWwzRzgZODfgJOAMuA14Etmtsrd5/VTlkrgETMDKAAuNLMWQkU1LTzhYRaQYmb17v7dri/g7vcSKjQqKyuH7H+vLzhhFBVFa9mycz95mSlBxxGRQaLXwnD35wgdfgLAzJKACYRK5DRC5zI+Mncf1+F7zAGedPf5wPwOy68CKrsrCzkkMcGYOaOCmx5Zwc6G5qDjiMggcSRTg7S4+5vu/pC7fzvS55nZXGAJcEz48tirzew6M7uurxmkd/960mjSkhOormvgf//+Nk+t/oDavU1BxxKRASySQ1L9wt2v6MO2V/WwfA6hQ2TSi8QEY1xBJpt2NPC7F6q4p3UDAKX5GUwqyWVSaS6TS3I5ZuQwEhMs4LQiMhDErDAk9oanJXP86Gzu/8qprH5/N8uqd/J69S5eWLedJ5aHphDJTElkYkkOk0tyOaU0l0nFuWRnJAecXETikQpjCEhLTmRyaR6TS/MAcHc279zP6xt3hkpk405+vXB9+81+FUVZTC7JZXJpLpNKcziqIIsEjUJEhjwVxhBkZhTnZVCcl8FFE8cAsK+phZWbd7F84y6WVe/kqbc+4NHwmzJlpyczqSSHSeESObk4h8xU/eqIDDX6Wy8AZKYm8bHyAj5WXgCERiEbtu9jWfVOlodHIgveDd3bkmBw7Mjh7SOQySV5FOelE74kWkQGKRWGdMvMKC/Morwwi8sqiwHYvb+ZFZt2hc+F7OSJ5Vt48OVqAAqyUplUksPk0tAo5IQx2aQlJwb5RxCRfqbCkIhlpydz9tGFnH10IQCtB5012/YeOhdSvZOn39oGQHKicfzo7PbDWJNKcxiVnR5kfBH5iFQYcsQSE4wJo4YzYdRwvnhaKQDb65vaz4O8vnEnD79SzewXQ2+dMjo7jUmlue0lctzo4SQnHukM+yISayoM6VcFWamcd9wIzjtuBAAHWg7y9tY97aOQ5Rt38eQbWwFIS07gpDE5oXtCSnOZVJJDflZqkPFF5DBUGBJVKUkJnFycw8nFOfzbGaHZX7bu3s/r1bvaS+S3L2zg7kWhS3rL8jM6jUKOHqEbC0XihQpDYm5UdjqfPCmdT540CoDG5lZWbdndfhhr8ZrtPB5+b/Ks1CQmFueweWcDuRmaSFEkSCoMCVxaciKVZXlUlh26sXDTjs43Fm7Z1ciWXY1c++BSvvHxo5kwanjAqUWGHhWGxB0zoyQ/g5L8DC4+JXRj4WfvfJEP9jTx0vo6nlr9PBeeOJKbzj2aY0YOCzityNChS1RkQEhKTGBsbjov3DyDG88dz+I12zn/tsXMnLucdTV7g44nMiSoMGRAyU5P5pvnHc0LN5/D9dPLefbtbZz3y8V845HlbKitDzqeyKCmwpABKScjhW//y7G8cPMMrj2rnKdWb+Pjv1jEf/xhJdV1+4KOJzIoqTBkQMvLTOG7FxzL8zefw9VnjuPJN95nxs8X8Z15K9m0oyHoeCKDigpDBoWCrFS+/8njeP4753Dl1FLmr3ifc362kO89/gabd6o4RPqDrpIaxB69dmrQEWKuaHgat3zqeK49q5y7Fq5j7qubmLdsM5dPKeaGcyo0n5XIR6ARhgxKI7PT+H8XncDCb0/nsspiHn1tE2f/ZCG3/GkV2/Y0Bh1PZEBSYcigNjonnR9dciILvjWdz04ew8OvbOSsnyzgv//yFjV7gymOy+9ZwuX3LAnke4t8FCoMGRLG5mbwv585ief+YzqfPnk09y+p4qyfLOBHf32L7fVNQccTGRBUGDKklORn8NPPncyz3zybC08cxW9feI9psxbw47+/w459B4KOJxLXVBgyJJUVZPKLyybyzDfP5hPHj+CexeuZNus5fvbUu+xqUHGIdEeFIUNaeWEWt33+FJ7+xlmcc2wRdyxYx7RZC/jFM2vYvb856HgicUWFIQKMHzGMO74wiX98Yxpnji/g9mfXcuas57j92bXsbVRxiIAKQ6STY0cO564vTeavN57J6Ufl84tn1nDmrAX8esE66ptago4nEigVhkg3jh+dzX1XVvKXr59JZWkuP33qXabNeo67F62n4YCKA3R58FCkwhA5jBPHZvPbq6Yw/4YzOGlsDj/++ztMm7WA+xZvYP+B1qDjicSUCkMkAhOLc7j/K6fy2Nc+xnGjh/Ojv73NtJ8sYPYL79HYrOKQoUGFIdIHk0tzefDq0/jjdVM5ekQW//3kW5z90wXc/1IVTS0qDhncNPmgDAjxNpHilLI8fv/V01myvo5f/nMNt/x5NXcvWs8N51RwWWUxKUn6v5gMPjH7rTaz2WZWY2aretluipm1mNml4a9Lzex1M1thZqvN7LrYJBbp3dTyfB695nQe/vfTGJ2Tzn/OX8U5P1vII69upLn1YNDxRPpVLP8bNAc4/3AbmFkiMAt4usPircBUd58InAZ818xGRymjSJ+ZGWdUFDDvuqk88JVTKRyWyncff5MZP1/IH5ZuokXFEfd0xVdkYlYY7r4Y2NHLZjOBx4CaDs874O5ts8OlovMuEqfMjLOOLuSJ6z/G766aQm5GCt+Z9wYf/8UiHn99s4pDBry4+cfXzMYAlwB3dbOu2MzeADYBs9z9/R5e4xozW2pmS2tra6MbWKQHZsY5xxbxpxvO4DdXVpKRksQ3/7CST/xyMX9asQV3DzqiyBGJm8IAbgVudvcP/TfM3Te5+0lABfBlMxvR3Qu4+73uXunulYWFhdFNK9ILM+Pjx43grzeeyd1fmkxKUgI3PbKCN7bsprpuH3Nf3chrVTvYqVlyZYCIp6ukKoFHzAygALjQzFrcfX7bBu7+fvik+TRgXiApRfrIzDj/hJF84rgR/GP1B3x73kq27W3ie4+/2b5NfmYK5UVZVBRlUVEY/lyUxajsNMJ/J0QCFzeF4e7j2h6b2RzgSXefb2ZjgTp3329mucCZwC8DiilyxBISjAtPHMX9L1Xh7vzscxNZV7uXdTX17R9PrnyfPY2Hph7JTEkMFUlh1qFCKcqiNC+DpMR4OkAgQ0HMCsPM5gLTgQIz2wzcAiQDuPvdh3nqBODnZuaAAT9z9zcPs71I3DMzSvIzKMnPYMaxh46wuju19U2sq6lnfVuR1Nbz0vo6Hl++pX275ESjLD+zvUAqirIoLwx9pKckBvFHkiEgZoXh7lf0YdurOjx+BjgpGplE4o2ZUTQsjaJhaXysvKDTur2Nzayv3ddpRPLOB3t5avUHHPS258OYnPQPHdqqKMoiJyMlgD+RDCZxc0hKRA5vWFoyE4tzmFic02l5U0srVdsbDhVJbejzkvV1NLUcuoakICuF8i4lUlGUxcjhOk8ikVFhiAxwqUmJHDNyGMeMHNZpeetBZ8vO/R86T/KXLudJslKTKC/M/NBJ9xKdJxmQ2m5AjMZ0OioMkUEqMaFv50leXLedx18/dJ4kJTGBsoKM9hIp73CuJC1Z50mGIhWGyBBzuPMkexqbO5XI+pp63np/D/9Y1fk8ydjcdPbubyEtOZGHXq6mLD+T0vwMRuekk5igw1uDlQpDRNoNT0vmlJJcTinJ7bS8sbmVqrrOJ9yfe6eG3Y3N/Of8Q/OJpiQmUJyXHi6QTMoKMkKf8zMYk5OuQ1wDnApDRHqVlpzIsSOHc+zI4e3LLr9nCe7ObVecQtX2Bqrr9vFe3T6qtzdQVbePl9bXsb/Dm0slJRjFeRmU5me0j0jaPo/NzdCU8AOACkNEjpiZMSo7nVHZ6Uwtz++0zt2p3dtEVV2oQKq276M6/Hhp1U7qmw6deE8wGJMbGpl0LJOyglCZ6JxJfFBhiEhUmBlFw9MoGp7GqePyOq1zd+r2HaC6bh9V4RFJVV1olDJ/xRb2driKywxGZ6eHSqQgdHirNFwsJXkZulExhlQYIhJzZkZBVioFWalMLv1wmexqaKaqLjQieW/7vlCx1DXw9ze3srOhudP2I4enUZqfwbiCzPbzJaXhUUpmqv6J60/amyISV8yM3MwUcjNTPnTyHWB3QzPVO8Ijku3h8yZ1Dfzz7W1sr+8882/hsFTK2g9vdT5vMiwtOVZ/pEFDhSEiA0p2RjInZeRw0ticD63b29hMdV1D+7mStkNei9bU8sdlmzttm5+Z0l4gm3fuJz05gbe37mFcQabOmfRAhSEig8awtGROGJPNCWOyP7Su4UBLuEwOnS95b/s+lmyoY+vuRgAuuO15EgxK8jKoKBpGRVEW44uyGD8idMPiUD/ENbT/9CIyZGSkJDFh1HAmjBr+oXWfu+slGltaueasctbW1LOuJjSdyqI1NTS3HnqHxPaJHcNFEvo8jOyMoXF4S4UhEmPRmONHPpqEBCMjJYlPnTy60/Lm1oNU17VN7Lg3XCb1vPJeHY3NhyZ2LByWSkVhaCRS0aFICrJSBtXEjioMEZEeJCcmtBcAjGxffvCgs2XXftbW7GXttlCJrK2p5/HXt3S6vyQ7PbnTIa3xI4YxfgC/k6IKQ0SkjxLCd60X5314Ysdte5pYGz6ktbamnnXb6vnHqg86XQ6cmZIYmsgxPBJpO7xVnJcR13NxqTBERPqJmTEyO42R2WlMG1/YaV1deIbgtR3m4/rQDMFJCe3vWTK+w3mS0vzMuJg6RYUhIhID+Vmp5GelctpRnadQ2b2/mfW1oZHIutp61m7by/KNO/nLyvfbt0lKMErzM0KjkRGd35Y3lpcAqzBERAKUnZ7MpJJcJnW5SbHhQAsbavcdOry1rZ41NXt55u1ttIbnmjeD4tyM9pFIRVEW9U0tpEepRFQYIiJxKCMlqdt7Sjq+Je/a8JVb62vqeX7tdg60hq7cKsxKjUomFYaIyADS+S15R7Uvb2k9yKad+7nmgaWkROl9R1QYIiKDQFJiAuMKMsnLTIna9wj+tLuIiAwIGmGIyBHRHetDj0YYIiISERWGiIhERIUhIiIRUWGIiEhEVBgiIhIRFYaIiEREhSEiIhGJWWGY2WwzqzGzVb1sN8XMWszs0vDXE81siZmtNrM3zOzy2CQWkaHi0Wun6r6SCMRyhDEHOP9wG5hZIjALeLrD4gbgSnc/Pvz8W80sJ0oZRUSkBzErDHdfDOzoZbOZwGNATYfnrXH3teHH74fXFXb/dBERiZa4OYdhZmOAS4C7DrPNqUAKsL6H9deY2VIzW1pbWxudoCIiQ1TcFAZwK3Czux/sbqWZjQIeBP6tp23c/V53r3T3ysJCDUJERPpTPE0+WAk8YmYABcCFZtbi7vPNbDjwV+D77v5ykCFFRIaquCkMdx/X9tjM5gBPhssiBXgCeMDd5wWVT0RkqItZYZjZXGA6UGBmm4FbgGQAd7/7ME+9DDgLyDezq8LLrnL3FVELKyIyQEXz8uCYFYa7X9GHba/q8Pgh4KFoZBIRkcjF00lvERGJYyoMERGJiApDREQiosIQEZGIqDBERCQiKgwREYmICkNERCKiwhARkYiYuwedISrMrBaoPsKnFwDb+zFONA2UrMrZ/wZKVuXsX9HOWeru3c7eOmgL46Mws6XuXhl0jkgMlKzK2f8GSlbl7F9B5tQhKRERiYgKQ0REIqLC6N69QQfog4GSVTn730DJqpz9K7CcOochIiIR0QhDREQiosIQEZGIqDC6MLPzzexdM1tnZt8NOk8bMys2swVm9paZrTazm8LLf2BmW8xsRfjjwqCzAphZlZm9Gc60NLwsz8yeMbO14c+5AWc8psN+W2Fme8zsG/GwT81stpnVmNmqDsu63X8Wcnv4d/YNM5sUcM6fmtk74SxPmFlOeHmZme3vsF8P906bscra48/azL4X3qfvmtm/BJzz0Q4Zq8xsRXh5bPepu+sj/AEkAuuBo4AUYCVwXNC5wtlGAZPCj4cBa4DjgB8A3wo6Xzd5q4CCLst+Anw3/Pi7wKygc3b52X8AlMbDPiX0tsSTgFW97T/gQuDvgAGnA68EnPMTQFL48awOOcs6bhcn+7Tbn3X479ZKIBUYF/53ITGonF3W/xz4ryD2qUYYnZ0KrHP3De5+AHgEuCjgTAC4+1Z3fz38eC/wNjAm2FR9dhFwf/jx/cDFwUX5kHOB9e5+pLMD9Ct3Xwzs6LK4p/13EfCAh7wM5JjZqKByuvvT7t4S/vJlYGwssvSmh33ak4uAR9y9yd3fA9YR+vch6g6X08wMuAyYG4ssXakwOhsDbOrw9Wbi8B9lMysDTgFeCS/6enj4PzvowzwdOPC0mS0zs2vCy0a4+9bw4w+AEcFE69bn6fyXMB73aU/7L55/b79CaPTTZpyZLTezRWY2LahQXXT3s47XfToN2Obuazssi9k+VWEMMGaWBTwGfMPd9wB3AeXARGAroeFqPDjT3ScBFwA3mNlZHVd6aDwdF9d0m1kK8Gngj+FF8bpP28XT/uuJmX0faAEeDi/aCpS4+ynAN4Hfm9nwoPKFxf3Puosr6Pwfm5juUxVGZ1uA4g5fjw0viwtmlkyoLB5298cB3H2bu7e6+0HgPmI0bO6Nu28Jf64BniCUa1vboZLw55rgEnZyAfC6u2+D+N2n9Lz/4u731syuAv4V+GK43Agf3qkLP15G6LzA0YGF5LA/63jcp0nAZ4BH25bFep+qMDp7DRhvZuPC/+v8PPDngDMB7ccufwu87e6/6LC847HqS4BVXZ8ba2aWaWbD2h4TOgm6itC+/HJ4sy8Dfwom4Yd0+l9bPO7TsJ7235+BK8NXS50O7O5w6CrmzOx84DvAp929ocPyQjNLDD8+ChgPbAgmZXumnn7WfwY+b2apZjaOUNZXY52vi48D77j75rYFMd+nsTq7PlA+CF1xsoZQU38/6Dwdcp1J6BDEG8CK8MeFwIPAm+HlfwZGxUHWowhdYbISWN22H4F84FlgLfBPIC8OsmYCdUB2h2WB71NCBbYVaCZ0/PzqnvYfoaujfh3+nX0TqAw45zpCx//bfk/vDm/72fDvwwrgdeBTcbBPe/xZA98P79N3gQuCzBlePge4rsu2Md2nmhpEREQiokNSIiISERWGiIhERIUhIiIRUWGIiEhEVBgiIhIRFYaIiEREhSHSCzPL7zB99AcdpsOuN7M7o/D95pjZe2Z2XQ/r68Ofy9ty9HcGke7oPgyRPjCzHwD17v6zKH6POcCT7j6vh/X17p7V09ci0aIRhsgRMrPpZvZk+PEPzOx+M3vezKrN7DNm9hMLvYnUP8LzgGFmk8Ozii4zs6cimYY8PFXNkvBr/U+0/1wiPVFhiPSfcmAGoZlvHwIWuPuJwH7gk+HS+BVwqbtPBmYDP4rgdW8D7gq/VmBzRIkkBR1AZBD5u7s3m9mbhN7B7x/h5W8Seme0Y4ATgGdCc0mSSGQFcAahOYMgNPfRrH7MLBIxFYZI/2kCcPeDZtbsh04QHiT0d82A1e4+9QheWycbJXA6JCUSO+8ChWY2FULvb2Jmx0fwvBcJTbUP8MVohRPpjQpDJEY89D7xlwKzzGwloSmpPxbBU28i9K6FbxIfbxMqQ5QuqxWJM71dVtvN9rqsVmJCIwyR+LMb+GFPN+61abtxD9gWk1Qy5GmEISIiEdEIQ0REIqLCEBGRiKgwREQkIioMERGJyP8H7Pi8HISODCsAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "pyplot.errorbar(time, k[:, 0], yerr=k[:, 1])\n",
    "pyplot.xlabel(\"Time [d]\")\n",
    "pyplot.ylabel(\"$k_{eff}\\pm \\sigma$\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Due to the low number of particles selected, the uncertainty on each value is rather high. However, we can still see the decline in `k` over time due to fuel consumption."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can then examine concentrations of atoms in each of our materials. This requires knowing the material ID, which can be obtained from the `materials.xml` file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "_, u235 = results.get_atoms(\"1\", \"U235\")\n",
    "_, xe135 = results.get_atoms(\"1\", \"Xe135\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAERCAYAAAB2CKBkAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAAqRElEQVR4nO3deXhU9dnG8e+ThH1HgiIIAUSQXUCRLahlV8ECWq17rYhaRXhbtdVaKm0t2oJYRXHftSqgqOxWCSCyyr4joCACCgXZlzzvHzOxEZMwkEzOTOb+XNdcmTnnzOT2JObhd87vPMfcHRERSVxJQQcQEZFgqRCIiCQ4FQIRkQSnQiAikuBUCEREEpwKgYhIgovLQmBmz5vZNjNbGsG2g8xsuZktNrOPzKxWtnXXm9ma8OP66KYWEYlNFo/XEZhZOrAHeNndGx9n2wuB2e6+z8xuBS5w91+YWWVgHtAKcGA+0NLdd0Y5vohITInLEYG7ZwA7si8zs7pmNtHM5pvZdDNrEN72Y3ffF97sM6BG+HlXYIq77wj/8Z8CdCuk/wQRkZiREnSAAvQ00N/d15hZa2AkcNEx29wETAg/rw58lW3dpvAyEZGEUiQKgZmVBdoCb5tZ1uISx2xzDaHDQB0LN52ISGwrEoWA0CGu/7p785xWmlkn4D6go7sfDC/eDFyQbbMawCfRiygiEpvi8hzBsdx9N7DezC4HsJBm4efnAKOAnu6+LdvbJgFdzKySmVUCuoSXiYgklLgsBGb2BjALqG9mm8zsJuBq4CYzWwQsA3qFN38EKEvosNFCMxsH4O47gCHA3PDjwfAyEZGEEpfTR0VEpODE5YhAREQKTtydLK5SpYqnpaUFHUNEJK7Mnz//W3dPzWld3BWCtLQ05s2bF3QMEZG4YmYbc1unQ0MiIglOhUBEJMGpEIiIJDgVAhGRBKdCICKS4FQIREQSnAqBiEiCS5hC8N2egzz4/nJ2HzgcdBQRkZiSMIVg5rrvePHT9XQZlsFHK7YGHUdEJGYkTCHo2ex0xtzWjgqlinHTS/O4843P+W7PweO/UUSkiEuYQgDQ/IyKvH9He+7qVI8JS7fQeXgG4xZ9jTqwikgiS6hCAFA8JYm7Op3FB3d04IxKpbjzjc+5+eV5fLPrQNDRREQCkXCFIEv908ox5rZ23NfjbGas/ZbOw6bxxpwvNToQkYSTsIUAIDnJuDm9DhMHpNOoenl+P2YJv3xmNhu/2xt0NBGRQpPQhSBLWpUyvP7r8/nbz5uwZPMuuj6awbPTv+BopkYHIlL0qRCEJSUZv2xdkymD0mlbtwp/+XAFfZ78lNVbvw86mohIVKkQHKNahVI8d30rRlzZnI3f7eXix6bz2EdrOHQkM+hoIiJRoUKQAzOjV/PqTBnUkW6NqzFsymp6Pj6DxZv+G3Q0EZECp0KQhyplS/Cvq87hmetasWPvIS57YiYPjV/BgcNHg44mIlJgVAgi0LnhqUwZ1JErWp3BqIwv6D5iOrO/+C7oWCIiBUKFIEIVShXj732a8tqvW3MkM5NfPP0Zf3x3Kd+riZ2IxDkVghPU7swqTLornV+1q82rszfSdXgGH6/aFnQsEZGTpkJwEkoXT+GBSxsy+ta2lC6Rwo0vzGXQvxeyc++hoKOJiJwwFYJ8aFGzEh/e2Z47LzqTcYu+pvPwaXy4eIvaVIhIXFEhyKcSKckM6lKfcb9pT7UKpbj99QX0f3U+23ariZ2IxAcVggLS8PTyjL2tLfd2b8Anq7bTadg03pr3lUYHIhLzVAgKUEpyEv071mXCgA40OK08d7+zmOuen8NXO/YFHU1EJFcqBFFQJ7Usb/Y7nyG9GrFg4066PprBizPXk6kmdiISg1QIoiQpybi2TRqTBqZzblplBr+/nMtHzWLttj1BRxMR+REVgiirUak0L954Lv+8vBlrt+2hx4jpPPHxWg4fVRM7EYkNKgSFwMzo07IGUwd1pFPDqjwyaRW9Hp/J0s27go4mIqJCUJhSy5Vg5NUteeqaFmzfc5BeT8zk4Ykr1cRORAKlQhCAbo2rMXVgR3qfU52Rn6yjx2PTmbdhR9CxRCRBqRAEpELpYjxyeTNe/tV5HDycyeWjZjF43DL2HjwSdDQRSTAqBAFLPyuVyQPTub5NGi/N2kCX4RlkrN4edCwRSSAqBDGgTIkUBvdsxNu3tKFEsSSue34Ov317Ebv2qcW1iESfCkEMaZVWmfF3duC2C+oy9vPNdBo+jYlLvwk6logUcVEtBGa2wcyWmNlCM5uXx3bnmtkRM+sbzTzxoGSxZO7u1oD3bm9HatkS9H91Pre/toDt3x8MOpqIFFGFMSK40N2bu3urnFaaWTIwFJhcCFniRuPqFXjvN+34Xdf6TFm+lc7DpzFmwSY1sRORAhcLh4buAEYDus3XMYolJ3H7hWcyfkAH6qaWZdBbi7jxxbls/u/+oKOJSBES7ULgwGQzm29m/Y5daWbVgZ8DT+b1IWbWz8zmmdm87dsTb0bNmVXL8tYtbRh8aUPmrN9Bl2HTeGXWBjWxE5ECEe1C0N7dWwDdgdvNLP2Y9Y8C97h7no133P1pd2/l7q1SU1OjFDW2JScZN7SrzaS70mlRqxJ/fG8ZVz79GV9sVxM7EcmfqBYCd98c/roNGAucd8wmrYA3zWwD0BcYaWaXRTNTvDujcmle/tV5PNK3KSu/2U33EdN5ato6jqiJnYicpKgVAjMrY2blsp4DXYCl2bdx99runubuacA7wG3u/m60MhUVZsblrc5g6qCOXFA/lb9PWMllI2ey/OvdQUcTkTgUzRHBqcAMM1sEzAE+dPeJZtbfzPpH8fsmjKrlS/LUNS0ZeXULvtl1gJ6Pz+Cfk1dx8Iia2IlI5CzepiO2atXK583L9ZKEhLVz7yGGfLicMQs2c2bVsgzt05SWtSoFHUtEYoSZzc9tGn8sTB+VAlCpTHGGXdGcF288l/2HjtL3qU/58/tqYicix3dChcDMKkcriBSMC+pXZdLAdK49vxYvzNxA10czmLHm26BjiUgMy7UQmNn92Z43NLPVwPxw24jWhZJOTkrZEik82Ksxb93ShuLJSVzz3GzufmcRu/ariZ2I/FReI4Le2Z4/Agxw99rAFcDwqKaSAnFe7cqMH9CBWy+oy+gFm+k8bBqTlqmJnYj8WKSHhk539wkA7j4HKBW9SFKQShZL5p5wE7sqZUtwyytqYiciP5ZXIahjZuPM7H2ghpmVzrauWJRzSQFTEzsRyU1KHut6HfM6CcDMTuU4vYEkNmU1seva6DTuGb2YQW8tYtyir/nrz5tQvaIGeSKJStcRJKijmc4rszbw8KRVGHBv9wZc3boWSUkWdDQRiYKTuo7AzE4zsyfN7AkzO8XMBodvMvOWmVWLXlwpDGpiJyJZ8jpH8CKwHPgK+BjYD/QApgNPRT2ZFIpjm9h1GzGdJz9REzuRRJLroSEz+9zdzwk//9Lda2Zbt9DdmxdOxB/ToaHo2bb7AH98bymTlm2lcfXyPNynGQ1PLx90LBEpACfbYiL7updP4H0Sp6qWL8moa1vx5NUt+GbXQTWxE0kQef1Bf8/MygK4e/arjM8EVkc7mASne5NqTB2UTq/m1fnXf9Zy8WMzmL9xZ9CxRCRKNGtI8jRt9Xb+MGYJX+/azw1t0/htl/qUKZHXrGMRiUUn3X3UzBqY2c+yRgbZlncryIASuzqelcqkgelcl62J3fQ1iXffaJGiLK/po3cC7wF3AEvNLPsFZn+LdjCJHWVLpPDnXo15u3+oid21z80JNbHbpyZ2IkVBXiOCm4GW7n4ZcAHwRzMbEF6nq44S0LlpP25i12n4NCYuVRM7kXiX56whd98D4O4bCBWD7mY2DBWChJW9iV1q2RL0f1VN7ETiXV6FYKuZNc96ES4KlwBVgCZRziUx7tgmdp2GTWP0fDWxE4lHeRWC64Afjfvd/Yi7XwekRzWVxIWsJnbjB3TgzKpl+b+3F3HDC3PZ/N/9QUcTkROQVyHYBxwys8rhRyUzMwB3n1k48SQenFm1LG/f0obBlzZk7oYddBk2jVdmbSAzU6MDkXiQVyGYD8wLf50PLAC2mdlUM0srhGwSR5LUxE4kbuVaCNy9trvXCX/NeqQCI1HTOcmFmtiJxJ8T7hnk7mOAqlHIIkWEmXF5qzOY+n8duah+VYZOXMllI2ey/OvdQUcTkRyccCEIX2WspnNyXFXLleSpa1v+qIndPyat4sBhNbETiSW5No0xs0E5LK4E9AQej1oiKXK6N6lGm7qnMOSDFTz+8VomLN3Cw32b0rJW5aCjiQh5/8u+3DGPsoSmk17j7s8UQjYpQiqWLs4/r2jGS786jwOHM+n71CwGj1vG3oNHgo4mkvDUfVQK3Z6DR3hk4kpemrWRGpVK8VDvJnSolxp0LJEi7aS7j4pEw4+a2KWEmtj97m01sRMJigqBBObctMqMv7MDt11QlzGfq4mdSFBUCCRQJYslc7ea2IkE6oQKgZktiFYQSWw/amK3Qk3sRArTiY4I1H5aouaHJnZ3qomdSGE60ULwYVRSiGST1cTuzz0b/dDE7mU1sROJmhMqBO5+f7SCiGSXlGRc3zbthyZ2D7y3jF88PYt1amInUuB0slhiWlYTu39c3ozVW/fQfcR0Rn6yVk3sRAqQCoHEPDOjb8saTBmUzkX1q/LwxFVcNnImy77eFXQ0kSLhRGcNVTKzptEKI5KXnzaxm8kjk1aqiZ1IPh23EJjZJ2ZW3swqE7o5zTPhG9iLBKJ7k2pMHZTOZc2r88TH67j4senM37gj6FgicSuSEUEFd98N9AZedvfWQKdIPtzMNpjZEjNbaGY/aRBkZleb2eLwNp+aWbMTiy+JSk3sRApOJIUgxcyqAVcAH5zE97jQ3Zvn0uxoPdDR3ZsAQ4CnT+LzJYF1PCuVSQPTue78Wrw0awNdhmeQsXp70LFE4kokheBBYBKw1t3nmlkdYE1BfHN3/9Tdd4ZffgbUKIjPlcSS1cTurVvaUKJYEtc9P4ffqomdSMSi2obazNYDOwEHRrl7rv/iN7PfAg3c/dc5rOsH9AOoWbNmy40bN0YpscS7A4eP8thHaxiV8QWVyxRnSK9GdGtcLehYIoHLqw31cQuBmdUG7gDSyHZHM3fvGcE3ru7um82sKjAFuMPdM3LY7kJgJNDe3b/L6zN1PwKJxNLNu7j7ncUs37Kb7o1P48+9GlG1XMmgY4kEJq9CkOutKrN5F3gOeB84oat43H1z+Os2MxsLnAf8qBCEp6M+C3Q/XhEQiVRWE7unM75gxEdr+HTddzxwSUN6t6iOmVpmiWQXyTmCA+7+mLt/7O7Tsh7He5OZlTGzclnPgS7A0mO2qQmMAa5199UnkV8kV9mb2NULN7G7/oW5bNq5L+hoIjElkkNDvwTqAZOBH5rEu3ueLanDJ5XHhl+mAK+7+1/NrH/4/U+Z2bNAHyDroP+R3IYuWXRoSE5GZqbzymcbGTpxJQbc070B17SuRVKSRgeSGPJ7juAh4FpgHf87NOTuflGBpoyQCoHkx6ad+/jD2KVkrN7OuWmV+HufptRNLRt0LJGoy28hWAs0dPdD0Qh3olQIJL/cndELNjPkg+XsP3yUuzrV4+YOdSiWrNZbUnTl9+b1S4GKBZpIJEDZm9h1OjvcxO6JmSzdrCZ2kpgiKQQVgZVmNsnMxmU9opxLJOqqlivJyKtb8tQ1Ldi6+yC9nlATO0lMkUwf/VPUU4gEqFvjarSpU4UhHy7niY/XMWHpNzzcpymt0ioHHU2kUBx3RBCeKroSKBd+rIhk+qhIPKlQuhj/uLwZL//qPA4ezuTyUWpiJ4kjkjbUVwBzgMsJNZ6bbWZ9ox1MJAjpZ6UyeWA617dJUxM7SRiRzBpaBHR2923h16nAVHcPpGW0Zg1JYZm3YQd3j17MF9v30rdlDe6/+Gwqli4edCyRk5LfWUNJWUUg7LsI3ycS11qlVWb8nR24/cK6jP18M52GZTBhyZagY4kUuEj+oE8Mzxi6wcxuAD4EJkQ3lkhsKFksmd91bcC437Tj1PIluPW1Bdz66ny2fX8g6GgiBSaiNtRm1htoH3453d3H5rV9NOnQkATl8NFMnpn+BY9OXUOpYsncf/HZ9G1ZQ03sJC7k69CQmQ119zHuPij8GGtmQws+pkhsK5acxG0XnMmEAR0469Sy/O6dxVz3/By+2qEmdhLfIjk01DmHZd0LOohIvKibWpZ/92vDg70asWDjTro+msFLn24gMzN6N3kSiaZcC4GZ3WpmS4D64RvMZz3WA4sLL6JI7ElKMq5rk8akgem0SqvMn8Yt44pRs1i7bU/Q0UROWK7nCMysAlAJeAi4N9uq7919RyFky5HOEUiscXfGLNjMgx8sZ/+howzoVI9+6WpiJ7HlpM4RuPsud9/g7le5+0ZgP6F7D5cN31BGRAg1sevTsgZTB3WkU8OqPDJpFb0eVxM7iR+RnCy+1MzWAOuBacAGNH1U5CdSy5UIN7FryfY9oSZ2QyeqiZ3EvkjGrn8BzgdWu3tt4GfAZ1FNJRLHujU+jakDO9KnRXWe/GQdPUZMZ+6GwI6mihxXJIXgcPim8klmluTuHwN53k5SJNFVKF2Mh/s249WbWnPoaCaXPzWLB95byh41sZMYFEkh+K+ZlQUygNfMbASwN7qxRIqG9vWqMOmudG5sl8Yrn22k6/AMPlm17fhvFClEkRSCXsA+YCAwkdC9iy+NZiiRoqRMiRT+dGkj3unfllLFk7nhhbkMemshO/fGxN1fRSJrMRFLNH1U4tnBI0d54j9rGfnJOiqWLsafezamR5PT1KZCoi6/3UdFpICUSElmUJf6jPtNe6pVKMXtry/gllfms223mthJcFQIRALQ8PTyjL2tLb/v3oBpq7fzs2HTeGvuV8TbCF2KhrxaTHwU/qoGcyJRkJKcxC0d6zLxrnTOrlaeu0cv5prnZvPld2piJ4UrrxFBNTNrC/Q0s3PMrEX2R2EFFCnqalcpw5s3n89fLmvMoq920fXRDJ6bsZ6jamInhSSvXkN9gZsI3Yfg2LOz7u4XRTlbjnSyWIqyr/+7n/vfXcp/Vm7jnJoVebhPU+qdWi7oWFIE5HWyOJJ7Fv/R3YdEJdlJUCGQos7dGbfoawaPW8aeg0e446J69O9Yl+IpOqUnJy9fhSD8AT2B9PDLT9z9gwLMd0JUCCRRfLfnIIPfX877i76mwWnlGNqnKc3OqBh0LIlT+b1D2UPAAGB5+DHAzP5WsBFF5FinlC3Bv646h2eua8XOfYf4+ciZPDR+BfsPqYmdFKxIDg0tBpq7e2b4dTLwubs3LYR8P6ERgSSi3QcO89D4lbwx50vSTinNQ72b0qbuKUHHkjhSEBeUVcz2vEK+E4nICSlfshgP9W7C6ze3xoGrnvmMP4xdwu4Dh4OOJkVAJIXgIeBzM3vRzF4C5gN/jW4sEclJ27pVmDggnZs71ObNOV/SZVgG/1m5NehYEuciPVlcDTg3/HKOu38T1VR50KEhkZCFX/2Xe95ZzKqt39Or+ek8cElDTilbIuhYEqPyPWsolqgQiPzPoSOZPPnJOh7/eA3lShbjT5c2pGez09XETn5CTedEiqjiKUkM6FSPD+7owBmVSzPgzYX8+qV5bNm1P+hoEkdUCESKgPqnlWPMrW25/+KzmbnuW7oMy+D12V+SqTYVEoE8C4GZJZvZysIKIyInLznJ+HWHOky+qyNNalTgD2OX8MtnP2PDt7qhoOQtz0Lg7keBVWZWs5DyiEg+1TylNK/9ujV/792EZZt30/XRDJ7OWMeRo5lBR5MYlRLBNpWAZWY2h2z3Knb3nlFLJSL5YmZceV5NLqhflfvfXcrfxq/kw8VbGNq3KQ1OKx90PIkxkVxZ3DGn5e4+7bgfbrYB+B44Chw59oy1haY2jAB6ELov8g3uviCvz9SsIZET4+58uGQLf3pvGbv2H+a2C8/k9gvrUiIlOehoUojymjV03BGBu08zs1pAPXefamalgRP5DbrQ3b/NZV13oF740Rp4MvxVRAqImXFJ09NpV7cKD36wnMc+WsOEJaHRQYualYKOJzEgkqZzNwPvAKPCi6oD7xbQ9+8FvOwhnwEVwxeviUgBq1SmOMN/0ZwXbjiXvQeP0OfJTxnywXL2HToSdDQJWCTTR28H2gG7Adx9DVA1ws93YLKZzTezfjmsrw58le31pvCyHzGzfmY2z8zmbd++PcJvLSI5ubBBVSYNTOea1rV4bsZ6uj6awcy1uQ3aJRFEUggOuvuhrBdmlkLoD3wk2rt7C0KHgG43s/TjvSEn7v60u7dy91apqakn8xEikk25ksUYcllj/t3vfFKSkrj62dnc885idu1XE7tEFEkhmGZmfwBKmVln4G3g/Ug+3N03h79uA8YC5x2zyWbgjGyva4SXiUghaF3nFCYM6ED/jnV5Z8EmOg+bxuRlgbUSk4BEUgjuBbYDS4BbgPHA/cd7k5mVMbNyWc+BLsDSYzYbB1xnIecDu9x9ywnkF5F8KlksmXu7N+Dd29pxStkS9HtlPre/voDt3x8MOpoUkkhmDWWG20/PJnRIaJVH1qnuVGBsuPlVCvC6u080s/7hz32KUFHpAawlNH30xpP6rxCRfGtSowLjftOOUdPW8dhHa5m59lseuKQhPz+nuprYFXGRXEdwMfAUsA4woDZwi7tPiH68n9J1BCLRt3bb99wzegnzN+6k41mp/K13E6pXLBV0LMmHfLWhDvcausTd14Zf1wU+dPcGBZ40AioEIoUjM9N55bONDJ24EgPu6d6Aa1rXIilJo4N4lN821N9nFYGwLwhdLSwiRVhSknF92zQm3ZVOi1qVeOC9Zfzi6Vms274n6GhSwHItBGbW28x6A/PMbLyZ3WBm1xOaMTS30BKKSKDOqFyal391Hv+4vBmrt+6h+4jpPPHxWg6riV2RkdfJ4kuzPd8KZPUc2g7oYKFIAjEz+rasQfpZVRg8bhmPTFrFh4u38HDfpjSuXiHoeJJPulWliJywiUu38Mf3lrFj7yH6pddhwM/qUbKYmtjFsnw1nTOz2sAdQFr27dWGWiRxdWtcjTZ1qvDX8ct58pN1TFr6DUP7NuXctMpBR5OTEMmsoUXAc4QuKPvhoGAkbaijQSMCkdgyY8233DtmMZt27ue6NrW4u1sDypaI5FYnUpjyO310trvHTGtoFQKR2LPv0BH+MWk1L3y6nmrlS/LX3k24sH6kvSmlMOR3+ugIM/uTmbUxsxZZjwLOKCJxrHTxFB64tCHv9G9LmRIp3PjCXAb9eyE79x46/pslcJGM35oA1wIX8b9DQx5+LSLyg5a1KvHBne154uN1jPx4LdNWb+fPvRpxcZNqalMRwyI5NLQWaJi9FXWQdGhIJD6s2LKbe0YvZvGmXXRueCp/uawxp5YvGXSshJXfQ0NLgYoFmkhEiryzq5VnzK1t+UOPBmSs3k6nYdN4c86XxNuU9UQQSSGoCKw0s0lmNi7rEeVcIlIEpCQn0S+9LpPuSqfR6eW5d8wSrn52Nhu/2xt0NMkmkkNDHXNarumjInIiMjOdN+d+xUPjV3A4M5PfdqnPje1qk6wmdoUiX9NHY40KgUh827JrP/ePXcpHK7fR7IyKPNynKfVPKxd0rCIvX+cIzOx7M9sdfhwws6NmtrvgY4pIIqhWoRTPXt+Kx646h6927OOSf01n+JTVHDqiJnZBOW4hcPdy7l7e3csTajbXBxgZ9WQiUmSZGT2bnc7UQR25uEk1Rny0hkv+NZ3Pv9wZdLSEFMnJ4h94yLtA1+jEEZFEUrlMcR698hyev6EV3x84Qu8nP2XIB8vZd+hI0NESSiRN53pne5kEtAIORC2RiCScixqcyuSBlRk6cSXPzVjP5OXf8PfeTWl3ZpWgoyWESEYEl2Z7dCV0d7Je0QwlIomnXMli/OWyJvy73/mkJCVx9bOzuXf0YnbtPxx0tCJPs4ZEJOYcOHyUR6eu4ZnpX3BKmeIMuawxXRudFnSsuHZS00fN7IE8PtPdfUhBhDtRKgQiiWPJpl3cPXoxK7bs5uIm1RjcsxGp5UoEHSsunez00b05PABuAu4p0IQiIjloUqMC437Tjt91rc+U5VvpPHwaYxZsUpuKAhbRoSEzKwcMIFQE3gL+6e7bopwtRxoRiCSmtdv2cM/oxczfuJOOZ6Xy1583pkal0kHHihsnfUGZmVU2s78AiwnNMGrh7vcEVQREJHGdWbUsb9/Shj/3bMTcDTvoOjyDl2dtIDNTo4P8yrUQmNkjwFxCs4SauPtgd9fVHiISmKQk4/q2aUwemE7LtMo88N4yfvH0LNZt3xN0tLiW18niTOAgcITQjWh+WEXoZHH56Mf7KR0aEhEAd2f0gs0M+WA5+w8fZcDP6tEvvQ7Fkk/oOtmEcVKHhtw9yd1LZW8xEX6UC6oIiIhkMTP6tqzBlEHpdDq7Ko9MWkWvx2eydPOuoKPFHZVOEYlrVcuVZOTVLXnqmpZs33OQXk/MZOjElRw4fDToaHFDhUBEioRujU9j6sCO9GlRnSc/WUePEdOZs35H0LHiggqBiBQZFUoX4+G+zXj1ptYcOprJFaNm8cd3l/L9AbWpyIsKgYgUOe3rVWHywHRubJfGq7M30nV4Bh+v0qz33KgQiEiRVLp4Cn+6tBGjb21LmRIp3PjCXAb+eyE79h4KOlrMUSEQkSKtRc1KfHBne+78WT3eX/Q1nYdN4/1FX6tNRTYqBCJS5JVISWZQ57N4/472VK9Uijve+JybX57PN7t0axVQIRCRBHJ2tfKMubUt9/U4mxlrt9N52DTemPNlwo8OVAhEJKGkJCdxc3odJg5Ip1H18vx+zBJ++cxsNny79/hvLqJUCEQkIaVVKcMbN5/PQ72bsHTzLrqNyOCZjC84moBN7FQIRCRhmRlXnVeTKYM60v7MVP46fgW9R85k5Te7g45WqFQIRCThnVahJM9c15J/XXUOm3bu55LHZjBsymoOHkmMNhVRLwRmlmxmn5vZBzmsq2lmH4fXLzazHtHOIyKSEzPj0manM2VQRy5tdjqPfbSGSx6bwYIvi373/cIYEQwAVuSy7n7gLXc/B7gSGFkIeUREclW5THGG/6I5L9xwLnsPHqHPk5/y4PvL2XfoSNDRoiaqhcDMagAXA8/msokDWS2tKwBfRzOPiEikLmxQlUkD07mmdS2en7meLsMzmLHm26BjRUW0RwSPAncDmbmsHwxcY2abgPHAHTltZGb9zGyemc3bvn17NHKKiPxEuZLFGHJZY966pQ3Fk5O45rnZ3P3OInbtK1pN7KJWCMzsEmCbu8/PY7OrgBfdvQbQA3jFzH6Syd2fdvdW7t4qNTU1SolFRHJ2Xu3KjB/QgVsvqMvoBZvpNHwaE5d+E3SsAhPNEUE7oKeZbQDeBC4ys1eP2eYm4C0Ad58FlASqRDGTiMhJKVksmXu6NeC929uRWrYE/V+dz22vzWfb9/HfpiJqhcDdf+/uNdw9jdCJ4P+4+zXHbPYl8DMAMzubUCHQsR8RiVmNq1fgvd+043dd6zN1xTY6D8vgnfmb4rpNRaFfR2BmD5pZz/DL/wNuNrNFwBvADR7Pe1NEEkKx5CRuv/BMxt/ZgXpVy/Lbtxdx3fNz+GrHvqCjnRSLt7+7rVq18nnz5gUdQ0QEgMxM57XZG/n7hJU4cHfX+lzXJo2kJAs62o+Y2Xx3b5XTOl1ZLCKSD0lJxrVt0pg8qCPnplVm8PvLuXzULNZu+z7oaBFTIRARKQDVK5bixRvPZdgVzVi3fQ89Rszg8f+s4fDR3GbPxw4VAhGRAmJm9G5RgykDO9K50an8Y/JqLv3XDJZs2hV0tDypEIiIFLDUciV44pctGHVtS3bsPUSvJ2bw0IQVHDgcm03sVAhERKKka6PTmDKoI7849wxGTfuCbo9m8NkX3wUd6ydUCEREoqhCqWI81Lspr/+6NZkOVz79GfeNXcL3B2KnTYUKgYhIIWh7ZhUm3ZXOzR1q88acL+kyPIP/rNwadCxAhUBEpNCUKp7MfRc3ZMxt7Shfshi/enEeA978nO/2HAw0lwqBiEgha35GRd6/oz13darH+CVb6Dw8g/cWbg6sTYUKgYhIAIqnJHFXp7P48M4O1KxcmgFvLuTXL81jy679hZ5FhUBEJEBnnVqO0be25Y+XNOTTdd/ReVgGr83eSGZm4Y0OVAhERAKWnGTc1L42k+5Kp9kZFbhv7FKueuYz1n+7t1C+vwqBiEiMqHlKaV69qTUP92nK8i276fZoBqOmreNIlNtUqBCIiMQQM+OKc89g6qCOdDwrlYcmrOTnIz9l+de7o/Y9VQhERGLQqeVLMuraloy8ugVbdu2n5+MzeG7G+qh8r5SofKqIiOSbmdGjSTXa1DmFv3y4glqVS0fl+6gQiIjEuEplivPPK5pF7fN1aEhEJMGpEIiIJDgVAhGRBKdCICKS4FQIREQSnAqBiEiCUyEQEUlwKgQiIgnOgroRwskys+3AxpN8exXg2wKME03xklU5C168ZFXOghXtnLXcPTWnFXFXCPLDzOa5e6ugc0QiXrIqZ8GLl6zKWbCCzKlDQyIiCU6FQEQkwSVaIXg66AAnIF6yKmfBi5esylmwAsuZUOcIRETkpxJtRCAiIsdQIRARSXAJUwjMrJuZrTKztWZ2b9B5spjZGWb2sZktN7NlZjYgvHywmW02s4XhR48YyLrBzJaE88wLL6tsZlPMbE34a6UYyFk/235baGa7zeyuWNinZva8mW0zs6XZluW4Dy3ksfDv7GIzaxFwzkfMbGU4y1gzqxhenmZm+7Pt16cKK2ceWXP9WZvZ78P7dJWZdQ0457+zZdxgZgvDywt3n7p7kX8AycA6oA5QHFgENAw6VzhbNaBF+Hk5YDXQEBgM/DbofMdk3QBUOWbZw8C94ef3AkODzpnDz/4boFYs7FMgHWgBLD3ePgR6ABMAA84HZgecswuQEn4+NFvOtOzbxcg+zfFnHf5/axFQAqgd/ruQHFTOY9b/E3ggiH2aKCOC84C17v6Fux8C3gR6BZwJAHff4u4Lws+/B1YA1YNNdUJ6AS+Fn78EXBZclBz9DFjn7id7NXqBcvcMYMcxi3Pbh72Alz3kM6CimVULKqe7T3b3I+GXnwE1CiPL8eSyT3PTC3jT3Q+6+3pgLaG/D1GXV04zM+AK4I3CyHKsRCkE1YGvsr3eRAz+sTWzNOAcYHZ40W/Cw/DnY+GQC+DAZDObb2b9wstOdfct4effAKcGEy1XV/Lj/7libZ9C7vswln9vf0VotJKltpl9bmbTzKxDUKGOkdPPOlb3aQdgq7uvybas0PZpohSCmGdmZYHRwF3uvht4EqgLNAe2EBo2Bq29u7cAugO3m1l69pUeGtPGzHxkMysO9ATeDi+KxX36I7G2D3NiZvcBR4DXwou2ADXd/RxgEPC6mZUPKl9YzP+sj3EVP/4HS6Hu00QpBJuBM7K9rhFeFhPMrBihIvCau48BcPet7n7U3TOBZyik4Wte3H1z+Os2YCyhTFuzDleEv24LLuFPdAcWuPtWiM19GpbbPoy531szuwG4BLg6XLQIH2b5Lvx8PqHj7mcFFpI8f9axuE9TgN7Av7OWFfY+TZRCMBeoZ2a1w/9KvBIYF3Am4Idjg88BK9x9WLbl2Y8F/xxYeux7C5OZlTGzclnPCZ04XEpoP14f3ux64L1gEuboR//KirV9mk1u+3AccF149tD5wK5sh5AKnZl1A+4Gerr7vmzLU80sOfy8DlAP+CKYlD9kyu1nPQ640sxKmFltQlnnFHa+Y3QCVrr7pqwFhb5PC+usdNAPQjMwVhOqrPcFnSdbrvaEDgUsBhaGHz2AV4Al4eXjgGoB56xDaLbFImBZ1j4ETgE+AtYAU4HKQe/TcK4ywHdAhWzLAt+nhArTFuAwoePTN+W2DwnNFnoi/Du7BGgVcM61hI6vZ/2ePhXetk/4d2IhsAC4NAb2aa4/a+C+8D5dBXQPMmd4+YtA/2O2LdR9qhYTIiIJLlEODYmISC5UCEREEpwKgYhIglMhEBFJcCoEIiIJToVARCTBqRBIwjKzU7K1+f0mW9viPWY2Mgrf70UzW29m/XNZvyf8tW5WjoLOIJITXUcgQqh/PbDH3f8Rxe/xIvCBu7+Ty/o97l42t9ci0aIRgcgxzOwCM/sg/Hywmb1kZtPNbKOZ9Tazhy10g56J4T5RmFnLcJfI+WY2KZJ20eGWJ7PCn/WXaP93ieRGhUDk+OoCFxHqZPoq8LG7NwH2AxeHi8G/gL7u3hJ4HvhrBJ87Angy/FmB9RASSQk6gEgcmODuh81sCaE7nk0ML19C6E5S9YHGwJRQD0GSiewPeztCPWUg1BtnaAFmFomYCoHI8R0EcPdMMzvs/zuxlkno/yEDlrl7m5P4bJ2kk8Dp0JBI/q0CUs2sDYTuL2FmjSJ430xCLdEBro5WOJHjUSEQyScP3Qe7LzDUzBYRah3cNoK3DiB0p7clxMbtEiVBafqoSCE53vTRHLbX9FEpFBoRiBSeXcCQ3C4oy5J1QRmwtVBSScLTiEBEJMFpRCAikuBUCEREEpwKgYhIglMhEBFJcP8PpwfWdmqPvHMAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "pyplot.plot(time, u235, label=\"U235\")\n",
    "pyplot.xlabel(\"Time [d]\")\n",
    "pyplot.ylabel(\"Number of atoms - U235\");"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAERCAYAAAB2CKBkAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAAmU0lEQVR4nO3de5xddXnv8c8z95lkMpNJhgAJuQJaEBUMkOClVmsVVPAoWmyPVXuhterRHu3R1iq+tK0v26M9WkVKDxywRbBStWjxWq/tTCAXLgHCJXsIkJDb7Lnfb8/5Y62Z7Ayzd9Yks2btvdf3/Xrt11577bXXflgT5pnf7/mt38/cHRERSa+KpAMQEZFkKRGIiKScEoGISMopEYiIpJwSgYhIyikRiIikXEkmAjO72cyOmNlDEY59hZntMrMJM7t61nuTZnZ/+LgrvohFRIpXSSYC4BbgdRGPfRp4F/C1Od4bdvcXh48rFyg2EZGSUpKJwN1/AXTl7jOzTWb2fTPbaWa/NLPnh8fuc/cHgakkYhURKXYlmQjyuBF4v7u/BPgwcH2Ez9SZ2Q4z22Zmb4o1OhGRIlWVdAALwcyWApcB3zCz6d21ET66zt0PmNlG4CdmttvdM3HFKSJSjMoiERC0bHrc/cXz+ZC7HwifO8zsZ8CFgBKBiKRKWXQNuXsf8KSZvRXAAi8q9BkzW25mteH2SuClwCOxBysiUmSsFGcfNbPbgVcCK4HDwHXAT4CvAGcA1cAd7v4pM7sY+BawHBgBDrn7+WZ2GfAPBEXkCuD/uPtNi/3fIiKStJJMBCIisnDKomtIREROXskVi1euXOnr169POgwRkZKyc+fOTndvneu9kksE69evZ8eOHUmHISJSUszsqXzvqWtIRCTllAhERFJOiUBEJOWUCEREUk6JQEQk5ZQIRERSTolARCTlSu4+AoHHD/fzo0cOs6y+muUN1SxvqKE5fF7eUEN9TWXSIYpICVEiKEF/8/1H+fGeI3nfr62qoGVJDc0NNc9JFDMJY8mxxLG8oYbGuioqKizvOUWkfMWWCMzsLOCrwCrAgRvd/Quzjnkl8G/Ak+Gub7r7p+KKqRxMTE5xT0cXv7n5LD70G+fSPTRO99AYPUNjM9vdg8H29L49h/roCV9P5ZljsMKg+biWRfWxRLKkZta+Y9s1VepdFCl1cbYIJoAPufsuM2sEdprZj9x99pz/v3T3N8QYR1l56Nk++kcneNk5KzltWR2nLauL/NmpKad/ZILuoTG6ppPH4HQiOf75QM8IDz/bR9fgGKMT+Zd7XlJTGSSH41oYxyeR2S2TpbVV5KwkJyIJiy0RuPtB4GC43W9me4DVaPGXU9KW6QRgy8YV8/5sRYXR1FBNU0M161kS+XPDY5NBSyMnUXQPjdOT0/LoCvc93TVE9+AYfSMTec9XXWkzySE3SQQtj2OtjpYlx7ab6qupVNeVSCwWpUZgZusJloG8Z463t5rZA8CzwIfd/eHFiKlUtWeyPG9VI62NUZZkXhj1NZXU19RzZnN95M9MTE7ROzx+XBfVc7utgu0nOwfZNdRD9+AYE3n6rsygqf5Yq6KloSZ/19WSY++r60rkxGJPBOHC8v8KfDBcUjLXLoIF5AfM7Arg28A5c5zjWuBagLVr18YbcBEbm5hi+74urrm4+K9BVWUFK5bWsmJp9ITl7gyMTsy0OroGj2+BBEkk2Heob4Q9B/voHhpneHwy7zmnu66C4vmx+sd08mhuqKYlZ3t5Qw0NNZXqupJUiTURmFk1QRK4zd2/Ofv93MTg7neb2fVmttLdO2cddyNwI8DmzZtTu6Ta/c/0MDI+xdZN8+8WKgVmRmNdNY111ZzV0hD5cyPjk/QMjYeJI0gaXUNjJ911VVNVcVxNY/YIrOlRVzNdWBp1JSUuzlFDBtwE7HH3z+c55nTgsLu7mV1CcINbNq6YSl1bphMz2LKhPBPByaqrruT0pkpOb4peOJ+YnKJn+Fi3VW4Smd2F9fjhAboHx+gZHmcyT9fVzKir+qAG01wfJIqm+mqa6qtpbggf9TUsm34dvldVqe4rSVacLYKXAu8AdpvZ/eG+PwfWArj7DcDVwHvMbAIYBq5xLaKcV3smywvObKKpoTrpUEpeVWUFK5fWsnIeXVdTU07/6ETQwsjXbTU8Tt/wOJ0DY+w9OkDP0Dj9BVofAEtrq45LFkHyqDkuWQT7a2beb26opr5aXViyMOIcNfSfQMF/pe7+JeBLccVQTobHJrnv6R7e9dL1SYeSWhUVNvMX/roV0UddTUxO0T8yMdMC6R0ep3d4PLy3Y5ye4XDf0Dg9w+Mc6u2nd3iC3uExxifz/11UU1kRjAKrn26BBMniuUmleqa10twQdL1pBJbk0p3FJWLnU92MTZZvfaCcVVVWBMXpJTUwj2G77s7Q2CQ9M0libCZZTCeS3uGxmYTybM8Iew720zM0xuBY/gK6GTTWVs3cQHhc99V0S6ShhlXLalm1rI5VjXUsq9e9H+VMiaBEtHd0UlVhXLy+JelQZJGYGUtqq1hSW8XqeQzdhWCEWd/I8cliphUyPE5v2I01vW9/9/BMa2WuMkhddUWQFMLH6WGSOG1ZHacvq5tJGnXVmueqFCkRlIi2TJYXrmliaa1+ZHJiNVXzr4FAUAcZGJuga2CMI/2jHOob4UjfCId6RzjcP8rh3hEe3N/DD3tH5rzjvKm++lhLYlaSmH6sXFqjAnmR0W+VEjAwOsGD+3t5z69uSjoUKXMVFcayumqW1VWzfmX+bix3p294gsP9YZLoGwkSR7h9uG+EJw4PcHRg9DkjrSoMVi6t5fSmOk5rDBLF6dOJoilMHI11NDdUqztqkSgRlIDtT3YxOeVcpvqAFAmzY9OVnLuqMe9xk1NOdmCUw31B6+LwcY9R9ncPsfOpLrqHxp/z2ZqqipkkMbsLKrfFoWnXT50SQQloy3RSU1nBReuWJx2KyLxUVtjM5IgX0JT3uJHxSY725yaL0ZmEcah3hEee7eMne47MeRd5Y13VTFI4Lbd1ESaN1sZaViypVcIoQImgBLRlsly0rlmFOClbddWVnNXSUPCOcvfgPo6gZhEkiukaxnSLoyMzwJH+0TnnrKqrrmDFktpgLqoltbQ0hM9LZj+HU5LUV6fmbnElgiLXMzTGIwf7+JNfPzfpUEQSZXasfnH2afm7o6amnOzg2ExrIjs4SnYwuFt8+rlrcIwnOwfoGsg/1Hb6bvGWJcE0Ii3hEOAVBZ5L9Y81JYIit62jC3d0/4BIRBUVRmtjLa2Ntbxgdf7uqGkj45MzkxwWemSODtD9VLCdb4GnhppKljfUsGJpMA/VdJJomesRTkFSDK0OJYIi157ppL66khetaU46FJGyVFddyRlN9ZzRFO1ejakpp29k/DmtjNzWRtfQseTRNTjGUIFWx/KGPIki5zGdXFqW1FBbtfCtDiWCIteWyXLxhhbNqy9SJCoqLLwruwZao31mZHzyuNZF99AY2YHwOSeR7D0yMPP+XK2OP3j5Bj72+vMW9j8IJYKidrR/lCeODPCWl6xJOhQROQV11ZWc2Rx9caepKad3eHymZTH9KDRU91QoERSx9o5gRu6tJ7EspYiUrooKm5mfalPEVscpfV/8XyEnqz3TSWNdFeefuSzpUESkjCkRFLH2TJZLN6zQvCwiEiv9hilSB3qG2Zcd0rBREYmdEkGRas8E9QHNLyQicVMiKFLtmSwtS2p4XkyjBEREpikRFCF3pz3TyZaNLUVx16GIlDclgiL0VHaIZ3tH2LppZdKhiEgKKBEUoTbVB0RkESkRFKH2jiyrltWyscAKUSIiC0WJoMhM1we2blyhZfpEZFEoERSZJ44M0DkwxmWqD4jIIlEiKDLT9w/oRjIRWSxKBEWmLdPJmuX1BZfsExFZSEoERWRqytnW0aXRQiKyqCInAjNbamYXmVlzjPGk2iMH++gdHld9QEQWVd5EYGbX52y/DHgE+Byw28yuWITYUkf1ARFJQqGFabbkbH8aeJO77zKzjcC/AHfHGlkKtWU62di6hFXL6pIORURSJGrX0DJ33wXg7h1RPmdmZ5nZT83sETN72Mw+MMcxZmZfNLO9ZvagmV00v/DLx/jkFPc+qfqAiCy+Qi2C55vZg4AB681subt3m1kFUBPh3BPAh8JWRCOw08x+5O6P5BxzOXBO+LgU+Er4nDq7D/QyODbJ1o2qD4jI4iqUCH5l1uvB8LkF+MSJTuzuB4GD4Xa/me0BVhPUGqZdBXzV3R3YZmbNZnZG+NlUma4PbNnYknAkIpI2eROBuz+VZ38n8M35fImZrQcuBO6Z9dZq4Jmc1/vDfcclAjO7FrgWYO3atfP56pLRlunk+ac3smJpbdKhiEjKFBo19Hwz+56Z/buZbTKzW8ysx8zuNbPZrYW8zGwp8K/AB92972SCdPcb3X2zu29ubW09mVMUtdGJSXbs69awURFJRKGi743A9cA/Az8Bvg8sJxhB9KUoJzezaoIkcJu7z9WKOACclfN6TbgvVe57uofRiSkNGxWRRBRKBI3u/h13vx0Yd/c7PPAdgoRQkAVTZ94E7HH3z+c57C7gd8LRQ1uA3jTWB9oyWSoMLtmg+oCILL5CxeLKnO3Zv8ijjBp6KfAOghvQ7g/3/TmwFsDdbyC4F+EKYC8wBLw7wnnLzrZMlgtWN9FUX510KCKSQoUSwZfNbKm7D7h77l3GZwM/PtGJ3f0/CYaeFjrGgfdGDbYcDY1NcN8z3fzuyzYkHYqIpFShUUP/kGf/XuCDcQWUNjv2dTM+6SoUi0hiTmr2UTM74X0EEk17R5aqCuPi9Scsu4iIxOJkp6H+/QWNIsXaMllefFYzDTWFeulEROKT97ePmeUb829AfTzhpEvfyDi79/fwvl87O+lQRCTFCv0Z2gNc7O6HZ79hZs8893CZr3s7uphy2Kr6gIgkqFDX0FeBdXne+1oMsaROe0eW2qoKLlzbnHQoIpJihUYN/UWB9z4STzjp0pbJ8pJ1y6mrrjzxwSIiMYmyroCZ2X+fHilkZmvN7JL4QytvXYNj7DnYp/UHRCRxUUYNXQ9sBd4evu4HvhxbRClxT8f0spSqD4hIsqKMWbzU3S8ys/sAwsVpokwxIQW0ZbI01FTywjVNSYciIikXpUUwbmaVgAOYWSswFWtUKdCW6eSSDS1UV57srRwiIgsjym+hLwLfAk4zs78C/hP4TKxRlbkjfSNkjg6qPiAiReGEXUPufpuZ7QReTXAz2ZuAJ2KOq6y1T9cHtD6xiBSBQiuUfcfM1gG4+6Pu/mV3/xLBUpL3L1J8Zaltb5ZldVWcd+aypEMRESnYNXQH8FMz+5iZVZvZmWb2L8BfAe9cnPDKU1tHJ1s2rqCyouAs3SIiiyJvInD32wgWnF8L7AHaCdYh2OLuOxcnvPLzTNcQz3QNa1lKESkaJyoWnwdcAtwLjAKriDbkVPKYrg9o/QERKRaFagQ3Edw49sfu/lsErYMm4AEz+41Fiq/stGeyrFhSw7mrliYdiogIULhF8BDB7KPtAO4+6O4fBn4T+PhiBFdu3J32TJatm1ZgpvqAiBSHQpPO/V2e/buBl8cWURl7snOQQ30jqg+ISFHRba2LqC2j+oCIFB8lgkXU3pHljKY61q9oSDoUEZEZSgSLZGrK2ZbJsnWj6gMiUlzmlQjMbFdcgZS7x4/0kx0cU31ARIrOfFsE+lP2JLVnptcfUCIQkeIy30Tw77FEkQJtmSxrWxpYs1z1AREpLvNKBIXWMZb8JqecbR1ZTTstIkVJxeJF8PCzvfSPTKhbSESKkhLBIlB9QESK2XxHDS03sxfGFUy5astkOfu0pZzWWJd0KCIiz3HCRGBmPzOzZWbWAuwC/tHMPh/hczeb2REzeyjP+680s14zuz98fGL+4Re/8ckptu/rUn1ARIpWlBZBk7v3AW8GvurulwK/HuFztwCvO8Exv3T3F4ePT0U4Z8l5cH8PQ2OTSgQiUrSiJIIqMzsDeBvw3agndvdfAF0nG1i5aNubxQwu3aBEICLFKUoi+BTwA2Cvu283s40s3OL1W83sATP7npmdn+8gM7vWzHaY2Y6jR48u0FcvjrZMll85fRnLl9QkHYqIyJxOmAjc/Rvu/kJ3/+PwdYe7v2UBvnsXsM7dXwT8PfDtAjHc6O6b3X1za2vrAnz14hgZn2Tn093qFhKRonbCZSfNbAPwfmB97vHufuWpfHFYd5jevtvMrjezle7eeSrnLSa7nu5mbGJKw0ZFpKhFWX/428BNwHeAqYX6YjM7HTjs7m5mlxC0TrILdf5i0J7JUllhXLKhJelQRETyipIIRtz9i/M9sZndDrwSWGlm+4HrgGoAd78BuBp4j5lNAMPANe7u8/2eYtaWyXLB6iYa66qTDkVEJK8oieALZnYd8ENgdHqnuxecktrd336C978EfClKkKVocHSCB57p4dpXbEw6FBGRgqIkgguAdwCv4ljXkIevJY/t+7qYmHLVB0Sk6EVJBG8FNrr7WNzBlJP2TJbqSmPzOtUHRKS4RbmP4CGgOeY4yk57R5YL1y6nvqYy6VBERAqK0iJoBh41s+0cXyM4peGj5ax3aJyHDvTy/ledk3QoIiInFCURXBd7FGXmniezTDm6kUxESsIJE4G7/9zMVgEXh7vudfcj8YZV2to7stRVV/Ditc1JhyIickJRpqF+G3AvQdH4bcA9ZnZ13IGVsvZMls3rWqitUn1ARIpflK6hjwEXT7cCzKwV+DFwZ5yBlarswCiPHurnT197ZtKhiIhEEmXUUMWsrqBsxM+l0raOYOZt1QdEpFREaRF838x+ANwevv5N4HvxhVTa2jKdLK2t4oLVTUmHIiISSZRi8Z+a2ZuBl4W7bnT3b8UbVulqz2S5ZEMLVZVqNIlIaYgyDfVn3f0jwDfn2Cc5DvWO0NE5yG9dujbpUEREIovyZ+tr5th3+UIHUg7aO4KlFDS/kIiUkrwtAjN7D/DHwEYzezDnrUbgv+IOrBS17c3S3FDNr5y+LOlQREQiK9Q19DWCovBngI/m7O9399QvSj+XtkyWLRtWUFFhSYciIhJZ3q4hd+91933u/nZ3f4pg8RgHlpqZOsFneaZriAM9w1x2trqFRKS0RLmz+I1m9gTwJPBzYB8aPvocbZmwPrBRiUBESkuUYvFfAluAx919A/BqYFusUZWgtkyWlUtrOfu0pUmHIiIyL1ESwbi7Z4EKM6tw958Cm2OOq6S4O+2ZLJdtWoGZ6gMiUlqi3FncY2ZLgV8At5nZEWAw3rBKS+boIEf6RzVsVERKUpQWwVXAEPAnwPeBDPDGOIMqNe1hfUDzC4lIKYoyxcT0X/9TwK3xhlOa2jJZVjfXs7alIelQRETmTRPinKKpKWdbR5atqg+ISIlSIjhFjx7qp3toXMNGRaRk5U0EZvYf4fNnFy+c0jNz/4DqAyJSogrVCM4ws8uAK83sDuC4fg933xVrZCViW0eWDSuXcGZzfdKhiIiclEKJ4BPAx4E1wOdnvefAq+IKqlRMTE5xT0cXb3iRlqUUkdKVNxG4+53AnWb2cXf/9CLGVDIeeraP/tEJDRsVkZIWZfjop83sSuAV4a6fuft34w2rNLRnsgBsUaFYREpYlEnnPgN8AHgkfHzAzP46wuduNrMjZvZQnvfNzL5oZnvN7EEzu2i+wSetLdPJuauW0tpYm3QoIiInLcrw0dcDr3H3m939ZuB1wBsifO6W8Nh8LgfOCR/XAl+JcM6iMTYxxY593Vy2aWXSoYiInJKo9xE052w3RfmAu/8CKLSAzVXAVz2wDWg2szMixpO4B/b3MDw+qWGjIlLyokw69xngPjP7KcEQ0ldw/IplJ2s18EzO6/3hvoMLcO7Yte3NYgZbNigRiEhpi1Isvt3MfgZcHO76iLsfijWqWczsWoLuI9auLY7F0doynZx/5jKaGqqTDkVE5JRE6hpy94Puflf4WKgkcAA4K+f1mnDfXN9/o7tvdvfNra2tC/T1J29kfJL7nu5RfUBEykKScw3dBfxOOHpoC9Dr7iXRLbTzqW7GJqdUHxCRshClRnBSzOx24JXASjPbD1wHVAO4+w3A3cAVwF6C9Q7eHVcsC60t00llhXHx+pakQxEROWUFE4GZVQIPu/vz53tid3/7Cd534L3zPW8xaMtkedGaJpbWxpZHRUQWTcGuIXefBB4zs+Ko0BaBgdEJHtzfq/qAiJSNKH/SLgceNrN7yVmr2N2vjC2qIrb9yS4mp1z1AREpG1ESwcdjj6KEtGU6qams4CXrlicdiojIgohyH8HPzWwdcI67/9jMGoDK+EMrTu0dWS5a10xddWovgYiUmSiTzv0BcCfwD+Gu1cC3Y4ypaPUMjfHws31s3aj6gIiUjyj3EbwXeCnQB+DuTwCnxRlUsdrW0YU7XHa26gMiUj6iJIJRdx+bfmFmVQQrlKVOe6aT+upKXrSmOelQREQWTJRE8HMz+3Og3sxeA3wD+E68YRWn9o4sF29ooaYqyRuyRUQWVpTfaB8FjgK7gT8kuCP4L+IMqhgd7R/l8cMDbNVqZCJSZqKMGpoys1uBewi6hB4L7wpOlfaOYFlKrU8sIuXmhInAzF4P3ABkCNYj2GBmf+ju34s7uGLSnsnSWFfF+WcuSzoUEZEFFeWGss8Bv+buewHMbBPw70DKEkEnl25ooapS9QERKS9Rfqv1TyeBUAfQH1M8RenZnmH2ZYfYqvmFRKQM5W0RmNmbw80dZnY38C8ENYK3AtsXIbai0Z5RfUBEylehrqE35mwfBn413D4K1McWURFqy2RZ3lDN81Y1Jh2KiMiCy5sI3L1kFoqJk7vTnulk66YVVFRY0uGIiCy4KKOGNgDvB9bnHp+Waaif7hri2d4R3qP6gIiUqSijhr4N3ERwN/FUrNEUobawPqAbyUSkXEVJBCPu/sXYIylSbZkspzXWsql1SdKhiIjEIkoi+IKZXQf8EBid3unuu2KLqkgE9YEsLzt7BWaqD4hIeYqSCC4A3gG8imNdQx6+Lmt7jwzQOTCq9YlFpKxFSQRvBTbmTkWdFjP1Ad0/ICJlLMqdxQ8BzTHHUZTaMp2sWV7PWS0NSYciIhKbKC2CZuBRM9vO8TWCsh4+OjXlbOvo4rXnr0o6FBGRWEVJBNfFHkUReuRgH73D4+oWEpGyF2U9gp8vRiDFpn3m/gEVikWkvEW5s7ifY2sU1wDVwKC7l/XE/O0dWTa2LuH0prqkQxERiVWUFsHMTGsWDKa/CtgSZ1BJG5+c4p6OLG+6cHXSoYiIxG5eq6x44NvAa+MJpzjsPtDL4Nik7h8QkVSI0jX05pyXFcBmYCS2iIrAdH1gy8aWhCMREYlflBbBG3MeryVYneyqKCc3s9eZ2WNmttfMPjrH++8ys6Nmdn/4+P35BB+X9kyW55/eyIqltUmHIiISuyg1gpNal8DMKoEvA68B9gPbzewud39k1qFfd/f3ncx3xGF0YpLt+7r4rUvXJh2KiMiiKLRU5ScKfM7d/dMnOPclwF537wjPdwdBS2J2Iigq9z3dw+jElOoDIpIahbqGBud4APwe8JEI514NPJPzen+4b7a3mNmDZnanmZ0114nM7Foz22FmO44ePRrhq09eeyZLhcElG1QfEJF0yJsI3P1z0w/gRoJ1it8N3AFsXKDv/w6w3t1fCPwIuDVPLDe6+2Z339za2rpAXz239kyWF6xuoqm+OtbvEREpFgWLxWbWYmZ/CTxI0I10kbt/xN2PRDj3ASD3L/w14b4Z7p519+n5i/4v8JLIkcdgeGyS+57p1rQSIpIqeROBmf0tsJ1glNAF7v5Jd++ex7m3A+eY2QYzqwGuAe6a9R1n5Ly8Etgzj/MvuB1PdTE+6aoPiEiqFBo19CGC2Ub/AvhYzgpdRlAsLjjFhLtPmNn7gB8AlcDN7v6wmX0K2OHudwH/w8yuBCaALuBdp/Ifc6raMlmqKozN65YnGYaIyKLKmwjcfV53Hec5x93A3bP2fSJn+8+APzvV71kobZksLz6rmSW1USZlFREpD6f8y75c9I2Ms3t/D5epPiAiKaNEENr+ZBdTDltVHxCRlFEiCLVlstRUVXDh2uakQxERWVRKBKG2TJbN65ZTV12ZdCgiIotKiQDoHhxjz8E+1QdEJJWUCIBtHeGylEoEIpJCSgQE3UINNZW8cE1z0qGIiCw6JQKC9Ykv2dBCdaUuh4ikT+p/8x3pG2HvkQG2blS3kIikU+oTQXtYH9D8QiKSVkoEmSzL6qo478yCUyeJiJSt1CeCtkyWSzeuoLLCTnywiEgZSnUieKZriKe7hnT/gIikWqoTgeoDIiIpTwTbMllWLKnh3FVLkw5FRCQxqU0E7k5bJsuWTSvIWXRHRCR1UpsInuwc5FDfiOoDIpJ6qU0Eqg+IiARSmwjaMllOX1bH+hUNSYciIpKoVCYCd2dbJstlqg+IiKQzETx+eIDs4JimnRYRIaWJoC3TCWj9ARERSG0iyLK2pYE1y1UfEBFJXSKYnHK2dWQ1bFREJJS6RPDIs330j0yoW0hEJJS6RDBTH9BCNCIiQCoTQZazT1vKacvqkg5FRKQopCoRjE9OsX1fl+oDIiI5UpUIHtzfw9DYpLqFRERypCoRtO0N5hfaokQgIjIj1kRgZq8zs8fMbK+ZfXSO92vN7Ovh+/eY2fo442nvyHLeGctYvqQmzq8RESkpsSUCM6sEvgxcDpwHvN3Mzpt12O8B3e5+NvB3wGfjimdkfJIdT3Vr2KiIyCxxtgguAfa6e4e7jwF3AFfNOuYq4NZw+07g1RbTLHC7nu5mbGJKhWIRkVniTASrgWdyXu8P9815jLtPAL3Ac35Tm9m1ZrbDzHYcPXr0pIKprqzg157XysUbWk7q8yIi5aokisXufqO7b3b3za2trSd1jovXt/D/3n0Jy+qqFzg6EZHSFmciOACclfN6TbhvzmPMrApoArIxxiQiIrPEmQi2A+eY2QYzqwGuAe6adcxdwDvD7auBn7i7xxiTiIjMUhXXid19wszeB/wAqARudveHzexTwA53vwu4CfgnM9sLdBEkCxERWUSxJQIAd78buHvWvk/kbI8Ab40zBhERKawkisUiIhIfJQIRkZRTIhARSTklAhGRlLNSG61pZkeBp07y4yuBzgUMJ06lEqviXHilEqviXFhxx7nO3ee8I7fkEsGpMLMd7r456TiiKJVYFefCK5VYFefCSjJOdQ2JiKScEoGISMqlLRHcmHQA81AqsSrOhVcqsSrOhZVYnKmqEYiIyHOlrUUgIiKzKBGIiKRcahKBmb3OzB4zs71m9tGk45lmZmeZ2U/N7BEze9jMPhDu/6SZHTCz+8PHFUUQ6z4z2x3GsyPc12JmPzKzJ8Ln5UUQ5/Nyrtv9ZtZnZh8shmtqZjeb2REzeyhn35zX0AJfDP/NPmhmFyUc59+a2aNhLN8ys+Zw/3ozG865rjcsVpwFYs37szazPwuv6WNm9tqE4/x6Toz7zOz+cP/iXlN3L/sHwTTYGWAjUAM8AJyXdFxhbGcAF4XbjcDjwHnAJ4EPJx3frFj3AStn7fsb4KPh9keBzyYd5xw/+0PAumK4psArgIuAh050DYErgO8BBmwB7kk4zt8AqsLtz+bEuT73uCK5pnP+rMP/tx4AaoEN4e+FyqTinPX+54BPJHFN09IiuATY6+4d7j4G3AFclXBMALj7QXffFW73A3t47trOxewq4NZw+1bgTcmFMqdXAxl3P9m70ReUu/+CYO2NXPmu4VXAVz2wDWg2szOSitPdf+jB2uIA2whWHUxcnmuaz1XAHe4+6u5PAnsJfj/ErlCcZmbA24DbFyOW2dKSCFYDz+S83k8R/rI1s/XAhcA94a73hc3wm4uhywVw4IdmttPMrg33rXL3g+H2IWBVMqHldQ3H/89VbNcU8l/DYv53+7sErZVpG8zsPjP7uZm9PKmgZpnrZ12s1/TlwGF3fyJn36Jd07QkgqJnZkuBfwU+6O59wFeATcCLgYMEzcakvczdLwIuB95rZq/IfdODNm3RjEcOl0i9EvhGuKsYr+lxiu0azsXMPgZMALeFuw4Ca939QuB/Al8zs2VJxRcq+p/1LG/n+D9YFvWapiURHADOynm9JtxXFMysmiAJ3Obu3wRw98PuPunuU8A/skjN10Lc/UD4fAT4FkFMh6e7K8LnI8lF+ByXA7vc/TAU5zUN5buGRffv1szeBbwB+O0waRF2s2TD7Z0E/e7nJhYkBX/WxXhNq4A3A1+f3rfY1zQtiWA7cI6ZbQj/SrwGuCvhmICZvsGbgD3u/vmc/bl9wf8NeGj2ZxeTmS0xs8bpbYLC4UME1/Gd4WHvBP4tmQjndNxfWcV2TXPku4Z3Ab8Tjh7aAvTmdCEtOjN7HfC/gCvdfShnf6uZVYbbG4FzgI5kopyJKd/P+i7gGjOrNbMNBLHeu9jxzfLrwKPuvn96x6Jf08WqSif9IBiB8ThBZv1Y0vHkxPUygq6AB4H7w8cVwD8Bu8P9dwFnJBznRoLRFg8AD09fQ2AF8B/AE8CPgZakr2kY1xIgCzTl7Ev8mhIkpoPAOEH/9O/lu4YEo4W+HP6b3Q1sTjjOvQT969P/Tm8Ij31L+G/ifmAX8MYiuKZ5f9bAx8Jr+hhweZJxhvtvAf5o1rGLek01xYSISMqlpWtIRETyUCIQEUk5JQIRkZRTIhARSTklAhGRlFMiEBFJOSUCSS0zW5Ezze+hnGmLB8zs+hi+7xYze9LM/ijP+wPh86bpOBY6BpG56D4CEYL564EBd//fMX7HLcB33f3OPO8PuPvSfK9F4qIWgcgsZvZKM/tuuP1JM7vVzH5pZk+Z2ZvN7G8sWKDn++E8UZjZS8JZInea2Q+iTBcdTnnSHp7rL+P+7xLJR4lA5MQ2Aa8imMn0n4GfuvsFwDDw+jAZ/D1wtbu/BLgZ+KsI5/0C8JXwXInNISRSlXQAIiXge+4+bma7CVY8+364fzfBSlLPA14A/CiYQ5BKov1ifynBnDIQzI3z2QWMWSQyJQKRExsFcPcpMxv3Y4W1KYL/hwx42N23nsS5VaSTxKlrSOTUPQa0mtlWCNaXMLPzI3zuvwimRAf47biCEzkRJQKRU+TBOthXA581swcIpg6+LMJHP0Cw0ttuimO5REkpDR8VWSQnGj46x/EaPiqLQi0CkcXTC3w63w1l06ZvKAMOL0pUknpqEYiIpJxaBCIiKadEICKSckoEIiIpp0QgIpJy/x90PU0BFQzRMAAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "pyplot.plot(time, xe135, label=\"Xe135\")\n",
    "pyplot.xlabel(\"Time [d]\")\n",
    "pyplot.ylabel(\"Number of atoms - Xe135\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also examine reaction rates over time using the `Results` object."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "_, u235_fission = results.get_reaction_rate(\"1\", \"U235\", \"fission\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAERCAYAAABsNEDqAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAAwVElEQVR4nO3dd3wVdfb/8de5KYReI9IREBCQEgKIIsquXVyauxYUFRX76qr7/bpf64q6uqv+1rYqiF1cXRUbirqWRVREelFERMSgFFGaSEg5vz/uxI0x5QZyM/cm7+fjcR+5mTsz99wh5J2Zz8wZc3dERERiFQm7ABERSS4KDhERqRQFh4iIVIqCQ0REKkXBISIilaLgEBGRSqk1wWFmD5rZBjNbGsO8Q81svpnlm9nxxab3NbMPzGyZmS02sxPiW7WISOKpNcEBPAwcFeO8a4DTgaklpu8Axrl7z2BdfzezJlVUn4hIUkgNu4Dq4u4zzaxj8Wlm1hm4B8gkGgpnu/tyd18dvF5YYh0rij3/2sw2BMtujmvxIiIJpNYERxkmAee6+2dmNgj4B/CrWBY0s4FAOvB5HOsTEUk4tTY4zKwBcCDwLzMrmlwnxmVbAY8Bp7l7YUXzi4jUJLU2OIiO72x2976VWcjMGgHTgSvdfXY8ChMRSWS1aXD8Z9x9K/CFmf0WwKL6lLeMmaUD04BH3f2ZaihTRCThWG3pjmtmTwKHAi2A9cC1wFvAvUArIA34p7tfb2YDiAZEU2AnsM7de5rZKcBDwLJiqz7d3RdW1+cQEQlbrQkOERGpGrX2UJWIiOyeWjE43qJFC+/YsWPYZYiIJJV58+Z96+6ZJafXiuDo2LEjc+fODbsMEZGkYmZfljZdh6pERKRSFBwiIlIpCg4REakUBYeIiFSKgkNERCpFwSEiIpWi4BARkUpRcJShsNB5+qOvmLF0XdiliIgklFpxAeDucOCx2V+ybutODuzSnEYZaWGXJCKSELTHUYaUiHHjqF58uz2X219fUfECIiK1hIKjHL3bNmHcAR149IPVLMnZEnY5IiIJQcFRgcuO7EbzBnW48vklFBSqBb2IiIKjAo0y0rh6eA8W52zhiQ9L7fclIlKrKDhicFzvVgzp0oK/zfiUDdt2hl2OiEioFBwxMDOuH9GT3PxCbnj5k7DLEREJlYIjRp0yG3DeoZ15cdHXzPrs27DLEREJjYKjEs47tDMdm9fj6heWsjOvIOxyRERCoeCohIy0FK4f0Ysvvv2B+/+zKuxyRERCoeCopKFdMxneuxX3vLOS1d/+EHY5IiLVTsGxG64e3oM6KRGufmEp7rq2Q0RqFwXHbmjZKIPLjujKu599y/Ql34RdjohItVJw7KZTB3dk/zaNuf6lj9m2My/sckREqo2CYzcVNUHcuD2X29QEUURqEQXHHujdtgmnBk0Ql65VE0QRqR0UHHvosiO60ax+Ha6cpiaIIlI7KDj2UOO6aVw9fD8W5WxhqpogikgtENfgMLPVZrbEzBaa2dxSXh9hZouLXjezISVeb2RmOWZ2d7FpJwXrXGxmM8ysRTw/Qyx+06c1B3Vpzl9fUxNEEan5qmOPY5i793X37FJeexPo4+59gfHAAyVenwjMLPrGzFKBO4J19gYWAxfGpepKMDMmjuhFbl4hN05XE0QRqdlCPVTl7tv9v1fQ1Sd6q28AzKw/0BJ4vdgiFjzqm5kBjYCvq6nccnXKbMC5h3bmhYVf895KNUEUkZor3sHhwOtmNs/MJpQ2g5mNMrPlwHSiex2YWQS4Dbj8ZytzzwPOA5YQDYwewJQy1jshOPw1d+PGjVX1ecp1/qGd6dC8Hlc/v5TcfDVBFJGaKd7BMcTds4CjgQvMbGjJGdx9mrt3B0YSPTQFcD7wirvnFJ/XzNKIBkc/oDXRQ1V/Ku2N3X2Su2e7e3ZmZmZVfZ5yZaSlMHFEL1apCaKI1GBxDQ53Xxt83QBMAwaWM+9MoFMw2D0YuNDMVgO3AuPM7GagbzDv58EhrqeBA+P5GSqrqAni3W+rCaKI1ExxCw4zq29mDYueA0cAS0vM0yUYq8DMsoA6wCZ3H+vu7d29I9HDVY+6+xXAWqCHmRXtQhwOJNxo9NXDe5CuJogiUkPFc4+jJTDLzBYBc4Dp7j7DzM41s3ODecYAS81sIXAPcIKX85vW3b8G/gzMNLPFRPdAborjZ9gtLRtlcLmaIIpIDWW14S/i7Oxsnzv3F5eRxFVBoTPinlls2JrLm5cdQsOMtGp9fxGRPWVm80q7lEJXjsdJSsS4ceT+aoIoIjWOgiOO+rRrwimD1ARRRGoWBUecXX6kmiCKSM2i4IiznzVBnLMm7HJERPaYgqMa/NQEccZyNm7LDbscEZE9ouCoBj9vgvhx2OWIiOwRBUc1KWqC+PzCr3lfTRBFJIkpOKpRURPEq9QEUUSSmIKjGmWkpXB90ARxkpogikiSUnBUs0O6ZnJs71bc9fZKvtykJogiknwUHCG4JmiCeM0Ly9QEUUSSjoIjBC0bZXDZEV35z4qNvLJkXdjliIhUioIjJKce0IGerRtx/cvL2LYzL+xyRERipuAISWpKhBtH7c+Gbbnc/oaaIIpI8lBwhKhvuyaMHdSeR95XE0QRSR4KjpD98cjuNKufzpXPL1UTRBFJCgqOkDWum8ZVx/Zg0VebeVJNEEUkCSg4EsCIvq05sHNzblETRBFJAgqOBGBmTBwZbYJ40yufhF2OiEi5FBwJonNmA849pBPTFqxVE0QRSWgKjgRy/rAutG9Wj6teUBNEEUlcCo4EEm2C2JNVG9UEUUQSl4IjwRzabS+O3b8Vd6sJoogkKAVHArp6eA/S1ARRRBJUXIPDzFab2RIzW2hmc0t5fYSZLS563cyGlHi9kZnlmNndxaalm9kkM1thZsvNbEw8P0MY9m6cwaWHR5sgvrpUTRBFJLFUKjjMLGJmjSr5HsPcva+7Z5fy2ptAH3fvC4wHHijx+kRgZolpVwIb3L0r0AP4TyXrSQrjBkebIP75pWVsz80PuxwRkZ9UGBxmNjX4y78+sBT42Mz+WBVv7u7b/b/HYuoDPx2XMbP+QEvg9RKLjQf+Eixf6O418tzVnzVBfF1NEEUkccSyx9HD3bcCI4FXgX2AU2NcvwOvm9k8M5tQ2gxmNsrMlgPTiYYCZhYBbgMuLzFvk+DpRDObb2b/MrOWZax3QnD4a+7GjRtjLDexFDVBfPj9L9QEUUQSRizBkWZmaUSD40V3z6PYnkEFhrh7FnA0cIGZDS05g7tPc/fuwfonBpPPB15x95wSs6cCbYH3g/V+ANxa2hu7+yR3z3b37MzMzBjLTTxFTRCven4phWqCKCIJIJbguB9YTfRQ0kwz6wBsjWXl7r42+LoBmAYMLGfemUAnM2sBDAYuNLPVRINhnJndDGwCdgDPBYv9C8iKpZZkVdQEceFXm3nyIzVBFJHwVRgc7n6nu7dx92OC8Yg1wLCKljOz+mbWsOg5cATRMZLi83QxMwueZwF1gE3uPtbd27t7R6KHqx519yuC938JODRYxa+Bj2P7qMnrpyaIr6oJooiEr9Kn43pULKf5tARmmdkiYA4w3d1nmNm5ZnZuMM8YYKmZLQTuAU4oNlhelv8FrjOzxUTHWi6r7GdINkVNEHfmFfIXNUEUkZBZbbjALDs72+fO/cVlJEnnttc/5a63VjL17EEc2LlF2OWISA1nZvNKu5SizD0OM2sd35Kksi4oaoL4vJogikh4yjtU9YCZzTazm83sUDNLrbaqpFTFmyBOnqkmiCISjjKDw92PIToI/Q4wCphtZs8F10e0r57ypKRDu+3FMfvvzV1vrWTNph1hlyMitVC5g+PuvtPdZ7j7xcFxrsuIXktxt5nNqZYK5ReuGd6T1IhxzYtLa00TxMJC572V3/KHpxZy1N9n8vHXMZ0RLiJxsNuD42aW7u67qrieuKgpg+PFTZn1BRNf/ph7x2Zx9P6twi4nblZt3M6z83OYNn8tX2/ZScOMVOqkRsgrcB47cyC92zYJu0SRGquswXGdVZWk8gsK+c3d7/HdD7v492WH0KBOzRmC2rIjj5cWf82z83NYsGYzEYOhXTMZk9WWw3u0ZOO2XE5+YDabf8jj4fED6N+hWdgli9RICo4aFhwAC9Z8z+h732f8Qftw9fAeYZezR/ILCnn3s295Zn4Ob3y8nl35hXRt2YDj+7dlZN827NUo42fzf7PlR06e/CHrt+5kymkDGNy5eUiVi9RcZQVHpf5MNbOmQDt3X1xllclu69e+KScPbM9D733B6Kw29GzdOOySKm35uq08Oy+H5xd+zcZtuTStl8bJA9szJqstvdo0Imgs8AutGtflqQkHMPaBDzn9oTlMHpfN0K7J25NMJJlUuMdhZu8AvyEaMvOADcB77n5p3KurIjV1jwOih3V+ffs7tG1aj+fOO5BIpPRftIlk0/ZcXlwUPRS1dO1WUiPGsO57cXz/tgzrthfpqbE3NNi0PZdTpszh8w3b+cfYLA7rUWqzZBHZDbt9qMrMFrh7PzM7i+jexrVmttjde8er2KpWk4MDYNqCHP7w1CJuHNWLsYM6hF1OqXblF/L2pxt4Zl4Oby/fQH6h06tNI8ZkteU3fVrTvEGd3V735h27OO3BOSz7eit3ndSvRp8sIFKd9uRQVaqZtQJ+R/Tue5JgRvZtw9Mf5XDLq8s5sufetNiDX8JVyd1ZunYrz87P4YWFa/l+Rx4tGtRh/JB9GJPVlm57N6yS92lSL53HzhrEGQ99xIVPLuD2gkJG9G1TJesWkV+KJTiuB14DZrn7R2bWCfgsvmVJZRQ1QTz6jpnc9Mon3P67vqHWs2HrTp5fuJZn5uWwYv120lMjHN6jJcdnteXgfVuQmlL1t7pvlJHGo+MHctYjc7nkqYXk5hXyuwHtqvx9RERnVdUot772KXe/vZInzz6g2s8y2plXwBsfr+fZ+TnMXLGRQod+7ZswJqstx/VuTeN6adVWx4TH5jFzxUYmjujJqYM7Vsv7itREezLGkQmcDXSk2B6Ku4+v4hrjprYEx868Ao74fzNJSzFevXhopQaZd4e7M3/NZp6Zl8PLi79m2858WjXOYHRWG0ZntaVzZoO4vn9ZcvMLuOCJ+fz7kw1cdex+nHVwp1DqEEl2ezLG8QLwLvBvQC1ZE1hGWgp/HtGTMx76iMnvruKCYV3i8j5rN//ItPk5PDt/LV98+wN101I4qtfejMlqy+DOzUkJ+cyuOqkp/GNsfy55agE3TP+E3PzCuG0LkdooluCo5+7/G/dKpEoMC5og3vnmZxzXuzXtm9erkvXu2JXPjKXreGZeDh+s2oQ7DNqnGecd2plj9m+VcFeup6dGuPPEftRJXczfXvuU3LwC/nB41zKvCxGR2MXyv/1lMzvG3V+JezVSJa4Z3pP/fLqRa19cyoOnD9jtX5aFhc6HX3zHs/NzeHXJN/ywq4D2zepx8a/3ZUxWW9o1q5pQipfUlAi3/rYP6SkR7nxrJTvzC/nT0d0VHiJ7KJbguBj4PzPbBeQF09zdG8WvLNkTezfO4NIjujHx5Y95bdk6jupVuesavtz0A8/OX8tz83PI+f5HGtRJZXjv1ozp35YBHZsm1S/elIjxl9H7k54aYdLMVeTmFXDtcT2T4kJJkURVYXC4e9WcbC/V6rTBHXh2Xg7XvfgxQ/bNrPBQ0tadebyy+BuenZ/DR6u/xwyGdGnB5Ud048iee1M3PaWaKq96kYhx/YieZKRFmPzuF+TmF3LTqP0VHiK7KaYD02b2G2Bo8O077v5y/EqSqpCaEuHGUb0Yfe/7/P2NFVxVShPEguAeF8/My+G1ZevIzS+kc2Z9/ueobozq14ZWjeuGUHl8mBn/d8x+ZKSlcNdbK9mVX8hfj+8dl2tKRGq6CoPDzG4GBgBPBJMuNrOD3P1Pca1M9li/9k05aWB7Hnp/NaOz2tKjdfTo4soN23hm3lqeX7CWdVt30rhuGr/NbsuYrLb0bdckqQ5FVYaZcdkR3aiTGuHW11eQm1/I30/sS5rCQ6RSYrmOYzHQ190Lg+9TgAXqVZUctuzI41e3vUP75vUY3a8Nz8xfy6KvNpMSMQ7tmsmY/m359X57USc1eQ9F7Y7JM1dx4yufcHiPltx9cr9a9/lFYrGnbdWbAN8Fz5Ovd3ct1rheGlceux+XPr2IBWs2033vhlx17H6M6NuGzIaJ0dMqDGcP7USdtAjXvLCMCY/O4/5T+5ORpvAQiUUswfEXYIGZvQ0Y0bGOK+JalVSpUf3akBIxuuzVICnv2REv4wZ3pE5qhCueW8IZD33ElNOzqZeeWNejiCSimHpVBd1xBwTfznH3dTGt3Gw1sI3oFef5JXd5zGwEMBEoBPKBS9x9VrHXGwEfA8+7+4Ulln0R6OTuvSqqozYfqpKKTVuQw2VPLyKrfVMeOmMADTOqp6+WSKIr61BVmaOCZtY9+JoFtAJygkfrYFqshrl739LeHHgT6OPufYHxwAMlXp8IzCylttHA9krUIFKmUf3acvfJWSz8ajOnTJnDlh15FS8kUouVt19+KTABuK2U1xz41Z6+ubsX/+VfP1gvAGbWH2gJzACyi01vUKy2p/e0BhGAY/ZvRXpKhPOfmM9Jk2fz2JkD9+jmUiI1WZl7HO4+IXh6tLsPK/4Ajolx/Q68bmbzzGxCaTOY2SgzWw5MJ7rXgZlFiAbW5aUsMjF4bUd5b2xmE8xsrpnN3bhxY4zlSm12WI+WTD4tm883buekybPZsG1n2CWJJKRYTmB/P8ZppRni7lnA0cAFZja05AzuPs3duwMjiYYCwPnAK+6eU3xeM+sLdHb3aRW9sbtPcvdsd8/OzMyMsVyp7Q7pmslDZwwg5/sfOfH+2Xyz5cewSxJJOOWNcewdHC6qa2b9zCwreBwKxNTdzt3XBl83ANOAgeXMOxPoZGYtgMHAhcHg+q3AuOBCxMFAdjB9FtDVzN6JpRaRWB3YuQWPjh/Ihm25/O7+D/jqu3J3bkVqnTLPqjKz04DTiY4vfET0VFyArcAj7v5cuSs2qw9E3H1b8PwN4Hp3n1Fsni7A5+7uwYD7S0BbL1aUmZ0OZJdyVlVH4GWdVSXxsuirzYx7cA7101N44uwD2KdF/bBLEqlWlT6ryt0fCcYzTnf3XxUb4xhRUWgEWgKzzGwRMAeY7u4zzOxcMzs3mGcMsNTMFgL3ACd4WUkmUs36tGvC1LMHsTO/kBPu/4DP1m8LuySRhBBLy5GbgL+6++bg+6bAZe5+VfzLqxra45A9sWL9NsY+8CGFhc7jZw1iv1a6o4DUDpXe4yjm6KLQAHD374n9rCqRpNe1ZUOemnAA6akRTpo8myU5W8IuSSRUsQRHipn9dEK7mdUFdIK71CqdMhvw9DmDaVAnlZMnz2bel9+HXZJIaGIJjieAN83sTDM7k+gg9yPxLUsk8bRrVo+nzxlM8wbpnDrlQ2av2hR2SSKhqDA43P0W4EZgv+Ax0d3/Gu/CRBJR6yZ1efqcwbRuUpfTH5rDzBW6uFRqn5juYOPur7r75cHjtXgXJZLI9mqUwT8nHEDH5vU565G5vPnJ+rBLEqlWFQaHmR1gZh+Z2XYz22VmBWa2tTqKE0lULRrU4Z8TDqB7q4ac+/g8Ziz9JuySRKpNLHscdwMnAZ8BdYGziF5zIVKrNamXzuNnDaJ32yZcMHUBLyxcG3ZJItUi1kNVK4EUdy9w94eAo+JblkhyaJSRxqPjBzKgY1MueWohT8/9KuySROIuluDYYWbpwEIz+6uZ/SHG5URqhfp1Unno9IEM6dKC/3lmMY/N/jLskkTiKpYAODWY70LgB6Ad0VYhIhKom57C5HHZHLbfXlz9/FIeeHdV2CWJxE2FN1h29y+Di/5aufufq6EmkaSUkZbCP8b25+J/LuCG6Z+Qm1/IBcO6hF2WSJWL5ayq44CFRO/Eh5n1De73LSIlpKdGuOukfozs25q/vfYpt7+xAvXtlJqmwj0O4Dqi99F4B8DdF5rZPnGsSSSppaZEuO13fUlPjXDnm5+Rm1fAFUd3x8wqXlgkCcQSHHnuvqXED73+hBIpR0rEuHl0b9JTI9w/cxW5+YVce1wPhYfUCLEExzIzO5los8N9gd8T+61jRWqtSMSYOKIXGakpPDDrC3LzC7lxZC8iEYWHJLdYzqq6COgJ5AJTgS3AJXGsSaTGMDOuPHY/LhzWhSfnrOHyZxaRX1AYdlkie6TcPQ4zSyF6575hwJXVU5JIzWJmXH5kN+qkRrjtjRXk5hfy9xP6kpaiy6EkOZUbHO5eYGaFZtbY3XX3GpE9cNGv96VOWoSbXlnOrvxC7j65H3VSU8IuS6TSYhnj2A4sMbM3iF4ACIC7/z5uVYnUUBOGdqZOagrXvriMcx6bx32n9CcjTeEhySWW4HgueIhIFTjtwI7USY3wp2lLOPORj5hy2gCFhySVWK4c193+RKrYiQPbk5YS4fJnFnHOY/OYNK6/DltJ0tDonEhIxvRvy82j9+c/KzZywRPz2ZWvs60kOSg4REJ0woD2TBzZi39/soGL/7lAp+pKUohrcJjZajNbYmYLzWxuKa+PMLPFRa+b2ZASrzcysxwzuzv4vp6ZTTez5Wa2zMxujmf9ItXh1AM6cM3wHry6dB1/eHoRBYVqzCCJrcIxDjPrCvwR6FB8fnf/VYzvMczdvy3jtTeBF93dzaw38DTQvdjrE4GZJZa51d3fDu4R8qaZHe3ur8ZYi0hCGj9kH3YVFHLzq8tJSzFuPb6PrjCXhBXLWVX/Au4DJgMFVfnm7r692Lf1KdYDy8z6Ay2JduXNDubfAbwdPN9lZvOBtlVZk0hYzj2kM7vyC7n9jRWkp0S4adT+Cg9JSLEER76737ub63fgdTNz4H53n1RyBjMbBfwF2As4NpgWAW4DTgEOK23FZtYEOA64YzdrE0k4v//1vtGLA99eSVpKhOtH9FRjREk4sQTHS2Z2PjCNaL8qANz9uxiWHeLua81sL+ANM1vu7j879OTu04BpZjaU6KGpw4DzgVfcPae0/zRmlgo8Cdzp7qXeas3MJgATANq3bx9DqSKJ4bIjurKroJBJM1eRnhrhqmP3U3hIQoklOE4Lvv6x2DQHOlW0oLuvDb5uMLNpRO/rUXLMomjemWbWycxaAIOBg4PAagCkm9l2d78imH0S8Jm7/72c954UzEd2drZGGyVpmBl/Oro7u/ILmTLrC9JTI/zPkd0UHpIwYrkAcLdu2mRm9YGIu28Lnh8BXF9ini7A58HgeBZQB9jk7mOLzXM6kF0UGmZ2A9AYOGt36hJJBmbGtcf1IK+gkHvf+Zz0lAh/OLxr2GWJALGdVZUGnAcMDSa9Q3S8Iq+CRVsSPQRV9D5T3X2GmZ0L4O73AWOAcWaWB/wInODl3GfTzNoS7dK7HJgfrPtud3+gos8hkmzMovfzyCso5I43PyM9NaJ7mEtCsIruh2xmDwBpQFHrkVOBAndPmr/4s7Ozfe7cX1xGIpIUCgqdy/+1iGkL1nLVsftx1sEVHiUWqRJmNs/ds0tOj2WMY4C79yn2/VtmtqjqShOR8qREjL8d35td+YXcMP0T0lIinHZgx7DLklosluAoMLPO7v45gJl1ooqv5xCR8qWmRPj7iX3JKyjk2heXkZYS4eRBOltQwhFLcPwReNvMVgFG9AryM+JalYj8QlpKhLtO7se5j83j/6YtIS3F+G12u7DLkloolrOq3jSzfYFuwaRP3T23vGVEJD7qpKZw7yn9OfvRufzPs4tJT40wom+bsMuSWqbMJodm9qvg62iiV3R3CR7HBtNEJAQZaSlMOjWbQfs049KnF/HKkm/CLklqmfL2OA4B3iLa1qMkR3cFFAlN3fQUppw2gNMenMPvn1xAasQ4oufeYZcltUSFp+PWBDodV2qqbTvzOHXKHJZ9vYVJ47IZ1m2vsEuSGqSs03ErvB+HmV0c3BfDzOwBM5tvZkfEp0wRqYyGGWk8Mn4g3fZuyDmPzWPWZ2XdwUCk6sRyI6fx7r6VaMuQ5kQvANQNlEQSROO6aTw2fhCdWtTnrEc/YvaqTWGXJDVcLMFR1FntGOBRd19WbJqIJICm9dN5/KxBtGtaj/EPf8Tc1bE0rxbZPbEExzwze51ocLxmZg0B3RhZJMG0aFCHJ84aRMtGGZz+0Ecs/Gpz2CVJDRVLcJwJXEG09cgOon2rdAGgSALaq1EGU88eRLP66Yyb8iFL124JuySpgWIJjsFEL/rbbGanAFcB+mkUSVCtGtdl6tmDaJiRxilTPmT5uq1hlyQ1TCzBcS+ww8z6AJcBnwOPxrUqEdkjbZvWY+rZg8hITWHs5A9ZuWFb2CVJDRJLcOQH98gYQfTeF/cADeNblojsqQ7N6/PE2YMwM06e/CFffPtD2CVJDRFLcGwzsz8BpwDTzSxCdJxDRBJc58wGPHn2IAoKnZMnz+ar73aEXZLUALEExwlALnCmu68D2gJ/i2tVIlJl9m3ZkMfPGsSPeQWcOGk2azf/GHZJkuQqDA53X+fut7v7u8H3a9xdYxwiSWS/Vo14/MxBbN2Zx8mTZ7Nuy86wS5IkVl533FnB121mtrXYY5uZ6TQNkSTTq01jHh0/kE3bd3Hy5Nls2KbwkN1T3h7HWAB3b+jujYo9Grp7o2qqT0SqUL/2TXnojAGs27qTsZM/ZNN23VpHKq+84JhW9MTMnq2GWkSkGgzo2IwHTstmzXc7OGXKHDbv2BV2SZJkyguO4v2oOsW7EBGpPgd2bsHkcdl8vmE7p06Zw5Yf88IuSZJIecHhZTwXkRpgaNdM7js1i+XrtnL6Q3PYnpsfdkmSJMoLjj5Fg+FAbw2Oi9Q8v+rekrtOymJxzhbOeGgOO3YpPKRiZQaHu6cUGwxP3Z3BcTNbbWZLzGyhmf3iFnxmNsLMFhe9bmZDSrzeyMxyzOzuYtP6B+tcaWZ3mplavIvsgaN67c0dJ/Zl3pffc9Yjc9mZVxB2SZLgYrkAcE8Nc/e+pd1+EHgT6OPufYHxwAMlXp8IzCwx7V7gbGDf4HFU1ZYrUvsM792a23/Xlw9WbeLsRxUeUr7qCI4yuft2/+9Nz+tTbCzFzPoDLYHXi01rBTRy99nBco8CI6uvYpGaa2S/NtwyujfvfvYt5z8xn135uu2OlC7eweHA62Y2z8wmlDaDmY0ys+XAdKJ7HQT9sG4DLi8xexsgp9j3OcG00tY7ITj8NXfjxo17+DFEaoffDWjHDSN78dbyDVz05HzyChQe8kvxDo4h7p4FHA1cYGZDS87g7tPcvTvRPYeJweTzgVfcPafk/LFy90nunu3u2ZmZmbu7GpFa55QDOnDtcT14bdl6/vDUQvIVHlJCajxX7u5rg68bzGwaMJBfjlkUzTvTzDqZWQuiN4862MzOBxoA6Wa2HbiDaJPFIm2BtfH8DCK10RkH7cOu/EL+8upy0lMi/O23fUiJ6DwUiYpbcJhZfSDi7tuC50cA15eYpwvwubu7mWUBdYBN7j622DynA9nufkXw/VYzOwD4EBgH3BWvzyBSm51zSGd25Rdy2xsrSEuJ8JfR+xNReAjx3eNoCUwLzpZNBaa6+wwzOxfA3e8DxgDjzCwP+BE4odhgeVnOBx4G6gKvBg8RiYOLfr0vuwoKueutlaSlGhNH9EJnwItV/Hs6+WVnZ/vcub+4jEREYuDu3Pzqcu6fuYrxB+3D1cP3q1Hh4e5s3J7LV9/tYM13O/j+hzxG9mtDs/rpYZcWOjObV9qlFHEd4xCR5GdmXHF0d3YVFPLge1+Qnhrhf4/qllThsTOvgJzvo8GwZtMO1nz3Y/T5dz/w1Xc/8mOJ61amzPqCe0/JonfbJuEUnOAUHCJSITPjmuE9yCso5L7/fE56aoRLD+8adlk/Kb7X8OWmICC+2/HTXsT6rT9vH18vPYX2zerRoXl9Dt43kw7N69GuWT3aN6vH5h15/P7JBRx/7wdcP6InJw5sH9KnSlwKDhGJiZlx/W96kZfv3PnmZ6SnGBf+at9qe/+deQU/BUHJYFjz3Q525v33tGEz2LtRBu2a1WPovpm0b1aP9sXCoXn99HL3mF66aAgX/3MBVzy3hAVrNvPnET3JSEupjo+ZFBQcIhKzSMS4afT+5BUUcuvrK0hPjTBhaOcqWbe7s3Fb7s+CoSgcvty0gw3bSt9r6Ni8fjQcigVDmyZ19+gXfbP66Tx8xkD+3xsruPvtlXz8zVb+MTaLds3q7enHrBEUHCJSKSkR46/H9ya3oJCbXllOWkqEMw7aJ6Zlf9xVbKwhCISivYavvv/lXkOrYK/hkK6V32vYUykR4/Iju9GnXRMufWohx909iztO7MchXXVBsc6qEpHdkldQyIVT5/PasvXcOKoXYwd1+Gmv4cufBqF/fjiprL2Gnx5VuNdQlb749gfOe3wen67fxqWHdeWCYV1qxTUtZZ1VpeAQkd22K7+Qcx+fx1vLN9A5sz5rN/9Y5l5DaeEQ772GqrRjVz5/em4JLyz8msP224vbfteXxnXTwi4rrhQcCg6RuNiZV8AN0z9mw9bcnwVDh2b1aNO0LnVSE2OvoSq4O49+8CUTX/6YNk3rct8p/dmvVUy3J0pKCg4Fh4hUkXlffsd5j89n6848bh7dm5H9Sm3SnfTKCo5Q78chIpKM+ndoxsu/H0Lvtk245KmFXPvC0lp1/xIFh4jIbtirYQZPnDWIs4bswyMffMlJk2ezfuvOsMuqFgoOEZHdlJYS4arhPbj75H588s1Wjr1zFrNXbQq7rLhTcIiI7KHhvVvzwgUH0ahuKmMf+JAH3l1FTR4/VnCIiFSBfVs25IULDuLw/Vpyw/RPuHDqArbn5oddVlwoOEREqkjDjDTuPSWLPx3dnVeXfsPIe95j5YbtYZdV5RQcIiJVyMw455DOPH7mIL7/YRcj7p7Fq0u+CbusKqXgEBGJgwO7tOCli4awb8uGnPfEfP7yyifkF9SMU3YVHCIicdK6SV2eOucATjmgPffPXMWpU+bw7fbcihdMcAoOEZE4qpOawg0j9+fW3/Zh/prvGX7nLOav+T7ssvaIgkNEpBoc378tz51/IGmpxgn3f8Bjs79M2lN2FRwiItWkZ+vGvHzhwQzp0oKrn1/KZf9axI+7CipeMMEoOEREqlHjemlMOW0Alxy2L9MWrGX0ve+zZtOOsMuqFAWHiEg1i0SMSw7ryoOnDWDt9zsYfte7vL18Q9hlxUzBISISkmHd9+Lliw6mbdN6nPHwR/y/N1ZQWJj44x5xDQ4zW21mS8xsoZn94oYYZjbCzBYXvW5mQ4LpHcxsfjB9mZmdW2yZk4J1LjazGWbWIp6fQUQknto3r8dz5x/ImKy23PHmZ4x/5CM279gVdlnliuuNnMxsNZDt7t+W8XoD4Ad3dzPrDTzt7t3NLD2oLTeYZylwILAB+Bro4e7fmtlfgR3ufl15dehGTiKS6NydqXPWcN2Ly2jZKIP7TulPrzaNQ60pIW/k5O7b/b/JVR/wYPoudy+6SqYO/63Tgkd9i96ouBHRIBERSWpmxthBHXj6nMEUFDpj7n2ff839KuyyShXv4HDgdTObZ2YTSpvBzEaZ2XJgOjC+2PR2ZrYY+Aq4xd2/dvc84DxgCcGeBzCljPVOCA5/zd24cWPVfioRkTjp174pL100hP4dmvLHZxbzf9OWkJufWKfsxvtQVRt3X2tmewFvABe5+8wy5h0KXOPuh5WY3hp4HjgO+A6YAUwAVgF3Aevc/Yby6tChKhFJNvkFhdz6+gru+8/n9GnXhHvHZtG6Sd1qrSGUQ1Xuvjb4ugGYBgwsZ96ZQKeSg93u/jXRMY6Dgb7BtM+DQ1xPEx37EBGpUVJTIlxxdHfuOyWLzzdsZ/hds3hvZanDxdUubsFhZvXNrGHRc+AIogFQfJ4uwVgFZpZFdDxjk5m1NbO6wfSmwBDgU2At0MPMMoNVHA58Eq/PICIStqN6teKFCw+ief10Tp3yIfe+83norUpS47julsC0IBdSganuPqPo1Fp3vw8YA4wzszzgR+CE4Ayr/YDbzMyJDobf6u5LAMzsz8DMYJkvgdPj+BlERELXObMBz19wEP/77GJumbGchV99z62/7UPDjLRQ6onrGEei0BiHiNQE7s6D763mplc+oUOzetx3an+6tmwYt/dLyNNxRUQkdmbGmUP2YepZg9i6M5+R97zHS4uq/4oEBYeISJIZ1Kk5038/hB6tGnHRkwuY+PLH5FXj3QUVHCIiSahlowymnn0Apx/YkSmzvmDs5A/ZsG1ntby3gkNEJEmlp0a47jc9uePEvixZu4Xhd85i7urv4v6+Cg4RkSQ3om8bpl1wIPXSUzhx0mweeu+LuJ6yq+AQEakBuu/diBcuHMKh3fbizy99zCVPLWTHrvy4vJeCQ0SkhmhcN41Jp/bnj0d248VFXzPqnvdZv7Xqxz3ieQGgiIhUs0jEuGBYF/Zv05gnPvySZvXTq/w9FBwiIjXQ0K6ZDO2aWfGMu0GHqkREpFIUHCIiUikKDhERqRQFh4iIVIqCQ0REKkXBISIilaLgEBGRSlFwiIhIpdSKOwCa2Uait5ndHS2AxLhDfPlUZ9VLllpVZ9VLllrjXWcHd//FVYS1Ijj2hJnNLe3WiYlGdVa9ZKlVdVa9ZKk1rDp1qEpERCpFwSEiIpWi4KjYpLALiJHqrHrJUqvqrHrJUmsodWqMQ0REKkV7HCIiUikKDhERqRQFRxnM7Cgz+9TMVprZFWHXU5yZtTOzt83sYzNbZmYXB9OvM7O1ZrYweByTALWuNrMlQT1zg2nNzOwNM/ss+No05Bq7FdtmC81sq5ldkijb08weNLMNZra02LRSt6FF3Rn83C42s6yQ6/ybmS0PaplmZk2C6R3N7Mdi2/a+kOss89/azP4UbM9PzezIkOt8qliNq81sYTC9erenu+tR4gGkAJ8DnYB0YBHQI+y6itXXCsgKnjcEVgA9gOuAy8Our0Stq4EWJab9FbgieH4FcEvYdZb4t18HdEiU7QkMBbKApRVtQ+AY4FXAgAOAD0Ou8wggNXh+S7E6OxafLwG2Z6n/1sH/q0VAHWCf4PdCSlh1lnj9NuCaMLan9jhKNxBY6e6r3H0X8E9gRMg1/cTdv3H3+cHzbcAnQJtwq6qUEcAjwfNHgJHhlfILvwY+d/fd7TRQ5dx9JvBdicllbcMRwKMeNRtoYmatwqrT3V939/zg29lA2+qopTxlbM+yjAD+6e657v4FsJLo74e4K69OMzPgd8CT1VFLSQqO0rUBvir2fQ4J+ovZzDoC/YAPg0kXBocFHgz7EFDAgdfNbJ6ZTQimtXT3b4Ln64CW4ZRWqhP5+X/GRNueRcrahon8szue6N5QkX3MbIGZ/cfMDg6rqGJK+7dO1O15MLDe3T8rNq3atqeCI4mZWQPgWeASd98K3At0BvoC3xDdlQ3bEHfPAo4GLjCzocVf9Oh+dkKcE25m6cBvgH8FkxJxe/5CIm3DspjZlUA+8EQw6Rugvbv3Ay4FpppZo7DqI0n+rYs5iZ//gVOt21PBUbq1QLti37cNpiUMM0sjGhpPuPtzAO6+3t0L3L0QmEw17VKXx93XBl83ANOI1rS+6PBJ8HVDeBX+zNHAfHdfD4m5PYspaxsm3M+umZ0ODAfGBiFHcOhnU/B8HtGxg65h1VjOv3Uibs9UYDTwVNG06t6eCo7SfQTsa2b7BH+Fngi8GHJNPwmOb04BPnH324tNL34sexSwtOSy1cnM6ptZw6LnRAdKlxLdlqcFs50GvBBOhb/ws7/iEm17llDWNnwRGBecXXUAsKXYIa1qZ2ZHAf8D/MbddxSbnmlmKcHzTsC+wKpwqiz33/pF4EQzq2Nm+xCtc05111fCYcByd88pmlDt27O6RuGT7UH07JQVRJP7yrDrKVHbEKKHJhYDC4PHMcBjwJJg+otAq5Dr7ET0jJRFwLKi7Qg0B94EPgP+DTRLgG1aH9gENC42LSG2J9Ew+wbII3qM/cyytiHRs6nuCX5ulwDZIde5kugYQdHP6X3BvGOCn4mFwHzguJDrLPPfGrgy2J6fAkeHWWcw/WHg3BLzVuv2VMsRERGpFB2qEhGRSlFwiIhIpSg4RESkUhQcIiJSKQoOERGpFAWHiIhUioJDJEZm1rxY2+p1xdpwbzezf8Th/R42sy/M7NwyXt8efO1cVEdV1yBSGl3HIbIbzOw6YLu73xrH93gYeNndnynj9e3u3qCs70XiRXscInvIzA41s5eD59eZ2SNm9q6ZfWlmo83srxa9mdWMoMcYZtY/6GI6z8xei6X1edAC54NgXTfE+3OJlEXBIVL1OgO/Itpp93HgbXffH/gRODYIj7uA4929P/AgcGMM670DuDdYV2j9p0RSwy5ApAZ61d3zzGwJ0TsKzgimLyF6p7ZuQC/gjWi/SlKILQgOItqTCKK9lW6pwppFYqbgEKl6uQDuXmhmef7fgcRCov/nDFjm7oN3Y90alJTQ6VCVSPX7FMg0s8EQvbeKmfWMYbn3iLb4Bxgbr+JEKqLgEKlmHr2P/fHALWa2iGgr7ANjWPRiondRXEJi3L5UaimdjiuSoCo6HbeU+XU6rlQL7XGIJK4twMSyLgAsUnQBILC+WqqSWk97HCIiUina4xARkUpRcIiISKUoOEREpFIUHCIiUin/HzcM+mXlWVCZAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "pyplot.plot(time, u235_fission)\n",
    "pyplot.xlabel(\"Time [d]\")\n",
    "pyplot.ylabel(\"Fission reactions / s\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Helpful tips\n",
    "\n",
    "Depletion is a tricky task to get correct. Use too short of a time step, and you will incur a steep cost in computational time due to excessive transport simulations. Use too long of a time step, and you may miss physics happening on shorter timescales, leading to incorrect answers. Consider the xenon plot from above. Xenon-135 is a fission product with a thermal absorption cross section on the order of millions of barns, but it has a half life of ~9 hours. Taking smaller time steps at the beginning of your simulation to build up some equilibrium in your fission products is highly recommended.\n",
    "\n",
    "When possible, differentiate materials that reappear in multiple places. If we had built an entire core with the single `fuel` material, every pin would be depleted using the same averaged spectrum and reaction rates, which is incorrect. Pins experiencing different flux will deplete at different rates. The `Operator` can differentiate these materials using the `diff_burnable_mats` argument, but note that the volumes will be copied from the original material.\n",
    "\n",
    "Using higher-order integrators, like the `CECMIntegrator`, `EPCRK4Integrator` with a fourth order Runge-Kutta, or the `LEQIIntegrator`, can improve the accuracy of a simulation, or at least allow you to take longer depletion steps between transport simulations with similar accuracy.\n",
    "\n",
    "Fuel pins with integrated burnable absorbers, like gadolinia, experience strong flux gradients until the absorbers are mostly burned away. This means that the spectrum and magnitude of the flux at the edge of the fuel pin can be vastly different than that in the interior. The helper `pin` function can be used to subdivide regions into equal volume segments, as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [],
   "source": [
    "div_surfs_1 = [openmc.ZCylinder(r=1)]\n",
    "div_1 = openmc.model.pin(div_surfs_1, [fuel, water], subdivisions={0: 10})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x7f8a8c4de9d0>"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQMAAADrCAYAAABtsfOQAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAAl+klEQVR4nO2de7BmVXXgf0tIA2IpL4XmIQ9t8UEyID0osTQqDRKLNGTio3FGIIPFWAODVZaOdDkkFuAISg2GChOlgIBkioZQZWwaEqp5VSyTNrQT5DlAA8NAp1sIIFUKgYhr/jjn3Hvud89rv87Z59z9o25x73de+7v97d9da+19zhZVJZFIJF43dAMSiUQcJBkkEgkgySCRSOQkGSQSCSDJIJFI5CQZJBIJwJMMROQqEXlGRO6v2S4icqmIbBGRe0XkvaVtp4rIo/nXqT7ak0gkzPEVGVwNHN+w/XeBFfnXGcCfAYjIHsAfA+8DjgL+WER299SmRCJhgBcZqOrfAs837HIi8D3N2ATsJiLLgY8BG1X1eVV9AdhIs1QSiUQg+qoZ7Ac8Vfr56fy1utcTiUTP7Dh0A7oiImeQpRjs8Bu7HvmGvQ4duEVToX46+rJdX/V6pVd/uaxhq3i91lLmxW3/+59V9c2mx/Ulg63AAaWf989f2wp8eOb1u6pOoKqXA5cD7LbvkfrB038Uop0TZ3HHf+vRTwzQjoWowlObDq7YkgRhw4YLdnnS5ri+ZLAeOEtE1pEVC19U1W0icivw30tFw+OAtT21aQmwsPPH0PGrEFnctmpBJDmExIsMROQ6sr/we4nI02QjBL8BoKrfAW4BPg5sAV4C/jDf9ryInA/cnZ/qPFVtKkQmGhlH5+/CrCCSHMLjRQaqenLLdgXOrNl2FXCVj3YsTeYFEK7zu97m7t5p2+WQxODKaAqIiVmyDupHAM2dffURh9mfWeGmeyrnopUw78hlOSwUQ5KCLUkGo8JHFLC447t09jZEms9fL4vunboQQ4oW3EgyGAUuUcDCzh+y49tQJYtqQbR37BQtuJFkEC0uUcD8sbF1/i7MCmKxHLqLIUUL3UkyiA7bKGDcAmiiLAdTMaRooTtJBtFgIwH/Aujr+bhi2Rd9iCFJoZokgyhQAwn4EUBdp7/hstutz2lwdT515qrKLSaSsBXDYikkIUCSwcCYRAPZvrYCmO38/XT6OqTm+osl0VUO9WJolsIB738iRQk5SQaD0I8EygIYtvN3ZVYSC+VgKoYuUkipwzxJBr3TNSWwk8D4BNBEWQ7mYnCTwtITQpJBb3SNBswlMC0B1GEvBlMpLNXUIckgOGEkMJQAjjtl39Z9VJWN124L2Ao7MXSVwlJNHWSMay2O53kGXVICOwmEEECXju6LMMKYF0OXNKJroXFsqcOGC3b5iaquND0uySAY3UQwlAT67Phd8ScIWylMQwhJBtHQJS3oHg34kkCMnb8Ndzl0l0KXKGEsaUOSQRT4iwZ8SGCMAqjDTQw2UhhvlDCoDETkeOBPgB2AK1T1wpntlwAfyX98PfAWVd0t3/YacF++7f+p6uq268UpgzYRdIsGXCUwJQHUYS+GblIwixLiE8JgMhCRHYBHgGPJHnV+N3Cyqj5Ys/9/AY5Q1f+Y//wLVX2DyTXjk0G7CLpGAzYS8C2A0NGi2N6YUIGdGDIpuEYJsQrBVgY+hhaPArao6uMA+UNPTwQqZQCcTPaMxAnQrT4QKhrwIYGqjv/ImnOdz9vEO9adv+g1W0GICMedsq+hFIQbLrutNUoQgd87/LBaIUxtToIPGVQthPK+qh1F5EDgYOCO0ss7i8hm4FfAhar6Vx7a1AP+0oI+JTDb+UN3/CqqrjkrCFM5mEuhmKvQHCUsFEJ23Oz2qcxc7HvS0RrgRlV9rfTagaq6VUQOAe4QkftU9bHZA8uLqOzyxgNmN/eMe1pgEw3YSqAsgCE6fxdm21WWg4kY7KTQHCUsnqzUFiWMUwg+agZHA19T1Y/lP68FUNVvVOz7j8CZqvp3Nee6Gtigqjc2XXPYmoEfEYSWwBgE0BVbMYBpTaG9ljCGOsKQBcQdyQqIx5CtkHQ38BlVfWBmv3cCfwMcnD86nXzxlJdU9RUR2Qv4e+DEuuJjwXAycBOBaTRgKoEpCaAOWzF0l0L7iEPsQhh6aPHjwLfJhhavUtWvi8h5wGZVXZ/v8zVgZ1U9p3TcbwPfBX5Ntgjst1X1yrbrDSMDdxGEigaKf8NeBND149JDPyjE0FUKPqOEmIWQJh0FpUkE7YVCExEMLoGWj8Prt5/TvAOgKC/vc1HzTh77SDgpdBUCVL2hoYSQZBCMZhH4SgsGkUDFP32Xzu5KrSwc+0wYKbilDUMIIckgCG4i8B0NeJFA6Z+7j47flUWCcOg7JlLwHyUML4QkA+/EIwJnCUQqgDp8iaGrFKYmhCQDr4QXgUk0YCWBkQmgDh9ieMe68z1GCfELwVYGrwvRmHEThwhU1U4Emn29fvs5c19jRpC597HL9q/MvT8THllz7tzvs/FaIhz72eWtLbrhsttqHzVfzFisamQxMSlWkgw6068IHllzrpkIZiQwRQox2Eih+H32K4Q64ozGkwwWYLKYSekoTyKwigYCSkA9/+cLVym0RQk+hJDRFB3EJ4RUM5jDLj3oIgKTaKAz+T+bLwFUddbtt5/g5dz5Bdhn1YZFL4uH/HlBXcHgdF1qCe11hPoawlD1g1RAdGJEIvAogbIAvHb87g1YIAhXMdhIYYpCSDKwZjgRGA8ZqrsEBhdAHR7FMCcFAyFA8xBkWCEc0q2hHUkysKZOBuFF0Fc0EK0A6vAghhBRQgghhIgOkgysGIcIXCUwCgHUURKDkxSiEcIDla/7FEKaZ2CM+chBryIojRKYUlTvt99+wrhFACCZzLbfdoLVqIQgRqMObcOP7SMN2ShDPfGOMCxRGdhPLGrCpwhshgonJYFZHKSwaCiyBXchUDnkGPuEpCUqgzrs0wPfIjChTwkUY/RtX8FwlEI/Qqifg9AkhIzhooMlWDPov04QSgShagJNn4lLn3y20znOPvDNtdt8PirdpqZgUkdoqiGEqx+4jS4M+aj0LouonAZ8i+yxaAB/qqpX5NtOBf5b/voFqnqNjzb5ImYR+JLAbOfv2uGbaDrHrCic5JBHCoUUugihiBC6COGRNefWCqGIDuqFkEUInz5rVc12pboBda+HpZdFVHIZrFTVs2aO3QPYDKwk+w38BDhSVV9ouqZ9ZOA/KuhbBL6igfK/u4/O70JZDk5iMIwS+okQ+h9uHHI0YW4RFVV9FSgWUenCx4CNqvp8LoCNwPEe2lTBNETgWhcocvpLn3x27mtoym1xqjnM1BPad/dTQ3CvHyx+fYhiog8ZVC2isl/Ffn8gIveKyI0iUix80PXY3olVBLbMSiBWvEkhQiHUE0cxsa/RhJuAg1T1t8j++hvXBUTkDBHZLCKbX/2l6YfZPCqwxbcIyiMFtu0ZgwRmcZZCSQhtUjAVgi0mowtDRAc+ZLAVKC9xtD/zhUIAVPU5VX0l//EK4Miux5bOcbmqrlTVlct2ra9UdyfcMGKHS3cWgW1aMFYJzOIkBYO0wUQIfaUL+dXaG+QJHzK4G1ghIgeLyDKyJdTWl3cQkfJvaDXwUP79rcBxIrJ7vqDKcflrHjGbaRg8PTAUgSlTkcAszlLwKIS+0oW+ZyY6y0BVfwWcRdaJHwJuUNUHROQ8EVmd73a2iDwgIj8FzgZOy499HjifTCh3A+flr3nCb3owBhFMTQKzlKVgRDRCaEsXFr/eV7rgpWagqreo6jtU9W2q+vX8tT/SfDUlVV2rqu9R1X+jqh9R1f9TOvYqVX17/vXnPtpjS9cnFi0+blgRlKOBXtCWrx6wihICCMGctickDVdMnPAMRLOoIGh6EFgEwSRQ89H48slbGg5RLr5uRfXGQPNozj7wzWbzEzpOUOoyD8Fl/kHVZCQfMxMHnYG4lBlKBIXEvYpgpvM3dfo6BKk8rlISnuRw6ZPPzk1c6iSFPEJoE0KXmYpNMxTbUKXmkevDzEycaGTQT1QwpAi8SaD0z2/T+e0vOyMHT59xoyjBJEJY3rx25KHXX1B97ADRQYoMAuE8jNjCYCLIJdCnAMqUI4gFYnCUQhEl+IwQ8kY2tk1Vre5fMI8OwrGEbmG2iwpqz+YhKhhEBHmR78snbxlMBLMUYvjSyY96KUIajTbkQmhrX1tBscuaDFVnrhpqHGrewQRlYLf2QRWh0wMTnEUQoQRm8SkF0+HHziMMDdR9JmyGGotWLTxP2GHGCcqgCr9RgQ9MogIfIohZArMskoIlnYVgcC9D2y5jjg4mJoNxRAW9iaAUDYwRQZyjBJ9C6JouVB47guhgYjKooseoIDIRjCkaqMNHlGAqhLb2tKUL/UQH/lkCMjDHOirwiA8RTIkFUYIFJjUE13TBf3Rg2ABLJiSD7ilC7FGBtQhGnha04Zo2dPqdmqQLDfiKDvKzLdwzUKowIRlUYX5DUqiowOczC2suMIm0oA3XtMFXupCdrH6TbXSweP/+CokTkYG/wqHl5f2thmwTFQwSDQx4lxL2aYOvdCFMdNA9VQgRHUxEBt1pShGGjgriFcHiTv+hm9cu+vrghrWV+4bCRQgdTj5AdNA9VQjBhKcj+3+kWc1lvE8wMrl2OBHMf/g+dPPaTkeILN5XFX54wjfKe/loXOlsmRAuvm6F0anrpg8v2g+tnao8dyNTw30LXa/TRJEqVN2v4HPa8gQig7ApQrRRQRARVP/ld6EQRHXk4Ie5CMGATulCwOig8XQDpQpeZCAix4vIwyKyRUQW/ZkUkS+KyIP505FvF5EDS9teE5F78q/1s8f6xCZFaD6hn6ig1weTVLcAnwKoo14MngiRLuBeO6g8LsJUwVkG+SIqlwG/C7wbOFlE3j2z2z+SLaLyW8CNwDdL215W1cPzr9V4oacUobUVgUYQvEYFGlQAdRRimJeC4/ks6wd9RAc+HhPQxwSkXhZRUdU7VfWl/MdNZE9Bjp7GFGGoqMCbCOajgSERwVuUYCOEPqKDkKlCvrfx+avocxGVgtOBvy79vHO+HsImETnJ7NJ+JhqFeGZBkKjAowiGiAbq8Bkl2NQPvEUHhvhIFXzWDXotIIrIfyBbV/FbpZcPzJ/K8hng2yLytppjOy6iEkeK0IVhagXDRwN1LIwSHOk5Omi75hhShV4WUQEQkVXAV4HVpQVVUNWt+f8fB+4Cjqi6iP9FVJpxSRHijQriFUGBDyHYpAuu0cHwqYI7fS2icgTwXTIRPFN6fXcR2Sn/fi/gA8CDeKbvFKELRlHBEhFBgVchdGSo0Ry7VKEKd3P0tYjKt4A3AH85M4T4LmBzvrjKncCFWlrKveXKw05BHhXjEUGBt5TBZ3TA8KlC1UV81Q36WkRllaruPTuEqKp/p6q/mS+u8puqeqVjS7zUC0KnCP1GBeMTQYGrELxHBwOkCouuEbBuMIEZiIl6xiuCAi8RgufooE/6bM7kZRBbvaC/qGD8IiiYF4LFsQFqBy6pQhX+6gZuTF4GXhliFCFRop/ooBHHVCFs5OF27pHKIFzxsM9HmzU3JEUFZVzSBdPoIGSH9ffZ8l9EHKkMqoh/slE/k4ymJ4ICl3TBBC+pQkBCFREnJINhSSlCn1h2xJ5ShbGyZGUw1GSjTlinCNONCgpsowObexac8FZE7G9EYdIy8PoU5CGfaJToHZe6gd8iYn8jCpOWgSkhi4fDP8BkaoT9cxny3yqKAnUFSQaTYfopQoFTIdHrg5XimqCUYd+mUcpg2a6vDt2EBXgtHvbytOOlide6QYRFRBGchtxHKYNEIiPGv8x94vf9T0QG8c8xCMvSSREK+ppzECsh5hpMRAaJRMKVJINEogbnackjy2KWpAyMJxw5zjFIw4qR0bGTuvyb2aynYLooq2/6WkRlJxG5Pt/+YxE5qLRtbf76wyLyMR/tGS1pJMECsz+/vc9EHBF9LaJyOvCCqr4duAS4KD/23WTPTHwPcDzwP/PzJTqz9IqHBUu9iOibXhZRyX++Jv/+RuAYyVajPBFYp6qvqOoTwJb8fIlEomf6WkRlbp/8AaovAnt2PDaRSPTAaAqI5UVUXn7x+aGbk0hMjr4WUZnbR0R2BN4EPNfxWGDhIiq7vGkPD81OJBJlellEJf/51Pz7TwB3aDaIux5Yk482HAysAP7BQ5uWFJE90DcxUvpaROVKYE8R2QJ8ETgnP/YB4AayVZT+BjhTVV9zbdOYMb8TTvjhCd8I0pbYUYUfnnDh0M2YDDv6OImq3gLcMvPaH5W+/xfgkzXHfh34uo92dEVVyQYzDI5BEcyOMUbg4utWpLkGAVGUi9etGLoZUTKaAqJPNl67zewAgZf3ucj6eiLC2QeGXyw24ZcvHPQW62MV5eXlZp8ZVTX/bHpkScogkUgsZjIyWOpFtKX+/hPuTEQGwk333D90IwZk6RURl3rxUBVuuucBr+eciAwiYPLP1kuMAZcIcZQyePWXy4ZuwgIEYZ9VG3ydjIuvs6t2p1ShGa8jCQrLV93s51yeUIWnNh1iffwoZZCoYumkCks9RQjFEpBB9z+XIsI71p3fcCaHhTXS8OKocBlWbOPQ6y8Idm4XJi0DEfjUmasqtxk/0spxroEptuJJqcLw2MwxaDrbp8+q/gz7ZtIyaCLI5A6PC3va1Q2mnyrYpgixzzwcesIRTEwGQ/5V9FpEdCRFB4GJsHjog9HKYPEHPv65BqZ1A7tUIYsOpiiEvgqHIesFPggxxwBGKwPhqU0HhzmzaxFx8FQhO3hq6YKLCGJKEUIVD12HFWG0MjCluof6LiKGSBVcRjCmGB0MTkuK4Ld42C+Tl0HTiMIQBRujVMFDdDAFIfQZFQyRItQXD/sbSYAJyqCPD39vqULX69UyfiGkCUb94SQDEdlDRDaKyKP5/3ev2OdwEfl7EXlARO4VkU+Xtl0tIk+IyD351+Em1w9VRGysG3hIFfqLDrITjFUIriLwHhU4pgg+6gWhiofgHhmcA9yuqiuA2/OfZ3kJOEVVi4VSvi0iu5W2f1lVD8+/7ul+adMioqe6wUC43bw0PiGkiKA7PoqH4C6D8uIo1wAnze6gqo+o6qP59/8EPAP0Oi83RN3ANVWwiQ6WihB8iGAMtQKIp14A7jLYW1WLd7Id2LtpZxE5ClgGPFZ6+et5+nCJiOzk2B7AT90gdKpg3iDXdCE7SexCGEIEHU86eIoQmlYZiMhtInJ/xdeCJdTyR5/XfsxEZDlwLfCHqvrr/OW1wDuBfwvsAdQuW1teROXVX86vjmtWN+g5VfAZHcyd0rWt80KISQpFe4ZIDWKfZFQQsl4AHWSgqqtU9bCKrx8AP8s7edHZn6k6h4i8EbgZ+Kqqbiqde5tmvAL8OQ3rLJYXUVm2a9GButcN+k4VQkYHfoRwYTRRQiEBHyIYIiqwPq2HFMFXvQDc04Ty4iinAj+Y3SFfWOX7wPdU9caZbYVIhKze4G0+8dCpwnxD2q9hFB14E0J2siGjBN/RgI0IfEQFU0gRwF0GFwLHisijwKr8Z0RkpYhcke/zKeBDwGkVQ4j/S0TuA+4D9gI8/db6SxV6jw6yE3uoH8yfrBwlhBZD+Rq+ogEIOOU4kpuSQqcI4LiIiqo+BxxT8fpm4HP5938B/EXN8R91uf78ebI0oI0iVbjhstsXbdt47TaOO2Vfswvn0cHrt1eNqJYbmO1b364sOrj0yWfrd6o8rc+FXaTUMZUPblhbap/bmctyiWm4sI+ooPa4iEYRCiYwA7G+bhA8VSiu4yE6GDZdWHzy4q/2bMTQFjlU7Tt/rjAiCJYeeIgKQqYIPusF4Gl5tTjJUoXVRxxWsa36T7XNsmu+ogMrciF86eRHAy79JhWdeGHkUKbvv/xD3pHY101JfaQIADKWGXhldtv3SP3g6T8qvaIc8P4nKsJZrZSBKpWpAlCbKqgqj6w5t7pBCrts/0pjh1SU7befULu9fB3TdAElsBDixFYEvqIC28KhaYpQJ4O6yGDDBbv8RFVX1re8mgmkCeAzVbCSo6eRBbB8cGopZVgKay4U7zOYCDq2YchblX2nCDAZGWT4mIBUN+fAW+0gVF8VuHjdioB1hDgoJBA0NQhYKzB91mFfKQJMSgZ+JiBBuOggWDFxph1TFYJrfaCv9MCOYSYalZmQDOqp79v9Rgctl110LR9CmIIUXNKCgj6nHPuKCvpmcjLomiqEjA58pQvOQphA2uAjLegsgsiigj5TBJicDGwKiZ6jA5N0IbQQ8vaMMUrwEQ1A/yIIHRWEShFgcjLIGDQ6KI5t6XgmU5W9CKEUJcQshbIEYhJBOOKICmCSMognOvBVPyiu6bxWY8RS8CkB8F8jmHpUAJOUQR1+o4O+04Ximl4Wb62QwhBiKF/b53ChkQgCpwdtF48lKoAJT0euu3mp/qam6vnCG6/dxrGfXW4+TXnurM03ExVC2H7bCZ2mK9ve1FRz8fkOmM9inG2bT2aFE2KugG8RuBL7CEKZiUxHrkJ569FPVL5uOkUZLKcpZ5drnaqc7aadhVC+thcpVDdokRwKuryXKkJOFDJOCzqKIFx64GfqcRW205EnGxmAv+ggO6b6JqYiXagVQp4utAnBNEIoru0tSljcoOrO2yCJgr5vHIpRBGNkwpEB2EUHt1HXG62jg+ySwSKEog3BooRIsSoSehQBxBcVwEA3KnVZRCXf77XSU47Wl14/WER+LCJbROT6/BFpXqlzXdXrwYqJ0HmEwfYeBm/FxZEQswjazj7Uw0va6GMRFYCXSwulrC69fhFwiaq+HXgBON2xPTNkw4zmKy/VDzXGLoQvHPSWSUvhCwe9JXoR2BQNfUUFLgRfRKWO/CGoHwWKh6QaHd8ds3kH89FBREJIUnCSQDwiMIsK+hQB9LeIys75mgebROSk/LU9gZ+r6q/yn58G9nNsTy110YFNuuA8VGQghOWrbra+9XkKUrCWAMxJwKcIGs9hKYKh5hXM0tciKgfmBY3PkK21+DbThtYtotLx6JrowC5dAMf6QXbpzrMUXZ+FMEYpOEkAjOYQmIjA9+hBDOlBQS+LqKjq1vz/jwN3AUcAzwG7iUgxvLk/sLWhHRWLqJhhV0wMlC6AvRAmLAUfEhhCBD7Tg6HoYxGV3Ys1FEVkL+ADwIN5JHEn8Imm4/3RXEy0TRf6FoJL2lBuW9HpYhBD0Rbn+wkM0oJs92FFEFNUAP0sovIuYLOI/JSs81+oqg/m274CfFFEtpDVEK50bE8L/tOFvoWQ7e4eJZTbWO6MfcihfD0vNxQZRgPZIX2IoOH6kYkAJj/pqIr6Jyn/3uGH1c5YbJqMBPUTkrLjO0xKyprQaWLSwkPySUo0Ns8LJhObenmyUP7RNb2/oD8RDBMV2E46WoIygKaZibZCaLuZKaQQssP6k8LgWEogO9RfsTDW9GCJPyrdnGoHhqkfZMebpwwmtxUvqieMz/HtlNIBm2ggZhHEwBKVQV0xMdsWon4AhkJYfpFRHWH+0AlKwUEC2eGZBPoRQdNx9SIYslZQsERlAM1CsBtuBI9CyJpoFSVkh1ZIYUxiKLXZRQKmk4ncRWA+jBiDCGDJ1gzK+K8fgMcawnxzrGoJC09RqitAfLWF0kfR9aEjNjMKQ4qgz9GDVEC0pm50IdsWUgjZeQykkAsB3J9CFI0YPAogO112Qp/RAIQRQbEtycABvzKAUQkha5I3KWSnmxFDgW9BVHzUfD12zEYCMD0RQJKBB8IKAfCbNmTN8iqFhaeuEYQDIZ43aCsB6JYWQNuNaXGJAJIMPNEshKqnI0E3IUCgKCFrWjApxEpICUDXEYP4RABpnoEnmhdvrfNml1EGaB9pyM5lMNowd9DCociY1kPwTXmEIFYRNJ87jpGDKpIMKjCdkAT+hXDo9Rd4kcIUxFB+Ly4S6EsEMU8saiLJYBHtE5J8CSFIlJA1c67TjFUMVQKwffhIVwmEFkHMUQGkmkEDdgVF6F5DgG51hOycFrWERSeZry1AfPWFsqx8rHTc9UEk3WcUjkMEqYAYBB9CAB+FxfnzepACLBJDQV+CqIpSfC11bvI0oq7RADAKEUCSQUDahQB1i7KYjTRk5+lZCgtOWi2IMl1l0ZaS+Or4ZUwlAF2eZ9keDUA8IoAkg8A0CSHb7jNtgIGlUHuxdlkUhOjsdYSRALgWCoeqEQwiAxHZA7geOAj4v8CnVPWFmX0+AlxSeumdwBpV/SsRuRr4HeDFfNtpqnpP23X7lwH4EwL4rCXMn79HKUSC6cNJTWoDYJ8WFNuHKhYOJYNvAs+r6oUicg6wu6rW/unI5bEF2F9VX8plsEFVb6w7pophZACuQoCwUcL8NaYrBpunE/uMBrLzxSsCGE4GDwMfVtVt+dOR71LVQxv2PwP4HVX99/nPVzMqGYBfIUBIKWTXGr8YbB9PbioBaI4GsnPGLQIYTgY/V9Xd8u8FeKH4uWb/O4D/oaob8p+vBo4GXiFfnk1VX2m77rAygK5CgPrCIphFCeAmhex64xGDy/oEZhKArtEAxC8CCCgDEbkN2Kdi01eBa8qdX0ReUNW6xVeXA/cC+6rqv5Ze2w4sAy4HHlPV82qOPwM4A2CXNx5w5DFnP9L8zoKTCQGaOrz/KAHcpTB/7Xjk4GNxEhsJgJ9oAOIQAYwgTRCRLwDvUdUzarZ/GPiSqrbeKjd8ZFDGPW0A8ygB/ElhYTvCC8L/qkSmEoCu9xaMIS2YZSgZfAt4rlRA3ENV/2vNvpuAtap6Z+m15blIhGzE4V9UtW4l5znikgH4TxvARgrZ+eOaVRiK8ufWVALQLRqA8YkAhpPBnsANwFuBJ8mGFp8XkZXA51X1c/l+BwE/Ag5Q1V+Xjr8DeDPZJ/+e/JhftF03PhlAuxCyfcyiBDB9wsiUxWAvAOgqgew67TcaxSoCSJOOIqF7HaF5n3xPBynANMTgJgAwlQCMqz5QRZJBVPiLEsBdCrBQDBCvHGY/j3YCABMJZNcddzRQJskgOroLAdqjBPAjhYJZOWRt6FcQVZ89+84/d1bATAIwHRFAkkGkdEkbsv3spQC+nlxaJYhZTO6ZaMO9489dbe67EBKA8YgAkgwiJ4wUIJwY6ugiDPDZ0eswFwBMWwIFSQajoEvqkO1nKgXoXwz9YycA6C6BYt8xSqDAVgY7hmhMoo7skWrtUYLkH1ozKYjMd5LpiMFeAGAuARi3CFxIMugdyT9sXVIHOykU+9WLITt3nCyMVG0EAEkCNiQZDMZCKTR38sVSADsxQJ0c5q/VD9XpqW3nh/lODd2fTjz2lMAnSQaD0zV1yPad/5DbiaHYt6rTNUvCLy6dvoyNAMrHJRHMk2QQBSapw/wxPsSw4Iw1kogNWwGUj00SWEySQVQslgIMI4bY8CEASBJoIskgSqT0oTWJFrJj68TQ/RzDMzvibbNCUYoCzEgyiB7baCE7dmEnWiyH7ucKR9VUF9vlyVIUYE+SwWhwiRbmz7G4k1ULYsFRjrJom9fmY13CFAW4k2QwSlyihcXnau6M7bJoI9QipCkK8EuSwaipjhbmtnoJ/9tk0R+zEUYSgF+cZCAinwS+BrwLOEpVN9fsdzzwJ8AOwBWqemH++sHAOmBP4CfAZ1X1VZc2LV1kpnOEkkN/pM7fL66Rwf3AvwO+W7eDiOwAXAYcCzwN3C0i61X1QeAi4BJVXSci3wFOB/7MsU0JoIscIB5BVNUVUufvFycZqOpD0HqP+1HAFlV9PN93HXCiiDwEfBT4TL7fNWRRRpJBEGblAHWCmDvCsyiaComp4w9PHzWD/YCnSj8/DbyPLDX4uar+qvT6fj20JzFHlSAKlLceXS8KU9I9APHTKoOmRVRU9Qf+m1TbjrlFVIBXNlywy/19XTsgewH/PHQjPJDeR1zULnHYRKsMVNV1svpW4IDSz/vnrz0H7CYiO+bRQfF6XTsuJ1t1CRHZbPPwhthI7yMupvQ+bI57ne+GVHA3sEJEDhaRZcAaYL1mj1i6E/hEvt+pQG+RRiKRWIiTDETk90XkabLFU28WkVvz1/cVkVsA8r/6ZwG3Ag8BN6hqMXD9FeCLIrKFrIZwpUt7EomEPa6jCd8Hvl/x+j8BHy/9fAtwS8V+j5ONNphyucUxMZLeR1ws6fcxygeiJhIJ//RRM0gkEiNgFDIQkU+KyAMi8ut8Ude6/Y4XkYdFZEu+KnRUiMgeIrJRRB7N/797zX6vicg9+df6vttZR9vvV0R2EpHr8+0/zhfcjY4O7+M0EXm29G/wuSHa2YSIXCUiz4hI5RC7ZFyav8d7ReS9rSdV1ei/yO59OBS4C1hZs88OwGPAIcAy4KfAu4du+0wbvwmck39/DnBRzX6/GLqtNr9f4D8D38m/XwNcP3S7Ld/HacCfDt3WlvfxIeC9wP012z8O/DXZE27fD/y47ZyjiAxU9SFVfbhlt7lpz5rd7LQOODF864w4kWzaNfn/TxquKcZ0+f2W39+NwDES3wqvY/ictKKqfws837DLicD3NGMT2ZyexuWwRiGDjlRNe45tevPeqlqsO7Yd2Ltmv51FZLOIbBKRk/ppWitdfr9z+2g2pPwi2ZBxTHT9nPxBHl7fKCIHVGyPHeP+EM3zDGKZ9uxK0/so/6CqKiJ1QzkHqupWETkEuENE7lPVx3y3NVHLTcB1qvqKiPwnsmjnowO3KTjRyEDDTXvulab3ISI/E5HlqrotD9meqTnH1vz/j4vIXcARZHnukHT5/Rb7PC0iOwJvIpt2HhOt70NVy22+gqzWMzaM+8OU0oTKac8Dt2mW9WTTrqFm+rWI7C4iO+Xf7wV8AHiwtxbW0+X3W35/nwDu0LyaFRGt72Mmt15NNnN2bKwHTslHFd4PvFhKUasZuirasXL6+2Q5zyvAz4Bb89f3BW6ZqaA+QvZX9KtDt7vifewJ3A48CtwG7JG/vpLsCVAAvw3cR1blvg84feh2N/1+gfOA1fn3OwN/CWwB/gE4ZOg2W76PbwAP5P8GdwLvHLrNFe/hOmAb8K953zgd+Dzw+Xy7kD1U6LH8c1Q5Clf+SjMQE4kEMK00IZFIOJBkkEgkgCSDRCKRk2SQSCSAJINEIpGTZJBIJIAkg0QikZNkkEgkAPj/rrZKQlWZccMAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 258.065x264.901 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "div_1.plot(width=(2.0, 2.0))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The innermost region has been divided into 10 equal volume regions. We can pass additional arguments to divide multiple regions, except for the region outside the last cylinder."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Register depletion chain\n",
    "\n",
    "The depletion chain we created can be registered into the OpenMC `cross_sections.xml` file, so we don't have to always pass the `chain_file` argument to the `Operator`. To do this, we create a `DataLibrary` using `openmc.data`. Without any arguments, the `from_xml` method will look for the file located at `OPENMC_CROSS_SECTIONS`. For this example, we will just create a bare library."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_lib = openmc.data.DataLibrary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_lib.register_file(\"./chain_simple.xml\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_lib.export_to_xml()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "<?xml version='1.0' encoding='utf-8'?>\n",
      "<cross_sections>\n",
      "  <depletion_chain path=\"chain_simple.xml\" type=\"depletion_chain\" />\n",
      "</cross_sections>\n"
     ]
    }
   ],
   "source": [
    "!cat cross_sections.xml"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This allows us to make an `Operator` simply with the geometry and settings arguments, provided we exported our library to `OPENMC_CROSS_SECTIONS`. For a problem where we built and registered a `Chain` using all the available nuclear data, we might see something like the following."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = openmc.Model(geometry=geometry, settings=settings)\n",
    "new_op = openmc.deplete.CoupledOperator(model, \"./chain_simple.xml\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "len(new_op.chain.nuclide_dict)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['I135', 'Xe135', 'Xe136', 'Cs135', 'Gd157', 'Gd156', 'U234', 'U235', 'U238']"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "[nuc.name for nuc in new_op.chain.nuclides[:10]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['I135', 'Xe135', 'Xe136', 'Cs135', 'Gd157', 'Gd156', 'U234', 'U235', 'U238']"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "[nuc.name for nuc in new_op.chain.nuclides[-10:]]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Choice of depletion step size\n",
    "\n",
    "A general rule of thumb is to use depletion step sizes around 2 MWd/kgHM, where kgHM is really the initial heavy metal mass in kg. If your problem includes integral burnable absorbers, these typically require shorter time steps at or below 1 MWd/kgHM. These are typically valid for the predictor scheme, as the point of recent schemes is to extend this step size. A good convergence study, where the step size is decreased until some convergence metric is satisfied, is a beneficial exercise.\n",
    "\n",
    "We can use the `Operator` to determine our maximum step size using this recommendation. The `heavy_metal` attribute returns the mass of initial heavy metal in g, which, using our power, can be used to compute this step size. $$\\frac{2\\,MWd}{kgHM} = \\frac{P\\times\\Delta}{hm_{op}}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5.080339195584719"
      ]
     },
     "execution_count": 38,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "operator.heavy_metal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [],
   "source": [
    "max_step = 2 * operator.heavy_metal / power * 1E3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\"Maximum\" depletion step:  58.4 [d]\n"
     ]
    }
   ],
   "source": [
    "print(\"\\\"Maximum\\\" depletion step: {:5.3} [d]\".format(max_step))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Alternatively, if we were provided the power density of our problem, we can provide this directly with `openmc.deplete.PredictorIntegrator(operator, time_steps, power_density=pdens)`. The values of `power` and `power_density` do not have to be scalars. For problems with variable power, we can provide an iterable with the same number of elements as `time_steps`."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}