{ "cells": [ { "cell_type": "markdown", "id": "5d7dd3f0", "metadata": {}, "source": [ "# Analytical benchmarks\n", "\n", "Two analytical solutions from `pytes.validation` check the packed-bed model: Schumann's\n", "solution for a bed of a sensible material, and the melting of a sphere (the spherical\n", "Stefan problem). Both models are built from the packed bed's own submodels with `pytes.ComposedModel` (see [composing models](composing_models.ipynb)), swapping only what each benchmark needs." ] }, { "cell_type": "code", "execution_count": 1, "id": "3ee69c76", "metadata": { "execution": { "iopub.execute_input": "2026-10-23T17:49:04.887664Z", "iopub.status.busy": "2026-10-23T17:49:04.887408Z", "iopub.status.idle": "2026-10-23T17:49:05.873797Z", "shell.execute_reply": "2026-10-23T17:49:05.872642Z" } }, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "import numpy as np\n", "import pybamm\n", "\n", "import pytes\n", "from pytes import validation\n", "from pytes.materials import SensibleMaterial\n", "from pytes.submodels import (\n", " BaseSubModel,\n", " DeviceEnergy,\n", " EnthalpyStorage,\n", " WallCoupling,\n", ")\n", "from pytes.submodels.storage import INITIAL_TEMPERATURE\n", "\n", "base = pytes.PackedBedModel()\n", "param = base.param\n", "geometry = base.submodels[\"storage\"].storage_geometry" ] }, { "cell_type": "markdown", "id": "820fed7d", "metadata": {}, "source": [ "## Schumann: a sensible packed bed\n", "\n", "The capsules are a sensible material, conductive enough to be lumped (Biot number 0.008),\n", "and the fluid enters at 70 °C into a bed at 32 °C. In Schumann's variables, `xi` is the\n", "number of transfer units from the inlet and `eta` the time since the first fluid reached\n", "the outlet, in units of the particles' time constant." ] }, { "cell_type": "code", "execution_count": 2, "id": "23fa3c8a", "metadata": { "execution": { "iopub.execute_input": "2026-10-23T17:49:05.875540Z", "iopub.status.busy": "2026-10-23T17:49:05.875274Z", "iopub.status.idle": "2026-10-23T17:49:06.489190Z", "shell.execute_reply": "2026-10-23T17:49:06.488405Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhgAAAFzCAYAAAB8X3AUAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjExLjIsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvgI3uAAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAX9dJREFUeJzt3XlcFPX/B/DXLvd9ityXXKJ4oOB95Z2aB95HpZmmZV9TK/ml5kGipqVldpmWmlYeeWWemYgHauKVIXKIHCIoyM0Cu/P7g9zaAF2WxeF4PR8PHjEzn5l97Ujsm5nPfD4SQRAEEBEREWmRVOwARERE1PCwwCAiIiKtY4FBREREWscCg4iIiLSOBQYRERFpHQsMIiIi0joWGERERKR1umIHEINCoUBaWhrMzMwgkUjEjkNERFRvCIKAvLw8ODo6Qiqt+jpFoyww0tLS4OLiInYMIiKieis5ORnOzs5Vbm+UBYaZmRmA8pNjbm4uchoiIqL6Izc3Fy4uLsrP0qo0ygLj8W0Rc3NzFhhEREQaeFoXA3byJCIiIq1jgUFERERaxwKDiIiItK5R9sFQh0KhQElJidgx6BnQ09ODjo6O2DGIiBoUFhiVKCkpQWJiIhQKhdhR6BmxtLSEvb09x0UhItISFhj/IQgC7t27Bx0dHbi4uDxxEBGq/wRBQGFhITIyMgAADg4OIiciImoY6kSBkZGRgYSEBPj7+6v92GhaWhoyMjLg5eUFU1NTrWUpKytDYWEhHB0dYWxsrLXjUt1lZGQEoPzn0M7OjrdLiIi0QNQ/z6OjozF27Fi0aNECnTp1wuXLl5+6j0wmw5gxY9CsWTOMGjUKTZs2xeeff661THK5HACgr6+vtWNS3fe4mCwtLRU5CRFRwyB6gTF06FBERUWpvc+yZcsQGRmJuLg43L59G1u2bMHrr7+OixcvajUb78U3Lvz3JiLSLlFvkUyZMgUAkJKSovY+33zzDaZNmwYnJycAQEhICPz9/bFp0yYEBQXVSk4iIiKqnjrRB0NdaWlpSE9Pr1BIdOjQAdHR0VXuJ5PJIJPJlMu5ubm1lpEqGjZsGJydnbF+/fpKl4kaEkEQoJCXQlaSh8KifKSkJSHzQQZy83OQV5CD/MJ8FBUXQlYiQ5m8FH4t3KGQKCCXlyL+9l1kpGdBLi9DmVwOeVkZyhRlkJfJIVfI0aaLJ3R0JVAICsT9mYp7SQ+hEP5+XQiAQij/L4DAHp4wMNYDACT+dR8pcQ/KtwjlGQVBUO7Xtoc7TCwMIUBAUkwmkm5m/tMWf7f9u33rnm4wtzEEAKTEZiHhWmaV5yKgmzOsmpoAAO4lPMLty/erbOvf2RG2juVzW9xPykHMxXtVtvULdoSda3l/vQepebh5LrXKtj7t7GHvYQEAyEovwI3Iqv+gbdbaDk7eVgCAnMxCXPk9ucq2Hi1s4drcBgCQl12MyyeSqmzr5mcN95ZNAACFuTJcOJJYZVtnbys0a2MHACguLEXULwlVtnXwtIBPO3sAQKlMjrP74yrmbGmLIb2H48Xn/6/K49SWelVgZGVlAQCsra1V1tvY2ODhw4dV7hceHo4lS5bUajaxhYWFYeHChSrrLCws8OjRI5V1W7ZsQVhYGJKSkuDt7Y3w8HAMGTLkGSYlqpsUcjnS7sUjLvEaYhNikJR2B/fupyEz6wHy8vPQfrArZIpiyIQSnP85CelxuSiRyVEmU6BUpoBcpoCiTIBQJqDFNy0gkZbfdru7/i5yL1X9R03zz5tDx6i8Y3HK7hQ8Ov2oyrbX/NOhZ1leNKRdSEPWiawq295pVwgDYwMAQPr1dDw49KDKtvfalcHQqrxoyLidgYyDGVW2zWwHGDct77P04O4DpB9Or7JtRisFTJ3KO+FnpWYh7XBalW3Tfctg/nfRkJ2ejdTDVRcN99xLYelpCQDIycxB8uGqC4FUBxmsvcs/M/Ky85B0uOpCIMWmGDbNCwAABbkFSHxCIZBsVogmrYoAAEWFRYh/Qtsk/Xzcalv+R66sRIbbR+9U2faOJA/xQeV9wUrKShD7hLZWPa1wp2N5v8EyWRliKmmbbFYIj/SzVR6jNtWrAkNPr/x/rH9fjQCAoqIi5bbKhIaGYs6cOcrlxzPBNTQdOnTA+fPnq9x++PBhvPLKK/jmm28wZMgQbN68GSNGjMD58+fRrl27Z5iU6NmSl5Xh+p+ncSLyF/wV9xfupqbgUUEuWgyzQzaK8Egqxx9rkpD/V0HlB5ACD17QgUTn76IhowC5t6toC0AoEyDRL28rNZBCqi+BRFcCqa4EUl1p+X+lgEQqgWspYCAFJJBAYWMAnWbGkEglkEolkEig/F4qlaB1mRGMi/UhhQRxzlZIDZIAEgkkACCB8nuJRIIOEksYy/QBSJDoUoa7XXSUfY0et3+83F7fFmYl5QVGkrMEd3rq/dNW8rh9eZ42Jg6wKDWGBECyswES+lbdIb6VuSOsSsuvYKQ5GOF236o/clpYOcG21BSABOl2Jrj1hLZ+tg6wLy0vRjJtzHCzX9XdCX3sm8K+1BIAkGVZgBv9qmyKZo52cC4tv4LxyKwQ1v2qHgvJw9UWrqW2AIA8o2JY9pNX2dbV0xoepX9fwdCXwbxfeQFRWc8vJy9reJX+fQVDtxTm/avueO7gYQGf0r+vYEAO0/6yCm083GwRYN+1ymPUpnpVYDwelyItTbUKTktLg5ubW5X7GRgYwMDAoLbj1XmrV6/GkCFD8OKLLwIA5syZgx07dmDt2rXYunVrpfsIgoC1a9fi888/R3JyMgICAvDRRx+ha9euyu0ff/wxNmzYgJSUFHh6emLOnDmYOnWq2rmOHDmC//u//0NMTAwcHR3x6quvYu7cuXxclKpNIZcjOe0izsccRmzGDfy8PQrx1x4gP0MGRYmg2lgKFI40+ftKgxRSk/KfN6mBFHomOjAw0YGRqR6MTfRgYmKI3jJnGOoaw0DXCPcGeqCkmwBzM0tYmFvBytIWVpa2MDe1hKmpJVyd3WBsZA4TI3OYjDWGgZ4upFI1OhJPq8abnV6NtqSZ0Gq0facabec8vYnSm9Vo+3o12j4Ddb7ASEhIQH5+Plq1agVjY2N07NgRBw4cwIQJEwAAhYWFOH78eIXbA9oiCAKKSquuTGuTkZ5OtZ5uuHr1KszMzGBoaIjg4GCsWLECAQEBAMrfx/nz57F8+XKVfXr37o2dO3dWecxFixbhiy++wMaNG9GrVy/cunULmzZtUhYYy5Ytw48//ojt27ejZcuWuHjxIkaNGgVTU1OMHTv2qZnz8/MxfPhwrFq1ClOmTEFGRga++uorXL9+HW3atFH7vVPjdPnyCXy38zOcj76CxIRMPMoohs/HvsrbE8lZechNKS5vrAMY2+jDzEYf5pZGsLI2R+/S1rC39oSdpSdM19jCxd4H7k4eMDPSU68gIKIqiVpgZGZmIj4+HpmZ5Z2Ebt68CUNDQzg7O8PZ2RkAsHz5cpw/fx43btwAUN7XoF+/fli4cCE6deqETz75BNbW1pg2rTqlv/qKSuXwX3SkVo79NDeX9oexvnr/RPb29ti0aRMGDBiAR48e4b333kPXrl1x48YNuLi4IC8vDwUFBWjSpInKfnZ2dkhPr/weamFhIdasWYOPP/4YQ4cOBQAEBQUpO9kWFxdj1apV2Lt3L4KDgwEAPXr0wBtvvIFvvvlGrQIjIyMDRUVFeOGFF2BsbAx3d/cKRRDRY3l56djw7WL8+PMh3PrzPgozKs4XJE8rho+tEZrIzdGqe1sYdLNB64Bu6NllKPxcXWCgyytjRM+CqAXGhQsXsGzZMgDl/Qe2bNmCLVu2YOrUqcpL7M2aNUNRUZFyn169euHEiRNYv349IiIiEBAQgE2bNqk9AmhD9e9bElZWVvj222/h4eGBr7/+GkuXLq1yP0EQqrxKEhMTg6KiInTv3r3S7X/99RcKCgowYMAA5bEef3l6eqqV28PDAwMHDkSnTp0wduxY9OrVC71791aOrkn02+8/4kr6L4jKv45buqVIO5aBzJN/P7kgBazcjODoZgkPj2YIatcHg/q/jFYeLtDT4TD/RGIStcAYNGgQBg0a9MQ2oaEVb4J17969yg89bTPS08HNpf2fyWtV9tqa0tfXh4+PD27fvg0AMDMzg4mJifJq0WOZmZlo2rRppcd4/BhbVQXI48ngrl27Bn9/f41ySiQS/PLLL4iIiMCxY8fw3nvvYfr06Th27JjGx6T678+b5/H+qtk4dvwqclOL4fo/V5i3NQcggWeAKZrk6aJ16/YY8cJbGNClG0wN6vzdXqJGh/9XPoVEIlH7NkVdUlJSgtjYWAQGBgIofx8dO3bEyZMn8eab//Qa+u2339C5c+dKj9G8eXMYGRkhIiICfn5+VW7/9ddfa1QMSCQS9OjRAz169EBYWBiCgoKwadMmrF69WuNjUv0jLyvDJ1++i7Vfbcbd69n4ezgHSHQksEwuRR+/ALTyHIvnhvSBywecJ4ioruM1xAZi1KhRiIiIQF5eHhITE/Hyyy/j0aNHKn1T5s2bh4MHD2LLli3Izs7Gxx9/jOjoaMyePbvSYxobG+Ott97CggULcODAAeTl5eHSpUuYMWOGcvs777yDpUuXYseOHcjNzUVSUhI+//xzhIWFqZX7/PnzmD59Oq5evYri4mJcuXIFKSkpaNasWY3PCdUPpaVF+HzHTFg4m2DOGx/h7rXy4sK6mTF6j2+Hr/buR8TWJKx5bTte6vcCXKxZXBDVB/XvT3Oq1OzZs7F06VJERUXB2NgYwcHBiIqKgq+vr7LNgAED8M0332DZsmV49dVX4e3tjd27dz9xDIywsDBYWlpi9uzZSEtLQ+vWrVWuLCxevBh2dnYICwvDyy+/DEdHRwwePFjtp3qCgoJw7do1vPTSS4iNjUXTpk0xbdo0TJ/OZ/AaOpksHz/99h62pJ3APR1AbiiFjpEUvh2a4oWh0zFj4jy42pqIHZOINCQRHt9ob0Ryc3NhYWGBnJycCp1Di4uLkZiYCA8PDxgaGoqUkJ41/rs/O2WlpZizYDS+/eEwnEM9oGOkA+syBZzv2KNP1/cw7rkeMKxB/yMiql1P+gz9N17BIKJn5rtty/G/+cuQk1o+NkXx8Yd4vnt3DOyyBAMnN+fYE0QNCAsMIqp1SXdjMPylfoj+e/IoHWMpgvp6463/bcPILu1YWBA1QOzkSUS1auXaN+DdsqWyuPDu0gTLv9mOkzv/xOhu7VlcEDVQvIJBRLWipKQAq/aMx6otJ1GaJ4eJowEGjQvB6v/7mk+CEDUCLDCISOtS0i7hzV+n4rauHE6TnWDYNB/z5v+M8d3bVWt+HSKqv3iLhIi0avmamQh+uT9idcpgKVegD/pg37ZbmNCjPYsLokaEVzCISCsUCgWGT+qA/dsvAQCcPU0wpN/HeOul4ZwXhKgRYoFBRDVWUiJD54G++OO3JABA8x5N8e5LBxHSqZXIyYhILCwwiKhGHj3KRLvnmiMh+iEgAbqFBGDDx6fQ0tlK7GhEJCJet6Q6YenSpdiwYUOVy1Q3ZWdnoHnHZkiIfgiJngQDJvfED19HsbggIl7BaChKSkrw/fff49SpUygtLUVQUBCmT58OIyMjlXa3b9/Gp59+iqSkJHh7e+Ott96Ck5NTtdto2+XLl+Hs7FzlMtU9JSUFGB/WA+mxeZAaSjFs+khs/GArrEz0xY5GRHUAr2A0EJ07d0ZkZCR69uyJAQMGYPPmzejSpQuKioqUbW7fvo3g4GA8fPgQY8eORUxMDIKDg5GRkVGtNkSlJYV4c0dfpLTShcerThgxcwK2rvqexQURKfEKRgNx+PBh2NraKpf79OkDJycnHDp0CCEhIQCAZcuWwdvbG9u2bYNEIsHIkSPh5eWFNWvWYOXKlWq3qcy1a9ewceNGJCcnIyAgAHPnzoWFhYVy+5UrV7Bx40akpKTA09MTM2bMgLe3t9rv7969e1i/fj1iYmLg6OiIyZMnIzAwsLqnibSgpESG/23uhzOGeTBQCOgZMB4rXl0OY33+OiGif/AKxtMIAlBSIM5XNSa6/XdxAQBWVlbQ1dVFfn6+ct2RI0cwdOhQ5VgEenp6GDJkCA4fPlytNv+1f/9+BAcHo7S0FOPGjYOenh7Gjx+v3H7o0CH07t0bdnZ2ePHFF2FgYIDAwEBcvXpVrfdWVlaGbt264datW5gwYQJatGiBN954AzExMWrtT9ojCAK6D2mOzUsvQZEpQ6f8flg+9QOYGLC4ICJV/K3wNKWFwHJHcV77/9IAfRONdv3kk08glUrRu3dvAEBhYSEyMjLg4uKi0s7Z2Rl37txRu81/KRQKzJgxA2+++SZWrVqlXJ+bmwug/ANpxowZ+PDDDzFlyhQAwIgRI5CTk4OwsDDs3Lnzqe/lzp07iI+Px7lz59CkSRMAwKuvvori4uKnnwjSqtGTOyPqaCIgAVxjHBH2yUqYGeqJHYuI6iAWGA3Q4cOH8d577+Gzzz5TdpQsKSkBgAqdPo2NjZXb1GnzXzExMUhLS8PYsWNV1pubmwMAbt26hbt372LLli3Ys2cPBEGAIAhISkqCoOYVGhcXFzg5OWHatGmYOXMmOnfuDBMTE5iYaFZ8kWYWh0/Bru/OAwA6D2uOrz88BAtjFhdEVDkWGE+jZ1x+JUGs166m3377DSNGjEBYWBimTZumXG9qagpdXV1kZWWptH/48CGsrKzUbvNfj2/BVLU9JycHAPDiiy/C0VH1StB/C5mqGBgYICoqCuvXr0doaChu3ryJoUOH4ssvv1QWMlS7jh3fhmVLvgMA+Pe0xxfrzsDZihOWEVHVWGA8jUSi8W2KZ+3kyZMYMmQIFi5ciHfffVdlm66uLlq0aIHo6GiV9dHR0WjdurXabf7L3d0dAHDz5k14eHhUud3MzAwDBgzQ5G0BAJycnBAeHo7w8HCkpKSga9euWLt2LRYtWqTxMUk9aWkJGD7hFShkCjTxNcGHq39HgAvHuSCiJ2Mnzwbi1KlTGDx4MBYsWIDQ0NBK27z00kv46aefkJiYCKB8rIkjR47gpZdeqlabf7Ozs8OQIUOwePFiPHz4EAAgl8vx/fffAwCaNm2KoUOHYtGiRUhJSVHud/v2bezZs0et93br1i38/PPPymV7e3uYm5tDLpertT9pTlAo0GdMZxRklMDARg9vvPc5nm/nK3YsIqoHWGA0EMOHD4dUKsWpU6cwYMAA5de2bduUbWbNmoWBAwciICAAHTp0QNeuXTF9+nSV/hPqtPmvTZs2wdzcHJ6enujSpQvc3NyQnJys3L5582Z4e3vD19cXnTt3RvPmzTFs2DBYW1ur9d5sbGywdetWODo6omfPnvDw8ICRkRFmzZqlwZmi6thx/B0oQixh2twEIdMm470JE8WORET1hERQt6ddA5KbmwsLCwvk5ORUuIdfXFyMxMREeHh4wNDQUKSE1Xfs2LFK/6L38vKCl5eXyrr4+HjcvXsXXl5eFZ4YqU6b/7p9+zbS0tLQokWLCo/NAsDdu3cRHx8PR0dHeHt7Qyr9p76Njo6GoaEhmjdvXukyAKSnpyM2NhZNmzaFr692/4qur//utSkhKQKjf5sJmVSCztluWDltDyyNOZAWUWP3pM/Qf2OB0UAKDKoZ/rurys15iOeXBCOntTGaF0rxdv+jCPJsKnYsIqoD1C0weIuEiCroPyYYZz5OwKPd6ejktIzFBRFVGwsMIlKx8dslOH8kAQAQ6NwRswYNFjkREdVHLDCISCk3Jwuz310OAPDu0gTrlv4AXR3+miCi6uNvDiJSCpnSHQUZJdC31MXrc7+Bu239GAOGiOoeFhhEBAA48MvXOL73TwBA35Hd8cZQ3hohIs2xwCAiyMvKMHnWbEABuLS1xKrFu6EjlYgdi4jqMRYYRIRfIsNh1N8ahg4GmPL6Svg7WYodiYjqOc5FQtTIFRU9wrr4H2HZ0RIDfZtj/ouviB2JiBoAXsGgOuHAgQP4/fffq1ym2vPp/unI0JXArlSBF7p9DEM9HbEjEVEDwCsYDdC9e/fw448/wtvbG4MGDVLZVlZWhqNHjyIpKQne3t7o3bs3JBJJtdto2zfffANnZ2f07Nmz0mWqHafP7MWCmT/Dbrgderd5HoMC/cWOREQNBAuMBkahUGDixImIiorCgAEDVAqM/Px89OnTB5mZmejatSs++OADBAQEYP/+/dDT01O7DTUck//3KkqzSlF2PgfTlqyq9UKSiBoPFhgNzPLly2Fubl7pX/4rVqxAamoqrl27BisrKyQnJ6NFixb46quv8Prrr6vdpjIlJSU4fvw4kpOTERAQgM6dO6tsl8lkOHbsGFJSUuDp6YnnnnsOurrV+/E7ffo0YmJi4OjoiOeeew5GRkbV2p9U7dyzDvF/PAAkwKiJM9mxk4i0in0wGpAzZ87gq6++wtdff13p9p9++gljxoyBlZUVAMDFxQWDBw/GTz/9VK02/xUfH48WLVpgzpw5iIqKwttvv43p06crtyckJKBly5ZYvnw5Ll++jNDQULRv3x5ZWVlqv7eQkBBMnDgRFy9exNdff41OnTohLS1N7f2potkLFgEAPINt8P4bS0VOQ0QNDa9gPIUgCCgqKxLltY10jdS+ZJ2dnY0JEyZg48aNlU6VXlJSgri4OPj5+ams9/Pzw7Fjx9RuU5mJEyeiWbNm2L9/P/T1y6fzPnv2rHL7yy+/jKFDh2L16tUAys/pwIEDsXjxYnzyySdPfW93797Fnj17kJCQAA8PDwDlRQsv52tuw9f/h7S/ciHRkWDcS/Nha2ogdiQiamBYYDxFUVkROmzvIMprR42PgrGesVptX3nlFQwbNgz9+vWrdHtBQQEEQYClpaXKeisrK+Tl5and5r/u3r2L8+fP4+TJk8riAoDyFklqaipOnz6NTp06Yf369RAEAYIgwMjICGfOnFHrvZmZmcHAwAAHDhzA9OnTYWBgAE9PT7X2pYoUCgUWhq8DAPh1bYp3Jv9P5ERE1BBpfIvkjz/+wKxZs1Q6EX733XfIz8/XSjBS34kTJ/DLL7/AwcEBa9euxdq1a5GYmIi4uDisXbsWubm5MDYuL1T+Wyj8e5s6bf7r3r17AAB3d/dKt6ekpAAAHjx4gLi4OMTHxyMhIQFubm4YOnSoWu/PysoKP/30EzZv3gxra2v06dMH3333HQRBUGt/UvXZxnnISiyE1ECKl6ethLkhO+8SkfZpdAXjl19+wahRozBs2DAcOnRIuT4pKQlr167FggULtBZQbEa6RogaHyXaa6vDzs4OM2bMUH7YA0BRURFKSkpw584dlJWVwdzcHK6uroiPj1fZNz4+Ht7e3gAAAwODp7b5LxsbGwBAenp6pUXG4+2TJk2q0SOnL7zwAl544QVkZmbi4MGDeOONNyCXyzFlyhSNj9kYCYKASJ3f4bnIE6539PHmqAliRyKihkrQQGBgoLBz505BKP8TUrk+JiZGcHd31+SQz1ROTo4AQMjJyamwraioSLh586ZQVFQkQjLtGTRokBASEqKybtasWYKPj49QXFwsCIIgZGdnC7a2tsLy5cur1ea//Pz8hEmTJqmsi4uLU37v7+8vDBkyRJDL5cp1ZWVlws2bN5XLQ4cOFV5//fVKl+/fvy/cv39f5fg9evQQ5s6d++STUA0N5d/9ac5e3iS0/LalELiphbDp6K9ixyGieuhJn6H/ptEVjL/++gvPP/88AKh0tHNyckJqamq1j5eWloaMjAx4eXnB1NT0qe0FQUBycjKys7Ph4uICa2vrar9mY7RgwQIcPHgQvXr1Qt++fbF37144OTlh1qxZ1WrzX99++y0GDBiAPn36oEePHrh16xaKi4uxa9cuAMDWrVsxYMAAdOnSBQMHDkR2djaOHz+OGTNmoHnz5k/NnZOTg4EDB6Jbt27w9fVFTEwMLl68iDVr1tT8pDQyH51aB1gBbfItMa5nX7HjEFEDplEfDEtLSyQnJwNQLTDOnj0LZ2dntY8jk8kwZswYNGvWDKNGjULTpk3x+eefP3GfGzduoFWrVmjXrh1efPFFODs7Y8yYMSguLtbkrTRYQ4cOxeDBqtNt29nZITo6GhMmTEBxcTFmzZqFc+fOqRR16rT5rw4dOuDWrVsYOnQoCgsLMXjwYPz444/K7YGBgYiNjcWUKVNQUFAAd3d37N69GzNnzlS2GTJkCHr16lXpsre3N65cuYIePXogNzcXbdu2RWxsLNq1a1fj89SYbPxuCXa/fQ33d9xDsOdbHBKciGqXJpdH5s6dK/Tp00dITk4WpFKpUFxcLBw8eFBwcXERFi5cqPZx3nvvPcHR0VFISUkRBEEQdu3aJUgkEuHChQtV7tOlSxehX79+gkwmEwRBEBISEgRzc3Nh9erVar9uY7hFQtXTGP7dHf3NBQCCT1c7Ib+4VOw4RFRPqXuLRKMrGGFhYbC2toaLiwsUCgVMTU0xePBgdO7cuVodPL/55htMnToVTk5OAMoHU/L398emTZuq3Of+/fsIDg5WPhLp4eEBV1dX3L9/X5O3QtQo7N73GdJu5gJSYNyEd2FiwCfUiah2afRbxtDQED/++COWLVuGP/74AwqFAoGBgWrdT38sLS0N6enpCAoKUlnfoUMHREdHV7nf4sWL8e6778LLywtubm44evQo8vLyMGPGjCr3kclkkMlkyuXc3Fy1cxI1BAtWLgMAuAda462Xqu5PQ0SkLRoVGA4ODrh37x58fHzg4+Oj0Qs/Hib6vx00bWxs8PDhwyr369+/P/bs2YM5c+bAyckJSUlJWLBgQZXjMABAeHg4lixZolFOovru+o1IxFwov8I3eNh4WBhx3Asiqn0a3SIpKSlBdnZ2jV748cyc/76yAJSP31DVrJ3C30NMy+VypKWl4dq1a7h69SrWrFmD8PDwKl8rNDQUOTk5yq/HHVSJGoP/LZwGyAEbL2O8O2OF2HGIqJHQqMAYP348PvzwwxqNpOji4gKpVFphwqq0tDS4ublVuk9aWhouXbqkHC4aKB9B8oUXXsDevXurfC0DAwOYm5urfBE1BgUF2Thz8jYAoHvfrnC2NhE5ERE1FhrdIomNjcXRo0fxww8/wM/PT2UOCgBP/LB/zNjYGB07dsSBAwcwYUL5aIKFhYU4fvw4Fi5cqGyXkJCA/Px8tGrVCtbW1pBKpRXG2khLS6t0gi+ixm5vZDg8ljRD6elsvDu78ll2iYhqg0YFhre3d5VDR1dHWFgY+vXrh4ULF6JTp0745JNPYG1tjWnTpinbLF++HOfPn8eNGzdgZGSEKVOm4P/+7/8AAJ6enjh69CgOHTqE/fv31zgPUUMiKBT4MeUw9Cz10K9XR3TwcRU7EhE1IhoVGOvXr9fKi/fq1QsnTpzA+vXrERERgYCAAGzatEnlFkazZs1QVPTPdOmff/45goODcfjwYTx8+BBubm44ffo0unTpopVMRA3F2egtSNQVYKAQ0Lv9e2LHIaJGRiLUpCNFPZWbmwsLCwvk5ORU6I9RXFyMxMREeHh4wNDQUKSE4khOToahoSGaNGkidpRnriH+u3sHN8G9giI8N8gTe8KvQFdH48mTiYiUnvQZ+m8aXcHo2rXrE7dHRkZqcljSguLiYmRkZMDJyQk6OtUbCnry5Mlo3749Vqzgkwb1XcytC4i7/ACQA36TB7K4IKJnTqMCo0+fPirLCoUCt2/fxs6dO1X6T9Czk5SUhBkzZuC3335D06ZNkZ2djSFDhmDp0qVo1qyZ2PHoGXt76QxADlh5GmHeq4vFjkNEjZBGBcbixYsrXT9w4ED8+uuvNclDGho7dixMTU1x7949WFlZobS0FHv27MG5c+cqFBhFRUXIyMhQPipcmezsbCgUCtjY2Kisv3PnDszNzVUGSEtLS4NUKoW9vT0A1Vst2dnZKCkpQdOmTZXtU1JSKn1c+P79+8rxVWxsbCq9VfPfY1eWUZ02DVlpaQlOHL0OAOjUrS3sLIxETkREjZI2J0DJzc0VmjRpos1D1oqGNtlZUVGRIJFIhK1btz6xXWFhoTBt2jTBwMBAcHBwEGxtbYUtW7Yot/fu3VsYPny4EBQUJDg5OQmGhoZC3759hYKCAmWb1q1bCx9++KHKcceMGSO88sorKsd5/vnnhRYtWghOTk6Cnp6eMGjQIOH8+fNC8+bNletmzZqlcpz3339f8PX1FXx9fQUrKyvB09NTOHXqlEobdTKq06ayc1jf/t2rsnLtGwIAQcdYKhy+8IfYcYioganVyc6qcv36dZXp2xuSgoKCKr/+O1X8k9r++4mYJ7WtDkNDQzg6OmLbtm0Vxgj5tzfeeAMnTpzA1atXkZaWhtu3b+PevXsqbQ4fPox169YhJSUFd+7cwZ9//okvvviiWnkA4Pfff8fmzZuRkpKCa9eu4ejRoxg4cCC2b9+OlJQUREVF4bPPPkNUVJRyn8WLFyMmJgYxMTF4+PAhZsyYgXHjxlU4Z+pk1Nb7qI82fPs9AKB5sD36tmsrchoiaqw0ukUyderUCuuys7Nx5MgRvPbaazUOVReZmppWue3555/HL7/8oly2s7NDYWFhpW179OiB33//Xbns7u6OBw8eVGgnVPPhnh07duCll16Cs7MzfHx80KVLF4wcORLPP/88AODRo0f47rvvsH37dvj6+gIALC0t8c4776gcJyQkBJ06dQIANG3aFH379n3i5HNVCQkJUU5k5+fnhxYtWiAoKAht2rQBALRt2xaenp6Ijo5Ghw4dVPbNysrCgwcPMGDAAMyfPx83b95Eu3btqpVRW++jvvkr5gKSrpXfZho+Yjqk0oZZ8BNR3afRFYz8/PwKX3Z2dvjqq6+watUqbWckNXTr1g3x8fH4448/MGvWLGRkZGDQoEF4++23AQC3b9+GXC5HYGDgE4/j5OSksmxqaoq8vLxq53F0dFRZNjExqXRdfn6+cvn48ePw8fGBm5sbBg4ciJEjR0KhUFS4KqNORm29j/rmyNUv4TTZCe49rDF78rtixyGiRkyjKxhjx47FsGHDKt22d+/eKrfVZ//+IPyv/z4OmpGRUWXb/3aqvHPnTo1y/ZtEIkFgYCACAwPxxhtvYP78+VizZg0++OAD5dwtVV1Zqc5r/FdZWVmNjgkAcrkco0aNwvz58zFv3jzo6OigtLQURkZGUCgUNT5+Y/F73iVYdbPCgPZdYW1qIHYcImrENCowhg8fXuUl/Cdtq89MTNSfJKq22lZFEASUlZVVmIXW1bV8aGi5XI7mzZvDysoKhw4dQqtWrZRtysrKoKur/o9BkyZNVPptKBQKXLt2Dd27d6/Re7h37x4ePXqEkJAQZcF26tQpyOXyGh23Mfkz9hDi9RTQFQT0aveW2HGIqJHTqMCoSmZmJiwsLLR5SFJDaWkp/Pz8MHXqVAQHB8PGxgZXr17FkiVLMH78eBgZlT+mGB4ejrfeegs6Ojro2bMnbt26hT179mDPnj1qv9bzzz+PZcuWoVu3bnBwcMAXX3yBxMTEGhcYjo6OcHNzw6JFizB37lwkJCRg7ty5DbbTcG343+J5eIBsdG3TFP0mBYkdh4gauWoVGAMGDKj0e6D8L9lbt27V+IOGqk9fXx+RkZH44osvsHLlSjx8+BAODg5YsmQJpkyZomw3ffp0ODg44Msvv8T27dvRunVrfPzxx8rtrq6usLOzUzn247EtHnvjjTdQUFCAlStXwsTEBMOHD4elpaVKJ9jKjuPm5lZhxlsPDw/l+BRSqRSHDh3CwoULMWXKFDg4OODzzz/HokWLYGZmVq2M6rRpaPLzH+HU/liUFchh6tSNI3cSkeiqNRfJvHnzAABr1qzB3LlzVbbp6enB3d0d48ePV/lAqIs4Fwn9V33/d38v7CUsX7gF+pa6uHQlCQFujk/fiYhIA7UyF8nq1asBALa2tpg/f37NEhKR1ny/+wAAwD/IkcUFEdUJGl1HZXFBVHfcvRuDu9cfAQCGvjBZ3DBERH/TuJPnw4cPceLECdy9e7fCY4osQIienSVrZkOQCzB1NMCsl/j/HhHVDRoVGOfOncOgQYNgbW2N+Ph4tGjRAnFxcZDJZGjdujULDKJn6NCxswCAgPZusDGrf/1HiKhh0ugWydtvv43Q0FDExcUBAG7cuIHU1FQMHDgQ/fv312pAIqpaRmYsikwkkOhKMHrkLLHjEBEpaXQF4+rVq8q5N6RSKWQyGWxsbLBhwwZ06dIFK1eu1GpIMTTEwcKoavX133vPmbVwecMVro/kmDq6Yc4DRET1k8ZzkTweUMvOzg5JSUkAAAsLC2RlZWkvnQgejyJZUlIichJ6lh4Pof7f0VDrupMPzgEA/A39YWqg1XHziIhqpMa/kZ577jnMmzcP06ZNw/fff4+2bev39NC6urowNjZGZmYm9PT0KswdQg2LIAgoLCxERkYGLC0tK8wrU5ddvf4boh8WwLCpProHzBA7DhGRCo0KjH+P/rhq1Sq8/PLLGDt2LPz8/LB582athRODRCKBg4MDEhMTlVdmqOGztLSsd6N9/t+Kt3F7+200626Dfid6ih2HiEhFtUbybCjUGYVMoVDwNkkjoaenV6+uXDxm7W6M7KQiDJgUjF+3RIkdh4gaiVoZyfMxBwcHlRk1GyKpVFovh4ymxuHCxV+RnVQESICXJ7wndhwiogo06mBQUlKC7OxsbWchIjWFf7YEAGDrZYIRfQaLnIaIqCKNCozx48fjww8/rLeP9hHVd5GRNwAAbdr5QI8zpxJRHaTRLZLY2FgcPXoUP/zwA/z8/KCvr6+yfe/evdrIRkSVuHHjNB4kFAAAJo6d+5TWRETi0KjA8Pb2hre3t7azEJEawj5dAAiAlYcRxg4aK3YcIqJKaVRgrF+/Xts5iEhdATlwecMFQXCGgW79e/qFiBqHGg+0VVpaWu9GPySqr7Kz7yDGtAwW7S0wweMdseMQEVVJo95hpaWlWLx4MVxcXFT6X8ycOVM5ARoRad/PkZ9CLpHATSbg+U5DxI5DRFQljQqM8PBw7Ny5E2vWrFFZ3717dyxdulQrwYioorWf7cD9n+/DPssOhnq8PUJEdZdGBca3336LHTt2YPTo0Srre/bsiYMHD2olGBGpevToPv48mY7MfZlwN+8qdhwioifSqMBITU1VPkUikUiU63V1dVFUVKSdZESk4qMv3oOiRICBtS5enxQqdhwioifSqMDw8fHBmTNnAKgWGNu2bUOrVq20k4yIVPx8+FcAgFdAE1iZGYmchojoyTR6iiQ0NBQTJ07EggULAAC7d+/G4cOHsXnzZuzcuVOrAYkIkJeVIfZKBgCge9fnRE5DRPR0GhUY48ePh0QiQVhYGBQKBUaOHAkfHx9s3boVw4cP13ZGokZv649rUJJTBqmhFLOmhokdh4joqWo8XXtRUREUCgVMTEy0lanWqTvVLFFd0X2oH07vvwW3QCvc+SNL7DhE1IjV6nTtj2VnZ+P27dsAyocPt7KyqsnhiKgK9+VZkBpLEdS+rdhRiIjUolEnz/z8fEydOhV2dnbo0KEDOnToADs7O0ybNg0FBQXazkjUqN1JuQz9UU3RYp0fZs1YLnYcIiK1aFRgzJgxA2fOnMHPP/+M9PR0pKen4+eff0ZERARmzpyp7YxEjdrBqI0AAG+5Drq1DhY5DRGRejS6RfLzzz/j/PnzaNmypXLd4MGD4ebmhk6dOuG7777TWkCixi4i/hwEOwHNpM1UHgsnIqrLNLqCYWlpCXt7+wrrHRwc2A+DSItycx9g98KriJ0XCw+z7mLHISJSm0YFxtChQ7Fw4ULIZDLlOplMhgULFmDo0KFaC0fU2H3y9WIoSgRIygRMeIG3H4mo/tDoFklcXByOHj2KXbt2oVWrVhAEAdevX8eDBw/Qr18/DBs2TNl27969WopK1PjsPfoLAMC9hTXMjQ1FTkNEpD6NCgxvb2/lXCSP+fv7axwiLS0NGRkZ8PLygqmpqVr7KBQKxMbGwtjYGK6urhq/NlFdFnMtDQDQMZidO4moftGowFi/fr1WXlwmk+HFF1/E/v374ezsjLS0NKxevRozZsx44n779u3D66+/Dh0dHRgbG8PNzQ3ff/89bGxstJKLqC44cWo3CtJLAB1g+suLxI5DRFQtNRpoq6aWLVuGyMhIxMXFwcnJCbt378aoUaPQvn17BAUFVbpPZGQkQkJCsGHDBkybNg0AcOLECdy7d48FBjUoG7Z8BACw9TRBp1btRU5DRFQ9Gg0VLggC9u/fj8jISGRnZ1fYvnHjRrWO4+DggGnTpmHJkiXKdS1btkS3bt3w+eefV7pPnz59IJVKcfTo0erGVuJQ4VQfOAdYIvVGDnqGtMTJXdfFjkNEBKCWhwp/++23sWHDBnTv3h2WlpYaBUxLS0N6enqFKxUdOnRAdHR0pfuUlJQgIiICa9asQUFBAeLi4uDo6IgmTZo88bVkMpnKEy+5ubkaZSZ6VmSyPBj3s4SNgy5GvDBJ7DhERNWmUYGxefNmnDx5Eh06dND4hbOyyidssra2VllvY2ODhw8fVrpPZmYmSktLcfPmTXh7e8PW1hYJCQno3r07tm3bVuFYj4WHh6tcJSGq645FbYFBKzO08DfB1DGzxY5DRFRtGo2DoVAo0Lx58xq9sJ6eHgCoXFkAymdnfbytqn1+/fVXXL58GdeuXUNiYiJu376Nt99+u8rXCg0NRU5OjvIrOTm5RtmJatuZhMMAAJ9SCxgZ6Iuchoio+jQqMEaMGKF2P4uquLi4QCqVIi0tTWV9Wloa3NzcKt3H1tYWxsbGCAkJUY4k2qRJE4wePRqnTp2q8rUMDAxgbm6u8kVUlx349Q/k38yHpzFnTyWi+kmjWyQffPAB/P39sX37djRrVnF+hB9++OGpxzA2NkbHjh1x4MABTJgwAQBQWFiI48ePY+HChcp2CQkJyM/PR6tWrSCVStG3b98KRUlqaipsbW01eStEdc7VG5H486d7AID2kRwZl4jqJ40KjNmzZ0Mul8Pd3R1mZmYav3hYWBj69euHhQsXolOnTvjkk09gbW2tfPwUAJYvX47z58/jxo0bAIAlS5aga9euWLZsGbp06YKoqCh8//332LFjh8Y5iOqSTzetAACYuxiiT+e+IqchItKMRgXGgQMHEBkZibZta3b5tlevXjhx4gTWr1+PiIgIBAQEYNOmTSq3MJo1a4aioiLlcuvWrREZGYmPPvoIJ06cgKurK06cOIHu3TkRFDUMp89fBAD4+Ntz9lQiqrc0GgfDyckJN2/ehIWFRW1kqnUcB4PqKoVcDiMbQ5TklOHtsKlY9d7XYkciIlKh7meoRp08+/fvr7XhwonoHz/u+xwlOWWQ6Evwxivvix2HiEhjGt0iuXv3LjZv3ozt27fDy8urwmVczqBKpJmtu78FANh7m8HV3lncMERENaBRgeHn5wc/Pz9tZyFq9G7ExAMAAlo1EzkJEVHNiDqbKhH9Q1acD7uZjtBPtcGLga+JHYeIqEY06oPxb6WlpdrIQdToHbu4DTJdKRydDDCy32Sx4xAR1YhGBUZpaSkWL14MFxcX6Ov/M4zxzJkzERcXp7VwRI3J+YQjAACvMjMY6Fc+XD4RUX2hUYERHh6OnTt3Ys2aNSrru3fvjqVLl2olGFFj892nJ5HyVQqs8iofKp+IqD7RqMD49ttvsWPHDowePVplfc+ePXHw4EGtBCNqTNLSEpByJQePzj5CkN9gseMQEdWYRgVGamoqvL29AUDlEVVdXV2VUTeJSD2fbA4D5IChjR5GDZwgdhwiohrTqMDw8fHBmTNnAKgWGNu2bUOrVq20k4yoETkRcRIA4O5nDR2dGve9JiISnUaPqYaGhmLixIlYsGABAGD37t04fPgwNm/ejJ07d2o1IFFjEPtXOgCgXQ3n9yEiqis0KjDGjx8PiUSCsLAwKBQKjBw5Ej4+Pti6dSuGDx+u7YxEDdqff11AbnIxAOClcXNETkNEpB0aFRh79+7FuHHjMG7cOBQVFUGhUMDExES5bdiwYdrMSNSgffZt+fTspo4G6NOpj8hpiIi0Q6MCY/jw4Xg8CauRkVGV24jo6e7kxMLQxRDNfJpwenYiajA0KjCqkpmZWW+ncCcSTXvAq5MXppgMETsJEZHWVKvAGDBgQKXfA4BCocCtW7fQvXt37SQjagTupv2JZD0JJIKAfkGTxI5DRKQ11SowWrZsCQA4cuSI8vvH9PT0EBISgvHjx2svHVEDt+vEV1CUKOAJKfzcmosdh4hIa6pVYKxevRoAYGtri/nz59dKIKLG5NtNu/FXZBpsR3gBr4qdhohIezQa0YfFBZF2JMY+gFAmwM89QOwoRERaxSEDiURy4cpp5KfJAABTx88VOQ0RkXaxwCASycbtHwMAzBwNENS6s8hpiIi0S6MCY+/evRptI6J/RP1xEQDg4WMrchIiIu3TqMB40nDgHCqcSD2JtzIBAJ2Cg0VOQkSkfVq9RcKBtojUc/7KGeSllve/mDbhLZHTEBFpHwfaIhJBxNVdsB9jD/30EgS26iZ2HCIireNAW0QiSJZdg+1AWzxXaid2FCKiWsGBtohEcEeSAQAIaNJB5CRERLWjxgNtlZaWai0MUWPwx40LiL6QjZIHJegbxCt+RNQwaVRglJaWYvHixXBxcYG+vr5y/cyZMxEXF6e1cEQN0cbtq5G6MRWpq+/AzbHl03cgIqqHNCowwsPDsXPnTqxZs0Zlfffu3bF06VKtBCNqqM5fjAIAePjYiJyEiKj2aFRgfPvtt9ixYwdGjx6tsr5nz544ePCgVoIRNVSJseX9L4LbtxM5CRFR7dGowEhNTYW3tzcAQCKRKNfr6uqiqKhIO8mIGqBrf/2BnORiAMCUsW+KnIaIqPZoVGD4+PjgzJkzAFQLjG3btqFVq1baSUbUAH29fTUgAMZ2+ugS1FfsOEREtaZaj6k+FhoaiokTJ2LBggUAgN27d+Pw4cPYvHkzdu7cqdWARA3J2QvnAADuvux/QUQNm0YFxvjx4yGRSBAWFgaFQoGRI0fCx8cHW7du5VwkRE+QcKu8/0X7tm1FTkJEVLs0KjAAYNy4cRg3bhyKioqgUChgYmKizVxEDU5q5h04zHGFWWwhpk6cLXYcIqJapXGB8ZiRkZE2chA1eMcubIWOhR6atzFHN/a/IKIGTu0Co0+fPmof9Pjx4xqFIWrIbqSfA3QBT4WV2FGIiGqd2gVG+/btazMHUYO3b8dF5JpLYdeX/y8RUcOndoGxYsUK5fd79+7FsGHDKm23d+/emmYianDuZaQg/tRDCGUCWk/tJXYcIqJap9E4GE96UoRPkRBV9NW2DyGUCdAz18HQfpzgjIgaPo0KjKpkZmbCwsJCm4ckahBOninvl+TkZQEd3Rr3rSYiqvOq9ZtuwIABlX4PAAqFArdu3UL37t21k4yoAYm9lQwAaNnCS+QkRETPRrUKjJYty6eWPnLkiPL7x/T09BASEoLx43n5l+jf8grykJFQAAAY0i9E5DRERM9GtQqM1atXAwBsbW0xf/58rYVITk7G/fv34ePjA3Nzc7X3y8nJwfXr1+Hk5AQPDw+t5SHSpm9/Wgt5kQJSQyleGvmG2HGIiJ4JjfpgaKu4KC4uRkhICHx9fTFp0iTY29vj008/VWtfQRAwYcIE9OjRA+vWrdNKHqLacObKcUgNpWjqaQoDQ2Ox4xARPRMa9Tbr2rXrE7dHRkaqdZwlS5bgwoULiI+Ph4ODA/bu3Yvhw4cjODgYHTp0eOK+H3/8MeRyeYVbNUR1jaF/EZp/1hwDC5uJHYWI6JnRqMD476ieCoUCt2/fxs6dOzFt2jS1j7N582bMmDEDDg4OAIBhw4ahZcuW2Lx58xMLjD/++AMfffQRLl26VKGzKVFdUlpahni9Akh0pOjRapDYcYiInhmNCozFixdXun7gwIH49ddf1TpGWloa7t+/j3bt2qmsDw4ORnR0dJX75eXlYezYsfjss89gb2+v1mvJZDLIZDLlcm5urlr7EdXU6Su/4pGOFPoKAb2Dxoodh4jomdHqOBjDhw/HiRMn1GqblZUFALCxsVFZb2Njo9xWmRkzZqB3794YOnSo2rnCw8NhYWGh/HJxcVF7X6Ka+HzjGsS+G4vSA1kwMlK/AzMRUX2n1QLj+vXrkEgkarXV09MDUN7R89+Kioqgr69f6T779+/HwYMHMWLECERGRiIyMhIFBQVIS0t7Yr+P0NBQ5OTkKL+Sk5PVfEdENXPtxl8ouV8CE5mJ2FGIiJ4pjW6RTJ06tcK67OxsHDlyBK+99ppax3BxcYFUKkVqaqrK+tTUVLi6ula6T1lZGVq2bImlS5cq1927dw+FhYVIS0vDqVOnoKOjU2E/AwMDGBgYqJWLSFsUCgEpcY8AAN07dRM3DBHRM6ZRgZGfn19hnZ2dHb766iuMHavefWZjY2N07twZ+/fvx8SJEwEABQUFOH78uEofj7i4OOTl5aFt27YYMWIERowYoXKcNm3aoGfPnli7dq0mb4Wo1hw5/QsKM0oACfDapHlixyEieqY0KjB++OEHrbx4WFgY+vbti9DQUHTq1Amffvop7OzsVJ5EWbFiBc6fP48bN25o5TWJnpUf934BADB3MoSHm7/IaYiIni2t9sGorh49euDkyZNISkrCunXr0KJFC0RGRsLU1FTZxtvbG4GBgVUeo23btvD09HwWcYmq5fLVKwAAD+8m4gYhIhKBRBAEQZMdf/vtN6xatQp//fUXBEGAv78/3nnnHTz33HPazqh1ubm5sLCwQE5OTrWGJidSlyAIsHQzRm5yMWa8MxwbVu4ROxIRkVao+xmq0RWMbdu2oX///rC0tMRbb72FuXPnwtLSEv3798e2bds0Dk3UUNxMvAVdbyMYOBpgyphZYschInrmNLqC4efnh3fffReTJ09WWb9582asXLkSMTExWgtYG3gFg2rbhj2L8XnebjiVCjg8lf2HiKjhqNUrGAkJCQgJqTjtdEhICBISEjQ5JFGDEvvgPADACzZPaUlE1DBpVGA4OjpWOrBVREQEHB0daxyKqL77MzURijIFWtq0FTsKEZEoNHpM9c0338S4ceMwa9YsBAcHAwCioqLw6aefVjlPCVFj8dfdOzix/DYgBfxOLBE7DhGRKDQqMObMmQNLS0usWrUKK1asAAB4eXlh7dq1mDJlilYDEtU3W3/4CEKZAD0zHXTvNFzsOEREotCowACAKVOmYMqUKSgtLQXwz9wiRI3dmajyCf8cvcwhrWToeiKixkDjAuMxFhZEquJiUwAArQN8RE5CRCQeUUfyJGpo7mU/wv2E8rl6Qp4fKXIaIiLxsMAg0qItP62DvFABqYEEY4a+LnYcIiLRsMAg0qLfIg4AAOw8TWFgaCRyGiIi8dS4DwYR/cOwWTGajmyKHvacgI+IGje1C4yuXbuqfdDKBuEiauge5cuQ7SagSbMmeM3/VbHjEBGJSu0CY8CAAcrv7927hw0bNqBHjx4ICgoCAFy8eBGnTp3CzJkztZ+SqB74PfoQsnWl0FcI6NxqlNhxiIhEpXaBsWDBAuX3w4cPx7p16/Dmm2+qtFm3bh1OnTqlvXRE9cjBo5vx6MEjtPW0gL6BqdhxiIhEpdFsqra2tkhISKgwi1pOTg6aNWuGBw8eaC1gbeBsqlQbfLs2ReyZDHR63gtnf7ktdhwiolpRq7OpAsCFCxcqXSeRSDQ9JFG9VSArQ0r8IwBA3249xA1DRFQHaPQUyWuvvYbRo0dj9uzZCAoKgiAIuHTpEtauXYtZs2ZpOyNRnbfvxF4UppcAEuDVCW+LHYeISHQaFRjLli2Ds7MzPv74YyxZUj5bpJeXF8LDwzFt2jStBiSqDw78+g0AwMLJEM4uviKnISISn0YFhkQiwWuvvYbXXnsNJSUlkEgknJOEGrVr168CALz8moqchIiobqjxSJ4sLqixKy6VIzk+CwDQq3NnkdMQEdUNGhUYpaWlWLx4MVxcXKCvr69cP3PmTMTFxWktHFF9cPz878hLlQEApk2cK3IaIqK6QaMCIzw8HDt37sSaNWtU1nfv3h1Lly7VSjCi+uLm3V/h+5Ev2r/uCm/vdmLHISKqEzQqML799lvs2LEDo0ePVlnfs2dPHDx4UCvBiOqLO7mXoWelhx7t2LmTiOgxjQqM1NRUeHt7A4DKuBe6urooKirSTjKieqBUrkCyTiYAoJ1jR5HTEBHVHRoVGD4+Pjhz5gwA1QJj27ZtaNWqlXaSEdUDEVcu4tRnd5B5MBPdAzj/CBHRYxo9phoaGoqJEycq5yfZvXs3Dh8+jM2bN2Pnzp1aDUhUl+3cvQ551/JRliqDk2OA2HGIiOoMjQqM8ePHQyKRICwsDAqFAiNHjoSPjw+2bt2K4cOHazsjUZ31R/Q5AICbr7XISYiI6haNCgwAGDduHMaNG4eioiIoFAqYmJhoMxdRnVcqVyAxLgMA0Kl9G3HDEBHVMRoXGI8ZGRlpIwdRvXP25p/IulPeqXnySA6RT0T0b2oXGH369FH7oMePH9coDFF98tOujyGUCdC30EXXrrw1SET0b2oXGO3bt6/NHET1zh+XTwMAXH0sVZ6mIiKiahQYK1asqM0cRPVKSZkC2SU5kOhLEBzYUuw4RER1To37YBA1Rn8k3IXpWDs0H90E/9clVOw4RER1jsYFxsOHD3HixAncvXsXZWVlKtvmz59f42BEddm5Kz+iRCqBrQD4+6jfP4mIqLHQqMA4d+4cBg0aBGtra8THx6NFixaIi4uDTCZD69atWWBQgxd/PwIwB5oLFpBINRoQl4ioQdPoN+Pbb7+N0NBQ5dTsN27cQGpqKgYOHIj+/ftrNSBRXVNcKsferyNx+/9uwyiVA2wREVVGowLj6tWrmDat/Ll/qVQKmUwGGxsbbNiwAdu2bdNqQKK6Jio2Gffj8iFLk6FLQF+x4xAR1UkaFRj5+fmwsLAAANjZ2SEpKQkAYGFhgaysLO2lI6qD9v3yGeRFCkgNpHhh0HSx4xAR1Uk1vnn83HPPYd68eTh48CBmzpyJtm3baiMXUZ116cKvAABHLzPo6xuInIaIqG7SqMD4+OOPld+vWrUKRUVFGDt2LG7fvo0vv/xSa+GI6pqiEjniY5MBAEHtfEVOQ0RUd2n0FMns2bOV3zs5OeHYsWPaykNUp52NScT9uHwAwPhhE0ROQ0RUd/H5OqJq2Lt/HeRFCugYSTF08GtixyEiqrM0KjB++uknTJhQ8a+3CRMmYNeuXTUORVRXPSi+CqvuVvDt0AR6evpixyEiqrM0ukWyePFi7Nmzp8L6hQsXYtSoURg5cmS1jpecnIz79+/Dx8cH5ubmT20vl8sRGxsLXV1deHh4QFeXI55T7SuQlaHMOQdOU5zwlk1PseMQEdVpGl3BSEhIgL29fYX1TZs2VQ6+pY7i4mKEhITA19cXkyZNgr29PT799NMn7vPBBx/AyckJISEh6NevH9zd3XHw4MFqvwei6jofm4gEAzkAoGdA9YpoIqLGRqMCw8/Pr9IrGLt27YK3t7fax1myZAkuXLiA+Ph4/PXXX9i+fTvefPNNREVFVdpeLpejqKgIN2/exM2bN5GYmIipU6dizJgxSE9P1+StEKnt6O9fICepGNalCni4dBM7DhFRnSYRBEGo7k7bt2/HK6+8ghkzZqB79+4QBAERERH4/PPPsWnTJowfP16t49jb22PGjBl4//33lesCAgLQpUsXfPHFF2od4/79+7C3t8cvv/yC559/Xq19cnNzYWFhgZycHLVuyRABQIch3rhwMA4+HZvg1rkMseMQEYlC3c9QjTovjB8/HnK5HGFhYcoxMXx8fPD111+rXVykpaXh/v37aNeuncr64OBgREdHq53l4sWLAIBmzZpV2UYmk0EmkymXc3Nz1T4+EQDkFpci8fY9AEDHwACR0xAR1X0a946cNGkSJk2ahIKCAkgkEhgbG1dr/8dDitvY2Kist7GxUXu48QcPHmDWrFkYPXo0fH2rHvQoPDwcS5YsqVY+on87cfkaHiQUAgBeGjVN5DRERHVfjcfBMDExwf3793H58mWUlJSovZ+enh6A8o6e/1ZUVAR9/ac//peTk4MBAwbA3t4eGzdufGLb0NBQ5OTkKL+Sk5PVzkkEAHv3fQihVIC+hS56dh8ldhwiojqvWgVGamoq5s2bp7Ju7ty58PT0RLt27eDv76+c+OxpXFxcIJVKkZqaWuE1XF1dn7hvbm4u+vXrBx0dHRw+fBhmZmZPbG9gYABzc3OVL6LquHL5LADAq6UtpFKOT0dE9DTV+k354YcfwsPDQ7l86dIlfPTRRwgPD8fp06fh6uqq9q0IY2NjdO7cGfv371euKygowPHjx9G37z9TYMfFxan0yXhcXADA0aNHlbO6EtWWpIcFuHMrEwDQv2dXkdMQEdUP1eqDcfDgQZw4cUK5fODAAbRv3x7z588HAHzyyScYPHiw2scLCwtD3759ERoaik6dOuHTTz+FnZ0dpk375x73ihUrcP78edy4cQOlpaUYOHAgEhMTsWnTJly/fl3ZztvbG02bNq3O2yFSy/4T+5CbXH4r77UX54uchoiofqhWgZGWloYmTZool8+dO4eePXsql728vHD//n21j9ejRw+cPHkSn332GS5cuICAgABs3boVpqamyjbe3t7Kvh2FhYWQSCTw9vZGeHi4yrHmz59freKGSF1JD47AI9QDZoky+Pi0e/oORERUvQLD2dkZ58+fx3PPPYe8vDycPXsWr7/+unJ7amoqnJycqhWgS5cu6NKlS5Xb3333XeX3FhYWiIyMrNbxiWqiTK5AhuIvmPiaYGwrFhdEROqqVoExceJETJgwAZMnT8bZs2dhYmKCPn36KLdHRETgueee03pIIrFEJ2XijlEBACmeaz5E7DhERPVGtQqM0NBQ5OXlYefOnbC2tsauXbtgYmKi3L5nzx6sWLFC6yGJxPLzL5/jz73psAkwReA4Pp5KRKQujYYKr+84VDipq8eoAETsugFHf3Ok/pkjdhwiItGp+xnKB/qJqpBVUILbMeXjunQMai5yGiKi+oUFBlEVjv4Rjfu3CwAAL4+aKnIaIqL6hQUGURX2710JhUwBfQtdDBo4Rew4RET1CgsMokqUyhW4evkcAKB5G3sOD05EVE38rUlUiYsJmUi8+RAAMPL5gSKnISKqf1hgEFXi6KnvIJdKINGVYPrk98WOQ0RU71RrHAyixiKj8BR8Vvqg5QM9NGlSvdFpiYiIVzCIKrj7sBDp+ncBAINa9BA5DRFR/cQCg+g/jkSdxm1dOQCgd9vJIqchIqqfWGAQ/cf+n1fh5qxbKPw+HQ5NW4kdh4ioXmKBQfQv2QUl+OvGVSiKFXDQsxQ7DhFRvcUCg+hffv3jJpJuPAIATAoZK24YIqJ6jAUG0b/s27MMZfly6JrqYPLEBWLHISKqt1hgEP0tt7gUV/+IBAC0bOcAfX0DkRMREdVfLDCI/nb0ajwSrmcBAF4aNVrkNERE9RsLDKK//bxrGUpzyqBjJMW0l5eIHYeIqF7jSJ5EAApkZSgxi4HdCDt4C6YwNjEVOxIRUb3GAoMIwG837+K+UzHs3OywouUbYschIqr3eIuECEDkxW+RrSuFmUJAh1YviR2HiKjeY4FBjV5WQQnORf6AnKgcBJXaQE/PUOxIRET1Hm+RUKP3c9QNRB1IRMmjMth9+LzYcYiIGgRewaBGb9+uhSh5VAZdEx28NWO12HGIiBoEFhjUqN2+n4drF88BAIK6uvHpESIiLWGBQY3ajuOHcffaIwDA2zPfFjcMEVEDwgKDGi25QsDxA8shlAgwczDAsCHTxY5ERNRgsMCgRut0bAb+/CMWADBgQCAkEonIiYiIGg4WGNRo7T7xHYolACTA+/M+FjsOEVGDwsdUqVFKfFCAvLK98HzPE+2zjdDCv4PYkYiIGhRewaBG6fvfI/GneQ4AYMZzHBqciEjbWGBQo5NbXIpLkUuRn6eAp1yKoIBJYkciImpwWGBQo/PjuVs49uMl3Jp7C60e+rNzJxFRLWCBQY2KXCFg7/a5KEwvgVRPglmTOXInEVFtYIFBjcrBq6n4I/IMAKBLH2/Y2TqJnIiIqGFigUGNRkmZAt9tWYj7MQWABPhw8XqxIxERNVgsMKjR+CEqAWeP7AQANO9kj+DAPiInIiJquFhgUKNQICvDjm9fx/2/CiDRkWDzJ1vEjkRE1KBxoC1qFL75/TqSiq5BoidB94G+6NCur9iRiIgaNBYY1ODdfViIqzfeh05fWwS1McfWsQfEjkRE1OCxwKAGrVSuwILtOxFnFQ9AirltXoCLs5fYsYiIGjz2waAGbe2RP3F+7zykxxWhuUIXI3qGiR2JiKhR4BUMarDOJzzEgW8n4PZvDyCJeIjtZ3dDV9dA7FhERI1CnbiCkZycjEuXLiE3N7dW96HGIyrhIRZ98j7O7LsBABjzUld0DBoucioiosZD1AKjuLgYISEh8PX1xaRJk2Bvb49PP/1U6/tQ4/LLtXt496stOL97IxQlApxaWGDblyfFjkVE1KiIeotkyZIluHDhAuLj4+Hg4IC9e/di+PDhCA4ORocOHbS2DzUOaY+K8FVEAo6fP4XkiCXIS5FBz1wXh3f+Ch0dHbHjERE1KqIWGJs3b8aMGTPg4OAAABg2bBhatmyJzZs3V1ksaLIPNUyCICD1URGuJD/CyZhMHLiSDH/JDZTcXYe753IACbB546do2byT2FGJiBod0QqMtLQ03L9/H+3atVNZHxwcjOjoaK3tAwAymQwymUy5rO1+G5HRB7H515W4eiq5yjYeAU3g6mcNAMjLKsblE0lVtnXzt4F7C1sAQEFuMS4dqbqts48VmrVuAgAoKijFhUOJ5RsEoUJbB08L+LRrCgAoKS7DuQMJVR63qZsZfIPtAQDyMgXO/BxfZVtbZ1P4d3JQLkfsvF1FSwFW9iYI6OqoXBP5czwUckWlrS2aGKF1D2flvmf3J6JEJocgCBAEQCEIUPz9Nk2tdNBygDHi9XVwN6IYADD7nWmYMOq1KnMTEVHtEa3AyMrKAgDY2NiorLexsVFu08Y+ABAeHo4lS5bUJO4TxaddwW+5KbhztOpCIMWyGE0Cyj/4CgsKkfCEtneNChHbpgRAeZ+TuCe0TZIWIL59GQCgpLQEsU9oa/2cNe50Kv9ELisuQ8yxu1W2texiibtdJQAAhaDAzSe0NW9vjtQe/3TnufGEtqatTJHe+58fuz9/uwuhpGIxBADGvsbI6KevXP4rIhnyPHmlbY08jFAwtBkAoPPUZuj74gt485UFVeYgIqLaJVqBoaenB6D8A/TfioqKoK+vX9kuGu0DAKGhoZgzZ45yOTc3Fy4uLhrlrkwzxzbobeGKq/0lVbbxdG0Cl7LywijXtAjW/Sv/UAUAN09buJf9fQXDUAbL/pX/hQ8ALj7W8JTbAQCK9EtgPqCsyraOzazgLS+/KlGiUwbTASVVtm3qZgk/efmVBjkUMB5QXGXbJs7m8Fc4K5cNBxRW2dba3hQBCjflsknfYsjllZ8LyybGaK1wBwBIJIBF71KUlsghBSCVSiCVSKAjlUAikcDKxhiv+k1A2+ajYW/uXOnxiIjo2RGtwHBxcYFUKkVqaqrK+tTUVLi6umptHwAwMDCAgUHtjX/Qte1gdG07GAitxk5vV6PtnKc3UXqzGm1nVqPttGq0nVxLbV+uRlsiIhKVaI+pGhsbo3Pnzti/f79yXUFBAY4fP46+ff+ZiCouLk7Zv0LdfYiIiEhcEkGopDfgM3Lq1Cn07dsXc+fORadOnfDpp58iMTERV65cgampKQBg6tSpOH/+PG7cuKH2Pk+Tm5sLCwsL5OTkwNzcvNbeHxERUUOj7meoqANt9ejRAydPnkRSUhLWrVuHFi1aIDIyUqVQ8Pb2RmBgYLX2ISIiInGJegVDLLyCQUREpJl6cQWDiIiIGiYWGERERKR1LDCIiIhI61hgEBERkdaxwCAiIiKtY4FBREREWifqdO1iefxkrrZnVSUiImroHn92Pm2Ui0ZZYOTl5QGAVic8IyIiakzy8vJgYWFR5fZGOdCWQqFAWloazMzMIJFUPQNqdTyeoTU5OZmDd2kBz6f28ZxqF8+n9vGcaldtnU9BEJCXlwdHR0dIpVX3tGiUVzCkUimcnWtnSm9zc3P+j6FFPJ/ax3OqXTyf2sdzql21cT6fdOXiMXbyJCIiIq1jgUFERERaxwJDSwwMDPD+++/DwMBA7CgNAs+n9vGcahfPp/bxnGqX2OezUXbyJCIiotrFKxhERESkdSwwiIiISOtYYBAREZHWNcpxMLQtLy8Pt27dgq2tLdzd3cWOU++UlZXh0qVLsLKygq+vb6VtCgoKEBMTAysrK3h6ej7jhPVLWVkZYmNjYWBgAHd3d+jo6FTa7s6dO3jw4AH8/Pxgamr6jFPWLwUFBYiNjYWNjQ1cXV0rbSOTyfDnn3/CxMSkyp9jUpWfn48rV67A3t4eXl5eFbanpqbi3r17aNasGaysrERIWD9ERkZWWNesWTM4ODiorCsrK8Off/4JXV1d+Pv7a22gySoJVCNff/21YGxsLPj6+gomJibCwIEDhfz8fLFj1Qv5+fnCwoULBRcXF8HMzEwYM2ZMpe2+//57wczMTPDx8RHMzMyEXr16CY8ePXrGaeu+srIyYdGiRUKTJk0Ef39/wcXFRXB3dxeOHTum0i4/P18YOHCgYGJiIvj6+grGxsbCxo0bRUpdt2VlZQlTp04VrK2thcDAQMHa2lpo06aNcPPmTZV2hw4dEmxsbARPT0/ByspKaNeunZCWliZS6vpj9OjRglQqFV555RWV9SUlJcL48eMFQ0NDoXnz5oKhoaGwYsUKkVLWfQCEli1bCl26dFF+7d69W6XN2bNnBScnJ8HFxUWws7MT/Pz8hFu3btVurlo9egN37do1QSqVCtu3bxcEQRAyMjIEd3d3YdasWSInqx8SExOFxYsXC6mpqcKgQYMqLTBu374t6OnpCV999ZUgCIKQnZ0t+Pr6CpMnT37Wceu8vLw8YcmSJcriS6FQCPPmzRPMzc2F7OxsZbs33nhD8PT0FDIzMwVBEIStW7cKUqlUuH79uhix67SbN28KO3bsEORyuSAIglBcXCz07dtX6NKli7JNZmamYGZmJixdulQQBEEoKioSOnbsKAwcOFCUzPXFl19+KXTr1k3o1KlThQIjLCxMsLOzE+7cuSMIgiAcPXpUkEgkwokTJ8SIWucBqPCHxL8VFhYKjo6OwsyZMwVBKP9j5PnnnxcCAwNrN1etHr2BmzNnjuDl5aWybsWKFYKFhYVQVlYmUqr6qaoCY9GiRYKDg4OgUCiU69avXy8YGhoKhYWFzzJivRQXFycAEE6dOiUIgiCUlpYK5ubmwurVq1XaeXh4CHPnzhUjYr2zaNEiwd3dXbm8YcMGwdjYWCgoKFCu27VrlyCRSHgVowo3btwQHBwchDt37ghdunSpUGB4enoK8+bNU1nXsWNHYcKECc8yZr0BQNiyZYtw6dIl4eHDhxW279mzp8LPY2RkpABAiI6OrrVc7ORZA9HR0WjXrp3KuuDgYOTk5CAhIUGkVA1LdHQ0AgMDVe4VBgcHo7i4GDExMSImqx8uXrwIiUSi7LeSkJCA3NzcCj+3QUFBiI6OFiNivXDt2jWcOnUKX3zxBb788kssWbJEuS06OhrNmzeHsbGxcl1wcDAEQcCVK1dESFu3FRUVYcyYMVi9ejXc3NwqbM/NzUVCQkKlv1v5M1q1t956C1OmTIGjoyOGDRuGzMxM5bbo6Gg4Ojqq9MkIDg5WbqstLDBqICsrCzY2NirrHi9nZWWJEanB4TnWXFpaGubMmYPJkycrJ/d7fM4qO6c8n1X77LPP8M477yA0NBTt2rVD7969ldv4M1o9s2fPRps2bTB+/PhKt/NntPq++uorZGZm4urVq4iNjcWtW7cwdepU5fbKfkb19PRgZmZWq+eUBUYN6Onpobi4WGVdUVERAEBfX1+MSA0Oz7FmHjx4gP79+8PX1xfr169XrtfT0wOASs8pz2fVvvzyS0RFRSE1NRUGBgbo168fFAoFAP6MVsdvv/2G77//HuPGjUNkZCQiIyORm5uL+/fvIzIyEnK5nD+jGnj11VeVV3ldXV3x3nvv4cCBAygoKABQ+c8oUH6Oa/Oc8jHVGnBzc0NqaqrKusfLVT3KRtXj5uaGuLg4lXU8x0/28OFD9OnTB9bW1jh48CCMjIyU2x5fkk5NTUVQUJByfWpqKs+nGoyNjfH666+jT58+uHv3Ltzd3eHm5oaLFy+qtOPPaOVkMhnatGmD8PBw5bqkpCRkZGRg/vz5OHz4MOzt7WFgYFDp71aeT/U0bdoUgiDg3r178PLygpubG9LT06FQKCCVll9XyMjIQGlpaa2eU17BqIG+ffsiIiICOTk5ynX79u1D27ZtK1yOIs307dsXUVFRyMjIUK7bt28fvL29K71/29hlZWWhT58+sLCwwKFDh2BiYqKy3dbWFm3atMH+/fuV67KzsxEREYG+ffs+67h13uO/AP8tLi4OOjo6sLS0BFD+MxofH4+bN28q2+zbtw/W1tYIDAx8VlHrhYEDByqvXDz+CggIwODBgxEZGQlTU1Po6OigV69eKj+jMpkMhw8f5s9oJSr7GT169CjMzMyUxUPfvn2Rm5uL33//Xdlm37590NfXR/fu3WsvXK11H20EioqKBH9/f6F79+7C3r17hUWLFgk6OjrCr7/+Kna0euPMmTPC6dOnhc6dOwu9e/cWTp8+LURFRSm3l5aWCoGBgULHjh2FPXv2CB988IGgo6Mj7Nq1S8TUdVNhYaEQGBgoODk5CYcPHxZOnz6t/Hr8SKoglI/ZoKOjI7z//vvC3r17hW7dugn+/v5CUVGRiOnrpmXLlgkvv/yysGPHDuHIkSPCBx98IJibm1d44qZfv36Cv7+/sHPnTmHdunWCgYGBsGHDBpFS1y+VPUVy4cIFwcDAQJg9e7awf/9+YdCgQYKzs3OlT0g0dhs3bhSGDRsmbNmyRTh06JDw1ltvCbq6usL69etV2k2aNElwdXUVvv/+e+Hrr78WLCwshEWLFtVqNs6mWkMPHjzAypUrER0dDRsbG7z22mvo1auX2LHqjd69e0Mmk6mss7Gxwb59+5TLjx49wsqVK3Hx4kVYWVlh6tSp6N+//7OOWuelp6dj5MiRlW5bsmSJSsfEkydP4osvvsDDhw/Rtm1bvPvuu7C1tX1WUeuV3bt34+eff0ZmZiZcXV0xevToCn9JFxYW4qOPPkJERASMjY0xYcIEjBo1SqTE9cvMmTPh5uaGd999V2X9xYsX8cknnyAtLQ3NmzfH/PnzlZ2VSdWRI0ewY8cOpKenw8PDA1OmTFG5BQoApaWlWL9+PY4cOQJdXV2EhITg5ZdfrtXRPFlgEBERkdaxDwYRERFpHQsMIiIi0joWGERERKR1LDCIiIhI61hgEBERkdaxwCAiIiKtY4FBREREWse5SIjoidLT0xEREYHRo0c/89dOSEjAhQsXAAAtWrRAQECAxsfS5H08ePAAx48fBwC4uLigS5cuGr8+UWPDAoOokXr48CGOHTv2xDbt2rVDfHw8xo8fL0qB8dtvv+F///sfhgwZAqlUWqMC48qVK9V+Hw8fPsTevXtx6dIltGnThgUGUTWwwCBqpLKzs7F3717l8qVLl5CVlYV+/fop11lbW8PBwQFjxowRIWG5Jk2a4IcffqjxcTR5H76+vvjhhx/w2muv4cGDBzXOQNSYsMAgaqS8vLxUPrhfe+01XLp0qcKHeXp6OoYOHapcTklJwYULFzBs2DBcv34dSUlJaNWqFdzd3SGXyxEVFYWsrCwEBQWhadOmFV733r17uHTpEszMzBAYGAhzc/Nq5db09Zs2bVrp+xg+fDiuX7+Ou3fvokWLFvDw8KhWHiKqHAsMInqi/95aOH/+PF566SW0bNkS+vr6kEqlOHPmDNauXYstW7bAyMgIcrkc169fx/Hjx1UmXVq2bBnWrFmDjh07orCwEDExMdi+fTv69Omjdh5NX7+q99G5c2cUFBTA1NQUp06dwhdffIHJkydr8QwSNU58ioSIqq2wsBCTJ0/G6dOncerUKbz88suYNWsWZs6ciVOnTiEyMhIDBw5EeHi4cp8DBw5g/fr1uHbtGg4fPoyIiAiEh4dj0qRJKC4urvXXr+o4w4YNw9mzZ3H06FEsX74c8+fP1+icEJEqFhhEVG06Ojp45ZVXlMudOnWCoaEhXnzxRZV1sbGxyuXNmzejRYsWuHDhAnbu3ImffvoJurq6SE9Px19//VXrr1/VcaZNm6Zc7tmzJzIyMvDo0aNq5SGiiniLhIiqzdTUFHp6esplAwMDWFhYQCqVqqz795WJO3fuQCaTYdeuXSrHGjNmjMp+tfX66h4HQLWvqBBRRSwwiOiZMDc3h6urK7Zs2SJ2FCJ6BniLhIieiQEDBmDfvn24d++eyvrU1FSREhFRbeIVDCJ6Jv73v//h4MGD6NChA2bOnAlra2tcvnwZERERuHnzptjxiEjLeAWDiAAAQUFB6N+/f4X1/x2gysXFBSEhISpt3N3dMWzYMJV1Xl5eGDx4sHLZyMgIJ0+exPLly5GYmIjLly+jffv2uHz58hNzFRQU4IcffsD169dr9PrqvA9LS0uMGTMGRkZGAMqHCv/hhx8QHx//xIxEVJFEEARB7BBERJU5efIkvvzySwDAiBEjnvlw5bdu3cL7778PAOjYsSNmz579TF+fqD5jgUFERERax1skREREpHUsMIiIiEjrWGAQERGR1rHAICIiIq1jgUFERERaxwKDiIiItI4FBhEREWkdCwwiIiLSOhYYREREpHUsMIiIiEjrWGAQERGR1v0/UzXJ9TNViTAAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fluid = base.submodels[\"fluid\"]\n", "storage = EnthalpyStorage(param, SensibleMaterial(\"Capsule\"), geometry, label=\"Capsule\")\n", "sensible = pytes.ComposedModel(\n", " {\n", " \"fluid\": fluid,\n", " \"storage\": storage,\n", " \"coupling\": WallCoupling(param, param.fluid, \"Capsule\", geometry),\n", " \"device energy\": DeviceEnergy(param, fluid, [storage]),\n", " },\n", " param,\n", " name=\"Sensible packed bed\",\n", " defaults_from=base,\n", ")\n", "pv = base.default_parameter_values\n", "pv.update(\n", " {\n", " \"Capsule conductivity [W.m-1.K-1]\": 1e3,\n", " \"Capsule density [kg.m-3]\": 1000,\n", " \"Capsule specific heat capacity [J.kg-1.K-1]\": 2000,\n", " \"Heat transfer coefficient [W.m-2.K-1]\": 300,\n", " },\n", " check_already_exists=False,\n", ")\n", "h, R, eps = 300, pv[\"Capsule radius [m]\"], pv[\"Bed porosity\"]\n", "Z, A = pv[\"Bed length [m]\"], pv[\"Bed cross-sectional area [m2]\"]\n", "rho_f = pv[\"Heat transfer fluid density [kg.m-3]\"]\n", "c_f = pv[\"Heat transfer fluid specific heat capacity [J.kg-1.K-1]\"]\n", "a = 3 * (1 - eps) / R\n", "u = pv[\"Mass flow rate [kg.s-1]\"] / (rho_f * eps * A)\n", "xi = h * a * Z / (eps * rho_f * c_f * u)\n", "rate = h * a / ((1 - eps) * 1000 * 2000)\n", "\n", "t = np.linspace(0, 3000, 151)\n", "exact = [validation.schumann(xi, rate * (ti - Z / u))[0] for ti in t]\n", "fig, ax = plt.subplots(figsize=(6, 4))\n", "for n_z in (50, 200, 400):\n", " solution = pybamm.Simulation(\n", " sensible, parameter_values=pv, var_pts={base.r: 5, base.z: n_z}\n", " ).solve([0, t[-1]], t_interp=t)\n", " theta = (solution[\"Outlet temperature [K]\"].entries - 305.15) / 38\n", " ax.plot(t / 60, theta, label=f\"{n_z} cells\")\n", "ax.plot(t / 60, exact, \"k--\", label=\"Schumann\")\n", "ax.set_xlabel(\"Time [min]\")\n", "ax.set_ylabel(\"Scaled outlet temperature\")\n", "ax.legend()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "d3ac1e2b", "metadata": {}, "source": [ "With the default, second-order scheme the model is within 1 % of Schumann's solution on 50\n", "cells. On finer meshes about 0.3 % is left: the difference between a capsule that conducts\n", "heat and Schumann's lumped one. (With first-order upwinding, numerical diffusion leaves\n", "1.3 % on 400 cells: see the [advection scheme](advection_scheme.ipynb) notebook.)" ] }, { "cell_type": "markdown", "id": "c25ab236", "metadata": {}, "source": [ "## Spherical Stefan: a melting capsule\n", "\n", "One capsule, solid at its melting temperature, with its surface held 1 K above it (Stefan\n", "number 0.011), against the quasi-steady solution, exact as the Stefan number tends to 0." ] }, { "cell_type": "code", "execution_count": 3, "id": "4d4965fc", "metadata": { "execution": { "iopub.execute_input": "2026-10-23T17:49:06.490812Z", "iopub.status.busy": "2026-10-23T17:49:06.490678Z", "iopub.status.idle": "2026-10-23T17:49:07.285803Z", "shell.execute_reply": "2026-10-23T17:49:07.284822Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhgAAAFzCAYAAAB8X3AUAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjExLjIsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvgI3uAAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAZilJREFUeJzt3Xd4FNXbxvHvbnoHAiGF0ENvoQTpCNJEpIs/BVFEpajYldcGCIKKWMCORrGigPQuiESFEAgSIAESQkiBdNKzKXvePyKrkQSSsMmkPJ/r4jI7O2fmmWFl78ycOUenlFIIIYQQQpiRXusChBBCCFH7SMAQQgghhNlJwBBCCCGE2UnAEEIIIYTZScAQQgghhNlJwBBCCCGE2UnAEEIIIYTZWWpdgBaMRiNxcXE4OTmh0+m0LkcIIYSoMZRSZGRk4OnpiV5f+nWKOhkw4uLi8Pb21roMIYQQosaKjo6mSZMmpb5fJwOGk5MTUHRynJ2dNa5GCCGEqDnS09Px9vY2fZeWpk4GjKu3RZydnSVgCCGEEBVwoy4G0slTCCGEEGYnAUMIIYQQZicBQwghhBBmVyf7YJSF0WgkLy9P6zJENWZlZYWFhYXWZQghRLUkAaMEeXl5REZGYjQatS5FVHP16tXD3d1dxlMRQoj/kIDxH0opLl26hIWFBd7e3tcdRETUXUopsrOzSUhIAMDDw0PjioQQonqpFgHj0qVLRERE0LlzZ1xcXMrUJjo6mvj4eNq0aWPWR00LCgrIzs7G09MTe3t7s21X1D52dnYAJCQk4ObmJrdLhBDiXzT99TwoKIhJkybRtWtXBgwYQHBw8A3b5ObmMnHiRNq2bcu0adNwd3dn5cqVZqupsLAQAGtra7NtU9ReV0Nofn6+xpUIIUT1omnAOHnyJFOmTOHw4cNlbrNw4UICAwOJiIggNDSU7777jscff7xc2ygLuacuykI+J0IIUTJNb5Hcf//9AMTExJS5jb+/P7Nnzzbd8x43bhydOnXC39+f3r17V0aZQgghhCinatEHo6zi4uKIj4+nR48exZb7+fld9/aKwWDAYDCYXqenp1dajeJa48aNo0mTJqxatarE10IIUR0ppcgrNJKTV0hWXiHZhoKi/+YVkG0oJCuvgJy8QvILjeQXKgqNisICA/q8TPR5megKsrDMz8bCaEBvNKDLz8FGr7BUBiyN+ZwOO0dGeho52dnkZGeTlZPL0B6taeHRAFDoUFcrQafU9Uq9Lrt2Q+l662SznJPyqFEBIyUlBQBXV9diy11dXU3vlWTp0qUsXLiwUmvT2rFjx3jjjTc4cOAA+fn59OrVi6VLl+Lr61tsvTVr1rB48WKioqLw8fFh6dKljBkzRqOqhRCichmNivTcfFKz87mSnceV7Hyu5OSRmpXPlZx/lqVm55GWnYcypGOZl4513hVs8tNxVJm4kIlD/hWs8lKxNKSjM2RiNGRRmJeLi1UBd7SxxIFc7MnluZ2ZXEw3kparSDOov/8LabmKjm56gh9xNNU2b3Um4SnFh0O4zeYkfXVWZj0Hf1ragASM67OyKjrpubm5xZbn5ORct1Pm/Pnzeeqpp0yvr84EV5u88sorTJ8+nffffx+9Xs///d//MXToUEJCQvDy8gJg586dPPjgg3z++eeMGTMGf39/JkyYwKFDh665KiSEENVRQaGR5Kw84tNziU83EJ+eS3JmXlFAyCkKClf+DhNp2QZ0OSnUy0ugvkrFyy6fRro0GurSCA8NxZiTjjEvm4LcXGzyDDjl5pFhUDSvr+edkbbw9/e8x9sZXM4s+QqCr7ueWW3/CQ3bzuUTkVryuil5VkQ5dadQb02BzhqfZsE41TPgYGdj+pPTqg2BjRsDoHQ6oKifl6Li/b0cWg+ocNubUaMCxtVxKWJjY4stj42NpWnTpqW2s7GxwcbGprLL09TWrVuLvV61ahVffvkle/bsMfV1Wb58OWPGjOG+++4D4KmnnuL777/n3Xff5euvvy5xu0op3n33XT766COio6Pp3LkzK1asoH///qb333nnHT788ENiYmJo2bIlTz31FDNnzixz7bt27eL//u//CAsLw9PTk4ceeoinn35aHvsUog4xGhWp2XlFoSEjl4R/BYj4dAMJGblcTsslKdOAtTGbBtlRuBYm0aqejka6KzTmCns2HiAzI4tcQx45hnwyDEbSDQqjgsHNLVg/3cG0v4Y7MkjOKSU0eFhS4OCOsq2Hzr4+Dg6/Q2Ymer0OFwdbXBztcXFyxNnJgXatW8DMBWDtCDaOvOS5mdwChUv9Bri4uFzz59/DKmx/opJPqsaqfcAIDw8nIyMDX19f7O3t6du3L5s3b2bq1KkAZGVlsXfvXhYsWFAp+1dKkZNfWCnbvhE7K4sKP6WQmppKYWGhaVwRpRSHDh3i9ddfL7be0KFD+emnn0rdziuvvMLHH3/M6tWrufXWWzlz5gxffPGFKWC89tprrF27lu+++45OnTpx5MgRJk+ejKOjI3ffffcN68zMzGT8+PG8+eabzJgxg4SEBD799FNCQkLo1q1bhY5dCFH95OYXEpOaQ3RKNhdTsolOySb2Sg6X03NJ+DtA5BcacSYbd10KHroU3Ejm/NE/yM9IIis9jcL0LAozDJxPL+SMEW5racG70/4JDa9eyCAxu+TQkGm0heYDwNENHBoxot8esvJ1uNSrh0u9BjjXb4BLAzdcGjbGq2kLLEePNrX9Y0wC9vb2ODg43PDf5PsfnmOeE1YLaBowEhISOHv2LImJiQCEhIRgaWlJ06ZNTVckli1bxqFDhzh58iQAixcvZtiwYcyfP58+ffqwcuVK3NzcePjhhyulxpz8Qjq8sqtStn0jpxeNwN66Yn9F8+bNw9vbmxEjRgCQkZFBVlYWjRo1Kraem5sbly9fLnEb2dnZvP3227zzzjuMHTsWgF69etGrVy+g6FbVm2++ycaNG/Hz8wNg0KBBPProo3z++edlChgJCQnk5ORw5513Ym9vT/Pmza8JQUKI6s9oVCRkGIhOzeZi8t8hIjXbFCji0w2AwpksXBL/wiH9AtaZsRjTE7DPSKVBeiaJ6bk0bahj6z3/DHLoFlByaNDrIE9vBy0GgWNjcHRjxVNRWDnWp76bFy5uTXBp3BQX9+a4NGhYNDDev8LBt6PeKPOxubm53dS5qas0DRhHjx5lyZIlAPTr14+1a9eydu1aZsyYwYwZMwDw8fEpNunYoEGD2L9/Px988AGBgYF07tyZr7/+GkdHxxL3URe9+OKLbN++nX379t1wNFKlVKmJPCwsjJycHAYOHFji+6GhoWRlZTFy5EjTtq7+admyZZlqbdGiBaNGjaJPnz7cfffd3HrrrQwdOtQ0SqYQovrINBRwMbl4cDD9NzWHvAIjFioft6xzuFwJwzo1koLUODrY5PHlAAea6BJx0uXg9kPpVxqs9Xqwqw/OXuDsyX3DIyiwsMfbuylNmrfEu3VHvNt0xaNFGywti3+FTR1RFWdBlJWmAWPUqFGMGjXquus8//zz1yzr168f/fr1q6yyirGzsuD0Im0+tXZW5e+DsGDBAt5//3127NhhutIA4OTkhIODg+lq0VWJiYk0/rtD0X+pvx+LKi2AXJ0M7sSJE3To0KHctV7d9rZt2/jtt9/Ys2cPL774Io888gh79uyp8DaFEBVjNCpir/xzG+Piv8JDdEo2KVlFv+xZqnzc86Lwscugue4yg3WXWbv1IJGXrhCVksf5guLb9Wmgp/2gf34J7OZlS1KuBd6N6+Pt2bgoPLTwwdunE83a9wCf9qZ1l99bJYcuKkG174OhNZ1OV+HbFFVt4cKFvP3222zfvt3UR+IqnU7HLbfcwv79+3n88cdNy/ft20ffvn1L3F779u2xs7Pjt99+o127dqW+v2PHjpsKAzqdjkGDBjFo0CAWL15Mr169+OKLL1i+fHmFtymEKN3VIHEuIYOz8Zmcjc/gXHwm4QmZ/+pzpnDOvEjDKyewSw3HIzUap9RkElMziErJw8pFT8Dcf0LDqsRMwhKKfumw0EPzBjb4eDWgdTMv2rdrA/feD/Wagos3uxfIPE91Qc345hQ3tGTJEpYvX8727dsZMKDkR5KeeeYZxowZw5o1axgzZgxffvklwcHBfPTRRyWub29vz5NPPslLL72Eh4cHgwcP5syZM3z++ed89NFH2Nvb89xzz7Fo0SI8PT0ZPXo0qampbN++neTkZF566aUb1n3o0CH8/f2ZM2cObdu2JSwsjJiYGFq1anVT50MI8U+QCE8oChFn4zM5l5BBeEIm2Xn/BIn6WVG4phynQ9ZFJnV2pKN1HC1UDP2/jGN/vLHEbV9MUxgbdUTfsBU0aMlrDbMwOjTGp1tfmnfshZXM51TnScCoJV577TUMBsM1/SVeffVV0xM2I0eO5PPPP+e1117joYcewsfHh/Xr1193DIzFixdTr149nnjiCeLi4ujatWuxKwsLFizAzc2NxYsXc//99+Pp6ckdd9zByy+/XKa6e/XqxYkTJ5g+fTpnz56lcePGPPzwwzzyyCPlPwlC1GHJmQbOXM4g7HJG0X/jMzgXn1EsSHiSTIP4QFpdCiE/OYaUpCQiE7KIyi4KEbaW8GMXJywKi26LdmlsQWaBJT5eDfBp4U3r1m3w6dgNn259aNahJ/p/Pf4/clhVH7Go7nRK3cT4ozVUeno6Li4upKWlXTPVe25uLpGRkbRo0QJbW1uNKhQ1hXxeRFXLzS/kXHwmYZfTOXM5gzPxRaEiMeOf6RAc8lNxTwrELvEkutSLvDnKhRYqFluVw8Qfs9kQWryThA5o2dCGjs0a8uVL91G/pS80aotybY3OWm5niOKu9x36b3IFQwghqqksQwEnYtI4Hn2FkNgrhF3K4EJyFkbTr4WKJrpEmsYfwvPiMXISY7h4KYWwpDxO/+tXxw9vccTWRQ96SwZ28MRgV0jHNi3p2LkrHXsOoP0tw7B3rn/N/mWuYHEzJGAIIUQ1UGhUnEvI4PjFKxyPLvpzNj7DFCZ0qhDPzFA6JxzBKukMz/a3p6tVNLbGLJ4Ny2X5n3nFtufmaEG3Fo3o3K4V+jH3Q+cB0KAl8yysmFf1hyfqIAkYQgihgfj0XIJNYSKVkJg0sv7uL6HDiFfOOTrEBGCReIaUy3GciUvjj6x/LkvMaWqPbQtLsLBmkG8TzqscfDt3xLf3AHwHj8ajZYcKjwQshDlIwBBCiEqWnVdAyN+3Oq7+uZR2ddJGRWNDDD4JvzOocRYjGibSqiCcL0KuMGd78YkdLfTQ3sMR37bNcBl5DwwZCw3bcIeFFXdU/WEJcV0SMIQQwowKjYrwhEz+ir5C8L9udRT+fa/DxphNs6Q/8Lh8lLz480TGJBCUmI9Rwe132NLBuejxzlua2dPV25JeHVrSs2cvfPsPp/OA27FzkFGLRc0gAUMIIW5CQkauqd9E8MUrhMSmkWm4+pSGohmXuFMfQT/7SKwuneB/q8M4W3Dtdpo2sKGwaT8Y+wB4dsO3YVuOfyj/RIuaSz69QghRRkajIvRyOkciUzh28QpHo1KJvZJjel+vCmiZfoyGl/4g79IZwi5cpn8bHUuG2oIRUp0VOQXgbKunV2s3evt2pvfA2/AbPgn3pmWbv0eImkIChhBCXEdMajYB55IICE/ij4hk03wcALYY8FNh2JzeTvLFc/wVmcQvGcVHvnS1tQTv3tCkF/Wb9OTc3Q1p2aUveovyzzUkRE0iAUMIIf4lLTufP88ncfBcEr+HJ3EhOdv0npMxg54pv9FJhTOrQxZumWHojPk0+iWT5JyiPhaWeujZ0pV+PTvTf8go+t1xD3g0MW2jdZUfkRDakIAhqoVFixbRsGFD5syZU+JrISqLoaCQo1GpBPwdKEJi00xjT9iqbDolHaRR/GFiz5/l8PkrnDTAYWcdrz7hWPQYqEsTHh9ji0WDpvQfPo5ewyZh7ygdMYWQgFGLXLhwgU8//ZRjx47x3HPPMWTIkGvWMRqN/PDDD2zfvh2dTsc999zDqFGjiq1z7tw5Vq5cSVRUFD4+Pjz55JN4eXlVau3Hjh2jSZMmpb4Wwlyu9qO4etvjyIUUcvOLbmtYk09PXTh3OJ/j993b+CkwjjM5xWdTcLGzwLdtU7JGvoJju1uhXlNekfEmhLiGBIxawt/fn9dee40ZM2awa9cupk6des06RqOR8ePH89dff/HCCy/QqFEjPv74YwoKChgzZgxQFC78/Py4/fbbufvuu/n222/x8/MjODgYNze3qj4sIcwiJjWb38OLbnv8ux+FDiOts07gGbefjIunWDOyEEcLA+TBsdxcUnMUjjZ6BnZqwuABA7h17L34DhiOhfSfEOKGJGDUEqNHj2b69Ono9fpSZzL99NNP2b17N6dOnaJly6Ie6xMnTiQ9Pd20zmuvvYaPjw/ffPMNOp2OSZMm0bp1a95++23eeOONUvd/4sQJVq9eTXR0NJ07d+bpp5/GxcXF9P7x48dZvXo1MTExtGzZktmzZ+Pj41Pm47t06RKrVq0iLCwMT09PHnjgAbp3717m9qJuudqPIiA8iYBzxftRuObH0yd+N/roo4ScjWLv5XzTe4fb2zO0sye0GMicLm34n2Mbetx6B5ZWVlochhA1mgSMG1EK8rNvvF5lsLKHMl56LcvVBX9/f+68805TuLjq37Ph7dq1i8cff9w0xLCVlRVjxoxh586dpQaMzZs3c9ddd/HAAw/wv//9jzNnznDPPfewbds2ALZv3860adOYN28eQ4YM4ciRI3Tv3p2AgAC6du16w7oLCgoYMGAA3bp149577yUhIYFHH32UL774gnbt2t2wvaj9DAWFHIu6QkB4IgHhyYTEXDH1o7Amn/4WZ5nkcoaEv/Yz78dIjv1rHAqdDro3b8DwgX40nTUPeo8AnY422hyKELWGBIwbyc+G1z212ff/xYG1g9k2FxISwsSJE3n//ffZv38/bm5uTJw4keHDhwOQnZ1NQkIC3t7exdo1adKECxculLhNo9HI7Nmzefzxx3nzzTdNy69eFVFKMXv2bN566y1mzJgBwIQJE0hLS2Px4sX89NNPN6z7woULRERE8Oeff9KoUSMAHnroIXJzc2/QUtRmmYYCfj2TwM6Tl/n1TOK/BrcCz4KLdE7cQ0b4Ue5umc6drYyQDcedCsktAK961gzv3ZHho8cxdPJMGrlr9P+4ELWYBIw6QilFbm4ub7/9NqNHj2bq1KmcPn2asWPHsnz5cubOnUteXtF9aTs7u2Jt7e3tTe/9V1hYGHFxcdx9993Fll+9KnLmzBkuXrzImjVr2LBhA0oplFJERUWhlCppk9fw9vbGy8uLhx9+mDlz5tC3b18cHBxwcDBf+BI1Q0pWHntD49l18jIHw5PIKyjqnKnDiG/eX3hG7+ZsyF8cDE/nz6J5w9BnWHFnl2bQeihdxt3KqZlNaN+jr0wEJkQlk4BxI1b2RVcStNq3meh0OlxcXGjRogVffPEFUNT/IjMzk2XLljF37lwcHR2xtLQkJSWlWNvk5GTq169f4nYzMzMBSn0/LS0NgPvuuw9Pz+K/Jf43yJTGxsaGw4cPs2rVKubPn28KRp988kmx2zuidrqUlsPuU/HsPHmZwAsppjk9rMnnLpcIxtodZ+bSH9mYUDwEN3W1ZfSA7kz83/0w6UHQ69EDHar+EISokyRg3IhOZ9bbFFry9fW9pq9GixYtSEpKQimFpaUlHTt2JDg4uNg6wcHBpfaVaN68OQCnT5+mRYsWpb7v5OTEyJEjK1y7l5cXS5cuZenSpcTExNC/f3/effddXnnllQpvU1Rf5xMz2XUqnp2nLvNX9BXTcsf8ZHok7aRJ5ikWdU/CwpAJBrDVF2Chh/7tGjN6xG3cfu8cOnTvI1cphNCQBIw65MEHH+SJJ54gNjYWLy8vDAYDa9eupX///qZ/iKdPn85rr73G888/T4sWLTh27Bi7du3iq6++KnGbbm5ujBkzhgULFnDLLbfg6upKYWEhP/zwA/feey+NGzdm7NixvPLKK/Tp08c0tsW5c+cICQlhwoQJN6z7zJkznD59mvHjxwPg7u6Os7MzhYWFZjozojo4F5/B9pDLbA+5xJn4DNPyxjnn6Ri/g7hTQfx2JpFT+WBjAfPbOuHYwAPajuLr1R1o0nss9Ro21vAIhBD/JgGjljh+/DgvvPCC6fVbb73FN998w8iRI3niiScAuPfeezl69CgdOnSga9euRERE4OHhwbp160ztHnvsMYKCgujcuTMdO3YkJCSERx555Jo+Fv/2xRdfMGXKFFq2bEmnTp2Iiori0UcfNb3v7+/P9OnTadu2LV27diU1NRW9Xs8HH3xQpmNzdXXl66+/Zu7cubRp04aIiAg8PT157LHHynmWRHWilOLMv0JFeELm1XfoYhFFu9gt/H7wDw5FZhD4r+46TetbM26wL9l3vYZjl6Gg19NJkyMQQlyPTpW1p10tkp6ejouLC2lpadfcw8/NzSUyMpIWLVpga2urUYXll5yczJEjR65Z7u3tTceOHYsti42N5dy5c3h4eNCmTZsSLyNHRERw8eJFWrdufc1TJaU5d+4ccXFxdOzYkYYNG17z/sWLF03hwMfHB71eb3ovODgYW1tb2rdvX+JrgMuXL3P27FkaN25M27Zty1RTZaupnxetKKUIvZTB9pBLbD95ifOJWQBYUkC3rENMd7/AEMtgHHIu8e4hA0/uMgDQrakT44b1Z+y0OXQdOFpufQihoet9h/6bBIxaEjCENuTzcmNKKU7FpbM95BI7Tl4mMqkoVFhQSNfMP3C/sIugv05zJCaPz++0ZYavNVjZE1v/Fn4878j4GU/RvIOvxkchhLiqrAFDbpEIIcxOKUVIbBrbQy6z4+Qlov4eSVOHkW7ZQXhGbSXw2El+jvnnyQ+dDsL07eF/S6DlYLys7HhSqwMQQtw0CRhCCLM5n5jJhmOxbPorluiUnL+XKnpZnucR179on/orzd86z9/DV6DXweDOTZg0aTLjH3wad8/KnVRPCFF1JGAIIW5KWnY+W07Esf5YDMEXr5iWtzBG0+PSBvRxwbw5IA/SAD0MamFDrlV9/nfXJCbPegE3DwkVQtRGEjCEEOWmlOLQ+RS+C7zIrpOXySssuiThTCa3pW8l+sgeth2P59cC0AFP9nLDo9cY6DiBbc/2x8ZBBkgToraTgCGEKLOUrDzWHY3m+8BoU2dNPUZGGA9hDN7AtsNneSfDaFq/dWMH7hk3CsvHloNXMwBsNKlcCFHVJGAIIa5LKcWf55P5PjC62NWKTtaXecb9GP0y9/JVwEUe2ls0+VwDewvuHnEL9z36An63yiOlQtRVEjCEECUq6WqFszGN25J+5mLQPu5sls3gHtYATO7ekF0pdtx7/yPcPu1RrG3kOoUQdZ0EDCGESUlXKywpoG/Gr9if2caOwAg++vsWSFyiBQ//byx0+x8ubUby0wIJFUKIf0jAEGa1ZcsWnJycGDx4sNallMmuXbuwsLDgtttu07oUTZV0taKVLpauF74h8HAg318wmNZt6GjJtNEDmPHUq+A3SKuShRDVnAQMYVaff/45TZo0qXDA2Lx5My4uLgwaVDVfXF9//TW2trZ1MmCUdLXCmnwmWR9ljtNBWmYFc2dwNr9fKECvg5E9WzHj4TmMue9RrK2ttS5fCFHNScAQZjVmzBjq1atX4faffvoprVu3rrKAUReVdLWiSWEMvWK+59jhw8wdY0kLvR50ep66awB+GU24/8mFNGnWQuPKhRA1iQSMWuTSpUvs2LEDJycn+vfvz4kTJ4pd/v/2229p37493bt3N7XZsWMHVlZWpnX27dvHiRMngKJZTHv16kW7du2u2dfBgwcJCwvD09OTIUOGYGdnBxRN3+7k5HTDWktqv3v3biIjI8nIyODdd98FiqaPr1+/PgaDgT179hATE0PLli0ZMmQIlpb/fHzLWvd/z9G/HT58mLNnzzJt2rRrao2Njb3ujLLVnVKKw5EpfHf4Ijv/vlphQSG3Zu/DKmQjPwde5PfsommJVp9yZcmLT0P3aQx2acJgbUsXQtRQEjBqicDAQG677TZ8fX1p2rQpzzzzDHZ2dvTv398UHt566y2mTp1aLGB89dVXODo6mtZJTk7mwoULQNEX7pw5c1i4cCFPPfWUqc3EiRMJCgpixIgRJCQk8OKLL7J9+3Y8PT3LdIuktPZJSUnk5OSQnp5uqiE/P5/z588zYsQIGjVqRKdOnfjss88oLCxk3759NGjQoMx1X+8cAeh0Ou6//34GDx5cbAbZOXPmMG7cuAr9vWjtSnYe64/F8t3hKCL+nrm0nkpjbPKPhP6xl69Cs7g622GTBnbMvGc8Dz6zBJo116xmIUQtoeqgtLQ0Bai0tLRr3svJyVGnT59WOTk5xZZnZmaW+qc862ZnZ5dp3fLy8/NT999/v+n1kSNHlE6nUw8++KBpWdeuXdVbb71VrN2UKVOKrfNfv//+u7KxsVGJiYlKKaWioqIUoM6fP29aJyIiQsXFxSmllBo7dqyaO3duqdu7UfvRo0erefPmFWszYMAA9fTTT5teG41GNWLECPXYY4+VuW6lynaOunTpohYtWmR6ffjwYaXT6VRERESJ+ynt86K16JQs9eqmk6rtS9tVs+e3qmbPb1UTXv5AHXl3irryoqtyskZB0Z9h3Vuqjd98qvLz87UuWwhRA1zvO/Tf5ApGGTk6Opb63u233862bdtMr93c3MjOzi5x3UGDBvHrr7+aXjdv3pykpKRr1lNKXbOsNJcvXyYwMJCVK1ealvXs2ZPevXuXeRv/du7cOY4cOUJSUhJGo5GCggJCQ0MZMGAATk5O2NjYsGXLFh555BFsbGxo2bJlmbdd3vaxsbEcPHiQPn36sGrVKpRSKKWws7Pj999/L3PdZT1HM2fO5N133+Wll15Cp9PxxRdfMGjQoHIdo5bOXM7gkwMRbP4rjgKjwpICJqoDOJzeyMLuSehSdWAJjw7xJrNeB+b83xu069xV67KFELWQBIxaIDY2FgAvr+KTRjVp0qTc21q4cCFvv/02Q4cOxdPTEysrK3Q6HSkpKQDUr1+fH3/8kVdffZX58+fTp08fpk2bxn333XfNiI0RERFs2bLF9Hr48OF06NChzO0BYmJiAEhKSsJg+OdRyWbNmuHr61vmust6jqZOncpzzz3Hr7/+yi233MIPP/xQLJRUV0ejUvjo1wj2hiYA4KQyGZ+5jr9+3cE7p4tugwx3c6L/qEnQ+xFef7VX0fzoQghRSSRglFFmZmap71lYWBR7nZCQUOq6er2+2Our/QZuRuPGjQFITEws9gWakJCAi4uL6bWlpSWFhYXF2mZlZZmuzmRnZ7No0SJ27tzJsGHDgKLjfu+994pdUbnzzju58847SUxMZOvWrTz66KMUFhYyY8aMYtvOzs4udnxXz2FZ20NRh02AadOmldqvoyx1l/Uc1a9fnwkTJvDFF1+YQsmkSZNK3K/WlFIEhCex8pdwAi8UBSlP4yUGxa5h1/4/ee1SgWndkX5tsJv+PgwcoVW5Qog6RgJGGTk4OGi+bmm8vLxo27Yt3333Hd26dQPg/Pnz/PHHH/j4+JjW8/b25tSpU6bXV65c4c8//zR1YMzIyMBoNNKoUSPTOl999VWxfV0NT25ubjRq1IgHHniAr776itOnT19TV+fOnU1Pg5S1vbOzc7HbS61bt6ZDhw6sWLGCgQMHmgJaYWEhZ8+epX379mWqu6znCOChhx7i9ttv5+zZs9x9992mJ2SqC6UUv55N5P1fzpmmR+9pGcH9xs3M+nAfyzKKQpWtlY77xgxm3qsr6NClm3YFCyHqpGoRMKKjo4mPj6dNmzY4O994GmelFJGRkaSmpuLt7Y2bm1sVVFl96XQ6VqxYwdixY4mPj6dZs2Z8++23pt/ar5o1axZjxozB3t4eDw8P1q9fj62tren9xo0bc9ttt3H33Xczbdo0zp8/z/bt24s9DpqWlsaoUaMYMGAAbdu2JSwsjCNHjvD222+XqdYbte/fvz+vvvoqTZs2xdHRkenTp/P1118zcuRI+vXrx6hRo0hNTWXv3r3Mnj2b9u3bl6nusp4jgMGDB+Pl5XVNnw2tKaX4JTSB9/ed40RMGjqMDNcH8VKD/TTN/ItCo8LWEtxdbHh0xv94ZP6bNPxX6BJCiKqkacDIzc3l3nvvZceOHTRr1oyoqCjeeOMNHnvssVLbnDhxgilTppCamoqXlxdhYWGMHDmSb775ptr9plmVbr/9do4cOcL69etxdHRk/fr1LF++vNg6I0aMICAggB07duDg4MDatWs5cuQINv+amGrr1q2sWbOG8PBwOnfuzOuvv86KFSto1aoVAD4+Phw/fpx169Zx9uxZfH19WbJkiem2w40G2rpR+1mzZuHm5kZwcDBJSUnk5+fTvXt3zp49y08//UR4eDjNmzdn/fr1tGnTpsx1l3aOTp48iZWV1TV1Dh8+HBsbG/z8/Mr/l2FmRqNi9+l43v/lHKcvpWNBIbelbyf3yE/sCU9l5VxHsLbGoutktm2+nea9hhf7OxVCCC3oVHkeVzCz+fPn88033xAYGIiHhwcbN25k/PjxHDp0qNQnIPr164eTkxNbt27F0tKSixcv0qVLF1566SWeeeaZMu03PT0dFxcX0tLSrrlikpubS2RkJC1atCj2231NNHXqVGxtbVm9erXWpdQoBQUFtGjRgueee+66YRcq9/NiNCp2nLzMyn3nCLucgQ15DElZz6U/N7P5ZIZp/IqfX72LcU+9C84eZt2/EEKU5Hrfof+m6RUMf39/Zs+ejYdH0T+M48aNo1OnTvj7+5caMBITE4uN4ti0aVOaNGlCYmJildUtaq/Vq1ezY8cOdDpdiZ1Oq0KhUbH1RByr9oVzLiETR7IZk/A9Z3/fycdnc0zrjenbgRcWr6DvrdJxUwhR/WgWMOLi4oiPj6dHjx7Flvv5+REcHFxqu0WLFvHMM8/QvHlzmjVrxu7du8nOzmbOnDmltjEYDMUecUxPT7/5A6gBRo4cWeLlf1G6mJgYOnXqxFtvvWWWDrjlUVBoZMuJOFbuC+d8YhbOZPGc7W76p26hq39R51i9DqYM68X8pavo3F372zdCCFEazQLG1fEJrj6GeJWrq6vpvZIMGTIEPz8/5s+fT5MmTTh//jwvvPACTZs2LbXN0qVLWbhwoXkKr0GmTp2qdQk1zoIFC6p8n/mFRjYGx/LB/nAuJGfjpDK5P2c98xv9jk1hJtSHOzvXo1Grrjz/+gf4tO9Y5TUKIUR5aRYwrv5mnZubW2x5Tk5OqVNBK6UYNWoUzZo1IyYmBmtra6Kjo/Hz86OgoICXXnqpxHbz588vNidFenp6sbkmhNBCXoGRDcdi+ODXcKJTcnBUWYy7/BVBB/ayLCaPhx53xKtVRxj0HD+/dCd6y2rx0JcQQpSJZv9ieXt7o9frTYMZXRUbG1vq1Yi4uDiOHTvG4sWLTSHE29ubO++8k82bN5caMGxsbKRXvag2DAWF/BQUw0e/RhB7JQcnspkQv4bA/bt5LyoPAGtLHX96z2LS7DdAr0d/g20KIUR1o1nAsLe3p2/fvmzevNl0KT8rK4u9e/cWu0wdHh5ORkYGvr6+NGjQAL1ebxo++qro6OhigyyZg4YP14gaxGg0lnnd3PxCfgi8yMcHznM5PRdbDExO/Y7gX7byTkRRHyFrSx0PTx7JC298gpdcZRNC1GCaXnNdvHgxw4YNM81JsXLlStzc3Hj44YdN6yxbtoxDhw5x8uRJ7OzseOihh/i///s/CgsLadmyJbt372bnzp3FJhu7GVfnsEhMTKRRo0Ylzo8hhFKKvLw8EhMT0ev1pd7Wu7ru1hOXWLo9lLi0XCwpYI5DAPfmrcPniwvkFIClXsdDk4fz4lufSbAQQtQKmgaMQYMGsX//fj744AMCAwPp3LkzX3/9dbGZS318fMjLyzO9/uCDD+jduze7du1i/fr1NGvWjD/++INbbrnFLDVZWFjQpEkTYmJizDJPiKjd7O3tadq06TVzzFx1Oi6dBVtOERiZgg4jkzjAq/V34JQTAxbw+MCGXHbowKvvfE6LVq2ruHohhKg8mg60pZWyDBJSWFhIfn5+FVcmahILCwssLS1LvMqVmpXH23vO8N3hixiVYlD+H9ge/pwvDyWwf7o9fdt5wqDnUL73obOS/kFCiJqjRgy0VZ1ZWFhcM0uqEDdSUGjk28MXWbHnLGk5+bTPD6XdyVV8/msUGX9fiPs5uxd9520DawfkBpwQoraSgCGEmfwZkcyCzac4E5+BuzGR4eHv8e2ev9iZWXSR0Le1O2+sWMWwMRM1rlQIISqfBAwhblJCRi6vbwtl4/E4HMjhJdut/Pz99yw/W3TJooWbE0tef50pD8wpta+GEELUNhIwhKigQqPim0NRLN91hmyDgXst9jPf7mccC1KxaKfn92hLXn7mMR59camMwyKEqHMkYAhRAcejr/DSxhBOxqbTIy+IBodX0cUjC8dOVuDamvtef5Uxq/vQ0MzjswghRE0hAUOIcriSncebu87wfeBF3AsTmBC2nG/3hZKQpdjhbMHYx9/Ctt/DWFhY0VDrYoUQQkMSMIQoA6UU647GsHRHGFlZmdxx6Qt+37OTdy4VANC2SQPeXfUxtgMna1ypEEJUDxIwhLiBsMvpvLzxJEcupNA/5wAFv33KquPpADjbWbLghad4dP5i0wR+QgghJGAIUaosQwHv/XKOzwMiaaIu8bXNVxRePsqtx7MBmDF+CEs/+g63xo01rlQIIaofCRhC/IdSip0nL7Nwy2mS0zO5L28d/+eyAyuVD63seXWaHyNmvkyfgUO0LlUIIaotCRhC/MuFpCxe3XyKA2cT8c0Lpv0f77DieAozZjvi7TsERq9ggWsrrcsUQohqTwKGEBRNpf7xgQg+/DUCu/wrTD6/gu92HSUuo2gUzp+tJ/H4tC9AZtcVQogykYAh6rz9ZxJYsPkUUclZDM/aTuwef5afyQWglbsLH6/257bR4zWuUgghahYJGKLOSsnK46WNIWwPuYwnSdwTuoQPdp0j3QBWFjqenzOd/3vjQ+zs7LQuVQghahwJGKJO+iM8iSfWHichI5eplvt42fp7Xk1OId0AvTs05fPvN9GxSzetyxRCiBpLAoaoU/ILjbyz5ywfHYjAU13mcz5jqGUoGOHVe/rRasIQZsybj4WFhdalCiFEjSYBQ9QZ0SnZPPZ9MH9Fp3BH2k8Eb/+BFwuNDH7YFYthr2LX+xEe0kuwEEIIc5CAIeqETcdjeennkzgbYpj012t8sDeKnAJwtLEgZOBquvUZo3WJQghRq0jAELVabn4hr246xdqgi4zI3MKZbf4sv5APwNBe7Vn9w1aat2ypcZVCCFH7SMAQtVZkUhazvzlK/KUYJoUt4bNdp0k3gL21BctfX8Csp15EJ+NaCCFEpZCAIWql7SGXeG7dCXrlH+FLm08ZfSSWdAP06dScNet30rpNW61LFEKIWk0ChqhV8gqMLN0Ryg+/h/Gi5bdMtf4FgG8e7MTGwsE8u+htLC3lYy+EEJVN/qUVtUbslRwe/e4YOReCGBCwkHDrLBhsC30epf2Ql2lvZat1iUIIUWdIwBC1wq9nEnjih2D6Xf6OQ5u+ZdelQvQ6mPrKF7Qecq/W5QkhRJ0jAUPUaIVGxbt7z7Jm/3FGhr3Gh9tPkZkHDRytWf3ZZxIuhBBCIxUKGAUFBXz11Vf8/vvvpKSkXPP+xo0bb7YuIW4oMcPAvB+CyTxzgA77lvBmcAYAA7u15tvN+2ji7a1xhUIIUXdVKGDMmzePb7/9ltGjR9OkSRNz1yTEDR0+n8zj3x1lbNZPbPjen6Nxheh08PK8h3j5rQ+lI6cQQmisQv8K//jjj+zdu5eePXuaux4hrstoVHx68Dyf7jrKcosPGGJzHLceVrycbc23369lyEgZkVMIIaqDCgUMvV5Pu3btzF2LENeVlp3P0z8dJ/rk77yX8w4DGl5BWdry4CvvMbn1OFzq1dO6RCGEEH/TV6TRsGHD2LBhg7lrEaJUJ2KuMHrlQTjiT+r3zzL1m2iSrLzQPbgbXc/7JVwIIUQ1U6ErGLa2tsyYMYNNmzbRunXra4ZbXrZsmVmKE0IpxTeHL7JsywnuvLicL34+SEKWor6DFef6vkNDj65alyiEEKIEFQoYFy5cYPDgwaSlpXH06FFz1yQEAFmGAuZvCOHQ8RBGHn+Rt/bEUKigayt3Nuz8jZatfbQuUQghRCkqFDD27t1r7jqEKCb0UjpzvzuG0+VAmu9eyIoT2QDcc8dgPlu7DXt7e40rFEIIcT0V6oMhRGX6IzyJiR/9QdfknbQIeJWfTmRjoYd3Fr/IN5v3SbgQQogaoMKDBURGRvLee+8RGhqKUooOHTowb948WrRoYc76RB2z53Q8j34XxDx+YI71ZlIHWnEoyYEl73/J0FF3aF2eEEKIMqrQFYx9+/bRvn17Dhw4QIsWLWjVqhUHDhygffv27Nu3z9w1ijpi0/FYnvrmd+bFL2CO5WYA6g9/hj/D4iVcCCFEDaNTSqnyNurVqxejR49mwYIFxZYvWLCAbdu2ceTIEXPVVynS09NxcXEhLS0NZ2dnrcsRwDeHovho4z56HHqBlQcSWHWHI3OXfA5d7tK6NCGEEP9S1u/QCgUMGxsb4uPjqfefsQdSU1Nxd3fHYDCUu+CqJAGjevno1wh2bPsZ612vsj6kqDPnc7Om8cZHazSuTAghxH+V9Tu0QrdI6tevz9mzZ69ZfubMGRo0aFCRTYo6SCnFmzvD+GPTZ6T+OJ/1IdlY6nV8+t4bEi6EEKKGq1Anz6lTp3LXXXexePFi/Pz8ADh8+DAvvvgiU6dONWuBonYyGhULtpwiedfbHPlxDWeSjTjbWbJ+3QZuu13mExFCiJquQgHj9ddfR6/XM3PmTNPtEBsbGx5//HEWL15s1gJF7VNoVDz703G8jizjs69/4nKmwruhI9v2/Ebnbr5alyeEEMIMKtQH46rc3FwiIiLQ6XS0bNkSW1tbc9ZWaaQPhnYKCo08szaIgaELmWARwIo/DXx9vgHbfj2Cp5eX1uUJIYS4gUrt5FnTScDQRl6Bkee/C+CO0OcZanMao84C/Z3vk9t+Uo0Jp0IIUdeV9Tu0zLdIHn30UQBWrVpl+rk0q1atKutmRR1hKCjk+TX7sV57H0//Fc8vMxviOv1b8LkNiRZCCFH7lDlgxMTElPizOURHRxMfH0+bNm3KfEWhsLCQ0NBQ7O3tadmypVnrEeaVm1/I8/47Kfh6Ou8FJAPwk+NMZvncpnFlQgghKkuZA8bGjRtNPy9btox27dqVuF5YWFiZd56bm8u9997Ljh07aNasGVFRUbzxxhs89thj1223YcMG5s6di729Pfb29ri7u/P999/TsGHDMu9bVI2cvELmr95Eqv8DfB2UDsDrLz7FI89LZ2AhhKjNKvQUSfv27Smt68b13vuvhQsXEhgYSEREBB4eHmzcuJHx48fj5+dH7969S2zz22+/MXnyZD777DNmzJgBwK+//kp8fLwEjGomy1DAS5+sJXb1w6wLyUangw+XL2bWUy9qXZoQQohKVqFOnjqdrsQQkZmZibu7O5mZmWXajru7O7Nnz+bVV181LevcuTP9+vXj448/LrHNkCFDsLKyYteuXeUt20Q6eVa+jNx8Fny4hnOrH2VLaC4Wevj6sw/534zZWpcmhBDiJpi9kyfAE088UeLPAEajkePHj9O9e/cybSsuLo74+Hh69OhRbLmfnx/BwcEltjEYDAQEBLBixQoyMjI4e/Ysnp6eeHh4XHdfBoOh2PDl6enpZapRVExaTj5LP/FnZvSLDIsxYG2p48dv1zD2LhmETQgh6opyBYzw8PASfwawsrKiR48ezJs3r0zbSklJAcDV1bXYcldXV9N7/5WUlER+fj4nTpygbdu2uLu7Ex4eTp8+ffjuu++u2dZVS5cuZeHChWWqS9ycK9l5LPnYn1fTXsHRJY/NT/clueezDBs9VuvShBBCVKFyBYytW7cCcP/99/Pll1/e1I6trKyAoo6e/5aTk4O1tfV12+zZs4fjx4/j5uZGcnIyffr04ZlnnsHf37/EdvPnz+epp54yvU5PT8fb2/um6hfXSs40sGDlp4w99yqOTQvJ8upP9+k/gbW91qUJIYSoYhWa7OxmwwWAt7c3er2e2NjYYstjY2Np2rRpiW0aNmyIg4MDEyZMwM3NDSi64jF58mQOHjxY6r5sbGxwdnYu9keYV0JGLi+9+yknP3uGO9aksjG1DQ4SLoQQos6qUMAICAhg9uxrO+vNnj2bgICAMm3D3t6evn37snnzZtOyrKws9u7dy7Bhw0zLwsPDTX0y9Ho9w4YNuyaUxMTE0KhRo4ocijCD+PRcXn5/NWf8n+HXyDysLC1oOG6phAshhKjLVAX06tVLHTly5JrlgYGBqnfv3mXezq+//qqsrKzUCy+8oDZt2qRuu+021apVK5WRkWFa58EHH1QdO3Y0vT5x4oRycnJSr7zyitq1a5datGiRsrS0VBs2bCjzftPS0hSg0tLSytxGlCw2NVvNfG2lGtbaWgHKwcZCHdy3W+uyhBBCVJKyfodW6ArG1U6W/9W2bVtOnDhR5u0MGjSI/fv3ExUVxXvvvUfHjh0JCAjA0dHRtI6Pj0+xJ1M6d+7MH3/8QUxMDG+++SaRkZEcOHCA8ePHV+RQxE2IvZLDS6u+JO7rZ9gTnoedlZ5tW7fS/9ZhN24shBCiVqvQOBg+Pj6sWLGCMWPGFFu+adMmnn766WueMKluZByMm5dpKODx977nypqH+flULrZWerZu2sTQUXdoXZoQQohKVNbv0ApdwZg1axYPPvggn376KWFhYYSGhvLJJ58wc+ZMZs2aVeGiRc1QaFQs+nonC9JfxlZXgJWFjvU/rZVwIYQQwqRCQ4U/9dRTpKSk8MQTT5CTkwOAnZ0dTz75JE8//bRZCxTVz6otfzDr4tM0tUzlswc682SvN+k1QCYuE0II8Y8K3SK5Kjs7m9DQUHQ6He3atcPevmY8NSC3SCru5z9Pk/TBaB5slUyuvRf2s/aCi5fWZQkhhKgilTJU+H/Z29tfM9S3qL2OhccStGQ0r2+7wJ4uDny/Z5OECyGEECWqcMAoKCjg9OnTXLx4kYKCgmLvjRs37mbrEtXMpdQMdr40kmXbLwDQauAULNzaaFuUEEKIaqtCAePcuXOMGzeOM2fOUFhYiJWVFfn5+QA4ODiUeTZVUTMY8gv4/v8m8Nr6kxgVPHDXHbz+/mqtyxJCCFGNVegpkieeeIJBgwaRlZUFFM0ncvToUbp168aSJUvMWqDQ3levz2HJV7+QWwDD+nXn029/RqfTaV2WEEKIaqxCnTxdXV0JDQ3Fzc0NnU5HXl4eVlZWnDhxgvHjxxMREVEZtZqNdPIsu61r3ubxec8SeUXRxcebP4JDcXBw0LosIYQQGqnUcTBSUlJMk401bNiQy5cvA9CqVatr5gkRNVf40b2oXQuJzVA0aejEnoNHJFwIIYQokwoFjH/z8/Nj2bJlhIeHs2TJElq3bm2OuoTGUmPCcN3yAGN8dHz8kB+7f/0Tt8aNtS5LCCFEDVGhTp6PPPKI6edly5Zxxx138OGHH1K/fn1+/PFHsxUntFGQkUTKp+NppU/njL4VE9/YhbNzPa3LEkIIUYPc1EBbVymliI2NpXHjxlhZWZmjrkolfTCuo7CAd+/ryjvbQlk9xZMWz/1K61ZyVUoIIUSRSu2DceuttxZ7rdPpaNKkSY0IF+L6fn7tHl5Ye5qLaYo1ObdKuBBCCFEhFQoYR48elbEuaqG/fl7B3LfXYSiEnr6d+erLr7QuSQghRA1VoYAxatQovv/+e3PXIjSUEhbAI48/z6VMhXfj+uzedxC9/qb7AAshhKijKtTJ08nJiUceeYT169fToUMHrK2ti72/bNkysxQnqobKiGfe1Ns5HFOAo60l63b+Rv16LlqXJYQQogarUMC4cOECQ4YMoaCggBMnTpi7JlGVCvP5YPZgvjmagQ5488PP8evWSeuqhBBC1HBlDhhBQUH07NkTgL1791ZaQaJqXV73LPd6xbKltTWNBkxj9gP3aV2SEEKIWqDMN9l79epl+lkG06od0o5twD3Un/p2OiY88xZfffaJ1iUJIYSoJcocMOrVq2caBry6zzUibiwv8Tx7334IpRQ/Wo9n6v1zsLCw0LosIYQQtUSZb5GMGzeOzp0706JFCwDT7ZKSBAUF3XxlovIU5LFg+q0s3ZHCxB4NWbZlBQ42FeqOI4QQQpSozN8qq1evZuzYsYSHh3Ps2DEmTZpUmXWJSrT+tXt5c+dFAJoNf4TWHg00rkgIIURtU6GhwidNmsS6desqo54qUZeHCj+7dw23Tbyf6HRFvz5+HPz9EDqdTuuyhBBC1BBl/Q41y1wkNU1dDRg5iZHcNbADW8Ny8WzoxImz0bjWl/EuhBBClF2lzkUiaiBjIe/MGs7WsFysLHR8sXaLhAshhBCVRgJGHXFy7QJe3xIOwJx5TzJiyCCNKxJCCFGbScCoA4yXQmh75kM23m3PwFt8WfHWW1qXJIQQopaTgFHbFRhI/24GVhSgb3ELP+36XSYxE0IIUenK/JjqzJkzy7zR1atXV6gYYX5HPn2MhhfDKKhXj7gByxjibKd1SUIIIeqAMgeM3Nxc08+ZmZls2rQJLy8vfH19AQgODiY2NpaxY8eav0pRIdlh+7n35c+5lGlk/PS7WT24h9YlCSGEqCPKHDC++eYb088PPfQQ8+bNY/ny5VhaFm2ioKCAZ555huzsbPNXKcrPkMn/PTKZcylGGjjZMuXh57G2lFsjQgghqkaFxofetm0bJ0+eNIULAEtLS15++WW6dOlituJExR18/yHe/y0ZgP4zFzC6p0xQJ4QQoupU6FfajIwMoqKirll+8eJFMjIybroocXNyzv3Gg8t+RAGtu/Zg5YuPaV2SEEKIOqZCAeOuu+5i8uTJ/Pjjj0RGRnL+/Hl+/PFHJk+ezJQpU8xdoyiPAgML507hXIqR+o42TH/lI5q62mtdlRBCiDqmQrdIVq1axfz587nvvvswGAwA2NraMmvWLJYuXWrWAkX5/PXlsyzfexmAZmOf5Kk7fDWuSAghRF10U3OR5OTkcP78eQBatmyJnV3NeASy1s5FkhBG9sp+vLY/k/UprXnj212M922idVVCCCFqkbJ+h1boCsZVdnZ2dOzY8WY2IcxFKdj+DPYWBQy5tS8nGy1kbFcvrasSQghRR1X4ucUffviBYcOG0apVK9OyxYsXk5CQYJbCRPlkBH6HMfI3cpUVCwum8+IdHdHrZRp2IYQQ2qhQwPjss894/PHHGTBggOkWCUCDBg14/fXXzVacKKO8LB6aNYc+n2fxclx/2rfvxC0tXbWuSgghRB1WoYCxYsUK1q1bxyuvvFJs+e23387atWvNUpgou1/em8va4+kExRnZYBzE/FHttC5JCCFEHVehPhiRkZH06tULAJ3un8vw9evXJyUlxTyViTLJuxTK3DeLRllt392PSWOH07KRo8ZVCSGEqOsqdAXDy8uLU6dOAcUDxrZt22jdWkaMrEorn/kfZ5IKqe9ghfG2Z3l8qI/WJQkhhBAVCxizZ8/mgQceYO/eveh0Ok6dOsXbb7/N7NmzmTt3rrlrFKWI/2Mti9b/BUDjgf9j3ihfGjhYa1yVEEIIUcFbJE8//TQZGRmMHTuWwsJCOnXqhJ2dHc899xxz5swxd42iJAUGXnx6DukGaObRgAYDpjG9b3OtqxJCCCGACl7BOHPmDAsXLiQ5OZljx44RFBREUlISCxYsICwsrNzbi46OJigoiPT09HK1S0tLIyAggMjIyHLvs6bL++1d/rp4BQD9kCd5/vYO2FpZaFuUEEII8bcKBYz27dsDRcOD+/r60qNHD+zt7Yu9Vxa5ublMnDiRtm3bMm3aNNzd3Vm5cmWZ2iqluPfeexk0aBDvvfde+Q+iJku/hPWf73LoQQcm3PM/uvS5lTFdPLSuSgghhDC5qZE8/yszMxMHB4cyr79w4UICAwOJiIjAw8ODjRs3Mn78ePz8/Ojdu/d1277zzjum2zN1zm9vQn42x2nDUe97WDe6fbHOtkIIIYTWyhUwnnjiiRJ/BjAajRw/fpzu3buXeXv+/v7Mnj0bD4+i377HjRtHp06d8Pf3v27AOHr0KCtWrCAoKIiRI0eW5xBqvMKEs3zy0Sfc382SZcb/MaqTBz2bN9C6LCGEEKKYcgWM8PDwEn8GsLKyokePHsybN69M24qLiyM+Pp4ePXoUW+7n50dwcHCp7TIyMrj77rv54IMPcHd3L9O+DAaDadZXoNx9PaqTNQsfZO72HJYfsUc/sy2/yKBaQgghqqFyBYytW7cCcP/99/Pll1/e1I6vDsjl6lp8SGtXV9frDtY1e/Zshg4dytixY8u8r6VLl7Jw4cKKFVqN5EYe4ZVvfgfAutMI7u3bkmauZb8lJYQQQlSVCnXyvNlwAUVXPKCoo+e/5eTkYG1d8lgOmzdvZuvWrUyYMIGAgAACAgLIysoiLi6OgICAUvc1f/580tLSTH+io6Nvun4tfPzyTGLSFa5Otlj1ncFjQ2RQMyGEENVTma9gPProo2Xe6KpVq264jre3N3q9ntjY2GLLY2Njadq0aYltCgoK6NSpE4sWLTItu3TpEtnZ2cTFxXHgwAEsLK59VNPGxgYbG5sy118dZZzczes/hwBQv99dzBvRgXr2MqiWEEKI6qnMASMmJsasO7a3t6dv375s3ryZqVOnApCVlcXevXtZsGCBab3w8HAyMjLw9fVlwoQJTJgwodh2unXrxuDBg3n33XfNWl+1ohQrX5lLYraicX0H3AdN5b4+zbWuSgghhChVmQPGxo0bzb7zxYsXM2zYMObPn0+fPn1YuXIlbm5uPPzww6Z1li1bxqFDhzh58qTZ919TpJ/YwfIdRZ1qHfpN47lRHbG2rNDdLSGEEKJKaPotNWjQIPbv309UVBTvvfceHTt2JCAgAEfHf2YD9fHxue6jr76+vrRs2bIqytVM5i9vM7SlJZ4NnWncZwK3d5ZBtYQQQlRvOqWUqkjDH374gc8//5zz588TEREBFF2RePjhh3FzczNrkeaWnp6Oi4sLaWlpODs7a13O9V0IgC9Hk48lfTOWMWPMEGYPbqV1VUIIIeqosn6HVugKxmeffcbjjz/OgAEDOH/+vGl5gwYNeP311yuySVGaA28A8EPBYNLsmjCll7fGBQkhhBA3VqGAsWLFCtatW8crr7xSbPntt9/O2rVrzVKYgOwzB3jio92cTdXxUcGdjOniKdOxCyGEqBEqFDAiIyPp1asXQLE5MOrXr3/dQbJE+Xyx7GneO5zHoG8KiFWu3O0nVy+EEELUDBUKGF5eXpw6dQooHjC2bdtG69Yy+JM55F8OY/nPxwBw6TmGlo0c6dmsvsZVCSGEEGVToYAxe/ZsHnjgAfbu3YtOp+PUqVO8/fbbzJ49m7lz55q7xjrpx+VPEpWmaOBgTXanyUzq0URmTBVCCFFjVGi69qeffpqMjAzGjh1rmjLdzs6O5557jjlz5pi7xjpHZSay/Nu9AHj2uJUsaxsmdm+icVVCCCFE2VUoYOh0OhYuXMj8+fMJDQ3FaDTSvn177O3tzV1fnbTvk//j+OUC7Kz0XOl2H0PbNMLdxVbrsoQQQogyq1DAuMrW1hZfX19z1SIA8nNYsfpbADr49iDJzoW7ekrnTiGEEDVLuQLGyJEjy7Tezp07K1SMAPXXD3RxLeSwvQUJ3R6kgb0VQ9tX74HLhBBCiP8qV8DYtWsXzZo1o3///pVVT92mFLrDn7D0Nlu6THqE+ZebMK6bFzaW184QK4QQQlRn5QoYzz77LGvWrCEwMJAZM2Ywffp0PDxkXgyzufgnJIairOxZnugHwOSe0rlTCCFEzVOux1TffPNNoqOjWbZsGb/99hvNmzfnzjvvZNOmTRQUFFRWjXXG9k8WsCeigLMNh5NcaEdHT2c6erpoXZYQQghRbuUeB8PKyooJEyawfft2IiIi6NWrF0888QRNmshv2jdDZSTw9Gf7GP5NNo8dKLqwNLmHnFMhhBA1001N156WlsaVK1fIysrCxUV+074Zv325kLCkQhxtLAhrPBJLvY4xXT21LksIIYSokHIHjPT0dD799FNuueUWevfuTUpKChs2bODMmTOVUV/dYCzkky+LHk0d0Lcneht7BrVphKujjcaFCSGEEBVTrk6e06dPZ926dXTu3JkHH3yQu+++Gycnp8qqrc5ICtrI+r/SAEjscC8A43y9tCxJCCGEuCnlChhr1qyhWbNm1KtXj/Xr17N+/foS15NxMMrn6w/fJK8QurRoRKJjSxxtLBnWobHWZQkhhBAVVu4rGMK8VHYKn+84CkDbgWMJBEZ1csfWSsa+EEIIUXOVK2B8+eWXlVRG3ZV48CvyCxW2VjpOuY8AYHx3uT0ihBCiZrupp0jEzXO7uJWwuQ5sevdJsrDDw8WWW1q4al2WEEIIcVMkYGgp8QzEBqGzsGKXXdE8L3d280Sv12lcmBBCCHFzJGBoKO6XT8ktUOS3vI3N4fkATPCVwbWEEELUfBIwtKIUjy9dTePlGSz5qx75hYr2Hs60dZfHfoUQQtR85erkKcznSuhvbDmVTl4hHLbsAQrG+8rInUIIIWoHuYKhkZ8/X05eIbTzcuG00ROdDsZ2k6dHhBBC1A4SMLSgFGu3/QpAj74D0el09GvVkMbOttrWJYQQQpiJBAwNJJ/8hb1nMwGIaj4WgPEyNLgQQohaRAKGBjZ+8Q6FCjo1rU+03h1bKz0jOrlrXZYQQghhNhIwqprRyPqdvwHQoVd/AIZ3cMfRRvrbCiGEqD3kW62qxQbx9hBFHy8ndrUcD8CdXeXpESGEELWLXMGoamd20L6RBc8+OJ4YvRs6Hfi1bKB1VUIIIYRZScCoaud2AxDm3BeAdu7OONtaaVmREEIIYXZyi6QK5SWe56FPjnB7a0vOd2oHZNKreX2tyxJCCCHMTgJGFTr404es+SufXZE6eg4omtCsZ3O5PSKEEKL2kVskVWjrls0AjOzXhdOXMwDkCoYQQohaSQJGVcnPYduR8wB0GjQWowKvenZ4uNhpXJgQQghhfhIwqkjEbz9yLrkQSz3ktxkFgF8LuT0ihBCidpKAUUV2rv8agH4dvAhJygegp9weEUIIUUtJwKgiuw4GATD8tiEEX7wCQC/p4CmEEKKWkoBRBVRaLClpRZObtR48hZz8QlzsrGjdyFHjyoQQQojKIY+pVgHdhYMEzHDgskNHNtdrC5yhZ7P66PU6rUsTQgghKoVcwagKkUWTm7l3u42gqCsA9JIOnkIIIWoxCRiVTSkKzv1a9GOLQQRFpQLS/0IIIUTtJgGjkuVeOkOjV0K57esc/spvRkpWHrZWejp7uWhdmhBCCFFpqkXAiI6OJigoiPT09DKtX1hYSGhoKOfOnaOgoKCSq7s5gTu/50ounErWcSqtaFKzbt71sLasFqdeCCGEqBSafsvl5uYyceJE2rZty7Rp03B3d2flypXXbbNkyRK8vLyYOHEiw4cPp3nz5mzdurWKKi6/3/bvAWBg11am2yN+cntECCFELadpwFi4cCGBgYFEREQQGhrKd999x+OPP87hw4dLXL+wsJCcnBxOnz7N6dOniYyMZObMmUyZMoXLly9XcfVl81vQaQAGDhhAYGQKIB08hRBC1H6aBgx/f39mzpyJh4cHAOPGjaNTp074+/uXuL6FhQWLFy+mQYN/vqBnz55NdnY2x44dq5Kay6MwK5VD59MAaDdoHLFXcrDQ6+jeVEbwFEIIUbtpNg5GXFwc8fHx9OjRo9hyPz8/goODy7ydI0eOANCqVatS1zEYDBgMBtPrsvb1uFkn9/1ERh442ehJa9ARCKGjpzMONjL8iBBCiNpNsysYKSlFtwtcXV2LLXd1dTW9dyNJSUk89thj3HXXXbRt27bU9ZYuXYqLi4vpj7e3d8ULL4c/920H4JZ2HhyKLOp/0VtujwghhKgDNAsYVlZFT1Tk5uYWW56Tk4O1tfUN26elpTFy5Ejc3d1ZvXr1ddedP38+aWlppj/R0dEVL7wcPNUlRvtYMnxQXwLOJQHQ36dRlexbCCGE0JJm1+q9vb3R6/XExsYWWx4bG0vTpk2v2zY9PZ3hw4djYWHBzp07cXJyuu76NjY22NjY3HTN5WI0cqdbLHfeY0/MlKdZ9VUS1hZ6eYJECCFEnaDZFQx7e3v69u3L5s2bTcuysrLYu3cvw4YNMy0LDw8v1ifjargA2L17Ny4u1XTAqpQIyMsASzt+TW0IQI9m9bGzttC4MCGEEKLyadrbcPHixQwbNoz58+fTp08fVq5ciZubGw8//LBpnWXLlnHo0CFOnjxJfn4+o0aNIjIyki+++IKQkBDTej4+PjRu3FiLwyhR0umD5GUY8WzXkd/Ci/pf9PdpqHFVQgghRNXQNGAMGjSI/fv388EHHxAYGEjnzp35+uuvcXT8ZxpzHx8f8vLyAMjOzkan0+Hj48PSpUuLbeuFF17gjjvuqNL6r8f/6+94bnUmD9x2mSN9kwHo31oChhBCiLpB8+cl+/XrR79+/Up9//nnnzf97OLiQkBAQFWUddOOhYQBUN+rJRmGAtycbGT+ESGEEHWGTIhRGZQi+HwCAHmNuwAwtH1j9HqdllUJIYQQVUYCRiXIuhTO2cR8AIIsOgIwrIObliUJIYQQVUoCRiU4+ft2FODmZMklowv21hb0bSX9L4QQQtQdEjAqQcjRPwFo6VE05sXwDo2xtZLHU4UQQtQdEjAqwYmTRTOoOrgW3RYZ372JluUIIYQQVU7zp0hqo1Et9VjeYs3Zpl1wc7KhXyvXGzcSQgghahEJGOamFKM8khk1wpZhhoGM7+6FpYVcKBJCCFG3yDefuWVcAkM6BUpPtM6D+/o017oiIYQQosrJFQwziw/9k/CLBdg08GRo9yZ41bPTuiQhhBCiyknAMLOf1q/nMf9sujVPZu2zPlqXI4QQQmhCbpGY2fGQoidIPDw8aNP4+tPICyGEELWVBAwzuxBzCYA2beTqhRBCiLpLAoYZKaW4lJQGQNdu3bQtRgghhNCQBAwzik7OJDqlaGr5Hn6lzxArhBBC1HYSMMwo+PgxMvIUAG269tG4GiGEEEI7EjDMKPyvQwC4O1li6+CocTVCCCGEduQxVTOqb5HFylG2JNs21boUIYQQQlMSMMzI1TKLGX7WnPXopXUpQgghhKbkFokZWWXFA6B38tC4EiGEEEJbEjDMKDL8HH9GF2CwaaB1KUIIIYSmJGCYiVKKj7efpO8X2RyLvKJ1OUIIIYSmJGCYSYahgKQMAwCt2rTXuBohhBBCWxIwzORKZg6JmYUANG3dUeNqhBBCCG1JwDCT2PNnKFSgA9xbddK6HCGEEEJTEjDMJO78KQBc7fVY2thpXI0QQgihLQkYZhJ/MRIAV0drjSsRQgghtCcBw0wS42MBqO9oo3ElQgghhPZkJE8zaeNVj/dG2mBo2ELrUoQQQgjNScAwk2YNrPlfbxuC6rfTuhQhhBBCc3KLxEx0uVcAKLB20bYQIYQQohqQKxhmci78PIUpBWQ0sNK6FCGEEEJzcgXDTD7acpQB/tn8GRavdSlCCCGE5iRgmElWdtEw4Q71XDWuRAghhNCeBAwzycrNA8CpvpvGlQghhBDak4BhJtmGAgAc6jfSuBIhhBBCexIwzCTDYATAuV59jSsRQgghtCcBwwyUUmTnKQCcXBpoXI0QQgihPQkYZpCTk4P6+2cXuYIhhBBCyDgY5qDX61k+3Ja8AoWzSz2tyxFCCCE0JwHDDGytLHm6T9EsqlF2DhpXI4QQQmhPbpGYQ6HB9KONrZ2GhQghhBDVg1zBMIMrKUmExRTgYKWjia291uUIIYQQmpMrGGZwIiSEPp9nc8f32dja2GhdjhBCCKE5CRhmkJWVCYClXoeNpZxSIYQQolp8G0ZHRxMUFER6enqltqksOX8HDCsLPXq9TuNqhBBCCO1pGjByc3OZOHEibdu2Zdq0abi7u7Ny5Uqzt6lsuYZcAAkXQgghxN807eS5cOFCAgMDiYiIwMPDg40bNzJ+/Hj8/Pzo3bu32dpUNkNODiABQwghhLhK0ysY/v7+zJw5Ew8PDwDGjRtHp06d8Pf3N2ubylaQV3QFQ6eTgCGEEEKAhlcw4uLiiI+Pp0ePHsWW+/n5ERwcbLY2AAaDAYPhn7EqzN1vIz+/aKp2dNWiS4sQQgihOc2+EVNSUgBwdXUtttzV1dX0njnaACxduhQXFxfTH29v75sp/RpejRuwYJANd/WUqdqFEEII0PAKhpWVFVDUafPfcnJysLa2NlsbgPnz5/PUU0+ZXqenp5s1ZHTzG0hjw0Pg5G62bQohhBA1mWYBw9vbG71eT2xsbLHlsbGxNG3a1GxtAGxsbLCpxAGwmrXvQbP2PW68ohBCCFFHaHaLxN7enr59+7J582bTsqysLPbu3cuwYcNMy8LDw039K8raRgghhBDa0imllFY7P3DgAMOGDePpp5+mT58+rFy5ksjISI4fP46joyMAM2fO5NChQ5w8ebLMbW4kPT0dFxcX0tLScHZ2rrTjE0IIIWqbsn6HavrYw6BBg9i/fz9RUVG89957dOzYkYCAgGJBwcfHh+7du5erjRBCCCG0pekVDK3IFQwhhBCiYmrEFQwhhBBC1E4SMIQQQghhdhIwhBBCCGF2EjCEEEIIYXYSMIQQQghhdhIwhBBCCGF2mg0VrqWrT+aae1ZVIYQQora7+t15o1Eu6mTAyMjIADD7rKpCCCFEXZGRkYGLi0up79fJgbaMRiNxcXE4OTmh0+nMss2rM7RGR0fL4F0VIOfv5sk5vDly/m6enMObU1POn1KKjIwMPD090etL72lRJ69g6PV6mjRpUinbdnZ2rtYfjOpOzt/Nk3N4c+T83Tw5hzenJpy/6125uEo6eQohhBDC7CRgCCGEEMLsJGCYiY2NDa+++io2NjZal1Ijyfm7eXIOb46cv5sn5/Dm1LbzVyc7eQohhBCicskVDCGEEEKYnQQMIYQQQpidBAwhhBBCmJ0EjHIIDw/n6NGj5OTkVGqb2iozM5OgoCAuXLhQ5jaXL1/m+PHjMqz7386fP8/Ro0fJysoqV7vExEQCAgKIiYmppMpqhuzsbI4ePUpERES52sXFxREcHIzBYKikymqOCxcuEBQURGZmZpnWLywsJDw8nGPHjpGcnFzJ1VV/hYWFBAYGcvr06TK3ycnJ4ejRo4SHh1diZZVAiRtKTExUffv2VS4uLqp169bKxcVFbdiwwextajN/f3/l4OCg2rZtqxwcHNSIESNURkZGqev/+uuvqnfv3qpx48aqa9euys7OTj322GOqsLCwCquuPq5cuaJuvfVW5eTkpNq0aaOcnJzUt99+W6a2eXl5qnfv3kqv16ulS5dWcqXV1w8//KCcnZ2Vj4+PcnJyUoMGDVIpKSnXbRMfH69GjBihnJycVM+ePVXz5s3Vxo0bq6ji6iUjI0ONGDGi2P/H/v7+123z559/qlatWilPT0/l6+ur7Ozs1LRp01ReXl7VFF2N5OTkqEWLFqlmzZopZ2dnNXr06DK127BhQ7Hvkb59+6rExMRKrtY8JGCUwaRJk1T37t1VZmamUkqpt956S9nZ2amYmBiztqmtTp06pSwsLNSaNWuUUkolJSWpVq1aqdmzZ5faZvXq1SowMND0+sSJE8rR0VG99957lV5vdXT//ferDh06qCtXriillProo4+UlZWVCg8Pv2HbZ599Vk2fPl01bty4zgaM8+fPK2tra/Xhhx8qpZRKS0tTHTp0UNOmTSu1TUFBgerZs6caPHiwSktLU0oplZqaqr755psqqbm6mTVrlvLx8VHJyclKqaJfGiwsLNTp06dLbdOhQwc1efJk0y8GYWFhysbGRn322WdVUnN1cvnyZfXiiy+qqKgoNWXKlDIFjOjoaGVnZ6dWrFihlCoKeV27dlWTJ0+u7HLNQgLGDaSkpCgLC4ti/6gYDAbl4uKi3nrrLbO1qc2ee+451bx582LLli9frhwdHcv1m8zIkSNrzP9Y5pSdna1sbW3Vxx9/bFpWWFioGjdurF599dXrtt25c6dq1aqVSk9Pr9MBY9GiRcrNza3YFbCPP/5Y2djYmH4J+K/169crQJ05c6aqyqy28vLylKOjo3r33XeLLW/atKl6/vnnS23XqFEjtXz58mLLvLy81Ouvv14pddYUZQ0Yb775pqpfv77Kz883Lfvyyy+VpaWlSk1NrcQKzUP6YNxASEgIhYWF9OjRw7TM2tqarl27EhwcbLY2tVlwcHCxcwHg5+dHZmZmme8pGgwGQkJCaN26dWWUWK2FhoaSm5tb7Bzq9Xp69Ohx3c9TfHw8M2bMYM2aNTg5OVVFqdVWcHAwvr6+xSZm8vPzw2AwlHov/JdffqFt27a0adOG06dPExoaSl5eXlWVXK2cO3eOzMzMa/4/7tmz53U/g0uWLGHlypWsWbOGPXv2MHfuXJydnZkxY0Zll1wrBAcH06VLFywt/5k2zM/Pj4KCAkJCQjSsrGzq5GRn5ZGSkgKAq6trseWurq6m98zRpjZLSUmhRYsWxZZdPTdlPR9PP/00OTk5zJ071+z1VXfX+zxFRkaW2EYpxbRp03jooYfo27dvpddY3aWkpODl5VVs2Y0+g3FxcTg7OzNw4ECSkpLIzc0lKyuLTz75hHHjxlV2ydXK9T6DFy9eLLXdyJEj2bBhA88//zzu7u5ERUWxcOFCGjduXKn11hYpKSklnvOr71V3cgXjBqysrADIzc0ttjwnJwdra2uztanNrKysSjwXQJnOx5IlS/D392fDhg3XfEnUBRX5PH322WecPn2awYMHExAQQEBAAPn5+URFRXH48OFKr7m6qchn0MrKiiNHjvDAAw9w+vRpzp8/z6xZs5g6dSoJCQmVXnN1UpHPYEFBAUOHDqVx48ZER0cTHBzM4cOHeeWVV1i1alWl11wb3Oy/nVqTgHEDzZo1AyA2NrbY8tjYWJo2bWq2NrVZs2bNSjwXwA3Px7Jly1iyZAlbtmxh0KBBlVZjdVaRz5NOp6N58+a89NJLvPDCC7zwwgtkZGSwY8cOXnvttUqvubqpyGewefPm2Nracv/995uWzZo1i6ysLI4dO1ZptVZHFfkMnjlzhnPnzvHII4+YLvH7+PgwfPhwNm/eXLkF1xI3829ndSAB4wY6dOiAp6dnsf8hzp8/T0hICMOGDTMtCwkJITQ0tFxt6ophw4YREBBAamqqadmmTZvo3Lmz6VJpeno6AQEBxZ6tf/PNN1m0aBFbtmxhyJAhVV53ddG8eXNat25d7PN06dIlAgMDi32eTp8+zcmTJwF46KGHTFcurv5p0KABs2bNYuvWrVV+DFobNmwYQUFBXLp0ybRs06ZNtGjRglatWgGQlZVFQEAAaWlpAIwYMQKDwUBiYqKpzdVxRBo1alSF1WvP3d2dTp06FfsMJicn8/vvvxf7DIaFhZn6Blw9R/8deyU6OrrOnb+yysnJISAggCtXrgBFn9sTJ04QFRVlWmfTpk14eXnRvn17jaosB617mdYEX331lbKyslLLly9X69evV926dVP9+/cv1iN90KBBauzYseVqU1fk5uaqTp06qf79+6uff/5ZLViwQFlYWKgtW7aY1jl48KAC1JEjR5RSSn3wwQcKUAsXLlQHDx40/Tlx4oRWh6GpdevWKUtLS7VkyRK1YcMG5efnp7p3716sd/no0aPV0KFDS91GXX6KpKCgQPXq1Uv16tVLbdiwQS1dulRZWlqqtWvXmtYJDg5WgNq/f79p2e2336769eunNm3apNauXas6dOighg0bpoxGowZHoa0tW7YoCwsLtXDhQvXzzz+rfv36qc6dOyuDwWBaZ+LEiapfv36m11OmTFGenp5q9erVaufOnWr27NnKwsJCHTx4UItD0NyhQ4fUwYMH1dChQ1Xfvn3VwYMH1R9//GF6PzQ0VAFqx44dSqmip8X69eunfH191fr169Vbb72lLC0t1VdffaXVIZSLdPIsg/vuu4969erx1VdfsXPnTsaMGcOzzz5brEd6586di/XUL0ubusLGxoYDBw6wbNkyVq5ciaurK7t27WLo0KGmdVxcXOjXr5/pHEZFRdGvXz92797N7t27Tet16dKFDz/8sMqPQWsTJ05k69atrF69mn379jFkyBCef/75Yr3LO3bseN2nHHr37o23t3dVlFvtWFhYsGfPHt58800++OAD6tWrx+bNmxk1apRpHUdHR/r164eLi4tp2fr163nvvff44IMPsLe3Z+bMmcyZMwedTqfFYWjqjjvuYNeuXXzyySccOHCAvn378sILLxTrC9C+fXvc3NxMr7/++ms+//xzdu3axZUrV2jZsiVBQUF069ZNgyPQ3muvvWa6OqHT6XjhhRewt7c3/Rtnb29Pv379qF+/PlD0tNiOHTt46623+Oijj3B2dmb9+vXceeedWh1Cuch07UIIIYQwu7r367QQQgghKp0EDCGEEEKYnQQMIYQQQpidBAwhhBBCmJ0EDCGEEEKYnQQMIYQQQpidBAwhhBBCmJ0EDCGq2OXLl/nxxx+1LuOGjh49yk8//cSRI0fKtP5/j0vL4yzLvqt7fULUdDLQlhBmkpyczJ49e667To8ePYiIiOCOO+6goKCgiiq7lsFg4Oeff2bChAklzso4d+5cNm3aRN++fbntttt4+OGHb7jNnTt3Fjuu/76uSv/dd1xcHL///juTJ08udZ3KouW+hdCSDBUuhJmkpqayceNG0+ugoCBSUlIYPny4aVmDBg3w8PBgypQpGlT4j3379jF//nzuvvvua94zGo18/vnnbNmypdZMznfs2DGmTZtW7Eu+qv4etNy3EFqSgCGEmbRu3ZoffvjB9HrWrFkEBQUVWwZFl8fHjh1reh0TE0NgYCDjxo0jJCSEqKgounTpQvPmzSksLOTw4cOkpKTQq1cv0+yz/3bp0iWCgoJwcnKie/fuODs737DWrVu3cscdd1yzPD4+nm3btmEwGDhy5AjJycnceuutZGZmEhYWxujRo03rJiQksG/fPqZMmVKhuTmq6riTk5M5ePAgRqPR9HfRtm1bPDw8Svx7GD9+PCEhIVy8eJGOHTvSokWLa7YZGBhIYmIiHTt2xNHRkb1795Z4Hsq776r8DAhR6bSda02I2uuRRx5RPXr0uGb5jh07lIWFhen1Tz/9pOzt7ZWfn5/q37+/GjhwoLKwsFArV65UvXr1UgMHDlT9+vVTzs7OKjAwsNi2Fi1apFxcXNSIESPUgAEDVKNGjdSePXtuWFuzZs3Uzp07r1keEhKiJkyYoAA1fPhwNWXKFBUSEqI++ugj1apVq2Lr7t+/XwGmGV3/e1z/ff1flXnc/973uXPn1IABA5Rer1dTpkxRU6ZMUV999VWpfw+33Xab6tOnjxo2bJiytrZWX3zxhWmdwsJCNX78eOXi4qJGjhypmjVrpkaOHFnsPPxbefddlZ8BISqbBAwhKkl5AgagPvroI9OyBx98UAHK39/ftGzKlClq/PjxptebN29Wbm5uKioqyrRs9erVyt3dXeXk5JRa14kTJ5SDg4PKzc0t8f3U1FQFqCNHjpiWVVbAqKzj/u++t2zZomxsbIrtv7S/h1WrVpmWLV++XLm5uZlef/fdd8rBwUGFh4crpZTKzMxUvr6+pQaM8u67qj4DQlQFeYpEiGrAwsKCBx980PS6T58+2Nract999xVbdvbsWdNrf39/OnbsSGBgID/99BM//vgjlpaWXL58mdDQ0FL3tXXrVoYNG4aNjU3lHEw5VOVxl7Wef3doHTx4MAkJCaYpttetW8f48eNp1aoVAA4ODsydO/em9vnvfVencyHEzZI+GEJUA46OjlhZWZle29jY4OLigl6vL7YsNzfX9PrChQsYDAbWrVtXbFtTpkwp1u6/tm7dyowZM8xYfcVV5XFXtB7AtP/o6GhGjBhRrE3z5s1vap/X27eW50KImyUBQ4gaytnZmaZNm7JmzZoyt0lOTiYwMJD169eXa196vR6j0Vhs2b+/6KpSRY7bXBo0aGC6mnFVampqlddxlZbnQogbkYgrRA01cuRINm3axKVLl4otj42NLbXN9u3b6datG+7u7uXal5eXF5cuXSIrK8u0bP/+/eUr2EwqctyOjo7k5+dfE5LKq1+/fmzbtq3Y+BWbNm26bhtz7bskFTkXQlQVuYIhRA01b948tm7dSu/evZkzZw4NGjTg2LFj/Pbbb5w+fbrENqU9nnojt956Kw0bNmTy5MlMmjSJ48ePX3NZvqpU5Lg7duyItbU1zz//PD169KBt27YV2vdjjz3GJ598wogRI5gyZQpBQUHs3LkToNRHdc2175JU5FwIUVXkCoYQlaRXr17X3K+HawdZ8vb2ZuLEicXWad68OePGjSu2rHXr1sXCgZ2dHfv37+f1118nMjKSY8eO0bNnT44dO1ZiPQUFBezevfuGAcPa2popU6bQoEED0zJ7e3sOHTpEt27dOHToEO3atWP79u3F7vX/97huNJhUZR73f/fdqFEj9u7di8FgYNOmTYSEhJTp76FevXpMmTIFOzs70+vAwED8/PwIDAykY8eOrFy5Ejs7OywsLEo8zoruuzI+A0JUJRkqXIg6Yv/+/dx7773ExsZWaGAsUSQlJaVY+Jo1axanTp3i4MGDGlYlRPUjt0iEqCNSUlJ44YUXJFzcpMmTJ9OnTx+aN2/OwYMH+emnn4oNES+EKCJXMIQQohxSUlJYvXo1Z86coUmTJtxzzz1m7VchRG0hAUMIIYQQZiedPIUQQghhdhIwhBBCCGF2EjCEEEIIYXYSMIQQQghhdhIwhBBCCGF2EjCEEEIIYXYSMIQQQghhdhIwhBBCCGF2EjCEEEIIYXb/DxDxvfVz5i0iAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "class SurfaceAt(BaseSubModel):\n", " \"\"\"Holds the PCM surface at a temperature, and sets the initial temperature\"\"\"\n", "\n", " def __init__(self, param, T_surface):\n", " self.T_surface = T_surface\n", " super().__init__(param, \"Surface\")\n", "\n", " def build(self):\n", " self.variables = {\n", " \"Phase-change material boundary temperature [K]\": pybamm.PrimaryBroadcast(\n", " self.T_surface, \"flow path\"\n", " ),\n", " INITIAL_TEMPERATURE: self.param.T_init,\n", " }\n", "\n", "\n", "capsule = EnthalpyStorage(param, param.pcm, geometry, inner_condition=\"temperature\")\n", "melting = pytes.ComposedModel(\n", " {\"storage\": capsule, \"surface\": SurfaceAt(param, param.T_m + 1)},\n", " param,\n", " name=\"Melting capsule\",\n", " defaults_from=base,\n", ")\n", "pv = base.default_parameter_values\n", "pv[\"PCM smoothing width [K]\"] = 0\n", "pv[\"Initial temperature [K]\"] = pv[\"PCM melting temperature [K]\"]\n", "tau = 861 * 213e3 * R**2 / (6 * 0.15 * 1)\n", "t = np.linspace(0, 1.05 * tau, 2001)\n", "fig, ax = plt.subplots(figsize=(6, 4))\n", "for n_r in (20, 160):\n", " solution = pybamm.Simulation(\n", " melting, parameter_values=pv, var_pts={base.r: n_r, base.z: 1}\n", " ).solve([0, t[-1]], t_interp=t)\n", " ax.plot(t / tau, solution[\"Average liquid fraction\"].entries, label=f\"{n_r} cells\")\n", "melted = np.linspace(0, 1, 101)\n", "ax.plot(validation.sphere_melting_time(melted), melted, \"k--\", label=\"quasi-steady\")\n", "ax.set_xlabel(\"Time / full melting time\")\n", "ax.set_ylabel(\"Melted fraction\")\n", "ax.legend()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "a40ee977", "metadata": {}, "source": [ "Both meshes follow the quasi-steady solution. Early on, a coarse mesh deviates slightly:\n", "the melting front crosses the large outer cells of the sphere one at a time (the enthalpy\n", "method with a sharp melting point). Refined, the model is within 1 % of the quasi-steady\n", "solution from half melted onwards." ] } ], "metadata": { "kernelspec": { "display_name": "pytes (3.14.4)", "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.14.4" } }, "nbformat": 4, "nbformat_minor": 5 }