pH calculation speedup.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import pHcalc.pHcalc as ph\n",
    "import time"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create a simple system consisting of one polyprotic acid."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = ph.Acid(pKa=[1., 3., 9.], charge=1, conc=0.001)\n",
    "s = ph.System(a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The pHsolver function uses the scipy.optimize.minimize function, which is probably overkill, and as such, it is a little slow."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 4.05 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "100 loops, best of 3: 1.88 ms per loop\n"
     ]
    }
   ],
   "source": [
    "%%timeit\n",
    "s.pHsolve()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3.2109146118164045"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "s.pH"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Under the hood, the minimize function is using the private function System.\\_s\\_diff\\_pos\\_neg as the objective function. In the case of the minimize function, it passes a single pH value to this function, and it then returns a difference between the positively and negatively charged species in solution. In a system at equilibrium, this difference should be zero. Now, you can also pass a Numpy array of pH values directly to the System.\\_s\\_diff\\_pos\\_neg function. The function then returns an array of charge differences. You can then use Numpy's argmin function to find the index of lowest value in that array, which is then the index of the optimal pH. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "phs = np.linspace(0,14,1400)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAD8CAYAAAB3u9PLAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xt0nNV97vHvb0Ya3W+W5AuS7xZ2HMxVgGm4hBAIEE5M\nmiaFpoEkHFxaaNOctClpu9LVdU5Oc9qcpKEhUBpI4CSNC7kUJ6EFAi0smhhsEzDY2Fj4JtmSLVv3\n+2hmnz/mlSOMpHkljfXO5fmspaWZ9yL9RkuaR3vvd+/XnHOIiIiEgi5ARETSgwJBREQABYKIiHgU\nCCIiAigQRETEo0AQERHAZyCY2bVmtsfMmszs7gn2m5nd4+3fYWbnj9v3kJkdM7PXTzlnnpk9bWZ7\nvc9Vs385IiIyU0kDwczCwL3AdcBa4GYzW3vKYdcBDd7HRuC+cfu+A1w7wZe+G3jGOdcAPOM9FxGR\ngPhpIVwENDnn9jnnRoBNwIZTjtkAPOIStgCVZrYIwDn3PNAxwdfdADzsPX4YuHEmL0BERFIjz8cx\ndUDzuOctwMU+jqkDWqf4ugucc2P724AFEx1kZhtJtDooKSm5YM2aNT5KzlyxuGNXaw+LKgqpKS0I\nuhyRrPLm0V4K8sIsrS4OupQ5tX379uPOudpkx/kJhNPOOefMbMI1NJxzDwAPADQ2Nrpt27bNaW1z\nbWQ0zpl/+W/86QdWc+eVq4IuRySrXPI3z3Dpqhr+7qPnBF3KnDKzg36O89NldBhYPO55vbdtusec\n6uhYt5L3+ZiPWrJeJC9EftjoHx4NuhSRrNM7NEpZYX7QZaQtP4GwFWgws+VmFgFuAjafcsxm4Bbv\naqP1QPe47qDJbAZu9R7fCjw+jbqzWnEkT4EgkmKxuKNveJSywrToGElLSQPBOTcK3AU8CbwBPOqc\n22lmd5jZHd5hTwD7gCbgn4A/GDvfzL4P/BJYbWYtZnabt+vLwNVmthd4v/dcgJJImP6RWNBliGSV\nPu+fLAXC5Hz9ZJxzT5B40x+/7f5xjx1w5yTn3jzJ9hPAVb4rzSElBXkMjKiFIJJKvUNRAMrVZTQp\nzVROQ8UFefQNq4Ugkkrdg14gFKmFMBkFQhoqiYQZ0BiCSEqNBUJFUSTgStKXAiENFUfyNIYgkmI9\nJwNBXUaTUSCkodKCsMYQRFLsZAuhWIEwGQVCGiou0GWnIqnWrRZCUgqENFQSCdOvQWWRlOoaiBIO\nGSWRcNClpC0FQhoqjuQxGI0Ri0+4moeIzED3YJSKonzMLOhS0pYCIQ2VFiQui9M4gkjqjAWCTE6B\nkIaKCxJN2gFdaSSSMt2DUcoVCFNSIKShkkiihaCBZZHU6VELISkFQhoqjqiFIJJq3YNRKhUIU1Ig\npKGxMYQ+tRBEUkZjCMkpENJQcYG6jERSKR53CgQfFAhpSC0EkdTqGxkl7jQpLRkFQhoaW69dgSCS\nGt0DmqXshwIhDZ1sIQwpEERS4ddLXysQpqJASEPFkTAhS9z/VURmTyud+qNASENmRmlBnrqMRFJE\nC9v5o0BIU2WF+WohiKSIlr72R4GQpkoL8k7eA1ZEZmcsEDQxbWoKhDRVVqguI5FU6R6Mkheyk6sA\nyMQUCGmqVIEgkjJdWvraFwVCmiotyNNlpyIpolnK/igQ0lRZYR49CgSRlOjsH2FeSSToMtKeAiFN\nlRXm0zesQWWRVOjoH6FKgZCUAiFNlRbkMRSNE43Fgy5FJON1Dowwr1iBkIwCIU2VasVTkZRwztHZ\nH1ULwQcFQpoaW+BOk9NEZqd/JMZILM68Eg0qJ6NASFMKBJHU6OgbAaBKXUZJKRDSVGlB4r8ZzUUQ\nmZ2OgUQg6Cqj5BQIaar05D0RdKWRyGx09nstBAVCUr4CwcyuNbM9ZtZkZndPsN/M7B5v/w4zOz/Z\nuWZ2rpltMbNXzGybmV2UmpeUHdRlJJIaHV4g6Cqj5JIGgpmFgXuB64C1wM1mtvaUw64DGryPjcB9\nPs79W+CvnXPnAl/0nounrECBIJIKnWNdRqUKhGT8tBAuApqcc/uccyPAJmDDKcdsAB5xCVuASjNb\nlORcB5R7jyuAI7N8LVmlVLfRFEmJjv4R8kJ28p8smZyfn1Ad0DzueQtwsY9j6pKc+8fAk2b2FRLB\n9BsTfXMz20ii1cGSJUt8lJsdivLDhEOm9YxEZqlzIDFLWQvbJRfkoPLvA591zi0GPgs8ONFBzrkH\nnHONzrnG2traOS0wSGN3TdM9EURmp6Nfs5T98hMIh4HF457Xe9v8HDPVubcCP/IeP0aie0nGKS3I\no1ddRiKzkpilrElpfvgJhK1Ag5ktN7MIcBOw+ZRjNgO3eFcbrQe6nXOtSc49AlzhPX4fsHeWryXr\nlBVqCWyR2TrRP6w5CD4lHUNwzo2a2V3Ak0AYeMg5t9PM7vD23w88AVwPNAEDwKemOtf70rcDXzez\nPGAIb5xAfk13TROZvc6BqGYp++Rr2N059wSJN/3x2+4f99gBd/o919v+AnDBdIrNNaUFeRz3pt2L\nyPTF4o6uAd0LwS/NVE5jZYX5GlQWmYWewShxp3WM/FIgpLGKony6BxUIIjOldYymR4GQxiqK8ukZ\nGiXRIyci03W8dxiAmtKCgCvJDAqENFZelEcs7jSwLDJD7X1eIJSpheCHAiGNVRQlrp3u0aWnIjMy\n1kKoVQvBFwVCGhsLhO4BjSOIzER73zDhkGlQ2ScFQhorL/QCQQPLIjNyvHeE6pIIoZDWMfJDgZDG\nyosUCCKz0d43TG2Zuov8UiCksV+PISgQRGaivXdYVxhNgwIhjVUUe4GgFoLIjBxXC2FaFAhprDSS\nh5m6jERmwjmnQJgmBUIaC4WM8sJ8tRBEZqB7MEo05tRlNA0KhDSn5StEZqZ9bA6CWgi+KRDSnAJB\nZGZOzlIu1RwEvxQIaa68KE+BIDIDYy2E+Woh+KZASHNjC9yJyPSc7DIqLQy4ksyhQEhz6jISmZnj\nfSNEwiHKi3zdB0xQIKS98kIFgshMJCalRTDTshV+KRDSXHlRPiOjcYaisaBLEckox3qHdIXRNCkQ\n0tzJ5SvUShCZlqM9Qyys0PjBdCgQ0lyFFrgTmZHW7iEWVRQFXUZGUSCkOa14KjJ9fcOj9A6NsqBc\nLYTpUCCkOa14KjJ9bd1DACxSl9G0KBDS3FggdOmuaSK+jQWCxhCmR4GQ5uZ5t/7rVCCI+NbWoxbC\nTCgQ0lxZYR4hg87+kaBLEckYbd2DABpDmCYFQpoLeTcI7xxQIIj41do9RFVxPoX54aBLySgKhAxQ\nWZyvQBCZhrbuIRbqktNpUyBkgHklETrUZSTiW1vPkMYPZkCBkAEqiyO6ykhkGtq6hzR+MAO+AsHM\nrjWzPWbWZGZ3T7DfzOweb/8OMzvfz7lm9odmttvMdprZ387+5WSnecVqIYj4NRSNcaJ/RC2EGUi6\nLqyZhYF7gauBFmCrmW12zu0ad9h1QIP3cTFwH3DxVOea2ZXABuAc59ywmc1P5QvLJlUliRaCc04r\nN4okcawncR8EzUGYPj8thIuAJufcPufcCLCJxBv5eBuAR1zCFqDSzBYlOff3gS8754YBnHPHUvB6\nslJVcT4jsTj9I1rxVCSZI94lp2ohTJ+fQKgDmsc9b/G2+TlmqnPPBC4zsxfN7Dkzu3Cib25mG81s\nm5lta29v91Fu9qkq8SanqdtIJKnmjgEA6quKA64k8wQ5qJwHzAPWA38KPGoT9Ic45x5wzjU65xpr\na2vnusa0UHVytrICQSSZls5BzOCMSrUQpsvPveUOA4vHPa/3tvk5Jn+Kc1uAHznnHPCSmcWBGiA3\nmwFTmFeSWM9Iy1eIJNfcOcDC8kIK8jQpbbr8tBC2Ag1mttzMIsBNwOZTjtkM3OJdbbQe6HbOtSY5\n91+BKwHM7EwgAhyf9SvKQpXF6jIS8aulc5D6Kk1Km4mkLQTn3KiZ3QU8CYSBh5xzO83sDm///cAT\nwPVAEzAAfGqqc70v/RDwkJm9DowAt3qtBTnF2AJ3uvRUJLmWjgHWr6gOuoyM5KfLCOfcEyTe9Mdv\nu3/cYwfc6fdcb/sI8LvTKTZXlRflYwZdGkMQmdLIaJzWniHq52lAeSY0UzkDhENGZVE+HQoEkSm1\ndg/iHOoymiEFQoZIrHiqQWWRqTR3JOYgLNYlpzOiQMgQVSURDSqLJNHSmZiDsHieWggzoUDIEFXF\n+WohiCTR3DlAOGQs1MJ2M6JAyBCJJbCHgy5DJK21dA5yRmUheWG9tc2EfmoZorq0gBN9I8TjujJX\nZDIHTwxo/GAWFAgZoqa0gNG4o3tQ3UYiE3HOsa+9jxW1JUGXkrEUCBmipjQxOe2Euo1EJtQ5EKVn\naJTlNaVBl5KxFAgZora0AID2Xl1pJDKR/cf7AFhRoxbCTCkQMkRNWSIQjvephSAykbfa+wFYrkCY\nMQVChqj27omgQBCZ2P7j/eSHTbOUZ0GBkCGqiiOEQ6ZAEJnE/vZ+lswr1iWns6CfXIYIhYx5JRFO\n9GkMQWQi+473aUB5lhQIGaSmtEAtBJEJxOKOAycGdMnpLCkQMkhNaYR2tRBE3uFI1yAjo3FdYTRL\nCoQMUlNawPFetRBETvVWe+KSU11hNDsKhAxSUxrhRP8wurGcyNvtaesFYPXCsoAryWwKhAxSU1rA\nUDRO/0gs6FJE0sqetl4WlBecvP+4zIwCIYPUeLOV1W0k8na723pZvbA86DIyngIhg2i2ssg7jcbi\nNLX3sUbdRbOmQMggYwvcKRBEfu3AiX5GRuOsXqBAmC0FQgb59QJ3CgSRMbs1oJwyCoQMUl1aQDhk\nHO1RIIiM2dPWSzhkrJqvWcqzpUDIIOGQUVtaQFvPUNCliKSN3W29LKsupjA/HHQpGU+BkGEWVBRy\nVIEgctKuIz28a5GuMEoFBUKGWVBWoEAQ8RzvG+Zw1yDn1FcGXUpWUCBkmIUVhRpDEPG81tINwLr6\nioAryQ4KhAyzoLyQ7sEoQ1HNVhZ5taULMzirToGQCgqEDLOgvBCAtm51G4nsaOlmVW0ppQV5QZeS\nFRQIGWahFwgaR5Bc55xjR0sXZ2v8IGUUCBlmQXlicpouPZVcd6R7iON9I5yzWN1FqeIrEMzsWjPb\nY2ZNZnb3BPvNzO7x9u8ws/Once7nzMyZWc3sXkpuWFChFoIIwCuHugDUQkihpIFgZmHgXuA6YC1w\ns5mtPeWw64AG72MjcJ+fc81sMXANcGjWryRHlBXkUZQf1pVGkvO2HuigKD/Mu8/QHIRU8dNCuAho\ncs7tc86NAJuADaccswF4xCVsASrNbJGPc78GfB7QHV98MjMWVhSqy0hy3ov7O7hgaRX5YfV8p4qf\nn2Qd0DzueYu3zc8xk55rZhuAw865V6f65ma20cy2mdm29vZ2H+VmvwXlBRzVVUaSw7oHouxu6+Gi\n5fOCLiWrBBKtZlYM/DnwxWTHOucecM41Oucaa2trT39xGWBBuVoIktu2HezAORQIKeYnEA4Di8c9\nr/e2+Tlmsu0rgeXAq2Z2wNv+spktnE7xuWpRRRFHe4aIx9XTJrnppf0dRMIhzl2sAeVU8hMIW4EG\nM1tuZhHgJmDzKcdsBm7xrjZaD3Q751onO9c595pzbr5zbplzbhmJrqTznXNtqXph2ayuspBozNGu\nG+VIjtqyv4Oz6yu0wmmKJQ0E59wocBfwJPAG8KhzbqeZ3WFmd3iHPQHsA5qAfwL+YKpzU/4qckxd\nVREALZ2DAVciMvc6+0fY0dLFe1bpSvVU8zXf2zn3BIk3/fHb7h/32AF3+j13gmOW+alDEuoqiwE4\n3DXIBUurAq5GZG49v7cd5+C9qzWmmGq6XisDjbUQDquFIDnouT3tVBXna0LaaaBAyEClBXlUFOVz\nuGsg6FJE5lQ87njuzXYuP7OWcMiCLifrKBAyVH1VkVoIknNeP9LNif4RdRedJgqEDFVXWcThLgWC\n5Jandh4lZHB5gwLhdFAgZKg6r4WQGM8XyX7OOZ54rZX1K6qpLi0IupyspEDIUHWVRfSPxOgejAZd\nisic2N3Wy77j/Xzw7EVBl5K1FAgZql5zESTH/GxHKyGDD7xbCxqcLgqEDDV+LoJItnPO8dMdR1i/\nopoadRedNgqEDKXZypJLXtrfwYETA/zm+fVBl5LVFAgZqqo4n7KCPJo7NBdBst+/bG2mrCCP69ep\nu+h0UiBkKDNjaU0xB070B12KyGnVPRDlZ6+18qFzz6A44mu1HZkhBUIGW1pdwsETaiFIdvvhyy0M\nj8a56cIlQZeS9RQIGWxZdTHNHQOMxuJBlyJyWozG4jz4wn4al1axrr4i6HKyngIhgy2tLmE07jjS\npbunSXb62WutHO4a5PeuWBl0KTlBgZDBllWXAGgcQbJSPO74x+f2sbK2hKvWzA+6nJygQMhgy6oT\ncxEOKhAkC/3stVZ2tfbwB+9dRUgrm84JBUIGqy0roCg/zAENLEuWicbifOWpPaxZWMaN59UFXU7O\nUCBkMDNjaXWxWgiSdb675SAHTwzw+WtX674Hc0iBkOGWVherhSBZpbV7kP/71Jtc1lDDlas1djCX\nFAgZbllNCYdO6NJTyQ7OOf7q8Z2MxuN86cZ1mKl1MJcUCBluVW0pI7E4zVrTSLLAo9uaeWrXUT77\n/jNZ4l00IXNHgZDhVs0vBWDv0d6AKxGZnTeP9vJXm3dy6aoabr9sRdDl5CQFQoYbC4Sm9r6AKxGZ\nueN9w9z28FZKC/L46sfO0WWmAdFKURmurDCfheWFNB1VIEhmGhyJcfsj22jvHWbTxkuYX14YdEk5\nSy2ELNCwoJS9xxQIknn6h0f55Ldf4tXmLv7+t8/j3MWVQZeU0xQIWWDV/FLeau8jHndBlyLi24m+\nYT7x4ItsO9jJ1377XK49S/c6CJoCIQusml/KwEiMI9260kgyw64jPXzoG//FziM9fOPm89hwrmYj\npwMFQhZomF8GQJO6jSTNxeOOB1/Yz4e/+V/E4o7H7riE69YtCros8WhQOQs0nLz0tI/3amanpKk9\nbb188fHXeXF/B1etmc/ffGQd88s0gJxOFAhZoKokwvyyAt5o6wm6FJF3ONY7xNd/vpfvv3SI0oI8\n/vYjZ/PRxnrNQk5DvgLBzK4Fvg6EgW855758yn7z9l8PDACfdM69PNW5ZvZ3wH8DRoC3gE8557pS\n8aJy0dozytl1RIEg6ePQiQH+8fm3eGx7C7G445ZLlvGZqxqoKokEXZpMImkgmFkYuBe4GmgBtprZ\nZufcrnGHXQc0eB8XA/cBFyc592ngC865UTP7P8AXgD9L3UvLLWsXlfPC3uMMj8YoyAsHXY7kqFjc\n8fyb7Wzaeoindx0lLxTiIxfU8XuXr2RZTUnQ5UkSfloIFwFNzrl9AGa2CdgAjA+EDcAjzjkHbDGz\nSjNbBCyb7Fzn3FPjzt8C/NZsX0wuW3tGOaNxx96jfZxVp3vPytxq7hjgsW3NPLa9hdbuIapLItx+\n+Qo+/Z7lLNBEs4zhJxDqgOZxz1tItAKSHVPn81yATwP/MtE3N7ONwEaAJUuW+Cg3N61dVA7ArtYe\nBYLMieHRGD/fdYxNWw/xQtNxAC5vqOWLN6zlqnctIJKnixgzTeCDymb2F8Ao8L2J9jvnHgAeAGhs\nbNTMq0ksrS6hOBLWOIKcdnuP9rJpazM/ermFzoEodZVFfOaqBj7auJi6yqKgy5NZ8BMIh4HF457X\ne9v8HJM/1blm9kngBuAqr7tJZigcMtYsLGNXqwJBUm9wJMZPdhxh00uHePlQF/lh4+q1C/hY42Iu\na6jVXc2yhJ9A2Ao0mNlyEm/mNwG/c8oxm4G7vDGCi4Fu51yrmbVPdq539dHngSucc7rlVwqsPaOc\nx391hHjcabVISYm27iEe/uUB/vnFQ3QPRllZW8JfXP8uPnx+HTWlBUGXJymWNBC8q4DuAp4kceno\nQ865nWZ2h7f/fuAJEpecNpG47PRTU53rfelvAAXA0971yFucc3ek8sXlmnV1FXx3yyH2n+hnZW1p\n0OVIBnujtYd/fO4tfrqjlbhzXLN2IZ++dDkXLqvS/IEs5msMwTn3BIk3/fHb7h/32AF3+j3X275q\nWpVKUucurgLglUNdCgSZkX3tfXzt53v5yatHKC3I45ZLlvGp9yxj8TzdvSwXBD6oLKmzan4ppQV5\n/Kq5k49cUB90OZJBjnQNcs8ze3lsewuRcIg7r1zJxstWUlGcH3RpMocUCFkkHDLOrq/glWZN+BZ/\nTvQN883/fIv/t+UgOPjE+qXceeUqass0PpCLFAhZ5tzFlTzw/D6GojEK8zVjWSbWMxTlW8/v48EX\n9jMYjfFbF9TzR1c1UF+lrqFcpkDIMuctqWI07nj9cDeNy+YFXY6kmcGRGA//8gD3/edbdA9G+eC6\nRXz26jNP3ptbcpsCIcuM3YLwV4e6FAhy0vBojEe3NnPPs0209w7z3tW1/Mk1qzWrXd5GgZBlassK\nWFpdzEsHOrj98hVBlyMBi8bi/HB7C//wbBOHuwa5cFkV9/7O+Vy0XP8syDspELLQ+uXV/PvONk1Q\ny2HRWJzNrxzhnmf3cvDEAOcsruRvfnMdlzXUaB6BTEqBkIUuWVnNv2xr5o22Ht59hroEcknPUJTv\nv3iI7/ziAK3dQ7z7jHIevLWR962ZryCQpBQIWejiFYnugF++dUKBkAOcc+xo6eax7c38+OXD9I/E\nuGRFNV/68Fm898z5aiWKbwqELLSooohl1cVs2dfBf79M4wjZal97H0/tOsqPXz7MnqO9FOSFuH7d\nIm67dLkGi2VGFAhZav2Kap54rZVY3GklyizRMxRl+8FOtuw7wc93HeWt9n4AzllcyZc+fBY3nH0G\nFUWaWSwzp0DIUpesrGbT1mbe0A1zMlL3YJQ9bb3sbuthd1svO1q62HWkh7iDvJBx8Yp5fGL9Ut6/\ndoEmk0nKKBCy1G+srAHguTfbFQhpIh539A6N0jkwQufACF0DUTr6E4+PdA1xuGuAw12DHO4cpHMg\nevK8iqJ81i4q5w/f18DFy+dx7pJKiiP605XU029VlqotK+Ds+gqe3X2MO6/MrIVlo7E4fUOjDEZj\njIzGGYnFGY7GGYnFGI7GGY7FGRmNMxpzxJzDOUcs7oi7xJtuzDnizhH3tiX2JT5icU7ui40/Zmxb\nfOxrcvLx2Pa4w/sab/++Y18zFndEY3EGozEGR2IMj8YZHIklnnuvZTJF+WHqqoqoqyzi7PpKlswr\nZvXCMtYsLGNheaGuEJI5oUDIYleuns8/PLuXzv4RqkoigdYyOBJj//F+DnUMcKx3iGM9wxztGeJY\n7zBdg1F6h6L0Do3SOxRlKDr5G+fpEA4ZIYOQGeGQETYj5G1L7LOTn0MhTu4P29g2I+xtzw+HKC3I\no6a0gKL8cOIjEqbQe1xSEKaqOMK8kgiVxflUFUeoKo5QXpSnN30JnAIhi125Zj5ff2Yvz73Zzo3n\n1c3Z9z3eN8wrh7p4pbmLV1u6eOtYH0e6h952TDhk1JRGmF9WSFVJhPqqIsoL8ygtyKOsMJ+ywjyK\n8sNE8kJE8kIU5IW9z4nnkXCI/HCIcCjxRj72pm3em3jYDDv5Ro73Bj+2nbe98YtIggIhi51dV0FN\naYRndx87rYEQjztePtTJM7uP8R+7j7G7rRdIvOmuXlDGRcvnsaK2lBW1JSyrLmFBeSHzSiK6+kkk\nzSgQslgoZFxx5nx+/sZRorE4+eFQSr/+ka5BHtvWwmPbm2npHCQvZDQuq+LPrl3DBUurWFdXQVFE\nS3CLZAoFQpa75t0L+OHLLfzirRNccWZtSr7m/uP9fPM/mvjxrw4zGndcuqqGP7lmNe9713zKC3Ud\nvEimUiBkuSvOrKWsII+fvnpk1oHQPRjla0+/ySO/PEB+OMTvrl/KbZcu1/12RbKEAiHLFeaHuXrt\nAp7c2caXPryOSN7Muo3+/fU2/vJfX+NE/wgfv3gJf3RVA/PLClNcrYgEKbWdypKWbjhnET1Do7zQ\n1D7tc4eiMf7yX1/jju9uZ2FFIZvvvJT/deM6hYFIFlILIQdcuqqWiqJ8fvJqK+9bs8D3eSf6hrn9\nkW28fKiLjZev4E+uWT3jFoaIpD/9deeASF6ID569iH97vZWeoWjyE4DmjgE+ct8v2Hmkh29+/Hz+\n/Pp3KQxEspz+wnPETRcuZiga5/FfHU567JGuQW7+py10DkT559vXc/26RXNQoYgETYGQI9bVVbB2\nUTnff6kZ59ykx7X3DvPxb71I90CU7952MRcsrZrDKkUkSAqEHGFm3HzRYna19vDa4e4JjxkejfH7\n391Oa/cg3/n0hayr1yqpIrlEgZBDNpxXR1F+mId/cfAd+5xz/NXjO9l2sJOvfPQcLlg6L4AKRSRI\nCoQcUl6Yz29fuJjNrx6m7ZTF5n6wvYVNW5u588qV3HD2GQFVKCJBUiDkmNsuXU4s7vj2f+0/ua25\nY4C//skuLlo+j/9x9eoAqxORICkQcsziecVcv24R//ziIXqGosTjjs89+ioAX/3YOVqBVCSH+QoE\nM7vWzPaYWZOZ3T3BfjOze7z9O8zs/GTnmtk8M3vazPZ6n3U5yxy544qV9A6P8u0XDvCD7S28dKCD\nL96wVvfmFclxSQPBzMLAvcB1wFrgZjNbe8ph1wEN3sdG4D4f594NPOOcawCe8Z7LHDirroIPvHsB\nX/v5m9z9ox1cuKyKjzbWB12WiATMTwvhIqDJObfPOTcCbAI2nHLMBuARl7AFqDSzRUnO3QA87D1+\nGLhxlq9FpmFsrCDu4K8/dJZu3ygivtYyqgOaxz1vAS72cUxdknMXOOdavcdtwISL7JjZRhKtDpYs\nWeKjXPFj9cIy/ueNZ2HA2jPKgy5HRNJAWixu55xzZjbh9Fnn3APAAwCNjY2TT7GVafvE+qVBlyAi\nacRPl9FhYPG45/XeNj/HTHXuUa9bCe/zMf9li4hIqvkJhK1Ag5ktN7MIcBOw+ZRjNgO3eFcbrQe6\nve6gqc7dDNzqPb4VeHyWr0VERGYhaZeRc27UzO4CngTCwEPOuZ1mdoe3/37gCeB6oAkYAD411bne\nl/4y8KgJIxLrAAAEyklEQVSZ3QYcBD6W0lcmIiLTYlOtfJluGhsb3bZt24IuQ0Qko5jZdudcY7Lj\nNFNZREQABYKIiHgUCCIiAigQRETEo0AQERFAgSAiIh4FgoiIAAoEERHxKBBERARQIIiIiEeBICIi\ngAJBREQ8CgQREQEUCCIi4smo5a/NrJ3EvRNmogY4nsJyTrdMqjeTaoXMqjeTaoXMqjeTaoXZ1bvU\nOVeb7KCMCoTZMLNtftYDTxeZVG8m1QqZVW8m1QqZVW8m1QpzU6+6jEREBFAgiIiIJ5cC4YGgC5im\nTKo3k2qFzKo3k2qFzKo3k2qFOag3Z8YQRERkarnUQhARkSkoEEREBMiRQDCza81sj5k1mdndQdcz\nGTNbbGb/YWa7zGynmX0m6JqSMbOwmf3KzH4adC3JmFmlmf3AzHab2RtmdknQNU3FzD7r/R68bmbf\nN7PCoGsaY2YPmdkxM3t93LZ5Zva0me31PlcFWeN4k9T7d97vwg4z+7GZVQZZ45iJah2373Nm5sys\n5nR876wPBDMLA/cC1wFrgZvNbG2wVU1qFPicc24tsB64M41rHfMZ4I2gi/Dp68C/O+fWAOeQxnWb\nWR3wR0Cjc+4sIAzcFGxVb/Md4NpTtt0NPOOcawCe8Z6ni+/wznqfBs5yzp0NvAl8Ya6LmsR3eGet\nmNli4Brg0On6xlkfCMBFQJNzbp9zbgTYBGwIuKYJOedanXMve497Sbxh1QVb1eTMrB74IPCtoGtJ\nxswqgMuBBwGccyPOua5gq0oqDygyszygGDgScD0nOeeeBzpO2bwBeNh7/DBw45wWNYWJ6nXOPeWc\nG/WebgHq57ywCUzyswX4GvB54LRdCZQLgVAHNI973kIav8mOMbNlwHnAi8FWMqW/J/ELGg+6EB+W\nA+3At70urm+ZWUnQRU3GOXcY+AqJ/wZbgW7n3FPBVpXUAudcq/e4DVgQZDHT9Gng34IuYjJmtgE4\n7Jx79XR+n1wIhIxjZqXAD4E/ds71BF3PRMzsBuCYc2570LX4lAecD9znnDsP6Ce9ujTexut/30Ai\nyM4ASszsd4Otyj+XuJ49I65pN7O/INFd+72ga5mImRUDfw588XR/r1wIhMPA4nHP671tacnM8kmE\nwfeccz8Kup4pvAf4kJkdINEN9z4z+26wJU2pBWhxzo21uH5AIiDS1fuB/c65dudcFPgR8BsB15TM\nUTNbBOB9PhZwPUmZ2SeBG4CPu/SdlLWSxD8Gr3p/b/XAy2a2MNXfKBcCYSvQYGbLzSxCYmBuc8A1\nTcjMjEQf9xvOua8GXc9UnHNfcM7VO+eWkfiZPuucS9v/YJ1zbUCzma32Nl0F7AqwpGQOAevNrNj7\nvbiKNB4E92wGbvUe3wo8HmAtSZnZtSS6PD/knBsIup7JOOdec87Nd84t8/7eWoDzvd/plMr6QPAG\nje4CniTxB/Woc25nsFVN6j3AJ0j8t/2K93F90EVlkT8EvmdmO4Bzgf8dcD2T8loyPwBeBl4j8bea\nNkstmNn3gV8Cq82sxcxuA74MXG1me0m0cL4cZI3jTVLvN4Ay4Gnvb+3+QIv0TFLr3Hzv9G0liYjI\nXMr6FoKIiPijQBAREUCBICIiHgWCiIgACgQREfEoEEREBFAgiIiI5/8DauT6U6qOn9MAAAAASUVO\nRK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x2256f7765c0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(phs, s._diff_pos_neg(phs))\n",
    "plt.ylim(-0.001, 0.01)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 loop, best of 3: 732 µs per loop\n"
     ]
    }
   ],
   "source": [
    "%%timeit -n 1\n",
    "s._diff_pos_neg(phs).argmin()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3.2122944960686208"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "idx = s._diff_pos_neg(phs).argmin()\n",
    "phs[idx]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is a little bit faster, and if you lower the resolution of the pH values, this method will be faster still. However, for subsequent pH calculations, you can minimize the pH window of interest, assuming the concentrations are not changing much for each step. \n",
    "\n",
    "First, decrease the concentration of the acid by a small amount."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "a.conc = 0.0008"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, we will narrow down our pH window, which leads to much faster minimization."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "phs = np.linspace(phs[idx]-0.5, phs[idx]+0.5, 100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 loop, best of 3: 104 µs per loop\n"
     ]
    }
   ],
   "source": [
    "%%timeit -n 1\n",
    "s._diff_pos_neg(phs).argmin()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The pH estimation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3.2779510617251866"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "idx = s._diff_pos_neg(phs).argmin()\n",
    "phs[idx]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And using the normal minimization method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1000 loops, best of 3: 1.76 ms per loop\n"
     ]
    }
   ],
   "source": [
    "%%timeit\n",
    "s.pHsolve()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3.2815544128417948"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "s.pH"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here we'll just set up three different scenarios to see how fast these different calculations can run."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The standard method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "def type0(n, a, s):\n",
    "    cons = np.linspace(0.001, 0.00001, n)\n",
    "    type0 = []\n",
    "    for i in cons:\n",
    "        a.conc = i\n",
    "        s.pHsolve()\n",
    "        type0.append(s.pH)\n",
    "    return type0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using the guess keyword argument."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "def type1(n, a, s):\n",
    "    cons = np.linspace(0.001, 0.00001, n)\n",
    "    type1 = []\n",
    "    pH = 7.\n",
    "    for i in cons:\n",
    "        a.conc = i\n",
    "        s.pHsolve(guess=pH)\n",
    "        pH = s.pH\n",
    "        type1.append(pH)\n",
    "    return type1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Only using the pHsolve method once and then using Numpy arrays for subsequent calcs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "def type2(n, a, s):\n",
    "    cons = np.linspace(0.001, 0.00001, n)\n",
    "    type2 = []\n",
    "\n",
    "    # Get an initial pH value using the minimizer\n",
    "    a.conc = cons[0]\n",
    "    s.pHsolve()\n",
    "    pH = s.pH\n",
    "    type2.append(pH)\n",
    "\n",
    "    # For subsequent pHs just use a narrower pH window\n",
    "    # Of course, this assumes that your pH will not change by more that 0.5 for any step\n",
    "    # in your calculation.\n",
    "    phs = np.linspace(pH-0.5, pH+0.5, 200)\n",
    "    for i in cons[1:]:\n",
    "        a.conc = i\n",
    "        idx = s._diff_pos_neg(phs).argmin()\n",
    "        newpH = phs[idx]\n",
    "        # Adjust the pH window by the difference\n",
    "        phs = phs - (pH - newpH)\n",
    "        pH = newpH\n",
    "        type2.append(pH)\n",
    "    return type2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The setup."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "a = ph.Acid(pKa=[1., 3., 9.], charge=1, conc=0)\n",
    "s = ph.System(a)\n",
    "n = 100"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Timing tests"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 loop, best of 3: 168 ms per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit type0(n, a, s)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "10 loops, best of 3: 124 ms per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit type1(n, a, s)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "100 loops, best of 3: 9.82 ms per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit type2(n, a, s)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The pH values from the second two methods are quite close, and the differences average out to be very close to zero... However, there is some nonrandomness to the differences that I can't explain at this time. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Ave =  -0.000120571102034\n",
      "Std =  0.0014133203688\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAD8CAYAAABkbJM/AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsvWmULFd1LvidGHKo+Y7SVYEGkLi2AOMyAtmI94wBc5Ft\nTJlu94NlbNx+DcbY7dVgLhLtYV08gZewn58NBg8829g8Y8wTBS1kC00gJEBjaQCE0JWEhrrzUHNm\nxnT6R8Q+cSLinIjIobKqRHxraelWVWRmZGbE2Wd/37f3ZpxzVKhQoUKFCkUwNvsEKlSoUKHC9kAV\nMCpUqFChQilUAaNChQoVKpRCFTAqVKhQoUIpVAGjQoUKFSqUQhUwKlSoUKFCKVQBo0KFChUqlEIV\nMCpUqFChQilUAaNChQoVKpSCtdknMEjs3r2bX3jhhZt9GhUqVKiwrXDvvfee4pzvKTruWRUwLrzw\nQtxzzz2bfRoVKlSosK3AGHuyzHEVJVWhQoUKFUqhChgVKlSoUKEUqoBRoUKFChVKoQoYFSpUqFCh\nFKqAUaFChQoVSuFZ5ZKqUKFChWcT5uYXcM0Nj+DIYgvnTTVx8MB+zM5Mb9r5VAHjWY6tdsFVqFCh\nHObmF/D+ax9Cy/UBAAuLLbz/2ocAYNPu4YqSehaDLriFxRY44gtubn5hs0+tQoUKBbjmhkdEsCC0\nXB/X3PDIJp1RFTCe1diKF1yFChXK4chiq6vfDwMVJfUsxla84CpUqJAPopG55u/nTTWHej4yqgzj\nWQzdhbWZF1yFChX0kGlkFZq2iYMH9g/5rGJUAWObYW5+AVd86BZcdPUXccWHbsnVIw4e2I+GnfyK\nN/uCq1Chgh4qGpkwPdXEB9/04solVaEcunVNzM5MY7nt4vc+/y0A4QVXuaQqVNi6yKOL77j61UM8\nEzWqgLGNkCdi64LAKy/eDQB4xfN34X++/Uc3/BwrG2+FCr3jvKmmko4aqZmbcDZZVJTUNkIvIvZq\nxwMAeL5OQhscKhtvhQr94eCB/ahbyWWZAXjBOWObc0IpVAFjG6EXEXu1HQYMNwg25JxkVDbeChX6\nw+zMNN7xn58nfp6eamLfZAO7x+qbeFYxqoCxjXDwwH7YJkv8rkjEXokyDNff+IBR2XgrVOgfL79o\nJwDg3975Y7jj6lfjnMkG2u7G379lUAWMbYTZmemEHlDGNUEZxjAoqcrGW6FC/3C8MDgQNdWwTLQ1\nzqlhoxK9txmevyfkMg+94VL88hUXFR6/0nYBDCfDOHhgP9732QfgSMGpsvFuPiojwvZCJwoYtShg\nNGsmVlbczTwlgSpgbDMstcILxwvKZQyrgpLa+AxjdmYa9zx5Bv/8jacAbH0b7/fDQroVG9hVyEec\nYYTOqIZtbBlKqgoY2wwUMMoGgBXhkhrOBXfpvkkAwAd+9oV42ysuHMpr9oLvl4W0Fyt2hc1Fxwu/\nr9oWpKQqDWObYZkyjJIBgDQMZwgZBgCxOA2DAusH3y+OrsqIsP3QSWkYdftZFjAYY69njD3CGDvM\nGLta8XfGGPuL6O8PMsZ+pOixjLFrGGPfiY7/HGNsahDnut0hMowuKSlvCLZaAOLCpot+q6LsQtpN\nK5atiMqIsP3gpDSMrURJ9R0wGGMmgI8CuBLApQDewhi7NHXYlQAuif57B4CPlXjsjQBexDn/IQDf\nBfD+fs+1H2yVhaPXDGMYLikAaDnbI8Mos5A+GwoRDx7Yj6adrBKujAhbG+kMo7GFMoxBaBgvB3CY\nc/44ADDGPg3gjQC+LR3zRgCf5JxzAN9gjE0xxvYBuFD3WM75l6THfwPA/z6Ac+0Jg+S7+xVauxW9\nVwQlNZwFfLtQUgcP7Mf7/teDYjcHZBfSZwP/PzszDc453v2ZBwDojQjfDwaAQWAYn5NwSZmRS8o2\n4QUcnh/AMjdXRRhEwJgG8LT08zMALi9xzHTJxwLArwD4177PtEcMauEYROBZbndXiDds0Xs9yjCc\nLU5Jzc5M47snVvBXtz4GQL2Q9sL/b8WF93UvPBfAA3jenlHc8luvyvz9+8UA0C+G9Tl1PB81ywBj\nYZEudZxuewHGNjlgbHnRmzH22wA8AJ/S/P0djLF7GGP3nDx5ckPOYVDCYb9CK+dcckmVpKQ64fEB\nB/ySWUk/aIsMYzgUWD+4/KJdAIC3/6eLcMfVr87c9N3y/1uVwlosuGa+XwwA/WJYn5PjBYl+Uo2I\nUtwKtNQgAsYCgOdKPz8n+l2ZY3Ifyxj7ZQA/A+AXIjorA87533DOL+OcX7Znz55e30MuBiUc9ht4\n1hxfLPplNQnSMIDh0ESkYQyLAusHbpQFUVaURsj/l58nslUX3qV10r3U10zlpCqHYX1OnXTAiOox\nWprrdJgYRMC4G8AljLGLGGM1AG8G8IXUMV8A8EuRW+pHASxxzo/mPZYx9noA7wPws5zz9QGcZ884\neGC/cCwQehEO+w08JHgD5XbwnHOsdjzRf6qs7tEPaMHc6pQUEAc13cCa2Zlp/NHsi8TPRa1YturC\nu9hyAOivmcpJVQ7D+pzCDCM2KtSjTQvVZ2wm+g4YnHMPwG8AuAHAwwA+wzn/FmPsnYyxd0aHXQ/g\ncQCHAfwtgHflPTZ6zEcAjAO4kTF2P2Ps4/2ea6+YnZnGm18WJ0K9Tr7q17GyJAWMMjbZjhfA9Tmm\nRmoA4h31RmK7iN5AfI55O7ef+qHzAAAvmp5Q0lYyturCu1xASVVOqnIY1ufU8YLEBjWmpDb/nhpI\npTfn/HqEQUH+3celf3MAv172sdHvLx7EuQ0Kz90xAgD4hcvPxx/93It7eg5abN77bw/ACzjOm2zg\nfa//gdKBJxEwSmQY5JDaOVLDyZWOtsX5IIXa9nbKMAooKSDOQtY6xbu7gwf2J0RRYGssvIuCklJ/\nJ/RdX33tg2i7AfZNNnBVF9fl9wvo83j3Z+4H58B5Uw2878DgP6eO6ycoqWaBhjFMo8WWF723Co4u\ntQH0LxzPzkyLHf9Nv/XjPVlqgXI7eCra2zFqA1AHmUELtevbpA4DiCmavAyDsrIVSQvSYXZmGh98\n04vFzb53vL7pM5iBcsWeszPT+NHnhSaAz73rik0/562KN/7weTAi99J1//d/2pDPyfHLZxjDNlpU\nAaMkji2HPHS/7h/OOc6u53PKOhC1MNm0S+kRJHjvIEpKsYgPWqiNRe+t75JyIk5Yp2EAse6z1ikO\nGEC48F524Q4AwH/7Lz+8JRbexZLFnq1tYoneTLTdQGwa5Q3cINFx0y6pyFaruE6HbbSomg+WBGUY\n/bbYWG55ktOpu+eiC3TXWK3UDn4lstTuGKWAkV3Eu2mRUSbtjSmpzRfoikCfx7qjDwa0eLbc0KFm\nGkx7bPp5T685AzjL/kHXDVmrde9BfHfbIDvcLKxKG4fFdQfA6MBfw/GDhFZCGYZqYzNso0WVYZTE\nMREw+ts5n17riH93+1zLLReMhRlDNxrGjpGIklIEu0G3yGhtozoMp4ToLQfm1ZJZBm0Ezqx2Co4c\nDshWC+RThdvJ4bZZSASMjcowosI9AtlqVRnGsI0WVcAoAT/gOLES3vz9VkyfkXad3fL8y20P43UL\nNdMolelkKCkvu4iXcX6UTXs559vTJVWCkgLK01L0mDNbLMMA8jcprSrDKIRc17S8QQEjW7gXV3qn\nMWyHW0VJlcCp1Y6gkfoVveVFpNvnWmq5mByxYZkM624JDSNa4HYSJaUIMmnnVj8tMjpeACqv3A67\n1FIuKel9lA0YW42SojoMIBLx6+rjKg2jGETzArH7bNDI2GprYUDoKDY2wrn1r/eDY+OHllUZRgmQ\nfgH0T7UkM4weAkbThm0apTId4ZIqqMOYnZnGORMNNG2zrxYZMrWzHXaplGF0vEAbvPuhpMjcsNlI\nuOtyMtMqYBRDzjA2KmBkMowcSgoI79+aZeBX//PzCmuF+kUVMErg2FK4k27YRt8Zhrzr7FZAX2q5\nmGjYsAxWWsOoWQZGoh1KHh2x2vH6LuySqZ3tsOjIAVt3M8qfWemAET3m9OrWCBiL6258DWiuG5lO\ndPz4s9gqbf23CtYkg8SGuaRSGYZtMhhMX7jHOc9YcTcKVcAoAcownrNjpG9uXs4wup1RsSxlGHk7\nRcJK28V43RItkXXnzjnHWseDF3CoWnZRfQEtOrr6Alpwmra5LTQMeciTjpZye6KkItF7C1BSfsCx\n0vawa0xvrQbCjDAQdGL4j63aTHEzQRlGwzYSVN8gERbuxRs0xljuTIzwvo3boW8kqoBRAseW2qiZ\nBvaO1weqYXS7qBIlZZnlMozVjoexhiUuJB0F5viB2BXrjpmdmcbrX3guAODv/8+XKdNeojQmmta2\ncEnJn7/OKeUmMoxyVmF/C4neK+1wF7x7LBQudN+Lik7cqs0UNxM0LmB6qplwnw0SqmyhYZtac0Z6\nQt9GogoYJXBsuY1zJxuwTKP0aFQd+ha9mzYso6SG0fYwVrdgUfNBzWPkthd5QYzScd2iQzugyaa9\nTSgpKcNw1dlDbxlG+PmcXXcQDKHhYx6IZ6eAoaNBVXTiVm2muJlYbXuwDIa9440NoaSCgMP1eULD\nAICGpR/TWgWMLYajS2HAsA0Gv8/CvTNrjrDJdbMLb7s+Ol6AiaYN22SlAtdKx8N4wxLdanWPkRfC\nvMylqO0H/X2yaW8r0RvIyTB6Eb2DAIyFhXIbxXOXBb2+yDAU1moglWFEC9BWbaa4mViLsvYdo/aG\n1GHQfaPKMNqaYlh6jF1RUlsDx5baOHeiAbOk2JyHM2sOzploAOhO9F6OqIWJblxSbQ9jdVtcSDqX\nlCzk5S30Yl635nlolzrRsOH6gVIP2UqQs6BylFRZlxQXC/RmW2sXRcDQW6uBdIYR/rvqYpvFSifM\n2ieb9oa4pDouzfNOfu4N21TaaoEqw9hS4Jzj2FIb+yYb4UI9gErvveMRPdBF8JH7SHWjYYw3YtFb\nF6DkDCOfksov7CJKaqJpg/PhzN/oB44ft8kYtOh9zkT4HW+2jpHOMHTXjSyoUuZLZoeJhhU9R21L\nNFPcTBDNO9msYbnlDnxT1IkcahlKys6hpHwKMlXA2HScWXPg+AHOnaQMo3eqpeX4aLsB9kYZRjei\nN934Ew2rK5fUWN2CHS2KuoaAspibF4hajpd7TEuipICtX+3teoFYDHWCIr0H02Bd2Wr3joff8Zm1\nzW0PshTVgsSidz6dCCQ3BLMz0/jlV1wIAPjwz7/k+zpYAJGRpG5haiSkXfO6BPQCyjCUlFRRhlFR\nUpsPstTum2yEO/s+ds3UR4oyjG5E7+VWuFiFojcr1D9o2t5YFGAAvei93ilHSRVpGC0pwwD0fPlW\ngesHIrgVUVJTTTtRtKUD5xx+wMV3vNmUlNywEsj57qT330lRjvTzVhjgs9mge4qum0HTUrpsYau4\npKrWIAU4vhwGjHMnm6UL5nQgeoI0jG4cV0sJSiosIOScgzF151GathdSUuSS0mUY5Sip9QJKiv5O\nu/YwvbZz3tXmwpEChq5jLVFSUyN2QuvRgQI5BYyzm61hrLto2iZGa+F3os0Oc4ouKWAMY0ToMIcB\n9YLVjofzd45gKrpullruQE0AsYahoqTyRe8qYGwBJDOM/jSM0yJgEJ/cAyXVtAXF5PocNUsdMCgI\njNfjDEO30JdxSXHOxaKaZ6tlDBiPAsZWr8VwvAB7ooW9pdk9UwDdMVIrVYdBOlGzZmG8bm2JDGMq\n6j8G6HWsdm7ACP/W2eAMgwoFKXhRoSCALRM0VtuhLjg5sjEZBn3WGdHbMrUZHm1qKpfUFsCxpTZM\ng2H3WD3MMPqw1Z6JWkUQv9276J0vYgNxRWqSktLYakv0gOp4cSVwHq3RtE2x09nqtRiuH2CiQZSU\nOnugDcLUSK2U6E1B0jYZdo7VNl30XhTdAfJ1LMoODZZsDQLEgUJn6xwUtkOh4GrHw2gtpqSWBlzt\nraOXGjVTm+F1hphhVAGjAEeX2tg7XodpMFiGAb+PXTM1o6MMo4xwTVhqhf2AbNOI6ypyzoVmYYzV\nbZgGA2P9uaRkUTRPw2jaJmqmmXvcVkGYoRlo2qbWJeVIlFQZDYOyRstg2DGy+QFDblgJ6LNaWX9K\na0+0IOkokUFhqxcK+gHHuuNjrGGJMcuDrrMh+i9buKfPMIYpeleUVAGOLbdw7mSYEVgm62qRT+P0\nmgPLYOJi60b0phsfCBcjIJ/SojbMY/XwK7ZNoy9KSub4dZlDy/XRsM14N7vFMwzHC2CbDCM1E+s5\nLinLYBirW6UyDMpILNPArtFaotPxZmBp3cUFu0birFRHJzohnThWtzLXicgwNpiSOm+qiQVFcNiM\nQkGVlvITP7AXQPgZTW2U6K3LMPI0DE2Q2QhUGUYO5uYXcOfjZzD/1CKu+NAtePzkan+i96qDHaO1\nuJCui+eiTrUAYFtESeV0n412w6Qn2DmC/VqZ7CFxjN5WO1Izxflt9WpvN+rZ06yZaGsyDC/gsE0D\n4w0Lq45X6Lun78Q2GXaO1ja9xTlpGLG1uiA7tAythrHRGcZWKRTUNV289r5nAIT31EjNhGWwgVd7\nxxlGtnDPC7jy/qxcUlsA4UXzoFgAFhZbOL7cFh1dde6kPJxec7BrtCZ24N2I3svtOMOwjfzus4Ak\nelPAsAzt8WUoqVJBxfXRrJmoKyrLt6L7xfED2GbY/j2PkrJNhtG6Bc5Dam60rr9tYkrKwM6xGk6v\nOT1fL4PAYsspTUmN1EzUTENrq03/ftCg6+E9n7kfAd/4YUA66LSUj335MQAhzcsYw9SIPXBKivQj\nVYYBhEE7LW4P0yVVZRgahBdN8gah4NFrx9qz6w52jtZgGcUZQhpLLU/UNxTZZIE4YBAlZRmGNjMI\nhTzSHYopqbwspGHHGQY911Ztk+36AWpmpGHkUFK2aYggUURL0Xu2TIZdozU4XpAItsNExwsLRadG\napJLSi96N2wTdStLXTqiDmPj38fszDQs08DusfqGDwPSQaeZnIzGNI/Ww3tlsmkPvGOt3lZLQ5Sy\n9x7dj5VLahORJ7T1aq09s0YBg0TrLjKMlouJZrT4l3BJrUguKSCkSLSFe44vdJUylFRea5CmbUo2\n3vAxW9X9EmYPBZSUH1JSY9EiUVTtTd+JbRpi0uGZTRqklLBiF1ir2wlKKuWSGmLhXsfz4XjBphom\ndJoJjTqmrH2yaQ98JkZe4R6gDtpVL6ktgDyhrdeAcXq1g12jNRhGOEGru0pvmZIq55KqmYbgQm0z\nn5IqaueRoKR0HU+FS4psteFxW9H94gccAUc0kdDStzf3A9gWw1g9/HzWCmoxKOuzDCaqq09vUnsQ\n2v1SdwD5/NJoOSGdWLOymajQMIZQuEfa22YaJnRayhtesg8AxLUwNVIbvEsqpzUIoC6e7AzRJVUF\nDA0OHtgvtAYC/dyLtdb1Ayy3PeyIdimWqaeI0vADjhVpUS+aoAcAqx1X7ITCx+hboq92POwYpYCh\nW1BKUFLEg1vJDGortsmW0/hcW60fwDYMQUOQ+6zM8+4cjaq9N0n4psVsqilZq3O+u0YU7DOid7SI\n6bqlDhKUwW1mhjE7M40/mn2R+Pm8qQY++KYX49J9kwDirH1qAzrW6gv3SMPIEb2rgLF5mJ2Zxqte\nsAcAwBAKcD/7kvMAdFc/QaAWEbuigGF30cjw03c9BQD485sexRUfugV3f+90eB55GkbbExd2+Hr6\nlujrjo+pZj4lRQtqLUc8bzk+GjWJkoou5IMH9mdS7M1ukx3PEGAlKSnSMAoyjCgomwYT3/VmzfZe\nlDIMxhhsw9AW7gmHmypgDEn0BmIq1Qv4pg6fuvLF+8S/r/21KzA7My2m7dG1MLEBGgZ99unNKmUY\nqn5SoXmDwTA23lhRBYwcjDdsnDfZwBMf+mnccfWr8dILdgLoTfQ+E+0yaddZts3I3PwCfv+6b4uf\nFxZb+Ic7ngSQ77KirpoE29I3LFzteKLVQT/DkdKV3vRcszPT+H9ee4k4bnqqueltsmXvelEdhm2x\nLkTvyCUV2WqBzWtxLjKMkdgskZdhCA0jI3oPx1YLxHNfgM21ZcvUD/WTI7qMDCJTIzZWOl5fHayz\nrxugbhkZV12ehuFGWtwwUNlqc3BkqYV9Em3Si1hNIOGTFpGw42zx81xzwyOZnR3dSHkBZ7mdDBih\nS0rtsHC8QBQi5bmkjKiwS3UM51wsOiqB9RXP3w3gEfzKFRfh995wqfa8h4VuKCnLMDAefZaFordo\nDWJE9JyxaQFjUWonQ+eku2bouwPTNx/ME70HZZuWq+ldPxAL5bAh33MiYHTCRo5ECdM9s9z2xH09\niNdVidexrVZBSSlmgG8UqgwjB0ejwUkEsib2kmFQEzoRMExW6nnyhGFdwJmbX8B9T57FnU+cwRUf\nugVz8wuoaUTv9YhimSqRYYzWrPB5FNSE44e9pkg4BVIT7VzqdLs5FtM0SLgnl5TjBcrvw/M5apKt\ntqxLyjIYGAtpqc1qQLjUcqNmkBQwWE52GKBRy9pqOedSwFB/d4O0Ta8kAsbmUVJyo8XjkZ12teMn\naN64AeHgvt8ww8gGyWaBS2oY+gVQBQwtgoDj6FIb01KGYZZwJ6kwN7+A35n7JgDgbf/jLszNL+TW\nRcjIdWspHk83r1xw+P5rH8LZ9Y7y+FUnnrMRPqcmYHRCF01IbWWPaTvh72SXlPz+6EIf9oyMufkF\nXPGhW3DR1V8UwROIA5cdUVKAmh8mSmqkZoKx8nUYtmlgbn4Bp1Y7+Oy9zyRee1hYWncwXrfEdWvl\n6Fht18eIQvSWg4fOJTVI2/SKREltpvAtU1InRIbhiUwTgND9BumUciJKKo0iW+22yjAYY69njD3C\nGDvMGLta8XfGGPuL6O8PMsZ+pOixjLGfZ4x9izEWMMYuG8R5doPTaw4cL0hkGES1dJNh0AJOF9Wx\n5Tbef+1D6Lheqc63KrdWXbQGyT5ed/M+daaldEnRAjhaD7MHbTdTVxJFVZlKZEttKkRvIL7Qh8lL\n5+18yfJbMw00o1kRqpkYbkRJMcYwWrMKMwy6Nm579CTef+1DIoAMu1hxbn4Bn7nnGSy3PRGsbEvd\nHoZa11N2KH9vMjWja28+SNu0/PluprVWft8nlqMMo+0mqvxFhjHAgNHxfGXAqBMlpfhMOtuJkmKM\nmQA+CuBKAJcCeAtjLE1SXwngkui/dwD4WInHfhPAmwDc1u859oKjS+HFvk+ZYZS/kHUL+GLLK9WX\nanZmGq+79BwAsVvr3ZGArMpQdDdpxwuUVBLdoKM1C7ap11VajoeRmqWt56DCvqZtwjQYzJRGQ9zr\nMBeBvJ0vnVvNYiLdV03dcyOXFIBSDQjpef/hju9tWrGibq5E2/GVmwaiExt2dkNA35dpMC0lNUjb\ntExJbRnReyXOMGRdULQ4H6BTSpctiDqMZwEl9XIAhznnj3POHQCfBvDG1DFvBPBJHuIbAKYYY/vy\nHss5f5hzvmmlwEcWw4vkvMn4ord70DB0C7iukZgK0ztG0LAN4dY68KJ90XNkH6+7SRu2oTyeNIzR\nuhU6t7TFfXGvIVWgosWJLux08GltQoaRt/OVRe9CSir63kfrZmnRm9pIlD2nQUIXKM+uu8pNQ4JO\ntMJpjnSN0057omEpd7fAYG3Ty+3iep9hgLKp8YaF41GGsZKyqstT9wb2ujpKysrb1GyjDAPANICn\npZ+fiX5X5pgyj80FY+wdjLF7GGP3nDx5spuH5iLOMGJKyjT0VJAOugXcLil6A6GoRnwpILu1so8/\neGC/KPIhNG0TP3juhFrDEJQU7S71lFRIN2k0jGiBatbiyvKOipIaYoaRt/ON/e6GOGeVU4p6SQFh\nhlE0dY+ujb3RzJOy5zRI5G1SVNduS/ru0oYF2tFONm04XqCsjZidmcav/vjzxM/92KYT44I3cSY8\n6TUX7BoRGsaak9QwNmKut6MRvW0z7A6h0pG2W4axqeCc/w3n/DLO+WV79uwZ2PMeXWqjZhmi+AqI\nW3J00+Jc12Zg32Sj9ExveRYGEGspqoV7dmYa73ndC8TPdPNetHtUubtfkyipWgElNRpRUqpFvyXt\nUoFQZ1FmGEMMGHntsuUOn2UpqdFSlFT4nf76q56/aa268zYpqk0GaTeUQQJxJkhBn64/XfHejz1v\nNwDgp1+8r6+mgSuJOozNc9RRhnH+zhGhZ6aLYa2ooHOQ/aQ6nq/MFhhjaNjqIUrbTfReAPBc6efn\nRL8rc0yZx24KFhZbOG+ykSigIQ2jm15SszPT+OCbXgwqwqQF/JyJRumCn8V1VwhsQHG32ssv2gUA\n+MTbLhM3b1i0pV8sylJSujbpFBBGpAwjKXqH/x4mzUCfPVEHO0dssfOlhbMmU1LaDCP8vMfqVuHU\nPfr8fuYl4WvThmP3WG1oxYoHD+wXvn1C0zbxnKlm7ndH3WoBKcMgSkoEDH2BI4C+F8/VtifuM0fK\nMHRut40Cve/n7hwBAJxYaYddnVOt7Sebg21xrqOkgPA7VLqkthkldTeASxhjFzHGagDeDOALqWO+\nAOCXIrfUjwJY4pwfLfnYTcHRxRb2TSZ3amV6OKnwhpecBw7gN199cbyAG0bpTCWTYRTMw1gRw5OS\nWYnqeKJYxuqWdgcKxLMudBoGBZ5YwzBSovfwNQwgWWX+u2+4VCzYKg1DVe2dpaSK6jDi9uazM9P4\n1NsvBwB84GdfNLTK9tmZaVx95Q+In2mTcu5kU3nNCTpRnsfuJ00KFDB0xXtEdfW7eK60PdHll76j\nzWiPT4Hx/ChgPH2mBdfnCdF7bn4BJ1bauPa+hYEFsbxsITfD2C6UFOfcA/AbAG4A8DCAz3DOv8UY\neydj7J3RYdcDeBzAYQB/C+BdeY8FAMbYzzHGngHwYwC+yBi7od9z7QZHl9oJ/QKItYNuC/fOrjvg\nHNg1FvPa3Yx7XWq5YpdMjwX0mQ6l9XLzQV3AWOuEFdwN28gd47ruhLurshpGuuvpZmgYBDoPos3k\n87BNFvfpUdhqvTQlpThG9Vp0rewdD6+hEyvDHdV6WdTG5q9/8aWJLFN1zQk6UaVhRAsnTXvUjwkN\n33e/fP48HMhcAAAgAElEQVRK2xVZGZ3DIOo8us1QKMOggPHYyVUA8T1FQWzQtum8DKNuG2oNww/E\nDJqNxkBag3DOr0cYFOTffVz6Nwfw62UfG/3+cwA+N4jz6xaeH+D4cjvhkALihbrbwr3TqbYgQLig\nlA08Og1DRx+lm6TR66kCzFqkTTDGwtYRiuf0A462G6Bpm1grYaulc5T5bvr7Zlgl6TXlxU7WMEZE\nHYa+sRsQdild63i5E/TkiXsAsGPEhm0ynNC4pjYKp1bD19s9Jmlwuu9OyjDSNTTE5dP1pyveExlG\nvwGj4+H8XeEiTefab52HzmYMQJv1yRoGADx+cg1AfE/lBbF+MsncDMNSN8l0vEBMudxobHvReyNw\nYqWDgEORYXRfuAfE8xB2STdv2fbmjhdEA47kgJEfuIiSmpApKY32sCbxsjpKStYntKK3mxS9011t\nyZK5GRkGBS75Bhd1GAW2WprpDcR9tPK6tpKRgb4jxhj2jNVF8dewQJuUXaNSVquZ6050YrMmV+mr\nRW8dJUXH99OMLwg4VjueOGcK6v3WefSSoVBmde5kA5bBRIZBAWOjZryEhXvq/lkNXYaxzUTvZx3o\nS09fkDEV1N0Nobp58ybgyVhKNZADwkXINJj2PIiSSrY3D4NBmOzFWOv4GKnnD1kSLpq6pRywA8jC\naXhJpR1XIsPYhIChGjHqCkrKiLqDZkVvqkcQlFQUWPKcUp4fRLMn4gxkz3h96JSUapOi2zSoNIxO\nmpKKpj3qZmLI18RygTFAhzXHA+dxJk7Pmed2K4NeFveOF8Bg4YZi73gdj59KBoyNmvHSk4axzUTv\nZx2OLGWL9gAUTi3TgbqVJjIMo1x786XIdTI5kuyGqdstAmGGMVozhdsE0Lc1WXPi6lVtwIiE8RFb\nX4fRcjw0bVMslOlMhBaezaCk3BxKyo5aSas61sqtygFILc71dk8v4OI6IewZb2gL+TYKp1cd1Cwj\n2eLe0GSQTjZgOKmMUGgYmoAvXxO9NuMjQ8HOlIZBbrdmtBk5d6LRleOsl8WdmgAyxrB3ooFnzobB\nhTZh/Qax/Nft0iW1nUTvZyOOLmaL9oBys7RVOL3aAWMQ7o/wucq1N1dlGAAt7pr5Fim/uHzu6ces\ndUINI3xOnU8/SUmphtuQi4qQpqRoUVJVGm80HCUlFdtqgfC9pWkLCuh0DAmeeU4pWSQn7J2oDz1g\nnFp1sHu0lsh0dLZpQSfWJFutlpLSD5oi9OqUIiqVNlby9TM7M40rLg5rPT719su70gkOHtifWVCL\nFveO64v+TedM1EGJOQVgCmI0iXEQM178gMMLeA4llRMwqgxj83B0qY2xupXQAACpcK9rDcPBjpFa\nYseflyHIINfJVCpgWGYOJdVxE5ZaQNI9Uo9Z7fjioteLojElJYoGU8/TcoLEjitdNd7exAwjpqTU\nLikgXCzTlJSbOqZMi3MvCERGQtg7XsfpNWeoNSin1zoJVx4Qfb+Ka7fleGAsLLasmeF3mHZJFQWM\nRIbRZ8CIKank5yW0KM3sEh1mZ6bx5pfH5V5lFve2G4h2HOdMxBtHeSM2OzONX7j8AjRts69iRYIY\ntZrnkkpRUkEUZKoBSpuII4utRJdagtkjJXV61UlUjAPlJ+7JYzYTj9cMRALCG288lWGIhT61wye7\nLB2jFsbjDENuXS7XMLVdP1EsFrbJjm9skWH4YXYyjHGSBHpPadHbNmOtIaSkPOXjKDsrM3XP9bkw\nRxD2jIcL96nVTqa2Z6NwetVJUKBAtr8XgYYnMcZgW1Q0p3ZJ6QR/OTPt1SlF2hsFjPRr0c+6YVd5\neNF0OI/7N19zCd7zky8oODoSn0WGEa8F4/XsRmxQG4F4nneehpF877LbbxioMgwFjiy1lPxmTEl1\nr2GkJ3LZOaK1jPSYTfH4nCK75baXyTB0tRtlXFLrEsctMpXUzaympOQ6jPj4Xmai94MO0SuyhpEa\na9msWYKaIbgpSmqsTIYh2XAJVIsxTFrq9GonYbIAoC0WFdP2AKVLymBxsNSL3vFn1y8lpcswHBEw\nuhfVacOiO/80ZC1hbxTwqV5Jho6i7QVFGUbD0gcMXZAZNKqAkcLc/AK+fWQZX/nuyUyBTyx6d7fg\nnVrrYHeKHgj55BIZRouK8BSUlDbDcPUZRuoxqx1PuH90dRhESY3WLVEglH6edcfDiC0XCjLlPAxg\n+E4ptYaR5H1HbDNTuCcoKStuDQIUuKQCrqSkAAzNWss5x6k1J1GDAUA7/KrlxKNQM6K3H4q/9ZwR\noUDynui1eI8C8WTTzrTHB/rLMOi7V1mnVZAn3+2NMoyxupWpv1GNI+4V9P70GYaRMR0UBZlBowoY\nEubmF3D1tQ+CNgvp6s2iCmsdVBlGWdF7ueViomEl9A8gbA+ia1642k521QTUtRtUkEe7R0vTrVam\npHQ3SMsNR3zGr5dtDULnsGkBQ1ponJQ43azluKSMJCWVl2G4fiBatxCIkhpW8d5qx4PjBVlKSkNj\ntqPhWEC88MhZWd02BJ+vrfT2uXBl9dpPKu5QYCtb0BDF2VPAiB5TVv+QBxmdE3UeTm/agGxG1g86\nBYt/0zbhp8YiiIBRuaSGj2tueCSzg5ILfGjh6EbDcP0Ai+tu5ubVVV6nsbjuJBoPEnTZAKDWMOJz\njx9DbS5o56zrVitslzX1+FUAaDu+sD0C4UVPQYVzjpbrCx582LOaY1ttkPhdLRUwMqK3NG4VAG7+\n9jEAwB9+8WFtewnPz2YYlF0OqxZDVfcDhJuUgCNDn9C0PQCoZ0TvkJrJa68NRFScwfpqxrfSDsV3\nsm+nNxax6N0DJdVthuEGIqu693tnAYQbyPT3XlRE2w1iDUPvkgJS9URST7RhoBK9JRQV+IQFWYDf\nBQd/lmow0hmGEQ6qyWszAVAfqVrm97rus64foOX6CpdUdqEnamWkFoveqiAkbLVS64h0YJF5cACJ\n2dCuzxHw0Mt/atVJLARz8wu45oZHcGQx1I0OHtg/8CZ9qsK9tBVxxM7aauWpfHPzC/jtuW+Jv+na\nS3hBIOamEGqWgZ2jtaFpGKqiPUC6BoIAdSP+rlqur6Wk5HqEhm1qx7S6UT+jqRG7D9E7rAkyDJbY\ncBCcfiip6DE6Si2Njhdgomljbn4Bf/jFb4vfp793HUXbC5wSlBQQvodIFqsoqc1EmQIfy2Cl51gA\noaUWgMLiWI7eWkz1kRLnobFIropOtTpKSsowxLQ9M35OTeuImmXAMg0trZRXh0EL8Xj0PmjOwbC6\nkKp6Sclty4GQbsujpMq2lwjnZ2Q3AHvH60OjpE5FGUZaNxPXXOo7brmB1AMseZ3I8xkatqnNMNwg\ndIdNNu2+bLVkZbdNI2OsoAxjrQ8NQ0eppUGUVPi961kHMT9kADRrkYZRV2QYRTTWoFEFDAkHD+zP\nVOmmC3woMygLVeNBQHJcFaSySy1XTUkZatF7VdF4MPF6gRww4uFJQEhJOX6QaR+y7sQct25HFVJS\nye64AQ91EnKmTERBjC7yQXQhLQOd6C2n8Y0CSqpsewkvCDLXEEDtQYaUYaxmOwsAMS2p/u7iTYPB\nkpXetIDVrWwdAMH1AtRMFmYYPQaM1Y4b06M5GUZPlJTTJSUVve+i7722ARlGXmsQQN9EcxioAoaE\n2Zlp/PBzJ2EwgEFd4GMp3Bt5IHog7VgRY1YL6K2ldV2GoaaklttqV1WcGWQpKbkOA8i2D1l3fCmo\nZKktznk0wjW+nOSup3STTqQ0jI1q4JaGo9AwOilb7YhtwfGDRBCWKamy7SVcn4vgLGPPeB0nl4el\nYYTXXMbKbWW/OwBYdz2xIQCSi7VsL9VVGgPkDjMw2az17JKStbe0aYJzLjj+vlxSJR/bjmi6ou9d\nRfX2ijjD0GgYVtapJmisSvTeHJiGgR85fwee+NBPK6s3rS5mcQNyhpH2xBcXAXLOM7MwCGHVbjbY\nxJ1q1bbaRIbhxMOTAH37EFkUVWkYrh826UtoGFKLiXaq+Isu8o1q4JaGnGFQ9pSx1So61sqUVNne\nQWGzQhUl1cDJ1U4me9sInF5zMN6wMgtP3KkgRbE4SYebrD913NhempdhUBv4yaaN5Zbb0/tckVra\nhL3I4ufwAi7ci91WesuPKU9JhYGy6HvPG5fcLYTobauXZboHZVqw0jA2GQuLLUzv0C9YplGuLTnh\nzJoDg6laexT3pVpzfHgBV2cYmkxHNW0PkGdoKETvepK/TlMBCUpKoYXIIz4JNUnrEBlGIxkwwpux\nux4/vSDZBDEW4mspSgpILkYyJUW9g2iXp2sv4flBptIbCDUM1+d9Dxgqg1Or2bofQNoQeMnrt502\nLFhmolutrGHoRrR6EcU3NWLD8YPS1I+M1U5ccJruRSZ/h0VDrFRY71bDcAPUbVN879NTTSXrQPdD\nXsv7siiyyFYuqS0Gzw9wdKmN5+QEjLJtyQmn1zrYOVrPtMLQCZAydFXegL7wTzVtD4gzGjkYpPUO\nWhjS76+VCBjZLCQ9bS95XCD+Ti2y6SKfnZmG5wd472cfBBDejBvhkkrXgzRsE64fJLKwEVufYdB3\nNTszjdsePYk7Hz+DO65+tea11BmGXIuxYzTrehskVK1oAHU/MbI8ywGjbhkZWy0QunT0LqmwroU2\nRovrrnDflYVccFpL2WrlBbkXSqrdhYZB9Be979mZae01uRF1GLoMI66FUdRhVBnG8HF8pQM/4Jie\nGtEeY3YxKQ/Q37xlajqoTfSkwlar62EjgoCOkpJej1osiMI9Q01JrTmeuPlVIh/tytM8OB2XzjDk\nm/8nX3guAOA/XbK7rwZueSM4HS8Q1ex0LunWIGKut5O/exuvWyIoqxCK3uoMAxhOLUbYeFB1zWSv\nAUEnJoJ9fG05XiDcOXUrxyXlh00XKRvuRfhelgpO0xqGHDx6oqS6qMMg+qtMu42NsNXqM4zw961K\n9N4aWIh63udRUrqOnzqcXss2gQOkca85lJSutTmdh8qSu1Jgq5UpsFVpzoV8TKbGwvG1Gsbc/AL+\nj7/+OgDgD677tlioZdGbXFJCw5CevyMt4L2iyJ7r+IEQ3Gl3lh4601AEDFpY5XnJ4w0bq9GYVhVU\nhXtA3F5iGLUYYeNBBSVlKOjE1GhdIBK9NRlGXrda2zSEo69b6q3j+XC8QCt6y1RYP6J32826ALPn\nki8+y9gYDaN84V6nIMgMGlXAkLCwuA4gpEZ0CDOM8heHqi0IUC7DoAIoJSWlaSS33HZRs4ys4Knw\ni691QncM0WU6i2DokkrpHF4gFmqyi55Zc8VCLYveGZdUosdUvID3ijx7Luccrh/rQHHX3GSlt6Ck\nnOzuzZboxLGGhYDrFy03CJR88t4htQfxA44z6+EsjDRUi1tLQScmXVJ+yiWl/p68iIqjItOlLtuD\nrKa0t/Ac4uubrlvTYH01HwSK9QbaxOioIRlxHcYAmw9qFv9bHzkOAHjfZx8UWXRRsd+gUQUMCSLD\nyAkYlmZqmQ56AbJ43Gt+hqEXvdMOqeTrhec+N7+Af7nrKaw7vrj4ylBSsoaRt1DL9ltaZOi85OAg\n5mT0kWHk2XPpteKJcTTIKdlL6q4nzgAA3vqJO8XnoaSkCoYoeX524h4Q0n4jNXPDGxCeXXfAebZQ\nFJCdcvH3KwKGpkq/48aZWMPSi97pDKNbSiqjp6Xa49MiP9W0u84wSKchuquI0ioqoJNRs9RZeS/o\neAEYg1IDm5tfwB9f/x3xM2XR9z0ZXrfDEr2rgCFhYbGF3WO1xG4rjW5stY4XYKXtKTOMMqL3Yq7o\nre5Ftdr2MkV74evFu0vKDOjGS198ZSmpvIU6UYfhJDMMVRfbfgJGnj1XjBglSopmi/uB6EI7N7+A\nj9x6WDyOPo/7n1oEkKSk6LPV6Ri6OgyAqr03VsPQFe0BEg2q0AMauZSUGR2TU7gXzQGRRe9ukKZS\n0632RcAYsbvWMEinIbNBkY6xWZQUjVpVtQrS9bm76eETACoNY1PwzNlWbnYB5A8uSkM1y1t+HqA4\nw7BNlvGB551H6DRRVYYXZwZz9x+JjkmKjV7AM5SU6we5C7VMb9GuPm2rBeKbsx9KKs8rHwcMS7xP\ngCqTw8dcc8MjGZqi5fq49ZHwZpQzBnoPtMCloav0nptfwJHFNq578Ki2ceEgQEV76caDgOySkjMM\n6ieWstVGFf9he/Oo0juncC+sa2FRR2PWdXsQKjgdK9AwdozUsC7V05QBBZiyAYPeY3r2hQqD1TD0\n87x1m7Pl6DqsAsYmYOFsfg0GEC4eZTOMU+Lmzdnt5bqkXEw2a8odh62p9FZ1qgXimQ5eTmZAfa9c\nhZOqmaKkHD/IXajlmg7a1Y+rKKkBZBjp+cp7x+vCK0/vRVBSkmZCn0nRzSin+7SgaQOGQvSmjI7e\n90b1zAKAU2vURyrPJSVnGOG/ExqGyeB6Qcbm2bAMdDy1aOxFvaQY661j7aooOJU0DC+5cQHCDMMP\neFd1DxQgdo4kdSwdeskwBtNLytcK3rrN2XjdAmNQblI2AlXAiMA5D4v2ijIMzUKtwhlN40GgnOi9\n3HIx2VR72cMRr2oNQxUw5D5CuouPFhl5QVlPWWbjUa9cLNR0jvsmG2Khlm+kthegYRvKRm20gPdb\n+DQ7M42fiyy5n3jby4Q9N01JiQxDEr3zbkbTYIlZJEUaRrpHFTC8nlmAlGHkXHNK0TtNSflBxoFD\ni5nqu5JtypPN7jvWZimpZC+pmJKKsoQuaCl6jzvE6NeCgEGidxkNQ2ygBtMaRCd46zZnL7twh5bG\n2ghUASPCqVUHHS8oDBi2ZqFWgfpIKV1SJUTvxZYjbpDMeUTie3q3p6WkpIxGd/H9wuUXAEhmAOmA\nQQuoXHz3qz/+fADAre99lVio63IdRtTczjBYOOdgwBkGQbQxl1sn+MnGh23Xh+cHCHgc/HQV5y+7\ncEdm51akYYQ77eRjhtUzCwg1DFVnAUA9u0FdpR/u7tM2TzpGVbznBXH336mRWtdDlOjzlItIVXUY\nO6IsYb2AVpJBWTJl+pRV6VBUQCfDHrDorXtN2pzRtUUV5+fvGh0aHQVUAUNgIbp5n7NDX7QHhAtm\nmcFHc/MLOPSFcH7CW//uzgz9YJex1WpamwP6+eIrHbXozRiLhjYF4uKjdY0uvtf+4DkA1JSUXLWb\ndmi1HR+MJXdkiQxDmreQblutCxh5hXg60I2enHsRvpdJUYfhi/dHN1r4efwQaJmnz+OC3aOZHd94\njobBeSiupkXvYfXMAvSdBQB1P7G2pugyrJ9JuoXEPAbFDp0qvQH0Rkl1UrbaVKv9bIZR3lpL14N4\n7ABFb+EGHNCGJ6+eYnZmGi+9YAdeftFOUeSatodvNKqAEaFM0R6gr3+QQZz1Uiu8qI8utTOcdakM\nY13deBBQV+0GAcdqR22rpdekm3B2Zhp1y8T/9cqLxMUn6xyEdIZBr53OQpq2mUiL5QpYufVEum01\nzSiW26r3OidDtoKK35GtVgoYor5CutFmZ6axd6KO/3LZc8Xn4fk84ZAC5Awju2BR8LZTi3XZxoX9\nYm5+AZ+7bwGnVjvKIKvSzYRGpaCk6HOK25vrx7RSpTcQZje9uKTqliGCuG2GYwRIL5Q1DCCe5VIG\nlFHsLO2SKk9J0VC1wWUY+UFqvGElrr30ILCNRhUwIjxzNiraKyF6F1FSZTjrMqMdl9ZdsdCloeoL\ntOZ44Fw9ezh8TJzm+0HoTZdbiBDH7SREUcUONOVgSfciomPC5wrrMOhGkD3+QMwXy6/bK+evpKSi\n343WQj2i5fpx2/KUOD1SsxJUh+tnHU+mETqBVBoGBe90hkEZ3blRtfdk01Y2LuwHFGQpAKuCrLpw\nTyV6m6Ho7SZ32vLEtzTkne5kDzMxVjpJ7S1N9cguKaC7am8KihQwihoQpt93Hhhj0QZqEIV7fmGb\n8rDTQPzZpjsWbDSqgBFhYbGF8YYlXBo66OofZJThrM0CW63nB1jpeMoaDEDdHl3XFoQgBwzVoCXV\nrIs1JSVlJDqeynUa6eciSoo0gnSGIQuotLjnfX65PaMUcy/k5mwNy0DLCcTv0uL0SM3EuhQIZJpF\nRrjLyy6IFLxVhVezM9P42tWvhm0yvOXl5w+8wWKZIJu+ZubmF/BXXw7rT17zp1+O27pYDB0/yOy0\nGzkZhuwOm2zaWGl7pZt0xpmRI77TeMMRZaCpDIPswGUgXFJlA0YXGgaQ3UD1ijwNgzBWV2QYFSU1\nfCyUqMEAogyjYDdRdtQroNcwyNJZqGGous/qKCnp3NPDkwBkKKm5+QX8ztw3AQC/9D/uTCwohRmG\nXIchz4xOZRjpOduA/vObbNq5VJVKw4iHIBlo1sLmefLvZKTHtKZnZhDG6lZ+hqGxOBoGw77JptDL\nBokymxSZJqSMhBafhcWYNq1H31E7o2GoXVKcc3hBHFyJRl3WWI9l0HnQok7f6cNHlwEkZ4sDvWUY\ndD3sKG2rLU9JAeEGYVCmjaLFf7xhYbUd9zLblpQUY+z1jLFHGGOHGWNXK/7OGGN/Ef39QcbYjxQ9\nljG2kzF2I2Ps0ej/OwZxrjosLLYKBW8gXKiL6jDKcNaqNg2EufkFXPnfbwMA/MXNjyq5e1UR1opm\n2p78mm5OwJBtl3QjExd9fLkjFhSVhjGSyjBsqao4o2HoAkb0nLrPj7Es/yzvojupxUX+d80M+2u1\nnThgpLOHZs1K9ClSUVJA+PkqNYzoeXWV3kAoqC9E9OcgUWaTUqZ485obHhELEH0Wog5DUFLJx8lz\nQ4BYXC5DS+nO4+aoglnummsaTGTE611oGBRcdmyA6A1kiwx7RakMo2HBC3iynmg7ZRiMMRPARwFc\nCeBSAG9hjF2aOuxKAJdE/70DwMdKPPZqADdzzi8BcHP084Zgbn4BjxxbwU0PHy905JQZ0To7M40/\n/rkXiZ+Vo17NrMBM5/L+ax/C8ajn0Nl1Vyn4xnUc8eOXCykppqCkktoEEOoORX2i0h1PG6kFnsRA\nJ51hpEVvBX1EnL+ZshDqhFTaRTuKDCNu/8xEhtHRUFKjmQwjj5LKBgwK3ipKijC9Y2MyjIMH9md2\nxJl59NI1l5eRUMCg9xhP3FNTUum5IZOiPUixtVZ3HlQpTtRnx/NRMw2xwemmASFdx6N1C7WcqYEE\nen/lM4zs/PFeQO8xD8KlF+kY2zHDeDmAw5zzxznnDoBPA3hj6pg3AvgkD/ENAFOMsX0Fj30jgH+M\n/v2PAGYHcK4ZzM0v4OprHwTt04scOWU0DAB47aWhRfX9V/6Acs6DrWn0V1bwVTledONZ48fENSTk\nMhmtqSmpoj5RaS9/OsNgjIX0kx9SGwlbraIOA0haa2dnprFrtIbpqab4/Ip20dSsThbSXZFhmGjY\nBlpObKvNLLAKSirtkgIiWkBJSUUZhmIeBmF6qokTK52BUBgyZmem8dbLzweQP48eCANbblsXkwJG\nuCjFE/fIVps895iKC//+wDNhD66f+6uvFW7AdOdB9JGsYdRtI55b0kUdhmz7bua0NyFQAZ3KmqxC\nWDMymG61RVnNeMqlJ7duGQYG8UrTAJ6Wfn4m+l2ZY/Ieew7n/Gj072MAzhnAuWaga+qlc+SEttri\nm/3UKrVoyFbcAoCpsdWWLfJSeeq7oaRWCyipvAUlXYexLnWzlVETXLgvFpuMS0rzb/pZ/l0R1Re3\n5M7aamtisQi0lNSogpJKW2QBEh4VordwSeVnGJwDx5YG34Tw/F2jAIA7f/s1yk1K6OgJv7u8z7IW\nLVorHcowUhpGasEVNmXLwNz8Aj725cfE34o2YAcP7M/0bGraJmZnzgufOzWXo24ZITXZZaU32b6b\ntlmsYbjdLcK1VG1RryhDSYlOAxQwKtE7Cx4qPMoQzhh7B2PsHsbYPSdPnuz6ubutwrVKFu5Ri4bd\n4+qAIUTv1HOVLfJSieZ0EakK94AkJbWmcEnF/Z/01eBhn6i0DhFkKCkgXECUdRgFGob8N7mNA1FV\ndDOfN9VI7KKVlJSgnxgathnaaqXfyUiL3p6WkrLFZy2DgndehvGc6Ht8ZnHwOsaRxRZqpoHdisaD\nBNrw0GdJi42ckdDnkqGkdBkGucMMpm3kqNuAzc5M492vfYH4mc7jlRfvAZDUMGpW2AJjxDa7qsOQ\nNbaGbZSqwyjrkAKyJpBeUWbxT9cBbUdKagHAc6WfnxP9rswxeY89HtFWiP5/QvXinPO/4Zxfxjm/\nbM+ePV2ffLdVuGZJSiqv8SCgLrwDyhd5yU0ACSttT9QJ6F7TS2UYsqOKdqDygkLnIi8o6bYNYYaR\nfc1khhFTUrJnXV585EASzlUOMtTN7Mw0fvi5UwCAG9/944lddBwwNLbaiI7oSDtiGSM1Cx0viIvF\nNJTUWN3CmuNnzA9eiQyDrisqFB0kjiy1ce5kI5dKsVPFm8/bM4rXXXpOIiOhBWi1XS7DIJ3BNo2e\n2qC8eHoSAPA/3355XESaqhmR26yP1K2ubbV07rRpyEOnBDUkY3AaRpkMI2QPqBZjO9Zh3A3gEsbY\nRYyxGoA3A/hC6pgvAPilyC31owCWIrop77FfAPC26N9vA/D5AZxrBt1W4dolKamTESW1R5NhkCic\nfq5051UVFw2oA85K28VY3dI2IrOMuI+Tqg4jPCYOBrMz07ji4l24dN9EYkFRaRiqGSK2xbDm+Ah4\nXBhWt5KDcdqun5jiR5DppXS/LJUbSv452UsqTUn5kq6RtdUCsaCqo6R0DQi9EqL3vqmweG8jhO+j\niy3sm2zkHpPWkNYdP0FLAnGAINqNFjH6fUb0pszKZD21QaEphHvH43OXJzYCycZ86UywCG0pw23W\nymkY3VBSg3BJeX64UaGW+zrQtUcGF9fbZi4pzrkH4DcA3ADgYQCf4Zx/izH2TsbYO6PDrgfwOIDD\nAP4WwLvyHhs95kMAfpIx9iiA10Y/Dxy0QE9PNbVioQzTYAh42IYjD6dW9I0HCbahng8+OzONA5ee\nizJCh+4AACAASURBVOfsaCq5aEDtstJ1qhWvZ8bBbq0TZiPpGyM9uGat44vglTwmfJ4gsvipZnbY\npoHlyO1Cr5Mp3HP93DkZQNYYEPPavvr3qgzDNAQdke4lRWiKgBE+r56SUjcgLCN61y0Te8frG5Nh\n9NBtWZUdpl1StFCHXVGzld5x5bzRUxsUGip1zkS8uUq3DZdpoqbdXcCQKalSorfrd7VrT+tyvUC0\nYelWwxhyhqFfXboA5/x6hEFB/t3HpX9zAL9e9rHR708DeM0gzq8IszPTpStvxaS8gKOWk/qfWu1g\nx4idG/1Ng2mzlcWWq63yls8jUYehaTwoP4Z2wWsdD6M1M5ONqOimyVTHXDkFp928KsOomYZYdOKJ\nfSxRJd52A0w0bZxec5RdbIHIbijdFB3hhoqPDwIu3ls6w7BNJgTPXNG7ngwYcn8kGWN1ogWSGUYZ\n0RsIhe8jS4MNGH7AcXylIzIYHWzTSLaT6WQzDNrlrkSbCqorYYwpx7TS92mZBq588T4AEMV401NN\nHDywP/f+Or7cQdM2c7sOyPz+aN3qzlbrJCmpZU2nYUKZnk4y0pusXlA0z5swmtIw8lqibwQGEjC2\nCk6fPo1Dhw5t6Gvc/b0zWDx8Ch/4wJ25weA/HjiCpXUHhw59TXvM6a8exo2HJ+DevTfzt6/f/RRs\n08Ch0zcrH3tsqY3Fu5/CJ1Zuxa17xsLH3Ps0OAcOLd6ifMwDDxzBcsvFoVM34UvfOoazZ9Zx6NDX\nk8/71cfx5QdHcOiBcwEAD339e9g5WsOhp+KYfu83j2FhqYVDK1/GuuNh8fbHccPJvTh261TiuZ64\n6yk4XoDFdQefP3MuHv33CXzjO8fx5PFVHHK+CgD4zh1PoGGbWFxu45Nnb8Zd+yYAhEVfi3c8AQD4\ng9//esKF9ejtT2C57eLD3u3Cheb5ARZvD9tc3POtJg49GZ7vbd89icWFJRw6dBduf/Qknnl6EZ96\nZi8WHz6Ov+RfS/TqOnxiFYsPHsE17lexd7yBw7c/gdZUE4e+e27ifT15eg2L8wv4s/ZtCaqFfv/3\n7dvwpZyd/sJDR3FiuY1Dz/yH9phusdJ2cfr2J/D1E3tx6OtT2uOe+dr34IzXcejb+xBwjqO3Porb\nntwF565d8fmdbWHx3qdx3/318Ho5dJf429nbH8MNj44nrtmjSy0s3v00PrX6ZdyxO3RqTT9xCvc9\neRZve/XFuP/zt+H+HEL5+oeOorXcxgc+8A3xu1OrHSx+40n8w9ItuOOccTxw11Oo2yYOHb0Bj973\nDBw/KP353X/XU2hGj33owaM4tdrR3lcAcC/dRydvLPX8Dz5wBEstF4dO31TqeBVWO+F9dN2pvXjq\nJv33BwCrXzuMzz81ibO378GRWx7FVx+dwqH7utdvewLn/Fnz30tf+lK+0fjb2x7jF1x1HV9qObnH\nvemv7uBv/uuv5x4z8/tf4r/zuYeUf/uJD9/K3/Wpe7WP/ebCIr/gquv4vz90VPzup/77bfxX/v4u\n7WPe9c/38tf86Zc555z/6ifv4a+N/i3jlX9yM3/3p+fFz6/44M38Pf96f+KY9/3bA/zyP7qJc875\nU6fX+AVXXcf/9e6nMs/1v/3VHfylf3Ajv+Cq6/gXHzzCOef8A1/4Fn/h7/2HOObyP7qJv/XvvsEv\nuOo6/um7nhS/f/T4Cr/gquv4BVddxxfOriee92V/GD7nA0+fFb9bXHfE8bMfvV38/nfnHuI//IEb\nOOec/9mXHuEXXHUd/+TXnuAXXHUdP77cSjzvbd89wS+46jp+1xOnxbm9798eyLyve588wy+46jp+\ny8PHE7+/5eHj/IKrruP3PXkm8xgZf3z9t/kl/+/13PeD3OO6wT3fO608pzRe+6df5u/8p3s455yv\ntF1+wVXX8b/5ymOJY+6L3t8r/+Rm8dkRfvSPb+Lv/UzyevjGY6f4BVddx29/9KT43T9qPmMVfv7j\nX+M//7GvJX732Inw+//cfc9wzjm/8s9v4//1H+7mnHP+9n+8mx/4b18pfF7Ca//0y/zX/jl8z+/5\n1/v5Kz54c+7xb/zI7fwXP3Fn6ef/9U/dy3/iw7eWPl6FJ0+F99FnFPdRGpf94Y38qs8+wIMg4Bde\nfR3/8A3f6eu1OeccwD28xBq7LWy1WwlUfewXpKCnVjtaSy0hr/PtUk5rcyBbhzE3v4DvHFvBzd85\noS2WsmRbreNlqAh6XpkaCo9LaRiSjZAcJyqXlG0aIv1vaJoPtj2dhqEu6Av/pmoBIgvpSQ2DPiui\nxZZT3DxhpFaOkqLCyJUMJaWmutKYnmrC8QPhpBsEjiyGOkApSiq6dqnR4khdr2Gk3UJ1y1DYapOt\nQQBg32SYYZWpNzm50sHeieS9knYByhpGt6K37JJq1oyBi96DaD747988AgA4+NkHC4sdxxsWVjoe\n/ICD82Iaa5CoAkaXID7XLWhxfmqlo5yrLCPtNiJwzgs1DLkOg9qJkM1TVywlz/LQzv42UhpGx88U\n5clBhW5cpegt1VzEzQfDRm08cj61XV+chxwA5EU/64bKahjalumST72REnMzvaRs6lMUC4qqxZ80\njIzoHZTUMEQtxuB0jKORJlI0lCnUsWjTkK32B2SXlJcRVBsK0ZiuBfl9k1uLAlkeTiy3Ew4pINm8\nEogWcXJJdalhJFxSpWy1fvcuKa93DWNufgF/duOj4ueiYsfxqGOtXDA5LFQBo0vQQp3XgLDl+Fhz\nfG2VN0Eneq9Gu4epZo7DSvKpl20nUrOShXvphQIIswcKKm40RCfjopF2VFQ1qxO9CXIvqfC5wxSX\nRG8gWVOiyzB4VJ+hO6ZmGpliQHpNOkdqilckens+VzpQ0k4VAn0meS4pIJ65MsgRrUcW2xirl2nP\nbyBdvJn9fsOf/YBnFs66bWYCOF0vtUSGEQaAowXi/mrHw5rjazMMskBTaxAAGOnDJUV1GLRhUSGs\n9O5C9O6zcK/bYsewcNQtLZQPElXA6BJFbcmBuGivKMPQ9aWiJnu5GYbk1ipbLCXXWKx1iikp1bQ9\nOoYyIyqgUrYGseIdZzMVMBw/bvtB/XH0LUOk6usoDU8fQ+c83rAyMzbohqKgRTRZul6imepTpOtW\nO1ILO+emGxDSNZFXhwHEGcYgrbVHStRgANQ8M6KkKMNIu6SkIJG2eTasLKXjKjKMnaM11CyjkJI6\nsZy11MrnoKvDKFr0CZzzRKeBhm2C82xXARllCuhk9Fu4122xI83EkItSh4UqYHQJeaHWIQ4Y+RmG\nrRn3GgeMchlG6XYikgd/teMlOtXK5+SKgJHtN0WvTeMzafylrg6DIFd6A0hMdGvWzLCoMEEraeim\nAp1jomlrMww6h5W2F9UUJBd2yrjWO15mxoMMxphyJoZoDVKi4+hEwxpo8d7RpXapGeE1S6rFcTQZ\nhhwwUjvthm1mK70V2g1jDPsmGzhSEDCoK3OakkpPpHQkq2uzZoFz9eS/NMLCT6Ah1WEAQNvJCxjd\nUVL91mF0W+xIzS87VcDY+lC1FU+jqPGgeC5TLXovtsLH59ZhSN1uy7S2BiIqKQj1gzVFhS+QpKSo\nX08mw5DGZ65rFh16Pfl8gHSGET5/wzYzPabkYNDR6BMdhc4x0bQToqw8OpTOYbnlKrMA+vu6oy/u\nI0w07IyfX9RhlOhyOr1jZOAZxnkFgjeQyjA66gxD/mzSdEdd0R5czMNIUXHnTjRwrICSUhXtyc8l\nF+7FdRjh97RWQsdopTQ2yiLzdIxeWoP0Q0kdPLAf6Usmr9hxrJHUMLZbt9rvK9jdZBhFLimN6C0y\njByXlFzpPTszjV971fPF33TV6lZUYET9knSUVCbDqKULu+LshnbzuuaDBLlbLUCjW2NBPFsBrs4q\nisa7TjQs+AEXAT0hekuUlEooNAyGhm1g3fEkPUK9+I/VrYyG4RU8hjA3v4AnTq7mOtq6Qdv1cXrN\nEc6kPCQ0jDIZRpqSss1EYSQQv2/bSr7v86aahaL3yagjwp5UhmEYYV8zxw/g+QECHi+MtPiX6Vib\ndvHRdagLGJzzqM14d5RUwPN1zTzMzkxjasRG0zZKdZsI53p74h4ZZmuQZ1Xh3jAgZnHnaRgr+Y0H\nCTpbLQ2PmSypYQDAzPnhQML/9Ws/hpdesFP5GKKSaGesqgq3jLg5YJ6GAYQ7S90xQL7o7fhxk7+G\nbSjanvvqf7vqY+KAEX5mbS/AWPScpK/EGYanFQrDFuf6qXwE1RCl2CWlv4HJ0UZZEDliAPQ85/to\nRPuUoqRMQ5wnucF0GwIgu3tt2EbiOwD0Yv++yQaOL7fhB1zY0dM4vtxGwzaUM1zsqG14mnqh77OM\n8N1KbWiKgk2387zl83L9AKZRPjMhtF0fZ9dd/OarL8G7f/IFhceT5kcDqirRewtDDC7KsdWeWu1g\nvGEpd92J59LMB1+KLgTdPG9ApqTC8yiahQHEi9/Suj5g1KzYuSXoJoWGQa9NN63qvdKNRMNrgHSG\nET3WUlFSxRqGyno70bSiv/ni9+kBQCttVxsImjUzMWRJZ1kcUwxRKtN8sKyjrRscjbSQ88qI3ola\nnCjYp7QsGn4FqDWMdIbhKlxSQBgwvICLVv8qnFjpYO94Q9kwk7Jd+v7rImAkm0TmIU1J0XWafg8E\nuo4aXVFSNBagN1rq8ZNr4By4eO9YqePJpXdqLQoYFSW1dVHGVntq1cGeAv0CiJoBalxSIzUzl0c1\nDAaDxZnOcoum7RXXblAGU0RJxVP51PO6adGvW4ZyB0nHNay4Z5XIMCRKqm4bqFmGaDsOpHtJFWce\ndLOKDENyOqUpKTmIpDFas7AmUVKqbrVATAvIKNN8sJf230U40kWGIdfirDseLIMpd6j0+WQzjGwd\nhsolBcTFe3nC9/HlNvZqqFtqhR9nGFF781SBZR5iSiqZZbZ1GQaNZ+0lw+hR+H7s5CoA4Pl7ygUM\nGklwJgrEVcDYwrBSO3sVTq52CgVvgFxLakoqT7+IHx83koszjPxutQBwNtqZ6CgpV1pQgGyGIafg\n6466tTkQ+/nliWq2Rv+oW2apDKPISUU1HfR4x8/aasPzUAcCGtNaREmppu65JWy1vbT/LgIFm3NL\nZBjJWhw/sghnz5e+4/RiRKK3bGlVVXrL55MnfJ9Y6eCcCfV516OsM5thdEFJiTqhZC2OTsMQlFSX\nojeQ7apcFodPrIIx4Hl7RksdTyzC6SrD2PqgXVRehnF6tYNdBTUYQHJxlrG47mY6xKpgS5TWctvN\nHZ4ExAtZnGGoaCSWyTBGbL2G0XL9zN/FcZEIKltu5QxDTv9VlBRt7nVjXFWZB3HhFIzkOgw5sOlu\nshERMPIpqQmlhhFEc070AaOX9t9FOLrUwq7RWiEFCiRrcdY17WEA5FJSgDqgpwMlBcE84fvEckc7\nM4ba6NN3my7ALEVJpTQM+r8+YEQZRpeiN5C/iczD4ZOreO6OkVLfHxBv9ETAqDSMrYtyhXtOuQxD\nI3ovtZzSGYYnNIyw1UfeYlVGw5ApKZHOZ+ZhxDdIKzfDyO7s6abvJDIMA/W06O364vzKiOHpDIPo\nLleq1m5Y2UwnjZG06J3jkkpPBPR8XuiQkuev0DnlOWKKMDe/gGvvW8DpNaeU40quxVmTKqDToGCf\npmZoIZUDhheEBY7pa2/HiI26ZWirvdcdD6sdT5th0LXYSWUY6Yr8PLQcDSWlqeEQNGlXASN83+lq\n7bJ47MRqaf0CiDdFZ1arDGPLQ9RhaDIMxwuw1HK7oKTUGUZeDQYhnG2Q3xsq+XoRJRWJ6qrdpcxx\nr3XUHLcs8umm7QHIaAeARvRW2GrbbuhusgymbROi2uXGfamkDCM6D8s0xLnrA4aJlqxh6CgpxdQ9\nVzNwKY3ZmWnccfWr8TM/tA97Jxp9BYv3X/uQ+ByKehABcS0OELqkijOMrIYBJOthXJ8r+2dR8d5R\njYZxQhTtqe8Vms0SO5ciDcPu3iWVFr0LKaku5mHU+sgw/IDj8VNreH5JOgqQNIwqw9j6UE26k3F6\njWowSlBKOtG7oPFg/PhYA1luuYV9hGixPLueI3pbLNEaRMVxixvECwv36AbOvl7SnQQk9Y+2uDkN\nZeFe3TYEjx3/PhDnoNQwos+A/iZ3qw3PxRSPV2G0bmKtBCVFPLKsY3iBurutDhfvHcPTZ9cLu6fq\n0IvjimpxgPwMgwRmefc6N7+AD38pfO43fOR2EZhcTZNGIBS+dQHjeNQWJN1HimCbRiKLS1OLrT4K\n97Sidw+UVLpRYjd45uw6HC/oKsOga+/UWiV6b3kISkqTYZxaKVflTc+VDjyccyytu5jMaTwoHi/d\n/GUyDLqpyb+ttNWaJfpNWXGW1XID0XZB93oJ7UDKMDpyhmEameruuhW5pxRCd9gzKumSMlgcBNtR\nr6H0CMt0PUgaTduKbLXFlBSQ7Cfl+ryw8aCMi/eOgfPYJdMtenFcWYYRtcXmoYah6AEGyC6p8POi\nbIaKSo8vd0Q2kx8wGtp+UjTLW0dJ1QQllXQu1SwDlsG6yzCocM/KL9xL019l0I+GcfhE+N13EzBG\n7LCX2ZlK9N76IFpHN8eibB8peq70TO+W68Pxg3IZhiRgLreLMwxhq1130bRNpRXWMuKq1XUN3STP\nKmg5nlb0jrUDtehNO2sKDOlMQumeIvts084I4DXLSPjsKZjKNz/tNHVOpjDDkCipHNEbSAYMLxoH\nWxaX7B0HEC8a3aIXx5XcLXi942cccIR6ipLKy2bC2efq971vqoFjUfFeGiLD0IneVrghUnVlLTsT\no+X4MKOqcSCmJbUBQ2gY3buknB5anNN3X9ZSC4SW+rG6JYJ3FTC2MIpE75NRwChTh6HKMOgi2FEi\nYMgaSJhhFFFSUYbRcrXcdaJPlKYFOn0GrlegYVAdhvT3pK02AGPhQpAOGFTfEWYYcqV37IZSdaUl\n+qvjBvG8AGkxiwOGvnCP87h1eTcahpdT0azChbtHYLBQ9OwFBw/sz3bcLXBcxRlyEA7HKtCfKGDk\nZTOOH2gzq3Mnm/ADLlqAyDi50kHNMrQFqlT9L9x0thwwys3EoE61Mq2qqiUhpLOZMqhJ90y3eOzk\nKnaP1XIbjaowLt2/lYaxhVHUrTbuI1WOUkoHnrOiyrucLZcyneW22xUlpepUCyQryHUct7xLzXNJ\nCQ1DkWF0vJBqoKK+TC+pqAFcXfF7gDKMJIVVs+Jix7brK3em8uQ/FShbopkZOtcTBefVTqxh5FEz\nKtQtExfuGsWjPQaM2ZlpvPT8HTAYSvUgAqQBYB5XDsciUCAi8Tcvm9HNDQHCOecA8GMfvDnh4pqb\nX8Anv/4kHC/AK//kVqVQn670proeoHyGse74GbtqMy9gSP3NyqJfSqqb7IIgbw6rgLGFYeZoGHPz\nC/jorYcBAD/5Z7cVWhxVovdSiVkY8ePDlD0IOFY7nrCU6kDBbnE9J8Mw4wxK59NPtwZRtTYHZN98\nfJnRjtWJMgy5KWG6rkJkGIoWIGN1KyN61604w2h7gbiBa1LAKhK9iaKhgKFbCFUahh8U22rTeP7e\nsZ4pKSCk6F524U488aGfxh1Xv7rQcSXahlOGodk4iMK96HPKqx/RzQ2Zm1/AP339SQAAR+zi+p25\nh/D+ax8StJDO3WVHmwVVfydq4VKEtpvd9DRs/WP7qcPopsX53PwCrvjQzbjvqUU8tLDUdQNKynAt\ng8Ho8prrB1XA6BK2pr05iYJU7FbG4qiqw6CiujIBwzLDDGPV8cA5lA3cEuduxoKfnpKSgkFHp2HE\nrUFaihsy/XpyhmFLO9y2NGs544ZyScNIZhiOpFWku9jWLUNkGB033pnKtE16Lkca9F6KMwyN6N3l\nbu/ivWN44tRaT7tTzjkOH1/FC84ZL/0Yet+rbQ8BVw++AuIgS4s01Y9QG/IdI7bIZnR2Yt0kuX+5\n8+lyEyJFhuGLnwnUwqUILcWGJm9Maz+id9leUrRWLEQFjeuOX7hWpEHX3zD1C6AKGF3D1FR692Zx\nNMSoUkLc2rwMJRVmGMvR4lbWVguoHVKANIcgyh5UHDfduCudMFAV1WHIfzcNBtNgcHwf7UjYpmMd\nX5r1rc0w/CgwZLWNWtTTyjYZ2p6vHDBTpGGQZqMb40qoW6F4mhC9g+5EbwC4ZO8YvIDjydPrXT0O\nCLvUrnQ8vOCc8pRGtp9YftGlvHDOzkzj//uNVwIAfut1+0U242rEfp3u4Wsm5aWPT2sYvWQY666f\ncfE1aqa2cK+X1iBxHUY50XsQDSjp/q0CxhaHPLhIRm8Wx2zwKTM8SZxLVOlNi1Zh4Z4kTGoDhhVT\nUmuOp9yB0iJKgUpHSd1x+BQA4C9vOZzgr2khIGGbfsd5TPV13CCqwzATTQnJbltP2W0Tk/UsM6Fh\nyItekYaRnvutO+7z9x+BH3B8/CuPifdWptI7DbJT9kJLfff4CgDgkh4yDHp/+gwjaasl0OIrawBh\noMx+Tjrdw9R0I0gfPwiXVNvJtq5p2kaOS8oHY8VjdmV0W4cxiAaUpGEMU78AqoDRNUyxyCcvjl4s\njioBfWndjbj44h2OHc0EFxlGgYYhz9jWaxipDEM1xtVKLzrZY+bmF/CRWw6Ln2WKzjZZHDDSczLE\nhLVQ9A4zjGSld90yUU8Jl3LPqPBvsoaR3JkCsYMrjTKUFFEK9LXRezu62OqakiLB8/CJla4eBwCP\nHg+DTC+UFL0/nUuqbmUzDEA9T8L11JXeOt3jLZc/t9SEOXkehmmwxGfbLGurVbj4ckXvaEOS12In\njXicbLmAMYgGlLQ5HObwJKAKGF1DWEpTGcbBA/sTtj+g2OJoK9qMlG0LAsSUVm8ZhnqhoGNWOx78\ngGsyjPAzoEFMquCm46+vueER1CxTtK2mQqpswIipp6x7ykAjyjCIwgpnPkcZhh3SVbGtNlu4l9dL\nCoizJ1Udho5SePpsq2tKarRuYXqq2XOGsXusjp0Fg7pk0MJOM1d0dRjCJZV6/3ZUx7AutwbRZBik\ne9CcjpGaiQ++6cV47+v2g/Mwy81zd1HLe7njMBAG7Bu+eQwLi63C/lnrjpcJWvmid3fjWYH4Gikr\neg+iASXZaoc5nhWoAkbXoDkUaQ1jdmYaV1/5A+LnchbHbJuRxZZTSr8A4m6eK53i4UlAcvHT9hCK\nspC87IEC3XIOrZGXdpPA3ZFEb3kSH+ccbTdA3TazGkakVdTtsF6CAre8qNSjx6httWbms5CRzjBU\nKb/uveXVI+gwN7+AU6sdzN1/pOtxrd89sdqVfgHEgSAvw5ibX8C/3PU0AOCtn7gzc07pBTfPTjw7\nM42vvf81eNPMNGqWgZ/5oX249ZET4AD+6b++PNfdJSq9XV9sBoS5xClnLmm7gTrD0A5Q8rtehLvV\nMCiQkkll32Sj6waUY5XovX0gz6GQ8dLzw9Gof/2LLy1lcbQUF1rY2rxkhmGEGkY8PKnAJWWUEL1T\n+oSqcM8wGCyDiddVaRh5aTfNak7baoFwlyYPt69bZtIl5ccaBqBuMkiFWY6KkirqJRW938UcSkr3\n3iyporgMemkeSAgdUitd0VFAsYZB50QFiXILEMJISnDOq/QmvO6F52Jx3cVdT5zBjd8+jj3jdbzk\nOVOF58p56CKi76tbwZgK92Q0aiZajkb0doOuivboPIHu6jBmZ6bxy6+4EIwBd1xVvFakITSMKmBs\nfehGq1LRnq6/v+p5gGSbkaWSw5OAuJdUmfGs4fHFGQbtkMmtlW5tTrBNQyw6KpdUXtpds0Juuu1l\nM4yOl2xnXU9pGB03Luij4+n/ZAUly21sq82K3roMg97LcssFY1BWbuve2+6xWlcZRj9umYXFFtYc\nH5d0mWGkv9+0RlXmnNK2VMcPCrWbH3/BHjRsA1944Ai+8shJvPYH9xbWDwgLcMcTi3i3gvG643Wt\nYXQznhUIrxGDdVeHEZ5bKMj3UkchXFKVhrH1YRlM2RuH2h+UaQtCzwMk24x0o2HYUaX3cttDI+r4\nmnu8WSx6pykpXXM622RxwFBkGPLchzRPTRbatuuLm7Mu8cCix5SiDiPuYpsc5JMQvS1D65KKMwz1\nTUqN7QIeLlgq8ZPeGy229N5GalZX3Wr7ccv0IngDWUoqnWGUOadmNDOE4Pm8cOFq1ky8YO8YPn33\n01hzfHzpW8cLMym6nlc7nnj+bgTjIAipzayGEbqkuMLe25bor24gz5EpizXHR1NzfxVhYpNE797O\n9vscuovjZJcZBn3ZXspWW7avDLUWWWm7hdmF/HpAsheN6hiRYeTUWJDorTtmdmZay03TTG/a8VMA\ncPwgbs8gdavlnIMxho4XYMeIIW5qyj46CdHbxHLbjTUMRbfavBttpGZiue1pO9XSe/ve6TX8+U2P\n4isHXwXLNPCnNz7Sla32vKkmFhQLdBm3DFlqX7C3P0oqnWGUOaembSR26LpKbxlz8wt4+FjsBDu9\n5uD91z4EAFo6hoL6StsT18fBA/sTVeLh+agFY9pMqDIMP+DRcC2WeUy3ojcQ1xF1g/WcSvsibEsN\ngzG2kzF2I2Ps0ej/OzTHvZ4x9ghj7DBj7OqixzPGdjHGbmWMrTLGPtLPOW4EzJwMY7xhle5Dkxa9\n266Pthtom7GlQYFruURrcyDJxxdRUkU+fds0xC5TV7iXd95O1LZaZauVh9jUrWR9Brmh0pSU48U8\nd8M2tLbahxaWAADv+cwDWpGZ3rOOtiLQ90QuNb/LSu9e3TJz8wv485seBQD81F98tSuh3JIyDMaQ\noV/KnFO68Z/r88LP6pobHsmIwkX0m0xJ0XdI2R1tynaO1rSCcXp4EkHuaJxGL6I3kBwLUBZrOb28\nirBdNYyrAdzMOb8EwM3RzwkwxkwAHwVwJYBLAbyFMXZpwePbAH4XwHv7PL8NgW2qZ3GfXNHPJ1bB\nShUB0iK9o2yGYcR1GEVV3gCiedPhv/U9hFKURY6GQeg2YFDW0HbVtlq5n09NoVXUpSaDImBED8X3\nXQAAIABJREFUYnj4OFNpq52bX8CnI/cPoBeZ6T0Xpfv0mVOm5QbF4q+M9LhWALjqyv25AiiJ0kV9\nmHSQMwwVf55HJRIatomW5Fxz/SA3GwN6o9/kNibpivMv/MYVAID3vk7/eVFQy2QYOUOUyLbdLcKa\nke7am7dcfbfgInz10ZMAgBu/fbxrd10/6JeSeiOAV0X//kcAXwZwVeqYlwM4zDl/HAAYY5+OHvdt\n3eM552sAbmeMXdzn+W0IwgxDQUmtdErNwSCkRW/qVFtaw7CMiJIql2EwxmAb4e6+yCW1FFWc52kY\nBF2ltw51y8Ba5MKpp1xLHc9H2411jbg3lB/N0Q4zCZFhRFx0J+GSCjOMNCV1zQ2PZGgD2uXKiw5R\nbEWLIBVKklvM68FWS7TdNxeW8DN/eXvhZiFPlC7jtKFrbimnxb2OSiSM1JKisVeiS28v9FtCw0gt\n4vQ50T2Txtz8Av74+ocBAB+8/mHUTEO8J1F8qBC+yVTRLcKq9O4zjDL3bRryewPiTQOgp/cGhX4z\njHM450ejfx8DcI7imGkAT0s/PxP9ruzjc8EYewdj7B7G2D0nT57s9uE9wTJZZvARELqkusowRIUo\nx9z8At76d3cCAP7/9s49Vo6rvuPf38zO3pd9/Ur8NokDgeIQQqhLAqElzaOBQEkEEiQSatSHSlVo\nAbWA06BUlVJhGlrRB1RKKSWiKCECZxNBaEIetChAXmyw4zgmwXESr+37sH2fe+8+Zk//mDmzZ2bn\neWb27ut8JOvenbuzntndmd/5/b6/x633PRdrxWBohFqjYQ1Pih3Gsv7PwJCUN+0ywsMwdEosvBm6\n5qTtttRhCB4Gn/UNNBu7OSEpo+l51BsMjIn9jzxptXqyDBs+cjYqzMKFR+5h1ANmW8fhzVvGsXo4\nh58fORX6vLRtJfhnFdQnLA4jht4SkooKxcmE3/ixLlbrLTfxYUPHWF7HqYVWg8G9MD7R70y55vLC\nwuZ686SKpPAwaxLCJh6GcfuDh1t6YSXtRSVL5DtDRA8T0XM+/64Tn8eslIPkI6dS7s8Yu4Mxtpsx\ntvvss8+W/e8TkdMIZlBIKoGHwS+IHx+exM37DmDa/vJPL1RjhRlydp76bLkWWYMh7gOE1WE0W6AD\niJyml9S74Ps2q8Q9ISmzNa0WaM4paPaSaoakvJ6EN602aYZN7JCU42HwkFSyeRgiukZ4x7nr8cSR\n06HPS9tWQjw+2fi52PiPMYZaoxGYdcaJE+rywsOjjPkPNFq/Ku/rYUSlBvPvrPemWyiW8OrpMu6T\nKKKU0TDKVTNwQRZGFr2oZIn8xjDGrgr6GxFNENEWxtgJItoCYNLnaSUAO4TH2+1tABBn/65DHFzE\nWaqamK/UE2oY1gVx15OvSoUZ+Gr2TLkaK0sKsAyCrlFgnDYveBj5nBa4cvSb1x2XvK45FysXXZsh\nqYZrTKbXw+ChJ7Fwz5s+61SS1xuueQFxM2x4SCoq88cxGKKHkWI2wSXnrccjL0xicm4ZGwPmXH/2\nmjfhM9951rWyStJWIudKrZb0MPLNOgzT9u7iiP1RoS4v4sCkIZ/XXz+ax6nFVoMRdUN1PAxBwwjq\nD8aPOwpD1+TqMCSunzTZdWlJG5K6H8BN9u83AbjP5zlPATifiHYSUR7ADfZ+cffvOnJ20z+RpEV7\n1uto9r7+cdioFQNv0dGIMQvD2UfXrB4+Ac3VcsKNOyxkwT0RmVWqGOoR52EAHtHbcM+3qJsNmA3m\nFr2FUaxi4R4AzC/XXLHvuKtcfk5RGShOSGqpDsYY6o3k8zBELj1vAwDgiZeDvYy3v24dGIA1I+F9\nmILIxMMwdNRMhprZcK6DdtQDiDqZr4cxlscZH4MR5YU9edQK+934Hz93PIm0LcfzOf9EmDAWA0Yg\nR5FFLypZ0oreewHcQ0R/DOAVAB8BACLaCuDrjLFrGWN1IvokgAcB6AC+wRg7GLa//RpHAYwDyBPR\n9QB+jzH2fMrjzQS/Sm8eL5XxMDaM+a+UolYM4moxroeR0ykwHAW4L9KwG4ozHEkmJCXcXJz25qLB\nEMZkNj0M090yxGj1SETRG7Bu5N4bWZxVriN6R9wEx/I5aGR5GM6NM4WH8dLEPAjAX9xVxN4fvoDP\nXtOaAfSjQxMAgPs/+W6cs2Es8f9hZOBh8PdHvMEmbboYB3Fh4VcYuG4sj19NtDZtDPMkC8USvvbY\nr53t3JMIanceXxuiRBqG2bASNWSMNv9O3P7gYRyfWcLWtSO+35V2kMpgMMZOAbjSZ/txANcKjx8A\n8EDc/e2/nZvm2NqJX0jK8TCSZEnZF9kH37YVdz/5WqxiJPf+zYtofCT6oywUSzgxs4x6g+GyvY/6\nfskMIcsn7IbCL2AZl3rIx8NwaxjNtFpRw2iGqprbLXHbdL0G9z68HkZc+EUcFV7SNMLqYQOzSzVn\nASHrYRSKJdxSOOiEmoJCIg8/P4E3blolZSwAd8diWQ/DqWOoms75tsPDcC0sfBYmG8byOO2z0OLv\n1y33Wk0Ktwk31Mv2PurbRVkn8h3slEQb4v234sCTBmSuHyB5eC8rVrbqo0/gFdYiUxIeBr/Idp+z\nHl/80IVOTULcMIO4ml09FO5h8BgtXwkH5e9r9kQ8AKFtCxwNQ8LD8OvtJDYfFEVvpw5DEMPznjqM\nikfc5q85v1yX6rXDL+I4xmZ8JIe5pZrTjFJ2pR0nJDJbruHJo6dx9a7EyYQOLg9D8mbF359y1XSE\nXtnssDDyMTyMpZrp26r8+ou34XfeeDbesHGVqxFo2BRA7+edJMyTVPTmRa8yoncnUQZDAl1r1TCm\n5isgQrLZBEIdxvUXb8Nbt6/FJTvXx+p0C3g9jHCDkSRG66TehmkYPEtKRvT28TByuuY0cHN6SYlN\nBoW6Cm+3Wme70BoEsEJFMkVYcUNSgFW8N7dcb3oYkiGpKKG2UCzh8i8/BrPBcM9Tx6QLtYjIOcag\nWRhRiHUMNU9xZJYYPqFLkfV2LcbpgFqM6YUKzlrlvh6DPIZta0dw429ZuTmy2lCSwj1ehySjYXQS\nZTAkMPTWkNTUQgXrR/OJLpxmW2Tri3amXI1d5W3tL2oY4V+8JKl4hhNuitY6ZFxqt8HQXNt5Wq1m\nj8kc8gtVGRo0jZDXNVf67JCPhyFzI4sbkgJsg7FUc9q76JI3zjChlnuHZ+xU56mF1pbjScjFWBCE\nMZIXDQYXvdugYQiv6eft8cWZn/ANWMkk3kLaMMGYh/mKt14de9HmHGtO0sOQ/Aw6hTIYEug+onfS\ntiD8dYBmL6kz5RrWjcUTrwFPM8EIg5Ekf5+/bhwNQ6oOw7Vy1F3beUhqKKeDiFyV3t7QkzMoydMz\nqqlhtFYIx8HxMOKGpDIQvf1uZIAV6/70d55NlcHjJc6CIAxxTGs7PQy/5AgRbjD8EkYAYNqn8wLP\nlOOvLXoSE/PLyOe02L3cRJKK3k2DoTyMvscISKtNajCcSu8GA2MMM+X4nWoB9wo4KiSVJBUvjveQ\npg7DL60WsLSJSt09Yc23oI+n4tqjWIOypPxaSsQhbmsQwGpAOLdUTy16e1N+R+1z4F6FH7KFWnEW\nBGHwm5xoMJK2RImDS8PwadcR5mEs14Lroq6/eBved+Fm7Fg/4vIkJucq2Lh6KNE8b+f4EtZhLHLR\nu8c0jN4yb12Cbk+6E5mar+DchJkrPCPJNBtYqNRRbzCsi9lHCmhe+ETAqoiVSpJUvHghKXkPQyzC\nEvPrecFdo0GuIjzAzpISsqesn3ZFd4CHYR1n8ovf6VYbW8NIL3oD7syXy/Y+gvLMcujzZQu1HA1D\nutLbel/KQkjK2yY8CyI1jBAPg2ctejUMzuY1w5iYrTht8wFgYm4ZmwIKJqPIJwxJcaG+1zSM3jra\nLsHwiN6MMamQlNPevMGcVhyJPAx7/1VDuVhTu+Km4jkr0FDROyMNQwxJ2RqG2QhKt22t6HZpGE7h\nXvjKNIqRRCEpA+Wq6Qj1Wa20j0cYizSFWmk9DJ49t1w1nYVTOzwM0WD4eYrjwwZ0jXw9DF4MG9QM\ndMv4MKpmA6cXq9hgP2dyvoLzNyabYCgea5LCPS56Kw1jAPBqGPOVOir1RqIaDMAtevOeOElEb36R\nxmltnuy4orNomims8mm1XNgWX7NaN10tpp2WIbXWAr28Pb7V20vKq4skhd9I44SkeLU3rwfIKr00\nzHtImsHjJU2VPtD0KsvVeksL+SxxVXr7GAxNI6wbNfw9jHnuYfhfk5vXWO/vidmmYU7jYSRtPqhE\n7wEip2suD0OmBgNwi948Vp0sJGXtL9MiOfx1Y3gYKYRTsUmgGC/O59yiN2C917pGnkpvrmHYmoen\nK61oxGTSah9/cRoAcOfPXolsQse1I24wssoWCtKcvvLRtyXO4PGScz7ftGm1DWfh1I4sKSJyPtMg\nLSqoPYgTkgq4JressQzDSdtgLFVNzC/XsXE82TXMyetWe3O/sa9+cIMR1DW6W+mto+0SrMFFzdXE\ntKTB4LHkWsMSvIGkIan2eBg5R9COoWHkk9+QRYPh3e4NO/HfLQ/Dq2Fw0dtT6S2EpJLeyArFEm5/\nqJl9FNWEjr/3vM12VqGZdrZ/aNZhyK1uechvqVpva5aU9bqEqun2GkXWjeZ96zC4wdgQUBfFDcaJ\nOctgTM5bPzeulvcwGLNafsTxMsvVOojkFjSdRBkMCXK6u7150lneHF5EZTYaziopiYfBL/ysPYx8\nnMI9+zkjhkTzQXvfYc/FwkVvBreYPiTUZ/DH/OdCpe7qMWW9rlsXSULYrAFfg+HxMLKseG5X+wf+\nnsh6GESEEUN31WG0o9IbsHWkqhnqYbw42dpPanqhGjouecOqIeQ0wslZK9NsYs66hjdJehhc77Jm\ng0Q/f7FiYiwf3AS0W+kt89YlGLrmZMUAzZBUkml7HN5mhIekkuSA81Vd3OFJSV83SMMoFEv410df\nAmD160laQDYU4mE4s75dKZWap9Jbd36K23n4wtAJXH5IajCSzhrgPbxOOSGp7r+k0noYgBV7Fyu9\nZbSiOIg1N34EhaSmFsJn0+gaYdP4sKNhTNieRhoNA0Ds1Npytd5z+gWgDIYU1ohWt4aR0whrZQp+\nNCu7YqZcxfhwLlEe///+yhofcm+xlOlc31yIhsGrjvlEvlOL8YY9ifA5B96buVO4V2t40m11Tx1G\nM/TEW4MYenPuBRE5xijpDTzpgCIekjq9aC0a0szDWCn45xs0HCsOw4aOctV0QrNp2rqHYcQxGOUq\nGt66qBjjkjeNDzkaBjcYGxNGCTjeuS1RyM7C6DTKYEhgaISaaRXbFYolfPOnR1FvMPz2PzyW+Kat\n65YeYlV5x9cvCsUS/vGhXzmPg5oJypAPyaJJOzcAaKbk+noY9QaW66Y704lrFXYdhqvS206r9a5w\nvW3T45J01kCr6N3dl1ShWMIvX5sBALzn9h9Lf1/4XG/eP6kdojfgn/kmsn4sjwZrjhTmTC9UcNbq\n8Otpy5oRnLQNxdR8RbrKG2heM3FrMSwPo/cUge7+dncpui1s3vsLa7XNMx5kbto528M4k7DK+/YH\nD/u2ac5iri8Xbv1WQFmMh/R2lXW22wagUmu4/sa1Da8nYc3ubjhT+EScwUwJb+BJR4mO5XVo1AxJ\n6V3sYXDvkH9v0iwyRvKWh9EsWGxvSCpMwwBai/f8+kh52bxmGCdnl8EYs1Nq5aq8ATFFPp7BWKyY\n0nUwnaT3TFwXwAW+Lz8UvNqO3+XSEr1nyjVsCKhK9aOdc325gOcX485iPGRQlhQXt6tCWi1/fsWT\nbsufz0NSQQZD5kaWRGwmIoyPGJmn1baDMO8wqbg+YlhzvWv1NmdJ5fjiINxgiLO9q/UGZpdqkQZj\ny5phlKsm5pbrmJirYJNkhhSQ3GCUa6ZUCLvTKA9DAh6nFot+RJLctLnofXoxWafapLH2uBSKJTz8\nvDXV7bp/e7xl9ZnFeEjHYHjCDM3mg2ZrWq3fdkNzWoN4QxayISkZxocNp1K/HRXPWZHlIoPP9eb1\nSG3LkuLhR58RrUCz0PWUMOb41GK8JJTNQi3G5PyydA2GeJzVmC3OyxUleg8MXODjXzgvSW7ahqY5\ndRhrE6TUtmOuLw9Z8FXoidnllpBF0pCNH+EhKRPLtYYnS0p3QlV5lyHRHQPjNQxDhr+w3g7EuHe7\nbpxZkOUiYzRveRjVNmdJ8Rtx0Ov7eRjT87wtSJSGYddizC7ZjQflPQzeSyu+hmH2pIbRe0fcBXAP\n4xO/+3rc9v1DWBa0hKQ3bV0jLFXrWKyaiTyMdhR2xQ1ZpK0PeOjgSQBA4dnjeOroGee48znNqYEY\n8tZhOJ5Ea0M6v8l6wyvpYQjjcbtZ9A6bdZ0UJ0sq5eCoKIZyVqV/UBYWNxjiqNaoKm8Obw/y66lF\nzFfq0im1gISGUa0rDWNQ4KvIq3dtxplyzclW2iZx087pGqZsdzpJlhSQfWFXO3URTqFYwq33H3Qe\ni5XUPN0WQGsdht3GXAw9cZ3Cb+5FGg0jKWKlfTen1Wa5yBgx7CwpswGi9oj9hWIJT758GmbIDPph\nQ8doXncZDKeQNiIkZbUyB/Yfm3Eey5K8DsOUGg3QaZTBkKA5WpXhgq3jAIB9f/4uvP116xK/lqGT\n01okSZV3O8hC0I4irJL6Y5ee42zz8zC8oSduVOaWWwVOb3v0duIyGF3sYQDZLTKahXsMhqZlXrEc\nlNEFtLZo8RbvNVubhxsAQ9dw1qoh7D82C0C+aA9IVodRsxM7eq21OaA0DCm4sFk3G5i0WwrIrk5y\nGjkroiQhqXbQDl3ES5gX42cM+O+ViJCU1zA47dFXwsNwhaS618PIEt4ahKc6Z02Sep/1Y3lXWu30\nfBVjeT3WCn7LmmG8PL0IQL4tCND8nsVpcd6rnWoBZTCkcCblmQyTko0Hm6/VnNSVRPRuB1kI2lGE\nCa/uWd+ih6H7V4A7Iamaj8HgbUJWOiQ1GJfUSD4Hxqy5Du3wqpKER3m1N2dqoRKpX3A2C15FGtE7\niYax1KPjWQEVkpKC3xTMBsPk/DLWjRqBlajRr9VcnXXawwDa1/COEya8ivFfXw2j3nD1zeLPWa61\n1mHwz2NlRO/mMQ2Oh9EMB7bDKMcNjxaKJTxx5DSWaqajc8RpC8LhmVJDOc3lKSbFSFDpzcez9qLo\nPRjLoYzhAl/NDkmlWZmIq7NuMBjtJsyLCQtJORpGwNjOliwpz0zwdsJvNLpGPdd9VBa+OrYMRvbn\nHCc86k0D5zrHkamFyJRaDs+U2jQ+nOqzSyJ6lyvKwxgo+AVieRiVdAU/WrOStRezJmQI8mK8NRbO\ndl1Dg1ntFLxNCf32BVpHvLYTHpLq5gyprBm2v6tzS/W2eBhxMrqCdI5K3YztYZRmygCAV0+XAzOx\n4sC/Z/E0DNvD6MHrXRkMCZxJeY0GpuYrOO/sMenX4nrIIHgXUbi8B6P1d69WEWQ8CsUS7vzpUQDA\nx7/1DG659s1tDbPxkNQgGQze6XZuuda2YsWo8GiQztFg8UYNFIol3PP0Medx1LCsMJJoGFz07sUF\nogpJSSDO4p6azyYk1WnBuxtwid4+c7nnK/XAsBXfzsMUc8vWKm5qvpJZF98gHA+jy1Nqs2TE8TBq\nK5KJ5kdYuvc3f3o08jO//cHDLSEk2QaeTlptjJBUU8PovfX64HzDM4R7GKcWqqiajXQFP5ryMDj5\nIO/BXs0y5vYk/EJSWbRfTwrXMAZF8AYEg7Fc71g7FD+dgzO7VItcKGRZqMo/+zh1GE0NQ3kYAwH/\ncvAvVhoNg7dKXzemPAy/GgsAgUK333NWolrdS1PDGJzLid+ozQbrWDsUMYHCj6iFQpa9tQytGZIq\n2APNdu75ge9gM65h9KLoPTjf8AzhN3me9peuB41lfJLMwuhXgkRv0dsI8kL49nZ18Q1jNK8jp1FX\nNx7MGnFlb3TQUF5/8TY8vucKBL3zYQuFLAtVNY2Q0wgHS7O4ed8BlGaWwOA/c2RxUAv3iGg9Ef2I\niF60f/r2xiCi9xLRYSJ6iYj2RO1PRFcT0TNEdMD+eUWa48waLm46HkaKkFRT9FYehugxiJ1s3R6G\nf0iK77sS1epe+EyMbm48mDXizY7PrOgkMguFrAtVDV3Dz+yaEBGvp1Ou1qFrtCJta7Im7RHvAfAI\nY+x8AI/Yj10QkQ7gqwDeB2AXgBuJaFfE/tMAfp8xdiGAmwB8K+VxZgq/yfN5GKlEbx6SUh6Gp/mg\n6GH4NyX0+30lqtW9FIolzC7V8PL0Yqaz1buZYcFgdEMoTnahwD2Ul/e+H4/vuSLV98TQCQuVuu/f\nRE+nXDUxaug9WbOTNoh2HYDL7d/vBPBjAJ/3POcdAF5ijB0BACK6297v+aD9GWNFYf+DAEaIaIgx\nVkl5vJnAL5DjM0tYPZRLlR6nQlJNeFiJyC0gB6bb+mRJAe2vVhfhWVmmPUgoTWpmL+EKSXWBZ9WO\ndv9Jyec0jNqja72Ink65YvpOs+wF0hqMTYyxE/bvJwFs8nnONgCvCY+PAbgkwf4fBvCLbjEWQPNm\ndmqxmqoGAxBEbxWScm76Qzl399OhgPAUESFvV4GvRIGeH1mOPe0lDF2DoZPVrbZLtJuVXCj4Yega\n3rJ1DPtLs66OzF5PZ7Fa78lOtUCMkBQRPUxEz/n8u058HmOMAYg3n9AHv/2J6AIAXwLw8ZDj+1Mi\nepqInp6ampL97xMh9v5Po18UiiV8++evAAA+/7321gr0AkGzvt0eRsAoVr0zK7ZOZGV1CyMrOHOk\nF8jnNGxdO4I/e8/rnW3DOa0lJLpU7WMPgzF2VdDfiGiCiLYwxk4Q0RYAkz5PKwHYITzebm8DgMD9\niWg7gHsB/AFj7Nchx3cHgDsAYPfu3dIGKwniBSKrX3j74EwvVAYilBEGNwxBnWf9/6b7DlBaKVZi\nhki3MpLXO1qH0W0YuoaayZyRve88bwOOTC+0XM+L1TpGjT71MCK4H5YoDfvnfT7PeQrA+US0k4jy\nAG6w9wvcn4jWAvgBgD2MscdTHmPmZOFhdKLArNvhoQ1v51/Re2jtSrtyg5L86ERWVrfA6wg6Vend\nbRi6hqrZwC9fm8Gm8SFcvWsTJuYqOGknx3DKPexhpP2k9wK4moheBHCV/RhEtJWIHgAAxlgdwCcB\nPAjgEIB7GGMHw/a3n/8GALcS0bP2v40pjzUzxLxz2RqMQQ5lBME1Ce/NP0joFh93ysPoRFZWt8BD\nh8rDsMjrhJrZwP5js7ho+1pctGMtAOCX9ghYzmKldzWMVEfNGDsF4Eqf7ccBXCs8fgDAAwn2vw3A\nbWmOrZ3owgUiW+U9yKGMMIZ0LVzDyHk1jJXrShtEp8XWTsFrMZSGYWHoGqYXKjgyvYgP/+Z2XLB1\nHDmNsP/YDK65YLPzvKUenecNqEpvKcSupLKT9gY5lBFEoVjCYrWOA6VZVz1DqIfB516om9aKo0Rv\nN4au4dCJeQDARdvXYtjQ8cZNq52Z4ZzFqtmTrc0B1d5cipxLw5ALSXVD3ng3wZMA7HIGVz3DB966\nxXleUEiqF6tme50Rx8NQISkAMHKaU49z4fY1AICLdqzFD/YfB2PMSRUvV+sY7cFOtYAyGFK4RO8U\njQcHNZThR1Q9g64RzAYLDEnJjshVyMM9jG6o9O4GuJd73lljTqbURdvX4K4nX8XRU2XsPGsM1XoD\nNZM580R6DfVJS0BkNRobNjSs7tGVQrcRlQTgpNwa3SV6DzIjKzjVsBfI2z21uNgNAG/dbv2+3xa+\nl3jjwR69b6hPWoJCsYR6g2G51sC7v/TYwBfcZUFU87ggrWJI3bQ6Bg9JDdKkwSAKxRIee8EqHH70\nhQnnnvDCCUu/+NTdz+KyvY9iX9Ga8NerGoa6yhLCY+0cv/bFiuREJQH4eRiFYgmPHJoAALz/X36i\nPoMVZkRlSQFoLcKdXarj5n0H8IXCAdxSOOg8rzSzhL0/fAFAb45nBZTBSIwquGsPUfUMXg+DX6S8\n0duJ2WVluFeYUUOJ3kDwPeGuJ15r2V6xR7gOZB3GIKIK7tpHWBJAXtega+TMzR7Upn/dhPIwLIKu\nfZMFdyras28/vlDZ1XPf1cH+pCXoxEQ3hZUFJabOKsPdeRwNY8ANRtC1r4fMu5heqPakRzzYn7QE\nquCuM3hbhijD3XlGVEgKQPA94cZLdrRsF+nFULYyGAkZ5N5BnaJQLOH543M4U645FeDKcHeeAyUr\nVZRnAPXaajkrgu4Jt11/obM9iF7ziImFxNl6jd27d7Onn36604ehyBBvBgpgGYYvfuhCAKpSvlMU\niiV87rv7UTXdg4LU4smfy/Y+6ts7btvaETy+54oOHJEbInqGMbY76nlK9FZ0NWHidtoZzAp5bn/w\nsMtYACrpIIzPXvMm34VPr3nEymAouholbncn6nNJRr/0jlMGQ9HVqDbw3Yn6XJLTD73jlOit6GqU\nuN2dqM9lMFEehqKr6RdXvt9Qn8tgorKkFAqFYsCJmyWlQlIKhUKhiIUyGAqFQqGIhTIYCoVCoYiF\nMhgKhUKhiIUyGAqFQqGIRV9lSRHRFIBXUrzEWQCmMzqcXmEQzxkYzPNW5zw4JD3vcxhjZ0c9qa8M\nRlqI6Ok4qWX9xCCeMzCY563OeXBo13mrkJRCoVAoYqEMhkKhUChioQyGmzs6fQAdYBDPGRjM81bn\nPDi05byVhqFQKBSKWCgPQ6FQKBSxUAYDABG9l4gOE9FLRLSn08fTDohoBxE9RkTPE9FBIvqUvX09\nEf2IiF60f67r9LG2AyLSiahIRN+3H/f1eRPRWiL6LhG9QESHiOid/X7OAEBEn7G/388R0V1ENNyP\n501E3yCiSSJ6TtgWeJ5EdLN9fztMRNfI/r8DbzCISAfwVQDvA7ALwI1EtKuzR9UW6gD+ijG2C8Cl\nAD5hn+ceAI8wxs4H8Ij9uB/5FIBDwuN+P+9/BvA/jLHfAHARrHPv63Mmom0A/hLAbsar6TLXAAAC\n6UlEQVTYWwDoAG5Af573NwG817PN9zzt6/wGABfY+3zNvu8lZuANBoB3AHiJMXaEMVYFcDeA6zp8\nTJnDGDvBGPuF/fs8rBvINljneqf9tDsBXN+ZI2wfRLQdwPsBfF3Y3LfnTURrAPwOgP8EAMZYlTE2\ngz4+Z4EcgBEiygEYBXAcfXjejLH/A3DasznoPK8DcDdjrMIYexnAS7Due4lRBsO6ab4mPD5mb+tb\niOhcABcDeALAJsbYCftPJwFs6tBhtZOvAPgcgIawrZ/PeyeAKQD/ZYfhvk5EY+jvcwZjrATgywBe\nBXACwCxj7CH0+XkLBJ1nZvc4ZTAGDCJaBeB7AD7NGJsT/8aslLm+Spsjog8AmGSMPRP0nD487xyA\ntwP4d8bYxQAW4QnD9OE5w47ZXwfLYG4FMEZEHxOf04/n7Ue7zlMZDKAEYIfweLu9re8gIgOWsfg2\nY2yfvXmCiLbYf98CYLJTx9cmLgPwQSI6CivceAUR/Tf6+7yPATjGGHvCfvxdWAakn88ZAK4C8DJj\nbIoxVgOwD8C70P/nzQk6z8zuccpgAE8BOJ+IdhJRHpY4dH+HjylziIhgxbQPMcb+SfjT/QBusn+/\nCcB9K31s7YQxdjNjbDtj7FxYn+2jjLGPoY/PmzF2EsBrRPQme9OVAJ5HH5+zzasALiWiUfv7fiUs\nra7fz5sTdJ73A7iBiIaIaCeA8wE8KfMfqMI9AER0Law4tw7gG4yxv+/wIWUOEb0bwE8AHEAzlv83\nsHSMewC8Dlan348wxrxiWl9ARJcD+GvG2AeIaAP6+LyJ6G2wRP48gCMA/hDWArFvzxkAiOjvAHwU\nVlZgEcCfAFiFPjtvIroLwOWwutJOAPhbAAUEnCcR3QLgj2C9L59mjP1Q6v9VBkOhUCgUcVAhKYVC\noVDEQhkMhUKhUMRCGQyFQqFQxEIZDIVCoVDEQhkMhUKhUMRCGQyFQqFQxEIZDIVCoVDEQhkMhUKh\nUMTi/wF2GGLOe+6hfwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x2256eb1c2e8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "diffs = np.subtract(type1(n, a, s), type2(n, a, s))\n",
    "print(\"Ave = \", diffs.mean())\n",
    "print(\"Std = \", diffs.std())\n",
    "plt.plot(diffs, 'o-')\n",
    "plt.axhline(diffs.mean(), c='k', lw=0.5)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}