Skip to content
Snippets Groups Projects
check_timeseries_stationarity.ipynb 29.2 KiB
Newer Older
many's avatar
many committed
{
 "cells": [
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Import the necessary modules. We need stationarity from ntrfc and some other modules for the definition of a signal and for rendering a plot"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
many's avatar
many committed
   "execution_count": 1,
many's avatar
many committed
   "outputs": [],
   "source": [
    "from ntrfc.timeseries.stationarity import stationarity\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Lets define a signal generator. We want to generate a signal with some noise, an initial transient and a constant deterministic fluctuation."
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
many's avatar
many committed
   "execution_count": 2,
many's avatar
many committed
   "outputs": [],
   "source": [
    "def signalgen_abatingsine(amplitude, noiseamplitude, frequency, mean, abate, time):\n",
    "    resolution = 2048\n",
    "    step = (resolution * frequency ** -1) ** -1\n",
    "\n",
    "    times = np.arange(0, time, step)\n",
    "    noise = np.random.normal(-1, 1, len(times)) * noiseamplitude\n",
    "\n",
    "    values = amplitude * np.sin(frequency * (2 * np.pi) * times) + mean + np.e ** -(times * abate) + noise\n",
    "    return times, values"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Lets define the input arguments for the signal generator and lets generate a signal."
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
many's avatar
many committed
   "execution_count": 3,
many's avatar
many committed
   "outputs": [],
   "source": [
    "test_amplitudes = 0.1\n",
    "test_noiseamplitude = 0.01\n",
    "test_frequencies = 6\n",
    "test_times = 40\n",
    "test_mean = -1\n",
    "test_abate = 1\n",
    "\n",
    "\n",
    "timesteps, values = signalgen_abatingsine(amplitude=test_amplitudes, noiseamplitude=test_noiseamplitude,\n",
    "                                          frequency=test_frequencies, mean=test_mean, time=test_times,\n",
    "                                          abate=test_abate)\n"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Lets compute the index of the stationary timestep."
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
   "execution_count": 4,
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/many/PycharmProjects/NTRfC/ntrfc/timeseries/stationarity.py:128: RuntimeWarning: invalid value encountered in true_divide\n",
      "  intersection = np.true_divide(np.sum(minima), np.sum(hist_2))\n"
     ]
    }
   ],
many's avatar
many committed
   "source": [
    "stationary_timestep ,_= stationarity(timesteps, values)\n",
many's avatar
many committed
    "\n",
many's avatar
many committed
    "well_computed_stationarity_limit = -np.log(0.05) / test_abate\n",
many's avatar
many committed
    "well_computed_stationary_time = timesteps[-1] - well_computed_stationarity_limit\n",
    "stationary_time = timesteps[-1] - stationary_timestep\n"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
many's avatar
many committed
     "name": "#%%\n"
many's avatar
many committed
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Lets plot the result"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
many's avatar
many committed
   "outputs": [
    {
     "data": {
      "text/plain": "<Figure size 640x480 with 1 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAi8AAAGdCAYAAADaPpOnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/av/WaAAAACXBIWXMAAA9hAAAPYQGoP6dpAABJgUlEQVR4nO3deVhU9f4H8PcMyyC7yK6sboi7qIi5JaikebXMsizTTMukm2mLVmrrxV97luXtesvqara5ZWmaayoKoriDGwqKAyrBsMg2c35/IAMDM8AMc+Yw8H49zzwPc+Z7Zj5nDjPnPd/zPefIBEEQQERERGQl5FIXQERERGQMhhciIiKyKgwvREREZFUYXoiIiMiqMLwQERGRVWF4ISIiIqvC8EJERERWheGFiIiIrIqt1AWYm0ajQVZWFlxcXCCTyaQuh4iIiBpBEAQUFBTA398fcnn9fSstLrxkZWUhICBA6jKIiIjIBJmZmejQoUO9bVpceHFxcQFQufCurq4SV2N5RWVF8P/AHwCQtSALTvZOEldERETUMJVKhYCAAO12vD4tLrxU7SpydXVtleHFpswGcKj829XVleGFiIisSmOGfHDALhEREVkVhhciIiKyKgwvREREZFUYXoiIiMiqMLwQERGRVWF4ISIiIqvC8EJERERWheGFiIiIrArDCxEREVkVhhciIiKyKgwvREREZFUYXoiIiMiqMLwYITE9F2sOX5G6DCIiolatxV1VWkwP/jsBANChrSOGd/GSuBoiIqLWiT0vJjh5NU/qEoiIiFothhcTlKsFqUsgIiJqtRheTKDWMLwQERFJheHFBGqB4YWIiEgqDC8m0LDnhYiISDIMLybgbiMiIiLpMLyYILugVOoSiIiIWi2GFxP8ejxL6hKIiIhaLYYXIiIisioMLyboF+gudQlEREStFsNLIwk1Do8O9nSSsBIiIqLWjeGlkSpqHGEkg0zCSoiIiFo3hpdGqnl4tJzZhYiISDIML41UrtZo/5bLmF6IiIikYpHwsmLFCgQHB8PBwQGRkZFITEyst/1PP/2EsLAwODg4oGfPnvj9998tUWa9dHpeGPmIiIgkI/pm+IcffsD8+fOxdOlSHD16FL1798aYMWOQk5Ojt/3Bgwfx8MMPY+bMmTh27BgmTpyIiRMn4tSpU2KXWq8KnbPqsueFiIhIKqKHlw8//BCzZs3CjBkzEB4ejpUrV8LR0RFfffWV3vaffPIJYmNj8eKLL6Jbt25466230K9fP3z22Wdil1qvmj0vAi/MSEREJBlRw0tZWRmSk5MRExNT/YJyOWJiYpCQkKB3noSEBJ32ADBmzBiD7UtLS6FSqXRuYqg55kXD8EJERCQZUcPLzZs3oVar4ePjozPdx8cHSqVS7zxKpdKo9vHx8XBzc9PeAgICzFN8LTV7XmrkGCIiIrIwqx96umjRIuTn52tvmZmZorxOBXcbERERNQu2Yj65p6cnbGxskJ2drTM9Ozsbvr6+eufx9fU1qr1CoYBCoTBPwfVwa2MHha0cpRUaqBleiIiIJCNqz4u9vT0iIiKwc+dO7TSNRoOdO3ciKipK7zxRUVE67QFgx44dBttbiqezAi/HhgEANMwuREREkhG15wUA5s+fj8cffxz9+/fHwIED8fHHH6OoqAgzZswAAEybNg3t27dHfHw8AOC5557D8OHD8cEHH2DcuHFYt24djhw5gi+//FLsUhtkc+fUuhywS0REJB3Rw8tDDz2EGzduYMmSJVAqlejTpw+2bdumHZSbkZEBeY2zvg0ePBhr167Fa6+9hldeeQWdO3fGxo0b0aNHD7FLbVDVZQE07HohIiKSjOjhBQDi4uIQFxen97E9e/bUmTZ58mRMnjxZ5KqMJ2fPCxERkeSs/mgjS6q6phEPlSYiIpIOw4sRbO6EFx4qTUREJB2GFyNUXUyah0oTERFJh+HFCNVHG0lcCBERUSvG8GKEqjEvPNqIiIhIOgwvRqjabcSjjYiIiKTD8GIEnqSOiIhIegwvRqjebSRxIURERK0Yw4sRtOGFPS9ERESSYXgxgpyHShMREUmO4cUIPFSaiIhIegwvRuCh0kRERNJjeDECD5UmIiKSHsOLEbjbiIiISHoML0bgbiMiIiLpMbwYgYdKExERSY/hxQhVh0qfzynE94kZEBhiiIiILI7hxQhVY14AYNH6k9h/4aaE1RAREbVODC9GkMlkOvdTMvKkKYSIiKgVY3gxglw3u+CDHeekKYSIiKgVY3gxgjK/ROoSiIiIWj2GFyPcLldLXQIREVGrx/BihL6BbaUugYiIqNVjeDFC7TEvREREZHkML0aQ1zraaEgnT4kqISIiar0YXoxQK7tAYcu3j4iIyNK49TVC7Z4XhR3fPiIiIkvj1tcIdcKLrY1ElRAREbVeDC9GqD1g186GI3iJiIgsjeHFCLUvD2Aj59tHRERkadz6GoE9L0RERNJjeDFC7TEv7ZwUElVCRETUejG8GKF2eLlZWCpRJURERK0Xw4sRZLXere8OXZGmECIiolaM4cUItXtefFy524iIiMjSGF6MUHvA7rie/tIUQkRE1IoxvBihds+LRhAkqoSIiKj1YngxQu1rGzG8EBERWR7DixHY80JERCQ9hhcj1A4vao1EhRAREbViDC9GqD1gV2DPCxERkcUxvBih9rWNuNuIiIjI8hhemoC7jYiIiCyP4aUJuNuIiIjI8kQNL7m5uZg6dSpcXV3h7u6OmTNnorCwsN72zz77LLp27Yo2bdogMDAQ//znP5Gfny9mmSbjbiMiIiLLEzW8TJ06FadPn8aOHTuwZcsW7Nu3D7NnzzbYPisrC1lZWXj//fdx6tQprF69Gtu2bcPMmTPFLNNkamYXIiIii7MV64nPnj2Lbdu2ISkpCf379wcAfPrppxg7dizef/99+PvXPbV+jx498Msvv2jvd+zYEe+88w4effRRVFRUwNZWtHJNwp4XIiIiyxOt5yUhIQHu7u7a4AIAMTExkMvlOHz4cKOfJz8/H66urgaDS2lpKVQqlc7NUjjmhYiIyPJECy9KpRLe3t4602xtbeHh4QGlUtmo57h58ybeeuutenc1xcfHw83NTXsLCAhoUt3GOJddiN1pORZ7PSIiIjIhvCxcuBAymazeW2pqapMLU6lUGDduHMLDw/H6668bbLdo0SLk5+drb5mZmU1+7ca6kFOIGV8n4cTVPIu9JhERUWtn9CCSBQsWYPr06fW2CQ0Nha+vL3JydHslKioqkJubC19f33rnLygoQGxsLFxcXLBhwwbY2dkZbKtQKKBQKBpdvxjOXlehVwd3SWsgIiJqLYwOL15eXvDy8mqwXVRUFPLy8pCcnIyIiAgAwK5du6DRaBAZGWlwPpVKhTFjxkChUGDz5s1wcHAwtkSLk0HWcCMiIiIyC9HGvHTr1g2xsbGYNWsWEhMTceDAAcTFxWHKlCnaI42uXbuGsLAwJCYmAqgMLqNHj0ZRURH++9//QqVSQalUQqlUQq1Wi1VqkwngwF0iIiJLEfXY4zVr1iAuLg7R0dGQy+WYNGkSli9frn28vLwcaWlpKC4uBgAcPXpUeyRSp06ddJ4rPT0dwcHBYpZLREREVkDU8OLh4YG1a9cafDw4OFjncOMRI0ZY5eHH3G1ERERkOby2EREREVkVhhciIiKyKgwvREREZFUYXsxAbYXjdIiIiKwVw4sZVGgYXoiIiCyF4cVI/zepZ51parVGgkqIiIhaJ4YXIz00IBBfzxigM83DWdrLExAREbUmDC8mkMt0z+tiK+d5XoiIiCyF4cUEtbOKmmNeiIiILIbhxQS1e140PNqIiIjIYhheTFAruzC8EBERWRDDiwlq97zwYCMiIiLLYXgxgU2tQS8ajnkhIiKyGIYXE9QesMvdRkRERJbD8GICWe3dRgwvREREFsPwYoI6RxtxtxEREZHFMLyYgOd5ISIikg7DiwnqnudFokKIiIhaIYYXE/AkdURERNJheDGBvNa7xt1GRERElsPwYgLuNiIiIpIOw4sJeJ4XIiIi6TC8mKT25QEYXoiIiCyF4cUE7HkhIiKSDsOLCXiSOiIiIukwvJigzlWl2fNCRERkMQwvJqiVXbBi90XMXXtUmmKIiIhaGYYXE9QOLwDw24nrKKvQWL4YIiKiVobhxQS1dxtV4cBdIiIi8TG8mMBQeOEh00REROJjeDFB7UOlq3DgLhERkfgYXkwgM9TzomZ4ISIiEhvDiwkMZBdUcLcRERGR6BheTMAxL0RERNJheDGBoTEv+bfLLVsIERFRK8TwYgJDY15WH7xs2UKIiIhaIYYXUxjYO1RaobZsHURERK0Qw4sJKjT6z6TLMS9ERETiY3gxgbujvd7pDC9ERETiY3gxgY2BEbsML0REROJjeDGjLj4uUpdARETU4jG8mNGNwlKpSyAiImrxGF5M9P7k3nWmHb50S4JKiIiIWheGFxM9ENEB7z7QS2farKGhElVDRETUejC8NIFNrZPVuTjYSVQJERFR6yFqeMnNzcXUqVPh6uoKd3d3zJw5E4WFhY2aVxAE3HPPPZDJZNi4caOYZZpMXuvdUws82oiIiEhsooaXqVOn4vTp09ixYwe2bNmCffv2Yfbs2Y2a9+OPPzZ4Gv7movYFGgWGFyIiItHZivXEZ8+exbZt25CUlIT+/fsDAD799FOMHTsW77//Pvz9/Q3Om5KSgg8++ABHjhyBn5+fWCU2We3zvfA8L0REROITreclISEB7u7u2uACADExMZDL5Th8+LDB+YqLi/HII49gxYoV8PX1bfB1SktLoVKpdG6WUrvnZf6PxxG39qjFXp+IiKg1Ei28KJVKeHt760yztbWFh4cHlEqlwfmef/55DB48GBMmTGjU68THx8PNzU17CwgIaFLdxqgdXgBgy4nrFnt9IiKi1sjo8LJw4ULIZLJ6b6mpqSYVs3nzZuzatQsff/xxo+dZtGgR8vPztbfMzEyTXtsU208bDmFEREQkDqPHvCxYsADTp0+vt01oaCh8fX2Rk5OjM72iogK5ubkGdwft2rULFy9ehLu7u870SZMmYejQodizZ0+deRQKBRQKhTGLYDaGLtBIRERE4jE6vHh5ecHLy6vBdlFRUcjLy0NycjIiIiIAVIYTjUaDyMhIvfMsXLgQTz75pM60nj174qOPPsL48eONLVV0o8J98NWBdKnLICIialVEO9qoW7duiI2NxaxZs7By5UqUl5cjLi4OU6ZM0R5pdO3aNURHR+Pbb7/FwIED4evrq7dXJjAwECEhIWKVajJDV5cWBKHZH+ZNRERkrUQ9z8uaNWsQFhaG6OhojB07FkOGDMGXX36pfby8vBxpaWkoLi4WswzR2Bh493jINBERkXhE63kBAA8PD6xdu9bg48HBwQ2e2K05n/jNUO+KWhDEfWOJiIhaMV7bqAlqX9uoSjPOW0RERFaP4aUJ9J3nBeBuIyIiIjExvDRB7QszVtGw64WIiEg0DC9NYKjnRaOxcCFEREStCMNLExjcbcSeFyIiItEwvDSBoUOluduIiIhIPAwvTWAoo2g4YJeIiEg0DC9NUFym1jud2YWIiEg8DC9NYOjyABzzQkREJB6GlyYwfLQRwwsREZFYGF6awFDPy83CUgtXQkRE1HowvDSBgeyCTSlZli2EiIioFWF4aQK5gfTS3d/VwpUQERG1HgwvTWDowowBHo4WroSIiKj1YHhpAgPZhSepIyIiEhHDSxO4Otjpnb7zbA6u3CqycDVEREStA8NLE7R1stc7/b/70zH8vT2WLYaIiKiVYHghIiIiq8LwQkRERFaF4YWIiIisCsMLERERWRWGFyIiIrIqDC9NtG72IAwIbit1GURERK0Gw0sTDQpth7E9/aQug4iIqNVgeDEDQ1eXJiIiIvNjeDEDmaHrBBAREZHZMbyYgaELNAq8xhEREZHZMbyYgaG9RswuRERE5sfwYgZyA+lFzfRCRERkdgwvZiA3sNtIrWF4ISIiMjeGFzOwMfAu/nvvJSjzSyxbDBERUQvH8GIGhnpePvrzHB75zyELV0NERNSyMbyYgaHwAgCXbhZZsBIiIqKWj+HFDOoLL0RERGReDC9mYGjMCxEREZkfN7tmwDPsEhERWQ7Dixlcz7stdQlEREStBsOLGRy4eEvqEoiIiFoNhhczMHRtIyIiIjI/hhczsDF0cSMiIiIyO4YXMzB0bSMiIiIyP4YXM3BW2EhdAhERUavB8GIGjw0KlroEIiKiVoPhxQycFbZSl0BERNRqiBZecnNzMXXqVLi6usLd3R0zZ85EYWFhg/MlJCRg5MiRcHJygqurK4YNG4bbt5v3eVQaOtjoj9NKlFVoLFMMERFRCydaeJk6dSpOnz6NHTt2YMuWLdi3bx9mz55d7zwJCQmIjY3F6NGjkZiYiKSkJMTFxUEub94dRA529Y95eeq7ZHz85zkLVUNERNSyyQRBEMz9pGfPnkV4eDiSkpLQv39/AMC2bdswduxYXL16Ff7+/nrnGzRoEEaNGoW33nrL5NdWqVRwc3NDfn4+XF1dTX4eY6366xLit6ZCrdH/dnq7KJD4aozodRSVFcE53hkAULioEE72TqK/JhERUVMZs/0WpUsjISEB7u7u2uACADExMZDL5Th8+LDeeXJycnD48GF4e3tj8ODB8PHxwfDhw7F//34xSjS7J4eG4qlhoQYfrzAQaoiIiMg4ooQXpVIJb29vnWm2trbw8PCAUqnUO8+lS5cAAK+//jpmzZqFbdu2oV+/foiOjsb58+cNvlZpaSlUKpXOTSr1nayOZ4IhIiIyD6PCy8KFCyGTyeq9paammlSIRlM5oPWpp57CjBkz0LdvX3z00Ufo2rUrvvrqK4PzxcfHw83NTXsLCAgw6fXNob6rS/PK00REROZh1DG+CxYswPTp0+ttExoaCl9fX+Tk5OhMr6ioQG5uLnx9ffXO5+fnBwAIDw/Xmd6tWzdkZGQYfL1FixZh/vz52vsqlUqyAFPfNY6YXYiIiMzDqPDi5eUFLy+vBttFRUUhLy8PycnJiIiIAADs2rULGo0GkZGReucJDg6Gv78/0tLSdKafO3cO99xzj8HXUigUUCgURiyFeGxtDCeUGwWlFqyEiIio5RJlzEu3bt0QGxuLWbNmITExEQcOHEBcXBymTJmiPdLo2rVrCAsLQ2JiIoDK3Sovvvgili9fjp9//hkXLlzA4sWLkZqaipkzZ4pRptnxAo1ERETiE+3UsGvWrEFcXByio6Mhl8sxadIkLF++XPt4eXk50tLSUFxcrJ02b948lJSU4Pnnn0dubi569+6NHTt2oGPHjmKVaVa2DC9ERESiEy28eHh4YO3atQYfDw4Ohr5TzCxcuBALFy4UqyxRseeFiIhIfM371LVWhj0vRERE4mN4MSMeDk1ERCQ+hhczkjO8EBERiY7hxYxs+G4SERGJjptbM7Jp4OrXPyZloqi0wkLVEBERtUwML2Z09e/ieh9/6ZcTWLLptIWqISIiapkYXsyoXK1psM22U9ctUAkREVHLxfBiRhWauuetqc3Olm85ERFRU3BLakYOtjYNttE0IuAQERGRYQwvFqYq4YBdIiKipmB4MSM1e1WIiIhEx/BiRrY2DZ+krrO3swUqISIiarkYXsxoxuAQdPNzrbcNz8JLRETUNAwvZuTmaIetzw3FpH4dDLZhdiEiImoahhcR1Hdx6VRlAT7fc6FR54QhIiKiuhheRNDQrqF3t6Xhu4QrFqqGiIioZWF4EUFjdg2lKlXiF0JERNQCMbyIoDHhReBR1URERCZheBFFw+mlMZcSICIioroYXkTQmJ6XDceuiV8IERFRC8TwIoL6jjYiIiKipmF4EYGsEbuNiIiIyDQMLyJo7InoTmfli1sIERFRC8TwIoLGHkk0bvl+ZOYWi1sMERFRC8PwIgJjjiQ6e53neyEiIjIGw4sIjBmwK+PFjoiIiIzC8CICY64cfS67QMRKiIiIWh6GFxGkGRFIdqfmiFgJERFRy8PwIoLE9NxGt+V5domIiIzD8EJERERWheFFBP/o7d/otnY2HLBLRERkDIYXEbx2bzeE+bo0qu3AkHYiV0NERNSyMLyIwNvFAdvmDUMXH+cG27LfhYiIyDgMLxLzcLKXugQiIiKrwvAiosZcoNHLRWGBSoiIiFoOhheJFZZUQG3E5QSIiIhaO4YXib30ywl0fOV3aBhgiIiIGoXhRUTGXLZo/4Wb4hVCRETUgjC8iChV2fjLBFy6UYirfxeLWA0REVHLwPDSTLz+6xkM+b/dKClXS10KERFRs8bw0swUlFRIXQIREVGzxvDSzNjIedo6IiKi+jC8EBERkVVheGlmNAIPmSYiIqoPw0szw/O9EBER1U+08JKbm4upU6fC1dUV7u7umDlzJgoLC+udR6lU4rHHHoOvry+cnJzQr18//PLLL2KV2CyVVmggsPeFiIjIINHCy9SpU3H69Gns2LEDW7Zswb59+zB79ux655k2bRrS0tKwefNmnDx5Evfffz8efPBBHDt2TKwym52h7+7GKxtOSl0GERFRsyVKeDl79iy2bduGVatWITIyEkOGDMGnn36KdevWISsry+B8Bw8exLPPPouBAwciNDQUr732Gtzd3ZGcnCxGmaKbPjjYpPm+T8w0byFEREQtiCjhJSEhAe7u7ujfv792WkxMDORyOQ4fPmxwvsGDB+OHH35Abm4uNBoN1q1bh5KSEowYMcLgPKWlpVCpVDq35mLxveFY+2SkSfNWqDVmroaIiKhlECW8KJVKeHt760yztbWFh4cHlEqlwfl+/PFHlJeXo127dlAoFHjqqaewYcMGdOrUyeA88fHxcHNz094CAgLMthxNZSOXYXAnT5Pm7f3Gdry95YyZKyIiIrJ+RoWXhQsXQiaT1XtLTU01uZjFixcjLy8Pf/75J44cOYL58+fjwQcfxMmThseALFq0CPn5+dpbZmbL2OVSVKbGqv3pUpdBRETU7Nga03jBggWYPn16vW1CQ0Ph6+uLnJwcnekVFRXIzc2Fr6+v3vkuXryIzz77DKdOnUL37t0BAL1798Zff/2FFStWYOXKlXrnUygUUCgUxiwGERERWTGjwouXlxe8vLwabBcVFYW8vDwkJycjIiICALBr1y5oNBpERuofA1JcXHlFZblctzPIxsYGGk3rHf+xJy0Hg0LbwcHORupSiIiImgVRxrx069YNsbGxmDVrFhITE3HgwAHExcVhypQp8Pf3BwBcu3YNYWFhSExMBACEhYWhU6dOeOqpp5CYmIiLFy/igw8+wI4dOzBx4kQxyrQK079OwmsbT0ldBhERUbMh2nle1qxZg7CwMERHR2Ps2LEYMmQIvvzyS+3j5eXlSEtL0/a42NnZ4ffff4eXlxfGjx+PXr164dtvv8U333yDsWPHilWmVfg5+arUJRARETUbRu02MoaHhwfWrl1r8PHg4OA6Z5Lt3LlzqzujLhERERmH1zYiIiIiq8LwQkRERFaF4cVKfJtwGWUVrfeoKyIioioML1ZiyabTmP3dEanLICIikhzDixXZk3ZD+3ftwc5EREStBcOLlSkpV+OzXecxKH4n0m8WobRCLXVJREREFiXaodIkjrDF27R/3/3+Hrg42CJlyWjYyGUSVkVERGQ57HmxADsb8YJFQUkFbhWWivb8REREzQ3DiwU8HhUs6vOX1jgKqaCkXNTXIiIikhrDiwW8FBuG72YOFO35Z6xOwuubTwMAPt99QbTXISIiag4YXizA3laOoZ0bvhq3qS7kFGL1wcvIKy5DVn6JaK9DRETUHDC8tCA3C8ukLoGIiEh0PNqoBVm88RQOXMoG2lTeV2t4LhgiImp52PPSgiRcuqVzv/87O7DlRBYA4MCFm5j/Ywryizmgl4ikt+/cDTyxOglK7uomEzC8tGAl5RrErT0GQRAwddVhrD96Df/3R6pOm/zicl4zicgEJeU8QWRTTPsqEbtSc/DKhpNSl0JWiOGlFRjwzp/avw9dvIXtp5UQBAE3CkrR+83tGPnBHumKI7JCqUoVwhZvw9JNp6QuxepZoudF3+VUjlzOxfhP9yP5Sq7or0/mx/DSCtQcyHvpZhFmf5eMDceu4feT1wEAV/++bZE6ytXs4WltBEEwuYfi4MWbeH3zaYv1cGw5kYUnVic1atfqJ3+eBwB8k3BF7LJaPLFH5t0qLEXkv3bind/O6Ex/YGUCTl7Lx6QvEjB55UHk3254vecWleG7hMvc/d4MMLxY0J4XRiDu7k5SlwEAmP/jcSy9c26YxtpxJhvzf0zBhZxCnM8uMGrehb+cQK/Xt9f5lZV85W+cupZf77x/nb+BeeuOmeULQxAEaGoMZC4pV+v8KsvKu90iBzprNILJyyUIAgpLK4yaJ1tVgpTMPMStPYawxdtwLU83IN8qLG0wlDzyn8NYffAyPt9zUWd6/u1ynXVYW3FZ42otLK3AhzvOIU1Z+b8ct/YYdqXm4L3tqbhyq6jeeQ29l7tSs/H65tMGd8WackHV4rIKg8G/qLQCBy/cRIURPww0GgFJl3PNfkLLhpatsLQC725Lxeks3c/72esqPLfumPZM4Sv3XsS721L1PUWjqTUCBEFAaYUaK/deRE5BKf7zV7rB9kmX/8bKvRcNPl7lyW+SsHjTacz/MUVn+qFLtxC/9ax2vfPCueJjeLGgYE8nPBYVJHUZBlWoNfjv/nSsP3oVdy3bhRd/Oo6B7/yJP89kAwBmfXsE649eQ8yHezHqo33IVpUYvDDkiat5uP/zAzhyubJLdl1SJm6XqzEofif+8dl+3C5TI6+4DJO+OIh7P90PjUZAcVkF3vj1NJIu63bjPvbfRGxMyaozXqfKqWv5jbpEgiAIePg/hzB2+V84dOkWfkm+irDF2zDzmyMAgJ1nszF42S7M+V9yg89VVFqB0R/txVtbzkCtEXC7TI11iRnIUVWHs4KSchy6dAsP/jsBRy7n4kyWqsHnNaTmxmlPWg5SMvMAVI5ZeuCLg/jukOEeAEEQMHb5X7j7/T3a56n6cm+MuWuPosfSP/QGVo1GwKFLt+psCCP/tRMTVxzAb3d6974/nKF9LFtVgoi3/8SQ/9vVqNevGSQybhWj9xvbMeU/h/S2HfbuboQv+QOL1p/Ax3+e0y7j30VldQLYu9tSsXzneYz5eJ/O9P8dysDw9/Zg68nrdX6N59+ufL+33/lMVClXa3D2ugpPrD6C1Qcv4/vEDJ3Hyyo0uFFQikHxO/Hmr7o9AEWlFdidlqM38JSUqzH8vT0Y9eFevetrxtdJeGTVYe2G9/eT1zHlywRk3/k/LCytqBNsvk/KwOSVCRjwzp84l12A5Cu5eOGn47hpwmVGDly4idiP9+G+zw9g2Hu7sXLvRcT/flZvre//kYbP91zEuOX7tfVV2ZSShbd/OwuNRsCyran4fM9FXLpR2OggWqW4rPJyKcPf242pqw6j35s76oSW/Nvl+CX5ap15i/QE9LIKDSZ8th9L7uwePJqRBwDYmZqj027Kl4fw772X8M3By1CVlCNk0e8Yt/wv7eOFpRWiBhpBEFCh1mDR+hOYuuoQNh/PwoZjdZexpq0nr+Nfv1e+58aE3+aCh0pbmLeLQuoS9BIEAd8duoK3tlR/sf505wP+5LdHkLJkVJ15Iv+1E64Otvh6xgCEejqjrZO99rGHvzyEojI1HliZgJ+fjtKZ78TVfGxMuYZ+gW21075PysBvJ67j4MVb+PrAZVxeNg6CIOjs8qravaXWCKjQaKCwtcHJq/kY/9l+AMDlZeOQU1CCz3dfxNTIQHT2cdHOq9EI+Lu4DIcuVQajKV9Wb/x23fki+nLfJQDA9jPZ2HryOrJVJdhy4jrevq8HlPklGN7FCzJZ5XWq1h+7hnPZhTiXXYgfkzJRcOeLL9DDEfteuhu9Xv8DqpLqL8MHViYAAF4ZG4Ye/m4Y3MkTAHCzsBRO9rZoY2+DGwWleGvLGUyNDEREUFvY2lT+tnh1w0n8dvI6tj8/DKXlGkz/OgkAMCrcBzvubESPXPkbjw0K0o5l8nZ10L52aYUGqXd6F67+fRteLgpEf7AXEcFtseKRfnXWa833rKCkAr+fVAIAvk24grcm9gBQ+f+SmJ6Lk9fy8fZvZxHu54rf/jlE+/7Ulne7DFduFSGonRO2n1beWXbd8xJtPXkdr2w4iedHdcHZ69VB7/idoHb5ZhFWH7wMAEhM1z9OISO3GADwfWImACDczxX9gz3Q760dlcvwxEB083OFl4tCGwAB/b01c9YcBQA8M6Ij+ga2xahwH3x78DKOXPm7TtvOr27Vub9082mk3yzC4nvDUVhSgbs/2IPcosrl/epAOpaMD9e2ffp/yfjr/E08NSwUi8Z2w8Zj1/BDUiZWPhaBm4WluFFQihsAVCUVcGtjp/M6iTV+HMSN7Ixn7tT85q9n8OaE7oh4u3K826R+HRDu74on7grGqxsqN8Ql5RqM/qg6uBWVVuDzqf0MrsOcghJ4ONpr/y9LytWYuuqwTptlWyt/YIzu7ouIoOrP942CUpy4mqe9H/mvnXWe/9rft6GusYEf+cFeAMDZN2OhsJUj6XIuurd3g7PCFltOZOGd387i86n90DewLdYcvoL//pWOSzerg66+3eH7z9/Eit0X6hyZaciu1Gwcv5qP41fz8eaEHnUeT76Sq/3+AID0W0X47s6uxNNZKpRVaHAuuwD3frof9/dtjw8f6qNtW67W4PLNInTydoZMJsOK3ReQrSrBG//orl0HGo2AQ+m30N3fDcVlFfB1dYBMJkNxWQVOXM3HgGAP2MhlmPnNEZy4mqf9TB24ULl8I7p463wv11T1//2fvy7Bw9Eeu18cAVcHO1zIKYSXiwIl5Wq4O9pBYWvTqPfK0hheLMzQF4PUQhb9Xu/jfd7coXe6qqQCk75IgIvCFovvDUdEcFt09HJGUVl1j0zVhrumRetPwrfGBrbqC7WKIAh4+ZcT+PHIVZ1pADD2k79wJbcIa56MxE9HdH9dPLE6CaeuqbD64GUcWzxK+8F94psk7Em7YXD5lu88r7NRqvpgA0Dsx5W/oF4Y3QUjw3wAAHtqfGEV1PjFVrXxrBlcavrX75Vf7k8OCUFnH2e8/MtJ2NnIcP6dsXht40n8cTobm49XHt7eoW0btLGzwfmcQgDAT0eu6gS+HbV+/QPA65tP45uEK/j04b4Y19MPq/ZfwvbT1e2u55fgWObfUKpK8NuJ61jxSOX0sgoNtp9RYtnWVMweForScg32nb+Bv87f1M5rW+MCo/sv3MRj/03U3j9zXYWQRb9j2f09MWVgYJ26/ncoA/87lIG/XrobH90ZL1Jb1Xu+ZJPu7szLt4qxaP0JbSCpKbeoDMt3nkdMNx+9u0Fmf6fbizbtq0S4ONji5OtjoKmxoQxf8ofemgBod1tdXjbOqN1nqw9e1oat+lS9x2sPZ2BeTBfM+yEFAPDpzvM6PbXnswvQydsZG45dQ0RQW9yqFewXb6z+DOXdLtP5f//l6FX8chT4LsFwPVtPKRGy6HesmtYffQLd4elc/UPr1LV83PvpfoT5usDeVo7hXbzw6S7DlyL5PjEDF3MK4e2qQFmFps560Cfxcq7e74qLNwpxLDMPizeeQu8Ad/xv5kDErT0GALjv84OY0Mcfm1KyGnx+AHj0v4cNPiYIlYFMVVKOolI19l+4iXPK+nePT/pCt97Scg3+LqpeL2VqjbZXbP2xa3Cwt0Fnb2esOZyBv4vKcKuoDE8OCUFbJ3u890caACDM1xUd2rbBsC5e+PFIJhaurz4a64m7QrBkfDie+q4y8L44pivm3t1JJ0DVVFyuhqy4DP9clwJPZ3tsSsnCqsf74+6u3jrLfauoDNtOKtE7wF2nJ7Lqu+nijUJ8tT8dTw/vCAAI8HCs932xBJnQwnbOqVQquLm5IT8/H66urlKXo1fwwt9Ee24NSpDZ5gEAQMDtnyGHQwNzmN/lZeOavIyezgqju7EHd2yHgxd1f1F9Nb0/Rob5iPqe1/bupF546ZcTRs3z6KBA/O9QRsMN61H7fXdR2OoEKwCYPSwUHdq20QaEAI826Orjgj/P6v/yq2nW0BC8Oq6yx2DWt0f0hicAiOnmgz/P6n9sRFcvnY3q5WXjcPJqPp5bd0znV3NjJL0ag/jfz2L9sWtGzVf1usb+T6S+FYuwxduMfi19/hndGZ29nbFk0yn8Xc9YLl9XByhVph2NM6yLF/adMxzYG+KisIWdrRzJr8VoQ7EUOno54eIN4/43xFD7f+aH2YPw0Jf6d19WmTkkBNvPKJGZa/xBEf0C3bW7qWp6dmQnbXD0dLbHzvkj0PvN7Xqf4+DCkfhy36U6IXp8b3/8elw38L33QC+oSip0et+ByuXu//YOnZ7SL6b2wz09/YxepoYYs/1meJFASw8vQzp5Yv+Fmw03tIBwP1f8PCeq3l/WLYU5QmNjX2fiigM6u11MNbGPPzY28ldzbTJZ5a9GU6THj22wt5EqzRwSgiu3ig0G0tZiXE8/7Rgua/HUsFDsO39TZzessd6a2EOnVw8A7G3kOPfOPU0trw6GF4YXANKFF5KGu6Md8ixwCOfyh/vin98fE/11iKj5urxsnNmf05jtN482ImohLBFcADC4EJHkGF6IiIjIqjC8EBERkVVheCEiIiKrwvBCREREVoXhhYiIiKwKw4sEhnfxgrx5nmiXiIio2WN4kcDqGQNwdHHdawURERFRwxheJCCTybQXNyMiIrImLg7SXxaRW1CJcK8RERFZo+aw/WJ4kUgzvbg0ERFRvWTNYAPG8CIRB1sbhHo6SV0GERGR1WF4kYhcLsPWeUOx5slIqUshIiJqtGbQ8cLwIiWFrQ3c2thJXQYREVGjNYPswvAitQqNIHUJREREjcYxL4QKtUbqEoiIiBpN+ujC8CK5NvY2UpdARETUaKUV0v/oZniRWEcvZ6lLICIiajR1MxjuwPAiMQc7G6S+FYt2TvZSl0JERNSgZjDkheGlOXCws4EdLxdARERWoBlkF/HCyzvvvIPBgwfD0dER7u7ujZpHEAQsWbIEfn5+aNOmDWJiYnD+/HmxSmxWlKoSqUsgIiKyCqKFl7KyMkyePBlz5sxp9Dzvvvsuli9fjpUrV+Lw4cNwcnLCmDFjUFLCDTsREVFz4OEs/TAH0cLLG2+8geeffx49e/ZsVHtBEPDxxx/jtddew4QJE9CrVy98++23yMrKwsaNG8Uqk4iIiIyQmXsbM1cnSVpDsxlokZ6eDqVSiZiYGO00Nzc3REZGIiEhQcLKLMORh0wTEZGVSL9VJOnrN5vwolQqAQA+Pj460318fLSP6VNaWgqVSqVzs0afPdJX6hKIiIgaRepBu0aFl4ULF0Imk9V7S01NFatWveLj4+Hm5qa9BQQEWPT1zUVhy54XIiKyDlJfIsDWmMYLFizA9OnT620TGhpqUiG+vr4AgOzsbPj5+WmnZ2dno0+fPgbnW7RoEebPn6+9r1KprDLADApth7u7emF32g2pSyEiIqqX1D0vRoUXLy8veHl5iVJISEgIfH19sXPnTm1YUalUOHz4cL1HLCkUCigUClFqsiQbuQxfzxiIsZ/8hTPXrXPXFxERtQ5Sn6hOtDEvGRkZSElJQUZGBtRqNVJSUpCSkoLCwkJtm7CwMGzYsAFAZRfUvHnz8Pbbb2Pz5s04efIkpk2bBn9/f0ycOFGsMpudcH9XqUsgIiJq1ozqeTHGkiVL8M0332jv9+1bOSB19+7dGDFiBAAgLS0N+fn52jYvvfQSioqKMHv2bOTl5WHIkCHYtm0bHBwcxCqz2VkwugsS03ORkVssdSlERER6ySTecSQTBEH6KyyZkUqlgpubG/Lz8+Hqar29GMELfzNpPg1KkNnmAQBAwO2fIUfrCX5ERGQZHb2csHPBCLM+pzHb72ZzqDQRERFZh4s3eJ4X0uOuTu2kLoGIiKhZYnhpppaO7y51CURERM0Sw0sz1cXHBWfeHCN1GURERM0Ow0sz5mgv2sFgREREJhsY4iHp6zO8EBERkVGkPsMuw0sz52DHVURERM1Liz3DLpnHtueG4cUxXaUug4iISEvqk9QxvDRzwZ5OmHt3J6nLICIi0mLPCxEREZERGF6IiIjIKOx5oUbxcLKXugQiIiIAwPnsQklfn+HFSuxeMAKb4+5CRFBbqUshIqJWLqegVNLXZ3ixEm6OdujVwR1rZ0Xiz/nDpC6HiIhIMgwvVkZha4NO3i5wa2MndSlERESSYHixUseXjkZ6/FipyyAiIrI4hhcrJpN6uDcREZEEGF6s3PpnBmPGXcFSl0FERGQxDC9Wrl9gWzxxV4j2/tk3YyWshoiISHy2UhdATRfg4Yj/TOsPd0cO4iUiopaP4aWFGBXuAwAoKiuSuBIiIiJxcbcRERERWRWGlxZs5WP90N69DeLu7oSYbj5IfCUaF//Fw6uJiMi6cbdRCza8szfGLgypM33mkBBsO6XEtbzbElRFRETUNOx5aYUW3xuO/S/fLXUZRKIa18sP658ZLHUZVMszIzpKXQK1AAwvrZRMJsO2eUPx41NR6BvojgHBbfHLnCj0CXCXurRm476+7fFr3BD8b2ak1KWY1cMDA6QuwSJs5TL4ujqI8txRoe3grKjbcf3syE6ivF5L8lJsmNQlSKrq4ApqGoaXVizM1xUDQzywfs5g/PhUFCKCPPDV9AF4ZkRHvD+5Nzp5O5vttV4Y3UX7t41chvbubTAvpjMeGxRkcJ57evg2+XV/fCrKpPmmDw7GRw/1Qc8ObhjS2bPJdUgtwKON9u/norvgg8m9zfr87o52WPFIv2Z1yYrgdk7wd2/TcEMTfD97EJIXx+DDB3XfxwWju2r/frB/h0Z9huLuNhx4wv1cTS+yAXteGKFzf+tzQ01+Lid7G2x9biiOLxlt9LyPDgpssI2h00BEhbYz+vX0cbK3QdKrMdr7vTu4YWgjPvdDOhn/3fCfaf0b3fbiv8Yippu30a/RWANDPLDwHtPCZERQWzNXYxyGF4JMJtNeasDDyR4vxYbhgYgO+HP+8EYHmF/mDIans0LvY4NCPRA3srP2/opH+uKvl+7GvJgueHNCdzw8MBBL7g3H1ueGoncHN227Lx6NQMKikVjViA/78SWj8ddLd8PBTo77+7UHALw6thsGhnjotHN1aHiY1/9mRtb7gQ5u51hn2uVl4/DV9PrrXHJvOI68FlNvGwD1vudvT+yB18Z1M/h40qsxOPn6aPwaN0Q7bWhnT/z+z+oNk0YQINSYZ8uzQxDbvTooHlw4Uvv3Y4OCcPqNMfXW6+2iwLHFozCulx9kMhleNtMv6wciOmDBqC6YGhmI8b394VKjp6O9extMiwrCnhdG4OGBdTd+M+4KxtPDdXdPdDMQBEaF+9TZUMlklT039VHY2uD+fh3w1PBQAJXvIwB8Nb0/Hh4YiDcn9MCTQ6rHnD09vKNOiASAL6b2wwtjuuLtiT0wuGP1hjjczxWrpvXH788N1dn1Vd/u3qB2jnhqWGid6f+bGYkx3ev+2g/2dMKg0MrPh7+bQ5N+rBxfOhrd/Fzh5miHY4tHIWXJqHrbVwWW+Pt74u2JPRt8/qOvjarzWQYqQ6QhG+fehcOvROOd+3rorN9AD0e8PbEHHo+q/uH0UmwYHOyqN4f/90AvfDNjICb28cdsPe9plf892XCv7Kpp/bHl2SGwt5Xj44f6GGz3vp4fFHIZIK9xGRh94XBsT1/8Mqfuj7Rt84bWCde1DQhuq/M5+WRKH7w4pms9c1SzaeDzITYO2KV6rXkyEr8ez4Krgx3+OK3EztQcve0igtqivbsDbhaWAgCOLh6F/NvlWHPoSp0Pv6O9LeR3/vFlMhni76/+8to49y7sv3ATA4Irv6j83NrAz83wr+c/5g2Dv7sDXBzs4OZoh9S37gEAfDC5t95rP61/5i7EfLgXAHBs8SjM+yEFe8/d0GlTX0/LMyM64sUxXVFYWgE7GzkOXbqF4HZOAIC7u3pj3exB6OjlDFVJOext5PjjtBId2joi6XIuHosKgp2NHGufjMQnO8/jcHqu3tf4c/5wnMlSYdX+S/j1eBbK1QL2vXg3ZLLKExKWlKvx9m9nAQCpb8VixtdJSLh0CwDg5VIZIHvWCIGRIR5wtK/+qHs42UMQquNLj/ZueHNid6hKyjE1Mgj+7m1wf7/2WH/0Gp4YEgJHexuD78cvc6LQyctF572eM6IjpgwIwLjlfyErvwSDQj3Q1tEeW08p9T7H7GGh+HLfJZ1pns72db7Mj2X8jfs+PwigciNe9ZpvTuiO+/u1xz+/P4br+SUY18sPS8d31843LSoIV24V4+vpA/DuH2lYufei9rH7+7XHhw/20d7XaATcKiqDk8IGhaUVGPjOToPLXmXRPd2wMDZMW8/IMB+MDKsMCz3aV6+HAcFt0cXHGfN/PA4XhS1O1giFjw4KwqODghC88DcAwLAuXoi5s3uhX2BbfPhgb/i6OaBDW0fE3d0Jn+2+gJFh3njirhCoBQEbj13D0vHhcHe0h7ujPbJVJVh98DIAoJufC4SaabWGTx/uh68PpOPhgYGws5FjckQHbD2lRGFpRZ22XX1c8OuzQ3A+pwAfbj+H0d198NaWs/h25kDY2lRv+Ns62QMA7G3lKKvQaKcvu78nenVwBwC8+Y8emD20IwL1/BAAKj+/Pydf1f5fy+UyaDT6F2Lrc0Px+Z6L8HFRYPuZbEyLCoKLg612F/jUyCCcupavbb/vpeoA+PSIjjiWkYfY7r4oLldrp8tlMsjlMnw8pS8A4JWx3XAt7zb+b2sqNh/P0v9m6uHWxk67Hs+9fY/eNp28nfGfaf0R4umEF346rvOYTCbDqHAfbD+TDQ8ne7g52uHJISFYtT+9ug1kiAjywJQBATiVlY9T11QAKnvWw3xd0SfAHdfybmPJptMoLK3AjYJS7bzDu1T26jweFYTD6bkY090XXx2ofu4+Ae4Y1sULyvzbGBnmgzK1BmlKFf699xIWjwtv9PsgBoYXqpePqwOeHFoZPh4cEKD9cu3m54qz11U6bZ8d2RlPfnsE9/VtDw8ne3g42eO1e6v/wV8c0xVpyoJ6u1plMhmGdvaqM333CyOQpixAW0c7nLiaD4WdHJ7OCnT1dTH4PPqnV//t4mALhW3jOh/3vDACu9Ny8PDAQMhkMrg4VHZjj+ha3aUrk8kw6E43dlWIqHrvYmvsAhvcyROH0nPrhJdHIgPxQEQHAEC4vys+fLCPzoa1ioOdDRJfjYaNTAYHOxu8+0AvvLbxVJ2QuPfFEfjr/E082D8ANnIZEl+JhloQ4GBXN4x4uzhg7azqX7EfTO6Nf93XEw52NjpBp7aIoLq/hoHKDdjW54bhXE4B+ge1hVojoPvSP1BaY2O28tF+6BfYFt6uDjiXXYA9adUhsioQ6jyno73275rr185GjgHBHvhlzmD8ejwLU2r1xLw5oYf272dHdsLJa3mI7eGHB/t3gL2N7vqXy2Xadedob4sfZg/C9fwSrNp/SbtR0MfQ/1vNX6caAZjYp/Kz0d3fTW/7xfeGY1PKNcyp1Wt0f78O2r+fH9UFd4d5o0d7VyhsK9fl8C7Vn5k5dwbEFpZWQCMIaOeswIiu3th+JhtyGTB7WEeMCq/8v/VyUeiMQXlvcm+8N7k3dqfl4GpuMRZvOg2gskfP3dEOdjZydPd3w3+nDwAAPDTA8C6flCWj8HPyVSy58xw114tcLtMJLomvRuNEZj4igtpCQGXALq3QaMMLAKgN/B9283PFpw9Xhoya3zc1GfoX9nNrA7+elT+OanYk6OtUaO/eBssf7guZDNiUkoX7+lb28PYNdMexjLw67VdN66/9EVbbqTfG4PClWxjS2VO7DoHKI0E3pWRpfwQCwKR+HeDpokDPO0H4tXvDsWhsN3R85XcAQNid78Blk3oh+crfmPTFQZ3XCvVyRqiXM3a/MAIl5Wocufw3Ono74Xp+CfoFVu76eaPGZ6SLd/V36sa5d9Utvrc/novuAvtGfneKheGFTPLtEwOx40w2XtlwUjstJtwHia9Gw8vA7qO59ezbb0iIpxNCPCs3aJFG7ueu+mWf/FqMzobbRi7T2X0yOtwHUwwMZg32dMIMz7qHnZtq5l0h2HIiC/f28se9vfzg4+oAtzaNv7yDt0v1QNQAD0d888TAOm2C2jkhqEYI8HbVnac+sjvBqOrvKqum9cfWU0r8cvRqgzW6Odppv7xtbWQ4vrSyyzu3qAwuDrbaAAgAnzzUF5uOX9Nu6PT9Ig/2dMLbE3vAw8m+zmMA4O/eBk8Nr/9IFieFLdY8aXhXQ21V/2sT+7bH+qNXMf/H4w3MoaurT/WGYGCIB+RymU7grW3mkBDMHFL//5mNXNao8QY1e64eGhAADyc79AloC1+3hgcx332nxlAvZ9jZyLWBzhiO9rbwbuR83i4OiAnXrat2gDDU89IYgzt5Yl1SZr1tZKh+QUNhFAD+b1IvTOzTHlF3dvN9P2sQLt4oRLlawPnsArz48wkAQL+gtnAzMFbHWWGL6G51d+UtvjccE/r44x+fHdBOk8tl2vVRxUYuw/bnh2Hn2RydC/M2tFvcwc5G27NsqEc7ups33p7YQ6fXsDapgwvA8EJGSnw1GkWlani5KPT+Oqm5UW0uXhnbDa+MrR4nsv/lu2FvI4dMJtP5RfalEQPpmsrN0Q67Foyw2OvVNii0Hd6a2AOdjRznENTOERP7+jcqvNRWFYb0DaJ1c7TDtKhgdPZ2wU9HMg12ST9azwBvsd3byx8bU7K040QaQy6X4eK/xqK0Qq2z687SbOQyxPbwM3q+u0wYkFqTWxv9QbMxxnT3xeJNp7S9mW9M6IFJXxyE+k6I6dG+8YOZx/fyg4OtvN4Ncs28YlNPeHGws8HdYd4696t60voEuMPFwQ6lFWqDIbsh+npG9eni44IuPro9z519XDB/VJdGh0Z9ZDKZpJ+zxmJ4IaN4uzgAdz4vTnoOFbUGHdpW/6qfFhWEP89m6wyWbC3qO9KrtkcHBSK3qAydvJ3RydsZn0zpU+eL0xyiOrbT/qJtbuxt5fhWTw9XQ2zkMkmDi5QGhXpgxl3BJg0Gbutkj1NvjNHu2usT4I7Ut2Irxy8dSMczRvTkymQyjO7e+KMXmzIYNbaJR0l29nbG41FBOj2lxvhndOeGG7UAMqG+HdpWSKVSwc3NDfn5+XB1Fe8ww+aqqKwIzvGVXxSFiwrhZF937IC5VKg1iFt7DBFBbTGrnhH5zd3Vv4vh6+qgM+iQiFqXcrUGnV/dCqDyqB5Du3xIPMZsv1vnzwEyC1sbOVY+FiF1GU1WsyeGiFonOxs5vps5EOVqDYOLFWB4ISIiAvQe6UjNE/vJiYiIyKowvBAREZFVYXghIiIiq8LwQkRERFaF4YWIiIisCsMLERERWRWGFyIiIrIqDC9ERERkVRheiIiIyKqIFl7eeecdDB48GI6OjnB3d2+wfXl5OV5++WX07NkTTk5O8Pf3x7Rp05CVlSVWiURERGSFRAsvZWVlmDx5MubMmdOo9sXFxTh69CgWL16Mo0ePYv369UhLS8M//vEPsUokIiIiKyTatY3eeOMNAMDq1asb1d7NzQ07duzQmfbZZ59h4MCByMjIQGBgoLlLJCIiIivUrC/MmJ+fD5lMVu9up9LSUpSWlmrvq1QqC1RGREREUmm24aWkpAQvv/wyHn74Ybi6uhpsFx8fr+3lqam1hpiisiKgpPJvlUoFtb1a2oKIiIgaoWq7LQhCw40FI7z88ssCgHpvZ8+e1Znn66+/Ftzc3Ix5GaGsrEwYP3680LdvXyE/P7/etiUlJUJ+fr72dubMmQZr5I033njjjTfemuctMzOzwZxgVM/LggULMH369HrbhIaGGvOUdZSXl+PBBx/ElStXsGvXrnp7XQBAoVBAoVBo7zs7OyMzMxMuLi6QyWRNqqU2lUqFgIAAZGZmNliXNWrpywe0/GXk8lm/lr6MXD7rJ9YyCoKAgoIC+Pv7N9jWqPDi5eUFLy8vkwtrSFVwOX/+PHbv3o127doZ/RxyuRwdOnQQobpqrq6uLfafEmj5ywe0/GXk8lm/lr6MXD7rJ8Yyurm5NaqdaIdKZ2RkICUlBRkZGVCr1UhJSUFKSgoKCwu1bcLCwrBhwwYAlcHlgQcewJEjR7BmzRqo1WoolUoolUqUlZWJVSYRERFZGdEG7C5ZsgTffPON9n7fvn0BALt378aIESMAAGlpacjPzwcAXLt2DZs3bwYA9OnTR+e5as5DRERErZto4WX16tUNnuNFqDGiODg4uHEjjCWkUCiwdOlSnTE2LUlLXz6g5S8jl8/6tfRl5PJZv+awjDKhuScGIiIiohp4YUYiIiKyKgwvREREZFUYXoiIiMiqMLwQERGRVWF4aaQVK1YgODgYDg4OiIyMRGJiotQlmc3rr78OmUymcwsLC5O6LJPt27cP48ePh7+/P2QyGTZu3KjzuCAIWLJkCfz8/NCmTRvExMTg/Pnz0hRrooaWcfr06XXWaWxsrDTFmiA+Ph4DBgyAi4sLvL29MXHiRKSlpem0KSkpwdy5c9GuXTs4Oztj0qRJyM7Olqhi4zRm+UaMGFFnHT799NMSVWycL774Ar169dKexCwqKgpbt27VPm7N665KQ8tozetPn2XLlkEmk2HevHnaaVKuR4aXRvjhhx8wf/58LF26FEePHkXv3r0xZswY5OTkSF2a2XTv3h3Xr1/X3vbv3y91SSYrKipC7969sWLFCr2Pv/vuu1i+fDlWrlyJw4cPw8nJCWPGjEFJSYmFKzVdQ8sIALGxsTrr9Pvvv7dghU2zd+9ezJ07F4cOHcKOHTtQXl6O0aNHo6ioSNvm+eefx6+//oqffvoJe/fuRVZWFu6//34Jq268xiwfAMyaNUtnHb777rsSVWycDh06YNmyZUhOTsaRI0cwcuRITJgwAadPnwZg3euuSkPLCFjv+qstKSkJ//73v9GrVy+d6ZKuR6OumNhKDRw4UJg7d672vlqtFvz9/YX4+HgJqzKfpUuXCr1795a6DFEAEDZs2KC9r9FoBF9fX+G9997TTsvLyxMUCoXw/fffS1Bh09VeRkEQhMcff1yYMGGCJPWIIScnRwAg7N27VxCEynVmZ2cn/PTTT9o2Z8+eFQAICQkJUpVpstrLJwiCMHz4cOG5556Trigza9u2rbBq1aoWt+5qqlpGQWg566+goEDo3LmzsGPHDp1lkno9suelAWVlZUhOTkZMTIx2mlwuR0xMDBISEiSszLzOnz8Pf39/hIaGYurUqcjIyJC6JFGkp6dDqVTqrE83NzdERka2qPUJAHv27IG3tze6du2KOXPm4NatW1KXZLKqM3F7eHgAAJKTk1FeXq6zHsPCwhAYGGiV67H28lVZs2YNPD090aNHDyxatAjFxcVSlNckarUa69atQ1FREaKiolrcugPqLmOVlrD+5s6di3HjxumsL0D6z6BoZ9htKW7evAm1Wg0fHx+d6T4+PkhNTZWoKvOKjIzE6tWr0bVrV1y/fh1vvPEGhg4dilOnTsHFxUXq8sxKqVQCgN71WfVYSxAbG4v7778fISEhuHjxIl555RXcc889SEhIgI2NjdTlGUWj0WDevHm466670KNHDwCV69He3h7u7u46ba1xPepbPgB45JFHEBQUBH9/f5w4cQIvv/wy0tLSsH79egmrbbyTJ08iKioKJSUlcHZ2xoYNGxAeHo6UlJQWs+4MLSNg/esPANatW4ejR48iKSmpzmNSfwYZXgj33HOP9u9evXohMjISQUFB+PHHHzFz5kwJKyNTTZkyRft3z5490atXL3Ts2BF79uxBdHS0hJUZb+7cuTh16pRVj8Oqj6Hlmz17tvbvnj17ws/PD9HR0bh48SI6duxo6TKN1rVrV6SkpCA/Px8///wzHn/8cezdu1fqsszK0DKGh4db/frLzMzEc889hx07dsDBwUHqcurgbqMGeHp6wsbGps4I6uzsbPj6+kpUlbjc3d3RpUsXXLhwQepSzK5qnbWm9QkAoaGh8PT0tLp1GhcXhy1btmD37t3o0KGDdrqvry/KysqQl5en097a1qOh5dMnMjISAKxmHdrb26NTp06IiIhAfHw8evfujU8++aTFrDvA8DLqY23rLzk5GTk5OejXrx9sbW1ha2uLvXv3Yvny5bC1tYWPj4+k65HhpQH29vaIiIjAzp07tdM0Gg127typs2+zJSksLMTFixfh5+cndSlmFxISAl9fX531qVKpcPjw4Ra7PgHg6tWruHXrltWsU0EQEBcXhw0bNmDXrl0ICQnReTwiIgJ2dnY66zEtLQ0ZGRlWsR4bWj59UlJSAMBq1mFtGo0GpaWlVr/u6lO1jPpY2/qLjo7GyZMnkZKSor31798fU6dO1f4t6XoUfUhwC7Bu3TpBoVAIq1evFs6cOSPMnj1bcHd3F5RKpdSlmcWCBQuEPXv2COnp6cKBAweEmJgYwdPTU8jJyZG6NJMUFBQIx44dE44dOyYAED788EPh2LFjwpUrVwRBEIRly5YJ7u7uwqZNm4QTJ04IEyZMEEJCQoTbt29LXHnj1beMBQUFwgsvvCAkJCQI6enpwp9//in069dP6Ny5s1BSUiJ16Y0yZ84cwc3NTdizZ49w/fp17a24uFjb5umnnxYCAwOFXbt2CUeOHBGioqKEqKgoCatuvIaW78KFC8Kbb74pHDlyREhPTxc2bdokhIaGCsOGDZO48sZZuHChsHfvXiE9PV04ceKEsHDhQkEmkwnbt28XBMG6112V+pbR2tefIbWPoJJyPTK8NNKnn34qBAYGCvb29sLAgQOFQ4cOSV2S2Tz00EOCn5+fYG9vL7Rv31546KGHhAsXLkhdlsl2794tAKhze/zxxwVBqDxcevHixYKPj4+gUCiE6OhoIS0tTdqijVTfMhYXFwujR48WvLy8BDs7OyEoKEiYNWuWVYVtfcsGQPj666+1bW7fvi0888wzQtu2bQVHR0fhvvvuE65fvy5d0UZoaPkyMjKEYcOGCR4eHoJCoRA6deokvPjii0J+fr60hTfSE088IQQFBQn29vaCl5eXEB0drQ0ugmDd665Kfcto7evPkNrhRcr1KBMEQRC/f4eIiIjIPDjmhYiIiKwKwwsRERFZFYYXIiIisioML0RERGRVGF6IiIjIqjC8EBERkVVheCEiIiKrwvBCREREVoXhhYiIiKwKwwsRERFZFYYXIiIisioML0RERGRV/h8ZD+BeU7f8cwAAAABJRU5ErkJggg==\n"
many's avatar
many committed
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
many's avatar
many committed
   "source": [
    "plt.figure()\n",
    "plt.plot(timesteps, values)\n",
    "plt.axvline(stationary_timestep, color=\"green\")\n",
    "plt.show()\n"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}