Check_SEI_Models_Temp.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pybamm\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "param = pybamm.ParameterValues(\"Chen2020\")\n",
    "# param = pybamm.ParameterValues(\"ORegan2022\")\n",
    "# param"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [],
   "source": [
    "# param['Initial inner SEI thickness [m]'] = 0.5e-12\n",
    "# param['Initial outer SEI thickness [m]'] = 0.5e-12\n",
    "# param['Inner SEI open-circuit potential [V]'] = 0.3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "cycle_no = 1\n",
    "exp_CCCV = pybamm.Experiment(\n",
    "    [(\"Charge at C/10 until 4.2 V\",\n",
    "      \"Rest for 10 hours\",\n",
    "      \"Discharge at C/2 until 2.5 V\",\n",
    "      )] * cycle_no\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [],
   "source": [
    "from os import name\n",
    "\n",
    "\n",
    "sei_models = [\"interstitial-diffusion limited\"]\n",
    "\n",
    "# sei_models = [\"reaction limited\", \"solvent-diffusion limited\", \"interstitial-diffusion limited\"]\n",
    "# sei_models = [\"ec reaction limited\", \"reaction limited\"]\n",
    "# Diffferent SEI models can be used, by setting the index of sei_models\n",
    "\n",
    "\n",
    "name_models = [\"SPM\"]\n",
    "models_SPM = [pybamm.lithium_ion.SPM(options={\n",
    "    \"SEI\": sei_models[index],\n",
    "    \"SEI porosity change\": \"true\",\n",
    "},\n",
    "    name=\"SPM \" + sei_models[index]) for index in range(len(sei_models))\n",
    "]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [],
   "source": [
    "# models_SPM[0].variable_names()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [],
   "source": [
    "# please note that the size of sims will be equal to the number of models times the number of experiments\n",
    "solver = pybamm.CasadiSolver(mode=\"safe\")\n",
    "sims_SPM = [pybamm.Simulation(model, parameter_values=param,\n",
    "                              experiment=exp_CCCV, solver=solver) for model in models_SPM]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Combine all solutions for the different SEI modesl in DFN and SPM and models in the list 'sols' Odd =DFN and Evens = SPM (unifrm diffusion in the particle)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [],
   "source": [
    "sols_SPM = [sim.solve(initial_soc=0) for sim in sims_SPM]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Vaiables to be postprocessed, you can add as many as output variables to the following list\n",
    "# you can find the list variables to be be output by typing model name dot varaible_name() e.g. here, models[0].variable_names()\n",
    "out_vars = [  # 'X-averaged negative particle surface concentration [mol.m-3]',\n",
    "    # \"Voltage [V]\",\n",
    "    'X-averaged SEI thickness [m]',\n",
    "    # 'Total capacity lost to side reactions [A.h]',\n",
    "    # 'Loss of lithium to SEI [mol]',\n",
    "    # 'X-averaged inner SEI thickness [m]',\n",
    "    # 'X-averaged outer SEI thickness [m]',\n",
    "    # 'X-averaged negative electrode porosity',\n",
    "]\n",
    "\n",
    "# Label for the output variables y-axis in the plots\n",
    "yaxis_labels = [\"SEI_thickness\",\n",
    "                'Total capacity lost to side reactions [A.h]', 'Total lithium in particles', 'X-averaged inner SEI thickness [m]',\n",
    "                'X-averaged outer SEI thickness [m]',]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "# from matplotlib import lines\n",
    "\n",
    "# for sol in sols_SPM:  # loop through the solutions of the DFN models\n",
    "#     LS_cycles = [sol.cycles[k].steps[2]['X-averaged SEI thickness [m]'].entries[-1]\n",
    "#                  for k in range(cycle_no)]\n",
    "#     plt.scatter(np.arange(1, cycle_no+1), LS_cycles,\n",
    "#                 label=models_SPM[sols_SPM.index(sol)].name)\n",
    "# plt.xlabel(\"Cycle number\")\n",
    "# plt.ylabel(\"SEI thickness [m]\")\n",
    "# plt.legend()\n",
    "# plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAHACAYAAABeV0mSAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABUAElEQVR4nO3deVhUZcMG8PvMsO8gu6wuuLGIqIj7Qi6Zr0uZpblnn1tlaiUtasv7kpampUVlilmWWi5luYXiigsooqnkgoLKIioMi2wz5/uDnERABxg4w8z9u665dM6cmbkZBubmnOc8RxBFUQQRERGRnpBJHYCIiIhIm1huiIiISK+w3BAREZFeYbkhIiIivcJyQ0RERHqF5YaIiIj0CssNERER6RWWGyIiItIrLDdERESkV1huiIiISK8YdLk5cOAAhgwZAnd3dwiCgK1bt0r+fKIoYv78+XBzc4O5uTnCw8Nx8eLFes1FRESkTwy63BQUFCAoKAgrV67UmedbvHgxPvvsM0RFReHYsWOwtLTEgAEDUFRU1CAZiYiIGjuBJ84sJwgCtmzZgmHDhqmXFRcX4+2338aPP/6InJwc+Pv7Y9GiRejdu3e9PJ8oinB3d8ecOXMwd+5cAEBubi5cXFwQHR2N5557rs7PS0REpO8MesvN48ycORNxcXH46aefkJSUhJEjR2LgwIH1tpsoJSUFGRkZCA8PVy+ztbVFaGgo4uLi6uU5iYiI9A3LTTVSU1OxZs0abNq0CT169EDz5s0xd+5cdO/eHWvWrKmX58zIyAAAuLi4VFju4uKivo2IiIgejeWmGmfOnIFSqYSfnx+srKzUl/379+Py5csAgAsXLkAQhEde5s2bJ/FXQkREZFiMpA6gq/Lz8yGXy5GQkAC5XF7hNisrKwBAs2bNcP78+Uc+TpMmTTR+TldXVwBAZmYm3Nzc1MszMzPRvn17jR+HiIjIkLHcVCM4OBhKpRJZWVno0aNHleuYmJigdevWWntOX19fuLq6IiYmRl1mFAoFjh07hmnTpmnteYiIiPSZQZeb/Px8XLp0SX09JSUFiYmJcHBwgJ+fH8aMGYNx48ZhyZIlCA4Oxq1btxATE4PAwEAMHjxYq8/n5eUFQRAwa9YsfPjhh2jZsiV8fX3x7rvvwt3dvcJRVURERFQ9gz4UPDY2Fn369Km0fPz48YiOjkZpaSk+/PBDfPfdd7hx4wYcHR3RpUsXvPfeewgICND68wHlh4MvWLAAX3/9NXJyctC9e3d88cUX8PPzq/HzERERGSKDLjdERESkf3i0FBEREekVlhsiIiLSKwY3oFilUuHmzZuwtraGIAhSxyEiIiINiKKIvLw8uLu7QyZ79LYZgys3N2/ehKenp9QxiIiIqBbS0tLg4eHxyHUMrtxYW1sDKH9xbGxsJE5DREREmlAoFPD09FR/jj+KwZWb+7uibGxsWG6IiIgaGU2GlHBAMREREekVlhsiIiLSKyw3REREpFcMbsyNppRKJUpLS6WOQUSNgLGxMeRyudQxiOgfLDcPEUURGRkZyMnJkToKETUidnZ2cHV15fxZRDqA5eYh94uNs7MzLCws+IuKiB5JFEUUFhYiKysLAODm5iZxIiJiuXmAUqlUF5smTZpIHYeIGglzc3MAQFZWFpydnbmLikhiHFD8gPtjbCwsLCROQkSNzf3fGxyrRyQ9lpsqcFcUEdUUf28Q6Q6WGyIiItIrLDekNbGxsRAEoVEdaaZpZh8fHyxbtkzjx42OjoadnV2N80yYMAHDhg1TX+/duzdmzZqlvl5YWIinn34aNjY26txVLasrQRCwdevWOj/OoyxcuBDt27dXX3/4a6+N+nwP1vZ7SkQNj+VGT9y6dQvTpk2Dl5cXTE1N4erqigEDBuDw4cPqdXx8fCAIAgRBgKWlJTp06IBNmzapb1+4cCEEQcDAgQMrPf7HH38MQRDQu3fvajN07doV6enpsLW11Ti3Nj7QNPVwUQAqZ67uA+zEiRN46aWXGiBlRZs3b8YHH3ygvr527VocPHgQR44cUeeualldpaenY9CgQXV+nJpYvnw5oqOj6/QYmn4/iUi/8WgpPfH000+jpKQEa9euRbNmzZCZmYmYmBjcvn27wnrvv/8+pkyZAoVCgSVLlmDUqFFo2rQpunbtCqD8MNZ9+/bh+vXrFU4pv3r1anh5eT0yg4mJCVxdXbX/xWmgpKQEJiYmNb6fppmdnJxqE6vOHBwcKly/fPky2rRpA39//0cuqyspvo/aKGVSvgeJ6pNKJaJUpYJSJaJMJUIUAYiACBEymQAbM2P1ujmFJep1RIjlC0VABCCXCXC0MlWvm5VXhFKlCFEUKz2nXCbAzdb833UVRShRqqrMJxMEuNuZV3mbFFhu9EBOTg4OHjyI2NhY9OrVCwDg7e2Nzp07V1rX2toarq6ucHV1xcqVK/H999/jt99+U5cbZ2dnhISEYO3atXj77bcBAEeOHEF2djZGjhyJc+fOVZsjNjYWffr0wd27d2FnZ4fo6GjMmjULGzZswKxZs5CWlobu3btjzZo1cHNzw8KFC7F27VoA/w7G3LdvH3r37o20tDTMmTMHu3fvhkwmQ48ePbB8+XL4+PgAKN/ik5OTg06dOmHlypUwNTVFSkoKvvjiC3z66adIS0uDra0tevTogZ9//hkTJkzA/v37sX//fixfvhwAkJKSgqtXr6ozJyYmYuLEiRXyLFiwAAsXLoSPjw9mzZql3vKzdOlSrFmzBleuXIGDgwOGDBmCxYsXw8rKSuPvm1KpxOuvv47Vq1dDLpdj8uTJlX7B9O7dG+3bt8eyZcvQu3dv7N+/X53v/vf64WX3d81s2bKlwlYxOzs7LFu2DBMmTEBJSQlmz56NX375BXfv3oWLiwumTp2KiIgI9WM9eP8zZ87g1VdfRVxcHCwsLPD0009j6dKl6q/3/veje/fuWLJkCUpKSvDcc89h2bJlMDb+95fuo9x/jPu7w3r37o2AgADI5XKsXbsWJiYm+PDDDzF69GjMnDkTP//8M1xcXPD555+rtzI9+B581PezuLgYb7/9Nn788Ufk5OTA398fixYtqrBlMjo6GvPnz0d2djYGDBiA7t27a/R1kGETRRF5xWXILSxFTmEpcu6VlP9bWAIna1MM9P93HqTpPyRAca8MhSVlKCxR4l6pEoUlSpQpVQj2ssfqCZ3U63b875+4U1BS5XP6N7XB9pd7qK8PWXEIaXfuVblucydLxMzprb7+wqpj+Dszv8p13W3NcCSin/r6S+sSkJiWU+W6dhbGSJzfv8rbpMByo6HCkrJqb5MJAsyM5Vpd18JE82+NlZUVrKyssHXrVnTp0gWmpqaPvxMAIyMjGBsbo6Sk4g/MpEmT8MYbb6jLzerVqzFmzBiN8zyosLAQn3zyCdatWweZTIYXXngBc+fOxQ8//IC5c+fi/PnzUCgUWLNmDYDyLRWlpaUYMGAAwsLCcPDgQRgZGeHDDz/EwIEDkZSUpN5CExMTAxsbG+zZswcAEB8fj1deeQXr1q1D165dcefOHRw8eBBA+S6Pv//+G/7+/nj//fcBlG+NuXr1qjpr165dsWzZMsyfPx/Jycnq17YqMpkMn332GXx9fXHlyhVMnz4db7zxBr744guNX5slS5YgOjoaq1evRps2bbBkyRJs2bIFffv2rXL9zZs3Y968eTh79iw2b96sfh2qWvY4n332GX799Vds3LgRXl5eSEtLQ1paWpXrFhQUqL8fJ06cQFZWFl588UXMnDmzwm6kffv2qbf8Xbp0CaNGjUL79u0xZcoUjV+Th61duxZvvPEGjh8/jg0bNmDatGnYsmULhg8fjrfeeguffvopxo4di9TU1EpTODzq+zlz5kycO3cOP/30E9zd3bFlyxYMHDgQZ86cQcuWLXHs2DFMnjwZkZGRGDZsGHbu3IkFCxbU+uugxkulEnEhIw93Ckpwp7AEOYX3y0p5YfFztcbUXs0BAEqViNbv7kCpsvJWEADo0dKxQrk5+Hc28oqr/gxQ3Ks4pYBcVvuj8QQBuH9v2UNH9ZkYyWBiJENVj25iJKt03dSo6tEs1S2XCsuNhtrO31XtbX1aOWHNxH+3koR88CfulSqrXDfU1wEb/i9Mfb37on1VtvGrHw3WOJuRkRGio6MxZcoUREVFoUOHDujVqxeee+45BAYGVnmfkpISLFmyBLm5uZU+TJ966ilMnToVBw4cQEhICDZu3IhDhw5h9erVGme6r7S0FFFRUWjevPyHf+bMmepyYWVlBXNzcxQXF1fYlfD9999DpVJh1apV6r+416xZAzs7O8TGxqJ///K/DiwtLbFq1Sr1B/rmzZthaWmJp556CtbW1vD29kZwcDCA8l0eJiYmsLCwqHa3hYmJCWxtbSEIwmN3bTw4dsfHxwcffvghpk6dWqNys2zZMkRERGDEiBEAgKioKOzaVf37zMHBARYWFpV2vVS17HFSU1PRsmVLdO/eHYIgwNvbu9p1169fj6KiInz33XewtLQEAKxYsQJDhgzBokWL4OLiAgCwt7fHihUrIJfL0bp1awwePBgxMTF1KjdBQUF45513AAARERH46KOP4OjoqH7M+fPn48svv0RSUhK6dOlS4b7VfT9TU1OxZs0apKamwt3dHQAwd+5c7Ny5E2vWrMH//vc/LF++HAMHDsQbb7wBAPDz88ORI0ewc+fOWn8tpDvKlCrsS76FOwXFuFNQqv73bmEJbheUIMTLHvOHtAUAqEQRT352sNrH6tHSUV1u5DIB5sZylCrLYGYsg525CewsjMsv5ibwb2pT4b4L/9MOMhlgbmwECxM5LE3lMDOWw9RIBvOH/sDd81pPCIIAI5kAuUxQl5QHi8t9B17vo/HUBA9u8XmcjQ98duk6lhs98fTTT2Pw4ME4ePAgjh49ih07dmDx4sVYtWoVJkyYoF7vzTffxDvvvIOioiJYWVnho48+wuDBFYuUsbExXnjhBfVuFz8/v2pL0uNYWFioiw1QPqbn/jT11Tl9+jQuXboEa2vrCsuLiopw+fJl9fWAgIAKWyqeeOIJeHt7o1mzZhg4cCAGDhyI4cOH18ukjH/++SciIyNx4cIFKBQKlJWVoaioCIWFhZWeLzU1FW3btlVff+uttzBjxgykp6cjNDRUvdzIyAgdO3asct+3tk2YMAFPPPEEWrVqhYEDB+Kpp55Sl8aHnT9/HkFBQepiAwDdunWDSqVCcnKyuty0a9euwsy8bm5uOHPmTJ1yPvi+k8vlaNKkCQICAtTL7j/3495TDzpz5gyUSiX8/PwqLC8uLlbPTH7+/HkMHz68wu1hYWEsNzqsqFSJP89nIktRjMy8ItxSFCMrrxi38opxu6AEPVo64tNR7QGU76Z8aV08qvtRszb996PRSC5DcydLGMlksLc0hr2FCews/ikt5sbwdbSscN+YOb1hbWZUYQt9dZ4O8XjsOvfZWWg+ppBzLrHcaOzc+wOqve3hzXwJ74ZrvO6hN/vULdgDzMzM8MQTT+CJJ57Au+++ixdffBELFiyoUG5ef/11TJgwAVZWVnBxcan2h2DSpEkIDQ3F2bNnMWnSpFpneni8hSAIj/3wzs/PR0hICH744YdKtz04sPfBD1ugfDzRyZMnERsbi927d2P+/PlYuHAhTpw4odUjZq5evYqnnnoK06ZNw3//+184ODjg0KFDmDx5MkpKSiqVG3d3dyQmJqqvPzxIuD5U9To/OHNuhw4dkJKSgh07duDPP//Es88+i/DwcPz888+1fs6qvtcqVdWDD+vymA8uu//+rcnz5OfnQy6XIyEhodJpEmoyZoqkcbegBPHX7uJU6l242ZljbJfyrY6CALz846lqC0tWXpH6/3KZgG7NHWEsF2BvaYImliZwsDSFg6UxHCxN0fShgbEPjlF5HCdrzYYFUP1iudFQTcbA1Ne6NdW2bdtKc5U4OjqiRYsWj71vu3bt0K5dOyQlJWH06NH1lLB814FSWXEXXocOHbBhwwY4OzvDxsammntWzcjICOHh4QgPD8eCBQtgZ2eHvXv3YsSIEVU+lyZ5HpaQkACVSoUlS5ZAJivfz7xx48ZHZqrqNXdzc8OxY8fQs2dPAEBZWRkSEhLQoUOHx32Zj+Xk5IT09HT19YsXL6KwsLDCOjY2Nhg1ahRGjRqFZ555BgMHDsSdO3cqla82bdogOjoaBQUF6kJ5+PBhyGQytGrVqs5Z61NV38/g4GAolUpkZWWhR4+qN8m3adMGx44dq7Ds6NGj9ZaTHu1SVj72nMtEzPlMnEy9C9U/BSbYy05dbkyN5Bjk7woBApxtTOFsbQZna1M4WZuiiZUJnK3NKjzm9y+GPvw0pEdYbvTA7du3MXLkSEyaNAmBgYGwtrZGfHw8Fi9ejKFDh9b6cffu3YvS0tJ6nSfEx8cHu3btQnJyMpo0aQJbW1uMGTMGH3/8MYYOHYr3338fHh4euHbtGjZv3ow33nijwiHqD9q+fTuuXLmCnj17wt7eHn/88QdUKpX6A9jHxwfHjh3D1atXYWVlVeUWFB8fH+Tn5yMmJgZBQUGwsLCotCWmRYsWKC0txeeff44hQ4bg8OHDiIqKqvHX/uqrr+Kjjz5Cy5Yt0bp1ayxdulRrk8/17dsXK1asQFhYGJRKJd58880KWzyWLl0KNzc3BAcHQyaTYdOmTXB1da3yez1mzBgsWLAA48ePx8KFC3Hr1i28/PLLGDt2rHq3kK6q6vvp5+eHMWPGYNy4cViyZAmCg4Nx69YtxMTEIDAwEIMHD8Yrr7yCbt264ZNPPsHQoUOxa9cu7pKSgCiKGPXVURy/eqfC8hbOVgjxskdn34o/w1+MCWnIeKTDdGt4M9WKlZUVQkND8emnn6Jnz57w9/fHu+++iylTpmDFihW1flxLS8t6nwBtypQpaNWqFTp27AgnJyccPnwYFhYWOHDgALy8vDBixAi0adMGkydPRlFR0SO35NjZ2WHz5s3o27cv2rRpg6ioKPz4449o164dgPJBo3K5HG3btoWTkxNSU1MrPUbXrl0xdepUjBo1Ck5OTli8eHGldYKCgrB06VIsWrQI/v7++OGHHxAZGVnjr33OnDkYO3Ysxo8fj7CwMFhbW1ca51FbS5YsgaenJ3r06IHRo0dj7ty5FUqatbU1Fi9ejI4dO6JTp064evUq/vjjD/WWqAdZWFhg165duHPnDjp16oRnnnkG/fr1q9N7q6FU9/1cs2YNxo0bhzlz5qBVq1YYNmwYTpw4oZ7LqUuXLvjmm2+wfPlyBAUFYffu3eqBzVR/ikqV2HshU31dEAR42JvDWC6gl58TPhjaDofn9cWfs3th0TOBNRqzQoZFEBti9KIOUSgUsLW1RW5ubqUPyqKiIqSkpMDX1xdmZmbVPAIRUWX8/VE7oigi/tpdbD55HduT0pFXVIY/XumBtu7lv59v5NyDubEcDpY1n6ST9MujPr8fxt1SRETU4NJz72HzyRv4OeE6UrIL1Mub2pkjM68IbWGjvk5UUyw3RETUoBKu3cHIqDj1wGBLEzmeDHDD0yEe6OzjAFkdJqwjAlhuiIioHomiiDM3cpGlKEZ42/IB6IEedmhiZYpmjpYY2dETg/xdYWnKjyPSHr6biIhI67Lzi7H11A1sir+O5Mw8uNmaoU9rZ8hlAozlMsTM6VXhZI9E2sRyUwUDG2NNRFrA3xtAqVKF2ORb2Bifhn0XslD2z34nEyMZOvk4IK+oVD3TLosN1SeWmwfcnweksLAQ5uYcxEZEmrs/SaKmZ0HXRx/tuIBvD6Worwd52OKZjp74T6A7bC0M93Whhsdy8wC5XA47Ozv1eWosLCx4jg4ieiRRFFFYWIisrCzY2dlVOqWDvsrKK8KviTcR6tsEAR62AICnAt2wLfEGhgc3xciOnvBzsX7MoxDVD5abh9w/e3BNTsRHRGRnZ1ejM7M3RvdKlNh9LgObT97AoUvZUKpEPN/ZE5Ee5Sc4be9ph7iIfjCWc35YkhbLzUMEQYCbmxucnZ0rnGiQiKg6xsbGervFRqUSEXflNjafvIGdZ9NRUPLvubqCvewQ4v3vKRAEQYCxnFu7SXosN9WQy+V6+8uKiKgm3vg5CTdy7gEAvBwsMCy4KYYHN4Wvo6XEyYiqxnJDREQAgOSMPGxPuokDF7Ox6f/CYGIkg0wmYFyYN1LvFGJEh6bo4GXPsYik81huiIgM2OVb+dh+Oh3bk27iYla+evmhS7fQt3X5pHv/16u5VPGIaoXlhojIAB29chvv/3YO59IV6mUmchl6+jniqUB3dPZtImE6orphuSEiMkAmRjKcz1DASCage8vyQvNEWxfYmnM+Gmr8WG6IiAxESnaBehBwBy97fPVCCDr5OMDe0kTiZETaxckIiIj0XJaiCFPXJWDAsgO4fOvfcTX927my2JBe4pYbIiI9JYoifk64jg+2n4OiqAxGMgEJ1+6iuZOV1NGI6hXLDRGRHrp+txBvbTmLA3/fAgD4N7XB4qeD0NbdRuJkRPWP5YaISM+sP5aKD38/h8ISJUyMZHgt3A9TevjCiKdFIAPBckNEpGfuFpagsESJTj72+OjpQO6GIoPDckNE1MiVKVXIyiuGu505AOClns3gbmeGoUFNIZNxNmEyPNxGSUTUiF3MzMOIL49g/OrjKC4rP6mlsVyG4cEeLDZksLjlhoioESpTqvDNwRR8uudvlChVsDEzwsXMfPg3tZU6GpHkWG6IiBqZS1l5mLMpCafTcgAAfVs7I3JEAFxszKQNRqQjWG6IiBoJpUrEt4eu4JPdf6OkTAVrMyPMf6otngnx4Jm6iR7AckNE1EgIAPZduIWSMhV6+Tnho6cD4GZrLnUsIp3DckNEpMOUKhGlShXMjOWQyQQsfiYQcZdvY2RHbq0hqg7LDRGRjkrJLsDrm06jjZsNPhjmDwDwdLCAp4OFxMmIdBsPBSci0jEqlYhvD6Vg0PIDiL92F1tP3UB2frHUsYgaDW65ISLSIdduF+D1TUk4fvUOAKBbiyZY9HQgHK1MJU5G1Hiw3BAR6QBRFLH+eCr++/t5FJYoYWkiR8STbTAm1Itja4hqiOWGiEgH3C0sxce7klFYokSXZg74+Jkgjq0hqiWWGyIiHeBgaYL/DQ/AzZx7mNTNl6dOIKoDlhsiIgnkFpZiwa9n8WSAG/q3cwUAPBngJnEqIv0g6dFSkZGR6NSpE6ytreHs7Ixhw4YhOTn5sffbtGkTWrduDTMzMwQEBOCPP/5ogLRERNpx6GI2Biw7gK2JN/HO1rMoKlVKHYlIr0habvbv348ZM2bg6NGj2LNnD0pLS9G/f38UFBRUe58jR47g+eefx+TJk3Hq1CkMGzYMw4YNw9mzZxswORFRzd0rUWLhr3/hhW+PIUNRBF9HS0SNDYGZsVzqaER6RRBFUZQ6xH23bt2Cs7Mz9u/fj549e1a5zqhRo1BQUIDt27erl3Xp0gXt27dHVFTUY59DoVDA1tYWubm5sLGx0Vp2IqJHOZ2Wg9c2JuLKrfI/3sZ28UbEk61hYcLRAUSaqMnnt05N4pebmwsAcHBwqHaduLg4hIeHV1g2YMAAxMXF1Ws2IqLaSskuwIgvj+DKrQI4W5ti7aTO+GCYP4sNUT3RmZ8slUqFWbNmoVu3bvD39692vYyMDLi4uFRY5uLigoyMjCrXLy4uRnHxvzN7KhQK7QQmItKQr6Mlnu7QFIUlSnw4zB92FiZSRyLSazpTbmbMmIGzZ8/i0KFDWn3cyMhIvPfee1p9TCKiRxFFERtOpKFva2c425gBAP47PADGcp3aWE6kt3TiJ23mzJnYvn079u3bBw8Pj0eu6+rqiszMzArLMjMz4erqWuX6ERERyM3NVV/S0tK0lpuI6GF3Ckrw0roEzNt8BnM2nYZKVT6skcWGqOFI+tMmiiJmzpyJLVu2YO/evfD19X3sfcLCwhATE1Nh2Z49exAWFlbl+qamprCxsalwISKqDwf+voUByw5gz7lMmMhl6OXnJHUkIoMk6W6pGTNmYP369di2bRusra3V42ZsbW1hbm4OABg3bhyaNm2KyMhIAMCrr76KXr16YcmSJRg8eDB++uknxMfH4+uvv5bs6yAiw1ZUqsTHu5Lx7aEUAEALZyssf6492rnbSpyMyDBJWm6+/PJLAEDv3r0rLF+zZg0mTJgAAEhNTYVM9u8Gpq5du2L9+vV455138NZbb6Fly5bYunXrIwchExHVl+t3C/Hi2nhcyMgDAIwL80bEoDYwN+HcNURS0al5bhoC57khIm0qLCnDU58dQu69Uix+JhD92rg8/k5EVGM1+fzWmaOliIgai7sFJbA1N4ZMJsDCxAhRY0NgZ2EMZ2szqaMREXTkaCkiosZi34UsPPHpfqw+nKJe5udizWJDpENYboiINFBUWn5eqInRJ5CdX4JtiTehVBnUXn2iRoO7pYiIHuPyrXy8vP4UzqWXz3A+sZsP3hzYGnKZIHEyIqoKyw0R0SP8knAd7247i8ISJRwsTbDk2SD0aeUsdSwiegSWGyKiaqTeLsSbvyShTCWia/Mm+HRUe7jYcGwNka5juSEiqoZXEwtEPNkG90rKMK13C+6GImokWG6IiP4hiiLWHL6KLs2aoK17+Twak7s//rQwRKRbWG6IiFB+wsvXN51GzIUsNHOyxO8v9+Asw0SNFMsNERm8o1du49WfTiFTUQwTIxkmdPWBmTFnyiBqrFhuiMhglSlV+HzvJXy+9yJUItDMyRIrnu+g3iVFRI0Tyw0RGaTce6X4v3XxOHrlDgBgZIgH3hvaDhYm/LVI1Njxp5iIDJKVqRFkggBLEzn+OzwAw4KbSh2JiLSE5YaIDIZKJaJMJcLESAa5TMDy54KhKCpFcycrqaMRkRZxxBwRGYTb+cUYv+Y43t/+l3qZk7Upiw2RHmK5ISK9d+LqHQz+7BAOXszGzwnXcf1uodSRiKgecbcUEektlUrE1wev4ONdyVCqRDR3ssQXY0LgYW8hdTQiqkcsN0Skl3IKSzBnY/mkfAAwtL07/jc8AJam/LVHpO/4U05EekelEjFm1TH8dVMBEyMZFg5ph+c7e0IQeG4oIkPAMTdEpHdkMgGvhfvB19ESW6Z3xehQLxYbIgPCLTdEpBcKS8pw5VYB/JvaAgDC27qgp58TTIz4NxyRoeFPPRE1einZBRi+8gjGfnuswpFQLDZEhok/+UTUqO05l4n/fH4IyZl5MJLLkJ1fInUkIpIYd0sRUaOkVIn4dM/fWLHvEgCgo7c9vhjTAc42ZhInIyKpsdwQUaNzp6AEr/50CgcvZgMAJnbzwVtPtoGxnBujiYjlhogaoa/2X8bBi9kwN5bjo6cDMLQ9T3pJRP9iuSGiRue1J/xw/e49vNyvBVq72kgdh4h0DLfhEpHOKypVYs3hFKhUIgDAzFiOlWM6sNgQUZW45YaIdFp67j3837oEJF3PRU5hKV57wk/qSESk41huiEhnxV+9g6nfn0R2fjHsLIwR4m0vdSQiagRYbohIJ60/looFv55FqVJEa1drfDOuIzwdeDZvIno8lhsi0iklZSq899tf+OFYKgBgcIAbPh4ZCAsT/roiIs3wtwUR6ZQr2fnYlHAdggDM7d8K03s350kviahGWG6ISKe0drXBx88EwtrMCH1bu0gdh4gaIZYbIpLc1lM30NzJCgEe5Wf05qR8RFQXLDdEJJkypQqLdl7ANwdT4G5rht9f6QF7SxOpYxFRI8dyQ0SSyC0sxcwfT6rPDzW8Q1PYmBtLnIqI9AHLDRE1uCu38vHi2nhcyS6AubEcn4wMwuBAN6ljEZGeYLkhogZ1+FI2pn2fAEVRGdxtzbBqfCe0dedpFIhIe1huiKhBrTl8FYqiMgR72eGrsSFwtjaTOhIR6RmWGyJqUEtHBSEq9jJe6dcSZsZyqeMQkR7iWcGJqF4pikqx9shViGL5Gb1tzIzxxsDWLDZEVG+45YaI6s212wWYvDYel7LyUaYSMbm7r9SRiMgAsNwQUb04duU2pn6fgLuFpXC1MUOor4PUkYjIQLDcEJHWbTyRhre3nkGpUkSghy2+GdcRLjYcOExEDYPlhoi0RqUSsXhXMqL2XwYADA50wyfPBMHchONriKjhsNwQkdacvZmLrw+UF5tX+7XErPCWPKM3ETU4lhsi0ppADzu8P9Qf5sZyPB3iIXUcIjJQLDdEVCdXswsglwnwdLAAALzQxVviRERk6DjPDRHVWvzVOxj+xWGMX3McOYUlUschIgLAckNEtbQ96SZGrzqGu4WlsDI1QqlSlDoSEREA7pYiohoSRRFR+69g0c4LAIAn2rpg+XPtYWHCXydEpBv424iINFamVGH+r39h/bFUAMDEbj54Z3BbyGU8IoqIdAfLDRFpbNHOC1h/LBWCALw7uC0m8XQKRKSDWG6ISGMv9WyOgxezMfsJP/Rv5yp1HCKiKrHcENEj5RWVwtrMGADgZG2K31/pwd1QRKTTeLQUEVXr7I1c9F2yHz8nXFcvY7EhIl3HckNEVTp0MRujvorDrbxifBd3FUoVD/UmosaBu6WIqJJtiTcwd9NplCpFdG3eBFFjQ7jFhogaDZYbIqpg1cEr+PD38wCApwLdsOTZIJga8azeRNR4sNwQEYDyyfn+98d5fHMwBUD5HDbvDm4LGbfYEFEjw3JDRAAAQRBgbly+hSZiUGu81LMZBIHFhogaH0kHFB84cABDhgyBu7s7BEHA1q1bH7l+bGwsBEGodMnIyGiYwER67rUn/PDLtDD8X6/mLDZE1GhJWm4KCgoQFBSElStX1uh+ycnJSE9PV1+cnZ3rKSGRflMUleL9387hXokSQPnWmxBvB4lTERHVjaS7pQYNGoRBgwbV+H7Ozs6ws7PTfiAiA5KdX4xx3x7HuXQFbuUX4/Png6WORESkFY1ynpv27dvDzc0NTzzxBA4fPvzIdYuLi6FQKCpciAzdjZx7eDYqDufSFXC0MsHUXs2kjkREpDWNqty4ubkhKioKv/zyC3755Rd4enqid+/eOHnyZLX3iYyMhK2trfri6enZgImJdM/lW/kY+eURXMkuQFM7c2z8vzC0c7eVOhYRkdYIoijqxLSjgiBgy5YtGDZsWI3u16tXL3h5eWHdunVV3l5cXIzi4mL1dYVCAU9PT+Tm5sLGxqYukYkanbM3cjFu9XHcKShBcydLfP9iKNxszaWORUT0WAqFAra2thp9fjf6Q8E7d+6MQ4cOVXu7qakpTE1NGzARkW4qU6ow/YeTuFNQAv+mNlg7sTOaWPFng4j0T6PaLVWVxMREuLm5SR2DSOcZyWVYMToY/Vo7Y/2ULiw2RKS3tLLlRqlU4syZM/D29oa9vb3G98vPz8elS5fU11NSUpCYmAgHBwd4eXkhIiICN27cwHfffQcAWLZsGXx9fdGuXTsUFRVh1apV2Lt3L3bv3q2NL4NIL+XeK4WtuTEAINDDDt9O6CRxIiKi+lWrLTezZs3Ct99+C6C82PTq1QsdOnSAp6cnYmNjNX6c+Ph4BAcHIzi4/BDU2bNnIzg4GPPnzwcApKenIzU1Vb1+SUkJ5syZg4CAAPTq1QunT5/Gn3/+iX79+tXmyyDSe9sSb6Dn4n04lXpX6ihERA2mVgOKPTw8sHXrVnTs2BFbt27FjBkzsG/fPqxbtw579+597OHZUqrJgCSixmxjfBre/CUJoghM6OqDhf9pJ3UkIqJaq8nnd6223GRnZ8PV1RUA8Mcff2DkyJHw8/PDpEmTcObMmdo8JBFp0Q/HruGNn8uLzehQL8x/qq3UkYiIGkytyo2LiwvOnTsHpVKJnTt34oknngAAFBYWQi6XazUgEdXM6kMpeHvLWQDlZ/b+7zB/ntmbiAxKrQYUT5w4Ec8++yzc3NwgCALCw8MBAMeOHUPr1q21GpCINPdl7GUs2nkBAPB/vZph3sDWPAEmERmcWpWbhQsXwt/fH2lpaRg5cqR6Hhm5XI558+ZpNSARaUapEnE85TYA4JV+LfFaeEsWGyIySFqboTgnJ6dRnMySA4pJnxWVKrHrrwwMbd9U6ihERFpV7wOKFy1ahA0bNqivP/vss2jSpAk8PDyQlJRUm4ckoloQRRH7krNw/28UM2M5iw0RGbxalZuoqCj1CSj37NmDPXv2YMeOHRg4cCDmzp2r1YBEVDVRFLFoZzImrjmBj3ZckDoOEZHOqNWYm4yMDHW52b59O5599ln0798fPj4+CA0N1WpAIqpMFEV8vCsZUfsvAwCa2vPkl0RE99Vqy429vT3S0tIAADt37lQfLSWKIpRKpfbSEVEloihiye6/8UVsebFZOKQtxoX5SBuKiEiH1GrLzYgRIzB69Gi0bNkSt2/fxqBBgwAAp06dQosWLbQakIgq+vTPi1ixr/ycbPOfaosJ3XwlTkREpFtqVW4+/fRT+Pj4IC0tDYsXL4aVlRWA8nNBTZ8+XasBiehfy/+8iM9iLgIA3hncBpO6s9gQET2sVuXG2Ni4yoHDr732Wp0DEVH1mtqbQyYAbz3ZBi/2aCZ1HCIinVSrMTcAsG7dOnTv3h3u7u64du0aAGDZsmXYtm2b1sIRUUXPhHhg16yeLDZERI9Qq3Lz5ZdfYvbs2Rg0aBBycnLUg4jt7OywbNkybeYjMng3c+6huOzfgfotXawlTENEpPtqVW4+//xzfPPNN3j77bcrnCizY8eOPCs4kZbN2pCInov34eiV21JHISJqFGo15iYlJQXBwcGVlpuamqKgoKDOoYio3NErt3E85Q5M5DL4NLGUOg4RUaNQqy03vr6+SExMrLR8586daNOmTV0zEdE/Pt9bfmTUs5084GprJnEaIqLGoVZbbmbPno0ZM2agqKgIoiji+PHj+PHHHxEZGYlVq1ZpOyORQYq/egeHL92GsVzAtN6cP4qISFO1KjcvvvgizM3N8c4776CwsBCjR4+Gu7s7li9fjueee07bGYkM0md7yyfqeybEA03teHoFIiJN1arcAMCYMWMwZswYFBYWIj8/H87OztrMRWTQEtNycODvW5DLBEzrxa02REQ1Uetyc5+FhQUsLCy0kYWI/nHw71sAgOHBTeHVhD9fREQ1Uatyk5mZiblz5yImJgZZWVkQRbHC7Tx5JlHdvNyvJXq3coadhbHUUYiIGp1alZsJEyYgNTUV7777Ltzc3CAIgrZzERm8AA9bqSMQETVKtSo3hw4dwsGDB9G+fXstxyEybGl3CmEsl/GwbyKiOqjVPDeenp6VdkURUd19tOMCei7eh43xaVJHISJqtGpVbpYtW4Z58+bh6tWrWo5DZLguZubhj7PpKFGqEMhdUkREtVar3VKjRo1CYWEhmjdvDgsLCxgbVxz0eOfOHa2EIzIkK/ZdgigCA9u5orWrjdRxiIgarVqVm08//ZSDiIm06PKtfPx2+iYA4OV+nNeGiKguan20VHXu3btX2yxEBmvlvktQiUB4G2e0c+cuKSKiuqjVmJtXXnmlyuUFBQV48skn6xSIyNCk3i7EtsR/ttr0bSlxGiKixq9W5eb333/HggULKiwrKCjAwIEDUVZWppVgRIbizI1cGMsF9PRzQpCnndRxiIgavVrtltq9ezd69OgBe3t7zJo1C3l5eRgwYACMjIywY8cObWck0muDA90Q2swBBcX8w4CISBtqVW6aN2+OnTt3ok+fPpDJZPjxxx9hamqK33//HZaWltrOSKT3HK1M4WhlKnUMIiK9UKvdUgAQGBiI7du346233oKFhQV27NjBYkNUA9n5xUi4xmkTiIi0TeMtN8HBwVUe/m1qaoqbN2+iW7du6mUnT57UTjoiPfbNgSv46sAVTOzmgwVD2kkdh4hIb2hcboYNG1aPMYgMy92CEnx/9BoAoEdLR4nTEBHpF43LzcNHRxFR7a05chUFJUq0dbNBn1bOUschItIrtRpzc+LECRw7dqzS8mPHjiE+Pr7OoYj0maKoFNGHUwAAL/dtwdm+iYi0rFblZsaMGUhLq3zW4hs3bmDGjBl1DkWkz9bFXYOiqAwtnK0woJ2r1HGIiPROrcrNuXPn0KFDh0rLg4ODce7cuTqHItJXhSVl+PZQ+VabmX1aQCbjVhsiIm2rVbkxNTVFZmZmpeXp6ekwMqrV1DlEBuHa7UJYmMjh3cQCTwW6SR2HiEgvCaIoijW90/PPP4/09HRs27YNtrblJ/nLycnBsGHD4OzsjI0bN2o9qLYoFArY2toiNzcXNjY2UschA1SqVOHG3XvwceS8UEREmqrJ53etys2NGzfQs2dP3L59G8HBwQCAxMREuLi4YM+ePfD09Kxd8gbAckNERNT41OTzu1b7kJo2bYqkpCT88MMPOH36NMzNzTFx4kQ8//zzMDY2rlVoIn1WqlRhe9JNPBXoDmN5rScGJyIiDdRqy01jxi03JIVfEq5jzqbTaO9phy3Tu/LwbyKiGqqXLTe//vorBg0aBGNjY/z666+PXPc///mPpg9LpPdUKhFfHbgMAOjfzoXFhoiontXo9AsZGRlwdnZ+5KkYBEGAUqnURjYivbD3Qhb+zsyHlakRxoR6Sx2HiEjvaVxuVCpVlf8nokeL2l++1WZMFy/YmnNMGhFRfePIRqJ6dOLqHcRfuwsTuQyTu/lKHYeIyCDUesa9mJgYxMTEICsrq9KWnNWrV9c5GJE+iIot32ozokNTONuYSZyGiMgw1KrcvPfee3j//ffRsWNHuLm5cYAkURVKylQoLlNBEICXejaTOg4RkcGo1aHgbm5uWLx4McaOHVsfmeoVDwWnhnY1u4CzERMR1VFNPr9rNeampKQEXbt2rVU4IkPDYkNE1LBqVW5efPFFrF+/XttZiPTG3guZuJVXLHUMIiKDpPGYm9mzZ6v/r1Kp8PXXX+PPP/9EYGBgpVMuLF26VHsJiRqZuwUlmPHDKahEEb+/0gMtnK2kjkREZFA0LjenTp2qcL19+/YAgLNnz2o1EFFjt+7oNdwrVaKduw2aO3GXFBFRQ9O43Ozbt68+cxDphaJSJb6Luwqg/AgpHklIRNTwajXmZtKkScjLy6u0vKCgAJMmTapzKKLG6tfTN5GdXwI3WzM8GeAmdRwiIoNUq3Kzdu1a3Lt3r9Lye/fu4bvvvqtzKKLGSBRFfHswBQAwvqsPjOWcAJyISAo1msRPoVBAFEWIooi8vDyYmf0746pSqcQff/wBZ2dnrYckagwOXcpGcmYeLEzkeL6Tl9RxiIgMVo3KjZ2dHQRBgCAI8PPzq3S7IAh47733tBaOqDG5cqsAJkYyPNvRE7YWPEEmEZFUalRu9u3bB1EU0bdvX/zyyy9wcHBQ32ZiYgJvb2+4u7trPSRRYzC+qw8GB7qh5nN+ExGRNtVoUECvXr3Qu3dvpKSkYNiwYejVq5f6EhYWVqnYTJ8+HdnZ2dU+3oEDBzBkyBC4u7tDEARs3br1sRliY2PRoUMHmJqaokWLFoiOjq7Jl0BUrxytTOFkbSp1DCIig1arEY/e3t4aHeL6/fffQ6FQVHt7QUEBgoKCsHLlSo2eNyUlBYMHD0afPn2QmJiIWbNm4cUXX8SuXbs0zk6kbXcLSnD2Rq7UMYiI6B+1Oiu4ph53Ts5BgwZh0KBBGj9eVFQUfH19sWTJEgBAmzZtcOjQIXz66acYMGBAnbIS1da6o9ewdM/fGBfmjfeH+ksdh4jI4DWqY1Xj4uIQHh5eYdmAAQMQFxdX7X2Ki4uhUCgqXIi05cFJ+0K87aUNQ0REABpZucnIyICLi0uFZS4uLlAoFFXOuwMAkZGRsLW1VV88PT0bIioZiN+T0jlpHxGRjmlU5aY2IiIikJubq76kpaVJHYn0yP2tNi908eakfUREOqJex9xom6urKzIzMyssy8zMhI2NDczNzau8j6mpKUxNefQKaV9iWg5OX8+FiVyG5zpxiyARka6o1z81X3jhBdjY2Gjt8cLCwhATE1Nh2Z49exAWFqa15yDS1HdHrgIAngpyQxMrFmgiIl1Roy03SUlJGq0XGBgIAPjyyy8fuV5+fj4uXbqkvp6SkoLExEQ4ODjAy8sLERERuHHjhvp8VVOnTsWKFSvwxhtvYNKkSdi7dy82btyI33//vSZfBlGdlSpVOHHtDgBgXJiPtGGIiKgCQXzc8doPkMlkEAShykO87y8XBAFKpVKjx4uNjUWfPn0qLR8/fjyio6MxYcIEXL16FbGxsRXu89prr+HcuXPw8PDAu+++iwkTJmj6JUChUMDW1ha5ubla3apEhqdUqcKhS9no04rnUyMiqm81+fyuUbm5du2aRut5e3tr+pANjuWGiIio8anJ53eNdkvpcmkhaihZeUVoYmkKuezxs3QTEVHDq9GA4sWLF1eYT+bw4cMoLi5WX8/Ly8P06dO1l45IB722IRE9F+9D3OXbUkchIqIq1KjcREREIC8vT3190KBBuHHjhvp6YWEhvvrqK+2lI9Ixl7LycPjSbaTn3oOnQ9XTDxARkbRqVG4eHp5Tg+E6RHrh+6OpAIB+bVzgYW8hcRoiIqoKp1Ql0lBRqRKbT14HAIztwvFnRES6iuWGSEN/nEmHoqgMHvbm6N7CUeo4RERUjRqffmHVqlWwsrICAJSVlSE6OhqOjuW/6B8cj0Okb346Xn5esuc6eULGI6WIiHRWjea58fHxgSA8/pd6SkpKnULVJ85zQ7VxNbsAvT+JhVwm4Mi8vnCxMZM6EhGRQam3eW6uXr1al1xEjZZ3Ewtsmd4Vp9NyWGyIiHRcjcbcxMXFYfv27RWWfffdd/D19YWzszNeeumlCvPeEOkLQRAQ7GWPCd18pY5CRESPUaNy89577+Gvv/5SXz9z5gwmT56M8PBwzJs3D7/99hsiIyO1HpJISpzygIiocalRuTl9+jT69eunvv7TTz8hNDQU33zzDWbPno3PPvsMGzdu1HpIIilN+S4BEZvP4EbOvcevTEREkqtRubl79y5cXFzU1/fv349Bgwapr3fq1AlpaWnaS0cksZTsAvx5PhMbTqRKHYWIiDRUo3Lj4uKiPhKqpKQEJ0+eRJcuXdS35+XlwdjYWLsJiSS0Mb68rPfyc0JTO55ugYioMahRuXnyyScxb948HDx4EBEREbCwsECPHj3UtyclJaF58+ZaD0kkBZVKxLZT5edOG9nRU+I0RESkqRodCv7BBx9gxIgR6NWrF6ysrLB27VqYmJiob1+9ejX69++v9ZBEUjiWcgc3c4tgbWaEvq2dpY5DREQaqlG5cXR0xIEDB5CbmwsrKyvI5fIKt2/atEk9ezFRY7flVPl5pAYHuMHMWP6YtYmISFfU+PQLAGBra1vlcgcHhzqFIdIVRaVK7DiTAQAYHtxU4jRERFQTtSo3RPquuEyFCd18cDzlDjr5sLQTETUmLDdEVbA1N8ac/q2kjkFERLVQo6OliIiIiHQdt9wQPSQ2OQulShG9/JxgYsT+T0TU2PA3N9FDPv3zIqZ8F48fj3NWYiKixojlhugBV27l43RaDuQyAU8GuEkdh4iIaoHlhugBW/+ZkbhHS0c4WZtKnIaIiGqD5YboH6IoYktiebnh3DZERI0Xyw3RP+Kv3UXanXuwNJGjf1tXqeMQEVEtsdwQ/WPLP7ukBvq7wdyEp1sgImqsWG6I/nE1uwAAMKIDd0kRETVmnOeG6B/rp3TB35l5aO7Ek78SETVmLDdED/BzsZY6AhER1RF3S5HBKywpQ35xmdQxiIhIS1huyOBtOXUDHT/cgyW7k6WOQkREWsByQwZvy8kbKCpVwcqUe2mJiPQByw0ZtNTbhYi/dheCAAxtz6OkiIj0AcsNGbSt/8xI3K25I1xtzSROQ0RE2sByQwZLFEX1uaR4ugUiIv3BckMG6/T1XFzJLoCZsQwD/Hm6BSIifcFyQwbr/labAe1cOZiYiEiP8Dc6GayXejaDi40ZOvrYSx2FiIi0iOWGDJa7nTmm9W4udQwiItIy7pYiIiIivcJyQwbnbkEJJqw5jm2JNyCKotRxiIhIy7hbigzOr6dvIjb5Fm7lFXPiPiIiPcQtN2Rwfjl5HQDwTIiHxEmIiKg+sNyQQUnOyEPS9VwYyQT8J8hd6jhERFQPWG7IoHx/9BoAILyNC5pYmUqchoiI6gPLDRmM/OIybP5nl9TYMG+J0xARUX1huSGDseXkdRSUKNHcyRJdmzeROg4REdUTHi1FBsPDwQIh3vYYEugGQRCkjkNERPWE5YYMRp9WzujTyhkqFee2ISLSZ9wtRQZHJuNWGyIifcZyQ3ov7vJtrNh7EYqiUqmjEBFRA+BuKdJrSpWID7afw7l0BXLvleLtwW2ljkRERPWMW25Ir/2ckIZz6QpYmxlhWu8WUschIqIGwHJDeutOQQk+3pUMAHi1X0s4WJpInIiIiBoCyw3prXe3nUV2fglaOlthXJiP1HGIiKiBsNyQXtp88jp+T0qHXCZg6bPtYWLEtzoRkaHgb3zSO8VlSizZ/TcA4OW+LRDgYStxIiIiakg8Wor0jqmRHJumhmHtkat4uW9LqeMQEVED45Yb0htFpUr1/93tzBHxZBvIOWEfEZHB0Ylys3LlSvj4+MDMzAyhoaE4fvx4tetGR0dDEIQKFzMzswZMS7ro8q189FuyH3+cSZc6ChERSUzycrNhwwbMnj0bCxYswMmTJxEUFIQBAwYgKyur2vvY2NggPT1dfbl27VoDJiZdk5yRh1FfHcWNnHtYsfcSypQqqSMREZGEJC83S5cuxZQpUzBx4kS0bdsWUVFRsLCwwOrVq6u9jyAIcHV1VV9cXFwaMDHpkrM3cvHc13HIzi9GWzcbrJvcGUZyyd/WREQkIUk/BUpKSpCQkIDw8HD1MplMhvDwcMTFxVV7v/z8fHh7e8PT0xNDhw7FX3/91RBxScecSr2L0d8cxd3CUgR52GL9lFA0sTKVOhYREUlM0nKTnZ0NpVJZacuLi4sLMjIyqrxPq1atsHr1amzbtg3ff/89VCoVunbtiuvXr1e5fnFxMRQKRYULNX4nrt7B2G+PQ1FUhhBve6x7MRR2FpyBmIiIdGC3VE2FhYVh3LhxaN++PXr16oXNmzfDyckJX331VZXrR0ZGwtbWVn3x9PRs4MRUH3b/lYH84jKENWuC7yZ1ho2ZsdSRiIhIR0g6z42joyPkcjkyMzMrLM/MzISrq6tGj2FsbIzg4GBcunSpytsjIiIwe/Zs9XWFQsGCowciBrVBUztzPNfZC2bGcqnjEBGRDpF0y42JiQlCQkIQExOjXqZSqRATE4OwsDCNHkOpVOLMmTNwc3Or8nZTU1PY2NhUuFDjFH/1DkrKyo+EkskETOjmy2JDRESVSL5bavbs2fjmm2+wdu1anD9/HtOmTUNBQQEmTpwIABg3bhwiIiLU67///vvYvXs3rly5gpMnT+KFF17AtWvX8OKLL0r1JVAD2J50E899fRSv/HiKh3oTEdEjSX76hVGjRuHWrVuYP38+MjIy0L59e+zcuVM9yDg1NRUy2b8d7O7du5gyZQoyMjJgb2+PkJAQHDlyBG3btpXqS6B69kvCdbz+82moRMDChFtqiIjo0QRRFEWpQzQkhUIBW1tb5ObmchdVI7DxRBre3JwEUQSe6+SJ/w4P4CkViIgMUE0+vyXfckNUnR+PpyJi8xkAwNgu3njvP+0gY7EhIqLHYLkhnfRgsZnQ1QcLhrSFILDYEBHR47HckE7ydbSEmbEMY0K98c7gNiw2RESkMZYb0kldmjXBjld7wqeJBYsNERHViOSHghPdty3xBi5m5qmv+zpastgQEVGNsdyQTtiedBOvbUjEs1/F4UbOPanjEBFRI8ZyQ5KLOZ+JWT8lQiUCA/1d4W5rJnUkIiJqxFhuSFJHLmVj2g8nUaYSMbS9Oz4cFsBdUUREVCcsNySZszdyMeW7eJSUqfBEWxd8MjKIE/QREVGdsdyQJK7fLcTE6BMoKFGiSzMHfP58MIzlfDsSEVHd8VBwkoSVqRF8mljAwcIEX43tyLN7ExGR1rDckCTsLEywbnIoFEWlsDU3ljoOERHpEe4HoAYjiiL2/31Lfd3MWA5nax4ZRURE2sVyQw3ms5hLGL/6OJbsTpY6ChER6TGWG2oQO8+m49M//wYAeNibS5yGiIj0GcsN1btzNxV4bcNpAMDEbj4Y1clL4kRERKTPWG6oXmXnF2PKd/G4V6pEj5aOePvJNlJHIiIiPcdyQ/WmTKnCy+tP4UbOPfg0scCK5zvAiHPZEBFRPeMnDdWbQ5eyEXflNixN5Fg1viNsLXjINxER1T/Oc0P1pncrZ3w9NgRKlYgWztZSxyEiIgPBckP1qn87V6kjEBGRgeFuKdKq4jIl3t5yBum596SOQkREBorlhrRq6Z6/8cOxVIz55hiUKlHqOEREZIBYbkhrEq7dwdcHrgAA3hzUGnKZIHEiIiIyRCw3pBWFJWWYs/E0RBF4uoMHBnCsDRERSYTlhrRi8c5kXL1dCFcbM8wf0lbqOEREZMBYbqjOjlzKRvSRqwCAxc8Ewtac89kQEZF0WG6ozqL+GWczOtQLPf2cJE5DRESGjvPcUJ19PTYEXx+4gkndfaWOQkRExHJDdWdmLMcr/VpKHYOIiAgAd0tRLeUUlmDN4RTOZUNERDqHW26oVhb++he2Jt7EmRu5WPpse6njEBERqXHLDdXYzrPp2Jp4EzIBeKGLt9RxiIiIKmC5oRrJzi/G21vOAgCm9mqODl72EiciIiKqiOWGNCaKIt7Zcha3C0rQ2tUar4ZzEDEREekelhvS2K+nb2LnXxkwkgn4ZGQQTI3kUkciIiKqhOWGNFJUqsSHv58HALzSryX8m9pKnIiIiKhqLDekETNjOdZN7oxnQjwwrXdzqeMQERFVi4eCk8Zau9rgk5FBUscgIiJ6JG65oUe6frcQSddzpI5BRESkMZYbqpZKJWLuptMY/sURbIpPkzoOERGRRlhuqFqrDl3B0St3YGokQycfB6njEBERaYTlhqp0Pl2BT3b9DQB496m28HG0lDgRERGRZlhuqJKiUiVe25CIEqUK4W2c8VwnT6kjERERaYzlhipZtPMCLmTkoYmlCSJHBEIQBKkjERERaYzlhipIuHYHaw5fBQAsejoQTtam0gYiIiKqIc5zQxV08LLHwiFtcbugBOFtXaSOQ0REVGMsN1SBIAiY0M1X6hhERES1xt1ShOIyJSJ3nEfuvVKpoxAREdUZy42BK1OqMGfjaXy1/womRZ+AKIpSRyIiIqoTlhsDVqZUYfbG09ielA4jmYBZ4S15ZBQRETV6HHNjoO7PZbPjbAaMZAJWjumAHi2dpI5FRERUZyw3BigrrwhTvkvA6bQcGMsFrBzdAf3buUodi4iISCtYbgxQxC9ncDotB3YWxoh6IQRdmjWROhIREZHWsNwYoIgn2+BWfjE+ey6Y54wiIiK9wwHFek4URWw5dR1L9/ytXtbC2QrbZnRjsSEiIr3ELTd6SqUSsed8Jr7Ydwmnr+dCJgAD27mirbsNAPCoKCIi0lssN3qmTKnCb0k38cW+y7iYlQ8AMDeWY2bfFmjuzC01RESk/1hu9MhfN3Mx9fsEpN25BwCwNjXCuK7emNjNF45WPAEmEREZBpabRk4URfUuJi8HC+QWlqKJpQkmdffF2DBv2JgZS5yQiIioYbHcNFJ3Ckqw9shVJFy7i3WTO0MQBFibGSN6Ume0cbWBuYlc6ohERESSYLlpZK5mF+DbQynYlJCGolIVAODI5dvo1sIRANDBy17KeERERJLTiUPBV65cCR8fH5iZmSE0NBTHjx9/5PqbNm1C69atYWZmhoCAAPzxxx8NlFQaSpWI2OQs/N+6ePRZEot1R6+hqFSFgKa2WDE6mJPwERERPUDyLTcbNmzA7NmzERUVhdDQUCxbtgwDBgxAcnIynJ2dK61/5MgRPP/884iMjMRTTz2F9evXY9iwYTh58iT8/f0l+Ar+pVKJuJFzD6IIiBAhioBKFCECEEXA0lQON1tz9frnbiqg+ucs3BXXFWFtZoQWztYAgPirdzBhzQn1/fq0csJLPZujSzMHHtJNRET0EEEU//l0lUhoaCg6deqEFStWAABUKhU8PT3x8ssvY968eZXWHzVqFAoKCrB9+3b1si5duqB9+/aIiop67PMpFArY2toiNzcXNjY22vtCABSWlKHt/F3V3j7I3xVfvhCivu4z7/dq1x3QzgVfje0IoLw0Df/iMIK97DE61At+LtbaC01ERNQI1OTzW9ItNyUlJUhISEBERIR6mUwmQ3h4OOLi4qq8T1xcHGbPnl1h2YABA7B169Yq1y8uLkZxcbH6ukKhqHvwasgEAWbGMggQIBPKJ8oTBEAAIJMJsDCp+HI7W5uq7ycI5f8CgCAAZsb/DgiWyQRsndGNW2mIiIg0IGm5yc7OhlKphIuLS4XlLi4uuHDhQpX3ycjIqHL9jIyMKtePjIzEe++9p53Aj2FmLMeFDwZpvP7xt8M1XpfFhoiISDM6MaC4PkVERCA3N1d9SUtLkzoSERER1SNJt9w4OjpCLpcjMzOzwvLMzEy4urpWeR9XV9carW9qagpTU87OS0REZCgk3XJjYmKCkJAQxMTEqJepVCrExMQgLCysyvuEhYVVWB8A9uzZU+36REREZFgkPxR89uzZGD9+PDp27IjOnTtj2bJlKCgowMSJEwEA48aNQ9OmTREZGQkAePXVV9GrVy8sWbIEgwcPxk8//YT4+Hh8/fXXUn4ZREREpCMkLzejRo3CrVu3MH/+fGRkZKB9+/bYuXOnetBwamoqZLJ/NzB17doV69evxzvvvIO33noLLVu2xNatWyWf44aIiIh0g+Tz3DS0+pznhoiIiOpHTT6/9f5oKSIiIjIsLDdERESkV1huiIiISK+w3BAREZFeYbkhIiIivcJyQ0RERHqF5YaIiIj0CssNERER6RXJZyhuaPfnLFQoFBInISIiIk3d/9zWZO5hgys3eXl5AABPT0+JkxAREVFN5eXlwdbW9pHrGNzpF1QqFW7evAlra2sIgqDVx1YoFPD09ERaWhpP7fAYfK00x9dKc3ytNMPXSXN8rTRX36+VKIrIy8uDu7t7hXNOVsXgttzIZDJ4eHjU63PY2Njwh0BDfK00x9dKc3ytNMPXSXN8rTRXn6/V47bY3McBxURERKRXWG6IiIhIr7DcaJGpqSkWLFgAU1NTqaPoPL5WmuNrpTm+Vprh66Q5vlaa06XXyuAGFBMREZF+45YbIiIi0issN0RERKRXWG6IiIhIr7DcaMnKlSvh4+MDMzMzhIaG4vjx41JH0kkLFy6EIAgVLq1bt5Y6luQOHDiAIUOGwN3dHYIgYOvWrRVuF0UR8+fPh5ubG8zNzREeHo6LFy9KE1Zij3utJkyYUOk9NnDgQGnCSiwyMhKdOnWCtbU1nJ2dMWzYMCQnJ1dYp6ioCDNmzECTJk1gZWWFp59+GpmZmRIllo4mr1Xv3r0rvbemTp0qUWLpfPnllwgMDFTPZxMWFoYdO3aob9eF9xTLjRZs2LABs2fPxoIFC3Dy5EkEBQVhwIAByMrKkjqaTmrXrh3S09PVl0OHDkkdSXIFBQUICgrCypUrq7x98eLF+OyzzxAVFYVjx47B0tISAwYMQFFRUQMnld7jXisAGDhwYIX32I8//tiACXXH/v37MWPGDBw9ehR79uxBaWkp+vfvj4KCAvU6r732Gn777Tds2rQJ+/fvx82bNzFixAgJU0tDk9cKAKZMmVLhvbV48WKJEkvHw8MDH330ERISEhAfH4++ffti6NCh+OuvvwDoyHtKpDrr3LmzOGPGDPV1pVIpuru7i5GRkRKm0k0LFiwQg4KCpI6h0wCIW7ZsUV9XqVSiq6ur+PHHH6uX5eTkiKampuKPP/4oQULd8fBrJYqiOH78eHHo0KGS5NF1WVlZIgBx//79oiiWv4+MjY3FTZs2qdc5f/68CECMi4uTKqZOePi1EkVR7NWrl/jqq69KF0qH2dvbi6tWrdKZ9xS33NRRSUkJEhISEB4erl4mk8kQHh6OuLg4CZPprosXL8Ld3R3NmjXDmDFjkJqaKnUknZaSkoKMjIwK7zFbW1uEhobyPVaN2NhYODs7o1WrVpg2bRpu374tdSSdkJubCwBwcHAAACQkJKC0tLTCe6t169bw8vIy+PfWw6/VfT/88AMcHR3h7++PiIgIFBYWShFPZyiVSvz0008oKChAWFiYzrynDO7cUtqWnZ0NpVIJFxeXCstdXFxw4cIFiVLprtDQUERHR6NVq1ZIT0/He++9hx49euDs2bOwtraWOp5OysjIAIAq32P3b6N/DRw4ECNGjICvry8uX76Mt956C4MGDUJcXBzkcrnU8SSjUqkwa9YsdOvWDf7+/gDK31smJiaws7OrsK6hv7eqeq0AYPTo0fD29oa7uzuSkpLw5ptvIjk5GZs3b5YwrTTOnDmDsLAwFBUVwcrKClu2bEHbtm2RmJioE+8plhtqUIMGDVL/PzAwEKGhofD29sbGjRsxefJkCZORvnjuuefU/w8ICEBgYCCaN2+O2NhY9OvXT8Jk0poxYwbOnj3LMW4aqO61eumll9T/DwgIgJubG/r164fLly+jefPmDR1TUq1atUJiYiJyc3Px888/Y/z48di/f7/UsdS4W6qOHB0dIZfLK40Ez8zMhKurq0SpGg87Ozv4+fnh0qVLUkfRWfffR3yP1U6zZs3g6Oho0O+xmTNnYvv27di3bx88PDzUy11dXVFSUoKcnJwK6xvye6u616oqoaGhAGCQ7y0TExO0aNECISEhiIyMRFBQEJYvX64z7ymWmzoyMTFBSEgIYmJi1MtUKhViYmIQFhYmYbLGIT8/H5cvX4abm5vUUXSWr68vXF1dK7zHFAoFjh07xveYBq5fv47bt28b5HtMFEXMnDkTW7Zswd69e+Hr61vh9pCQEBgbG1d4byUnJyM1NdXg3luPe62qkpiYCAAG+d56mEqlQnFxse68pxps6LIe++mnn0RTU1MxOjpaPHfunPjSSy+JdnZ2YkZGhtTRdM6cOXPE2NhYMSUlRTx8+LAYHh4uOjo6illZWVJHk1ReXp546tQp8dSpUyIAcenSpeKpU6fEa9euiaIoih999JFoZ2cnbtu2TUxKShKHDh0q+vr6ivfu3ZM4ecN71GuVl5cnzp07V4yLixNTUlLEP//8U+zQoYPYsmVLsaioSOroDW7atGmira2tGBsbK6anp6svhYWF6nWmTp0qenl5iXv37hXj4+PFsLAwMSwsTMLU0njca3Xp0iXx/fffF+Pj48WUlBRx27ZtYrNmzcSePXtKnLzhzZs3T9y/f7+YkpIiJiUlifPmzRMFQRB3794tiqJuvKdYbrTk888/F728vEQTExOxc+fO4tGjR6WOpJNGjRolurm5iSYmJmLTpk3FUaNGiZcuXZI6luT27dsnAqh0GT9+vCiK5YeDv/vuu6KLi4toamoq9uvXT0xOTpY2tEQe9VoVFhaK/fv3F52cnERjY2PR29tbnDJlisH+oVHV6wRAXLNmjXqde/fuidOnTxft7e1FCwsLcfjw4WJ6erp0oSXyuNcqNTVV7Nmzp+jg4CCampqKLVq0EF9//XUxNzdX2uASmDRpkujt7S2amJiITk5OYr9+/dTFRhR14z3Fs4ITERGRXuGYGyIiItIrLDdERESkV1huiIiISK+w3BAREZFeYbkhIiIivcJyQ0RERHqF5YaIiIj0CssNERER6RWWGyKS1IQJEzBs2LAGf97o6GgIggBBEDBr1iz1ch8fHyxbtuyR971/Pzs7u3rNSES1YyR1ACLSX4IgPPL2BQsWYPny5ZBqonQbGxskJyfD0tKyRvdLT0/Hhg0bsGDBgnpKRkR1wXJDRPUmPT1d/f8NGzZg/vz5SE5OVi+zsrKClZWVFNEAlJcvV1fXGt/P1dUVtra29ZCIiLSBu6WIqN64urqqL7a2tuoycf9iZWVVabdU79698fLLL2PWrFmwt7eHi4sLvvnmGxQUFGDixImwtrZGixYtsGPHjgrPdfbsWQwaNAhWVlZwcXHB2LFjkZ2dXavchYWFmDRpEqytreHl5YWvv/66Li8DETUwlhsi0jlr166Fo6Mjjh8/jpdffhnTpk3DyJEj0bVrV5w8eRL9+/fH2LFjUVhYCADIyclB3759ERwcjPj4eOzcuROZmZl49tlna/X8S5YsQceOHXHq1ClMnz4d06ZNq7DFiYh0G8sNEemcoKAgvPPOO2jZsiUiIiJgZmYGR0dHTJkyBS1btsT8+fNx+/ZtJCUlAQBWrFiB4OBg/O9//0Pr1q0RHByM1atXY9++ffj7779r/PxPPvkkpk+fjhYtWuDNN9+Eo6Mj9u3bp+0vk4jqCcfcEJHOCQwMVP9fLpejSZMmCAgIUC9zcXEBAGRlZQEATp8+jX379lU5fufy5cvw8/Or9fPf35V2/7mISPex3BCRzjE2Nq5wXRCECsvuH4WlUqkAAPn5+RgyZAgWLVpU6bHc3Ny08vz3n4uIdB/LDRE1eh06dMAvv/wCHx8fGBnx1xqRoeOYGyJq9GbMmIE7d+7g+eefx4kTJ3D58mXs2rULEydOhFKplDoeETUwlhsiavTc3d1x+PBhKJVK9O/fHwEBAZg1axbs7Owgk/HXHJGhEUSppgYlIpJQdHQ0Zs2ahZycHEnuT0T1h3/SEJHBys3NhZWVFd58880a3c/KygpTp06tp1REVFfcckNEBikvLw+ZmZkAADs7Ozg6Omp830uXLgEoP0zd19e3XvIRUe2x3BAREZFe4W4pIiIi0issN0RERKRXWG6IiIhIr7DcEBERkV5huSEiIiK9wnJDREREeoXlhoiIiPQKyw0RERHpFZYbIiIi0iv/Dy0aR28knNE0AAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from matplotlib import lines\n",
    "\n",
    "\n",
    "for out_var in out_vars:  # Loop through the output variables\n",
    "    for sol in sols_SPM:  # loop through the solutions of the SPM models\n",
    "        plt.plot(sol['Time [h]'].entries, sol[out_var].entries,\n",
    "                 linestyle='--',\n",
    "                 label=models_SPM[sols_SPM.index(sol)].name)\n",
    "\n",
    "    plt.xlabel(\"Time [h]\")\n",
    "    plt.ylabel(yaxis_labels[out_vars.index(out_var)])\n",
    "    plt.legend()\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "env",
   "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.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}