diff --git a/.gitignore b/.gitignore index b130bb9287..a2ae0d0425 100644 --- a/.gitignore +++ b/.gitignore @@ -59,6 +59,7 @@ htmlcov/ pyproject.toml # virtual enviroment +env/ venv/ venv3.5/ PyBaMM-env/ diff --git a/CHANGELOG.md b/CHANGELOG.md index f54e2a4d5a..9600e4e69c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,4 +1,7 @@ # [Unreleased](https://github.com/pybamm-team/PyBaMM) +## Features + +- Add ambient temperature as a function of time ([#872](https://github.com/pybamm-team/PyBaMM/pull/872)) ## Bug fixes diff --git a/examples/notebooks/create-model.ipynb b/examples/notebooks/create-model.ipynb index 0ae441220b..cb1f50f01e 100644 --- a/examples/notebooks/create-model.ipynb +++ b/examples/notebooks/create-model.ipynb @@ -2,7 +2,10 @@ "cells": [ { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "# Creating a Simple Model\n", "Before adding a new model, please read the [contribution guidelines](https://github.com/pybamm-team/PyBaMM/blob/master/CONTRIBUTING.md)" @@ -10,61 +13,80 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "In this notebook, we will run through the steps involved in creating a new model within pybamm. We will then solve and plot the outputs of the model. We have choosen to implement a very simple model of SEI growth. We first give a brief derivation of the model and discuss how to nondimensionalise the model so that we can show the full process of model conception to solution within a single notebook. " ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "Note: if you run the entire notebook and then try to evaluate the ealier cells, you will likely recieve an error. This is becuase the state of objects is mutated as it is passed through various processing. In this case, we recommend that you restart the Kernal and then evaluate cells in turn through the notebook. " ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "## A Simple Model of Solid Electrolyte Interphase (SEI) Growth" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "The SEI is a porous layer that forms on the surfaces of negative electrode particles from the products of electrochemical reactions which consume lithium and electrolyte solvents. In the first few cycles of use, a lithium-ion battery loses a large amount of capacity; this is generally attributed to lithium being consumed to produce SEI. However, after a few cycles, the rate of capacity loss slows at a rate often (but not always) reported to scale with the square root of time. SEI growth is therefore often considered to be limited in some way by a diffusion process." ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "### Dimensional Model" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "We shall first state our model in dimensional form but to enter the model in pybamm, we strongly recommend converting models into dimensionless form. The main reason for this is that dimensionless models are typically better conditioned than dimensional models and so several digits of accuracy can be gained. To distinguish between the dimensional and dimensionless models, we shall always employ a superscript $*$ on dimensional variables. " ] }, { - "attachments": { - "SEI.png": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAACGQAAAGRCAYAAAApCfcPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAABM5QAATOUBdc7wlQAAABl0RVh0U29mdHdhcmUAd3d3Lmlua3NjYXBlLm9yZ5vuPBoAACAASURBVHic7N13YFvlvf/x75Esyxree8R2YmfhTGeQQBZhhU0ZDZRLaYEW2tsftNALLR20hV5ob0tvS8u6UKClYTWhpGkLBEoYCRkEEmeHJI7jvWTLQ7KscX5/0NAQfI5k50iWnffrP/IcneeRJcviPJ/z/SryL5dffvlvVVX9TwEAAIhzD/7i0MtJ1qBtuNcBAAAAAACi77ebbq3+8+7PjR3udQAAAIRlUm/dcVfW3k/+czjXAgAAAAAAAAAAAAAAMBoRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwGIEMAAAAAAAAAAAAAAAAgxHIAAAAAAAAAAAAAAAAMBiBDAAAAAAAAAAAAAAAAIMRyAAAAAAAAAAAAAAAADAYgQwAAAAAAAAAAAAAAACDEcgAAAAAAAAAAAAAAAAwWMJwLwAjz/z588XpdA445na7ZfPmzTFeEYB4l5ycLPPmzdMc37t3r9TW1sZwRQAAAAAAAAAAAEB0EcjAoF155ZVSUlIy4Nj+/fsJZAD4jMzMTLn55ps1xx9//HECGQAAAAAAAAAAABhVCGQAAAAAAAAAAADghKh+r6n+b98p0xp3jl3QkTb9yrZYrglA5JpeuXtsoK9zwL1ja+Y4T/bCW+tjvSZgNCCQEUVWq1Xuueee4V7GgPx+v3zve98b7mUAAAAAAAAAAIBRINTvNXW8/8xY7QOCCoEMIH51bHuxKNDTbBtozF40q51ABjA0BDKiyGQySVmZZhh0WPl8vuFeAgAAAAAAAAAAAAAAoxaBDAADOu+88yQzM3PAse7ubnn55ZdjvCIAAAAAAAAAGB3cO1/OrFv1/6YN9zoGYi+a6Rp7/cvbh3sdAACMBgQyAAxo8eLFMn78+AHHGhsbCWQAAAAAAAAAwFCF/Eqov9cy3MsYSKjfy94RAAAG4Y8qAAAAAAAAAAAAgLDUoF85/IflmtVdHKXzXDln3FEbyzUBsdK5/cWsjg+eLdQazz//no+Scis8sVwT4h+BDAAAAAAAAAAAAADhqUGl5+BbuVrD5kR7IJbLAWLJ17rPoff+D/S0HhbNUZysTMO9AAAAAAAAAAAAAAAAgNGGChnD7Nlnn5X9+/fHfN5QKBTzOQEAAAAA8eFgtS/pjX/2pO/e2+c8VO1zNLcEkjpcAavHq5q9faEEERGLRQklJipBp8MUsFiUUFaWxVdUaPEUj7F4x5ZYvZMnJfXOnmnrsSQq6nA/HwAAAAAAgHhEIGOYHT58WKqqqoZ7GQAAAACAUW7rh17nM8+68t94szu3vsFvD3e836+a/H7V1NsbsoiI1NX7Hdu2S8axx1itSnB8eVLX1Iok9+xKu3vZ2cmuvDyL38h1X3LloRkfbvemG3lOo/zuf4u2X3R+qmuwj9v6odf5ueWH5ugd85cXxm2unGHrHfrqAAAAAADAcCOQAQAAAADAKPba691pv/x1S9nWDz2ZqsG1LHw+1bxzlzd95y5v+rMvdMh/3SXqpIk29+KFjtaLzkttO3WuvftE5/B6VfPRUEi88ftFGcrjQkFVwj2nYEAd0rkBAMDoUXbjmvfMtrRArOc1We2U2AYAwCAEMgAAAAAAGIXqG/yJt367btK6t3vyjA5iaAmFRNm9x5u2e4837eHH2sbfcF3mwZ//d8HB2MwOAAAwulizx3vNjqyYBzIAAIBxCGQAAAAAADDKrPlHV8at366b2tkZtA7nOvr9qmk45wcAAAAAABhOXBgBAAAAAGAUefzJ9rwbbj5SOdxhDAAAAAAAgJMdFTIAAAAAABglVjzfkf3dHzZMDYVEGexjFUUkJdncb7MpQafT7DeZRPV6Q+aenpClpzdk8VPtAgAAAAAAYFAIZGBUUBRFcnNzJScnR7KyssRut4vFYhERkf7+funt7ZWWlhZpbGyUjo6OYV5tZJKTkyUvL0+ys7PF4XCI0+kUEZFQKCQej0c8Ho+43W6pra0Vt9s9zKuNb2azWXJzc6WoqEhSU1PFZrOJ2WyWvr4+cbvd0tTUJHV1ddLf32/IfFarVfLz8yUnJ0dSU1M/ee1ERDwej3R1dUlzc7M0NjaK1+s1ZM5oMJlMUlBQIGPGjJG0tDRJSkoSERGv1ysej0fq6+ulrq5OfD7fMK906JxO5yevldPpFLvd/slYT0+PdHR0SHNzszQ0NEgwGBzGlQIAAIS3vcrr+PZ366dFGsYwmxW1cqat/ewzUlpnzrR1zZxu601PM2v2KHe1BxM2bO5N2fqBJ2XHLm/K7j19qc0tAZtxzwAAAAAAAGB0IZCBESsrK0tOP/10mTp1qkycOPFTG6l6WltbZdeuXfLee+/Jtm3bJBDQvN4YUw6HQ2bPni2VlZUyYcIEycnJifix3d3dUlNTIzt37pSqqio5cOBAxJvH+fn5MmXKlM/8e2pqquZj7Ha7nH322RGv73jV1dVy4MCBsMdVVlZKZmbmgGN+v1/WrVun+djk5GRZuHChVFZWypQpUz4J6Gh54IEHZP369WHXpKWkpEQWLlwoU6ZMkXHjxonZbA77mFAoJDU1NbJjxw5Zv359RD+TaFMURWbNmiWLFy+WadOmfSpMMpBQKCT79u2TLVu2yLp16+I+HJSQkCAzZ86UuXPnyqRJk6SgoCCix/l8Pjlw4IBs2bJFNmzYIO3t7VFeKQAAwOAEg6LcfEvtFJ9PDftF1GpVgtf9R2b1t76RU5uTY/ZHOkdGpjlw4XkprgvPS3Ed/bct73ucf36pM3ftP7tza4706395jJLFC5zNP/1xwf7hmLt4TOLITScDAAAAAICoI5CBEWf69Oly6aWXytSpU0VRBl2FV7Kzs2XJkiWyZMkS6ezslL/+9a/yyiuvSF9fXxRWG9748ePloosuklNPPVUSEob2K5mcnCxTpkyRKVOmyFVXXSUej0fWr18vr7/+ethN/okTJ8rNN988qPlSU1MH/ZhjrVy5MqLwwYUXXijTp08fcKyrq2vAQEZaWppcccUVsnTpUrFaI2+ZPZT3kqIocuqpp8pll10mZWVlg368yWSSsWPHytixY+Xiiy+W2tpaWbVqlaxfvz7m1RgURZH58+fLNddcI3l5eRE/zmQyyeTJk2Xy5Mly1VVXyRtvvCHPPfec9PT0RHG1g5eeni4XXnihnHXWWWFDJgOxWq1SUVEhFRUVct1118nGjRtl5cqVUl1dHYXVAgAADN4TT7Xn7f/Ip52q/pdZM+3tj/2ueGdpicWQIMGc2faeObPtPT8TObhjV5/98d+3Fb68xl3U3RPST0QbyOE0ByZPssZv6TkAAAAMK1/LXpunYbvT315t83c3WUP9HrMa9JtERExWZ8BsT/Nb00u9SfnTeuxFlb1iMqvDvWY9QW9HgufIFqev/aDd33EkKej3mEPeLouIiGKxhkwWR8Cc5AwkphX3JeVP6bUVTO9VEpJCw73u+KSKt77K0ddU5fC1HnAE+joTQn3dCaKqismaHEhMze+zZI7zpk29rE0xW074faEG/Yq3douzr3mPw9dx2Bb0uC0hX3eCiIhitoRMVkfQ4sz1JWaVeRxj5nRbMkoJgMMwqr/X1FuzKdnX8pG9v+OwLejrSQj5ev71/ksMmazOQGJqYV9i1livs+T0LnNyTsQ3cEAfgQyMGCUlJXLTTTfJxIkTDTtnWlqaXHvttXL++efL448/Lps3bzbs3OEUFBTI9ddfLzNnzjT83EcrWJx99tlSXV0td95550nRbmHp0qVy3XXXDWnDfbCBjNLSUvna174m5eXlg55Ly5gxY+TWW2+VSy+9VB599FHZt2+fYefWk5aWJrfccotm+CVSiYmJct5558lpp50mjz32mGzcuNGgFZ7Ymi6//HK55JJLwlZJidTR8Mq8efPk1VdflT/96U/i8XgMOTcAAMBQPfx/bWPDHbN0SXLTH58o2ZmUpETlYuzUiiTPr39Z9NF//6Tw4ONPtuX/8VlXcTTmAQAAAPSofq+ps+rPWZ07Xs7z1m9LD3o7Ir5zz5RoD9gKK11pUy5uSq+8qkWxOIY/yBAKKN3716Z17nw511OzKaO/s9YhqhrxBW3FZFYtaUW9jpJ57c7yM9pTJi3rMFmTI9owcL3/h1xP3daUT68nqDu3t3l3at1fbh0f6fqOlznn+gZb4fTecMfVr/52uRryD7iWpJxJvVmnfa1B67HuXX/NcG97Pr/n8HtZQW9n2PeHvWj229assiHd1RvsbUtwffCn3K7df8/ra9qZFvL3hS+v/S8Jzhyvo2Ree/qMzzcmTzynM97DQke1rnugyNdZo9neMnPOdQ22wsqwr7FhQkGlYc13x4VCPpPWIfnn/uSQ2ZYaPPYxdatv+8wGkK9xp+6NEG0bHh3TueOlyEvgHyPBlubPO/fHh4fyWC19jVX2jg9W5HcffDurv+1Aihrm9/ffFEnMKO5xli9tyai8uslWNCu+7sIdYQhkIO4piiJXXnmlXHHFFRG1gRiKzMxMufPOO+Xvf/+7PPXUU1ENLyiKIpdeeqlcffXVUXs+xxo7dqyYTKZRHcgwmUxy/fXXy3nnnTfkcwwmkHHRRRfJtddeG7XXr6SkRH7605/KCy+8IC+++KKoavS+Y5WVlcldd90laWlphp0zNTVVvv3tb8uKFStk1apVhp13sMaNGye33377oCp+DIaiKLJs2TKZPn26/OxnP5Pa2tqozAMAABDOm2/1pB6p1W8XUlKc2BPNMMaxHA4ldOs3sutv/UZ2fUODPzHa8wEAAAAiIn53Q2LLWw8Ud25/cUzI1zOku7NC/Z6E3up3c3qr381pWvvT/vTKL9TknnnnkUgDDEYKejsSWt/5TVHHB8+PCfQ0a25uh6OGgkq/q8bZ76pxdnz4fIlitoSc5UuaM+dcV5c86bwOvcf2fPRmpnvX6qLBzHd0rqGuN3ncwo5IAhkdH6woUQO+AS/SO8YuaBkokNGx7fnsljd/Ud7ffih5cKsa/DV6X9vBpJZ1/1Pq3rm6UGud4QR6WmzuXauL3LtWF1nSCnuz5n+1OnPeTY1GVOuIpqDfY+54/xnNmwZC3k5L8dVP747Vety7V2e0b35Cs8y5NWeS22xL/ejYf1PVkOg9By3d+16LrE/6ABKcuV5jAhmquKtWZrW++1Cpt2F7xlDP0e+qcbo2P+l0bX5qnH3M7LbcM+846CxfGt996+MUgQzENZvNJrfddptUVlbGZL7zzz9f8vLy5Oc//7n4/cZX4on18zkZKIoi3/jGN2Tx4sUnfJ5Ijvna174mZ5555gnNFel6li9fLsXFxfKrX/0qKoGayZMny1133SV2u93wcyuKItdcc42YzWZ58cUXDT9/OIsWLZKvf/3rhlXF0JOfny/33Xef3HvvvbJ3796ozwcAAHC81X9zh7375hf3F+6KRRjjeAUFlv5YzwkAAICTixr0Ky1v/qy4bf0jZSG/17B9r2CfO7Ftw8PjO6tWFRWc/5M9qdOuaDPq3LpCQaX1nV8Xtr7zYHmwr8vwgLMa9Ju6963N7963Nt9ROq913I1/+9DoOeJNsLfVcuSFm07pOfhW7pBOEHFVAREJBZSmN+4vblv/UPlQgxgD8XfWOxr/cfeU9s1PlxRe/MAuZ9nCLqPObbTs025qaHv3d+PVYP+AFSm697+eH/R27Dfb0gOxWI9ryx91Q0XpM6+ui8U6YsFbv81R//K3JnsbqoYYxBiIKp7aLVnVT30+K7XigrrCi3+13+zIiMlrN1polmYBhpvNZpMf/OAHQw4vqKo6pE3syspKueOOOwyvfpCcnCx33333CYUxAoHAqK50MRTLly8/4TCGSPhAhqIo8vWvf33IYQyfb2it3ubPny+33Xab4e/HwsJC+c53vjPoMIbX65W2tjbp6OiIKLS0fPlyWbBgwVCXOSRnn3223HLLLUMOYwQCg/8eYbPZ5Pvf/76MHz/kSnwAAABDtvUDb7re+NjSxO6li53cxQIAAIBRx9d2IOmj3y2e27LuVxONDGMcK9DTbDvyws2V9au/XT6ojfkh8LUdTDrw8Jmzm9b+9JRohDGOF+h1jfqKdr6WvbaPHj577pDDGCISaYUMf1eT5cAj58xqfetXE40MYxyrv/1Q8uGnLz+15fX/jtsWkWZHtj9lwpmNWuMhf5/ZtenJ6JS2Po6/szaxt/pdzZsYTImOQObc65pisZZoa3nj/uKDj54739gwxrFUce9aU/TRw0tO9dZ/4IjOHKMTFTIQlxISEuR73/ueTJw4MaLja2trpaqqSnbv3i1HjhyR7u5u6e7uFpGPgxCpqakyceJEqayslJkzZ4rVqt8SrLKyUq677jr5/e9/f8LPRUQkMTFRvve970W8WauqquzZs0e2bdsm+/fvl8bGRmlvb/+kdYXNZhOn0ykZGRlSXl4u5eXlMmHChCG1Zjhy5Ii89NJLn/n3JUuWSHr6wNd1e3p6ZO3atYOe66jdu42pRDV58mS54oorNMf7+vpk9+7dUl1dLW63W7q6usTj8YjNZpPs7GwpLS2VadOmSUpKSthAxhVXXCFLly6NaF09PT3yzjvvyLZt2+Tw4cPS1vbv4HR6erqUlpbKjBkzZNGiRZKSkqJzpo/NmzfP0PdjUlKS3HXXXeJ0hq8a19/fL++++65s3rxZ9u3bJ11d/w7dKooiubm5MnXqVFm4cKFUVFR85vGKoshNN90kv/71rw1ZeziLFi2Sm266KaKKJx6P55PPjT179khnZ6d0d3eL3++XxMRESUlJkcLCQpk5c6bMmTMn7O+XzWaTO+64Q+68805xuVxGPSUAAICwqmt8ul/szjojpSVWawEAAABipXv/62m1L3xlRiyCCyKquDY/Oc7vrk8qveaZXWIyG94yonv/a2m1L3x1ZrCvO/plf08Sftdh66EnLp0T6G1NivZcvpa9tkO/v2z2ibSXiZQaCirN6345ydd+0DFm+eN7RKKaExqSjHk31rn3/KNQa9y17fmi7CW3Rb0yRdvGxwpVnSBVyqRlDcPRkshIatCv1Kz4YsWJtEsZDH9nvePQE5fOLbl2xVbn2AVxW6klnhDIQFy66aabZPLkyWGP2759u7z00kuyY8cOzWOOhjPq6urkjTfekIyMDFm+fLksXbpUTCbtIjEXXHCB7Ny5UzZv3jyk53CUoihyyy23RBTGCAQCsnbtWlm9erW0tGhfM/V6veL1eqW1tVX27dv3yb+XlJTI6aefLqeddprk5+dHtL5Dhw7JoUOHPvPvU6ZM0QxkdHd3yzPPPBPR+aPFbDbL17/+9QE33hsaGmTVqlWyYcOGsJUpFEWR2bNnS2+vdku66dOny/Lly8OuyefzycqVK2XNmjWa83Z0dEhHR4d8+OGH8swzz8h5550nn//858Vm0/+OdsEFF8i+fftk/fr1YdcRzg033BBReOfNN9+UZ555Rjo7OwccV1VVmpqapKmpSdauXSuTJk2Sr3zlK1JaWvqp4+x2u3z5y18+4XWHM378eM33xLE6OjpkzZo18uqrr4rX6x3wmP7+fmlra5O2tjbZvn27PP3007J48WK5+uqrJSsrS/PcGRkZcuutt8qPfvSjTwJUAAAA0dTUErR4PCHd/7cvL0sM238ZAAAAGEm69/4j/cjzN1aG/H2DrEKgiDnJ6TdZU/xqoN8U7HMnarVVGHDefa8VHHn+y2rx1X/YNdg163HvfCmz9s//OWMoVRWUBGvQbEvvV8wJITXgN4X6uy2hfs9Jv/8X8nWbDz39+Uojwhiqql8Zxdd2MGkoYQyTxRYwJaX6FbMlpAZ8pqCnPUkvPHC8zh1/GSOiyMehjPjiLFvitmaP7/K1fjTgXan9bQdSeg6+kxLd1iuquKte0gyFiIhkzvvKiG5Xogb9yuGnr5jWc+jdQVeAURKsQbM12a8kWEMhX0/Cx+G2yPY1Qv29lpo/fmHW2OtXbbEXze4Z9MJPMif9BzLiz4IFC8JWIvB6vfLQQw/Jhg0bBn1+l8slDz/8sLzzzjty55136rZs+OpXvyq7du3S3awP56yzzpL58+eHPe7AgQPy4IMPSl3d0D/7a2pqpKamRp599lmprKyUyy+/POIqIyONw+EQh+PTFZFCoZC8+OKLsnLlyohbu6iqKlu2bNEct1qtcvPNN4fd5G9qapL7779famtrI5pXRMTv98vq1avl/ffflzvvvFOKinTbmMmNN94oVVVVn1R/GYqKioqwv1+BQEB+97vfydtvvz2oc+/du1e+853vyM033yxLliz51NhQqrcMRlJSktx2221h25Rs2LBBHnroIc0ghhZVVWXdunWyceNGuf3223VbD02ZMkXOOuusE6oiAwAAEKn2Nn/Yu+fy8yz9sVgLAAAAEAueI1ucR57/ysxIwhgmS1LQMW5Rs6P09I7k8jM6kvJP8RxfTSDY22px7/pbRs9Hr2d1ffTPvHChCPeuvxU2v35vT+5Z3685waciIiI9B9elDiaMkeDI6kuecFazo2yRy1k6v8uSVvzZuwNDQaWvea+t9/D6VE/9B6meug/T+tsOhC/XfJz0mcsbbQVTP7VhroYCSvMbP9O8o9aWe0pn6rTPNQx2rqOSCmYYsrnbsOa7Zf3tB5O1xq3ZE9yOklNd1pyJvQmOrP4EZ16/yZygBvu7zX3N+x19DdtSeqo3ZAV6mm0S0p5H9XtNNX+6ZkZkYQxF7IXTXSlTPteYMvFslzVn4qcuVKuBPlNvzcZk96412e6dqwuDnnb9cu8i0rnjpTHWrLLenDO/eyT8/LGVXnl1bdOrP/lsae1/cW1+vMhZttCYkuoDcO9cnenvatTcBLTlntJpL54z4PtNUUySd/b3PrO2nkPrM3oOrtPc8MiY88XqxLQxg9uA+BeTNWXQ/dRrVnyxItIwhjWrvMtZvqTVOW5hh7NssfszlUFCAcVzZIvTvftvWV37XsvT+/0R+TiUcWTFl2aM/9o/N5mTc8L3uD+JEchAXElJSZEbbrhB95jm5ma59957paFhyH/PRURk586d8v3vf1/uvvtuSU1NHfCY9PR0ueSSS2TFihVDmiMrK0uuu+66sMetX79eHnzwQfH7OShqKwAAIABJREFUjfm8UlVVtm7dKlu3bpUpU6ZIKKTzbWGUCAQC8stf/vKEK5oc7/LLL5ecHM32YiLycUWOH/7wh9LR0TGkORoaGuQHP/iB/PjHP5biYu22bykpKXLNNdfII488MqR5FEUJW6kiFAqd0M/R7/fLb3/7WwmFQhG3eDHCNddcE/Z1euqpp+Svf/3rCc3T19cn999/v3zjG9+QRYsWaR531VVXydtvvx22QgsAAMCJ8vaFIr6bDwAAABjp/O6GxJoV180M+b26+1smq9OfUfmFmpxF36oNt1FodmT7M+Z+qTlj7pea/a7DBxpf/WG5e/ffC0RVNe/Sa337NxOSy8902UvnD/3uOfm4pcaR526IKIxhSRvTm7PomwcyZl/bErZlismsJuVXeJLyKzyZIo0iIr7mPTbXhyvyunb/I6/fVa270XpU8qRlHcmTln3qwrca6DPpBTISM0p6sxdHvxWFnv72g87ewxuyPzOgmNTUScsaspd8+7CtcLrmnbjJE87pFJF6EVU6t/85y2xP19wor135nxO1qkAcKzGzrLvwovt3O8uXurWOURKSQs6yJW5n2RJ3wfn3Hmp+4/7itg2Ploer4tLy1q8m2kvmuZ3lZ2ieezhkzr2+qeXNX04M9fcO+Pvavf/1/KC3Y7/Zpv3zPRGuLU/r3gWbNusL2u9Tk1kd6H0c8nvMeoGM1CmXtDjLlsTkdWh69UelkbQpsRfNas9edGt1yikX6PdaNyWo9tL53fbS+d35599T7dr8ZF7LugfK9UIt/q5G+5E/f23y2C+vrBrCUzhpEMgYZkuWLJFJkybFdM4jR44M+s73WLnyyislJUX771ZPT48hYYyjampq5H//93/lhz/8oWYFhAsuuED++te/DqkqwfLly8O2oti4caP8+te/jriiw2Dt3LkzKueNNw899JDhYQyn0ynnn3++7jEej0fuu+++IYcxjurq6pL77rtPfvazn+n+DixdulRWrVql29JGy+zZs2Xs2LG6xzz77LMn/HNUVVUeffRRKSwsjEmFloKCAjn33HN1j1m5cuUJhzGOCgaD8tBDD0lhYaGUlZUNeExaWposW7ZMXn75ZUPmBAAA0JKWag574aq1LUAPagAAAIwKdau+MTlcJQJb7imdY65+qsqaVdY32PNbMkp9xVf/YVfH1j+21q++Y5rWRrgaCip1f/lmxYRbNmwKG47QpMqRF2+uCHo7w1RBUCRjzrXVBRfcd1BJSBry3ZfW3Mne/GX3VOcvu6faXbUqs/Wd34wLBXyjMuA90AZygjPHO+by3+5wjj9z4B7dA1IkbfqVbVqj3fteS3PvXD0m3FlSJ19QP+bzj+5RLLaIXz8lISmUd+6PDqdMPLe95tkvz9RrvaKGgkr96tsrJty66T3FbImbXtoma3Iw5ZQL6zu3PV8y0HjI32d2bXoyL3uJ8QEef2dtYm/1+s+Gco6uLdHhz5zzxSaj540Vz+ENyW3rHxqvd4xitoRyz/ivfdlLbo+8tPu/Hy0Zc69vSp36ubbDT39+uqfug0ytI3sOrsvr+HBFY/rML7QOfp6TA4GMYRZJKwujbdy4MS4DGVlZWXL22WdrjquqKr/4xS8MC2McVVVVJS+99JJcdtllA44nJSXJ4sWLZc2aNYM6b0FBwWfaNhyvpqZGfvOb30QtjHGy+Oc//ylvvfWW4eddtmxZ2EDNE088Ydh7sqWlRf7v//5Pbr/9ds1jzGazXHTRRfLEE08M+vwXX3yx7vj+/fvlL3/5y6DPO5BAICAPPvigPPDAA5KYmGjIObVcddVVYjZrB8g3b94szz77rKFz+v1+eeCBB+SBBx4Qq3Xg/18699xzZfXq1aKqcfP9FwAAjEL5eYlh25F8uM2T8qVrM5pjsR4AAAAgWlzv/yG35+BbuqX5kyctayi5+qldJ7opnT7r2pYEZ+77NSuum60VyvC1HUhp3/REXub8rzYOZY629b8r8NRuydI7RjGZ1YILf16VMfdLhn6fT512WXvqtMva+xqrtHu6jyKW9DE9ZTes3jpge5chU6Xh79+bLKL/VkubemntmOWP7zm+VU6k7KXzu8fd8JctBx+/eK5eC5N+V42z5c2fFRvVSscoWfO/Wte57YUSrZ+Ta9vzRdEIZLRvfLxQDQU0A0cpk85rVCyOkVlePhRQal+6dYoaCmq+qUyJjkDJtSved45d0KV1TCTMtvTAuBv/+uHhPyyfqtcapfmfPx+fPv3zbWJKYENkAAQyEDeWLVsmFov2jVvvvfee7NixIypzr1q1Ss4991xxOBwDji9dunTQgYxzzz1XTCbtcGkoFJIHH3yQlgYnyOPxyB//+EfDz6soipxxxhm6x+zdu1fWrVtn6LwbNmyQc845R6ZOnap5zKJFi+QPf/jDoFrc5OXlyeTJmlXkROTjlh5GtrdpbGyUv/3tb/K5z33OsHMeLzMzUzfY5vf75cknn4xKKKKpqUlee+01ueiiiwYcz83NlYqKipOmSg0AABgeDocSslqVoM+naiZU33q3J1sNyUfKqLz3DQAAYPSqfenWSSZzYkw3DFOnXtqcOuWS9ljOGQk16Fda3vwf3bvB7WNmtxsRxjgqeeI5nblL79zbtPaeU7SOaV3/u3GZ87/SONjNdtXfa2p9+8Fy/aMUKbjgvh1GhzGOlZQ/zROtc8cLc1JK/9jr/vyBsWEMkY6tf8rpbz+k2/rFXjjDVXTFI3uHGsY4ypozyVu8/LFth59ePkcvZNC+6anSnMW31w6mEke02Qpn9NoLp7s89dsyBhrvbzuQ0nPwnRRn2cITCg58miqdVSsL9Y7IOu2rw9pW50S0vfdIvt57TzElhIou/922Ew1jfHK+hKRQ8dVP7frot4uT/e76AUNc/o5ap2vL07kZp94wYquORBOBDMQFs9msu/kdCATkT3/6U9Tm93q98uqrr2pWySgpKZGcnJyI20RYLJaw1THeeOMNqa6uHuxScZy33npLuroM/Dv9LxMmTJC8PM02YCIi8uKLLxo+r4jI888/rxvIcDqdMmvWLNm4cWPE5zz99NM12/KIiOzevVv27ds3qHVGYvXq1XLRRRdJQkJ0/twsXbpUN/j0yiuvDKm9S6TWrFkj559/vmaFjtmzZxPIAAAAUTemyNJ74GC/Zt+72jq/8/d/cOXd8KUMLowAAACMIN17XymI9ZyJWeN64zGQ0bbhoQK/u0GzmkOCM8db+sXnthvdriF78Tfreg6+mal1Z7i/s87RtftvGSmnXOgazHlb3/5NUaC3TbMFhYhIxuz/qGZz88TlX3DfLmtW+aDb14TTuv6RUr1xU6I9UHTlozuNek86y5a4M0+9/lDbe49pBnmC3g5r23uP5WcvurXeiDmNkj7r2lqtQIaIiGvz74ucZQt3GzWfe+fqzIHa1hxly6/otBXN6jFqvlhSg36l9Z2HBu6l/i85i7+5P7XiokF9JoVjtqUHiq98ZPvBJy6eJ6o64GaTa+szRXxmDYz7YxAXKioqJC0tTXN827Zt0tQU3d/hcC0vpk+fHvG5KioqxOl0ao6HQiF56aWXIj4ftL3xxhtROW9lZaXuuMvlkqqqqqjMvXfv3rDv93DrO96sWbN0x6PVxqirq0u2bdsWlXOLfBw00fPKK69EbW4Rkba2Nt3AxcyZM6M6PwAAgIjIjKn2sD2Q77mvcfK27d6BSwICAAAAcc71/jPFeuM5i755wGxLD0Rj7tyzvn9Qr8KB64MVgwzOqOL68LkxekdY0op6Cy68/+Dgzovj2cfMbk+feVWr0efta6yy+1r2aG9qiUjG3C9VGx0EyTvnh4cTnDlevWPc2/+sWxliOKTPuqbFbM/UrFDSvf+1/KC3w7C7Ol1bni7SXc/Ma2qNmivWOrc9lx3oabZpjVvSxvRmL/l2VJ6fvfS0bue4RZp3wHobdmT4WvZqru1kRiADcSHc5vJ7770X9TXU19dLZ6f2dcwJEyZEfK7Zs2frju/cuVOam2nhfKJqamqiVmVk2rRpuuPvvPOOoe09jqWqathWKOHWd6ykpCQpL9eufhcIBKL6OxatsEd2draMGaP9/y3V1dVRD3KJiOzatUtzrLCwUGw2vn8AAIDouvCClLAlwbp7QpbLv1A9+w8rXLo9twEAAIB403NwXapeef7EzLLuzHk3Ru1CoL14To+j9FTNTX3P4Y1ZEgpG3JOie++r6f7OOt2wdO6Zd+5XEpLipu3ESJVx6vVHonFe19Y/5euNmyy2QO4Zdxg+t2KxhTLnXFejd4y3eU+ar3lPXF2UVswWNW3KpZotQkL+PrNr85P6Jcsj5HfXJ/ZWr8/WGjdZk/0Zc744YjfoOrau0A2b5J5xx0dGVwo6VvaiWw5rj6ri3rU6K1pzj2QEMhAX9DaXVVWV999/P+prUFVVdu/WrohUWloa8blOOUWzpZyIyKBaTUDb3r17o3LehIQEGTdunO4xGzZsiMrckZ4/Oztb0tPTIzpXeXm5ZksNEZE9e/ZIT0/0qnO9//77EgwGDT9vuFDKpk2bDJ9zIHqfG4qiDOqzAwAAYCjOPze1Iz/PErb3c2dn0Pqt/6qffsayA7P/uKIjJxA4wUbGAAAAQAx0bl+lGypOn3ZZvZjMUduAFBFJOeVCzQ3cYJ870VPznmZg5HjuXWs0N4tFRBIzSnrSZ15teFWHk43ZnulLn3Z5WzTO3XPoXd1N5+TyJc0ma7LxF8VFJHPulxoVU4JOWEcV9+41cbcpnnnaTfWKzu9px4fP6wYNItW+8bECNRTQ3P9OPeX8BsViG5Fhp2CvK8FTtzVTa9xsS+tPn/n5qH52OMuWuBOcuZpVWnoOvaO5vpOZYeVfMDTvvfdeTO7gPtaRI1EJBA5ZYmKiFBVpf862trZGdbP4+Lm0FBREVnXMYrHoPh8Rkf379w9qXRjYoUOHonLewsJCSUjQ/ngMBAJSU6MbQj1hDQ0N4vF4xG7XbHMmxcXF0tHREfZcxcW61fyi9nM8yufzSUNDg241i6EoK9NtkyaHDx82dD4tLS36N6Tm5+fLnj17YrIWAABwcjKbRb35xszqu+9tqojk+Kod3oxv/lddxt33NPafNs/Zes7ZzrZzzkp15eWY/dFeKwAAADBYvYd1Nr8VRU2feVXU73ZPnnC2q/Hv39cc7z2yOcU+dkFXJOfqqX5bN5CRWnFJg16LFEQmZeLZTWJKMDyoE+xtS/C17k/ROyZt5lWNRs97lDk5z28vmdfWW/1ujtYxvdXvZsgZ/xVXbTmsWWV99pJ5rb3V6wdct6/tQEpP9bspzgh/jwamSmfVSp0NOkUy531Vs1JHvHPvWp2p6lTjcZYtbo7Ge/54jjGzXO49fx+wNU5f877UaM8/EhHIGGbr1q2LSfWHeFZSUqJ7935dXew+G91ut+ZYUlKSOBwO6e3t1T1HcXGx7vPp6+uLu1DMSBXNQIae2tpa8fuje61aVVWprq6Wigrta+pjxoyR7du3hz1XuIBQLIILhw8fNjyQMXbsWN3xWH12dHXpfz/MzCQQCgAAou9rX81ueGGlu3DXHq9uH+NjubuCif94zV34j9fchSL1Ulhg8VScktQ5fZq9a+5su3v+XEe3zabE/Z1D26s8aV/+yhH9MoVGU0SefKxYu1QaAAAADOF31yf2u2qcWuPWrPFdloxSX7TXYc0q70tw5ngDPS0DtoLwNu6MqEKGv6vJ4u+s121Xkj5z+YhtpxBP7GPmaPeoPwG9NRtTRA1pboqbEh3+lEnnhb+T8gSknnJBs14gI143xTNP/XKtViBDRMS18Yki59gFQ/7/LPeuNRl+d6PmXa62gikdtsIZ+pt8cazn0NsZeuMpk8+LSWUdW/GpnVqBjKCn3ep31ydaUgv7Y7GWkYJABoZdbq5+C2OTySSnnXZa1OY3mUxis338HaqkpET32NTU1LCBjKws/UpQBw8ejEr7hpORy+WKynnDbaDHqvJCuEBGpBv94d6T0a72cXSOhQsXGnrOvDztlnKqqkpZWVnY0MaJSEpK+iR8FQgENKuqRNpaBgAA4ESYzaL+/pExVedcdHCeuyuYOJRz1Df47fUNfvtrr3cXfHxORS0tsfRMnmRzz6m0uc9YnNxRcUpS2NYosVbf4HfUN7h1L2obTfn48iuBDAAAEFU5S761z2Sxx/RismPswqhsYg+Vp+Y93UoESVnlsSmvLSIWR7ZPK5Dh76zTLnV8DM8R/edjSc33WHMmabYDQOQcJaeeQKUFbX0NVZoBIRGRpOzx3dFuoWMvnqv73AK9rUnB3laL2ZEdV1UQU6dc2m5J/YFHKzTRvX9tXtDbsd9sSw8M5fyuzU/p3p2aPvMLI7Y6hohIX8te3c8PW+HMmHweJqYX9umN+1r22glkfBqBDAy7cJvFM2bMkBkzZsRoNfqsVmvYY8I9n7a2qLQsOymFC8cMVUaGbshQ2tvbozLvYOeJdKM/3POJpO3JiTJ6DovFIsnJ2qFzRVHkW9/6lqFzDlVi4pD2QwAAAAatvNza96enSrd+4UuHZ3UNMZRxrGBQVQ4e6k8+eKg/ec3f3UV339skmelm37Rpto6FC5zt1y7PaM7INA/pQhkAAADCy5p/c73ZkXVSf9/yNGzXrTyRmFUWs7vdTfZ0zQ1Gf29L+M0DEfE2VOk+n6TcyVEJEZxsTImOgDU3OsEWX/tB3TC4NSf6r6Etf6rHZEkKhvx9muXavQ1VDuf4M+MqYCWiSNq0K+ta3/nNhIFGQ35vgmvzk3nZi28bdHDC765P1KsaYk5K9mfMvnbkVp8JBZT+9mrNzw+TxRawZpXpBiWMYnHm6YYt+jtqI/o8PJkQyMCwS0uLuKLusItkY1Vvk1hExOOJuxvKRqRAIBC1tiFOp27ANWavoder/30x3DqPcjj0bxYMN48RjP6ZpaWliaKMjD6KkQS5AAAAjDL/VHv3318q2/TFGw/POFTdH1HZ5MFo7wha33yrJ+/Nt3ry7v+f5skLTnO2XP359PqLL0h1JSRI1HvVAgAA4OTi76gbsCLFUX3Ne5KbXrk7emVyjxHoaU7SGgt5OyO6COjvbNA8h4iINWtCzCp+jGZme5pPJDrXj/3dLbqvYVJ+RfRfQ5NZtaSP7fa17NHcYOt318Xlhens026qb9vwcLka9JsGGu/48PmioQQy2jc+VqCGAgOeU0Qk5ZQL6xWLLe5bcmrxd9YlqsF+zeenJFhDTa/8KCafhcE+t26+INDbwl2qxyGQgWGXlKT7tyuuHG1PoCdcaCNaVR1ONtEMEVgslmGb+1jh3iuRVl7QOy6awZZjGR3IGEkhB5NJ8zsSAABAVEyeZPW+vXbCph/d2zj2mec6Svv6QuH/R2YI+vtV0z/Xdef9c1133r3FiT133ZG3/4rPpVISEAAAAIYJ9OhXnuje/3p+9/7XY7UcTaGA9kbpsQI9TbrPJ8GZ4zNmRSc3k8URtcoygd5W3dfQkpIXk9fQ4sz0+Vq0xwPuhri8iG5OzvM7y89o6t73WsFA4762Ayk91e+mOMcuGESlEVU6q1bqtCtRJGv+V0d0uxJfx2HdzdSgtzOx9d3fjo/VevSE+r1RuQYxkrFLhGEXbvN7pAn3fKiQYQxVjd7Nd/HyGoYLfkQayEhI0M7exSpcYvQ8IymQAQAAMBxsNiX0s58WHNz8zoR3ll+eVmOxKFG9E6jmSL/zpm8cqTznooOzPtjm1S/RBgAAAEQo0Ns2Ii4E6t25fqxgX7fuxWdzUupJ3aLGKGarM2o/x1B/r+7N7uaktJi8hmZrsu48wb6uuL0pP2v+V2r1xl2bntAJV3yWe9eaDL+70a41biuY5krKnzaiN+f8XU0jpuqEGugjf3AcfiAYdqPt7vFQSP86p97mOOJDvLyG4SqyhFvnUXrhlUiqvhjB6J/ZaPvcAAAAiJbCAkv/Q78Zs2/rholvf/vW3L0TJ1jd0Zxv6weezAsvOzjv4cfaBrzbCQAAABiMUNA3Mi4EqqoiEXTwCwV8uhdkzXYCGUZQLLZgtM6tBgZutXGU2Rab11CxpujOEwr0xW2VAmf5UndiZlm31nj3vrV5Qa874vW7Nj+lG+DImHXNiK6OISKi9vfG7et5PFUNjYx+8zHEzjCGXX9/v+74tm3bpLW1NUar0dfR0RH2mHDPx+HgZrF4F66Fh82m27bQMHa7ZqBTRMK/147Sez42m00URYlqxZGj8xgp3HN3uVyydetWQ+ccqkOHDg33EgAAAKSwwNL/3Ttyjnz3jpwjO3b12f/8UmfO++970rfv9KZ7vSFDrw34fKr5Bz9pnFLf6Lfee3d+tZHnHkjFZFvnpZekNkR7HgAAAMSe6u8fMZuQkQhXScNsSyeQYYjo7QerIf1AhsmaErUwyKfmSbTrzqNG2EZnuGRUXlXbtPanpww0FvJ7E1ybn8zPXvzNsEEKv7s+sbf63RytcXNSan/6rGt0mruMDKHACAmnYUAEMjDswm2sbty4UdauXRuj1Zy4vr4+3fFwm+wYfj6ffou3WL2G4eYJt86j9N6TiqJIUlJS1FuXGP0zC/fce3t75ZFHHjF0TgAAgNFiakWSZ2pF3mERORwIiLJxU0/yund703fu8qbs3deXXN/gd4RCJ3YFU1VFHn6sbbzDZg5+946cI8asfGAlJYm9t/2/7BF/xxMAAABGP8WUoFv2WO33suka7xSz7t2Nqj82G+dquA16s/57bbhlzL2hseWtX00I9XsG3Kvu2PZcUSSBjPaN/1eghgKaP4uUUy5sUBKS4vpnEQk1RNWJkYxABoad261fMTcnRzPYFpfCVdGgQkb8C/eejNVrGG6ezs7OiM7jdruloEC7YrTNZhtxgYxIPjdiUfkDAABgpEtIEHXB6c6uBac7u47+W3d3yLz5fU/ylg88Ke9v7U3btt2b3tEZHFLv7gcebJ44c4ate9k5yeHLDQIAAADHURISdTdSM+Z++ZA50RGTigThhd8vNSVYdZ9PwNNhMWw5iApTQmJQ7w0X6uuIyd5ryNejO48pzkMIZltqMGXSeQ2dVSuLBxr3tX6U4ql+N8U+dkHXQOMfU6WzaqVOuxJFsubfNCrC++ZE/TY8ttxTOp3jl8ZFuwHH2NP5///jEMjAsGtra9MdH2mBjPb2dt3x4uIB/7YgjoQL1eiFG4xUVKTb9kxcLldE5wn3fPLz8yM+11AZ/TPzer3S29urGVqxWq2SkpISNrgBAACAz0pONgXPPMPZeeYZzk4ROSIism271/Hci515/3jNnVdX7484oRwKiXL7d+qmLF44ab3NpsT1BUEAAADEH1OCVXcTMn3G8iZ78ZyeWK3nRJks+puqwb5O9u3iXLiQUMDrjslrGPR16wcyLPotTeJB1mk31WkFMkRE2jY9UVQ8dsFurXH3rjUZfneD5t2g9qKZ7Un5FZ4TXWc8UMIEMsyOLF/esh9HvWUohobSRxh2zc3NuuMjLcDQ2NioO15YWEiVjDgX7j05duzYmKyjtLRUd7ylJbK2Z+GeT7h5jFBSUmL4OZuamnTHR9pnBwAAQDybMd3We/+9+Qe3b560/vePFm+tmGyLrFybiDQ1B2w//1UzX84AAAAwaAmOLN3exb72j2yxWosREpJzdHueB70uKmTEuQS7/nsy6HXFpkJGX5fue8WSnBdZz/NhZCua1WMrmKp5R2n3vrV5Qa/brDXu2vK07l2t6ZXXjIrqGCIilpSCfr3x/q6GEfVZeLIhkIFhd+jQIQnqFHgqLi6W7OzsGK7oxLS0tEhXl3YFJUVRpLy8PIYrwmDV1NTojmdlZUlKSkpU15CYmCiFhYW6x4Rb51FHjui37I5FICMacxw6dEh3vLKy0vA5AQAAIHLJhantb79evvmeH+bvtFqViO66WvFcR7G/X6XnLQAAAAbFkpKnG2Dwtx8eUZuQlpR83U3y/vZqY3s/w3Bmp34go695rzMW6/C5DuvOk5BaEPeBDBGRjFnX1mqNhfzeBNfmJ/MHGvO7GxJ7D72jWWLfbEvrT591TWR3tY4A1swy3b7zga7GEfVZeLIhkIFh5/P5pK5OP6Q20jZWDx48qDteUVERo5VEh6KM7uuoHR0dYVtdTJo0KaprmDRpkpjNmsFPCQaDYX9vjjp8+LDu+MSJEweztEHLzc2VtLQ0w8/70Ucf6Y7PmjXL8DkBAADwb1+/Kavh5RfGbXI6Tf5wx7a1B5Je/ltXZizWBQAAgNHDklGq227A07QzunfOGcyaPaFXb9zbvC85VmvB0FjTS3Tfk33Ne6L+Gvrd9YlBT7tV75ik7Am6G/jxIn3WNS1mW5pm9YeObc8NWAWjfeNjBWoooLnPnTrl4nrFbFGNWGM8sKQV9ZsstoDWeKi/1+Jr2UcoI04RyEBc2LZtm+743LlzY7QSY3z44Ye64wsXLoz7UIPPpx2eTEgY/W3sdu3apTu+cOHCqM6/ePFi3fH9+/dLf79uhapP1NbW6lZtKSwsjGqVjAULFkTlvNu3bxdV1f4+VVhYGLbKCAAAAE7MnNn2nkcfLN5mMknYC13/eLVr5JQ+BAAAQFxwjJmjfWFTRDyHN2ZJKBjfF9uPYS85Vff59LftTw16O0b/BfgRLCmvokdv3Ne2P+qBDE/Ne7pBJFOiI2DNGRmBDCUhKZQ65eJ6rXFf60cpnup3j3u+qnRWrdJpV6JI5twbNc85MilizZ6g+/nh3vVyVqxWg8EhkIG4sGnTJt3x6dOnj6iN1a1bt+qO5+TkyOzZs2O0mqHxerX/VicmJsZwJcOjqqpKd3z27Nlit0enepzVapV58+bpHhNufcdSVVV27type0y0QhMi0QuvtLW1hW1bcv7550dlbgAAAPzbsnOSO5adk9oQ7rjtOzzGl00DAADAqOYYe3qXYjJrhn+Dfe7Enup3RkyVjKTcUzzmpBTNO+3UUFBx73g5vivLKdqvh4hISB3drQpZR24QAAAgAElEQVTtY+Z064373U32aFcq6D6wLkNv3Jo1rltk5LwMWfNvrhPFpPm+atv0+09tELp3rcnwu+s1N2jsRZXtSfkVupVMhkoJ8/6PJlvBNN3S7j0H3iKQEacIZCAu7N+/X5qamjTHFUWRq6++OoYrOjFNTU2yd+9e3WOuvPLKuK6SoRfIsNvtYjKN7o+PLVu2SDCo3Q47MTFRli5dGpW5Fy9eLElJSbrHbNy4cVDnDHf8okWLolL5ZMKECTJmzBjDz3vUW2+9pTt+5plnSm5ubtTmBwAAwMe+c3tOdbhjjtT6ncHgCLoqCAAAgGFnsiYHbflTO/SO6dy+agRdAFTEXjKvTe+Izh0r82O1mqFQzBZVMVtCWuNqwDeqNw+sORO9Ccl5OtUnVHG9/8e8qC0gFFC69r6m+x5xlMxzRW3+KLDmTPQ6Sk7V/L3o3vdaftDr/qTHu2vL0zrVMUTSZ//H/2fvvsOkKu/+j3+mbZntHRbYpS5IlaogFhSUGI09scQSY0wsj0YTE+NPY4yJxhhb9DGxRIzRRBNLjI9GRVAQEVFEWVBAlrq9t5nZ3Wm/P1aUNrPtzJzZ3ffruua62J17zvnO7uww5z6f8727t957L1jinKFPHEkKej0Re/2nTlgc9r3DXfpxlt9V44jU/tF7A/pNEf1HMBjUa6+9FnbM3LlzNWPGjChV1Hdvvvlm2PvHjBmjk046KUrV9JzbHTo8aLfblZ09sIN2jY2N+vTTT8OO+fa3v620tDRD9xsfH69zzjkn7JiSkhLt3r27R9v96KOP1NbWFvL+rKwsLVq0qEfb7I5IB6mWL18eNjzkcDh02WWXxXT4CQAAYCCYNDHBPSzfEfYKJL8/aKmr89N+GQAAAD2SUrSoOtz9TZteHrbvydpYlzr+pJpw97t2rslpqyiOTHtmg1gdCSFPSgd97f3md9FbSQVz6sLd3/z5qxEL1TR99kqm310XH25M8riFYeuLRZmzL94T6r6A12Ov//DJIZLkbSqPc+1YlRtqrC0xoz1j+nlh3zP6whafEjaQEeiIXCAjpWhhoy0hLXSHHX+HtfrdP4YNq8AcBDIQM5YvX67GxsawY66++mrl5PSPZYffe+89VVeHf8+/+OKLNXLkyOgU1ENlZeGX1yooKIhSJeZ54403wt6flJSkCy64wNB9XnjhhcrMDNttTK+//nqPt9ve3q7ly5eHHXPuuecqPd24LtJHHXWUpk6datj2DsXj8ei///1v2DEzZszQ6aefHtE6AAAAIE2elBD+gE5SZZWXq3UAAADQIxmzL6mw2OJCdmQItLc6KpfePiqaNfVFxvTv1NgSUrwhBwQDlsqlvxkTxZJ6zGJ3+kLdF/A0DfjP/GlTTq8Kd39Hw+7kpo0vRWTpmbo1j4c9OWNPymlLGbsg7NIWsSh9ypl14TqPNKx/drgk1a15LD/o94Y8v502+bQyi80RsWVFbAmpIV/7kuRz1cZFat+y2oOpE08Ou1xo/YdPjfS3VA74v8H+hkAGYobH49Gzzz4bdkxaWpp++ctfdnnC2kijR4/WpEmTevw4n8/X5fOJi4vTTTfdpCFDIte9qrd27AjfcXjWrFlRqsQ869at6/LncPzxx2v+/PmG7G/+/PlavHhx2DG1tbVdLtMRyquvvhp2GZbk5GT9+Mc/ls3W9wDz0KFDdfnll/d5O93x4osvdhnmuuCCC7Rw4cKo1CNJiYmJOuGEE6K2PwAAgFiQlWkPeaXOXtXV3shNTgEAAGBAcqQO8aYWLaoIN6Zh3dMj3bs/TI5WTX1hcSQGUieeEvaKyJaty4Y2bXw5Iif0jWCNTw55UrqjuSwxmrWYIW3iKfX2pOzQLaklVbx+24Sg32to6+bmz17NdO1cE/aq5dSJ3yyX1RaxQELEWG3BjGlnh1xqpL1ma5p7x6rUxg0vhO4AYbEEs+f+IPzVxn1kc2aGDlNJ6qjfHtHuNjnzr94jizXk7zfQ3urY8+/rxkeyBvQcgQzElGXLlmnz5s1hx+Tn5+t3v/udxo+P7PvJtGnT9Ktf/Up33323ioqKerWNlStXatOmTWHHZGVl6Te/+U2v9xHOiBEjdO2118pu73lX4J07dyoYDP1/9nHHHaf8/Py+lBfzgsGgnn766bBjLBaLrr76as2ZM6dP+5ozZ46uvvrqLpfWeOaZZ8KGKsKprKzUW2+9FXbMlClTdPXVV8tq7f1/Dzk5ObrllluUnByd4x+Px6MlS5aEHWOxWPSjH/1I3/3udw0JnISSmZmpCy+8UI8++qjOPffciO0HAAAgFiU5bV1+ULXZLf1vYhAAAACmyz3+hh0Wqz1kl4yg32vd/c8fTPO7avrFleF5x9+wy2KPD/P5Oaiyl6+f4q3fGXZpCrPEZ49pCXWf39MUH6t1G8ZqC6ZP/3bIJTYkydu4J6nqrTsKjdploL3FVv7qjRPDjbFY7YHso64IGWqIdVlzf1hmsTlC/p3v+ff1k7xNZSEDD87hs+ricyeEXuPcAInDpruk0OdxPJWbUiK5//jcCZ6U8eEDai1b3syveedeli6JIQQyTJaenq68vDxTbkYuTWCUQCCgBx54QC6XK+y4vSGG733ve0pNTTVs/0OHDtWZZ56pe++9V7/85S81ZcqUPm0vGAzqoYcekscT/v0/IyNDv/nNb3T++ecrPr5vn1MsFosOO+ww/exnP9N9992nY445psuT/Ifi8Xi0bdu2kPc7HA7dfvvtmj9/vhyOfvEZt1c++eQTrV69OuwYh8OhG264Qeecc06PT/bbbDadddZZuuGGG7r8ORYXF+vdd9/t0fYP9Oyzz6q1tTXsmGOOOUY333yz0tLSerz9iRMn6o477lBeXl5vS+yVVatW6Z133gk7xmKx6IwzztBdd92lyZMnG7bv+Ph4zZs3Tz/5yU/0pz/9Saeffrqczphe4hEAACAiaut9XXa/SEmx9i5dDAAAgEEtYegUd9rUs7o8AV7yl9Nn+FuqozJhHWhvsVW/dUeBp3Rdj69Mc6QXtmcc/u3d4cb4PY1x25ecMdPbVBZzXeYShkwMGciQgmr49J+50avGHLnHXr/HGpcUtltC7XsPj23Z/HpGn3cW8Ft2PfPdyd6m8rATz6mHLS6Pzx4btnNHLHOk5Xckjzk25HIwHXUlYcMOmbO+G/Ewii0xw+dIG+IOdb97z7osv6cpcleFShp28h1fhA90SZXL7jys7v1Hh0ayjn21lryTVv3WHWGX0xnMen7ZPAx1xRVXmLbvdevW6Y477jBt/6FUV1frnnvu0U033RS2s4PVatUpp5yihQsX6t1339Xq1au1ZcsWtbe3d3tfKSkpGjVqlCZMmKAjjzxShYWGhRW/Ul1drfvvv1833nhj2GDE3hPzixYt0htvvKGVK1eqvDzsUlBfsdvtGj16tI488kgdddRRys7ONqT2d999V+PGjQt5f3p6uq677jq5XC5t375dpaWl8ng8crtD/l8kSdqyZYs+++wzQ2qMhscee0xFRUVhf65Wq1Xnnnuu5s+fr3/9619as2aNfL7QS4nZbDYdccQROvvss7v1umtubtZDDz0UtmtJdzQ3N+uRRx7RT37yk7Djpk2bpj/+8Y96/vnntWzZsi5/p/n5+TrjjDO0YMGCQ77OV65cqWOOOaZPtXfl0Ucf1dChQ7vsnjNq1Cjddttt2rJli5YvX64NGzaourq62/ux2WwaMWKERo8erZkzZ2r69Ol9DlIBAAAMBDU13i4/FBUMj+/+ARsAAACwj/xv/raktWRFrq+lMuSSGO3Vm9O++PPCIwrOfniDc9T85kjU4W0qi6tecV9B46f/GhFob3UkFsxa25vtDP3G7dtbvngrz9tUEfIke0fD7uRtfz7xiMLznlzvLJgd/kq7Hgh6Pdba9/6U37rjvcxR33thQ08f7xw2PUwgQ6pd/ejorFkXVdpShoQNLPRntsQMX868H5VUvXPPhFBjggGfdfdzl00v+M7j61MmLG7ozX6Cfq9l198vmtS6fVXYqyCtjkTf0JNuK+nNPmJJ1pHf39Oy9a0eBwlszqz29MO/UxOJmg4Un13UEurvNuhrt5W/+ouxI85+eEuk9u/IHNmee/T/fFH19h9CvvYUDFjKX71pSlv158n537yzxGJPCNl5pPeCatrwQnbNqodHeso/zXSOmFWXu/CmsEGzwYpABmLSp59+qgcffFDXXHNNlx0HEhIStGjRIi1atEh+v1/btm3Tzp071dzcrNbWVrW0tMhutyspKemrW1ZWlkaNGqWcnLBLbRnmo48+0pIlS3TppZd2OTY1NVXnnHOOzjnnHFVXV2vr1q2qqKhQQ0OD2traZLValZycrKSkJKWnp2vMmDEqLCyMyDIMy5Yt09lnn91lF5KkpCRNmTKl2x1FXnjhhX4VyGhubtYf/vAH3XbbbV2eeB8+fLiuu+46ud1ubdy4Ubt27VJdXZ3cbrecTqcyMjJUWFioKVOmKCkpqVv79/v9uu+++1RbW2vE09Hq1as1depULVq0KOy45ORkXXLJJTrvvPNUXFysrVu3qqqqSi0tLXI4HEpJSdHw4cM1ZcoUjR49OmTgaOvWrXrllVciHshob2/XnXfeqV//+tcqKOg6iDl+/Pivwhu1tbXavHmz6urq1NLSotbWVrW3tyspKUlOp1NJSUlKTU1VQUGBCgoKBnRXGAAAgN7wdgQtGza2hW3DmOS0+nJzbQN2QhYAAMBozZtfz7TGp5jSYSwhd7w70q3/e8qWmOEbcdaDG3Y+de6cYMAf8upHb1OZc/uSM49In3b27ryFN+9wpOV39HXffk+DvfGTf+U0f/5qrmvXBzlBv7fPHeit8Sn+Yd+6d+OuZ747O9zz8bVUJm7/y6lHZs6+eEfeCTfttCWm9fo10V5bktCw7qkh9ev+Ueh318XH5xQ19WY7KeNOaLQ6En0Br+eQ5xn9nob4rQ8dNzfn2Gu3ZUw/v6ovNcey3ON/tqdh48v5HbXbQp5ECXg99l1/v3hm9vwrt+WdcNMui83R7asu26s+T9zzzx9O8VRt6rLlfc78q7c5Mkf2+wB8StGJjXGZo1o66nf0aOmP9Mmnl/bkZ9sXyWOOrW0tWREyINP4yXOFvpbKhJyj/2dn8tgFvfob60ruCTfubt2xKsu1c02YE51B1X/41KjWkhU5ucf//IsMgwIr7p3vpzR8+nxey9alQ8ItIYOvEchAzFq1apXa2tp0/fXXd/vqc5vNtt9J1ljy6quvKhgM6tJLL+32EiK5ubnKzTWvs1dbW5uWLFmia6+91rQaYsUXX3yhu+66S7/4xS+6dTLe6XRqzpw5mjNnTp/26/f7de+992rDhh6HlMN67LHHlJOTo8MPP7zLsfHx8Zo1a5ZmzZrV4/3U19frnnvuUXJyj7v29UpLS4tuueUW3XjjjTrssMO6/bjs7GzNnz8/gpUBAAAMbK/8tzmzudkftpXy+KKEiExEAQAADFSlL13b9eRdhGTPu+KLoSf/ZodZ+w8leezxTUMW37ap4rVbJkuhz70GA35Lw/rnChs3vDgieeyCytTxi2pTihY2ONILunXC2u+qt7t3f5DiKl2X6t71QaZ7z4dZRoQwDpQy/sTG3AU/3Vy17K6wk5lBv9dat+bxMQ3rny1Im3RqWfq071Qmj5rXIqst7AnooN9r8ZR9nNxasiK9+bPXhngqN6YrGOz5GucHsDgSAyljj69q+vzVYaHG+Fw1CRWv3Ty58vVbJznSh7scacPd1jin32pPDBvOyD7qyt1GdgOJKKs9WPidxzaUPHrykaHCKVJnp4yalX8sair+d37mnO/vypx1YWW4kIqndF1y7epHhzdtenl4d153SSPn1eSe8PMB05kgY/p3SquW/a77E/wWazBr3uVlESxpP5nTz6uuWn7XhKCvPeSV0q0lK/JaS1bkWeOTvfGZo1ptzswOW1yyL9zfrM2Z7h122n1fdK8KiwrP+1vxtkdPmt1Rtz1seKWjfldy6fNXTq9admdr+qTTKpLHHlefNGp+c7cCLAG/xVOx0ekp/SjFXfpxWuv2FTnhuvrg0AhkIKZ99NFHuvnmm3XdddcpPz/f7HL67LXXXlN9fb2uuuoqOZ394/1q5cqVGjlypE477TSzSzHdp59+qttvv10//elPu+waYgS32637779f69atM3zbfr9fv//973XDDTdo+vTphm9fklpbW3XnnXeqtrY2aoGMvfu97bbbdPHFF2vx4sXdDkABAAD0Vy+93JR1+qlpdRbDp4e7x++X5f6HasZ0NW7WTGevWvQCAAAA+8qed0W5z1XnqFlx//hwoQypM8jQsuXN/JYtb+ZLkiN9mCsuvdBtTUj12hPTvbbEjA5/R4s90Nbs8Le12APtrXZvU6nT21zp7GrbRsld8LM9HQ27Exs+/sfIrsYG2lsdDR//Y2TDx/8YaY1z+hJyipodqcM81qSMDltcii/gddkDbc12n6fJ4W2uSOyo254S9HdE5Egh57if7Gze8np+uO4eUmc4pqN+V3JH/a5uTRKnTfpmtfpLIENSwtCp7vxTfr+x9N/XTlMwEPZn0dGwO7nyjVsnVS29/bCEnPHNcdljWm3OdK/VkeTztzU5vE1liW3VW1LDLctzIEf6CFfheUuKpYEzD5515GUVNSsfGBcu5LKvpMI5tfHZY9siXddetpRcb/qUM/c0rO/e36ynojijO9u1J+d5hp2mbgYyJFtSpm/0957/uOSxU2d3p1OFt2FPcs2qh8bVrHpIVkeCPz5nXLPNmd1uT0zzWhMzvFarLehra7H725ocgfYWu89VH+dt2JES8LYZ36J/kCGQgZi3fft2/fSnP9UFF1ygxYsXR2RpjnDq6+u1e7dxwcI1a9Zo165duuKKKzRp0iTDthtJTz31lMrKynTRRRdF9cR6LNq0aZN+/vOf65prrulRB4aeKikp0QMPPKCyssiFOtvb2/W73/1O3//+93XiiScauu3q6mr99re/VWlpqaHb7S6v16vHH39c69ev1/e//33l5YVdXs9wfr9fGzdujOo+AQDA4HXtDaVT77i7su2aK3NLvntuRnW0gxl33VM9YtNnni5b6J59Rnp1NOoBAADAwDdk0c277M5Mb+Ubv5rUVSBgX97GsiRvY1n31pKOouFnPrTVYnME6z98alR3HxPocNvdZZ9kquyTSJYWUuKwaa7MWZdsr1v7ly7D2QNdxszzqwMdrRvLX/t/k7sKZUidHTM8VZvSu7MUSTiO9GGuMZe9/JEtKdvXl+3EGltihi91wkkVjcX/HtGd8RkzL4z6iYihJ9++veWL5Xm+1qpuh2ciwZFe2D7mh29+sOPJs2a0V29O6+7jAt42m6e8e0ER9J1J188APdPe3q4nnnhCP/7xj/Xee+/J74/scmPt7e1auXKlfv3rX+vyyy83vENBRUWFbr31Vt1///0qLy83dNt7BYNBFRcX65577pHP1/f/i5ctW6arrrpKzzzzjGkn2WNFdXW1brnlFv35z39WfX29odtubm7WkiVLdOONN0Y0jLGXz+fTI488orvvvlt1dXWGbPPtt9/WDTfcEBOvk3Xr1umaa67Rk08+adjzC6ekpERPPPGELrvsMj3wwAMR3x8AAMBe23d0pPz4htLDZx+9Ze5fn67P8/ujc3nUgw/XDLvvwaou14wcXxTfNHN6Yr+5yg0AAACxL/uoK8sLv/v0h/bkPI/ZtRhh2Gn3fTH0G7dvtNjjI3sCxED5p9y5PWX8iZE5ydHPZM29vGL4mQ9+YnUkRiUckTh0SsOYH/z3Q0d6YbeW4elvsub+aE93xtmTctoypp1dG+l6DmRLzPCNvPDpj23OLNN//o7UId5xV7z5YfrUs3YPpE4pAwkdMiIoGAyqtTV255va2nrXvcftdod8Xh5PZD/3lJeX695771VWVpYWLlyoI444QoWFhX3ebjAY1K5du1RcXKwNGzZo06ZNam+P7HtoMBjUu+++q1WrVmn27NlasGCBpk+fLofD0ettBgIBbdu2TZ988olWrFihyspKAyvuXArixRdf1IsvvqicnByNGjVKBQUFKiwsVHJyspKSkpSQkKCEhATFx8eH3E5HR0e39ufxeEK+1lwuV6+eg1GCwaCWLl2qd955R8cdd5wWLFigoqKiXi+PsX37dr3zzjt66623Iv7aO5Q1a9Zo/fr1Ovnkk7V48WJlZ2f36PHBYFAffPCB/vOf/2jLli0H3R8IBMK+H3q93h7X3F0+n0+vvPKKXnvtNR1xxBE6+uijNW3atLCv0e5qbm7Wxo0bv3rvMPpvDgAAoKd27OxIuf7nZdPuvq/afcZpaWU/uCS7vKDAYfgHzJaWgO3mX5WPfua5hlHBbnRyvu5/ckuMrgEAAABIKTqxseh/Vr1f/uqNY5s2/ntET7pl9JXFag84C+fUJmQXGXZiJPuoK8udhXOay1++fqKnom/dE6LCaguOvPDvG6vfvrup5t2HxgU6XIP6vGPG9HNrEnInfLDnhasm96RbQU9YrLZg+ozzdg479Q/bLDZHdNbVMYGzYHZr4tBJjV39HaRPOa1MVrspP4fEYTNc465cuqb0xR9PaN2+Mrptug9gcSQFRnz70c2pExbXVLx52wRvY2lUOwHZU4Z4UiecVBXNffYng/qNMdLa2tp08cUXm12G4W6++WazS1BdXZ2ee+45Pffcc8rNzdWECRM0duxYFRYWKjs7W1lZWQcFG4LBoFwulzwej6qrq1VWVqaysjKVlpZqx44dampqMuW5BINBrV27VmvXrlViYqImTpyo8ePHa8yYMcrNzVVOTk7I51JXV/fV89i5c6c2btwYtRBQTU2NampqtHbt2ojt4+67747Yto3i9Xq1dOlSLV26VDk5OZo8ebImTpyo4cOHKy8vT2lpB3/mamlpUVVVlcrKyvT555+ruLg4Jk7kt7e366WXXtLLL7+sSZMmafr06SoqKtLw4cOVkpJy0NjS0lKVlpaquLhYH3/8cdi/od27d5v+fuj3+7V69WqtXr1acXFxmjBhgsaNG6fRo0crLy9P2dnZBz1PqfO5ejweNTU17fe+sWfPHu3evVvB7pyBAAAAiLKKSq/z4Udqx/35sdqxU6c4G765OLXyzG+l14wc2bdwRkOD3/7Xp+uHPPxIzZi6Bn+3Eq6zZjjrzjkzPWJXLDU2+uJWvdeaGqntd+WoucnN0V4mBgAAAF+zJWX6Rnz70c3ZR121p2r570a1bnt7SNDvjdgntPicoqb0yadVZMy+pMKROsTwK82cw2e1jr1qxdq69x8bWrP64VHehj2GryNuc2a1p006tSx77g8MaNNsUe6Cn+3JmnNZRe37Dw9r+uy1Ie3VW9OkwTlvmjjscFfR1SvX1qx8YFjt+4+O9rlqEozatnPE7NqhJ936hXPk3BajthnLMmZ8d4/n1V+EDGRYrLZg5pGXR77VeBiO9ML2UZe+9Kl75+qU2jWPD28tWZHn9zTGmVVP2tQz61Innbq69t0Hh9V99GRhJJdossYl+VLGHVeVMeO7ZSnjFzXSnSM0Ahno96qrq1VdXa2VK1fu932Hw/HVFfBer9eUrgM95fF4tG7duoOWSNn3uXR0dHS7wwSiq6amRm+//bbefvvt/b6fnPz152WXyxXzJ/ADgYCKi4tVXFz81ffsdrsSEjo/N3o8nogvGxRpHR0d2rBhgzZs2HDQfUlJSV91OukPvy8AAIBwAgFZPvnUnfnJp+7M395VqWH5DtfMw50NU6YkNk+flthy2Ph495AhjpCTyA0NfvuHH7uTi4vbUla+15L1wYfubK832O3J7ZRkq/eRBws2GvNsDm31GlfOad/ekRPJfYRTVjL5rYQES8Cs/QMAgP7JmpDui88e22x2HaE40ob0fELdaguGe06O1F5sswcSh01zjbzwHxu9jXu21n3wRH7rtmW5nsrP0xQM9Oksoc2Z1e4smFWXPPro+rSJJ9dFZ4kIi7LmXl6RNfcHFU0bXshu+PT5oa7tq3ICXk+vz+vZU4Z4kkceWZtcdGJt+tQza43urmBLyvTlLbx5V97Cm3f5PU02T9m6ZE/pJylt1VuSg16Xzd/htgfaW+3+DrdNgdCBGWtCereW/EjIGtsS8LcfcjtxGSPcvX0ehrDagjnHXV+aPf/K8rq1S4Y0bXhpqLtsfWZvXou2hBRv8tjjKzOPuLQsedT8iL5nxGeParUlpBzy+NCMn6nFkRD2OMtZcERtfPaY3i0HYDDnyHktBSPnfS7p8/babQnuPetS2iqKk71NZYkBr9sW6HDbg163zd/eGvJv2J6cY8h7i8XmCOYcd31pznHXlTZt/HdW06cvDWndsSrH39bUp6CIxWoLJuQd1pg0cm598rjj65LHLGgayF1ajPTVH/5ZZ531UDAYvMrMYgAAALrjwT9sfzkh3p9odh0AAMSCgqJNC1yuQK/XPoyLswQy0m3tCQlWv8NhCXR0BK1tbQGbyx2w93W7f3288KMTT0hp7O02JOnEU0pmrlvvzurLNiKpN4GMDz9yJy8+rWReuDGvvzxm9exZzthdBxUAgCh66INrdzz/2RmjzK4D/Y/fVeNoLXk31b3nw7T2+h1JvuaKBJ+7Lj7Q7rIHA36rggFZ7Al+a5zTb3Uk+myJ6d649OHu+Oyxrvi8Ca6EIVNc8dljY+KEb9DXZnXteC+1defqtLaqz1N8zeWJPldtfKC91RH0+yyyWL98Lgk+a1yyz5Ga70nIHuOKzx3vco44ojlh6CRzQwqDnL+l0tG8bXm6a9cH6R11O5K8zRWJgbZmR8DrsSnot8piC1jjU3yO5Kx2e+owd8KQiS3Jo+Y3Jo2a3zxYT3qXPLJ4hnvPhyHXVx9+9sPrMw7/Tk00a+q3An6Lu/TjJNeuNWltFcUp3uayRG9zVUKgvdUR8Lr3ew1aHYk+m8PptyVnt+gqRyIAACAASURBVMdljXYl5E5oTcyb4E7In+ayxqf076t1o8UavLb4puzNe7+kQwYAAAAAAINYR0fQWlXtMzTo6HRafQ/fN/yTvoYxAAAAgL6wJeV406aeWZc29cw6s2vpK4s9IZA87oTG5HEn8Bm7H7KlDPFmTD+/JmP6+QQIuqG9tiTBXbouZDDfnpTTljH1rIgtjTngWG1BZ8HsVmfBbEL/JiCQAQAAAAAADFNYENf6yIMjNtDdAQAAAADQG3VrHh0WbomX9KlnlMpqH5SdQ9D/dHvdVwAAAAAAEHueeaJw3QXnZuwckmf3mFlHQoLVf9klWSXvLS9aQxgDAAAAANAbQb/X0rjhpeGh7rdYbcHsuT8si2ZNQF/QIQMAAAAAgH7s6PnJzUfPT26WtHX1GnfKv15sGLLi3dac3Xs6koNRuF4oNdXWcerJqeU/vTZvV0GBoz3yewQAAAAADFQN6/+R63fXxYe6P2nkUdWOzJEce6LfIJABAAAAAMAAMe9IZ8u8I50tkr6orPQ6li5vzXjv/dbMdevdGTt2diQHgwrZ8rUnMtJt7UfOSa49eXFKzRnfSq9NTLQEjNjuoYwYEeduavY7IrX9vrJZ1ePYi9NpDYwdE9/c1ZjeVwUAAAAA/VP9B48Xhrs/e97lu6NVC2AEAhkAAAAAAAxAQ4Y4vBeen1F94fkZ1ZLU0hKwFW9qc36+uS1p6xftSSU72pKqa/wJbrff5nYH7K2ugN3jCdhtNkswLs7idzgsgbRUmzctzdYxJM/eNrIw3j1+XLxr7pFJzePHxUdteZS//GnE59HaV7RMmpjg/mBl0Rqz6wAAAACAWNK06ZVMT8Wm9FD3x2UWtqZMWNwQzZqAviKQAQAAAADAIJCSYvXv00EDAAAAAIDYEfBZqpb+pijckMyZF+2WMY0fgaixml0AAAAAAAAAAAAAAGDwqnzz9pHttdtSQ91vT8ppy573w/Jo1gQYgUAGAAAAAAAAAAAAAMAUNSvuH17z3v+ODTcme96PtlsciYFo1QQYhSVLAAAAAAAAAAAAAABREfR6rB31O+NbSt7JaFz/3DBPRXFGuPGOtKHu7PlXlUWrPsBIBDIAAAAAAAAAAAAAAIbb8cRp01q3r8rr/RYsGnrSrzZbbI6gcVUB0cOSJQAAAAAAAAAAAACAmJM26dTStKln15pdB9BbBDIAAAAAAAAAAAAAADElLnts87DT791qdh1AXxDIAAAAAAAAAAAAAADEjLjMwtbRl768zpaY4TO7FqAv7GYXAAAAAAAAAAAAAACALJZg+pSz9gw77Q9fWONT/GaXA/QVgQwAAAAAAAAAAAAAgGkcacPcyWOOrc456so98XmHecyuBzAKgQwAAAAAAAAAAAAAgOEy51xSmjzm2LoDv2+xJwRsznSvPTHT5xx2eKstZYjXjPqASCOQAQAAAAAAAAAAAAAwXNrkMw4KYwCDidXsAgAAAAAAAAAAAAAAAAYaAhkAAAAAAAAAAAAAAAAGI5ABAAAAAAAAAAAAAABgMAIZAAAAAAAAAAAAAAAABiOQAQAAAAAAAAAAAAAAYDACGQAAAAAAAAAAAAAAAAYjkAEAAAAAAAAAAAAAAGAwAhkAAAAAAAAAAAAAAAAGI5ABAAAAAAAAAAAAAABgMAIZAAAAAAAAAAAAAAAABiOQAQAAAAAAAAAAAAAAYDACGQAAAAAAAAAAAAAAAAYjkAEAAAAAAAAAAAAAAGAwAhkAAAAAAAAAAAAAAAAGI5ABAAAAAAAAAAAAAABgMAIZAAAAAAAAAAAAAAAABiOQAQAAAAAAAAAAAAAAYDACGQAAAAAAAAAAAAAAAAYjkAEAAAAAAAAAAAAAAGAwAhkAAAAAAAAAAAAAAAAGI5ABAAAAAAAAAAAAAABgMAIZAAAAAAAAAAAAAAAABiOQAQAAAAAAAAAAAAAAYDACGQAAAAAAAAAAAAAAAAYjkAEAAAAAAAAAAAAAAGAwAhkAAAAAAAAAAAAAAAAGI5ABAAAAAAAAAAAAAABgMAIZAAAAAAAAAAAAAAAABiOQAQAAAAAAAAAAAAAAYDACGQAAAAAAAAAAAAAAAAYjkAEAAAAAAAAAAAAAAGAwAhkAAAAAAAAAAAAAAAAGI5ABAAAAAAAAAAAAAABgMAIZAAAAAAAAAAAAAAAABiOQAQAAAAAAAAAAAAAAYDACGQAAAAAAAAAAAAAAAAYjkAEAAAAAAAAAAAAAAGAwAhkAAAAAAAAAAAAAAAAGI5ABAAAAAAAAAAAAAABgMAIZAAAAAAAAAAAAAAAABiOQAQAAAAAAAAAAAAAAYDACGQAAAAAAAAAAAAAAAAYjkAEAAAAAAAAAAAAAAGAwAhkAAAAAAAAAAAAAAAAGI5ABAAAAAAAAAAAAAABgMAIZAAAAAAAAAAAAAAAABiOQAQAAAAAAAAAAAAAAYDACGQAAAAAAAAAAAAAAAAYjkAEAAAAAAAAAAAAAAGAwAhkAAAAAAAAAAAAAAAAGI5ABAAAAAAAAAAAAAABgMLvZBQAAAAAAAAAAAAAAgPDq1z4xpGr5H8b5XDUJ8TnjmtMnfasiY/YlFY7UIV6za8OhEcgAAAAAAAAAAAAAACCGtW57O63slZ9NUTBokaT26i1pVdV3p1WvuK/IWTinNuPwb5enTzun1mJPCJhdK75GIAMAAAAAAAAAAAAAgBjW8Mm/huwNY+wrGPBZXTtW57p2rM4tf+1mb2rRosqMmReUJ49d0GRGndgfgQwAAAAAAAAAAAAAAGKYxWoPdjUm0N7qaCx+aURj8UsjHOnDXGkTT63InPO98vjssW3RqBEHs5pdAAAAAAAAAAAAAAAACC1z1oXlFquty1DGXt7GsqTa1X8eu/WBI4/e9vDxs2rfeyjf72mgYUOUEcgAAAAAAAAAAAAAACCGOQtmt2bMvGBnjx8YDFo85Z9mVvz31smf3zXp2J1/O3dK08aXshTwHbT8CYxHAgYAAAAAAAAAAAAAgBiXf/IdJa3b3snpaNid3JvHB33ttpYtS4e2bFk61ObMak8Zu6Aqa87FZc6R81qMrhWd6JABAAAAAAAAAAAAAECMszgSA8NOu3eTLNZuL10Sit9dF9+44fmCksdPnbv5nsPnVb5+6yhv4544I+rE1whkAAAAAAAAAAAAAADQDySPXdCUMf07u4zcprdhT3LNqofGbbl35rHbHz9lev3aJ/OCXg9ZAgPwQwQAAAAAAAAAAAAAoJ8Y+o3fltgSM9qN3m4w4Le4dr6fU/afn0z77M6i43b+7bzJTRtfzpL63JBj0CKQAQAAAAAAAAAAAABAP2FLTPPnHHVlSST3Eehw21u2vJm/+9lLZ26+e+r8ild/Mbq9enNiJPc5EBHIAAAAAAAAAAAAAACgH8k55pqyuIyC1mjsy9tU7qx9/9GxW/941NFb7z/iyOq37ijwt1Q6orHv/o5ABgAAAAAAAAAAAAAA/YnVHsw9/sYvor3b9tptqVXv3DPh8z8cfuyOJWdOa1j/bE7Q77VEu47+gkAGAAAAAAAAAAAAAAD9TMb079TE505oMmPfQb/X2lqyIq/0haumf/67w44pfeHqIvfuD5PNqCWWEcgAAAAAAAAAAAAAAKAfyp73w51m1+D3NMQ3rP/HyJJHF8/bfM/h8ypfv3WUt3F3vNl1xQICGQAAAAAAAAAAAAAA9EOZM86vdqQNdZtdx17ehj3JNaseGrfl3lnHbH/sm9PrP/jLkKDXM2hzCYP2iQMAAAAAAAAAAAAA0K9Z7cHMWRftMruMAwUDfotr15qcsld+NvWzO4uO2/m38yY3bXw5SwqaXVpU2c0uAAAAAAAAAAAAAAAA9E72vCvKa1Y9PDbQ3uIwu5ZDCXS47S1b3sxv2fJmviN1qDtt4jcrMmdfUh6fd5jH7NoijQ4ZAAAAAAAAAAAAAAD0U9b4FH9q0cJKs+voDm9zhbN2zeNjtj44/+it9x9xZPVbdxR4mytjMkhiBAIZAAAAAAAAAAAAAAD0Yxkzzy83u4aeaq/dllr1zj0Tttxz+LHbH//m9Pq1T+YF/V6L2XUZiSVLAAAAAAAAAAAAAADox5LHHt/kSMt3e5vKnWbX0lNBv9fq2rkmx7VzTU7l0ts7UiecVJE588Jy58i5LWbX1ld0yAAAAAAAAAAAAAAAoJ9Lm/jNCrNr6Cu/pzGuYf1zhSWPnzJ38z2Hzyv/v5+Paa/ekmh2Xb1FIAMAAAAAAAAAAAAAgH4u/fDzKs2uwUjehj3JdWseH7P1wfnzSx5ZPKNuzeNDgl5Xv8o4sGQJAAAAAAAAAAAAAAD9XOKwaS57cp7H11rVbztKHFIwYHHv+TDbvefD7Mo3fuVLGn10dcb0cyvSJn+rTrKYXV1YBDIAAAAAAAAAAAAAABgAkgqPqGva9J/hZtcRKQGvx96y5c38li1v5tuTsttSDltckX3E98sShk51m13bofSrdh4AAAAAAAAAAAAAAODQksctqDO7hmjxuWoTGj56etQX/7tg/rYHj55T8849I/yuGofZde2LQAYAAAAAAAAAAAAAAANA2mGn1FmstqDZdUSbp+qz9Mq37jjs899POXbHkjOnNax/Nifo95q+nglLlgAAAAAAAAAAAAAAMADYkjJ9jozC1o667Slm12KGoN9rbS1ZkddasiKv/JWfeVPGnVCVMeO88pTxJzaaUQ+BDAAAAAAAAAAAAAAABojE3AnNgzWQsa9Ah8vRtOk/w5s2/We4I32EK3XCSZVZc75XHp87wROtGliyBAAAAAAAAAAAAACAASIhf0qz2TXEGm/jnqS6NY+P2frg/PnbHjxmTs2Ke4f7PU22SO+XDhkAAAAAAAAAAAAAAAwQzuEzW8yuIWYFgxZP1aZ0z9JN6dXv3FuUNObY6swZ55WnTvhGg6y2oNG7I5ABAAAAAAAAAAAAAMAA4SyY0yKLJahg0GJ2LbEs4PXYWza/nt+y+fV8e1JOW+rEb5ZnzrqoInHYNJdR+yCQAQAAAAAAAAAAAADAAGGNT/HbEjM7/O66eLNr6S98rpqE+g+fHF3/4ZOj47PHNqdPPq08Y/YllY60/I6+bNciyfC2G0Css9lsslqtZpcBAOilOIffq0GY601Osvosg/B5AwD6L4fDEoiPt/rNrgMAgJ5Iclp9NpuFefMYUtaQ2dbiz00wuw4AiFkWW9Aa5+TYC/txl36cGWhvcZhdR39msdoDzoJZdelTzynPmHFujcWeEOjyQdbgtcU3ZW/+ahsikAEAAAAAAAAAAAAAAHBI1rgkb8q4E6oyZpxXnjJ+UaNCXTVKIAMAAAAAAAAAAAAAAKDnHOnDXakTFldmzfleeXzuBM9+dxLIAAAAAAAAAAAAAAAA6AOLJejMn9aQNu3s8ozp51fZEtP8BDIAAAAAAAAAAAAAAAAMYrE5Akkj59XIottat614TJJPIpABAAAAAAAAAAAAAABglFJJ90m6j0AGAAAAAAAAAAAAAACAsW4lkAEAAAAAAAAAAAAAAGCsKqvZFQAAAAAAAAAAAAAAAAwwiQQyAAAAAAAAAAAAAAAAjPVPliwBAAAAAAAAAAAAAAAwhl/S3yVdaTe7EgAAAAAAAAAAAAAAgH7uM0n/krRE0i5JIpABAAAAAAAAAAAAAADQc/WSnpf0N0mrDryTQAYAAAAAAAAAAAAAAED3eCS9rM4QxpuSfKEGEsgAAAAAAAAAAAAAAAAIb506QxhPS6rrzgPsXz5osMgwuwAAAHooXpLT7CJMlibJ2o1xTZICEa4FANA93X3vBgAAAAAAAGLZbkn/kPS4pG09fbDF8HIAAACMtVNSYTfGFUn6IrKlAAAAABigHJKSzS4CgCTpCkm/7ca4xyX9LMK1ALEuRXTDBxDaG5LGml1EP9Uk6T+SnpK0TFKwtxviTRoAAMS67n7Q4UpsAAAAAL3lldRgdhEAJEnubo7rEH+3AH8DAMLJM7uAfsYr6XV1LknyiqQ2IzZKIAMAAMS67i5DQucvAAAAAAD6v+5emME8AAAAoY1UZxcddO1jdXbC+IekaqM3TiADAADEOiZiAAAAAAAYPJgHAACg72aaXUCMK5f0vKQnJa2P5I4IZAAAgFjHRAwAAAAAAIMH8wAAAPTddLMLiEFt6lyK5G+S/ivJF42dEsgAAACxjokYAAAAAAAGD+YBAADouxlmFxAjApLe19dLkrREuwACGQAAINYxEQMAAAAAwODBPAAAAH1jEUuWbFZnJ4ynJe02sxACGQAAINYxEQMAAAAAwODBPAAAAH0zTVKu2UWYoE7Ss+rshrHW5Fq+QiADAADEOiZiAAAAAAAYPLo7D2CNaBUAAPRfJ5ldQBR1SHpTnSGMl7/8OqYQyAAAALGOQAYAAAAAAINHoJvjmAcAAODQBkMg4zN1hjCWSKo2uZawCGQAAIBYRyADAAAAAIDBg3kAAAB6L0nSPLOLiJAySS9IekLSpybX0m0EMgAAQKxjIgYAAAAAgMGDeQAAAHrvZEnxZhdhoBZ1hjCekrRC3e+kFTMIZAAAgFjHRAwAAAAAAIMH8wAAAPTehWYXYAC/pGXqDGG8JMltbjl9QyADAADEOiZiAAAAAAAYPJgHAACgd3IkLTa7iD74XNI/JT0paaeplRiIQAYAAIh13W1BZo1oFQAAAAAAIBoIZAAA0DvnS3KYXUQP1Ut6XtLfJK0yuZaIIJABAABiHRMxAAAAAAAMHswDAADQOxebXUA3tUl6RZ1LkrwuyWduOZFFIAMAAMQ6JmIAAAAAABg8mAcAAKDnjpM03ewiurBOnZ0wnpFUa3ItUUMgAwAAxDomYgAAAAAAGDxYuhQAgJ77qdkFhLBdnSGMv0kqMbkWUxDIAAAAsY5ABgAAAAAAgwfzAAAA9MxkSSebXcQ+miT9R51LkixT9/9vH5AIZAAAgFjHRAwAAAAAAIMH8wAAAPTMDTL//0W/pLfV2QnjBUkuc8uJHQQyAABArGMiBgAAAACAwYN5AAAAum+8pPNM3P8n6uyE8XdJVSbWEbMIZAAAgFjHRAwAAAAAAIMH8wAAAHTf3ZIcUd5nhaR/SfqrpI+jvO9+h0AGAACIdUzEAAAAAAAweDAPAABA9yyQdGqU9uWR9JI6lyRZqs4lStANfGABAACxbraktG6MWyepIcK1AAAAAACAyMqTNOUQ3/dKat3n60ZJJVGpCACA2GOVtFbSzAjvZ52kRyU9K6k5wvsCAAAAAAAAAAAAAAAw1Y/U2VUqErctkm6WVBi1ZwMAAAAAAAAAAAAAAGCykersVmFkCKNB0lOSFopVNgAAAAAAAAAAAAAAwCBjkfSGjAlh+CQtlXSRJGc0nwQAAAAAAAAAAAAAAEAsuVJ9D2KslfQ/krKjXDsAAAAAAAAAAAAAAEDMGSmpRb0LYeyWdKekw6JdNAAAAAAAAAAAAAAAQCy7XT0LYbRI+qukEyRZTagXkuxmFwAAAAAAAAAAAAAAAMLqzhIjAUnLJT0l6UVJrohWBAAAAAAAAAAAAAAA0M+dotDdMDZJulHScNOqAwAAAAAAAAAAAAAA6KdukdSmzhBGtaQ/SpplakUAAAAAAAAAAAAAAAADQIqkkZJsJtcBAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAEBEWcwuAAAAxJQ4SZ/v8/UkSW0m1QIAAAAAAMx1q6T5X/77Hkmvm1gL+qfhkpZ8+e9GSeeYWAsAAAAAAKaKlxTc55ZobjkAAAAAAMBEL+jrOYJLTa4F/VORvn4NVZlcCwAAUWc1uwAAAAAAAAAAAAAAAICBhkAGAAAAAAAAAAAAAACAwQhkAAAweMR9eYv1bQIAAAAAgMiLxPkBzjkMLryGAADoAv+xAQAweFwvaaOkkwza3rGSPv5yuwAAAAAAoP/Ik/SZpIskWQzYnlPSryS9adD20D/8WtIrkgoN2t4MSaslfcOg7QEAAAAAEBVDJDVJCn55e0VSwSHGxe8zJigp8RBjhkp6SlLgyzEtkoYZXzIAAAAAAIiQp/T1sf87kiaFGPfCPuMuDTHmVEk79hl3jpGFImaNluRR5+/cpc5ATvwhxhXp69dGVYhtZUh6QJLvy3FbQ2wLAAAAAICYdIakVu0ftmiW9BNJjn3GhQtkONTZDaP5gDGtks6MbPkAAAAAAMAg2ZJqtf+xfYekuyQlHTA2XCBjrKT/HrCdoKT/i1ThiCm/0MG/+02SFhwwLlwgwyrpBzr06/H4SBUOAAAAAEAkZEv6naR27X+Qu0XSoi/HhApkHC1pwwH3eSU9Iik/OuUDAAAAAACDpKizo4Fb+x/rl6lzGZO9DhXISPzysZ4DHrtH0uWSbJEuHjFjrqSVOjiY8YqkEV+OCRXIOFydy5Ps+7iApH9KGheF2gEAAAAAiIgidR7c7l1yZN+D5TEHfG+4OkMXB45dKmlKtAsHAAAAAACGGqbO436v9j/uXy5pog4OZJwqafsBY+sk/VyHXvYUg8NCScU6uKPqryRN1v6BjHTtvzzJ3tsqSUdFuW4AAAAAACLmSEkrdPAyJuG+XvHl4wAAAAAAwMAxSdJ/tP8cQJukyn2+3nXA/W51duLMMKFexB6HpCskVWj/18m+rxuPOkMZ+97/iaTFJtQLAAAAAEBULJS0Xge3l9z39pmkc8wqEAAAAAAARMWhLt448OZXZ+fNkeaUiBjnVGfHlEaFfx3tFkvcAAAAAAAGCZukG3XoA+RfiIPjWJSkzlaxf5D0X0nbJNWrs81s4Mt/75L0lqQHJX1HXLUEAAAAAOieMyTV6OA5gvfU2U0D6EqOpCd18GvIJ+kGSfGmVQYAAAAAQBQlqfPKhSYdOpDhUmcL0jSzCsR+5kj6uzrbfIa70uRQN6+k/5O0SJIl2oUDAAAAAPqFSZJe0aGPKzskPSIpz7TqBofR6gy/hLuNMq26rjnU2f3iwOVJ9t42q7MTK3MTAAAAAIABa+/B8YFre4a61aozuMEVDOYYo84wRU9DGKFuqyVNj+ozAAAAAADEshHqDFv41PUxZYs6L95IMaXSgW+Kuv4dTDatutAs6gxabFP35ibWSDrGlEoBAAAAAIiQvQfHW7X/QXClpMsO+N4PdHBgY6u4iiHaLlfvOmJ0dfNJukX8LgEAAABgMMuSdK+kNu1/zPh/6gzz7/36SR08l1CmzrkEljo1Vn8MZCyUtE7719igzqVW935dL+l+dXbw3Hfci5ImRL9kAAAAAACMNVfSu9r/oNcv6Sl1TsDEH3BfojqXKnlABx8sr5W0ILrlDzpWSf8r44MYB96eF51PAAAAAGCwSVRnJ8wG7X+MWCrpoi/HvLDP9y9VZ7fNa9XZIWPfx7AEhbH6UyBjkqR/av/aAuqca8qTVLTP96u+fMxUSasOeMze+amhUawdAAAAAABDHKZDr/964LIVhwpk7HW49r8yZt8rZg6LbPmDkkXSnxT5MMbe28vqnFgDAAAAAAxsVnWGJ3Zo/+PCDnVekJG8z9gDAxl7DVfnyfMDjy3fF0tQGKE/BDJCLXGzWZ3dMvY6VCBD6pz3uEidHVv3fbxLncvhpEW2fAAAAAAAjHO89j+4rVPnFS3WA8aFC2RIXx8sVx0w7qRIFT6I3aaeBSrKJL2qzgmxJyX9W9L2Hm7j8Wg8MQAAAACAqTIkVWv/48HlkiYeYmyoQMZep+rgY88XjS950OkPgYwbtH89rZJ+JSnugHGhAhl7pf//9u49WI6yzOP49xxySEgCgcASFJIgiQExkCB35Kpmragocr+JCxgsLFxhiQuCFrisSFGLiTdYFY1ylQC1WNxkUagNJshlI4IKLELIpojcEiAk5JiTnNk/OqeS3c3pt2fm7e6Znu+nagoq8/b7/IaZKurtPO/bJI1AGzd29AOH5BFakiRJkqS83MqGIyP/ZpAxoYaMARsvlm+NG1PANJKjOkM3X/pIdqNMTZlrAnA5yY2RLE0Zp216GkmSJElShcxgQ3N/2jow1JAByb2DS4HVwDvAe6Kl7Fzt0JDRAzy9PsudJCdmbEqoIWPAxiezzokXU5IkSZKkYowF9g2MydqQMWAfYFzz0bSRUcBfCN94eZrkBk1WOwL/kWHetxn8JookSZIkqRq6gXOBEYFxWRoyBkwETmo+mmiPhgyAg4EjAmOyNmRA8rv8HDCm+WiSJEmSJLWeehsyFN8swjddHqWxZ6kOAW7LMP/cpj6BJEmSJKkq6mnIUDzt0pCRRT0NGZIkVc7/fW68JEmSyvNe4JzAmMXAdOCtBuZfC5wMPBwYdxzJLhdJkiRJkiRJktQgGzIkSZJax0ySUywGUwNOAZY1UWMNyRGyKwPjvtJEDUmSJEmSJEmSOp4NGZIkSa1hDHBaYMz1wPwItRYDVwTGTCc5IlWSJEmSJEmSJEmSJDWpC9h7o5fNm8WZSfqzYdcCEyLWGw68Hqg5K2I9SZIkSVL7OYjksZbHAe8pOUsn2YP09XoNmFxauvpsyYbf0CdLziJJkiRJkjrUQtJvtNyWQ81vBmq+TPojVCRJkiRJUnxVasiQJEmSJEkq1a6Eb7QcWVLdj+RQV5IkSZIkDc6GDEmSKsIdj5IkZbM9cHpgzFLg+hwzDAPOATZLGfMO8H2gP8ccim9a4P0VwH051H0WeAKYmjLmo8CvcqgtSZIkSe3s4ySPgkxzJ9CbY4b9gXGBMQuAl3LMIEmSJEmSFMV1hHcnnJRj/asz1D8jx/rKzx2kf6/35lh7dqD2EznWliRJkqR2dSbhNfr3cqw/GVgVqP84MDTHDMqPJ2RIkiRJkjrOCOBPpC+G3wR2yaH2sYG6NeDnOdRV/rqA5aR/t1/LsX7ot9UPjM6xviRJkiS1q7I2bmS5P/E2yWMq1Z5syJAkSZIkdaTJJI8FSVsQPwZsHrHmeGBZoOZzwFYRa6o4EwjfZDkix/o7ZKj/4RzrS5Ik9AFHygAADQNJREFUSVK7ytIY8QbxN278LFCzBpwcuaaKZUOGJEmSJKljfZ7wovhbkWr1AA8HavUCe0Wqp+IdQ/r320dyky9PzwUyzMy5viRJkiS1q6I3bpweqFUDrolUS+WxIUOSJEmS1NFuIH1R3A8cFaHO7ECdGnB2hDoqz9dJ/36fLyDDLwIZrisggyRJkiS1q6I2bkwCVgTqPAVsEaGWymVDhiRJkiSpo40EniF9Ybwc2LmJGh8naexIq3FrE/OrNYSeOfxAARm+F8jwmwIySJIkSVI7y3vjxjDgiUCNlcD7mqih1mFDhiRJkiSp4+1J+FjSR2jsWNKxwOuBuZ8HRjX1CdQKHiL9e/5pARn+MZDhpQIySJIkSVI7y3vjxrWBuWvAaQ2nV6uxIUOSpIoYUnYASZLa2JMkf5H93ZQx+wGXARfUMe8Q4GZg25QxfcApwFt1zNsqLgB2KTtEZJcDixu8dnzg/SUNzluPUPZ3kezG6i0giyRJkiS1o5XA8cBvGfyRIdsAtwCHAGvqmPsE4MzAmDm05+Mm9wc+VnaIyBYA95UdQpIkSZKkqriJ8LGkn6xjvisD89WAcyNlL8MCwp+v3V77NPHfozcw94wm5s7qwECGGrBTATkkSZIkqd2dQ3h9dUUd800k2YyRNt+zwJZx4hfu7yl/TR/7dVWE/y6ekCFJkiRJ0npbkzw+JG2RvAwYl2Gu6SQNHGlz3QV0Rf0ExbIhY4MRGeY+usG56zExQ44pBeSQJEmSpCqItXFjKLAwMNdq2nu9ZkPGptmQIUlSRXSXHUCSpAp4k+T40LTjRkcDN5D+uLAdSY4XTWu2WAJ8lmThrfaX9liaAatzT5GtxujcU0iSJElSNXwBWJTyfhfJI0ZCGzdmAXsFxnwR+H32aJIkSSqSDRmSJMXxOHBhYMwhwCWDvLcZSTPGdinXrwVOIjltQ9WwVYYxvbmnyNaQMSr3FJIkSZJUDVk3blzP4Bs3jgXODtSZC1xbdzpJkiQVxoYMSZLimQ3cERhzETBtE3/+deBDgWsvBuY3kEuta/MMY1rlhIwsWSVJkiRJiccIb9w4lE1v3JhAuNHiz8CMBnJJkiSpQDZkSJIUTw04A3gxZUw3yaNL3r3Rnx1B+CbNL4F/aSacWlKrNGT0En4Mjg0ZkiRJklSfRjZu9JDcN0g7pfCvJCdwrGgqnSRJknKX9hx7SZJUvzdIboo8xOB/gb09cCPwEZJHlNxI8siSwbwCnA70x4tZqq8C25YdIrLnG7yuJ8OYIh5ZUiO5oTcsZYwNGZIkSZJUn4GNG1OBnQcZM7BxYy9gKclmjAMC854PLIwTsXTzgPPKDhHZf5YdQJIkSZKkqvsyyY2XtNelwK8DY9YChxUbXQU6lPDvZFJBWVYHcpxZUA5JkiRJqpr9SJrg09ZcDwKfItmMkTbutoKzqxx7EL5fMLm0dJIkSZIklawL+AXhxXPo9bWig6tQBxD+DUwpIEcX4Zt+pxaQQ5IkSZKqKsvGjbWB918Ati46uEphQ4YkSZIkSQGjgRdpvBnjQdIfZaL29wHCv4PQUbUxbJEhx/EF5JAkSZKkqmp248Ya4MDCU6ssNmRIklQR3WUHkCSpwpYDJwJ9DVz7KnAKsC5qIrWaNRnGbJF7imw1smSVJEmSJG1aDTgdWNzg9RcCD8eLI0mSpCLYkCFJUr5+C1xa5zX9wGeApdHTqNW8lWHMsNxTZGvIeDP3FJIkSZJUbY1u3LgHmBU/jiRJkiRJUvvrBn5J9mNILysnpkowgvDv4egCckzMkGNKATkkSZIkqRNcRPZ7BEuA7cqJqRL5yBJJkiRJkuownuTxI6HF9DxgSEkZVY7VpP8mziogw4GBDDVgxwJySJIkSVIn6CY5UTO0DusDDi4po8plQ4YkSRXhI0skSSrGRWT7/+5KksYNdY7XAu+PLSDDuMD7NWBZATkkSZIkqRPsSbZTCJcBz+ScRZIkSTmyIUOSpPwdR/ZTDqYD5+eYRa1nceD9ULNEDOMD7/8F6C0ghyRJkiRV3UjgZmBYhrFjgOuArlwTSZIkKTceiS5JUr4mAD+q85pvAgvWv6poDsnRm1XyGeDpBq9dRPoRtEU0ZIRO4VhUQAZJkiRJ6gTXALvVMX468A/AVfnEKd3fAp8tO0RkdwM3lR1CkiS1BhsyJEnKz1DgFmBUndcNAX4OTAWWxw7VAnYF9i47RGQjmrg21OzQCo8seaGADJIkSZJUdZ8DTm3guoGNGw/HjdMSdgNOLjtEZC9jQ4YkSVrPR5ZIkpSfK2m88WAs8DM8lrQTPBl4fzzNNXxksXvg/VBGSZIkSVK69wPfbvDaHpKNG6PjxZEkSVIRbMiQJCkfRwFfDIx5KfD+J4AvxImjFva7wPtDgP1yrL8DMDEwZmGO9SVJkiSp6oaTnKA5PGXMKqCW8v444IcxQ0mSJCl/NmRIkhTfWOBa0k+3eA7Yk/Bxo1cBH4iUS61pEfBGYMwHc6x/cOD9GvBEjvUlSZIkqeq+S3JCRpoZwOzAmGNw44YkSVJbsSFDkqS4hpAcI7ptyphe4HhgOXDi+n8OZijJLpqtYgVUy6kB8wJjymzIeJL036gkSZIkaXAnAGcExlwN3AxcQHjjxreAvSLkkiRJUgFsyJAkKa7LgYMCY85lw4kD/w2cFRg/EY8lrbpfBd4/ENg8p9qHBd6/P6e6kiRJklR1WdbzTwEz1/97H9k2bszFjRuSJEltIe0odUmSVJ/pwN2k///1NuC4Tfz51cDZgfnPBH7SWLSWsi8wquwQkT0KrGji+l2BZwJjjgTuaqLGpkwCng2MmUa4YUSSJEmS9L8NJTntIu00i1XAPvz/9eCxwK2B+W8had5od7sDB5QdIrKngMeanGMPkhMrQ2P+0GQdSZIkSZLawg7AyySPnxjs9TyDNyIMBRYGrl8N7JnbJ1DZfkf69z83h5rfCNR8meQxPJIkSZKk+nyf9PVWDTgt5fprMlz/d/lEVwvYg/D3P7m0dJIkSZIkFagb+DXpi+Q1hHd8TATeCszzB2B49E+gVjCT9O9+LbBLxHrDgdcCNWdFrCdJkiRJneIYwn+ZHjoBcxjhjRsrSU6YUPXYkCFJkiRJ0nr/RHiRfF7GuU7IMFfo+bNqTzsAvaR/9z+OWO/CQK1+PJFFkiRJkuo1DlhG+nrrv4AtM8z1XpLHY7pxo/PYkCFJkiRJEnA4yckFaQvku4GuOuacE5ivBpwaJb1azQ9J/97XAQdFqDMOeDtQ654IdSRJkiSpk/QAC0hfa60GptYx54mB+WrAD+LEVwuxIUOSJEmS1PG2B5aSvjheAmxX57wjgD8G5n0b2LXpT6BWMwnoI/27fxEY3USNHmB+oEYNOLSJGpIkSZLUia4ivNY6q4F5f5ph3lOai64WY0OGJEmSJKmjdQP3kb4w7gMObnD+9wOrAvP/nuSZsqqW7xC+6fIoMKqBuXuA2zPMf3tTn0CSJEmSOs/HSB79mLbWmtvg3COAPwXmduNGtdiQIUmSJEnqaBcTXhhf1GSNszLU+HaTNdR6tgFeIfzdPwNMqWPenYB5GeZdBezc/MeQJEmSpI6xE/Aa6WutP9NYY/2AyYQ3bjwODG2ihlpHloaMw4F3lfDqye9jS5IkSZIEBwBrSF8UPwBsFqHWDYE6NeDTEeqotUwH1hH+7tcCPwb2S5lrN+BKwjfuBl5nRP80kiRJklRdQ4CHSF9nrSG5l9Cszwfq1IDZEeqofFkaMsp67Z/j55YkSZIkdbhtgBdJX5i+QrJjIIaRJCchpNV7A080qKJ/pr4bIq8C9wM3AdcD9wBL6pxjTiGfTJIkSZKq4xuE11pfiljvxkCtfuCoiPVUDhsyJEmSJEkdpwu4g/RF6TpgWuS6ewO9gbqPAJtHrqtydQH/SnE3VO7Fo20lSZIkqR5HkJxcmLbWuotkfRfLSODZQM3luHGj3dmQIUmSJEnqOOcTXpRellPtczPUviKn2ipPN3A1+d9MuR0YVtBnkiRJkqQqGAMsJX2ttQTYNofaewN/DdR+BOjJobaKYUOGJEmSJKmj7Ev4Zsc8kmfH5qEL+LdA/X7gyJzqq1wzgHeIfxOlD/gqcXdrSZIkSVLVdQP/Tni99cEcM5wXqF8DLs+xvvJlQ4YkSZIkqWNsDbxA+mJ0OTA+5xzbAIsCOV4F3p1zDpVjF+BO4t1AmQ9MKfQTSJIkSVI1XEJ4zfWVnDNk2bixDvhozjmUDxsyJEmqCHdDSpIUdjHw6cCYS4C7C8hyEPCdwJgHgS8XkEXl2Af4EnA0MLzOa/uAe0l+Qw+Q3EiRJEmSJGW3L3APsFnKmPnAp0hOsszT6PW1xqSMeZXkpI5lOWdRXO8DflN2iEFMAxaWHUKSpHZhQ4YkSVJ72gI4DDgc2B2YBGwHjCQ5Pncl8DrJqSp/JLlJdz+wooSskiRJkiRJkiR1nP8BL5lddiQu6vYAAAAASUVORK5CYII=" - } - }, "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "![SEI.png](attachment:SEI.png)" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "In our simple SEI model, we consider a one-dimensional SEI which extends from the surface of a planar negative electrode at $x^*=0$ until $x^*=L^*$, where $L^*$ is the thickness of the SEI. Since the SEI is porous, there is some electrolyte within the region $x^*\\in[0, L^*]$ and therefore some concentration of solvent, $c^*$. Within the porous SEI, the solvent is transported via a diffusion process according to:\n", "$$\n", @@ -95,14 +117,20 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "### Non-dimensionalisation" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "To convert the model into dimensionless form, we scale the dimensional variables and dimensional functions. For this model, we choose to scale $x^*$ by the current SEI thickness, the current SEI thickness by the initial SEI thickness, solvent concentration with the bulk electrolyte solvent concentration, and the solvent diffusion with the solvent diffusion in the electrolyte. We then use these scalings to infer the scaling for the solvent flux. Therefore, we have\n", "$$\n", @@ -123,21 +151,30 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "### Dimensionless Model" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "After substituting in the scalings from the previous section, we obtain the dimensionless form of the model given by:" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "Solvent diffusion through SEI:\n", "\\begin{align}\n", @@ -165,14 +202,20 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "## Converting the Model into PyBaMM" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "As always, we begin by importing pybamm and changing our working directory to the root of the pybamm folder." ] @@ -180,7 +223,11 @@ { "cell_type": "code", "execution_count": 1, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "import pybamm\n", @@ -191,7 +238,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "A model is defined in six steps:\n", "1. Initialise model\n", @@ -206,14 +256,20 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "#### 1. Initialise model" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "We first initialise the model using the `BaseModel` class. This sets up the required structure for our model. " ] @@ -221,7 +277,11 @@ { "cell_type": "code", "execution_count": 2, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "model = pybamm.BaseModel()" @@ -229,14 +289,20 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "#### 2. Define parameters and variables" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "In our SEI model, we have two dimensionless parameters: $k$ and $\\hat{V}$ aswell as one dimensionless function $D(c)$ which are given in terms of the dimensional parameters, see (5). In pybamm, inputs are dimensional so we first state all the dimensional parameters and then define the dimensionless parameters in terms of them. To define the dimensional parameters, we use the `Parameter` to create parameter symbols. Parameters which are functions are defined using `FunctionParameter` and should be defined within a python function as shown. " ] @@ -244,7 +310,11 @@ { "cell_type": "code", "execution_count": 3, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "# dimensional parameters\n", @@ -267,7 +337,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "We now define the dimensionless variables in our model. Since we solve for these, we do not need to write them in terms of the dimensional variables. We simply use `SpatialVariable` and `Variable` to create the required symbols: " ] @@ -275,7 +348,11 @@ { "cell_type": "code", "execution_count": 4, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "x = pybamm.SpatialVariable(\"x\", domain=\"SEI layer\", coord_sys=\"cartesian\")\n", @@ -285,14 +362,20 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "#### 3. State governing equations" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "We can now use the symbols we have created for our parameters and variables to write out our governing equations. Note that before we use the reaction flux and solvent flux, we must derive new symbols for them from the defined parameter and variable symbols. Each governing equation must also be stated in the form `d/dt = rhs` since pybamm only stores the right hand side (rhs) and assumes that the left hand side is the time derivative. The govenering equations are then simply" ] @@ -300,7 +383,11 @@ { "cell_type": "code", "execution_count": 5, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "# SEI reaction flux\n", @@ -316,7 +403,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "Once we have stated the equations, we can add them to the `model.rhs` dictionary. This is a dictionary which stores the right had sides of the governing equations. Each key in the dictionary corresponds to the variable which the equation is being solved for." ] @@ -324,7 +414,11 @@ { "cell_type": "code", "execution_count": 6, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "model.rhs = {c: dcdt, L: dLdt}" @@ -332,14 +426,20 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "#### 4. State boundary conditions" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "We only have boundary conditions on the solvent concentration equation. We must state where a condition is Neumann (on the gradient) or Dirichlet (on the variable itself). \n", "\n", @@ -357,7 +457,11 @@ { "cell_type": "code", "execution_count": 7, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "# electrode-SEI boundary condition (x=0) (lbc = left boundary condition)\n", @@ -367,7 +471,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "On the SEI-electrolyte boundary (x=1), we have the boundary condition\n", "$$\n", @@ -377,7 +484,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "which is a Dirichlet condition and is just entered as" ] @@ -385,7 +495,11 @@ { "cell_type": "code", "execution_count": 8, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "c_right = pybamm.Scalar(1)" @@ -393,7 +507,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "We now load these boundary conditions into the `model.boundary_conditions` dictionary in the following way, being careful to state the type of boundary condition: " ] @@ -401,7 +518,11 @@ { "cell_type": "code", "execution_count": 9, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "model.boundary_conditions = {c: {\"left\": (grad_c_left, \"Neumann\"), \"right\": (c_right, \"Dirichlet\")}}" @@ -409,14 +530,20 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "#### 5. State initial conditions" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "There are two initial conditions in our model:\n", "$$\n", @@ -426,7 +553,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "which are simply written in pybamm as" ] @@ -434,7 +564,11 @@ { "cell_type": "code", "execution_count": 10, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "c_init = pybamm.Scalar(1)\n", @@ -443,7 +577,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "and then included into the `model.initial_conditions` dictionary:" ] @@ -451,7 +588,11 @@ { "cell_type": "code", "execution_count": 11, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "model.initial_conditions = {c: c_init, L: L_init}" @@ -459,14 +600,20 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "#### 6. State output variables" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "We already have everything required in model for the model to be used and solved but we have not yet stated what we actually want to output from the model. PyBaMM allows users to output any combination of symbols as an output variable therefore allowing the user the flexibility to output important quanities without further tedious postprocessing steps. \n", "\n", @@ -481,7 +628,11 @@ { "cell_type": "code", "execution_count": 12, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "model.variables = {\"SEI thickness\": L, \"SEI growth rate\": dLdt, \"Solvent concentration\": c}" @@ -489,7 +640,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "We can also output the dimensional versions of these variables by multiplying by the scalings in the Non-dimensionalisation section. By convention, we recommend include in the units in the output variables name so that they do not overwrite the dimensionless output variables. We also `.update` the disctionary so that we add to the previous output variables." ] @@ -497,7 +651,11 @@ { "cell_type": "code", "execution_count": 13, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "L_dim = L_0_dim * L\n", @@ -514,21 +672,30 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "And thats it, the model is fully defined and ready to be used. If you plan on reusing the model several times, you can additionally set model defaults which include: a default geometry to run the model on, a default set of parameter values, a default solver, etc." ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "## Using the Model" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "The model will now behave in the same way as any of the inbuilt PyBaMM models. However, to demonstrate that the model works we display the steps involved in solving the model but we will not go into details within this notebook." ] @@ -536,7 +703,11 @@ { "cell_type": "code", "execution_count": 14, - "metadata": {}, + "metadata": { + "collapsed": true, + "deletable": true, + "editable": true + }, "outputs": [], "source": [ "# define geometry\n", @@ -583,7 +754,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "Using these outputs, we can now plot the SEI thickness as a function of time and also the solvent concentration profile within the SEI. We use a slider to plot the concentration profile at different times." ] @@ -591,7 +765,11 @@ { "cell_type": "code", "execution_count": 15, - "metadata": {}, + "metadata": { + "collapsed": false, + "deletable": true, + "editable": true + }, "outputs": [ { "data": { @@ -629,14 +807,20 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "## Formally adding your model" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "deletable": true, + "editable": true + }, "source": [ "The purpose of this notebook has been to go through the steps involved in getting a simple model working within PyBaMM. However, if you plan on reusing your model and want greater flexibility then we recommend that you create a new class for your model. We have set out instructions on how to do this in the \"Adding a Model\" tutorial in the documentation. " ] @@ -658,7 +842,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.3" + "version": "3.5.2" } }, "nbformat": 4, diff --git a/examples/scripts/DFN_ambient_temperature.py b/examples/scripts/DFN_ambient_temperature.py new file mode 100644 index 0000000000..a2dfea865f --- /dev/null +++ b/examples/scripts/DFN_ambient_temperature.py @@ -0,0 +1,52 @@ +# +# Example showing how to load and solve the DFN +# + +import pybamm +import numpy as np + +pybamm.set_logging_level("DEBUG") + + +# load model +options = {"thermal": "x-lumped"} +model = pybamm.lithium_ion.DFN(options) + +# create geometry +geometry = model.default_geometry + +# load parameter values and process model and geometry + + +def ambient_temperature(t): + return 300 + t * 100 / 3600 + + +param = model.default_parameter_values +param.update( + {"Ambient temperature [K]": ambient_temperature}, check_already_exists=False +) +param.process_model(model) +param.process_geometry(geometry) + +# set mesh +var = pybamm.standard_spatial_vars +var_pts = {var.x_n: 30, var.x_s: 30, var.x_p: 30, var.r_n: 10, var.r_p: 10} +mesh = pybamm.Mesh(geometry, model.default_submesh_types, var_pts) + +# discretise model +disc = pybamm.Discretisation(mesh, model.default_spatial_methods) +disc.process_model(model) + +# solve model +t_eval = np.linspace(0, 3600 / 2, 100) +solver = pybamm.CasadiSolver(mode="fast") +solver.rtol = 1e-3 +solver.atol = 1e-6 +solution = solver.solve(model, t_eval) + +# plot +plot = pybamm.QuickPlot( + solution, ["X-averaged cell temperature [K]", "Ambient temperature [K]"] +) +plot.dynamic_plot() diff --git a/examples/scripts/rate_capability.py b/examples/scripts/rate_capability.py index abfa578445..692264313e 100644 --- a/examples/scripts/rate_capability.py +++ b/examples/scripts/rate_capability.py @@ -16,13 +16,9 @@ for i, C_rate in enumerate(C_rates): experiment = pybamm.Experiment( ["Discharge at {:.4f}C until 3.2V".format(C_rate)], - period="{:.4f} seconds".format(10 / C_rate) - ) - sim = pybamm.Simulation( - model, - experiment=experiment, - solver=pybamm.CasadiSolver() + period="{:.4f} seconds".format(10 / C_rate), ) + sim = pybamm.Simulation(model, experiment=experiment, solver=pybamm.CasadiSolver()) sim.solve() capacity = sim.solution["Discharge capacity [A.h]"] @@ -35,12 +31,12 @@ plt.figure(1) plt.scatter(C_rates, capacities) -plt.xlabel('C-rate') -plt.ylabel('Capacity [Ah]') +plt.xlabel("C-rate") +plt.ylabel("Capacity [Ah]") plt.figure(2) plt.scatter(currents * voltage_av, capacities * voltage_av) -plt.xlabel('Power [W]') -plt.ylabel('Energy [Wh]') +plt.xlabel("Power [W]") +plt.ylabel("Energy [Wh]") plt.show() diff --git a/pybamm/input/parameters/lead-acid/experiments/1C_discharge_from_full/parameters.csv b/pybamm/input/parameters/lead-acid/experiments/1C_discharge_from_full/parameters.csv index c1b130f248..40bbce3e54 100644 --- a/pybamm/input/parameters/lead-acid/experiments/1C_discharge_from_full/parameters.csv +++ b/pybamm/input/parameters/lead-acid/experiments/1C_discharge_from_full/parameters.csv @@ -4,6 +4,10 @@ Name [units],Value,Reference,Notes # Temperature,,, Reference temperature [K],294.85,Room temperature, Maximum temperature [K],333.15,, +Ambient temperature [K], 298.15,, +Heat transfer coefficient [W.m-2.K-1],10,, + + ,,, # Electrical Number of electrodes connected in parallel to make a cell,8,Manufacturer, diff --git a/pybamm/input/parameters/lithium-ion/anodes/graphite_Chen2020/graphite_LGM50_electrolyte_reaction_rate_Chen2020.py b/pybamm/input/parameters/lithium-ion/anodes/graphite_Chen2020/graphite_LGM50_electrolyte_reaction_rate_Chen2020.py index 2216a5aa0d..9b2a6c92be 100644 --- a/pybamm/input/parameters/lithium-ion/anodes/graphite_Chen2020/graphite_LGM50_electrolyte_reaction_rate_Chen2020.py +++ b/pybamm/input/parameters/lithium-ion/anodes/graphite_Chen2020/graphite_LGM50_electrolyte_reaction_rate_Chen2020.py @@ -26,7 +26,7 @@ def graphite_LGM50_electrolyte_reaction_rate_Chen2020(T, T_inf, E_r, R_g): Reaction rate """ - m_ref = 6.48E-7 + m_ref = 6.48e-7 arrhenius = exp(E_r / R_g * (1 / T_inf - 1 / T)) return m_ref * arrhenius diff --git a/pybamm/input/parameters/lithium-ion/cathodes/nmc_Chen2020/nmc_LGM50_electrolyte_reaction_rate_Chen2020.py b/pybamm/input/parameters/lithium-ion/cathodes/nmc_Chen2020/nmc_LGM50_electrolyte_reaction_rate_Chen2020.py index 065f126c91..9868f6cb6f 100644 --- a/pybamm/input/parameters/lithium-ion/cathodes/nmc_Chen2020/nmc_LGM50_electrolyte_reaction_rate_Chen2020.py +++ b/pybamm/input/parameters/lithium-ion/cathodes/nmc_Chen2020/nmc_LGM50_electrolyte_reaction_rate_Chen2020.py @@ -25,7 +25,7 @@ def nmc_LGM50_electrolyte_reaction_rate_Chen2020(T, T_inf, E_r, R_g): : double Reaction rate """ - m_ref = 3.59E-6 + m_ref = 3.59e-6 arrhenius = exp(E_r / R_g * (1 / T_inf - 1 / T)) return m_ref * arrhenius diff --git a/pybamm/input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/electrolyte_conductivity_Nyman2008.py b/pybamm/input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/electrolyte_conductivity_Nyman2008.py index 8963ebe575..168b35e6e3 100644 --- a/pybamm/input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/electrolyte_conductivity_Nyman2008.py +++ b/pybamm/input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/electrolyte_conductivity_Nyman2008.py @@ -29,9 +29,7 @@ def electrolyte_conductivity_Nyman2008(c_e, T, T_inf, E_k_e, R_g): """ sigma_e = ( - 0.1297 * (c_e / 1000) ** 3 - - 2.51 * (c_e / 1000) ** 1.5 - + 3.329 * (c_e / 1000) + 0.1297 * (c_e / 1000) ** 3 - 2.51 * (c_e / 1000) ** 1.5 + 3.329 * (c_e / 1000) ) arrhenius = exp(E_k_e / R_g * (1 / T_inf - 1 / T)) diff --git a/pybamm/input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/electrolyte_diffusivity_Nyman2008.py b/pybamm/input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/electrolyte_diffusivity_Nyman2008.py index af26eb3549..a93baf8183 100644 --- a/pybamm/input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/electrolyte_diffusivity_Nyman2008.py +++ b/pybamm/input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/electrolyte_diffusivity_Nyman2008.py @@ -28,11 +28,7 @@ def electrolyte_diffusivity_Nyman2008(c_e, T, T_inf, E_D_e, R_g): Solid diffusivity """ - D_c_e = ( - 8.794E-11 * (c_e / 1000) ** 2 - - 3.972E-10 * (c_e / 1000) - + 4.862E-10 - ) + D_c_e = 8.794e-11 * (c_e / 1000) ** 2 - 3.972e-10 * (c_e / 1000) + 4.862e-10 arrhenius = exp(E_D_e / R_g * (1 / T_inf - 1 / T)) return D_c_e * arrhenius diff --git a/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Chen2020/parameters.csv b/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Chen2020/parameters.csv index f35d901edd..a826e8176f 100644 --- a/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Chen2020/parameters.csv +++ b/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Chen2020/parameters.csv @@ -4,6 +4,8 @@ Name [units],Value,Reference,Notes # Temperature Reference temperature [K],298.15,25C, Heat transfer coefficient [W.m-2.K-1],10,, +Ambient temperature [K], 298.15,, + ,,, # Electrical Number of electrodes connected in parallel to make a cell,1,, diff --git a/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Kim2011/parameters.csv b/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Kim2011/parameters.csv index 1ed5821c2f..335b774a1d 100644 --- a/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Kim2011/parameters.csv +++ b/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Kim2011/parameters.csv @@ -4,6 +4,8 @@ Name [units],Value,Reference,Notes # Temperature Reference temperature [K],298.15,25C, Heat transfer coefficient [W.m-2.K-1],25,, +Ambient temperature [K], 298.15,, + ,,, # Electrical Number of electrodes connected in parallel to make a cell,1,, diff --git a/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Marquis2019/parameters.csv b/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Marquis2019/parameters.csv index 80a2746b6a..37f1a11a3b 100644 --- a/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Marquis2019/parameters.csv +++ b/pybamm/input/parameters/lithium-ion/experiments/1C_discharge_from_full_Marquis2019/parameters.csv @@ -3,6 +3,7 @@ Name [units],Value,Reference,Notes ,,, # Temperature Reference temperature [K],298.15,25C, +Ambient temperature [K], 298.15,, Heat transfer coefficient [W.m-2.K-1],10,, ,,, # Electrical diff --git a/pybamm/models/event.py b/pybamm/models/event.py index 5a9cafb159..0f1764a7cb 100644 --- a/pybamm/models/event.py +++ b/pybamm/models/event.py @@ -13,6 +13,7 @@ class EventType(Enum): to the discontinuity and then restart just after the discontinuity. """ + TERMINATION = 0 DISCONTINUITY = 1 diff --git a/pybamm/models/full_battery_models/lead_acid/full.py b/pybamm/models/full_battery_models/lead_acid/full.py index 955cc15d50..02d6afc23b 100644 --- a/pybamm/models/full_battery_models/lead_acid/full.py +++ b/pybamm/models/full_battery_models/lead_acid/full.py @@ -126,4 +126,3 @@ def set_side_reaction_submodels(self): self.submodels["negative oxygen interface"] = pybamm.interface.NoReaction( self.param, "Negative", "lead-acid oxygen" ) - diff --git a/pybamm/models/submodels/electrode/base_electrode.py b/pybamm/models/submodels/electrode/base_electrode.py index a4e786f32a..6b2368f0d6 100644 --- a/pybamm/models/submodels/electrode/base_electrode.py +++ b/pybamm/models/submodels/electrode/base_electrode.py @@ -190,4 +190,3 @@ def _get_standard_whole_cell_variables(self, variables): ) return variables - diff --git a/pybamm/models/submodels/electrolyte/base_electrolyte_diffusion.py b/pybamm/models/submodels/electrolyte/base_electrolyte_diffusion.py index 0f4df2a527..cf99fa04fe 100644 --- a/pybamm/models/submodels/electrolyte/base_electrolyte_diffusion.py +++ b/pybamm/models/submodels/electrolyte/base_electrolyte_diffusion.py @@ -105,8 +105,10 @@ def _get_standard_flux_variables(self, N_e): def set_events(self, variables): c_e = variables["Electrolyte concentration"] - self.events.append(pybamm.Event( - "Zero electrolyte concentration cut-off", - pybamm.min(c_e) - 0.002, - pybamm.EventType.TERMINATION - )) + self.events.append( + pybamm.Event( + "Zero electrolyte concentration cut-off", + pybamm.min(c_e) - 0.002, + pybamm.EventType.TERMINATION, + ) + ) diff --git a/pybamm/models/submodels/external_circuit/__init__.py b/pybamm/models/submodels/external_circuit/__init__.py index 21966d10b5..f181a9cca9 100644 --- a/pybamm/models/submodels/external_circuit/__init__.py +++ b/pybamm/models/submodels/external_circuit/__init__.py @@ -8,4 +8,3 @@ LeadingOrderVoltageFunctionControl, LeadingOrderPowerFunctionControl, ) - diff --git a/pybamm/models/submodels/external_circuit/base_external_circuit.py b/pybamm/models/submodels/external_circuit/base_external_circuit.py index 9c69e00eab..908a9463c3 100644 --- a/pybamm/models/submodels/external_circuit/base_external_circuit.py +++ b/pybamm/models/submodels/external_circuit/base_external_circuit.py @@ -50,4 +50,3 @@ def get_fundamental_variables(self): Q = pybamm.Variable("Leading-order discharge capacity [A.h]") variables = {"Discharge capacity [A.h]": Q} return variables - diff --git a/pybamm/models/submodels/external_circuit/current_control_external_circuit.py b/pybamm/models/submodels/external_circuit/current_control_external_circuit.py index 0368852226..952a7f87bc 100644 --- a/pybamm/models/submodels/external_circuit/current_control_external_circuit.py +++ b/pybamm/models/submodels/external_circuit/current_control_external_circuit.py @@ -34,4 +34,3 @@ class LeadingOrderCurrentControl(CurrentControl, LeadingOrderBaseModel): def __init__(self, param): super().__init__(param) - diff --git a/pybamm/models/submodels/external_circuit/function_control_external_circuit.py b/pybamm/models/submodels/external_circuit/function_control_external_circuit.py index 6e3812f99a..93d1d7653f 100644 --- a/pybamm/models/submodels/external_circuit/function_control_external_circuit.py +++ b/pybamm/models/submodels/external_circuit/function_control_external_circuit.py @@ -91,4 +91,3 @@ class LeadingOrderPowerFunctionControl(LeadingOrderFunctionControl): def __init__(self, param): super().__init__(param, constant_power) - diff --git a/pybamm/models/submodels/particle/base_particle.py b/pybamm/models/submodels/particle/base_particle.py index f4234d6458..98ed2a050a 100644 --- a/pybamm/models/submodels/particle/base_particle.py +++ b/pybamm/models/submodels/particle/base_particle.py @@ -89,4 +89,3 @@ def set_events(self, variables): pybamm.EventType.TERMINATION, ) ) - diff --git a/pybamm/models/submodels/porosity/base_porosity.py b/pybamm/models/submodels/porosity/base_porosity.py index 1aebacf3df..8d9157635e 100644 --- a/pybamm/models/submodels/porosity/base_porosity.py +++ b/pybamm/models/submodels/porosity/base_porosity.py @@ -96,25 +96,33 @@ def _get_standard_porosity_change_variables(self, deps_dt, set_leading_order=Fal def set_events(self, variables): eps_n = variables["Negative electrode porosity"] eps_p = variables["Positive electrode porosity"] - self.events.append(pybamm.Event( - "Zero negative electrode porosity cut-off", - pybamm.min(eps_n), - pybamm.EventType.TERMINATION - )) - self.events.append(pybamm.Event( - "Max negative electrode porosity cut-off", - pybamm.max(eps_n) - 1, - pybamm.EventType.TERMINATION - )) - - self.events.append(pybamm.Event( - "Zero positive electrode porosity cut-off", - pybamm.min(eps_p), - pybamm.EventType.TERMINATION - )) - - self.events.append(pybamm.Event( - "Max positive electrode porosity cut-off", - pybamm.max(eps_p) - 1, - pybamm.EventType.TERMINATION - )) + self.events.append( + pybamm.Event( + "Zero negative electrode porosity cut-off", + pybamm.min(eps_n), + pybamm.EventType.TERMINATION, + ) + ) + self.events.append( + pybamm.Event( + "Max negative electrode porosity cut-off", + pybamm.max(eps_n) - 1, + pybamm.EventType.TERMINATION, + ) + ) + + self.events.append( + pybamm.Event( + "Zero positive electrode porosity cut-off", + pybamm.min(eps_p), + pybamm.EventType.TERMINATION, + ) + ) + + self.events.append( + pybamm.Event( + "Max positive electrode porosity cut-off", + pybamm.max(eps_p) - 1, + pybamm.EventType.TERMINATION, + ) + ) diff --git a/pybamm/models/submodels/thermal/base_thermal.py b/pybamm/models/submodels/thermal/base_thermal.py index 030a631297..9f532cadae 100644 --- a/pybamm/models/submodels/thermal/base_thermal.py +++ b/pybamm/models/submodels/thermal/base_thermal.py @@ -30,6 +30,9 @@ def _get_standard_fundamental_variables(self, T, T_cn, T_cp): T_x_av = self._x_average(T, T_cn, T_cp) T_vol_av = self._yz_average(T_x_av) + T_amb_dim = param.T_amb_dim(pybamm.t * param.timescale) + T_amb = param.T_amb(pybamm.t * param.timescale) + q = self._flux_law(T) variables = { @@ -64,6 +67,8 @@ def _get_standard_fundamental_variables(self, T, T_cn, T_cp): + param.T_ref, "Heat flux": q, "Heat flux [W.m-2]": q, + "Ambient temperature [K]": T_amb_dim, + "Ambient temperature": T_amb, } return variables diff --git a/pybamm/models/submodels/thermal/isothermal/isothermal.py b/pybamm/models/submodels/thermal/isothermal/isothermal.py index f4d6a7f446..a10fc80e9e 100644 --- a/pybamm/models/submodels/thermal/isothermal/isothermal.py +++ b/pybamm/models/submodels/thermal/isothermal/isothermal.py @@ -22,8 +22,8 @@ def __init__(self, param): super().__init__(param) def get_fundamental_variables(self): - - T_x_av = pybamm.PrimaryBroadcast(self.param.T_init, "current collector") + T_amb = self.param.T_amb(pybamm.t * self.param.timescale) + T_x_av = pybamm.PrimaryBroadcast(T_amb, "current collector") T_n = pybamm.PrimaryBroadcast(T_x_av, "negative electrode") T_s = pybamm.PrimaryBroadcast(T_x_av, "separator") T_p = pybamm.PrimaryBroadcast(T_x_av, "positive electrode") @@ -33,6 +33,7 @@ def get_fundamental_variables(self): T_cp = T_x_av variables = self._get_standard_fundamental_variables(T, T_cn, T_cp) + return variables def get_coupled_variables(self, variables): diff --git a/pybamm/models/submodels/thermal/x_full/base_x_full.py b/pybamm/models/submodels/thermal/x_full/base_x_full.py index 22ef947378..3b939d8a72 100644 --- a/pybamm/models/submodels/thermal/x_full/base_x_full.py +++ b/pybamm/models/submodels/thermal/x_full/base_x_full.py @@ -25,6 +25,7 @@ def get_fundamental_variables(self): T = pybamm.standard_variables.T T_cn = pybamm.BoundaryValue(T, "left") T_cp = pybamm.BoundaryValue(T, "right") + variables = self._get_standard_fundamental_variables(T, T_cn, T_cp) return variables diff --git a/pybamm/models/submodels/thermal/x_full/x_full_no_current_collector.py b/pybamm/models/submodels/thermal/x_full/x_full_no_current_collector.py index a919e6003b..f70acd4824 100644 --- a/pybamm/models/submodels/thermal/x_full/x_full_no_current_collector.py +++ b/pybamm/models/submodels/thermal/x_full/x_full_no_current_collector.py @@ -35,11 +35,18 @@ def set_boundary_conditions(self, variables): T = variables["Cell temperature"] T_n_left = pybamm.boundary_value(T, "left") T_p_right = pybamm.boundary_value(T, "right") + T_amb = variables["Ambient temperature"] self.boundary_conditions = { T: { - "left": (self.param.h * T_n_left / self.param.lambda_n, "Neumann"), - "right": (-self.param.h * T_p_right / self.param.lambda_p, "Neumann"), + "left": ( + self.param.h * (T_n_left - T_amb) / self.param.lambda_n, + "Neumann", + ), + "right": ( + -self.param.h * (T_p_right - T_amb) / self.param.lambda_p, + "Neumann", + ), } } diff --git a/pybamm/models/submodels/thermal/x_lumped/base_x_lumped.py b/pybamm/models/submodels/thermal/x_lumped/base_x_lumped.py index fa93592e60..88225fdf04 100644 --- a/pybamm/models/submodels/thermal/x_lumped/base_x_lumped.py +++ b/pybamm/models/submodels/thermal/x_lumped/base_x_lumped.py @@ -33,6 +33,7 @@ def get_fundamental_variables(self): T_cp = T_x_av variables = self._get_standard_fundamental_variables(T, T_cn, T_cp) + return variables def get_coupled_variables(self, variables): diff --git a/pybamm/models/submodels/thermal/x_lumped/x_lumped_0D_current_collectors.py b/pybamm/models/submodels/thermal/x_lumped/x_lumped_0D_current_collectors.py index 1d2025a16a..d8a75059be 100644 --- a/pybamm/models/submodels/thermal/x_lumped/x_lumped_0D_current_collectors.py +++ b/pybamm/models/submodels/thermal/x_lumped/x_lumped_0D_current_collectors.py @@ -13,11 +13,12 @@ def __init__(self, param): def set_rhs(self, variables): T_av = variables["X-averaged cell temperature"] Q_av = variables["X-averaged total heating"] + T_amb = variables["Ambient temperature"] cooling_coeff = self._surface_cooling_coefficient() self.rhs = { - T_av: (self.param.B * Q_av + cooling_coeff * T_av) / self.param.C_th + T_av: self.param.B * Q_av + cooling_coeff * (T_av - T_amb) / self.param.C_th } def _current_collector_heating(self, variables): diff --git a/pybamm/models/submodels/thermal/x_lumped/x_lumped_1D_current_collectors.py b/pybamm/models/submodels/thermal/x_lumped/x_lumped_1D_current_collectors.py index fa96d98e70..a897f92348 100644 --- a/pybamm/models/submodels/thermal/x_lumped/x_lumped_1D_current_collectors.py +++ b/pybamm/models/submodels/thermal/x_lumped/x_lumped_1D_current_collectors.py @@ -15,15 +15,21 @@ def __init__(self, param): def set_rhs(self, variables): T_av = variables["X-averaged cell temperature"] Q_av = variables["X-averaged total heating"] + T_amb = variables["Ambient temperature"] cooling_coeff = self._surface_cooling_coefficient() self.rhs = { - T_av: (pybamm.laplacian(T_av) + self.param.B * Q_av + cooling_coeff * T_av) + T_av: ( + pybamm.laplacian(T_av) + + self.param.B * Q_av + + cooling_coeff * (T_av - T_amb) + ) / self.param.C_th } def set_boundary_conditions(self, variables): + T_amb = variables["Ambient temperature"] T_av = variables["X-averaged cell temperature"] T_av_left = pybamm.boundary_value(T_av, "negative tab") T_av_right = pybamm.boundary_value(T_av, "positive tab") @@ -35,11 +41,11 @@ def set_boundary_conditions(self, variables): self.boundary_conditions = { T_av: { "negative tab": ( - self.param.h * T_av_left / self.param.delta, + self.param.h * (T_av_left - T_amb) / self.param.delta, "Neumann", ), "positive tab": ( - -self.param.h * T_av_right / self.param.delta, + -self.param.h * (T_av_right - T_amb) / self.param.delta, "Neumann", ), "no tab": (pybamm.Scalar(0), "Neumann"), diff --git a/pybamm/models/submodels/thermal/x_lumped/x_lumped_2D_current_collectors.py b/pybamm/models/submodels/thermal/x_lumped/x_lumped_2D_current_collectors.py index 811e20c9e1..fa3cc51af7 100644 --- a/pybamm/models/submodels/thermal/x_lumped/x_lumped_2D_current_collectors.py +++ b/pybamm/models/submodels/thermal/x_lumped/x_lumped_2D_current_collectors.py @@ -15,6 +15,7 @@ def __init__(self, param): def set_rhs(self, variables): T_av = variables["X-averaged cell temperature"] Q_av = variables["X-averaged total heating"] + T_amb = variables["Ambient temperature"] cooling_coeff = self._surface_cooling_coefficient() @@ -25,9 +26,9 @@ def set_rhs(self, variables): T_av: ( pybamm.laplacian(T_av) + self.param.B * pybamm.source(Q_av, T_av) - + cooling_coeff * pybamm.source(T_av, T_av) + + cooling_coeff * pybamm.source(T_av - T_amb, T_av) - (self.param.h / self.param.delta) - * pybamm.source(T_av, T_av, boundary=True) + * pybamm.source(T_av - T_amb, T_av, boundary=True) ) / self.param.C_th } diff --git a/pybamm/models/submodels/thermal/x_lumped/x_lumped_no_current_collectors.py b/pybamm/models/submodels/thermal/x_lumped/x_lumped_no_current_collectors.py index db66a8f36f..d67f0b47a8 100644 --- a/pybamm/models/submodels/thermal/x_lumped/x_lumped_no_current_collectors.py +++ b/pybamm/models/submodels/thermal/x_lumped/x_lumped_no_current_collectors.py @@ -27,13 +27,14 @@ def __init__(self, param): def set_rhs(self, variables): T_av = variables["X-averaged cell temperature"] Q_av = variables["X-averaged total heating"] + T_amb = variables["Ambient temperature"] # Get effective properties rho_eff, _ = self._effective_properties() cooling_coeff = self._surface_cooling_coefficient() self.rhs = { - T_av: (self.param.B * Q_av + cooling_coeff * T_av) + T_av: (self.param.B * Q_av + cooling_coeff * (T_av - T_amb)) / (self.param.C_th * rho_eff) } diff --git a/pybamm/models/submodels/thermal/xyz_lumped/base_xyz_lumped.py b/pybamm/models/submodels/thermal/xyz_lumped/base_xyz_lumped.py index f331416774..e4cbb17696 100644 --- a/pybamm/models/submodels/thermal/xyz_lumped/base_xyz_lumped.py +++ b/pybamm/models/submodels/thermal/xyz_lumped/base_xyz_lumped.py @@ -45,11 +45,12 @@ def get_coupled_variables(self, variables): def set_rhs(self, variables): T_vol_av = variables["Volume-averaged cell temperature"] Q_vol_av = variables["Volume-averaged total heating"] + T_amb = variables["Ambient temperature"] cooling_coeff = self._surface_cooling_coefficient() self.rhs = { - T_vol_av: (self.param.B * Q_vol_av + cooling_coeff * T_vol_av) + T_vol_av: (self.param.B * Q_vol_av + cooling_coeff * (T_vol_av - T_amb)) / self.param.C_th } diff --git a/pybamm/parameters/standard_parameters_lead_acid.py b/pybamm/parameters/standard_parameters_lead_acid.py index e48faefcb1..be73bc36e0 100644 --- a/pybamm/parameters/standard_parameters_lead_acid.py +++ b/pybamm/parameters/standard_parameters_lead_acid.py @@ -165,7 +165,7 @@ # Fake thermal -Delta_T = pybamm.Scalar(0) +Delta_T = pybamm.thermal_parameters.Delta_T # -------------------------------------------------------------------------------------- @@ -431,7 +431,8 @@ def c_p_init(x): # Thermal effects not implemented for lead-acid, but parameters needed for consistency T_init = pybamm.Scalar(0) Theta = pybamm.Scalar(0) # ratio of typical temperature change to ambient temperature - +T_amb_dim = pybamm.thermal_parameters.T_amb_dim +T_amb = pybamm.thermal_parameters.T_amb # -------------------------------------------------------------------------------------- "5. Dimensionless Functions" @@ -499,4 +500,3 @@ def U_p(c_e_p, T): current_with_time = ( dimensional_current_with_time / I_typ * pybamm.Function(np.sign, I_typ) ) - diff --git a/pybamm/parameters/standard_parameters_lithium_ion.py b/pybamm/parameters/standard_parameters_lithium_ion.py index a841dfe713..969ca016c7 100644 --- a/pybamm/parameters/standard_parameters_lithium_ion.py +++ b/pybamm/parameters/standard_parameters_lithium_ion.py @@ -379,6 +379,9 @@ def chi(c_e): / (pybamm.thermal_parameters.rho_eff_dim * F * Delta_T * L_x) ) +T_amb_dim = pybamm.thermal_parameters.T_amb_dim +T_amb = pybamm.thermal_parameters.T_amb + # Initial conditions T_init = pybamm.thermal_parameters.T_init c_e_init = c_e_init_dimensional / c_e_typ diff --git a/pybamm/parameters/thermal_parameters.py b/pybamm/parameters/thermal_parameters.py index 4041c4f5b1..5273d369f3 100644 --- a/pybamm/parameters/thermal_parameters.py +++ b/pybamm/parameters/thermal_parameters.py @@ -110,3 +110,14 @@ h = h_dim * pybamm.geometric_parameters.L_x / lambda_eff_dim T_init = (T_init_dim - T_ref) / Delta_T + +# -------------------------------------------------------------------------------------- +# Ambient temperature + + +def T_amb_dim(t): + return pybamm.FunctionParameter("Ambient temperature [K]", t) + + +def T_amb(t): + return (T_amb_dim(t) - T_ref) / Delta_T # dimensionless T_amb diff --git a/pybamm/solvers/solution.py b/pybamm/solvers/solution.py index 9e4855ca66..0adec1a630 100644 --- a/pybamm/solvers/solution.py +++ b/pybamm/solvers/solution.py @@ -331,4 +331,3 @@ def append(self, solution, start_index=1, create_sub_solutions=False): copy_this=solution, ) ) - diff --git a/tests/integration/test_solvers/test_idaklu.py b/tests/integration/test_solvers/test_idaklu.py index 69aabb5861..7f0a1bd1fc 100644 --- a/tests/integration/test_solvers/test_idaklu.py +++ b/tests/integration/test_solvers/test_idaklu.py @@ -72,4 +72,3 @@ def test_changing_grid(self): debug = True pybamm.settings.debug_mode = True unittest.main() - diff --git a/tests/unit/test_parameters/test_parameter_sets/test_LGM50_Chen2020.py b/tests/unit/test_parameters/test_parameter_sets/test_LGM50_Chen2020.py index 7fdd260e51..915db76e96 100644 --- a/tests/unit/test_parameters/test_parameter_sets/test_LGM50_Chen2020.py +++ b/tests/unit/test_parameters/test_parameter_sets/test_LGM50_Chen2020.py @@ -23,8 +23,8 @@ def test_load_params(self): electrolyte = pybamm.ParameterValues({}).read_parameters_csv( pybamm.get_parameters_filepath( - "input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/" + - "parameters.csv", + "input/parameters/lithium-ion/electrolytes/lipf6_Nyman2008/" + + "parameters.csv", ) ) self.assertEqual(electrolyte["Reference temperature [K]"], "298.15") @@ -34,9 +34,7 @@ def test_load_params(self): "input/parameters/lithium-ion/cells/LGM50_Chen2020/parameters.csv", ) ) - self.assertAlmostEqual( - cell["Negative current collector thickness [m]"], 12E-6 - ) + self.assertAlmostEqual(cell["Negative current collector thickness [m]"], 12e-6) def test_standard_lithium_parameters(self): diff --git a/tests/unit/test_parameters/test_parameter_values.py b/tests/unit/test_parameters/test_parameter_values.py index ba2ef688d7..75b1e4109f 100644 --- a/tests/unit/test_parameters/test_parameter_values.py +++ b/tests/unit/test_parameters/test_parameter_values.py @@ -14,9 +14,12 @@ def test_read_parameters_csv(self): data = pybamm.ParameterValues({}).read_parameters_csv( pybamm.get_parameters_filepath( os.path.join( - "input", "parameters", - "lithium-ion", "cathodes", - "lico2_Marquis2019", "parameters.csv" + "input", + "parameters", + "lithium-ion", + "cathodes", + "lico2_Marquis2019", + "parameters.csv", ) ) )