{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Tutorial 3 applying our knowledge to solve Chem Reactor problems\n", "\n", "So we already know how to do the basic operations such as variable operations, loops, conditions, numpy operations, plotting\n", "\n", "Now it is time to use all of those together with other cool python modules to solve engineering problems.\n", "\n", "If you don't know how to do those problems it's ok, start small and make your way up.\n", "\n", "The easiest example we could do would be to analyze the chemical reaction kinetics in the simplest case:\n", "\n", "## Problem 1:\n", "\n", "*A -> B *\n", "\n", "We assume it is a first order reaction and we use the design equation from the Batch reactor:\n", "\n", "So that means that $-r_A = k C_A$.\n", "\n", "Lets assume $k = 10^{-3} 1/s$ and $C_A(t=0) = 10$ moles\n", "\n", "So we can expect that the half-time of A should be around 10000 seconds which is ~3 hours.\n", "\n", "Ok lets see that through Python.\n", "\n", "We want to track the concentration which we will call $y_A$.\n", "\n", "$dy_A/dt = - k C_A = -k y_A$\n", "\n", "Main thing that we need to use here is the ODEINT module:\n", "\n", "`from scipy.integrate import odeint`\n", "\n", "In this module if you have an equation $\\frac{dy}{dt} = f(y,t)$ you need:\n", "\n", "- To define at least the Python function f which implements the mathematical function $f$. It must take a vector y and a time t, and return the new vector f(y,t).\n", "\n", "\n", "Some additional commands that you may need in your script:\n", "\n", "```python\n", "#lets import the modules # plotting everything inline\n", "# %matplotlib inline \n", "import matplotlib.pyplot as plt # plotting modules\n", "plt.style.use('presentation') # just have in your script for prettier plotting\n", "import numpy as np # our matlab-like module\n", "from scipy.integrate import odeint # integration of ODEs so you don't have to write your finite difference yourself\n", "```\n", "\n", "To plot your data use these commands:\n", "\n", "```python\n", "plt.plot(times, y_A_exact,'k--', label='analytical')\n", "plt.xlabel('$t [s]$')\n", "plt.ylabel('$C_A [moles]$')\n", "plt.legend()\n", "plt.savefig('odeint.pdf')\n", "plt.show()\n", "\n", "```\n", "\n", "\n", "\n" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/Users/bazilevs/anaconda3/lib/python3.6/site-packages/matplotlib/figure.py:2022: UserWarning: This figure includes Axes that are not compatible with tight_layout, so results might be incorrect.\n", " warnings.warn(\"This figure includes Axes that are not compatible \"\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAboAAAEYCAYAAAAqIzNgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3XmYFOW1x/HvGWg2AVc0Kq6oMagwAkYUF3AF9eIGcUOjXsU9iSGrZsEsJjEachU14h6MSwSNGsU1IEYligKixoVRjKhRRAWGRYaZc/94q2eKpmemp6dnevt9nqee7n7rrepTxQxnqupdzN0REREpVRX5DkBERKQtKdGJiEhJU6ITEZGSpkQnIiIlTYlORERKmhKdiIiUNCU6EREpaUp0IiJS0pToRESkpHXMdwD5sNlmm/n222+f7zBERCQDL7300qfu3ivb7csy0W2//fbMnj0732GIiEgGzOy91myvW5ciIlLSlOhERKSkKdGJiEhJK6hndGbWG/ghMAjoD3QFdnD3hSn1ugC/BMYAGwFzgR+6+8w2C66qCq66Cu64A6qroXt3GDMGxo2DPn3a7GtFpLDU1NSwaNEiVq9ene9QSk6XLl3o3bs3iUQip/u1QpqPzsyGAvcALwEdgMNIn+j+AhwJfB94B7gAGAHs4+5zm/ueQYMGeYsao0ybBqNGQU1NWJISibBMmQIjRmS+PxEpWu+++y49evRg0003xczyHU7JcHeWLFnC8uXL2WGHHdZZZ2YvufugbPddaLcuZ7r7Fu5+BHBvugpm1h84GbjY3W9096eAbwD/AX6R84iqqkKSW7ly3SQH4fPKlWF9VVXOv1pECs/q1auV5NqAmbHpppu2yZVyQSU6d6/LoNpIoIZw5Zfcbi1wN3C4mXXOaVBXXbV+gktVUwMTJuT0a0WkcGWc5Kqq4PzzoWdPqKgIr+efrz+MG9FWfzwUVKLL0G7Au+6+MqX8NaATsFNOv+2OO9ZJdB8ADwO/Aeqzck0NTJ6c068VkSI3bRr06wc33QTLl4N7eL3pplA+bVq+IywbxZjoNgE+T1P+WWz9esxsrJnNNrPZixcvzvzbqqvX+TgIOAq4BKhqop6IlLE2fuSx7777Nlvnj3/8IytXpl4P5N7ChQu58847W7zd6aefzpQpU9ogovUVY6IzIF0Lmiaved19krsPcvdBvXq1YCSZ7t3X+VgZez+3iXoiUsba+JHHc88912ydXCa6tWvXNrou20TXnoox0X1G+qu2jWPrc2fMmNCyMpI20SUScOqpOf1aESkCZumX66/PLNFde23j+2hC9+gP6xkzZjB06FBGjRrFrrvuyimnnIK7c/XVV/Phhx8ybNgwhg0bBsDjjz/OPvvsw4ABAxg9ejTV0V2oRx55hF133ZX99tuPb33rWxx11FEAjB8/nrFjx3LYYYdx2mmnsXDhQvbff38GDBjAgAED6pPtj370I5555hkqKyuZMGECtbW1fP/732evvfaiX79+3HDDDUBoVXnhhRfSt29fjjzySD755JOsT3uLuXtBLsBZhCu37VPKfwasAbqllI8HvgQ6N7fvgQMHesYWLHDv1s093GH3u0JMDvgRUZl36xbqiUjJe/311xs+JP8PaIulCRtssIG7u0+fPt179uzp77//vtfW1vrgwYP9mWeecXf37bbbzhcvXuzu7osXL/b999/fq6ur3d39t7/9rV922WW+atUq7927t7/zzjvu7n7iiSf6kUce6e7uP//5z33AgAG+cuVKd3dfsWKFr1q1yt3d33rrLU/+Pzp9+vT6bdzdb7jhBv/lL3/p7u6rV6/2gQMH+jvvvONTp071Qw45xNeuXesffPCBb7jhhn7vvfc2fX7rTzOzvRX5pBiv6B4EEsDoZIGZdQROAB539y9z+m19+oR+ct26QSKx/hVdt25hvTqNi0gefP3rX6d3795UVFRQWVnJwoUL16sza9YsXn/9dYYMGUJlZSW333477733Hm+88QY77rhjfb+1k046aZ3tRo4cSdeuXYHQUf7ss89mjz32YPTo0bz++utp43n88cf585//TGVlJXvvvTdLlizh7bffZubMmZx00kl06NCBrbbaioMOOii3J6IJBTUyCoCZjYreDoxeR5jZYmCxuz/t7nPN7B7gj2aWAN4FzgN2AE5pk6BGjIBXXoEJE9j5z3+m6/LlrAI+BD6ZOpXNhw9vk68VkQLX2IAb558fWlc2dfsykYCxY2HixFaF0LlzQ4+qDh06pH2e5u4ceuih3HXXXeuUz5kzp8l9b7DBBvXvJ0yYwBZbbMG8efOoq6ujS5cuabdxd6655hoOP/zwdcofeeSRvPU9LMQrunuj5dzo83XR58tidc4AbgV+RWjtvw0w3N1fbrOo+vSBiRPpsGwZ/WL/+PPmNjsQi4iUm3Hj1nm2n1YiARdf3GYh9OjRg+XLlwMwePBgnn32WRYsWADAypUreeutt9h1111555136q8C77nnnsZ2x9KlS9lyyy2pqKhg8uTJ1NbWrvc9AIcffjjXX389NVGSf+utt1ixYgUHHHAAd999N7W1tXz00UdMnz69LQ47rYJLdO5ujSxDY3VWuft33f0r7t7F3fd29xntFWPlttvWv5/7zDPt9bUiUixSHnmsI5Fol0ceY8eOZcSIEQwbNoxevXpx2223cdJJJ9GvXz8GDx7MG2+8QdeuXbnuuusYPnw4++23H1tssQUbbrhh2v2df/753H777QwePJi33nqr/mqvX79+dOzYkf79+zNhwgTOOuss+vbty4ABA9h9990555xzWLt2Lcceeyw777wze+yxB+eddx4HHnhgmx17qoIa67K9tHisyxR/OvlkzotuAZy83Xb8Jc09cREpTf/+97/52te+llnlqqrQhWDy5IbB4E89NVzJFchz/erqarp37467c8EFF7DzzjtzcRteaTYn3fkttbEui0Jl1FwXYO5//5vHSESkoEWPPFi6FGprw+vEiQWT5ABuvPFGKisr2W233Vi6dCnnnHNOvkPKuYJrjFIM9hg5ksqxY+kPDFyzBlatgqhlkohIMbn44ovzegXXHpTosrDBFlswZ5dd4K23QqurV1+FvfbKd1giIpKGbl1mqzLWo04tL0VECpYSXbbiiW7evPzFISIiTVKiy5au6EREioISXZbW7r47vwCOA/Z4/nnqmhjdW0Sk0Nx2221ceOGFzdb58MMP6z+fddZZjQ791ZQZM2bUDxadD0p0Weqw9dZca8b9wKt1dbzz9NP5DklEJKdSE91NN91E37598xhRdpTosmQVFVRu0jBb0NxHHsljNCKSL+PHj8fMMlrGjh273vZjx45dp8748eOb/c5jjjmGgQMHsttuuzFp0iQgTN1z6aWX0r9/fwYPHszHH38MwEMPPcTee+/NnnvuySGHHFJfnrR8+XJ22GGH+iG7li1bxvbbb8+9997L7NmzOeWUU6isrGTVqlUMHTqU5GAbjz76KAMGDKB///4cfPDBALzwwgvsu+++7Lnnnuy77768+eabWZ/XXFKia4X+O+5Y/37uv/6Vx0hEpJzccsstvPTSS8yePZurr76aJUuWsGLFCgYPHsy8efM44IADuPHGGwHYb7/9mDVrFnPmzOHEE0/kiiuuWGdfPXr0YOjQoTz88MMA3H333Rx//PGMHj2aQYMG8Ze//IW5c+fWz2IAsHjxYs4++2ymTp3KvHnzuPfeewHYddddmTlzJnPmzOEXv/gFl1xySTudkaapH10rVO61F7z4IgBz3347z9GISLm4+uqruf/++wF4//33efvtt+nUqVP9c7CBAwfyxBNPALBo0SJOOOEEPvroI9asWVM/JU/cWWedxRVXXMExxxzDrbfeWp8kGzNr1iwOOOCA+n1tEt3dWrp0Kd/85jd5++23MbP6q8R80xVdK1TGpqGY++mneYxERPJl/PjxGU8AmrzNGDdp0qR16jR363LGjBk8+eSTPP/888ybN48999yT1atXk0gk6qfBiU/Xc9FFF3HhhRcyf/58brjhBlavXr3ePocMGcLChQt5+umnqa2tZffdd28yBndPO+XOT3/6U4YNG8arr77KQw89lPa78kGJrhV2OewwkjMyfVBXx6cFcj9aRErX0qVL2XjjjenWrRtvvPEGs2bNarb+1ltvDcDtt9/eaL3TTjuNk046iTPOOKO+LHUKnqR99tmHp59+mnfffReAzz77bL3vuu2221p0XG1Jia4VOnbpwh7xuekefDCP0YhIORg+fDhr166lX79+/PSnP2Xw4MFN1h8/fjyjR49m//33Z7PNNmu03imnnMLnn3++zizjp59+Oueee259Y5SkXr16MWnSJI477jj69+/PCSecAMAPfvADfvzjHzNkyJD6+eoKgabpaaWxX/saN77xBgBXHnUU4x56KCf7FZHC1KJpeorIlClTeOCBB5g8eXJe42iLaXrUGKWV+vfrB1Gim/vaa3mORkSk5S666CKmTZvGIyXaTUqJrpUqhw2Dv/4VgLmxjpUiIsXimmuuyXcIbUrP6Fqp/9FHczUwE/hnTQ0USCsjEWk75fjIpz201XlVomul7ltuyUU77cT+wIZ1daDblyIlrUuXLixZskTJLsfcnSVLltClS5fmK7eQbl3mQmUlLFgQ3s+dCwMH5jceEWkzvXv3ZtGiRSxevDjfoZScLl260Lt375zvV4kuFyorYcqU8F5T9oiUtEQikXZ0ESlcunWZC9HcdJ8DLz/7bH5jERGRdSjR5cDnFRVsB2wCDJ0zh7oePeD886GqKt+hiYiUPSW61po2jY2OP56V0cflwMLqarjpJujXD6ZNy2d0IiJlT4muNaqqYNQobNUq+seK5wLU1MDKlTBqlK7sRETyqGgTnZkNMbPHzewTM1tmZi+b2ZntGsRVV4WEBlTGitdpjlJTAxMmtGdUIiISU5SJzsz6AU8CCeBs4HjgReBmMzuv3QK54460iW5OvE5NDeR57DgRkXJWrN0LTgQ6AP/j7tVR2RNm1h84Dbi+XaKorq5/G+85NwtwwNLUExGR9lWUV3RAJ6AGWJVS/gXteUzdu9e//SqwcfT+U+DtRuqJiEj7KtZEd1v0erWZbWVmG5nZ2cDBQPs9EBszBhIJIJzIfWOrnku+SSTg1FPbLSQREVlXUSY6d38VGAocDXxA6Kt9LXCuu9+dbhszG2tms81sds6G7hk3rj7RAQyJrarvNp5IwMUX5+b7RESkxYoy0ZnZzsBU4DXgf4BDgD8BfzKzU9Jt4+6T3H2Quw/q1atXbgLp0ycM/dWtGyQS61/RdesW1vfpk5vvExGRFivKRAdcTnhGd5S7/93dn3L3bwF/Bf7PzNrvuEaMgFdegbFj2atHDxLATsDeQO306WG9iIjkTbEmuj2Aee5ek1L+ArApsHm7RtOnD0ycSLdly/h4wADeBm4BOnzySbuGISIi6yvWRPdfoNLMOqWU7w2sBj5r/5CCjQ84oOHDc881XlFERNpFsSa6icAOwENmdrSZHWZmE4GTgOvdfU3eIhsSa5KimQxERPKuKBOdu08BjgA6AzcRGqbsB1wAfD+Poa2b6F54AdbkL+eKiEiRJjoAd5/m7kPdvZe793D3Sne/zt1r8xrYllsyZ6utuBI4bvVq5txzT17DEREpd0Wb6ArZVZ068X3gfmDG1Kn5DkdEpKwp0bWBIYMH179/bvbsPEYiIiJKdG1gyPHH179/9qOP8Lq6PEYjIlLelOjawG4jR9Izev9RXR0L//nPvMYjIlLOlOjaQIdOnRi86ab1n5+76648RiMiUt6aTXRmVpuD5WftcTCFZEi/fvXvn505M4+RiIiUt0wmXjXgPWBhFvs34IBma5WgIcOHw/TpADz7zjt5jkZEpHxlOsP4re7+i2y+wMzKsiXG18eMoeKHP6QOmL96Ncvef5+e22yT77BERMqOntG1kR5bbUX/rl0BcGDW5Mn5DUhEpExlkuh6Ab9vxXe0dvuiNWSnnerfP/voo3mMRESkfDV769Ldl7TmC1q7fTEbdcwxfGX+fIYAe3XokO9wRETKUotvXZrZr81stpltklKeMLPeuQut+B141llcCgwFNnjxRVi7Ns8RiYiUn2ye0R0CdHD3+jnfzGxP4APgPTN7x8yG5ii+4rbtttA7yv0rVoSZyEVEpF1lk+i2B1IHcPw1sBnwLLAR8Hcz69O60EpErD8dAwdCz55w/vlQVZW/mEREykg2ia4n8HHyg5ltBBwK/NXdDwC+Tug/94OcRFjMpk2DJ58E4HNgLsDy5XDTTSEBTpuWz+hERMpCNonuA8LVW9Jh0X4mAbj7AuBBQvIrX1VVMGoUi9asYXdgE2Bkcl1NDaxcCaNG6cpORKSNZZPoXgEON7NkM8JTgJVAfJyrd4GtWhlbcbvqKqip4SuEkwHwfrTUq6mBCRPaPTQRkXKSTaKbAGwHPGVmNwJHAY+4e7xJ4ebA6hzEV7zuuANqaugI7B0rfi5ep6YG1JFcRKRNtTjRufszwE+AIcD/AsuAy1Kq7QZ81Oroill1df3bIbHi6U3UExGR3MtqCDB3vxzYmtBFbCd3fz25zsx2JDRIeSEXARat7t3r3x4cK36MMCRYunoiIpJ7WY916e6fuPvMNCOf9ARuB6a2KrJiN2YMJBIA7Ask09lC4K1knUQCTj213UMTESknWSc6M9vVzL5pZqeb2UHJcnef6+5nuvuDuQmxSI0bV5/oOrHuVV39qJeJBFx8cTsHJiJSXrIZAqzCzG4FXgNuiZbHY+std+EVsT59YMoU6NYNEgmGx1Y9CqF8ypRQT0RE2kw2V3Q/BL4JPA+cS7hFGU9u+5rZQjMbloP4ituIEWHYr7FjObxbt/riGcCqmTPDehERaVPZJLozgDeBA939RuDV+Ep3fxZYC3yj9eGVgD59YOJEdlixgq926gSEfhczNW2PiEi7yCbRbQdMc/faJuq8RGiDITHD+/ZlK+BMYPPXX2+uuoiI5EA2iW4Z0LmZOh8AW2ax7xYxsyPMbKaZVZvZsmj6oIOa3zI/Lv/5z1kE3Azs+fLL+Q5HRKQsZJPoXgQOMbOmtq0DNswupMyY2TnAA4Srx2OB0cC9QLemtsunbocfjnXpEj688QYsXJjXeEREykE2ie5mYGfgF03U6Qe02cziZrY98Efg++5+sbs/4e6Pufvv3P3vbfW9rda1Kxx4YMPnxx7LXywiImUimyHApgL3AD82s6nALvH1ZnYsYXLWf+YkwvTOJFw1/qkNv6NtDG/oaFCnaXpERNpcth3GTwGuB44BTgQws+lmNh+YAtQAv8tJhOntB7wBnGhmVWa21swWmNkFbfidObH2kEOYBBwHbPfgg9SsWJHvkERESlq2Y13WufuFhJaVdwGfAgcSBnOeDxzl7i/lLMr1bUW4ffp74LeEOfGeACaa2bfTbWBmY6PGKrMXL17chqE1rUPfvlzeoQP3A4vcef7mm/MWi4hIOch6CDAAd/+Xu49x9y2ArkA3d6909ydyE16jKoAewDnufqO7/8PdzyMMOvLjdKOzuPskdx/k7oN69erVxuE1zioqGL5Lw93eR+++O2+xiIiUg2YTnZmdaGa7NFfP3b909/aagy7Z0CU1oT4ObEE7dG1ojeFHH13//tG5c/MYiYhI6cvkiu5O4N9mttTMnjazP5jZGDPrm8dxLV9rpDwZT117BZKNgy64gI7R+zmrVvHfefPyGo+ISCnLJNF9j5DsFhHmEP0OYRqe+cByM3vWzK6OZjHoZ2Yd2i7cevdHr4enlB8OLHL3/7ZDDFnr2bs3QzZs6Gb4+LXX5jEaEZHS1rG5Cu7+h+R7M+sGVAIDgQHRshewDw3ziX4Ztb58yd3Pz3nEwSOEybpvMLPNgHeAUYRGKWe00Xfm1PDBg3k66kf36GOPcVqe4xERKVXm7s3XamoHZp2B/qyb/HYHOrp7m13dmVlP4DeEBLcxobvBb939zua2HTRokM+ePbutQsvI3HvuYc8TTwRgUzM+Xr2aDtGgzyIi0sDMXnL3QVlv39pEl3anZglgD3cvyAEdCyHReV0dW3XsyH+j8z8L2LtHjzAz+bhxmqdORCTS2kTXmhnGNzSzYWZ2SGqrTHevKdQkVyjsscfWecD4KMDy5XDTTdCvH2jUFBGRnMgq0ZnZj4CPgCeBxwitMj8ys1+bWY9cBliSqqpg1CiGR1dznYEvkutqamDlShg1KtQTEZFWaXGiM7PTgMuBlcBkwuDK9xCa9P8YmGNmO+QyyJJz1VVQU8NwQquaz4AJqXVqamDCeqUiItJCLX5GZ2YvEzpl93P3JbHyCkKLxz8CHwP93b0gB3LM+zO6nj3DbcpM6i1d2vbxiIgUsHw8o9sVuD+e5KB+/MubCQM97wCMyzaoklddndt6IiLSqGwS3Uqg0aG+3P0pwnO747MNquR1757beiIi0qhsEt2rwMHN1JlHuKqTdMaMgUSi/uNawiCdZwJ/SRYmEnDqqe0fm4hIickm0d0G9DezHzZRZ6vswikT48atk+huIIxddith+nYgrL/44vaPTUSkxGQzw/hthFkDLjezu8ysMr7ezA4CvgG8kJMIS1GfPjBlCnTrBokEx9IwGvUM4KMuXcJ6dRoXEWm1bDuMHw08CJwAvGRmn5vZHDNbSEiCHYHLchNiiRoxAl55BcaOZauePRkaFTtw7wEHhPUiItJq2c4wvtrdjwVGAg8TLkj6A70JV3KHu/szOYuyVPXpAxMnwtKlnHjKKfXFdz33XB6DEhEpLa2dYfzv7j7S3TcCegJd3X0fd/9HbsIrH8f/7Gf1U0nMqq7m3Zkz8xqPiEipaFWii3P3anevydX+ys2mu+zCYb161X/+629+k8doRERKR2sGde5kZhvnMphyd+Jxx9W/v2vGjPwFIiJSQrIZ63IbM/sHoeP4p2b2hZnNMLM/mNkYM+trZtbcfmR9R//kJ3SJ3s9bvZp///3veY1HRKQUZHNFdx0wFPgAeIowruV+wHeA24H5wHIzezZHMZaNnr17c+RWDV0Q7/n97/MYjYhIacgm0e0PvAjs5O6HuftXgQ2BAwjJ7s/AO8BeOYuyjJwYzTpuwPtz50IbTIwrIlJOskl0XwIz4g1P3H2Fu//T3a9x9zPcvR+ggRqzcOQll/DHTp1YBNy8bBnMnZvvkEREilo2ie5J4KvNVXL3NVnsu+x13XRTvj16dMMYaoMHQ0VFmLLn/PM1GauISAtlk+h+DRxsZl/PdTAS2Xnnhvdr1oTbl8uXw003Qb9+MG1a/mITESky2Yx1+TpwEnC/mZ1kZh1yH1YZq6qCK65Iv66mBlauhFGjdGUnIpKhbLoXbAGcC2wO3AH818zuNbMfmtkhZrZJroMsK1ddBTU1OPAKcAmwXtfxmhqYMKHdQxMRKUbmLWzVZ2b3EWYR/wJYCGwHJDuOJ3f2HjDb3b+RmzBza9CgQT579ux8h5Fez56wfDn/JDRvBdgCWAT1Q4TV11u6tL2jExFpd2b2krsPynb7js1XWc9BhL5yQ9y9OgpiO2BgbBmAZhjPTnU1APsAXwH+S+io+AApJzSqJyIiTcumMUot8GgyyQG4+3vufp+7X+ruw919czTDeHa6h14ZHYCzYsUTG6knIiJNyybRPQPs2Fwld/9PFvuWMWPqZx8/h5DwIEzIOj9ZJ5GAU09t99BERIpRNonuMuAIM9st18G0hpk9amZuZr/KdyytMm5cfaLrDRwXW3Vt8k0iARdf3M6BiYgUp2wS3UmETuNPmNmhOY4nK2Z2EmHi1+LXpw9MmQLdukEiwYWxVZOBL7p2Dev79MlXhCIiRSWbRPc94ChCW4lHzWyBmV1vZv9rZpVmlk0Dl6yZ2UbABOC77fm9bWrECHjlFRg7lv179GCPqHglcNuQIWG9iIhkJJtEdxDwfeBO4E1ge8LjpEnAS4SZC14ws+tyFWQzrgBec/e72un72kefPjBxIrZsGRfFnsdd+/TT1K1dm8fARESKSzYjo8xw9z+4+xh37wv0BIYA3wJuA94g3EY8J5eBpmNm+wGnAee39Xfl08lXXslG0RR/C2pqeOzyy/MckYhI8ch6hvEkd1/p7s+7+7Xu/r/uvifQgzaepsfMEsANwJXu/mYG9cea2Wwzm7148eK2DC3nNth8c84cMAAIHRQTDz+c34BERIpIqxNdOu6+xt1fbot9x/wQ6EoYZLpZ7j7J3Qe5+6BevXq1bWRt4DsTJvAcMBs45MUXYcGCfIckIlIUmk10Zva6mWV9a7C12zeyz22BS4GfAp3NbKOoUQqxzyU12PQ2++/PPkceiUGYzeC69noEKiJS3DK5otsV2KwV39Ha7dPZEehCGFT689gCoVXo51DfWLF0XBjrbDBhguapExHJQKZdAYZa1BgiCy0bNTozc4FhacqnE5LfzUDp3durrQUzcMeBj935SnKeuttvD/3r1PVARGQdGSe6aCkI7v4FYVSsdUTJ+D13X29d0auqgm98g5Xu3EoY+7IjYSofq6kJU/eMGhX636kzuYhIvUwSXborp5ZamIN9lLdonro1wA8InccBZgIHJusk56mbuN4Q0CIiZavF89GVgoKej64x0Tx1AOcBf4qKjwXuS62neepEpIS0dj66NuleIG0gNv9cfPzL+wkPLNPVExERJbriEZt/bjfCFO9JP2uknoiItHGiM7PtzexYM7usLb+nLMTmqQP4JZBsB/sQMAs0T52ISBo5SXRm1tHM+pnZN81sgplNN7PPgSpgKutOli3ZiM1TB7A7Yb6kpJ+A5qkTEUmj1YnOzF4CqgmPin4N9AVeJox3+Q1gQ3ffurXfU/ZS5qkDGE/DDORPAdNPPlldC0REUuTiiq4vcCWwsbv3dvfD3X0coaP4v919eQ6+Q2Cdeero2ZOdKyo4PdaR/9K778br6vIYoIhI4clFohtAmKlgipn1y8H+pCnRPHUsXQq1tfzs2WfpFK2aX11N1S235DU8EZFCk4tpev7t7ocD1wP3m9mtZqZble1k23324TsDBjAOeAfY6brrQFd1IiL1ctbq0t3vI7R8/w9hZKoKINHkRpITv3v4Ya7s2pVeAHPmhC4GGvBZRARb3s43AAAZAUlEQVTIcfcCd1/t7j8HBgIPAk+Z2Y/MbINcfo+k+MpX4MgjGz6vWhWm8kkO+NyvH0yblr/4RETyKKNEZ2ZdzexpM7vZzDo3V9/dF7r7scDJwGnAu62MU5pSVQUps45/mXxTUwMrV4YBn3VlJyJlKNMrutOB/YEF7v5lM3XrufvjQD/gdy0PTTJ21VWwdi0Q7hufRZiMb028TnLAZxGRMpPRoM5m9gRQCfRuKtFZmCfnbsIFxbfd/fPG6uZTUQ7q3JRowOcaYHvgw6h4InBBaj0N+CwiRaa9BnXuDzze3NWch6x5G3AKcGRTdSWHooGcE8C3Y8WXAh+lqSciUk4yTXQbAe9lUtHdpwEfAEdlG5S0UGwg54uAnaL3S4FvNVJPRKRcZJrovgA2bMF+/0noaiDtITbgc1dgUmzVFELzVw34LCLlKtNE9x9Co5JMvQ9s2fJwJCspAz4PA86IrT4fWNaxowZ8FpGylGmiewLY18x2z7B+AtB9svaSZsDnK4HNo9UfAJdsvrkGfBaRspRporsZqAHuzLDz9y7A4qyjkpZLGfB5k4oK/i92lXfde+/x/KRJTexARKQ0ZZTo3H0BcAVhGrRZZva1xuqa2VeBw4AXchKhZC5lwOcTVq/miF69gDCVxNnnnceaHj00PJiIlJWMhwBz958BtxMamcw1sxvN7EAz6wKhD52ZDQP+Rpgm7ca2CFgyZxUVXP/AAyQvwXesq2N5dbWGBxORstKxJZXd/Qwzmw/8Cvhf4EzAzWwZocFfJ8CAm9390VwHKy237eabc02HDnSvrWUU4R+nXk1NWEaNCrc99QxPREpQiwd1dvc/ADsCvwfejvaxEdCZ0NfuO+5+di6DlFa46irOqKhgNClJLk7Dg4lICctoCLAmdxAap/QCVrh7UTRAKbkhwJoSDQ+WUT0NDyYiBai9hgBrlLuviGYrKIokV3bSDPtVDZwKTG+mnohIKcjpfHRSgFKG/aoCBgN3ACcAixqpJyJSKpToSl1seDCALjR0cFwMjAbWdOyo4cFEpGQVZaIzs1FmNtXM3jOzVWb2ppn9xsx65Du2gpMyPNjWwD2E/h8As4Dv1tZqeDARKVlFmeiA7wG1wCXAcOB64DzgCTMr1mNqG2mGBxsK/DZW5Vp37hg5MjRIUWdyESkxrW51mQ9m1iu18YuZnUbo0H6wu/+jqe3LqtVlUlVV6EIweTJUV+MbbMDoujqmrlgBhE6Qs4iN3J1IhGXKlDC8mIhInuS91WU+NNLC88Xodev2jKVopAwPZsuWcevjj7NrtHoVcBxhPiYg9K1buTJ0JteVnYgUsaJMdI04MHr9d16jKCI97riD+zp0qJ9moorQ7aA2XkmdyUWkyBXlrctUZrY1MAeY5+6HNlJnLDAWYNtttx343nsZTZhe2qLO5PcC34gV/wb4UWo9dSYXkTwpy1uXcWbWHXgAWMu6842uw90nufsgdx/UKxrRv+xFncRHAz+Iig4GLmqknohIMWrRoM6FJpo54UHC2JsHuvuiZjaRuO7d64cH+w2wLWGk7i7p6omIFKmivaIzswQwFfg6cIS7z89zSMUn1pm8AriANEkOYNUqdTsQkaJVlIku6iv3F8KdtqPdfVaeQypOKZ3J05kIXFlToznsRKRoFWWiA64lPFq6ElhhZoNjS+88x1Y80nQmj/sd4Xnd94HrkoXqdiAiRaZYE12yB/OlwPMpy1n5CqoojRgRJl0dO7ZhZJREghogfs12AXBTfDt1OxCRIlES3QtaqixHRmmJqNvBcuBQ4F+xVZcCvyD6C0ndDkSkHZR99wJpA1F3gh7Ao8CesVW/Bk4mjKTCsmUaH1NECp4Snawv1p1gI2AGYeTspHuAg4CPITRQUUMVESlgSnSyvpQ57HoCDxGe0yXNAvYGXotvp4YqIlKAlOhkfWm6HXQkdDX4Pxp+aN4DRgI1qduroYqIFBAlOllfE90OvkUYb20DoBNwG7Bex4SamjAdkIhIAVCik/TSdTvo2ROAo4BnCT32929sezVUEZECoe4F0jJR14PG/JXQIvM0wOIrNJGriGRJ3QukfaU0VIn7D2EepNOBUcCn8ZVqqCIieaJEJy3TxPiY44Fk9/H7gN2Bh1MrqaGKiLQzJTppmSYaqlwDnBf7/DHhed6xwIJkYU0NTJqk53ci0m6U6KTlGmmosgFh8OdHgK/Eqv8N6AuMA76AkOzU0VxE2okSnWSnTx+YODGMdVlbG1579ADCiNvzgVNj1WuAPwA7EaaeWBvfl57fiUgbUqKT3Ik1VNkM+DNhQOh9Y1WWEG5xpm3ru2YNHHOMbmuKSE4p0UnupGmo8nXgn4TxMbeLyq4kTSdzgLVr4dVXdVtTRHJKiU5yp5GGKgZ8A3gDuAM4MmWzWuAQ4I9Adeo+dVtTRFpJiU5yq7ERVRIJugCnkNKRnNBY5SngYmBb4CfAR6n7XbUKvvY13dIUkRZTopPcS9dQ5ayzGu1/d3Ps/eeEOe96E678pgBfQriVWVOjW5oi0mJKdNI+muhoPhW4HugTK6sjdFMYDWwJXAi8RKwRS/KW5nHHhUYwasAiIo1QopP20URH867AucCbhLEyh6Zs+jmhS8Ig4K7U/a5eDXfeqQYsItIoJTppP409v9t9d+jYkQ6EK7jpQBXwcxpaagJ0Bo5I2WUN8KT7unPi6WpPRGI0e4HkX1VVuAJbuXK9VXXA04R57zqy7vM8gCeBQ4ENgeHAQcAwQsd0AzALV3pJyVkUrrkGZs+GO+6A6mro3j0kxXHjwtWniBSM1s5eoEQnhWHatNCFoKYmLBn6FqEDeqqtCLdAk0t94ovr2DH03UtKJKBDBzj4YJg5UwlQpEBomh4pDeluayYS4YqsCVsQuiSk+hC4kzBt0C6EsTevTa20du26n2tqwjO/hx9O/8zvllvCrU/dChUpKrqik8LVxC3NOAdeIfTFmwHMpGG6oLhbgDNSyn5BuD3aF9gN2Bno1NSXpbsKbO5WKMBVV+k2qUiWdOsyC0p0RaSxW5qpz95iaoF5hKQ3g4bE9wqwR6yeE8bk/CxW1pGQ7JKJbydCt4c+wOakuf0Zly4JVkQ3Terq1o2/uQQ5alRoparkKNLqRIe7l90ycOBAlyKyYIH7BRe49+zpXlERXseMce/a1T2kuyaXWvA3wNemlFeFXJfxsjBl++Xg94LPAl+UZv8ZLx07rvu5Q4d1X5NLIuHerZv7zTe7n3eee48e7mbh9bzzwnlasKDl6556Krf7W7Ag3z8xUmKA2d6K//PznnSyDhy2IQycsRRYRpjUettMtlWiKxGPPBL+408k1k0IZhklmM/AbwX/HvgI8O2aSHKJNIns5ZQ6HcC3Ad8bfCT42eA/Ab8G/IFsk2CmyTGRcO/cOSyp56OpdU0l1Wz215pk3BYJtxC+q1DiKJRjzuIPodYmuqK8dWlm3Qh3p74kDI3owK+AbkA/d1/R1Pa6dVlCqqpgwgSYPLnhFt/IkTB1ahgfs4WWEwaffg34N/AOoU9fB+DFlLoPASMz3O9OwNspZfcB3wM2jZaNomXDaIm/3wbo37JDyb+W3Mrt0CEMF5d8zWSbbNe153cVShyFcszJW/ZTpoQGaBkqy2d0ZvZtwjyeX3X3BVHZDoT/S37g7n9oanslujLQ2LO91P98W+EZwpRDi4D3gcVN1B1CmK4obiJwUYbfNRJ4IKVsPPB7wl93yWWD2Puu0dIFOAw4MWX7JwhJvHO0dIq9JlJeewO9UrZfSvgLM0F4ttmR0Iy76XayIoQRkl55JePnza1NdB2z3TDPRgKzkkkOwN3fNbNngaMJSVDKWbK7QurV3qmnwoABcNFFLe6zl2r/aElaDXxAmHnh45RlxzTbL2nBd/VIU7YMWBktzdmI9RPdLcDdGX7/7wlXn3EjCQ19UnWgIfF1iJabgONS6h0K/IeQHDukebXofQVwHetf0f4PYWSceD2LLfHPEwldTJJWAOdE762RJb7uOtZtjfsh8MtYndTX+PtuwG9TYn8T+FPsc+o+4u+3AL6fsv1LhDkeU7dJZYS7CWellD8NPNbMdkmVhBGL4v4OzGpi+7hhwMGphTU14Xdz4sQM99I6xZrodmP9P3Ah3HFK/TeRcpWcRSHdL9OBB66fBA84AJ56KlzxZXEV2IWGFpqZ+C5wMiHhLSFcISWXL1I+75lm+/Xm7mtC1zRlX7Zg+3RdLho7I7XR8mUzdauAdzP8/uVpyp4g82O4IuXzl8BfMtwWQqKM+4x1E1VTNmb9RPcfwvyLmejL+onuVcIfH5k4hPUT3fPAbzLc/lTW/0/1cdIP1JBOBY0kusmTleiasQlhrN9UnxF+rtZjZmMJ/YfZdtt0XYylrDSWBNM982vqKrAVt0J7RMvOWR7Cnwi3LlbGlhWx19WxZY802x9OuB35JbAmev2ScJW0Jva6hnBVkWoDoGdUby0hudU1EmuHNGW1acoak25ki5Y8dEm96mnpAxvdjm0D1S35U611ijXRQfqf1UZ/Ht19EjAJwjO6tgpKilxLrwKbSoJNPZTPwbPCCqB7tGTjnOarNOnxNGV1NFzRraUhAaaLcQYhidbGtotv79FnB3ZPs/3fY9t4bBtP2dYJfSDjugOTU+rHF1I+p/5HuSXhdma8bnybeFnnNLHvQsPzldT6qe83TbP9QOB3NJ2wk+u2S7PuAMK8j01tl5Tuj6QjWf+cNma/xlZ0z/Ynt+WKtTHKx8Df3P2clPLrgNHunvrcfB1qjCI519iV4MUXh/UtSZCNtZDLYUMakbxKJMJwfxneuizXVpf/ADq5+34p5TMIx3RgU9sr0UnBaCxBHncc3Hdfbq4e27PZuZKxZKKdW13mveN3NgvwHcJdkR1jZdsTHheMa257dRiXopZupJgLLmjopNvSdU89ldv93Xxz+o787d1xvRC+q1DiKJRjTg4o8MgjLfqRpxxHRiE8B18AzCd0JxhJ6ED+DtC9ue2V6ETaWCEk3EL5rkKJo1COWSOjZM7MtgUmELrjGGHw+u+4+8LmttWtSxGR4lGuHcZx9/8Ax+c7DhERKWyaeFVEREqaEp2IiJS0on1G1xpmthh4rxW72Az4NEfhFDudi0DnoYHORaDz0KC152I7b6Z/dFPKMtG1lpnNbs2D0VKicxHoPDTQuQh0Hhrk+1zo1qWIiJQ0JToRESlpSnTZmZTvAAqIzkWg89BA5yLQeWiQ13OhZ3QiIlLSdEUnIiIlTYlORERKmhJdhsxsGzObYmZLzWyZmd0XjbdZdMxslJlNNbP3zGyVmb1pZr8xsx4p9TY2s5vM7FMzW2FmT5rZevMwmlkXM/u9mX0U7e95MzsgTb0KM/uxmS00s9VmNs/MCmoYNzN71MzczH6VUl4W58LMjjCzmWZWHf2czzazg2LrS/48mNkQM3vczD6JzsHLZnZmSp2cH5+ZnW1mb5jZl9Hv5LltdYxpvru3mV0THcfK6Hdg+zT18nbcZnaMmc2J9veemf3EzNJNXr++1owIXS4L0A14G3gVOIYwY8J8oArYIN/xZXE8s4C/AqcABxKmPfoiKq+I6hjwDLAIOAkYDjxN6PTZO2V/f4m2Pxs4GLgPWAVUptT7NfAl8D1gGHADYSLoI/J9TqL4TgI+Ikyy/KtYeVmcC8Kk4zU0DJZ+OPBD4KhyOQ9AvyjO6dHv+aFRTA6c11bHF+2nLqo/DPhV9Pm8tj7m6PuHAh8DjwCPRce7fZp6eTnu6GexltCoZRjwXWA18LuMji9fv1TFtADfjk7yTrGyHQhz4n033/FlcTy90pSdFv1wHxR9Pjr6PCxWZ0PgM+DqWFn/qN4ZsbKOwJvAg7GyzaMf/MtSvvcp4JUCOCcbAf8l/AeemuhK/lwQ5nNcRZgBpLE65XAeLgfWkDLdF+GPwOfb4viibT8Bbk+pdwvhj4hEOxx3Rez9WaRJdPk8bmAO8HRKvZ9F/1Zfae74dOsyMyOBWe6+IFng7u8CzxJ++YuKuy9OU/xi9Lp19DoS+NDdp8e2Wwo8xLrHPJJwFXBPrN5a4G7gcDPrHBUfDnQC7kj53juAPcxsh+yOJmeuAF5z97vSrCuHc3Em4S/pPzVRpxzOQydC7KtSyr+g4VFPro9vH6BXmnqTgU2B/bI9mEy5e10G1fJy3Ga2DVDZSL0EMKK5wJXoMrMb4bZlqteAvu0cS1s5MHr9d/Ta1DFva2bdY/XedfeVaep1AnaK1fuSMGFuaj3I43k0s/0IV7TnN1KlHM7FfsAbwIlmVmVma81sgZldEKtTDufhtuj1ajPbysw2MrPkbboJ0bpcH99u0Wvquc3770aKfB132nrRxcZKMjg/SnSZ2QT4PE35Z8DG7RxLzpnZ1sAvgCfdPTkjbVPHDA3H3Vy9TWKvX3h0z6GJeu3KzBKE5wdXuvubjVQrh3OxFbAz8Hvgt8BhwBPARDP7diyukj4P7v4q4XnV0cAHhOO4FjjX3e+OxZXL40u+pu4zr78baeTruBurlyxr9vwU7cSreZCuZ721exQ5Fv0V/gDheeMZ8VVkdsy5rtfefgh0JTwMb0w5nIsKoAdwurvfF5X9I2p592Mzu5oyOA9mtjMwlXBVcS7hFubRwJ/MbLW7/4W2OQ80UreQ5Ou4m6qX0c+KEl1mGvurYWPS/5VRFMysC/AgsCNwoLsviq3+jMaPGRqO+zMgXTeLjWPrk68bm5ml/KWXWq/dWOgecinh4Xvn2DMGos8bAcspg3MBLCFc0T2RUv44oXXllpTHebic8BzqKHevicqeMrNNgf8zs7vI/fHFr2A+itXbJGV9vuXruJu6st2IDM6Pbl1m5jUa7hPH9QVeb+dYciK6ZTcV+Dqhye/8lCpNHfN/3L06Vm8HM+uWpt4aGu7TvwZ0BvqkqQf5OY87Al0ID7k/jy0QmkV/DuxBeZyL1xopT/7FXEd5nIc9gHmxJJf0AqGBxObk/viS5z713ObzPKSTr+NOWy+629CNTM5PWzdbLYWF0M9sLbBjrGx7wl9+4/IdXxbHU0HoR7caOLiROscQbhUcGCvrSfjL/5pYWWVU75uxso6ERi0PxcqSTY5/nvI9TwLz83QeNiI8j0ldnNCiayjQvUzOxZFR7KNSyh8D3i+jn4kZwDtAp5TyOwm3MTvl+vgILQcXA7em1LspOredsj2eLM9BY90L8nbcwFxgekq9n5Bh94J2/0EqxgXYgPDXynzC/fqRwLzoF6J7vuPL4niuj35gfwUMTll6R3UqgOeA94ETCU2GZxBuE2yTsr+7CVc/ZxFap00hJNEBKfV+G5V/l5BEridcKfxPvs9JSpzOuv3oSv5cEK7c/hH9B3MuoTHKpOhcnF5G52FUdMyPRb/rhwETo7I/tNXxRee8LvqdHEpoHFYHXNDOxz6Khv8fzos+H5jv4waOiMpviOpdHO3/9xkdWz5+mIpxIdybngosIzy3+RtpRg4ohgVYGP0gp1vGx+ptQui8+RmhGe9TQP80++sK/IHQ4Xo18C9gaJp6HQh/hb1H+IvvFVKuIAphISXRlcu5IFydXUsYIWNNFNPJZXgeRhAS+OLod30uoetJh7Y8PsLING9F9d4Gzs/Dz326ZUYhHDdwHOEC40vgP4QO4x0yOTZN0yMiIiVNjVFERKSkKdGJiEhJU6ITEZGSpkQnIiIlTYlORERKmhKdiIiUNCU6EREpaUp0IiJS0pToRIqEmQ01M48tb7RiX5ul7EsjR0jJUqITKUBm9t0oAZ2cZvXTwGWEMRiztTLax2WEYZpESpbmoxMpTAOi15fSrJvh7uNbs3N3XwmMh3ClCGzXmv2JFDJd0YkUpoFANWGQWxFpBSU6kQJiZr+LnpftSpgLrzb2DK13hvvY38zuM7MqM1ttZp+Y2Qtmdnlbxi5SqHTrUqSwvATcDnyTMPfbE7F1lnaLeAWzS4BfE6YxeQz4FNgCGAQMBy7JcbwiBU+JTqSAuPtfzWwjQqK73d0nJddFz9IaZWZbECau/Cdh5vg1Kes3y33EIoVPty5FCk+yIcrLLdxuV8KEl2+lJjkAd/+0tYGJFCMlOpHCMwCoAea3cLvXgKXAmWb2oJmdYGYb5zw6kSKjRCdSQMysI7AH8Lq7f9mSbaMrtv2AKcDBwN3AYjObZmYDmtxYpIQp0YkUlr5AF1p+2xIAd3/V3UcDGwOHAvcRGqE8bmadcxalSBFRohMpLJXR65zW7MTd17j7k+7+DULjlE0JrS9Fyo4SnUhh2TR6XdaSjcxsTzPrk6Z8J2B3QneDRa0PT6T4qHuBSGFJDvn1azPbHVgBvObu9zaz3beAb5rZC4RGKZ8AOwAjo/VnuntdWwQsUuiU6EQKiLvPNLOLgOTSGbgcaC7RPUD4ff46MJrwnO9D4E7gd+6uocSkbCnRiRQYd59IC2cmcPe/AX9rm4hEipu5axoqkWIQjYwyPVb0prvvmuW+NgMWx8vcvdkhxkSKka7oRIrHQsL8cUmtGekkOR+dSMnTFZ2IiJQ0dS8QEZGSpkQnIiIlTYlORERKmhKdiIiUNCU6EREpaUp0IiJS0v4fsdOels8YVaYAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#lets import the modules # plotting everything inline\n", "%matplotlib inline \n", "import matplotlib.pyplot as plt # plotting modules\n", "plt.style.use('presentation') # just have in your script for prettier plotting\n", "import numpy as np # our matlab-like module\n", "from scipy.integrate import odeint # integration of ODEs so you don't have to write your finite difference yourself\n", "\n", "#reaction order\n", "def f(y, t):\n", " return -k*y\n", "\n", "k = 0.001 # L/(mol*sec)\n", "time_start = 0 #s\n", "time_finish = 10000 #s\n", "N_points = 50\n", "\n", "times = np.linspace(time_start,time_finish,N_points)\n", "y_A0 = 10. # moles\n", "\n", "y_A_calc = odeint(f, y_A0, times)\n", "plt.plot(times, y_A_calc, 'ro-', label='integrated')\n", "\n", "\n", "y_A_exact = y_A0*np.exp(-k*times)\n", "plt.plot(times, y_A_exact,'k--', label='analytical')\n", "\n", "\n", "plt.xlabel('$t [s]$')\n", "plt.ylabel('$C_A [moles]$')\n", "plt.legend()\n", "plt.savefig('odeint.pdf')\n", "plt.show()\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now lets do the same calculation as the case number 1 but we will use the conversion factor $X_A$ as our variable:\n", "\n", "$X_A = \\dfrac{N_{A0} - N_A}{N_{A0}} = \\dfrac{C_{A0} - C_A}{C_{A0}}$ if the volume $V=const$\n", "\n", "So the design equation looks like:\n", "\n", "$\\dfrac{d X_A}{dt} = (-r_A/N_{A0}) V $\n" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Once deleted, variables cannot be recovered. Proceed (y/[n])? y\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAboAAAEYCAYAAAAqIzNgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XucVXW9//HXZxhkO5j3g+coIl7ChHKSo6JcdEAFLz+J\nA2o3qcTStEPlSfGcvHCRrETLorSbUF5+FQIqpYQXmBRR5PzSKdEKL2khNhkIxTjizHx+f6y1me1m\n32Zm71lr7/1+Ph77MXt/11p7fdd3ZvZnf7/rezF3R0REpFLVRJ0BERGRUlKgExGRiqZAJyIiFU2B\nTkREKpoCnYiIVDQFOhERqWixCnRmdpCZzTezNWa23cw6zGxQgcf2M7N5ZvaambWE7zGm1HkWEZF4\ni1WgA44AzgE2A48CXRnktwC4ELgaOAvYBKwws6OLnUkRESkfFtcB42Z2IfAD4FB3fzXPvvXA08Cn\n3P32MK0PsB74vbtPKnV+RUQknuJWo+uuicAOYFEywd3bgZ8BE8ysb1QZExGRaFVKoBsKvOzurWnp\n64HdCJpERUSkClVKoNsX2JIhfXPKdhERqUKVEuhEREQyqo06A0WyBcg0DCFZk9ucYRtmFs+eOCIi\nFcbdLapzV0qNbj1wqJkl0tKHEXRSeSHbge6uR8pj5syZkeehJ4+2tjbuX7SIL595Jtc2NPDlM8/k\ngbvvZtOmTVw6YgQrEwk6CMatdAAP9+vHGf37888wLfXxAPBI+Hxmlu0PZUh34Mvh+3d127VZ0qPY\nlu8ash3b3Wvv7rYozplt28wY5SVOv4uoVUqg+wVBp5Nzkwnh8ILzgBXu/k5UGZPe09zczOdHjWL3\nT3yCuQ88wOzGRuY+8AC7nX8+0444ghvWrmVsayvJr5UG7Hj7bS7fvp3+Gd5vNTA2x/lWA6dk2VYb\nvn9Xt7WR/YOht7eNBhpzbPtTN48r9rbk9kdikp9q3pbcnu13EZXYBTozm2JmU4BjCT4PzgzTTgq3\nDzKzNjO7OnmMuz8D/By42cwuNLNx4evBBF+ypIK0t7fzwN13c9VZZzFz7FiuOuss7l+0iFlnn13U\nYJYrIOXbXqrg0pvbJgCLge1Ztq3Pse3nvbgNguv4Rv/+schPNW+D3L+LqMQu0AF3E4yHu4jg8+C7\n4etZ4XZLeaT6FLAQuA74JXAQMMHdm0qe4wrS0NAQdRZyylZrW//xj/OhdeuKGsxSA1JDnu3pShVc\nevNDqwa4HDi3f38eSSR2XqsDjYkE+7/3vVw2fDgrM2x7e/jwXtu2MpHgv0eMYN7atcwIm6ejys/J\nEVx/nLbl+l1EKbYzo/QGM/Nqvv44a29vZ8XSpTz+4x9T29JCW10dIz/5Se6/6SbmPfXULgHtKmAu\nmYPWTGB2lvPkOm45kCB7kFwO9AVOzbCtA5gO3AC75LUDuBS4KcM2gJeBz/Xvz5fa2xkX1k4dWJVI\ncMfQofQFPvLccztrrqXctqS+nmvuvZenH3+c1QsX7vxdjJk2jfGTggmHVtxzTyy21dTU0NHREZv8\nVPO29N/F9fffH2lnFAW6Kr7+uGpubmb2xImc09REQ8oH77zaWurb25mQ4XdWimCWK1gB/AP4cP/+\n3J2lWbS7AStuwaWmJo4NP1JOzEyBLioKdNEqZq2tVMEsV7BaUl/Ppbfdxi0XXsiUpqaiBiwFF6kk\nCnQRUqCLTrFrbaUMZrmCVb7mMgUsEQW6SCnQRaOjo4PpI0dyw9q1Rau1lTqYiUj3KdBFSIGu9DI1\nT9a9732ceMstjGtNn4O7+7U2gGbgv8z4ZG0tp77zjoKZSEwo0EVIga60sjVPXmLGre5Fr7VtB644\n/njOuvxy1vzkJwpmIjGhQBchBbrSydU8Wcpa28xlyxgwYEBRrkFEiiPqQFcpkzpLhLI1T05paspY\n80oOtM70Vz+BoNZ2PJlrbf2BPY87jh2XX87VabW2+aq1iUgGqtFV8fUXQ3eaJ1VrE6kuqtFJ2ero\n6GD2xIm7NE8acECWIAeqtYlI71Kgk25bsXQp53SjebKG4D7dFcDEmhomdHTsUmubFdbazjr33Azv\nICJSODVdVvH1d0Wm+3DNmzbxg6ef7lbzJARrwa2dPp2W559XD0mRChZ106UCXRVff6Gy3Yf7ghnf\nzlJ+hQwFmDFiBPPXrFFQE6lwUQc6NV1KTrnuw73HvUfNkzOXLVOQE5GSU42uiq+/EMsXLyYxdSpj\nM8xiouZJESmEanQSa6sXLmRuhiAH+XtPbgfu+eAHmf/1ryuoiUhkFOhkp0wdTra+8krWYQLJ5skZ\nwNlqnhSRmFLTZRVff6qsA7+BW8l8Hy6pHfjUMccw6MAD1TwpIruIuulSga6Krz8p17yUy4F+wLgc\nx69MJNhx112cPnly6TIpImUr6kCnr9uSc+D3BGAJwf22TLYDS+rrd66YLSISNwp0wuqFC2nI0uEk\ndZjAr2pqSNZ/naAmN2PECN2HE5FYU2cUobalJec9uAHAd4ALjjqKxwYP1tyTIlJWFOiEtj59sg78\nTjJg4ODBfOWXv+ylXImIFIcCXZXZZQhBezt1Tz/NKnJ3OFmVSDBm2rTeyqaISNGo12UVXX+2IQQr\ngZuAu9G8lCJSfFH3ulSgq5LrzzWEAOBl4HOJBF8CxqUEQS12KiI9pUAXoWoKdLnmrEzSvJQiUgpR\nBzrdo6sSueasTDrl7bdZ9fzz6nAiIhVFX9OrRL4hBBD0rKxtaemN7IiI9JrYBTozG2hmi83sTTPb\namZLzOzgAo892Mx+YmavmFmLmf3BzK4zs7pS5zvu2lpayNdI60BbXdUXlYhUmFg1XZrZ7sAq4C1g\napj8FWClmR3t7m/lOLYOeAToA1wF/Bk4DpgDHAF8tIRZj5VdhhDs2EHdb36jIQQiUpViFeiAi4DB\nwBB3fxnAzH4HbAAuBm7Ocewo4HBggrs/HKb92sz2A75kZgl3z32TqgKkDiGYmzaE4EZgBNmHECyp\nr2e+5qwUkQoTt6bLs4Enk0EOwN3/BDwOfCjPsbuFP7empW8luM7Ievz0lo6ODmZPnMgNa9cyNgxy\nEFz4KcB3gXN3351HEgnNWSkiVSNuNbphwL0Z0tcD5+Q59mGCmt8NZnYp8CpBBebzwK25mj0rRa5V\nCAAOBf6ro4Mnp09nZdoQAs1ZKSKVKm6Bbl9gS4b0zcA+uQ5097fNbAzBqjLrk8nAj9x9elFzGVMa\nQiAisqu4BbpuM7N+wCKCyfY/TtAZ5Xhgppm1u/ulUeavN2gIgYjIruIW6LaQueaWraaX6tPAScAR\nKff4VpvZNuD7Znaru/8u/aBZs2btfN7Q0EBDQ0M3sh0PbXV1eVch0BACESm1xsZGGhsbo87GTrGa\nAszMHgH6uvtJaemrANx9bI5jbwXOc/f90tKPBp4BPuLui9K2VdQUYMu/9S0SX/wiWQuJoOPJjrvu\n4vTJk3stXyJS3aKeAixuvQ+WASeY2eBkQvh8FHBfnmNfB/Y2s8PS0k8gqMhsLFYmY6mlhQk/+hGL\nCYYKZJIcQjBeQwhEpIrErUZXR1D7egu4JkyeQzD0q97dW8L9BgEvAbPcfW6YdgjQRBDwrifodXkc\ncDXwe3cfkeF8ZVmj22VAeF0do996iwmrVvEGMNuMKX37MnbHDq1CICKRi7pGF6tAB8EUYMA3gdMI\nbjc9DFzm7q+m7HMInYHuupT09wGzgBOB/Qk6pNwHXO/u6ePryjLQZVtTrhFYDMwE9v/e91ix//6s\nXrhQqxCISOQU6CJUboEu35py24EZ++3H/L/+lZo+fXo7eyIiGUUd6PT1vozkGxDeH5iyfTsP3pfv\ndqaISPVQoCsjqxcupCHPgPCxra08tmBBL+VIRCT+FOjKiAaEi4h0nQJdGUkOCM9FA8JFRN5Nga6M\njL7gAhoTiZz7aE05EZF3U6ArIxMmT2bxkCEaEC4i0gUKdGWkxoyZ/fszg2AhVa0pJyKSn8bRldP1\n3347fPKTdAArzFg9ejS1tbUaEC4isRb1ODoFunK5/s2b4cgj4Y03gteXXw7z5kWbJxGRAijQRSjO\ngW6X+SxffJHRf/4zE4Cagw+G556DPfaIOpsiInlFHejith6d8O75LOemzWc5HZh53XUMUJATESmI\nanQxu/6C5rMcMYL5a9bofpyIlIWoa3T6pIyZguazbGriwXvv7c1siYiULQW6mNF8liIixaVAFzOa\nz1JEpLgU6GJG81mKiBSXAl3MaD5LEZHiUqCLmQmTJ7P4qKM0n6WISJEo0MVMTU0NM487TvNZiogU\nicbRxe36N22Cww6jo7WVFcDqY4+l9j3v0XyWIlK2oh5Hp5lR4uZrX4PWVmqAM449ljOeegossr8P\nEZGypxpdnK5/40Y4/HB4++3g9S9/CWedFW2eRER6KOoandrA4uRrX+sMcscfD2eeGW1+REQqgJou\nI/SuFQq2bKHtyScZDcEKBbNmqclSRKQI1HQZ0fWnrlDQkLZCweL+/Zn54osMOOCASPImIlJMUTdd\nKtBFcP1aoUBEqknUgU6fohHQCgUiIr1HgS4CWqFARKT3KNBFQCsUiIj0HgW6CGiFAhGR3hO7QGdm\nA81ssZm9aWZbzWyJmR3cheOPMrNFZvY3M2sxs9+b2fRS5rmrtEKBiEjviVWgM7PdgVXAEGAqcD7w\nXmBluC3f8ccCTwK7ARcCZwA3An1KlefumDB5Movr67VCgYhIL4jV8AIz+wJBYBri7i+HaYOBDcAV\n7n5zjmMNeBZ43t3PKfB80Y2jW7uW2SecwBRgLOwcR7cqkWBJfT0zly1jwIABkeRNRKSYoh5eELdA\n9zDQz93HpKU3Au7uY3McOw54CBjj7msKPF90c11ecQUdN94YrFAwYAC1w4ZphQIRqUhRB7q4TQE2\nDMg0eGw9kK+WNir8WWdmTwD/DmwBfgZc6e65+/P3ptZWWLgwWKEAOGPBAk3eLCJSInGrNuxLEJzS\nbQb2yXPsgQQtgD8DfgWcCnwd+DRwVxHz2HN33w1//3vw/JBD4PTTo82PiEgFi1uNridqCG5z3eHu\ns8O0R82sFviqmR3p7n9IP2jWrFk7nzc0NNDQ0FD6nN56a+fziy+GPrHqKyMi0iONjY00NjZGnY2d\n4naP7nXgHne/JC39u8A57p51lmMzux64Epjo7venpH8Q+A3wUXf/edoxvX+PrqkJPvjB4HnfvvCX\nv4A6nYhIBYv6Hl3cmi7XE9ynSzcUeK6AY+MvtTY3ZYqCnIhIicUt0C0DTgiHFAA7hxeMAu7Lc+xy\nYAfBcm6pziBo0lxXpDx237ZtcOedna8vuST7viIiUhRxa7qsA54B3gKuCZPnEEzoX+/uLeF+g4CX\ngFnuPjfl+GuBq4F5wErgOOBa4KfufmGG85W86fJdi6tu2EDbhg3B4qpDh1Lz7LNaXFVEKl7UTZex\n6ozi7i3heLhvArcT9KJ8GLgsGeRClvJIPX6OmW0DLgW+BGwi6Hk5lwikLq46N21x1enbtzPzb3/T\noHARkRKLVY2ut5WyRqfFVUVEAlHX6PQJWyJaXFVEJB4U6EpEi6uKiMSDAl2JaHFVEZF4KHmgM7O+\npT5HHGlxVRGReOiNGt1VvXCO2NHiqiIi8dDjXpdmthdwEMGkygdleH6ku+/Vw3yWhHpdioiUXtS9\nLosR6J4Gjib47H4N+AuwN8HacM3AFHcf3cN8lkSpB4w3Nzcz+4gjmPKPf2hxVRGpWpUQ6PoCnwHW\nufu6MO1Kd/96+Pxcd7+7xzktgZLPjPL663QceCAr3FkN1I4cSds++2hxVRGpKlEHuh7PjOLu7wC3\nmNnRZjbR3ZdBZz+MuAa5XrF0KTXuweKqJ58MMVq2QkSkWhRtCjB3/62ZPW9mpxA0XcqiRZ3Pzzsv\nunyIiFSxkkwBZmZHAu8DVgDj3P2Bop+kCEradPn663DggeAONTWwcSP867+W5lwiIjEWddNlSW4S\nhSt5309w7+6OUpwj9pYsCYIcwEknKciJiESkZL0h3L3N3ecDj5fqHLGmZksRkVjI23RpZsvd/Yxu\nn8DsBHd/srvHl1LJmi43bYKDDupstnztNTjggOKfR0SkDJRD0+UxPTlBXINcSaU2W558soKciEiE\nCgl0A8zsoq68qZld3M38VAY1W4qIxEYhTZftwB8JVvn+VUFvavaUux9fhPyVVEmaLl97DQYO7Gy2\n3LQJNPuJiFSxcmi6/KG7HwWcambn59rRzA43sxuB4UXJXTlKbbZsaFCQExGJWJfG0ZnZfwPvuPtN\nKWk1wCTgs8A4guDp7t6nyHktumLV6Nrb21mxdCmP//jH1K5eTdu2bYwGJtxyCzWXXNLzjIqIlLGo\na3RdHjBuZtMIBoPPBz4NXAj8G8Gcxc8APwHmxHXFglTFCHTNzc3MnjiRc5qaaGht3TlxcyOwePhw\nZi5frombRaSqxT7QmdlUd78j5fVZwE3Ae8OkVuDnwPfc/alwnx+7+6dKkuMi6mmg01I8IiL5lUOg\nexU4Cfg4QQ1uEEHt7Q9ABzDe3TeWOJ8l0dNAt3zxYhJTpzK2tTXrPisTCXbcdRenT57c7fOIiJSz\nqANdIdWMgcCLwHUETZSLCOavPAqYAnzTzPYtXRbja/XChTTkCHIAY1tbeWzBgl7KkYiIpCu0Pe0l\n4EpgoLt/1N0bAdz998BlwA/N7JDSZDG+altayPcVxcL9REQkGoUEuo3Ake4+z93fSN8YNlt+GrjB\nzD4AYGajipvNeGqrqyNfw6eH+4mISDQKCXSr3b0j1w7uvgX4JHC5mX0RWFKMzMXd6AsuoDGRyLnP\nqkSCMdOm9VKOREQkXVHXozOzPsADwKnVMI5OvS5FRPKLujNK0RdeNbODgRfdfbeivnEJFG0c3fDh\nTNm4kbGwcxzdqkSCJfX1zFy2TOPoRKSqVVygAzCzde5+XNHfuMiKNTNKx6hRrFizhtVA7dChtB16\nKGOmTWP8pEmqyYlI1avUQHeWu9/fzWMHAjcDpxJUkB4Gvujuf+7i+/w3cD3BPcaTsuzT80C3bRvs\nuy+0t4MZNDfD/vv37D1FRCpI1IGuJNWNHgS53YFVwBBgKnA+wQwsK8Nthb7PYcBVwF+7k48uWbUq\nCHIAxxyjICciEjO1UWcgzUXAYGCIu78MYGa/AzYAFxPU9ApxC3AnwZycpe0U8+CDnc/Hjy/pqURE\npOvidgPpbODJZJADcPc/AY8DHyrkDczsYwSrov9PKTK4CwU6EZFYi1ugGwY8myF9PTA038Fmtjfw\nDeAKd3+zyHnb1UsvwQsvBM/r6mDkyJKfUkREuiZugW5fYEuG9M3APgUcfyPwB3e/vai5yuahhzqf\nNzRAv369cloRESlc3O7RdZuZjSHovHJMr51UzZYiIrEXt0C3hcw1t2w1vVTfA24DXjOzvQjnUwZq\nwtdvufuO9INmzZq183lDQwMNDQ2F5bStDR55pPO1Ap2ICACNjY00NjZGnY2dSjKOrrvM7BGgb/q4\nNzNbBeDuY3Mc20EwKUmmsRoOXObu3047pvvj6J54ovOe3MCB8OqrwTg6ERF5l6jH0cWtRrcMmGdm\ng8PelpjZYGAUMCPPsQ0Z0r5FcB/yPwnW1Cue1Ptz48cryImIxFTcAt0Pgc8B95nZNWHaHOAV4AfJ\nncxsEMEaebPcfS6Auz+a/mZm9ibQx90fK3pOdX9ORKQsxKrXpbu3AOOAPwK3A3cQ1MROCbclWcoj\n79sWO59s3QpPPhnmxOCUU4p+ChERKY641ehw978A5+bZ5xUKmPEk1z29Hkmd9mv4cE37JSISY7Gq\n0ZUNNVuKiJQNBbruUKATESkbsWu6jKv29nZWLF3K47fcQu2LL9IGjO7XjwkjRujbgohIjMVqHF1v\nK3QcXXNzM7MnTuScpiYaWlt3riLeWFPD4uOO0yriIiI5RD2OToEuz/V3dHQwfeRIbli7lv4Ztm8H\nZowYwfw1a7SauIhIBlEHOn0y57Fi6VLOaWrKGOQA+gNTmpp48N57ezNbIiJSIAW6PFYvXEhDa2vO\nfca2tvLYggW9lCMREekKBbo8alta8o5Kt3A/ERGJHwW6PNrq6vJOreLhfiIiEj8KdHmMvuACGhOJ\nnPusSiQYM21aL+VIRES6QoEujwmTJ7O4vp7tWbZvB5bU1zN+0qTezJaIiBRIgS6PmpoaZi5bxozD\nD2clnTNEO7AykWDGiBHMXLZMQwtERGJK4+gKvP6OSy9lxa23shqoHTyYtmHDGDNtGuMnTVKQExHJ\nIepxdAp0hV7/8OHw9NPB8wcfhNNOK13GREQqiAJdhAoOdNu3w157BUvzmMGWLcFrERHJK+pApza3\nQqxb17n+3NChCnIiImVEga4QTzzR+fzEE6PLh4iIdJkCXSEU6EREypYCXT7uCnQiImVMgS6fF1+E\nN94Inu+9Nxx5ZLT5ERGRLlGgyye1NnfCCaAxcyIiZUWf2vmo2VJEpKwp0OWjQCciUtY0YDzX9f/z\nn8GYuY6OYKD4m2/Cnnv2XgZFRCqABozH2bp1QZADGDZMQU5EpAwp0OWiZksRkbKnQJeLAp2ISNlT\noMvGHZ58svO1Ap2ISFlSoMvmhRc6B4rvsw8MGRJtfkREpFtiF+jMbKCZLTazN81sq5ktMbODCzju\nWDP7kZn90cy2m9krZnanmQ3uVkY0UFxEpCLE6tPbzHYHVgFDgKnA+cB7gZXhtlw+DAwFbgbOAK4E\nhgP/a2YHdTkzuj8nIlIRaqPOQJqLgMHAEHd/GcDMfgdsAC4mCGLZfN3d30hNMLM1wMvAZ4BZXcqJ\nAp2ISEWI1YBxM3sY6OfuY9LSGwF397HdeM/XgV+4+2cybMs8YPwf/wgmcNZAcRGRHtOA8XcbBjyb\nIX09QbNkl5jZUcAA4LkuHZg6UPz971eQExEpY3FrutwX2JIhfTOwT1feyMz6AN8DmoEFhRzT3t7O\niqVLefzqq6kF2oDRBxzAhI4OatQZRUSkLFXyp/d3gROAj7v71nw7Nzc38/lRo9j9E59g7h//yGxg\nLpD49a+ZPnIkzc3Npc6viIiUQNxqdFvIXHPLVtPLyMy+Bnwa+IS7P5Jr31mzZuHuPHDbbczZuJHU\nm4AGjH3nHY5fu5YZEycyf80a1exERPJobGyksbEx6mzsFLfOKI8Afd39pLT0VQCFdEYxs6uAOcB/\nuvutefZ1d2f54sUkpk5lbGtr1n1XJhLsuOsuTp88uZBLERGRkDqjvNsy4ITUQd7h81HAffkONrPP\nA9cBX84X5FKtXriQhhxBDmBsayuPLSjoVp+IiMRI3ALdD4E/AfeZ2UQzmwjcC7wC/CC5k5kNMrM2\nM7s6Je0jwDeB5UCjmY1IeRyV66S1LS3k+6ph4X4iIlJeYnWPzt1bzGwcQcC6nSC+PAxc5u6pUcZS\nHkkTwp+nh49UvwbGZTtvW10dnvZmu+Qt3E9ERMpLrO7R9TbdoxMRKT3do4uBCZMns7i+nu1Ztm8H\nltTXM37SpN7MloiIFIECHVBTU8PMZcuYceSRrCRopiT8uTKRYMaIEcxctkxDC0REypCaLlOuv+P6\n61lx1VWsBmoPPJC2Y45hzLRpjJ80SUFORKSbom66jFVnlKjVPPssZxCs8cO118LFF0ecIxER6SlV\nU1I1NXU+r6+PLh8iIlI0arpMXn9rK+yxB7S3B0vzbNsWvBYRkR6JuulSNbqk554LghzA4YcryImI\nVAgFuqTf/rbz+dFHR5cPEREpKgW6JN2fExGpSAp0SarRiYhUJAU6AHfV6EREKpQCHcCmTfD3vwfP\n3/MeOOSQaPMjIiJFo0AH7262/MAHQLOgiIhUDH2ig5otRUQqmAIdqCOKiEgFU6AD1ehERCqYpgBr\nbYX+/TtnRdm2LeiQIiIiRaEpwKKWPvWXgpyISEVRoNP9ORGRiqZAp/tzIiIVTYFONToRkYqmQKca\nnYhIRVOvy+SLPfaArVs1K4qISJGp12VcHH20gpyISAXSJ3uS7s+JiFQkBbokBToRkYqkQJekjigi\nIhVJnVGSLzT1l4hISagzShwcdpiCnIhIhYpdoDOzgWa22MzeNLOtZrbEzA4u8Nh+ZjbPzF4zsxYz\nW2NmY/IeqPtzIiIVK1aBzsx2B1YBQ4CpwPnAe4GV4bZ8FgAXAlcDZwGbgBVmljuS6f6ciEjFqo06\nA2kuAgYDQ9z9ZQAz+x2wAbgYuDnbgWZWD3wU+JS73x6mPQqsB+YAkzIddxUw+p//ZEJHBzUaRyci\nUnHi9sl+NvBkMsgBuPufgMeBD+U5diKwA1iUcmw78DNggpn1zXTQXCDxne8wfeRImpube5b7CtDY\n2Bh1FmJJ5bIrlcmuVCbxFLdANwx4NkP6emBonmOHAi+7e2uGY3cDjsh0kAFj336bG9auZfbEiXR0\ndHQxy5VF/6iZqVx2pTLZlcoknuIW6PYFtmRI3wzs04Njk9uz6g9MaWriwXvvzZdHEREpI3ELdJEa\n29rKYwsWRJ0NEREpolgNGDez14F73P2StPTvAue4+wE5jv0ZUO/uR6Wln0twn+797v582rb4XLyI\nSAWLcsB43Hpdrie4T5duKPBcAcdOMrNE2n26YQSdVF5IPyDKghcRkd4Rt6bLZcAJZjY4mRA+HwXc\nl+fYXxB0Ojk35dg+wHnACnd/p7hZFRGRchC3pss64BngLeCaMHkOQV+RendvCfcbBLwEzHL3uSnH\n/xQYD8wAXgYuBc4ETnT3lKXERUSkWsSqRhcGsnHAH4HbgTuAF4FTkkEuZCmPVJ8CFgLXAb8EDgIm\nKMiJiFQxd6+qBzAQWAy8CWwFlgAHR52vHl7TOcA9wKtAC/B74Hpgj7T99gZ+BPwN+CfwEEEnnfT3\n6wfMA14L328NMCbDfgb8D0Ht+S2C2vjkqMsjRzn9CugA5lR7uRC0dPwa+Ef4f/AU0FCtZUJwe2QF\n8FdgG/D/gAuq5e+EoFIwP8zr9vD/ZFCG/SIrA+AzwPNAK8Fn3MUFX1/Uf2C9/MvcnWA6sd8SzMJy\ndvh8A7B71PnrwXU9AdwNfAw4Cfg8wZjCNWn7rSYIhucRNPE2hn+wB6btdxfB+MNpwFiCLwMtwNFp\n+30l/OMmLI1PAAAHhUlEQVS8DDgZuBVoB06PukwylNFHw3+6dnYNdFVVLgTT6e0AbgROAU4DrgDO\nrMYyAT4Q5vmR8DPhlDB/HakfppVcJuH5NxG0hC0P85Ep0EVSBgRBrp3gVtbJ4c92Cgx2kf/T9fIv\n8wvAO8ChKWmDw7QvRp2/HlzXfhnSpoZ/CA3h6w+Fr09K2WdP4O/AzSlp9eE/+CdS0voQfIO6NyXt\nXwi+WV2bdt6HgWeiLpO0PO0T/hN/mLQaXbWVC3BI+IEzPcc+1VYm14f52z0tfQ3weLWVCcHE+LsE\nuqjKIDz2r8CCtP1uA5qBPvmuKVb36HpBT+bSjC13/3uG5HUEzQIHha/PBl5z90dTjttG0Fs19doL\nnTP0dKAvwTe3VHcCHzCzQ7p9QcX3deC37v7zDNuqrVySH2Lfz7FPtZVJX2CHu7+Vlr6Vzn4ME6mu\nMskkqr+LE4H9M+x3B7AfMDpfxqst0PVkLs1y0wA4neMPc137oLDHKxQ+Z+hQ4G13fzHDfkZMytPM\nRhMs9/S5LLtUW7mMIvhm/VEze8HM3jGzDWZ2aco+1VYmPwbMzL5tZv9mZnuZ2WcIOsZ9I9xnKNVV\nJplE9XeRHFudfu6Cy6raAl1P5tIsG2Z2EDAbeMjdnw6T880Fuk+B++2b8vPNAvaLTPjN8XvAPHff\nZcKAULWVy4EE6z3eQNBkdxrwIPAdM5se7lNVZeLu6wnuI/0HsJHgmuYDn3X3u8PdqqpMsoiqDJI/\n09+z4LKK28wo0kNm1p9gcP0OghvB1exKIEHwgS6BGmAPgvsnyUkYGs3sUILeb/Mjy1lEzOwIgs4S\nvyNYE7OVoCnu+2bW6u4/jTJ/0nPVFui2kLnmlu0bSFkxswRBr6nBBDeMX0vZnOvak9uTPwfl2G9z\nyn57F7BfJMzsYODLBPekEmHZJMdd9jOzvQi61ldVuRB0HDiC4IZ/qgcJ7p8cQPWVyVcJvhhOdPe2\nMG2Vme0PfAv4KdVXJplEVQbJ992HoFNKtv2yqramy57MpRlrZlZL8K10OHCGu6dfT65rf9U7B+Sv\nBw4NA0Oq9DlD1xMEjMMy7Jd6bzAqhxGM5bmT4B9lC8E/hBN0pd8MvJ/qK5f1Be5TTWXyfoLOSm1p\n6U8B+5nZAKqvTDKJqgyS9+LSz528N5e/rKLsxtrbD4LhBTuAwSlpg8O0ch5eYAQ9nLaTMug3bZ9k\n1+AxKWl7Am/w7q7BHyToGjw1Ja1P+MeU3jX4beCatPM8DDTFoEz2JBhTmP7oAH4CjAHqqrBczgyv\nd3Ja+grglSr9W1lFMJa2Ni39/4b/U7XVVCbkH17Qq2UQln8zcFvafsmB67V5rynqP7Je/gXWEUwv\n1kTQBXYiwUj8DUBd1PnrwXUlB7fOAUakPQ4K9zGCYRSvEIwnm0Aw2PON5D4p7/dTgiauCwl6ni0m\nGHtVn7bfV8P01MGebQQ1ysjLJUtZpY+jq7pyIRgY/TeCgeOnAT8MP8CmVmOZAFPC6/9V+JlwGvCd\nMG1etZRJWA5T6Pw8+Wz4+qSoyyD8W20jmN4xOWC8jaDDUP5ri/qPLIJf5kCCWURSpwDbZQaAcnoQ\nTJ/TnuVxbcp+yel73iCYvudBsk/fcyOd0/c8Qfbpe77Mu6fv+Y+oyyNPWbUDs9PSqqpcCDqjzCcY\nRN8a5u/DVV4mE4CVBPeAtgK/CT9crVrKhCC4ZfoMWRmHMiCYHeX34X5/oAtTgMVq9QIREZFiq7bO\nKCIiUmUU6EREpKIp0ImISEVToBMRkYqmQCciIhVNgU5ERCqaAp2IiFQ0BToREaloCnQiIlLRFOhE\nYs7MbjKzp8xsTcoqzvmO+YKZPWFmz5vZgaXOo0icKdCJxN++wDnuPtI7l0LJyd2/5e4nEswxWG3r\nToq8iwKdSEyY2WFmtjFcNLZob1vE9xIpSwp0IvExkV1XURaRHlKgE4mP0cBad98RdUZEKokCnUh8\njAYejToTIpVGN6lFImRm5wEXEDRZDgDGmtlxwC/c/dY8xx4MXANsJ1hAtRW4wd3fKm2uRcqLAp1I\nhNx9EbDIzC4C6oFTC2m6NLME8DDwf9x9g5kdBqwFngaWlTLPIuVGTZci8TAWWNeF+3OnAf8C/Dl8\n/SbwTeChEuRNpKwp0InEQwPw6y7s/yawN/CsmX0bGOru16vZUmRXCnQiETOzocABdCHQuftjwJXA\nbsDngEfN7JLS5FCkvCnQiURvHPAOsAbAzPYys4H5DnL3ee4+CHg/wf25z5U0lyJlSoFOJHqjgWdS\npvf6AkHgy8jMlprZM8nX7v48cCfwUklzKVKmFOhEotcH+BOAmR0LtLh7rtlR/p2gxyXhMQOAjwFz\nSphHkbKl4QUi0bsO+K6Z3QD81d1vzLP/+cDJZvYVYE+gDviCu/9vifMpUpYU6EQi5u6/BcZ0Yf/H\ngMdKlyORyqKmSxERqWgKdCIiUtEU6ETKg9aVE+kmBTqR+NsM3GNmT5lZXSEHmNllZrYOGAm0lTR3\nIjFn7h51HkREREpGNToREaloCnQiIlLRFOhERKSiKdCJiEhFU6ATEZGKpkAnIiIVTYFOREQq2v8H\n6pGTJQodX/gAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "%reset\n", "\n", "%matplotlib inline\n", "import matplotlib.pyplot as plt # plotting modules\n", "plt.style.use('presentation') # just have in your script for prettier plotting\n", "import numpy as np # our matlab-like module\n", "from scipy.integrate import odeint # integration of ODEs so you don't have to write your finite difference yourself\n", "\n", "\n", "plt.clf()\n", "k=0.001\n", "y_A0 = 10.\n", "time_start = 0 #s\n", "time_finish = 10000 #s\n", "N_points = 50\n", "\n", "times = np.linspace(time_start,time_finish,N_points)\n", "\n", "y_A_exact = y_A0*np.exp(-k*times)\n", "\n", "x_A = ((y_A0 - y_A_exact)/y_A0)\n", "# print(y_A_exact)\n", "# print(x_A)\n", "plt.plot(times, x_A, 'ro-')\n", "\n", "plt.xlabel('$t [s]$')\n", "plt.ylabel('$X_A$')\n", "plt.legend()\n", "# plt.savefig('conversion factor.pdf')\n", "plt.show()\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Problem 2:\n", "\n", "*2A ->B*\n", "\n", "We assume it is a first order reaction and we use the design equation from the Batch reactor:\n", "\n", "So that means that $-r_A = k [C_A]^2$.\n", "\n", "Lets assume $k = 10^{-3} 1/s$ and $C_A(t=0) = 10$ moles\n", "\n", "So we can expect that the half-time of A should be around 10000 seconds which is ~3 hours.\n", "\n", "Ok lets see that through Python.\n", "\n", "We want to track the concentration which we will call $y_A$.\n", "\n", "$dy_A/dt = - k C_A = -k y^2_A$\n", "\n", "The exact solution here would be\n", "\n", "$ 1/y_A = kt + C$\n", "$y_A(0) = 10$ => $C = 1/10$\n", "\n", "\n", "$1/y_A = kt + 1/y_A_0$\n", "\n", "So now it is easier to track $1/y_A$ rather than the function itself\n" ] }, { "cell_type": "code", "execution_count": 39, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Once deleted, variables cannot be recovered. Proceed (y/[n])? y\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "/Users/bazilevs/Downloads/yes/envs/mypython3/lib/python3.6/site-packages/matplotlib/figure.py:1999: UserWarning: This figure includes Axes that are not compatible with tight_layout, so results might be incorrect.\n", " warnings.warn(\"This figure includes Axes that are not compatible \"\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAboAAAEYCAYAAAAqIzNgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3XmYFNXVwOHfAUcBEQUhRkXFHVFB\nYWJQUMHdaFwiRIkL6CeoLFFjjBqX4G6MSqJGFDBiMC4RdwMGF3CLRAcVFUQRREGJIigMewPn++NU\nDzU11T3dMz3TPT3nfZ5+BqpvVd/GsU/fOveeK6qKc845V6ya5LsDzjnnXF3yQOecc66oeaBzzjlX\n1DzQOeecK2oe6JxzzhU1D3TOOeeKmgc655xzRc0DnXPOuaLmgc4551xR2yTfHShWbdu21Q4dOuS7\nG845V5SmTZv2naq2y6StB7o60qFDB8rKyvLdDeecK0oi8kWmbf3WpXPOuaLmgc4551xR80DnnHOu\nqBVUjk5EdgBGAEcCArwEXKSqX2Zw7k1AKdANaAOcrapjY9r1B04I2u4IPKiqA1JcsydwK7A/sBR4\nGLhSVVdl+94AEokECxYsYPXq1TU53VWjWbNmtG/fnpKSknx3xTmXypw5cPvt8NBDsHw5tGwJZ5wB\nl1wCu+5aJy9ZMIFORFoArwBrgP6AAjcAk0Wks6quqOYSw4D3geeBs9K0OwNoB7wI9E3Tn85Bm38D\nxwM7A38CtgdOzeAtVbFgwQK22GILOnTogIjU5BIuBVVl8eLFLFiwgJ133jnf3XHOxZk4Efr0gUTC\nHgDl5TBmDDz4IIwfD8cem/OXLZhABwwEdgH2VNXPAETkA2A2cB5wRzXnb6mqG0RkN9IHuqNVdUNw\n/WPStLsWWAD0VdVE0H4t8KCI/FFV383kTYWtXr3ag1wdERG23nprFi1alO+uOOfizJljQW7lyqrP\nJQNfnz7wwQc5H9kVUo7uBGBqMsgBqOrnwJvAidWdnAxeuWgnIiXAMcA/k0Eu8E9gbSb9SXPtzBrO\nmQODB0OrVtCkif0cPNiOu1j+BcK5Anb77RWjuOeB8rg2iQSMGJHzly6kQLc38FHM8RlAp3ruy65A\nMyL9UdXVwJw678/EidC5sw3ny8tBdePwvnNne9455wpV3Bf1MWOYn0hwMvBz4Oq48xIJGDcu590p\npEDXBvg+5vgSoHUe+gKp+9Mm5jgiMkhEykSkrMa30MLD+0Si8nOJhB3v06fGI7uDDjqo2jZ//vOf\nWRl3eyHH5s2bx8MPP5z1eQMGDGD8+PF10CPnXK2l+KL+j0SCTsDTQbO7gNiSGsuX57xLhRTowCag\nROXjflTyNbPqj6qOUtVSVS1t1y6jyjRVhYb3KdVieP+f//yn2ja5DHTr1q1L+VxNA51zrkCl+aK+\nPRAOYf+HTcqoomXLnHerkALd98SPlFoTP7KqS0uCn6n6syTmeHZE4h8jR2YW6P7619TXSKNl8Es0\nZcoUevXqRZ8+fejYsSOnn346qsqdd97J119/Te/evenduzcAkyZN4sADD6Rr16707duX5cE3rgkT\nJtCxY0d69uzJr3/9a44//ngAhg8fzqBBgzjqqKM466yzmDdvHgcffDBdu3ala9euFcH28ssv5/XX\nX2e//fZjxIgRrF+/nksvvZSf/OQndO7cmfvuuw+wGZVDhw6lU6dOHHfccXz77bc1/md3ztWhNF/U\newEDsLzP68AoYj5gS0rgzDNz3y9VLYgHtrTgjZjjU4BXs7jObthIbEAGbRcAY2OObwqsBm6IHG8W\nHL+2umt369ZNo2bOnLnxLzagr5tHGptvvrmqqk6ePFlbtWql8+fP1/Xr12v37t319ddfV1XVnXba\nSRctWqSqqosWLdKDDz5Yly9frqqqt9xyi1577bW6atUqbd++vc6dO1dVVU877TQ97rjjVFX1D3/4\ng3bt2lVXrlypqqorVqzQVatWqarqp59+qsl/m8mTJ1eco6p633336fXXX6+qqqtXr9Zu3brp3Llz\n9YknntAjjjhC161bp1999ZVuueWW+vjjj8e+v0r/xs65uvXZZ6oXXKC6xRaqIhWfQc8Gj+hn01LQ\nNek+u1q0sGtmACjTDONCIS0veBa4TUR2UdW5ACLSAegBXF6fHVHVtSLyAvBLERmuqsn7b32AzYK+\nNngHHHAA7du3B2C//fZj3rx59OzZs1KbqVOnMnPmTHr06AHA2rVrOfDAA5k1axa77LJLxZq1fv36\nMWrUqIrzTjjhBJo3bw7YQvmhQ4fy/vvv07RpUz799NPY/kyaNIkPPvigIv+2dOlSZs+ezWuvvUa/\nfv1o2rQp2223HYcddlhu/yGcc9mLWRP3FfBr4ElgG2AWsFXolFaprlVSYo/x4+tk0XghBbrRwFDg\nGRG5ChuVXQ/MB+5LNhKRnbCZj9ep6nWh44diC8F/HBwqFZHlAKo6PtSuExtnTTYHdhKRPsHfX1XV\n5CyS4cBbwD9F5K9AB2zB+HhVnVbrd6tx6T9sptKYMelvX5aUwKBBcPfdterCZpttVvHnpk2bxubT\nVJUjjzySRx55pNLx9957L+21N99884o/jxgxgm222Ybp06ezYcMGmjVrFnuOqnLXXXdx9NFHVzo+\nYcIEXzrgXCGJrIlbD9wDXMnGZQPfYB/gt8edX1ICzZtvrIxy5plw8cV1VhmlYHJ0apVPDgM+BcYB\n/wA+Bw5T1XAOU4CmVO37tcDj2GQegCHB3x+PtPtl6Hgb7NZx8u97h/rzPnA0sC3wL+Am4O9Y1Za6\nc8kl9kuQTkmJ/VLUkS222ILycvt17d69O2+++SaffWbLG1euXMmnn35Kx44dmTt3LvPmzQPgscce\nS3m9pUuXsu2229KkSRPGjRvH+vXrq7wOwNFHH83IkSNJBEH+008/ZcWKFRxyyCE8+uijrF+/noUL\nFzJ58uS6eNvOuUyFcnHvAwdiI7nw2rj/A34fd27yi/rSpbB+vf28++46C3JQWCM61GpanlJNm3nE\nzHxU1V4ZvsZwbLSWSdvXsP+G9WfXXW34Hi2TA3U+vE8aNGgQxx57LNtuuy2TJ09m7Nix9OvXjzVr\n1gBwww03sMcee3DPPfdwzDHH0LZtWw444ICU1xs8eDCnnHIKjz/+OL17964Y7XXu3JlNNtmELl26\nMGDAAC688ELmzZtH165dUVXatWvH008/zcknn8wrr7zCvvvuyx577MGhhx5aZ+/dORcRV5ty9WqW\nJxIMB/6MjeiSOmK34A5Jdb06/qIeRzTVLTRXK6WlpRrdePXjjz9mr732yuwCc+bYEoJx4+pteJ+t\n5cuX07JlS1SVIUOGsPvuu3NxPf8CR2X1b+ycSy+uNiVW2WQIEK62vxl26/J3wZ+rCH9Rz0E9SxGZ\npqqlmbQtqBGdC9l1VxvO1zIPV5dGjx7Ngw8+yNq1a9l///0577zz8t0l51yupKhNuQxbJrA4dOww\nYCSwR9x1mjTJ+xd1D3Suxi6++OK8j+Ccc3UkxZq4VtgEkwHA1li1/TOJySflaNJcLnigc845VzUX\nF6S1vga2izQ9C5tVeQ7QNtX18pCLS8UDnXPONXYxubgV2Ky9v2BVO8JVcgXLxcWqp0lz2SiY5QXO\nOefyIKY+5QRsrdVtQALbEDTlyt6Sksq7FAwaZHvK1cEGqjXlIzrnnGvMQrm4r4GLqLr4uC1WcPhH\n0XMLKA+Xjo/oXI2MHTuWoUOHVtvm66+/rvj7ueeey8yZM7N+rSlTplQUjHbO1VJ0r7iRI1mfSHAP\nsBeVg9zWwFisEHGVIAcFlYdLx0d0rs6MHTuWffbZh+22s1T2mDFj8twj5xq5mFzcdOzW5H8jTftj\nty5jJ5sUYB4uHR/R5dHw4cMRkYwegwYNqnL+oEGDKrUZPnx4Rq970kkn0a1bN/bee++KQswtW7bk\nyiuvpEuXLnTv3p1vvvkGgOeee46f/vSn7L///hxxxBEVx5PKy8vZeeedK8p2LVu2jA4dOvD4449T\nVlbG6aefzn777ceqVavo1asXyUX0L7zwAl27dqVLly4cfvjhALz99tscdNBB7L///hx00EF88skn\nNfp3dc7FiMnF3Q90o3KQ2x14GRvJVQlyBZyHS8cDXSP0t7/9jWnTplFWVsadd97J4sWLWbFiBd27\nd2f69OkccsghjB49GoCePXsydepU3nvvPU477TRuvfXWStfaYost6NWrF//6178AePTRRznllFPo\n27cvpaWl/OMf/+D999+v2MkAYNGiRQwcOJAnnniC6dOn8/jjdrOkY8eOvPbaa7z33ntcd911/P73\nsZXynHM1EbMuridWOBigBLgG+ABbAF5JSQkMGVJvtSlzzW9dNkJ33nknTz31FADz589n9uzZbLrp\nphV5sG7duvHiiy8CsGDBAk499VQWLlzI2rVrK7blCTv33HO59dZbOemkk3jggQcqgmQqU6dO5ZBD\nDqm4Vps2tv3i0qVL6d+/P7Nnz0ZEKkaJzrkspahPGQ10ewJXAJOx+pQdU12vgeTiUvERXR4NHz48\n441pw3u9JY0aNapSm0xuXU6ZMoWXXnqJt956i+nTp7P//vuzevVqSkpKKrbCCW/ZM2zYMIYOHcqH\nH37Ifffdx+rVq6tcs0ePHsybN49XX32V9evXs88++6Ttg6rGbrtz9dVX07t3bz766COee+652Ndy\nzlVj4kTo3Nm2+yovB1U2lJdzbyJRsbVL2JXYOrnYIFdSAi1aNJhcXCoe6BqZpUuX0rp1a1q0aMGs\nWbOYOnVqte233357AB588MGU7c466yz69evH2WefXXEsug1P0oEHHsirr77K559/DsCSJUuqvNbY\nsWOzel/OOWLzcB9iu1dfgC3ynhs5pYSY8l0NNBeXige6RuaYY45h3bp1dO7cmauvvpru3bunbT98\n+HD69u3LwQcfTNu2KYv9cPrpp/P999/Tr1+/imMDBgzg/PPPr5iMktSuXTtGjRrFL37xC7p06cKp\np54KwO9+9zuuuOIKevToUbFnnXMuC6E83ErstmRXIPl1djVwc7rzG3guLhXfpqeO1HqbngZm/Pjx\nPPPMM4wbNy6v/Sjmf2PnqkhRn/IFYDC2c3VSCXA5thlqs1TXa9HCRnANILj5Nj2uXg0bNoyJEycy\nYcKEfHfFucYjZk3c/4CLgUcjTQ/GJpuk/ArYwNbFZcsDnau1u+6KS3E75+pMZK+4DcBo4DJgaahZ\na+BPwNlE8lQlJdC8ecFu6pxrHujqWaoZh672/Da8azQia+IS2L5w4SB3BrZvXEOtT5lLPhmlHjVr\n1ozFixf7B3IdUFUWL15Ms2Ypsw/ONUzR2pStWtnSgVCg2wy4N/jzrsCLwDgadn3KXPIRXT1q3749\nCxYsYNGiRfnuSlFq1qwZ7du3z3c3nMudmDwc5eW8A5RSeVlAb2A88DOgefQ6UPR5uHQ80NWjkpKS\n2MoizjlXRSQPB7ar92+Ah7EJJ6dGTjkl7jpNmjSKPFw6Huicc64QhfJwG7ACzL8DfgievhA4Gtgq\n1fmNMBeXigc655wrBCnWxM3AttF5M9L8cCBtWYVGmItLxQOdc87lW0wubhVwA3ArsC7UdBdgJHBU\nqms14lxcKj7r0jnn8immPuWLwL7ATWwMcptgJb0+JBLkSkoqz8gskvqUuVRQgU5EdhCR8SKyVESW\niciTIrJjhufeJCKTRGSxiKiIDEjTdqCIzBKRNSLyiYicH9OmqYhcLCIficgKEVkoIk+JSOdavEXn\nnKsssiZuNBbI5oSaHAS8hwW+FuFzk3m4pUuLrj5lLhVMoBORFsAr2G4R/YEzsc1uJ4vI5hlcYhg2\nq/b5al5nIFYN5wngGOBx4B4RuSDS9HpsJ/mngZ9jud9dg/74HHbnXM1E18WNHFkp0J0MbB38eSvs\nw+p1IHbzK8/DZaSQcnQDsdvPe6rqZwAi8gEwG8vF3lHN+Vuq6gYR2Q04K66BiGwC3AiMU9Urg8OT\nRWQ74HoRGaOqyd+4AcBjqnpV6PwPgI+B47DfP+ecy1zcuriIttg37EnYh96P4xp5Hi4rBTOiA04A\npiaDHICqfo5NNjqxupNVdUMGr3Eg0A54KHJ8HPYlqmfo2KbAski75MzeQvp3c841BJFc3GrgGuxW\nUdQAbK1clSDnebgaKaQP7L2Bj2KOzwA65fA1iHmdGcHP8OvcA5whIieKSCsR2SU4tgB4LEf9cc41\nFqFc3MvYZJPrgTuB/1R3bpHuE1dfCunWZRvg+5jjS7Ai3Ll6DWJeZ0nkeVT1GhFZAzzJxi8EnwK9\nVHUJMURkEDAIYMcdM5pD45wrRtE1cS1bwurVLEok+A1VbymNwSacpOS5uFoppBEdQFy141yW+k9e\nq9qqysHklKuwpSy9gb5AOTApyOlVoaqjVLVUVUvbtWuXoy475xqUiROhc2crvFxeDqpoeTl/SyTo\nSOUgtyVWjHlMqmuVlNhmqJ6Lq5W0gU5E1ufgcU2Gffme0IgqpDXxI72aqDJyi/x9CYCItAFGALep\n6h9UdYqqjsdm/bYDLs1Rf5xzxSRmTdwsoBfwf2z8AAI4LXjuPGI+iD0Xl1PV3boU4AtgXg2uLcAh\nWbSfwcYcWlgnYGYNXj/VaxC8zsLIaxB6nT2wnS/eCZ+sqktEZA5pNup1zjVioTxcArsddHPw56QO\nWLI/NnR5fco6kUmO7gFVva4mFxeRTGZCJj0L3CYiu6jq3OD8DkAP4PKavH6Mt4DvgNOBl0LHz8C+\nbCXLyf0v+HlA0C+C/rQBdgPezVF/nHMNWYr6lGAfrq+yMcg1BX6LzbRsUeVCAc/F1YlCmowyGhgK\nPCMiV2F5tOuB+YTWrInITljRgOvCAVhEDsVuKyZn5JaKyHKA4LYjqpoQkauxBeJfYcHuMOAcYJiq\nrg3azROR54FLg2D9Krb84HfYSG9k3fwTOOcajGrWxAn2wdUZ6Br6cyxfF1enqgt07YCV1bTJyfmq\nukJEDsNyY+Ow35OXgYtUdXmoqWBfjqK3ta8FDg39fUjwSJ6TfJ17RUSBS7Bc25fAUFW9J3K9U4M2\n/YKfy7CRXE9VLcvkPTnnilRkrzjFNj09CSgJNdsTu03UlcgHVkkJNG++cUZmI94rrj6IarUTEF0N\nlJaWalmZx0PnitLgwTarMpHgE+B8YAqWj6s2z+J5uJwQkWmqWppJ20JbXuCcc4UlWpuyVSsYM4Y1\niQTDsduRU4Km1wJzq7ue5+HqXVY5OhG5EdvU9qjwomkRKQG2UdUFOe6fc87lT1werrycKdiygE9D\nTZtileW3SXUtz8PlTbYjuiOAppEgtz/wFfCFiMwVkV457J9zzuVHzJq474CzsQoS4SB3ADAN2yS1\nylYrviYu77INdB2AaOLpRqzg9pvYrhLPi4h/XXHONWyhNXEK/B3bQ2xsqMkWwN1Yrcou0fO9PmXB\nyDbQtQK+Sf5FRLYCjgT+qaqHYF9sBJuG75xzDUeafeIewDbJXBxq3gerbDIEu21ZhefiCka26+i+\nwkZvSUdhwXIUgKp+JiLPYsHPOecahmrWxP0K+CN2u3JH4K/A8amu5bm4gpPtiO4D4GgRSX6BOR1b\nJ/daqM3nQGzRY+ecKzgxubhoqGuGFV++BKsjWCnIlZRUnpHpubiCk+2IbgRWJeRlEZmN/fcer6rr\nQm1+BKzOUf+cc65uhXJxi7EqEguBCVTeOqV38KjE18Q1CFmN6FT1dWzrmh5YMe5l2NKRsGjBZOec\nKxwxuThNJBiHTTZ5AHiBDHdX9jxcg5D1gnFVvQnYHtt5YjdVrdhZINiF+wDg7Vx10DnnciZmr7jZ\n2KSCs7DlA0mvpruO7xPXoNSoqLOqfgt8G/NUK+BB4OnadMo553IuUp9yLbbu7QZgTajZDthkk5/H\nXaNJE69N2QDVKNCJSEfgp9gt7C9V9RUAVX0f2wnAOecKSygX9zpW2eTj0NNNgIuwXEzL6Lmei2vQ\nsi0B1gS4HxvlgwW6DcnriIioV4l2zuVbdJ+4li1h9WrWJxKcD4yJNC/FttHpmup6notr0LLN0V2G\nrZt8CyvY/QSVJyYdJCLzRKTK5CTnnKsXMXk4ysshkaAplZcOtAT+AkwlRZDzXFxRyDbQnQ18Ahyq\nqqOBj8JPquqbwDrgl7npnnPOZSFmTVzUbdguyidjty5/TUxlE18TV1SyzdHtBPxVVdenaTMNOKjm\nXXLOuRoK5eHWAncC52JFeJPaAtOxqeNVeC6uKGU7olsGbFZNm6+AbWvWHeecy0KK+pRvAPtji7+v\niDktNsiB5+KKVLYjuneAI0SkiapuSNFmA7Bl7brlnHPViKlP+T02kWB0qNm9WHWLtFtRe33Kopbt\niO5+YHfgujRtOlO5yLdzzuVWJBenwCNYZZNwkNscq1u4X/R8r0/ZqGQ1olPVJ0TkMeAKEdmLSE1L\nETkZ25x1fO666JxzEaFc3BxgMDAp0uRE4C5sAXglnodrdGqyYPx0bMR2QfKAiEzGcrydsNm7f8xJ\n75xzLsWauLWJBLdjt5fC37i3xzZDPSnV9TwP1+hkHeiC3NxQERkHDMPKxB0aPP0BcKmqTstdF51z\njVbcPnHl5QA8A/w+1LQJ9oF0PbbzdxWeh2u00gY6ETkNeFdVP40+p6r/Bf4btNsMEFX17Xmcc7kR\nqU0Z1Qf7hv0qNsNyFCkmnHh9ykavuhHdw4CKyHLgfWyN3LvB4+NkuS9VXZP6Es45VwOhPJwCS7CF\n3kmCle2agI3kqnyYeS7OBaoLdL/Fvix1xfagOxj7nQNYJSLTqRz8ZlSzmNw55+JFc3FB2dy52GST\nr7APmZLQKXsGj1iei3OBtIFOVe9I/llEWmCzdLthga8r8BPgQDYGvzUi8iEwTVUH10mPnXPFJyYX\nlwBux3YTSOZEbgcur+5anotzERmvo1PVlar6H1W9S1XPVtUuWM63OzAU+BtWOm4/bAeMGhGRHURk\nvIgsFZFlIvKkiOyY4bk3icgkEVksIioiA9K0HSgis0RkjYh8IiLnp2jXXESGi8jsoO03IvK8iGxa\nw7fonAuLqU/5FvZN+go2BjnBbl9W4WviXDVqtB9dUpCbe5vQjuIiUgLsW5PrBaPGV7B9EPtjI8Ub\ngMki0llVV1RziWFYLvF5Nm4lFPc6A7Hb+zcDLwGHA/cE2wyNjLyXicDOQduZQDtspmmVOrDOuRoI\n5eJ+wILbvZEm+2H/wx4QPdfzcC4DtQp0cVQ1gd1Kr4mBwC7Anqr6GYCIfADMxkaJd6Q5F2BLVd0g\nIruRItCJyCbAjcA4Vb0yODxZRLYDrheRMcF7ALgE+2K5t6rOD13miRq8N+ccxObiFPgncCHwTahp\nC2yd3IWk+LDyPJzLQLYlwAAQkS1FpLeIHCEie+SwPycAU5NBDkBVPwfexAodpJWm/mbYgdio7KHI\n8XHYpK6eoWODgccjQc45V1Mxe8Up8AvgNCoHueOwWyiXkGJGpe8T5zKUdaATkcuBhdgtv38DH4vI\nQhG5UURi12lmYW8ie9wFZmBVV3Jh7+Bn9HVmBD87AQR5wR2AuSIyOsgXrhaRl0WkSuk851w1UuwV\nJ0CXULNtsRqCz2H7glXieThXA1kFOhE5C7gJWImNgP4MPIbtWHAF8J6I7FyL/rTBCpBHLQFa1+K6\n0dcg5nWWRJ7fLvh5GXY79TSgHzYanBI3QUZEBolImYiULVq0KEfdda5IhHJxUZcDewFDsBltp2AB\nsEJJCQwZAuvXw9KllpPzkZzLULY5uouw0VxnVa3YoUBEmmC7j/8ZeFFEumQwcSQVjTkmMcdqKnmt\nuNcJS34JWAn8XFVXAohIGfAZ9v/kZeETVHUUVqCB0tLS6q7vXPFKUZ9yaSLB1dgHyS6h5s2wBbnN\nU13Pc3GuFrINdB2Bv4WDHFTkxu4XkXlYEfFLSL+VTyrfs3FEFdaa+JFeTYRHbgtDx9tEnk++xzeT\nQQ5AVeeLyCxsIb1zLipmTZyWlzMe+DXwP+BTbDpz+BtsbJDzNXEuB7LN0a0ksjVPmKq+jOXtTqlh\nf2awMYcW1gnLS+dCMhcXfZ1kDjD5OnOBVaQeYWYy8cW5xiUmDzcPOB74JRbkwD4kpqS7jufiXA5l\nG+g+wtacpTMdW3dWE88C3UWk4q6GiHTAyo89W8NrRr0FfIdtNxR2BjaaexMqlkn8CzhYRDYP9WdH\nrOrQOznqj3PFI5SHSwC3Yd8oJ4SabAs8DvSKO99zca4OZBvoxgJdROSyNG22S/NcdUZjXwCfEZET\nReQEbDeO+dh6UQBEZCcRWSci14RPFpFDRaQPcExwqFRE+gTHgIoAdjXQX0RuEJFeInIdcA5wjaqu\nDV3yD9gmxf8SkZ+LSF/s/9kfsC2vnGvc5syBwYM3ViYZORISCf6L1Qe8FLsNBHYbZDA22aQPKRLv\nnotzdUFVs3pgdx3WYzvX7xd57jDsdt9L2V43dI0dsQXZy4By4GmgQ6RNB+yW4vDI8SnB8SqPmNc5\nD0sVrMEWpA9O0Z8DgMnY/69Lg/7sVt376NatmzpX1CZMUG3RQrWkRNVKMOtS0CGgEvn/rzPoW0Gb\n2EdJiV1rwoR8vyvXQABlmmFcEdXsJgeKSLMgyJ0Y/BIvw0ZhrbF1ZxuAw1T19awuXGRKS0u1rKws\n391wrm7MmWMLvyN7xU0Beof+3hwrynwRlXcdoKQEmjffOCPT94pzWRKRaaoauwVhVNYLxlV1taqe\njFUx+Rcb13u2x2peHt3Yg5xzRS/FmrheWJFasPzBDOz2ZZUgN2iQ5eA8F+fqQY1rXarq81jxZESk\nJbBGN9aIdM4VixRr4tYlEnyGrTkKuw34GdAXz8O5wlCjQBdsUbO5qn4PoKrLc9or51xhiFkTR3k5\nb2NJ7oXY5JJw2aK22FKCKnxNnMuTbEuA7SAir2ATM74TkR9EZIqI3CEiZ4hIJxHJZRUT51y+xKyJ\nW4Yt+u6O7Yf1DVb7Ly1fE+fyLNsR3T3Ybfj5wCdYzdWewCFsXFi9SkSmq2qPXHXSOZcHoTycAk9h\nGz5+HWrSHCvlpcTcpvS94lyByDbQHYwtlO6ZzMcFi6n3Dx5dgW7YEhrnXEMSs08cwJfAUGw3gbCj\nsW++u5CC5+Jcgcg20K0BpoS60vMaAAAgAElEQVQnnagVb34jeAAVOTznXEMRk4tbB9yFVVcIV2j/\nEfAX4FTSTDbxXJwrINkGupew8ldpaeXqIs65QhbOxYUcj1WHCBsE3EJkzyxfE+cKXLaB7kbgvyJy\ngKq+XRcdcs7VsxRr4vqxMdDtjdXgq5J49zycawCymnWpqjOx3/+nRKSfiDStm2455+pMivqUUWdh\n6+FuAt4lJsiB5+Fcg5DViE5EtgHOx27TPwTcKSJTgDJs38R3VXVJ6is45/IqJhc3H7gQ+C1wUKip\nYBUhPA/nGrpsb12OxL7k/YDVt9wJ23vuFILlBSLyBVZsM3bNqHMuTyK5uPXYFhxXAcuxyubvUrlc\nV5Ug16SJ5+Fcg5NtoDsM+BDokayGIiI7YUsKko+u1HzjVedcXQnl4t7FJpZMCz39ETbbLHY5t+fi\nXAOWbaBbD7wQLvmlql8AXwBPJo8Fm5M65/IlRX3K5YkEVwN3YtuMJHXCJpv0THU9z8W5BizbQPc6\nadaHJqnqlzXrjnOu1lLUp3wWW/g9P9R0M2yd3KVA7OJXz8W5IpDtNj3XAj8Tkb3rojPOuVqKqU/5\nFfALbAPJcJA7AstDXElMkPP6lK6IZDui64fdxn9RRPqr6ot10CfnXE3FrIkrxzaOTGoHjAB+hden\ndI1DtoHut2ys3/qCiHwOvMjG5QUfqeq63HbROZdSivqUYR2xHQauBc4F/gi0SXU9z8W5IlSTWZdd\nQ489sG2pBgXPrxWRD7HlBYNz1kvnXFUxubjlwNvY/6hhlwNHUXmdXCWei3NFTDTmG2DGJ4u0ALpQ\nOfh1AjZR1UZdNaW0tFTLysry3Q1XrObMgc6dK9WnfA4YAnyHLRVIO2vM61O6Bk5EpqlqaSZta7TD\neJKqrgTeCh7JF98U2Kc213XOVSOUi/sK2wz1ydDTg4GJpKlq4nk414hkO+uyWqq6VlXfzfV1nWu0\norUpW7WCMWNYn0hwN7AXlYNcW+D0dNfzPJxrZGo1onPO1bEUa+LexxLj70Sanw38Cdg67lqeh3ON\nVNoRnYjMFJEaTyqp7fnONWoxa+JWYFOfS6kc5PYEpgB/IybI+Zo418hVN6LriN0Jqananu9c4xWz\nJu4XwKTQ3zfFFnxfhlU5qcRzcc4Bmd267CUSm9LORM2ndDrX2GSwJu4KNga63sC92BqfWJ6Lcw7I\nMNAFj3ohIjtghRuOxCaNvQRclEn9TBG5Cbur0w1bE3u2qo5N0XYgcAmwM7bl0AhVvTfNtXfBZm03\nB3ZX1c8yf1fOVSMmF7cB+x8g/DWzF/AbbE3Pmfhecc5lorpA1zsHrzEv04bBurxXgDVAf2xEeAMw\nWUQ6q+qKai4xDHgf2y/yrDSvMxAr1n4zFkgPB+4REVHVkSlOuwdYigU653Insk8cwHSsEsPFwKmR\n5rdHz/c1cc6llTbQqeqr9dWRwEBsneueyRGTiHyA7Ql5HnBHNedvqaobRGQ3UgQ6EdkEuBEYp6pX\nBocni8h2wPUiMkZVE5FzfgXsjwXGETV7a86lEMrFrcBKdd2B7Yl1IXA0sFWqcz0P51y1cr6OrpZO\nAKaGbwuq6ufAm1jx9bRUdUN1bYADsbq2D0WOj8MmrFXakktEWmOfO7/FdlZ3rnai6+JGjoREgglY\npYU/YUEO4Hvslz8lz8M5V606C3Qi0kFEThaRa7M4bW8sDxY1AystlgvJLYairzMj+Bl9nVuBWao6\nLkev7xqziROtdNeYMVBeDqosBH4JHEfl+/y9gA+C41WUlECLFp6Hcy4DtV4wHtwK7ITd2tsv9GiF\n5coXAn/I8HJtsC+xUUuA1rXta+g1iHmdJZHnEZGe2C3Q/TO5sIgMIihwveOOvsm6i4jk4jZgieLL\ngWWhZm2A24ABxEw2adLE83DOZalWgU5EpmEjpE2Br7FR0bvAwdiX1EmqWp7lZeOWJNR4fUOaa6Vd\n+hDU7LwPm405M5MLq+ooYBRYUefadNIVoVAubi5WpmtqpEl/7NZlu+i5notzrsZqe+uyE/bls7Wq\ntlfVo1X1EiyIfFyDIPc98VtltSZ+pFcTVUZukb8nn78oOHaniGwlIlsBLYLnthCRLXLUH1eMUtSn\nTAa6VtgMq6TdgZeBscQEOfBcnHO1UNtA1xX4CTBeRDrnoD8z2JhDC+sEZDSqyvA1iHmdZG5uZujv\nP8aKw38fPP4aPPcu8HqO+uOKTUwejvLySlVO2mLfEEuAa7BcXHQPOcBzcc7lQK0Cnap+rKpHAyOB\np0TkARHZvhaXfBboHizOBmxSC9AjeC4X3sK27IoWeD8DG80lJ7ndgq0jDD/+GGp7bo7644pJTH3K\nhVgNyqj+wCxsOUGz6JNen9K5nMnJ7gWq+qSITMAqFH2ABdCSGlxqNDAUeEZErsJugV4PzMfyZQCI\nyE7AHOA6Vb0udPxQ7M7Pj4NDpSKyPOjj+OBnQkSuxhaIf4UtGD8MOAcYpqprg3azsM8hQtfvEPzx\nv14ZxcUK5eE2YAnby7FKA3tha1uShJjNUT0X51zO5WybHlVdDfxBRB7AFlW/LCK3AXdlUNEkeY0V\nInJYcP447LPgZawE2PJQUwGaUnVEei1waOjvQ4JH8pzk69wrIoqVALsU+BIYqqr3ZPRmnUtKUZ/y\nQ6zCwVuhpoOwsj1N013Pc3HO5ZxoTOHYnFxY5Cjgz0BbVf1RnbxIASstLdWysrJ8d8PVpZj6lCux\nWxC3AetCTXfD7u8fkepa4fqUfpvSuWqJyDRVLc2kbbU5OhFpLiKvisj9IlJlJ5BUVHUS0JmNeS3n\nikdMLu7fWGWTW9gY5EqwbXQ+IBLkSkoqz8j0XJxzdSaTW5cDsHVxL6jqmmwurqrriKlB61yDF8rF\nfYMVX34k0qQnlliuUtLH83DO1atMZl3+AlhMNQWVxTwmIn8P6kM6VxyqWRN3HpWD3FbYrKpXSVG3\nzvNwztWrTEZ0XbAKJ2lHc6qqIjIW2yJnElWLJjvX8MTk4SivXAfhFmAisBZbs3I7sE3ctXyfOOfy\nIpMR3VbAF5lcTFUnYgusj69Np5wrCDF5uFVAItKsIzbrKvntrkqQ8zycc3mVyYjuB2DLLK75BrBv\nzbrjXAEJ5eEAXgQuwDZNvCzS9IK48z0X51xByGRE9yU2ezJT84Fta9Yd5/IoxT5x32K3JI/CqhRc\nixVlrpbn4pwrCJkEuheBg0RknwyvWQK0rHmXnMuDmPqUG4Ax2K3Jh0NNNwM+SXctr0/pXEHJJNDd\nj6UlHhaRzTNovwewqFa9cq4+xeTiZmIldgZSeduMX2F14Spl2XxNnHMFrdpAF9R0vBVbCztVRPZK\n1VZE9sTu8Lydsx46V9dCubhVwNXYzsFvhJrsgi0I/weRySbJPNzSpbB+vf28+24fyTlXQDLavUBV\nrwEexLa2eV9ERovIoSLSDCrW0PUGnsZK+Y2uqw47V2spcnGzsWT0DWycWbkJVqn8I+wbXBWeh3Ou\n4GVc1FlVzxaRD7HPgf/Dqv2riCwDmmO7jAtwv6q+UBedda7W4tbFBXag8je/g7DKJrHJaV8T51yD\nkdV+dKp6B3YX50/YBslNsHV2m2Fr7S5S1YG57qRzORGTiwtrhgW2rYB7sZ11qwQ5z8M51+BkvU2P\nqv4PW0Z0WTA5pR2wQlV9AoorbKFc3MfYPlA3Etq/CeiFfWNrFT3X18Q512BlsnvBn0XkEBGR6HOq\nukJV53mQcwUnRX3K1YkE12B17W4GHos5tUqQA8/FOdeAZTKiGwoMAxaLyLPAk8BLyZ24nSs4KepT\nvoxVMJkdanopVrV801TX8lyccw1eJjm67bDPh2nAGcBzwHci8qiInCoiW9RlB53LSkwebhHQH9sP\nLhzkugP/IkWQ81ycc0Wj2hGdqn4LjAJGiUgrrGDzycDPgF8Ca0TkZeAp4Fm/jenyKpSHU2As8Ftg\nSajJltiOA4OI+abnuTjnik62sy6XqerDqtoXm4RyEvAo8FNs7dzXwW7kF4rITrnvrnMRKdbEzQJ6\nY2tgwkHuVGwiyvmk+OX3XJxzRSfrWZdJwf50zwLPikgTrGLSL4ATgRHAHSLyvqp2y0lPnYtKsybu\nNmzj06QOwD1ESneFeS7OuaKV1YguFVXdoKqTVXWYqu6IjfD+iC0kdy73qlkTdwuwNVam51KssonX\np3SucarxiC4dVX0HeAf4fV1c37lwLu47bMuM8KaJbYG/A+2J2WPK83DONSo5GdGFicifRGROrq/r\nGrEUa+I0kWAsto3O5TGn/YwUGyl6Hs65RqUuRnRtsZSIc7WXYk3cJ9iEkilBs3uBs4AD013L83DO\nNUo5H9E5lzMxebg12A7fndkY5AB2Cp6L5Xk45xq1akd0IvL3LK95UA374lxloTwc2CzK86i8u3dT\n4DfAH4AquwJ7Ls45R2a3Ls/A1t5WqXWZhtasO65RmzPHgttDD8Hy5aD2a7QYmzn5QKT5AVglgy6p\nrue5OOccmQW6cmABMDjDa15Oij0qMyEiO2Dr8I7EgutL2PY/X2Zw7k1AKdANaAOcrapjU7QdCFwC\n7AzMA0ao6r2h51sBFwHHAHtig4eZwK2q+nQN355LJcWauLnYWpXvQk23wAoyn4/9R6nCc3HOuZBM\nAt10oIuqvlptS0BEBtS0MyLSAngFS7f0x0aGNwCTRaSzqq6o5hLDgPeB57G5CaleZyC29djNWCA9\nHLhHRERVRwbNdsSC+wPA9cAGoB/wlIgMVdW/1uxduirCubiInbFt7ZO/fKcAd2IFWCuUlEDz5jYK\nbNkSzjzTRnIe5JxzZBbo3gd6iMiuqlrXywYGYhu77qmqnwGIyAdYLd7zgDuqOX9LVd0gIruRItCJ\nyCbYNmTjVPXK4PBkEdkOuF5ExqhqAvgc2EVVw5++/w5GnJcBHuhyJZKLCxPsG8nPgduDn5V4Hs45\nV41MZl2+CnyArb3NxNPAdTXszwnA1GSQA1DVz4E3sdJiaanqhgxe40CsTudDkePjsGIaPYNrrYgE\nuaQyIgMKl6UU9SlfB/oC0ZC3J1afskqQA8/DOeeqlcnuBU8AT2R6QVV9Bnimhv3ZO8W5M7DPwFzY\nO/j5UcxrAHQCJqc5/xBgVo760vjE5OKWAL8D7g+alGJD5rAquTjPwznnMlRo6+jaAN/HHF8CtM7h\naxDzOksiz1chIoOwbcxuTvW8iJSJSNmiRb5bURWRdXGKDas7sjHIgd2fXpXqGr4mzjmXpUILdBC/\nNCGbpQ3VSV4rqyUQItILmwcxTlX/EddGVUepaqmqlrZr1652vSxGoVzcbGxa7ZnYxqhJJ2M7/Fap\nBl5SAkOGwPr1sHSp5eR8JOecy0ChBbrviR9RtSZ+pFcTqUZubSLPVxCRn2BbEr0C/F+O+lHcUtSn\nXJtIcAOwL/ByqPkO2D3rJ0mRDPZcnHOuhupk94JamMHGHFpYJ2wNW65eg+B1FkZeg+jriMi+wL+x\n2aenBDMyXTop6lO+ge3q/XGoaRPgQmz2Usu4a3kuzjlXS4U2onsW6C4iuyQPiEgHoEfwXC68ha0/\nPj1y/AxsNPdm6LV3B17E1i0fr6opU0cukGafuJeoHOS6YXs53UFMkPNcnHMuRwptRDcaGAo8IyJX\nYXm064H52HIqAERkJ2AOcJ2qXhc6fii2dODHwaFSEVkOoKrjg58JEbkaWyD+Ffb5exhwDjBMVdcG\n1/oRFuQ2xUopdhKplCp8L9hl3YWlWRN3OfAI8DW2kHEIKWZT+ro451wOFVSgU9UVInIYVgJsHDZx\n5GWsBNjyUFPBPiOjI9JrgUNDfx8SPJLnJF/nXhFRrATYpcCXwFBVvSd0biesKD5YpZWoZOmwxi1F\nfco52C/XTqGmzYBHsX2cdkh1Pc/FOedyTFS9/nJdKC0t1bKysnx3o27F5OLWArdhw/BDgYlkOGU2\nnIvz25TOuWqIyDRVLc2kbaHl6FxDEZOLewPYH7gSWI3N4Hks1fklJZVnZHouzjlXRwrq1qVrQEK5\nuO+x/NuoSJOuwO5x53oezjlXj3xE56qXYk2cJhI8glU2CQe5zbEk63+xmZVVeB7OOVePfETn0kux\nJm4OtofRpEjznwN3Y3scVeFr4pxzeeCBzqWWYp+4r7DKJuFFhdtj9dFOJmbySZMmvk+ccy5vPNC5\n1FKsidse+CXwIBbUhmGzLFtFG3ouzjlXADzQuY1SrIlTqo7SbsNGdjcBP0l1Pc/FOecKgAc6Z2Jy\ncQr8E7gF26Bvq1DztljZmFiei3POFRCfdeli18R9DvwMOA2rZn1FuvN9TZxzroB5oHOVcnEJ4FZs\na4cXQk2eB5bFnZvMwy1d6nvFOecKkge6xii6Lm7kSEgkeAtb93YZG2dUClZlewYxk03A83DOuYLn\nObrGJiYX9wPwe+BeKm+73gVbCH5A3HU8D+ecayB8RNeYxOTingH2AkayMci1wGZVlhET5DwP55xr\nYHxE15jErIv7Afhf6O/HAX+l8vY6gK+Jc841WD6iK1Yp6lNGA91ZQC9gW2A88BwxQQ48F+eca7B8\nRFeMUtSn/C+2W214AyfBdrjdAtgy7lqei3PONXA+ois2MXm4pdjMyQOBs7ElBGHtiQlynotzzhUJ\nH9EVm1AeToEngF8DC4OnPwJux/aPi+W5OOdckfFA19ClqE/5BTAE+Fek+bFYtZOUPBfnnCsyHuga\nsphc3DrgL8A1QHhznR8Hx/sSs40OeC7OOVe0PEfXUMXk4t7GdhL4LRuDnAAXAB9jW+tUBDmvT+mc\nayR8RNdQRdbEfYctEwhvhroPVtnkwOi5nodzzjUiPqJrCDJYE9cW+F3w5+bY1jrvEhPkwPNwzrlG\nxUd0hS7FmrhVWEALuxybXXkZsEvctTwP55xrhHxEV8hi8nDrgDuw6iVzI82bAfcRE+Q8D+eca8R8\nRFfIInm4MmAQ8F7w98HARFLMogTPxTnnHAU4ohORHURkvIgsFZFlIvKkiOyY4bk3icgkEVksIioi\nA9K0HSgis0RkjYh8IiLnp2h3koi8JyKrReQLEblKRJrW8O2ll2KfuGXAhcBP2RjkABYAi9Ndz3Nx\nzjlXWIFORFoArwAdgf7AmcDuwGQR2TyDSwzDUlfPV/M6A7G7fE8AxwCPA/eIyAWRdkcHbd7B1lr/\nBbgKuCnzd5WhiROhc2ebZFJeXrHw+ymgE3AnsCFo2gy4GZts0jbuWiUl0KKF5+Kcc47Cu3U5EEsx\n7amqnwGIyAfAbOA8LD2VzpaqukFEdsMK81chIpsANwLjVPXK4PBkEdkOuF5Exqhq8n7hLcAbqjoo\n1K4lcJWIjFDV/5EL4VxcYD4WtZ+JND0K2zuuUh6upASaN7fKKC1bwpln2kjOg5xzzhXWiA44AZia\nDHIAqvo58CZwYnUnq+qG6tpgM+7bAQ9Fjo8DtgZ6gt1CBfZL0a4EG+HlRiQX9xC2GWo4yP0IeAR4\ngZggN2gQLF0K69fbz7vv9iDnnHOBQgt0e2N1h6NmYHfwcvUaxLzOjOBnp3TtgsC7Mof9sTqVoUC3\nPbAi9PQgYBZWo7LKxBPPwznnXFqFFujaAN/HHF8CtM7haxDzOksiz6dqlzzWJnpQRAaJSJmIlC1a\ntCjzHi1fXumvvbEEZSfgDSyZWOXNex7OOecyUmiBDmx3maiUM+hrIHmtuNfJtF1sf1R1lKqWqmpp\nu3btMu9Ry5ZVDt2JzbDsEX3C18Q551xWCi3QxY6UsAFN3MiqJqIjNyJ/X1JNO4CtQs/X3hln2Agt\npBWwafhASQkMGeJ5OOecy1KhBboZbMyNhXUCZubwNYh5nWTObWa6diLSAWiRw/7AJZdUCXRVeC7O\nOedqpNAC3bNAdxGpmFgYBJYewXO58BZW7P/0yPEzsFHamwCq+iUwPUW7BFaUJDd23dVybS1aVA14\nnotzzrlaKbRANxqYBzwjIieKyAnYLPv52JwMAERkJxFZJyLXhE8WkUNFpA+2CBygVET6BMcACNbI\nXQ30F5EbRKSXiFwHnANco6prQ5f8PXCoiNwXtLsYWzD+l5ytoUs69ljLuQ0a5PvEOedcDolqdXMy\n6ldQ7msEcCQ26eNl4CJVnRdq0wH4HLhWVYeHjk8BDo27rqpWmkAiIucBl2D1kb8ERqjqPTH9+QXw\nB6xayzfAGOBGVV2f7n2UlpZqWVlZuibOOedqSESmqWppRm0LLdAVCw90zjlXd7IJdIV269I555zL\nKQ90zjnniprfuqwjIrII+KIWl2iLzQ5tjPy9Nz6N9X2Dv/eavvedVDWjyhwe6AqUiJRlev+52Ph7\nb3zvvbG+b/D3Xh/v3W9dOuecK2oe6JxzzhU1D3SFa1S+O5BH/t4bn8b6vsHfe53zHJ1zzrmi5iM6\n55xzRc0DnXPOuaLmga6AiMgOIjJeRJaKyDIReTKo/VnUgsLbT4jIFyKySkQ+EZGbRWSLfPetvonI\nCyKiInJDvvtSH0TkZyLymogsD37ny0TksHz3q66JSA8RmSQi3wbv+10ROSff/colEWkvIneJyFsi\nsjL4ve4Q066ZiPxJRBYG//+/JSKH5LIvHugKhIi0AF7Bikf3B84Edgcmi8jm+exbPfgtsB7bLeIY\nYCRwAfCiiDSa31ER6Qd0yXc/6ktQWP0ZYBpwMtAXeBzb77FoiUhn4CWgBBgInAK8A9wvIhfks285\nthvwS2zT7NfTtLsf+3e4BjgeWAj8W0T2y1lPVNUfBfAALsQ+7HcLHdsZWAf8Jt/9q+P33i7m2FmA\nAoflu3/19G+wFfA/oF/wvm/Id5/q+P12AFZhO5PkvT/1/N5vAtYCLSPHpwJv5bt/OXyfTUJ/Pjf4\nve4QadMlOH526NgmwCfAs7nqS6P5ttwAnABMVdXPkgdU9XNsI9gT89areqCqi2IOvxP83L4++5JH\ntwIzVPWRfHeknpwDbADuzXdH8mBTbPPmVZHjP1BEd9lUdUMGzU7A/i0eC523DngUOFpENstFX4rm\nH7UI7A18FHN8BtCpnvtSCJL7Cn6c117UAxHpiY1gB+e7L/WoJzALOE1E5gQbKX8mIkPy3bF6MDb4\neaeIbCciW4nIQOBwbC/OxmRv4HNVXRk5PgP7QrBbLl5kk1xcxOVEG+xedtQSoHU99yWvRGR74Drg\nJVUt6k39RKQEuA+4TVU/yXd/6tF2weNPWG52Dpaju1tENlHVv+Szc3VJVT8SkV7AU2z8cpMAzlfV\nR/PWsfxI97mXfL7WPNAVlrjV+xJzrGiJSEtsgsI64Ow8d6c+XAY0B27Md0fqWRNgC2CAqj4ZHHsl\nmJV3hYjcqUHCptiIyO7AE9io5XzsFuaJwL0islpV/5HP/tUzoR4+9zzQFY7vif/20pr4bzxFR0Sa\nAc8CuwCHquqCPHepTgVLR67EEvWbRfIRm4nIVkC5qq7PSwfr1mJsVvGLkeOTsJm32wJf13en6slN\n2AjueFVNBMdeFpGtgb+IyCMZ5reKwRIgbglV69DzteY5usIxA7tfHdUJmFnPfal3wS28J4ADgJ+p\n6od57lJ92AVoBjyEfZlJPsCWXHwP7JufrtW5GSmOJ7/JF/MH/b7A9FCQS3ob2Br4Uf13KW9mADsH\ny6vCOmEzUz+rekr2PNAVjmeB7iKyS/JAcBunR/Bc0QrWyv0DS8afqKpT89yl+vI+0DvmARb8epOj\n/9EL0FPBz6Mjx48GFqjq/+q5P/Xpf8B+IrJp5PhPgdXkaBTTQDyLrSfsmzwgIpsApwKTVHVNLl7E\nb10WjtHAUOAZEbkKu299PTAfm6xQzP6K/aLfCKwQke6h5xYU6y1MVf0BmBI9LiIAX6hqleeKyARg\nMnCfiLQF5gJ9gKMo/tzs3djC+OdE5B4sR3cCtoZyhKquzWfncklE+gR/7Bb8PFZEFgGLVPVVVX1f\nRB4D/hzc1fkcKxaxM3B6zvpRpPneBinI2YwAjsRu4byMLaidl89+1TURmQfslOLpa1V1eP31Jv9E\nRIEbVfWqfPelLolIK+BmLMC1xpYb3KKqD+e1Y/VARI7FJiLtjd2+noNtWXNfMeVkg9/lOK+qaq+g\nTXIy1q+wwgnTgcty+UXPA51zzrmi5jk655xzRc0DnXPOuaLmgc4551xR80DnnHOuqHmgc845V9Q8\n0DnnnCtqHuicc84VNQ90zjnnipoHOucaCRHpJSIaesyqxbXaRq7llSdcwfJA51wREpHfBAHoVzFP\nvwpci9VcrKmVwTWuBb6oxXWcq3Ne1Nm54tQ1+Dkt5rkpta0fqqorgeFgI0VS1yp1Lu98ROdcceoG\nLAdm57sjzuWbBzrnioiI/DHIl3UEWgLrQzm09hle42AReVJE5ojIahH5VkTeFpGb6rLvztUVv3Xp\nXHGZBjwI9Af+A7wYek5izwg3EPk9tmXKl8C/ge+AbYBS4Bjg9znur3N1zgOdc0VEVf8pIlthge5B\nVR2VfC7IpaUkItsA1wFvAIdHNwANNkh1rsHxW5fOFZ/kRJR3szyvI9AU+DRul2tV/a62HXMuHzzQ\nOVd8ugIJ4MMsz5sBLAXOEZFnReRUEWmd8945V8880DlXRERkE2BfYKaqrsnm3GDE1hMYDxwOPAos\nEpGJItI17cnOFTAPdM4Vl05AM7K/bQmAqn6kqn2B1sCRwJPYJJRJIrJZznrpXD3yQOdccdkv+Ple\nbS6iqmtV9SVV/SU2OWVrbPalcw2OBzrnisvWwc9l2ZwkIvuLyK4xx3cD9sGWGyyoffecq3++vMC5\n4pIs+XWjiOwDrABmqOrj1Zz3a6C/iLyNTUr5FtgZOCF4/hxV3VAXHXaurnmgc66IqOprIjIMSD42\nA24Cqgt0z2CfBwcAfbE839fAw8AfVdVLibkGywOdc0VGVe8my50JVPVp4Om66ZFz+SWqvo2Uc41B\nUBllcujQJ6rasYbXagssCh9T1WpLjDmXDz6ic67xmIftH5dUm0onyf3onCt4PqJzzjlX1Hx5gXPO\nuaLmgc4551xR80DnnCudSLAAAAAeSURBVHOuqHmgc845V9Q80DnnnCtqHuicc84Vtf8HMCCPKS0r\np1kAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# %matplotlib inline \n", "%reset\n", "import matplotlib.pyplot as plt # plotting modules\n", "plt.style.use('presentation') # just have in your script for prettier plotting\n", "import numpy as np # our matlab-like module\n", "from scipy.integrate import odeint # integration of ODEs so you don't have to write your finite difference yourself\n", "\n", "#reaction order\n", "def f(y, t):\n", " return -k*y**2\n", "\n", "k = 0.001 # 1/sec\n", "time_start = 0 #s\n", "time_finish = 10 #s\n", "N_points = 50\n", "\n", "times = np.linspace(time_start,time_finish,N_points)\n", "y_A0 = 10. # moles\n", "\n", "y_A_calc = odeint(f, y_A0, times)\n", "plt.plot(times, 1./y_A_calc, 'ro-', label='integrated')\n", "\n", "\n", "inverse_y_A_exact = (1/y_A0 + k*times)\n", "plt.plot(times, inverse_y_A_exact,'k--', label='analytical')\n", "\n", "\n", "plt.xlabel('$t [s]$')\n", "plt.ylabel('$1/C_A [moles]$')\n", "plt.legend()\n", "# plt.savefig('odeint.pdf')\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "# Solving Algebraic equations in Python:\n", "\n", "\n", "Let's say we have a nonlinear equation \n", "\n", "$x^2 - log_{10}(x) = 1$\n", "\n", "The equation can be easily solved using the method of \"thorough observation\"\n", "\n", "x = 1 is the root\n", "\n", "How do we find it numerically?\n", "\n", "\n", "```python\n", "%matplotlib inline\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "from scipy.optimize import fsolve\n", "```\n", "\n", "```python\n", "def func(x):\n", " return x**2 - np.log10(x) - 1\n", "fsolve(func, 0.5)\n", "```\n", "\n", "And that is it! Nothing more! Imagine how useful that would be for the Thermo 346 course =)\n" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.0\n" ] } ], "source": [ "%matplotlib inline\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "from scipy.optimize import fsolve\n", "\n", "def func(x):\n", " return x**2 - np.log10(x) - 1\n", "answer = fsolve(func, 0.5)\n", "print(answer[0])\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Problem 3: CSTR \n", "\n", "From solving ODE now we are going to solve \"Simple\" algebraic equations.\n", "\n", "Given a continuously stirred tank reactor with a volume of $66 m^3$ where the reaction $A \\to B$ occurs, at a rate of $−r_A= k C^2_A$ ($k=3 L/mol/h$), with an entering molar flow of $F_{A0} = 5 mol/h$ and a volumetric flowrate of $\\upsilon = 10 L/h$, what is the exit concentration of A?\n", "\n", "From a mole balance we know that at steady state $0=F_{A0}−F_A+ V \\cdot r_A$. That equation simply states the sum of the molar flow of A in in minus the molar flow of A out plus the molar rate A is generated is equal to zero at steady state. This is directly the equation we need to solve. We need the following relationship:" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The exit concentration is 0.005 mol/L\n" ] } ], "source": [ "#Copyright (C) 2013 by John Kitchin\n", "\n", "from scipy.optimize import fsolve\n", "\n", "Fa0 = 5.0\n", "v0 = 10.\n", "\n", "V = 66000.0 # reactor volume L^3\n", "k = 3.0 # rate constant L/mol/h\n", "\n", "def func(Ca):\n", " \"Mole balance for a CSTR. Solve this equation for func(Ca)=0\"\n", " Fa = v0 * Ca # exit molar flow of A\n", " ra = -k * Ca**2 # rate of reaction of A L/mol/h\n", " return Fa0 - Fa + V * ra\n", "\n", "# CA guess that that 90 % is reacted away\n", "CA_guess = 0.1 * Fa0 / v0\n", "CA_sol, = fsolve(func, CA_guess)\n", "\n", "print('The exit concentration is {0} mol/L'.format(CA_sol))" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] } ], "metadata": { "celltoolbar": "Slideshow", "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.4" } }, "nbformat": 4, "nbformat_minor": 2 }