{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "### Forward Euler for single-variable model\n", "We will develop code to use the forward Euler method to solve a single variable ODE, then test it with the \"test\" equation (the degradation model).\n", "\n", "First, make a single-variable degradation model." ] }, { "cell_type": "code", "execution_count": 114, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import scipy.integrate\n", "\n", "def degradation_model( t, y, params ):\n", " dydt = -params[0] * y;\n", " return dydt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Write the method\n", "\n", "The method should take in a function handle, the array of time steps, the initial value of y0 (a single state), and an array of params (though it will have only one entry).\n", "\n", "The outline is this\n", "\n", "- Make a 2D array to hold Y (1 x NT) \n", "- Set the value for the first time step\n", "- for each of the remaining time steps\n", " - compute dydt at current time step\n", " - set value of Y at next time step using Euler method\n", " \n", "The Euler method as we wrote in the notes is\n", "\\begin{equation*}\n", "\\widehat{y}(t + \\Delta t) = \\widehat{y}(t) + \\Delta t \\cdot f(t,\\widehat{y}(t))\n", "\\end{equation*}\n", "\n", "and since this is the code version, we really should add a vector of parameters to the function, i.e.\n", "\n", "\\begin{equation*}\n", "\\widehat{y}(t + \\Delta t) = \\widehat{y}(t) + \\Delta t \\cdot f(t,\\widehat{y}(t),\\vec{p})\n", "\\end{equation*}\n", "\n" ] }, { "cell_type": "code", "execution_count": 115, "metadata": {}, "outputs": [], "source": [ "def forwardEuler1D( fn, t, y0, params ):\n", " \"\"\" [t,y] = forwardEuler1D( fn, t, y0, params )\n", " Given a function handle fn for dydt, find the value of\n", " y at timepoints t, given y(0) = y0.\n", " fn takes t, y, and params as input.\n", " t is 1-dimensional (NT,)\n", " y is 2-dimensional (1,NT) (NOTE WE ASSUME WE HAVE JUST ONE VARIABLE)\n", " \"\"\"\n", " y = np.zeros( (1,t.size) )\n", " y[0,0] = y0;\n", " for i in range(1,t.size):\n", " h = t[i]-t[i-1];\n", " yprime = fn(t[i-1],y[0,i-1],params);\n", " y[0,i] = y[0,i-1] + h*yprime;\n", " \n", " return (t,y)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Compare the method's solution to the true solution\n", "\n", "Make an array of time steps starting at t=0, ending at t=4, with a time step of 0.1. Use y0 = 1 and a degradation rate constant of 2. \n", "\n", "Evaluate the true solution at those time steps, then use forward Euler to approximate the true solution at those time steps.\n", "\n", "Plot both solutions in an upper subplot and the difference between them in a lower subplot." ] }, { "cell_type": "code", "execution_count": 129, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjgAAAEYCAYAAABRMYxdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy86wFpkAAAACXBIWXMAAAsTAAALEwEAmpwYAABUHElEQVR4nO3deXwV1dnA8d9zl+whIRtLAiRA2HfCUhFUEMWqgIIVFQWrVYpW0VZr31qr2L7W6qu1xWpRXFAEKipSa+sCorIIhH0TDBAgrCGBkD259573j7mEEAJcIMlNbp7v5zOfO3PmzNznZMLw5MyZGTHGoJRSSikVSGz+DkAppZRSqqZpgqOUUkqpgKMJjlJKKaUCjiY4SimllAo4muAopZRSKuA4/B1AVXFxcSY5OdnfYSillFKqAVi9evURY0x81fJ6l+AkJyeTnp7u7zCUUkop1QCIyO7qys95iUpE3hCRwyKy6QzrRUT+KiIZIrJBRPpUWjdBRH7wThMuPHyllFJKKd/5MgbnLWDEWdZfA6R6p3uAVwBEJAb4PTAA6A/8XkSaXkywSimllFK+OGeCY4z5Bsg9S5VRwExj+Q6IFpEWwNXAF8aYXGPMUeALzp4o1bqC40dJ/9c//BmCUkoppepATYzBSQT2VlrO8padqfw0InIPVu8PrVu3roGQqrfxwz/zo8y/sz46kZ6Dr6u171FKKVX/lZeXk5WVRUlJib9DUT4ICQkhKSkJp9PpU/16McjYGDMdmA6QlpZWay/H6v2T/+HAc/OIXPQYx/sMpUl4WG19lVJKqXouKyuLyMhIkpOTERF/h6POwhhDTk4OWVlZpKSk+LRNTTwHZx/QqtJykrfsTOV+ExIWSfGVf6Kt2cu3M6f6MxSllFJ+VlJSQmxsrCY3DYCIEBsbe169bTWR4CwA7vDeTTUQyDPGHAA+A64SkabewcVXecv8qu2gMfzQdDBXHHyDZWvW+zscpZRSfqTJTcNxvsfKl9vEZwPLgY4ikiUid4nIJBGZ5K3yKbATyABeAyYDGGNygaeBVd5pqrfM71rf9jdsYij516PkFZX7OxyllFJK1bBzjsExxtxyjvUGuO8M694A3riw0GpPcFwKB9IeZGj6c0yf8yb3/PQef4eklFKqkcnJyWHYsGEAHDx4ELvdTny89UDelStXEhQUdFH7/+1vf4vL5eLZZ58FYPfu3VxxxRWsWbOG6Ojoi9p3Q9Bo30XVYsQj5Ia24arM5/lq0x5/h6OUUqqRiY2NZd26daxbt45Jkybx0EMPVSwHBQXhcrkuav+PP/448+fPZ+vWrQA8+OCDPP30040iuYFGnODgCCbyxpdIth3ihw//oJeqlFJK+d3EiROZNGkSAwYM4NFHH+XJJ5/k+eefr1jfrVs3MjMzAXj33Xfp378/vXr14t5778Xtdp+yr9DQUF588UXuu+8+Pv30U/Lz87ntttvqsjl+VS9uE/cXZ+oVHGs3kgkZH/L8B6P47e3X+jskpZRSfvDUvzazZf/xGt1nl5ZN+P31Xc97u6ysLJYtW4bdbufJJ5+sts7WrVuZO3cuS5cuxel0MnnyZGbNmsUdd9xxSr0f//jHzJgxgwkTJrBkyZILaUaD1agTHIDoUX+m9KWFXLL9WT7f1IerurXwd0hKKaUasZtuugm73X7WOgsXLmT16tX069cPgOLiYhISEqqte99991FcXEzHjh1rPNb6rNEnODRpgX3ob7nii//hVx++Qb+UR2gafnEDu5RSSjUsF9LTUlvCw8Mr5h0OBx6Pp2L5xHNgjDFMmDCBZ5555pz7s9ls2GyNb0RK42txNRwD76UkpjMPu2fwvx+v8nc4SimlFADJycmsWbMGgDVr1rBr1y4Ahg0bxrx58zh8+DAAubm57N69229x1kea4ADYHYSM/gstJYd2W17hv5sO+jsipZRSijFjxpCbm0vXrl2ZNm0aHTp0AKBLly784Q9/4KqrrqJHjx4MHz6cAwcO+Dna+kWsx9jUH2lpaSY9Pd0v3+35aDKe9XO4xf48/3h4PDF6qUoppQLW1q1b6dy5s7/DUOehumMmIquNMWlV62oPTiW2q6YiwRH8qnw6T8zf6O9wlFJKKXWBNMGpLDwO+/AnGWDbin3zPD7dqN19SimlVEOkCU5VfSZgWvblyeD3+NNHKzhSUOrviJRSSil1njTBqcpmR677P6I5zt3ls/nd/E3Ut3FKSimllDo7TXCq07I3knYX4+2fs2fzct5Pz/J3REoppZQ6D5rgnMnQx5HwWF6KfIffzd/A6t1H/R2RUkoppXykCc6ZhEYjV/2B9mXf80D459z7zmoO5BX7OyqllFIBZv78+YgI33//vd9i2L9/P2PHjq2Rfc2fP58tW7b4XH/z5s106NCB4uKT/8dee+21zJ49+6Li0ATnbHrcDJ2uY7LrXTqVb+GemaspLnOfezullFLKR7Nnz+bSSy+96P/QK3O5XOdVv2XLlsybN69Gvvt8E5yuXbty44038sc//rFi+/Lycm655ZaLisOnBEdERojINhHJEJHHqln/oois807bReRYpXXuSusWXFS0dU0ERr2MRCXxetg09u/fy6MfbNBBx0oppWpEQUEBS5YsYcaMGcyZM6eifPHixQwZMoRrr72Wjh07MmnSpIp3UkVERPDQQw/RtWtXhg0bRnZ2NgCXX345U6ZMIS0tjZdeeomFCxfSu3dvunfvzk9/+lNKS0tZtWoVPXr0oKSkhMLCQrp27cqmTZvIzMykW7duALz11luMHj2a4cOHk5yczLRp03jhhRfo3bs3AwcOJDc3F4DXXnuNfv360bNnT8aMGUNRURHLli1jwYIFPPLII/Tq1YsdO3awY8cORowYQd++fRk8eHC1PVVPPPEE77//PuvWreOxxx7j5Zdfvuif7TlftikiduBlYDiQBawSkQXGmIr0zBjzUKX6vwB6V9pFsTGm10VH6i+h0XDT24TMuIqPmr/NZevvo3OLSCZf3t7fkSmllKop/3kMDtbwA16bd4dr/nTWKh9//DEjRoygQ4cOxMbGsnr1avr27QvAypUr2bJlC23atGHEiBF8+OGHjB07lsLCQtLS0njxxReZOnUqTz31FNOmTQOgrKyM9PR0SkpKSE1NZeHChXTo0IE77riDV155hSlTpjBy5Egef/xxiouLGT9+PN26dSMzM/OUuDZt2sTatWspKSmhffv2PPvss6xdu5aHHnqImTNnMmXKFG688UZ+9rOfAfD4448zY8YMfvGLXzBy5Eiuu+66iktew4YN49VXXyU1NZUVK1YwefJkFi1adMr3hYWF8fzzzzNkyBAefvhhUlNTL/rH70sPTn8gwxiz0xhTBswBRp2l/i1AzfWz1Qcte8E1f6L10eX8LWkRz322jS+3HPJ3VEoppRq42bNnM27cOADGjRt3ymWq/v3707ZtW+x2O7fccgtLliwBrLeD33zzzQCMHz++ohyoKN+2bRspKSkV766aMGEC33zzDWD1lnzxxRekp6fz6KOPVhvXFVdcQWRkJPHx8URFRXH99dcD0L1794pkaNOmTQwePJju3bsza9YsNm/efNp+CgoKWLZsGTfddBO9evXi3nvvPeM7s66//nqio6OZPHmybz+8czhnDw6QCOyttJwFDKiuooi0AVKAyqlZiIikAy7gT8aY+RcWqp/1vRN2L+faTW+xPL49U+Y6+GjyJaQ2i/R3ZEoppS7WOXpaakNubi6LFi1i48aNiAhutxsR4bnnngNARE6pX3W5uvLw8PBzfm9OTg4FBQWUl5dTUlJS7TbBwcEV8zabrWLZZrNVjO+ZOHEi8+fPp2fPnrz11lssXrz4tP14PB6io6NZt27dOeM6sX+brWaGB9f0IONxwDxjTOWRuG28L8G6FfiLiLSrupGI3CMi6SKSfuJaYr0jAte9iMSmMtX1IomO49w9M51jRWX+jkwppVQDNG/ePG6//XZ2795NZmYme/fuJSUlhW+//RawLlHt2rULj8fD3LlzufTSSwEraTgxIPi9996rKK+sY8eOZGZmkpGRAcA777zDZZddBsC9997L008/zW233cavf/3rC44/Pz+fFi1aUF5ezqxZsyrKIyMjyc/PB6BJkyakpKTw/vvvA2CMYf369Rf8nefDlwRnH9Cq0nKSt6w646hyecoYs8/7uRNYzKnjc07UmW6MSTPGpMXHx/sQkp8ER8BPZmIvL2Re/GscPlbIfe+tweX2+DsypZRSDczs2bO54YYbTikbM2ZMxWWqfv36cf/999O5c2dSUlIq6oaHh7Ny5Uq6devGokWLeOKJJ07bd0hICG+++SY33XQT3bt3x2azMWnSJGbOnInT6eTWW2/lscceY9WqVaeNh/HV008/zYABAxg0aBCdOnWqKB83bhzPPfccvXv3ZseOHcyaNYsZM2bQs2dPunbtyscff3xB33e+5Fx3BImIA9gODMNKbFYBtxpjNlep1wn4L5BivDsVkaZAkTGmVETigOXAqMoDlKtKS0sz6enpF9GkOrB+Dnx0L1vb3801m4Yy8ZJknhzZ1d9RKaWUOg9bt26lc+fO/g6jWosXL+b555/nk08+OW1dREQEBQUFfojK/6o7ZiKy2nul6BTn7MExxriA+4HPgK3AP40xm0VkqoiMrFR1HDDHnJoxdQbSRWQ98BXWGBzfb46vr3qOgz4T6JzxOs90O8BbyzKZu2qPv6NSSimllNc5e3DqWoPowQEoL4EZV2Lysngo+m/8e4+d2T8bSFpyjL8jU0op5YP63IOjqlejPTjqDJwhcNPbiMfN87xImygHk95dzb5j+joHpZRqKOrbH/nqzM73WGmCczFi28GoaTgOrGZe+/9SWu7hnpnp+joHpZRqAEJCQsjJydEkpwEwxpCTk0NISIjP2/jyHBx1Nl1GwYCfE73iFWZf2ovrF8Vy33treGV8H4Iddn9Hp5RS6gySkpLIysqi3j6eRJ0iJCSEpKQkn+vrGJya4CqDN6+BI9tZMGA2D3yex9BOCZrkKKWUUrVMx+DUJkcQ3PQW2OyM3P4YfxrZnkXfH+bn766h1KWXq5RSSqm6pglOTYluBTdMh4MbGbf/Wf53dGdNcpRSSik/0QSnJnW4Cq58CjZ9wK2H/8IfR3fVJEcppZTyA01watqlU2DwL2HN29x2bLomOUoppZQfaIJTG4b+DgZMgu9e5rbi9/jjDd00yVFKKaXqkCY4tUEErn4Geo+Hr5/lNtfHmuQopZRSdUgTnNpis8H1f4WuN8IXv+M225ea5CillFJ1RBOc2mSzw43TocMI+PcvuS3kO01ylFJKqTqgCU5tszvhprchZTDM/zm3NdmgSY5SSilVyzTBqQvOEBg3GxL7wvt3clvsD5rkKKWUUrVIE5y6EhwBt70PCZ1gznhua5ZVkeTc/XY6x0vK/R2hUkopFTA0walLodFw+3zrqcfv3cxtSUf485geLN+Rww0vL2VndoG/I1RKKaUCgiY4dS08Du74GMJi4N0x/KR1Pu/ePYDcwjJGv7yUb7brW22VUkqpi+VTgiMiI0Rkm4hkiMhj1ayfKCLZIrLOO91dad0EEfnBO02oyeAbrCYtYcICcITAzFEMjDrGgvsvpWV0KBPfXMnr3+6kvr3lXSmllGpIzpngiIgdeBm4BugC3CIiXaqpOtcY08s7ve7dNgb4PTAA6A/8XkSa1lj0DVnTZKsnx7jhzWtoVbyVD35+CVd1ac4f/r2VX72/gZJyHXyslFJKXQhfenD6AxnGmJ3GmDJgDjDKx/1fDXxhjMk1xhwFvgBGXFioASi+I0z8NziC4c0fE57xL/5+Wx+mXJnKB2uyGDf9Ow4fL/F3lEoppVSD40uCkwjsrbSc5S2raoyIbBCReSLS6ny2FZF7RCRdRNKzsxvZGJSEznD3ImjRE96fiG3J80wZlsort/Vh28F8rp+2hPV7j/k7SqWUUqpBqalBxv8Cko0xPbB6ad4+n42NMdONMWnGmLT4+PgaCqkBiYiHOxZAj5th0R/go3u5pnMMH/z8Ehw2Gz/5x3I+XrfP31EqpZRSDYYvCc4+oFWl5SRvWQVjTI4xptS7+DrQ19dtlZczBG74Bwx9HDbMhbevp0uTUhbcP4heraJ5cM46nvnPVtweHXyslFJKnYsvCc4qIFVEUkQkCBgHLKhcQURaVFocCWz1zn8GXCUiTb2Di6/ylqnqiMCQR6xXOxzYAK8PJbZoJ+/ePYDxA1vzj693cvfbq/ShgEoppdQ5nDPBMca4gPuxEpOtwD+NMZtFZKqIjPRWe0BENovIeuABYKJ321zgaawkaRUw1VumzqbraLjzU3CVwevDce5cxB9Gd+cPo7vx7Q9HGD1tqY7LUUoppc5C6tvzVtLS0kx6erq/w6gf8vbB7Jvh0GYY8Sfofw/f7cplypx1ZBeUMumytjwwLJVgh93fkSqllFJ+ISKrjTFpVcv1Scb1WVQi3Plf6HAN/OdR+PcvGdimCZ89NIQbeyfy8lc7GPm3pWzMyvN3pEoppVS9oglOfRccATe/C4MehPQZMOsmoijkuZt68ubEfhwrLmP035fyf59vo8zl8Xe0SimlVL2gCU5DYLPB8Kkw6mXIXAKvXwn713JFpwQ+n3IZo3q15G+LMhg5bQmb9mlvjlJKKaUJTkPSe7z1eoeyAnhtGCz6I1FBhhd+0ovX70gjx/vCzhe/2K69OUoppRo1TXAamuRBMHk59PgJfPNneG0oHNzIlV2a8cVDQ7iuRwteWvgDo19eypb9x/0drVJKKeUXmuA0RKFN4YZXYdxsKDgE0y+Hr/9MdLDwl3G9+cftfTmcX8Kol5fw14U/UO7W3hyllFKNiyY4DVmnH8N9K6DLaPjqj9bYnMNbubprcz5/6DJGdGvBC19s54a/L2X17qP+jlYppZSqM5rgNHRhMTB2BvxkJuTthX8MgSUvEhNq52+39ObV8X04dLyUMa8sY/Ks1WQeKfR3xEoppVSt0wf9BZKCbPj3w7B1AST1g9GvQFwqhaUuXvt2J9O/2UmZy8P4gW34xdD2xEYE+ztipZRS6qKc6UF/muAEGmNg0wfw6a+gvBiGPQEDfg42G4ePl/CXhT8wd9Vewpx2Jl3ejrsuTSHEqU9CVkop1TBpgtPY5B+CT6bAtk+h9SXw4+egeTcAMg7n86f/bOPLrYdoERXCw8M7cGOfJOw28W/MSiml1HnSBKcxMgY2zIX//BpK8qDnOLjifyC6NQDf7czhmU+3sj4rj07NI/nNjztzWYd4PwetlFJK+U4TnMas+Cgs+QuseBWMB/rfA4N/CWExGGP498YD/Pm/29iTW8Tg1Dgeu6YTXVtG+TtqpZRS6pw0wVHW28kXPwPrZkFQBFw6xRqfExRGqcvNrO/28NdFP5BXXM6Irs25e3Bb+rZp6u+olVJKqTPSBEeddPh7WDgVtv0bIprDFb+BXuPB7iCvuJzXvtnJzOWZHC9x0ad1ND8b3JarujbXMTpKKaXqHU1w1Ol2L4cvfw97V0BsKlz5e+h0HYhQWOri/fS9zFi6i725xbSOCeOng5K5Ka0V4cEOf0eulFJKAReZ4IjICOAlwA68boz5U5X1DwN3Ay4gG/ipMWa3d50b2OituscYM/Js36UJTh0zxrrT6sun4Mg26/k5Vz5lvfMKcHsMn28+yGvf7mTNnmNEhTq5dUBrJl6STLMmIX4OXimlVGN3wQmOiNiB7cBwIAtYBdxijNlSqc4VwApjTJGI/By43Bhzs3ddgTEmwtdANcHxE7cL1r8HXz0D+fsheTAMmAQdrwGb9Zyc1buP8vq3O/ls80HsNuH6ni352eC2dG7RxM/BK6WUaqwuJsH5EfCkMeZq7/JvAIwxz5yhfm9gmjFmkHdZE5yGpKwI0mfAin9Yr36Ibm3dddV7vPWST2B3TiFvLs3kn+l7KSpzc2n7OCZeksxlHeNx2vXtH0opperOxSQ4Y4ERxpi7vcu3AwOMMfefof404KAx5g/eZRewDuvy1Z+MMfOr2eYe4B6A1q1b9929e7fvLVO1w+2yBiGv+AfsXgrOMOh5Cwy4F+I7ApBXVM6slbt5a2kmh/NLiYsIYmTPRG7sk0jXlk0Q0UHJSimlaledJDgiMh64H7jMGFPqLUs0xuwTkbbAImCYMWbHmb5Pe3DqoQPrYcV02Pg+uEuh3VDr8lX74WCzUeby8PX2bD5ck8XCrYcpc3vo2CySMX0TGdUrUcfqKKWUqjW1folKRK4E/oaV3Bw+w77eAj4xxsw70/dpglOPFR6B1W/CqhmQfwBi2kL/e6HXrRBijcM5VlTGvzYc4MM1WazdcwybwKWp8Yzpk8hVXZoTGqTvvVJKKVVzLibBcWANMh4G7MMaZHyrMWZzpTq9gXlYPT0/VCpvChQZY0pFJA5YDoyqPEC5Kk1wGgB3OWz52HoyctYqCIqEHj+B7mOh1UCwWeNwdmYX8NHafXy4Zh/7jhUTEezgx92bc2OfJPonx2DT5+oopZS6SBd7m/iPgb9g3Sb+hjHmjyIyFUg3xiwQkS+B7sAB7yZ7jDEjReQS4B+AB7ABfzHGzDjbd2mC08BkrYaV/4AtC8BVDE0SoesN0G0MtOwNIng8hpWZuXy4JotPNx6koNRFy6gQruzSjGGdmzGwbQzBDu3ZUUopdf70QX+qdpXmw7b/wqYPIONL8JRbl7C6jbGmhM4AFJe5+XzLQT7ZcIAlPxyhuNxNeJCdwanxDOucwNBOCcRGBPu5MUoppRoKTXBU3SnKhe8/sZKdXd9YL/hM6ArdbrSSnZgUAErK3SzbcYQvtx5m0dbDHDxeggj0ad2UYZ0TuLJzM1ITIvRuLKWUUmekCY7yj/xD1nidTR/A3u+sssS+1mWs9sOtW85FMMawef9xvtx6iIVbD7NxXx4ArWJCGdapGVd2bkZaclNCnHopSyml1Ema4Cj/O7YHNn9kJTsH1ltlTRKh3RXQbhi0vRzCYgA4mFfCou8Ps3DrIZZkHKHU5SHIbqNnqygGpMQyoG0Mfds0JSxI34ullFKNmSY4qn45tgd2LLKmnYuhJA8QSOxjJTvthlrvxbI7KC5zs3znEVbszOW7Xbls2peH22Nw2ITuSVEMbBvLgJQY0pJjiNAXgSqlVKOiCY6qv9wu2L/GSnYyFsK+dGvcTnATSBliJTvthkLTZBChoNRFemYuK3blsmJnDhuy8nB5DHab0K1lEwZ4E57erZsSEx7k79YppZSqRZrgqIaj+Kg1ODljoZX05O21yiOaQ6t+kNTf6t1p2QucoRSVuViz+xgrduWwYmcu6/Yeo8ztASAxOpQeSVF0T4qiR2I03ROjiApz+q9tSimlapQmOKphMgZyMqzLWHtXQtZKOJpprbM5oHl3K+Fp1R+S0iC6DSUuD2v3HGND1jE27Mtj0748ducUVeyydUyYN+GJontiFF0To4gK1aRHKaUaIk1wVOAoyLaeoHxi2rcayr0JTHiC1bvTqh+06Gndnh6RwLHicjbtO87GfXls3HeMDVl5ZB0trthlSlw4XVo0IbVZBKkJkaQ2iyA5Npwgh74dXSml6jNNcFTgcrvg8JaTCc/elZBb6X2uYbGQ0AWadbWmhK6Q0Imj5U5vwpPHhqxjbDuYz+7cIk78k3DYhOS4cFITIqypmZX4pMSF65OXlVKqntAERzUuRblwaBMc2gKHN3s/t5zs6UGsQcsVSU8XiO9ESWQrdhx1kXG4gB8OFbD9UD4ZhwvIzCnE4/2nYrcJbWLDaBcfQZuYMNrEhtE6NpzWMWEkRodqr49SStUhTXCU8njgWObJZOfQZmvK3WHdtQWAWM/miUmB2HbW6yZi2lLaJJld7gS2H3WTcSifHw4XsCO7gD25RZSUeyq+wibQIiqUNrFW4tMqJow2MVby0zo2jCYhDn0ys1JK1SBNcJQ6k/JiyN4GR36A3J2nTkVHTq0b2QJi2lkJUNNkTJNEjjoTyHLHkFESRWaeiz05hezOLWJvbhFHCspO2Tw8yE7zqBBaRofSvEkILaJDaREVYpVFhdI8KkSTIKWUOg9nSnD0qWhKOUOtW85b9jp9XfExOLqrUtKzC3J2wPbPoPAwAsR4px4A4fHeHqAkSE6kNLwF2RJHlieGXaVR7CgJZ1++hwN5JWw/lM3h/FKq/o1xIglqERVKQmQwcZHBxEUEER8ZTFzEySkmPAi7TRMhpZSqjiY4Sp1NaDSE9oaWvU9fV1YIx/dDXhYc3wd5++B4llWWswN2fUNw6XGSgCRg4IntQqKsu71aJOBpn0BRUAx5tmiOmCgOupuwtyySncUethfCiiOFZBeUUubynPb1NoGY8FOTn9jwIJqGBxEV6qRpWBBNw5xEhwXRNNxa1nd5KaUaC01wlLpQQeEQl2pNZ1JyvFLysw8KDkPhYe9nNrZDG4koyCaiNI9EoOdp3xGBiW2KJySaUmc0xY4mFEgkeUSSa8LIdoVzsDyUrPxQ9hwKZlVxENnlwZRS/ROcQ5w2moYFWUlPmJX0RIY4aBLqpEmIg8gQJ01CHUQGO62yUG9ZiIPwIAc27TFSSjUQmuAoVZtCmlhTQuez1ysvgcJsb/KTDQWHrPmio0jxUezFuYQVHyUsfzuxxUetpz0b9+n7sVuTsQfhdkZQ7oygzB5BsS2CIgmjgFDyTBjHPCHk5oWQcySYo64g9pY5OOYOpsgEU0gIRSbE+iSEcu9pQgQigx1EBDsIC3YQHuwgIthOeJA1Hx5st8qCrPUR3uWwIDshTjthQQ5CnfaK5dAgO6FOu15mU0rVCp8SHBEZAbyEdfp83Rjzpyrrg4GZQF8gB7jZGJPpXfcb4C7ADTxgjPmsxqJXKlA4QyC6lTX5wuOB0uNWolOca30WHYWSY1B6HCk5jqP0OI6S44SW5hNVehxKsq1tSo5bn1Qa/ONNjKrjFgfl9jDKbSGUSCilEkxJWRDFZUEUHw+iyBNEocdBgSeIfLeDQo+TIyaYLJwUE0ypcVKKkzK8n95lq8yBxx6MzRGCPSgEmzMEhzMYZ5CTEIedYKftlM8Qp41gp50Qh/UZfOLTbiPI4Z0qz1derlTutFvLTrtgt4kO6lYqAJ0zwRERO/AyMBzIAlaJyAJjzJZK1e4Cjhpj2ovIOOBZ4GYR6QKMA7oCLYEvRaSDMdX96amU8pnN5h0fFA2knP/2Hg+UFUBpvvVsoLICa0zRaVMB9rJC7GWFhJQXEllaAK4S686z8mJwFUL5EasHqrzIu67o3N9fnTLvBLix4cJJuThw4aAMB+XGThkOyoyDUmOnHIc1GTsu7LhwUI6dIuzkYcftLS+vtM7ar7WNBxsubLjFDuLAiB3sDusVIDYH2OxgcyA2B2Kzg83unbchNjtid2Cz2a15mwOb/cTniXmbd70Du82G2O3YbXbEbsdms2PzLtscDmxix2YXROzY7Dar3G7HLjZrWaxEzGYT7CLYbFSUyYl1AraKecFuw1pXUc9ab/NuXzEvVj1blfXCyTLx1lGqIfGlB6c/kGGM2QkgInOAUUDlBGcU8KR3fh4wTax/DaOAOcaYUmCXiGR497e8ZsJXSl0Qm+3k5bOaZgy4SsF1IgkqAVeZ97MU3KXe9aXVlJWAuxy7uwy7u4xgdzm4y7x1vPPuMoy7HOMqxeMqw7hKMW4Xxl0MHpdVz+MCj/UpHhfiKcfmcSHGjc24zhI7Vl9zPfoTzGMEDycmGx5sGKj4NN7yU8tOXfZgPXzSZQTDielkufHu35w2UfFpObkOsdadKAPBCJWWrXkAI6fWs1ad3L946wicUsd4653cV+X9ndw/3niosr+T6yvXOVl2si7eOes7pfL6SvutvO2p5aduczHrqZJImsrfWXnFWRPOsyejp29afX1zxu+QsyydXmhvcwn9R006a0y1wZcEJxHYW2k5CxhwpjrGGJeI5AGx3vLvqmybWPULROQe4B6A1q1b+xq7Uqo+ErEuuTlDILRp7XyFd7qgZ0YbAx63NwlyWWOZKi9XTO4q8+6TdU/59Jy6D+O2Hhzp8VifxlpvPG48HjcetwfjceH2eDBuFx6PB4/HjTEejMeaPN79ejwnylwY7z4r6hkDxlrGGIzHjcFb5jHeh1cazIk4MFY9b32MNw0y1r5sGMSc3O+J+laWYbzLAJ6TZcac+olBvPMGvPOcVgdzMn06sfrEvGAwhirrT02zTmx06v4r7aPKsxfkRDyn7MPa5vT9nrrNifmT1U5N987mtPXVxHVe25+h3Jxl3blI1RDP87t98UNQ+HlvUxPqxSBjY8x0YDpYD/rzczhKqUAmYl2Kstft6U8461AnpQJW1R6RuuLLH0D7gMojH5O8ZdXWEREHEIU12NiXbZVSSimlapQvCc4qIFVEUkQkCGvQ8IIqdRYAE7zzY4FFxnoHxAJgnIgEi0gKkAqsrJnQlVJKKaWqd84+Wu+YmvuBz7B6V98wxmwWkalAujFmATADeMc7iDgXKwnCW++fWAOSXcB9egeVUkoppWpbvXvZpohkA7tr8SvigCPnrBV4GmO7tc2Ng7a58WiM7dY2n1sbY0x81cJ6l+DUNhFJr+6to4GuMbZb29w4aJsbj8bYbm3zhbuguyyVUkoppeozTXCUUkopFXAaY4Iz3d8B+EljbLe2uXHQNjcejbHd2uYL1OjG4CillFIq8DXGHhyllFJKBThNcJRSSikVcAI2wRGRESKyTUQyROSxatYHi8hc7/oVIpLshzBrlA9tnigi2SKyzjvd7Y84a5KIvCEih0Vk0xnWi4j81fsz2SAifeo6xprmQ5svF5G8Ssf5ibqOsaaJSCsR+UpEtojIZhF5sJo6AXWsfWxzIB7rEBFZKSLrve1+qpo6AXX+9rHNAXf+BhARu4isFZFPqll3ccfZGBNwE9YTl3cAbYEgYD3QpUqdycCr3vlxwFx/x10HbZ4ITPN3rDXc7iFAH2DTGdb/GPgP1rsOBwIr/B1zHbT5cuATf8dZw21uAfTxzkcC26v5/Q6oY+1jmwPxWAsQ4Z13AiuAgVXqBNr525c2B9z529uuh4H3qvs9vtjjHKg9OP2BDGPMTmNMGTAHGFWlzijgbe/8PGCYiEgdxljTfGlzwDHGfIP1epAzGQXMNJbvgGgRaVE30dUOH9occIwxB4wxa7zz+cBWILFKtYA61j62OeB4j1+Bd9HpnareDRNQ528f2xxwRCQJuBZ4/QxVLuo4B2qCkwjsrbScxeknhoo6xhgXkAfE1kl0tcOXNgOM8XbfzxORVtWsDzS+/lwCzY+83d3/EZGu/g6mJnm7qXtj/ZVbWcAe67O0GQLwWHsvW6wDDgNfGGPOeKwD5PztS5sh8M7ffwEeBTxnWH9RxzlQExxVvX8BycaYHsAXnMyMVWBZg/Vulp7A34D5/g2n5ohIBPABMMUYc9zf8dSFc7Q5II+1McZtjOkFJAH9RaSbn0OqdT60OaDO3yJyHXDYGLO6tr4jUBOcfUDl7DbJW1ZtHRFxAFFATp1EVzvO2WZjTI4xptS7+DrQt45i8ydffhcCijHm+InubmPMp4BTROL8HNZFExEn1n/0s4wxH1ZTJeCO9bnaHKjH+gRjzDHgK2BElVWBdv6ucKY2B+D5exAwUkQysYZUDBWRd6vUuajjHKgJziogVURSRCQIa3DSgip1FgATvPNjgUXGO5KpgTpnm6uMRxiJdU0/0C0A7vDeYTMQyDPGHPB3ULVJRJqfuE4tIv2x/p036JO/tz0zgK3GmBfOUC2gjrUvbQ7QYx0vItHe+VBgOPB9lWoBdf72pc2Bdv42xvzGGJNkjEnG+v9qkTFmfJVqF3WcHTUSaT1jjHGJyP3AZ1h3F71hjNksIlOBdGPMAqwTxzsikoE1YHOc/yK+eD62+QERGQm4sNo80W8B1xARmY11J0mciGQBv8caoIcx5lXgU6y7azKAIuBO/0Rac3xo81jg5yLiAoqBcQ355O81CLgd2OgdpwDwP0BrCNhj7UubA/FYtwDeFhE7VsL2T2PMJ4F8/sa3Ngfc+bs6NXmc9VUNSimllAo4gXqJSimllFKNmCY4SimllAo4muAopZRSKuBogqOUUkqpgKMJjlJKKaUCjiY4SimllAo4muAopZRSKuBogqOUUkqpgKMJjlJKKaUCjiY4SimllAo4muAopZRSKuDUu5dtxsXFmeTkZH+HoZRSSqkGYPXq1UeMMfFVy+tdgpOcnEx6erq/w1BKKaVUAyAiu6sr9+kSlYiMEJFtIpIhIo9Vsz5YROZ6168QkeQq61uLSIGI/OqColdKKaWUOg/nTHBExA68DFwDdAFuEZEuVardBRw1xrQHXgSerbL+BeA/Fx+uUkoppdS5+XKJqj+QYYzZCSAic4BRwJZKdUYBT3rn5wHTRESMMUZERgO7gMKaClqdrqjMRXZ+KYfzSzl8vJTs/BJrPr+UvOJyOjaLpF9KDH1aRxMZ4vR3uEoppVSt8iXBSQT2VlrOAgacqY4xxiUieUCsiJQAvwaGA3p5qgZs2pfHR2v3cei4lcBke6eCUtdpdR02IT4ymPBgB4u+P8y0rzKwCXRu0YR+yTGkJTelX3IMzZqE+KElSinlm/LycrKysigpKfF3KMqPQkJCSEpKwun07Y/02h5k/CTwojGmQETOWElE7gHuAWjdunUth9QwHSsq4/nPtzFrxR6C7DZaRocSHxFMl5ZNSIgMJiEyhPjIYGu+STDxEcE0DQvCZrN+7oWlLtbtPcbKXbmk785l7qq9vLUsE4BWMaH0S47xTk1pFx/B2Y6XUkrVpaysLCIjI0lOTtZzUyNljCEnJ4esrCxSUlJ82saXBGcf0KrScpK3rLo6WSLiAKKAHKyenrEi8mcgGvCISIkxZlqVwKcD0wHS0tKMT5E3Eh6PYW76Xv783+/JKy5nwo+SeWh4B6JCz+8yU3iwg0Ht4xjUPg6AcreHrQeOWwlP5lG+2Z7Nh2uswxoTHsQ13Zrz00tTaBcfUeNtUkqp81FSUqLJTSMnIsTGxpKdne3zNr4kOKuAVBFJwUpkxgG3VqmzAJgALAfGAouMMQYYXCm4J4GCqsmNOrP1e4/xxMebWJ+VR7/kpkwd1Y3OLZrUyL6ddhs9kqLpkRTN3YOt7Dgzp4hVmbksyzjC+6uzmLViD0M7JXD3pSn8qF2snlyUUn6j5x91vr8D50xwvGNq7gc+A+zAG8aYzSIyFUg3xiwAZgDviEgGkIuVBKkLlFtYxnOffc+cVXuJiwjmxZt7MrpXYq3+AxcRUuLCSYkL5ydprXi8oJRZ3+3hne8yufX1FXRqHsndg9tyfc8WBDvstRaHUkopVRPE6mipP9LS0kxjfdCf22OYvXIPz3++jfwSF3deksyDV6b69a6nknI3C9bvZ8a3u9h2KJ/4yGDuGNiG2wa2ISY8yG9xKaUaj61bt9K5c2e/fX9OTg7Dhg0D4ODBg9jtduLjrQfnrly5kqCg2j8Xzp8/nxtuuIGtW7fSqVOnWv++6uzfv58HHniAefPmXfS+5s+fT4cOHejSpepTZ86uut8FEVltjEmrWlcTnHpizZ6jPPHxJjbtO87AtjFMHdWNDs0i/R1WBWMMSzKO8Pq3u/h6ezbBDhs39knirkuTaZ9Qf+JUSgUefyc4lT355JNERETwq1+dvDHY5XLhcNTuPTs333wz+/fvZ+jQoTz11FM1ss+6iPtMJk6cyHXXXcfYsWPPa7vzSXD0ZZt+lltYxqPz1nPj35eRnV/KX2/pzeyfDaxXyQ1Yl7AGp8bz9k/788VDQ7ixTyIfrsniyhe+YeKbK/luZ46/Q1RKqTozceJEJk2axIABA3j00Ud58sknef755yvWd+vWjczMTADeffdd+vfvT69evbj33ntxu93n9V0FBQUsWbKEGTNmMGfOnIryxYsXM2TIEK699lo6duzIpEmT8Hg8AERERPDQQw/RtWtXhg0bVjE49/LLL2fKlCmkpaXx0ksvsXDhQnr37k337t356U9/SmlpKatWraJHjx6UlJRQWFhI165d2bRpE5mZmXTr1g2At956i9GjRzN8+HCSk5OZNm0aL7zwAr1792bgwIHk5uYC8Nprr9GvXz969uzJmDFjKCoqYtmyZSxYsIBHHnmEXr16sWPHDnbs2MGIESPo27cvgwcP5vvvv7/gY3NCvXsXVWNSXObmlunfsSO7gHuHtOUXw1KJCK7/hyS1WSTP3NiDX13VkVkr9jBz+W7GTf+O4V2a8T8/7kxKXLi/Q1RKBain/rWZLfuP1+g+u7Rswu+v73re22VlZbFs2TLsdjtPPvlktXW2bt3K3LlzWbp0KU6nk8mTJzNr1izuuOMOn7/n448/ZsSIEXTo0IHY2FhWr15N3759AesS2ZYtW2jTpg0jRozgww8/ZOzYsRQWFpKWlsaLL77I1KlTeeqpp5g2zbrHp6ysjPT0dEpKSkhNTWXhwoV06NCBO+64g1deeYUpU6YwcuRIHn/8cYqLixk/fvwpCdsJmzZtYu3atZSUlNC+fXueffZZ1q5dy0MPPcTMmTOZMmUKN954Iz/72c8AePzxx5kxYwa/+MUvGDly5Ck9OMOGDePVV18lNTWVFStWMHnyZBYtWnSeR+RU9f9/0wD2xMeb2H44n7fu7M9lHU57EWq9FxsRzAPDUrlnSFtmLNnF37/KYPgLX3PHj5J5cFgqUWH6xGSlVOC66aabsNvPftPFwoULWb16Nf369QOguLiYhISE8/qe2bNn8+CDDwIwbtw4Zs+eXZHg9O/fn7Zt2wJwyy23sGTJEsaOHYvNZuPmm28GYPz48dx4440V+ztRvm3bNlJSUujQoQMAEyZM4OWXX2bKlCk88cQT9OvXj5CQEP76179WG9cVV1xBZGQkkZGRREVFcf311wPQvXt3NmzYAFhJ0OOPP86xY8coKCjg6quvPm0/BQUFLFu2jJtuuqmirLS09Lx+RtXRBMdP5q3O4v3VWfxiaPsGmdxUFuK0c98V7flJWite+GIbby3bxYdrs3hwWCrjB7bBadcroUqpmnEhPS21JTz8ZG+1w+GouDwEVDx12RjDhAkTeOaZZ864nxUrVnDvvfcCMHXqVEaOHFmxLjc3l0WLFrFx40ZEBLfbjYjw3HPPAaffOn2mu20rl1eO+0xycnIoKCigvLyckpKSarcJDg6umLfZbBXLNpsNl8t6uv7EiROZP38+PXv25K233mLx4sWn7cfj8RAdHc26devOGdf50P95/OCHQ/n8bv4mBraNYcqVHfwdTo2JjwzmmRt78O8HBtOtZRRP/WsLV7/4DV9uOUR9G8yulFI1KTk5mTVr1gCwZs0adu3aBViXXubNm8fhw4cBK2HZvXv3KdsOGDCAdevWsW7dulOSG4B58+Zx++23s3v3bjIzM9m7dy8pKSl8++23gHWJateuXXg8HubOncull14KWEnDibud3nvvvYryyjp27EhmZiYZGRkAvPPOO1x22WUA3HvvvTz99NPcdttt/PrXv77gn0t+fj4tWrSgvLycWbNmVZRHRkaSn58PQJMmTUhJSeH9998HrKRw/fr1F/ydJ2iCU8eKylxMnrWG8GA7fx3XG7st8B5e1blFE965qz9vTExDBO6emc5tr6+o8evmSilVX4wZM4bc3Fy6du3KtGnTKi77dOnShT/84Q9cddVV9OjRg+HDh3PgwAGf9zt79mxuuOGG075r9uzZAPTr14/777+fzp07k5KSUlE3PDyclStX0q1bNxYtWsQTTzxx2r5DQkJ48803uemmm+jevTs2m41JkyYxc+ZMnE4nt956K4899hirVq264PEwTz/9NAMGDGDQoEGn3N4+btw4nnvuOXr37s2OHTuYNWsWM2bMoGfPnnTt2pWPP/74gr6vMr1NvA4ZY/jl++v5aO0+3r1rQMVrEwJZudvDeyv28OKX28krLucnfVvxy6s7kBCpL/hUSvmmPt0mXp8sXryY559/nk8++eS0dRERERQUFPghqtqlt4nXU++vzuLDNft4YGhqo0huwHolxIRLkvn6V1dw16AUPlybxeXPLeblrzIodZ3frZJKKaWUrzTBqSPbDubzxMebGNQ+lgeGpfo7nDoXFebk8eu68MVDlzE4NY7nPtvGNS99q8/PUUqpC3T55ZdX23sDBGTvzfnSBKcOFJa6mDxrNZEhTv5yc2COu/FVclw4/7g9jbd/2p9yt4dx07/jV++vJ7ewzN+hKaXqsfo2nELVvfP9HdAEp5YZY/jtRxvZdaSQl8b1Ij4y+NwbNQKXdYjn8ymXMfnydsxfu49h/7eY99P36klMKXWakJAQcnJy9PzQiBljyMnJISTE9/GbOsi4ls1ZuYfHPtzIw8M7NMpLU77YdjCf3360kfTdRxmQEsMfb+hO+4QIf4ellKonysvLycrKqni2jGqcQkJCSEpKwuk89SGy+rJNP9h64DijX15K/5QY3rqzf6O+NHUuHo9hbvpenvl0K8Xlbn5+WTsmX9GeEOfZnxKqlFKqcdO7qOpYQamL+2atISrUyYs399Lk5hxsNuGW/q1Z+MvLubZ7C/66KINrXvqWpRlH/B2aUkqpBkgTnFpgjOE3H24kM6eQv97Sm7gIHXfjq/jIYP4yrjfv3jUAYwy3vb6Ch+au40jBxb+XRCmlVOPhU4IjIiNEZJuIZIjIY9WsDxaRud71K0Qk2VveX0TWeaf1InLDaTsPQO+t3MO/1u/nl1d1ZGDbWH+H0yBdmhrHf6cM4YGh7flkw36G/d/X/HOVDkJWSinlm3MmOCJiB14GrgG6ALeISJcq1e4Cjhpj2gMvAs96yzcBacaYXsAI4B8iEtAv+Ny0L4+n/rWFIR3i+fll7fwdToMW4rTz8FUd+c+Dg+nYLJJHP9jAba+vIPNIob9DU0opVc/50oPTH8gwxuw0xpQBc4BRVeqMAt72zs8DhomIGGOKjDEub3kIENB/frs9hgfnrCUmLIgXf9ITm467qRHtEyKZc89A/veG7mzcl8fVf/mGvy/OoNztOffGSimlGiVfEpxEYG+l5SxvWbV1vAlNHhALICIDRGQzsBGYVCnhCTiLtx1mR3Yhj1/XmVgdd1OjbDbh1gGt+fLhy7iiYwJ//u82Rk5byoasY/4OTSmlVD1U64OMjTErjDFdgX7Ab0TktKf0iMg9IpIuIunZ2dm1HVKteXv5bpo3CeHqrs39HUrAatYkhFdv78ur4/uSU1DK6JeX8vQnWygqC9i8WSml1AXwJcHZB7SqtJzkLau2jneMTRRwykuGjDFbgQKgW9UvMMZMN8akGWPS4uPjfY++HtmZXcA327O5bUBrnHa9Oa22jejWnC9/eRm39G/NjCW7uOrFb/h6e8NNjpVSStUsX/4nXgWkikiKiAQB44AFVeosACZ458cCi4wxxruNA0BE2gCdgMwaibyeeee73Tjtwrj+rf0dSqPRJMTJH2/ozj/v/RHBDhsT3ljJlDlrydFbypVSqtE7Z4LjHTNzP/AZsBX4pzFms4hMFZGR3mozgFgRyQAeBk7cSn4psF5E1gEfAZONMQH35LbCUhfz0rO4tnsLfdeUH/RPieHTBwfzwLBU/r3xAFe+8DUfrsnSW8qVUqoR01c11IB3v9vN4/M38eHkS+jTuqm/w2nUth/K59cfbGDtnmMMTo3jD6O70SY23N9hKaWUqiX6qoZaYoxh5vJMuidG0btVtL/DafQ6NItk3qRLmDqqK2v3HOOqF7/h5a8yKHPpLeVKKdWYaIJzkb7bmcv2QwXc8aM2iOhzb+oDu02440fJfPnwZQztlMBzn23jur99S3pmrr9DU0opVUc0wblIM5dn0jTMyfU9W/o7FFVF86gQXhnfl9fvSKOw1M3YV5fzmw83kldU7u/QlFJK1TJNcC7C/mPFfL7lEDf3a02I0+7vcNQZXNmlGZ8/NISfDU5h7qo9DHthMR+v26eDkJVSKoBpgnMR3luxx3rj9QC9Nby+Cw928Ntru7Dg/ktJjA7lwTnrmPDmKvbkFPk7NKWUUrVAE5wLVOpyM3vlHoZ1bkarmDB/h6N81C0xig8nD+LJ67uwZvdRhr/4tb7XSimlApAmOBfo040HyCksY8KPkv0dijpPdpswcVAKXzw8pOK9Vtf9dYkOQlZKqQCiCc4FenvZbtrGhzOofay/Q1EXqEVUKK/e3pfX7kgjv6Scsa8u55f/XE92vj4JWSmlGjpNcC7A+r3HWLf3GBN+lKy3hgeA4V2a8cXDl/Hzy9uxYP0+hj6/mBlLdullK6WUasA0wbkAM5fvJjzIzo19Ev0diqoh4cEOfj2iE59NGULvNk15+pMtXPfXJSzfkXPujZVSStU7muCcp5yCUv61YT9j+iYRGeL0dziqhrWNj+DtO/sx/fa+FJa5uOW17/jF7LUcyCv2d2hKKaXOgyY452lu+l7KXB5uH9jG36GoWiIiXNW1OV8+fBkPDkvls80HGfZ/X/PK4h36ygellGogNME5Dy63h1nf7eGSdrGkNov0dziqloU47Tw0vANfPnQZg9rH8ex/v2fEX77hm+3Z/g5NKaXUOWiCcx4Wfn+YfceKuUNvDW9UWseG8dodabx5Zz8McMcbK7n3nXT25upDApVSqr7SBOc8zFyeScuoEK7snODvUJQfXNExgf9OGcwjV3fkm+1HGPbC1zz73+85XqLvtlJKqfpGExwfZRzOZ2lGDrcNbIPDrj+2xirYYee+K9qz8JeXcV33FryyeAeXP7eYt5dl6m3lSilVj/j0P7WIjBCRbSKSISKPVbM+WETmetevEJFkb/lwEVktIhu9n0NrOP46M3P5boLsNsb1a+XvUFQ90DI6lBdu7sUnv7iUjs0i+f2CzVz14jf8d9NBfYmnUkrVA+dMcETEDrwMXAN0AW4RkS5Vqt0FHDXGtAdeBJ71lh8BrjfGdAcmAO/UVOB1Kb+knA9WZ3FdzxbERgT7OxxVj3RLjOK9nw3gjYlpOGzCpHdXc9Ory1m756i/Q1NKqUbNlx6c/kCGMWanMaYMmAOMqlJnFPC2d34eMExExBiz1hiz31u+GQgVkQaXIXy0dh+FZW5975SqlogwtFMz/vPgYJ65sTuZOUXc8Pdl3PfeGn1buVJK+YkvCU4isLfScpa3rNo6xhgXkAdUfUnTGGCNMaZBvejHGMPbyzLp2Sqanq2i/R2Oqsccdhu39G/N149czgPDUlm09TDDXljM059s4VhRmb/DU0qpRqVORsuKSFesy1b3nmH9PSKSLiLp2dn16xkjy3bksCO7kAk/0gf7Kd+EBzt4eHgHFj9yOWP6JPHm0l0M+fNXTP9mByXlbn+Hp5RSjYIvCc4+oPLI2iRvWbV1RMQBRAE53uUk4CPgDmPMjuq+wBgz3RiTZoxJi4+PP78W1LK3l2USGx7Ej7u38HcoqoFp1iSEP43pwX8eHEKfNk3530+/Z8ifv+LtZZma6CilVC3zJcFZBaSKSIqIBAHjgAVV6izAGkQMMBZYZIwxIhIN/Bt4zBiztIZirjMH80r4cushxvVvRYjT7u9wVAPVsXkkb93Zn9k/G0hyXDi/X7CZy59bzDvLMyl1aaKjlFK14ZwJjndMzf3AZ8BW4J/GmM0iMlVERnqrzQBiRSQDeBg4cSv5/UB74AkRWeedGsxT8r7Zno3HwKhe+tZwdfF+1C6WufcM5L27B5DUNJTffbyZK55bzKwVu/UdV0opVcOkvj2zIy0tzaSnp/s7DAAenLOWpRk5rPrtMETE3+GoAGKMYUnGEV78Yjtr9hwjMTqUXwxtz5i+STj1QZJKKeUzEVltjEmrWq5n0jMwxrBsRw6XtIvV5EbVOBFhcGo8H/z8Et66sx9xkcE89uFGhv7fYv65aq8+FVkppS6SJjhn8MPhArLzSxnUvurd7krVHBHh8o4JzJ98CW9O7Ed0aBCPfrCBK1/4mvfT9+LSREcppS6IJjhnsCzjCACXtIvzcySqMRARruiUwIL7B/H6HWlEBDt4ZN4GLntuMW8s2UVhqcvfISqlVIOiCc4ZLN2RQ+uYMFrFhPk7FNWIiAhXdmnGJ7+4lNfuSKNldAhTP9nCJX9axHOffc/h/BJ/h6iUUg2Cw98B1Ecut4fvduZwrT77RvmJiDC8SzOGd2nGmj1Hee2bnfx98Q5e+2YXN/RO5GdDUmifEOnvMJVSqt7SBKcam/YfJ7/ExSXt9fKU8r8+rZvyyvi+ZB4p5PUlO3k/PYu56XsZ1imBe4a0pX9KjA6EV0qpKvQSVTWWVoy/0QHGqv5IjgvnD6O7s+yxoUy5MpW1e49x8/TvGP3yUv694QBuT/165INSSvmTJjjVWL4jh07NI4mLaHAvPleNQGxEMFOu7MDSXw/l6dHdyCsu57731nDF89aA5Lzicn+HqJRSfqcJThUl5W5WZebq3VOq3gsNsnP7wDYs/OXlvDq+D3ERQUz9ZAsD/vdLHnl/Pev2HqO+PchTKaXqio7BqWLNnqOUujx6eUo1GHabMKJbC0Z0a8GmfXnMWrGHj9ft4/3VWXRLbMJtA9owsmdLwoP1n7tSqvHQHpwqlmXkYLcJA9rG+DsUpc5bt8QonrmxOyv+ZxhPj+6Gy234zYcbGfC/C/nd/E18f/C4v0NUSqk6oX/SVbFsxxF6JEURGeL0dyhKXbDIECe3D2zD+AGtWbPnKLO+28Pc9L28891u+rZpym0DWvPj7i0Icdr9HapSStUK7cGpJL+knPVZeQzS8TcqQIgIfdvE8MLNvVjxm2E8fm1ncgvLePif6xn4zEKe/mQLm/fn6VgdpVTA0R6cSlbuysXtMTr+RgWkpuFB3D24LXddmsLyHTm8u2I3M5dnMmPJLjo0i2B070RG9UokMTrU36EqpdRF0wSnkqUZOQQ7bPRp09TfoShVa0SES9rHcUn7OI4WlvHvjQf4aO0+/vzfbfz5v9sYkBLDDb0TuaZ7C6JC9VKtUqphkvrWNZ2WlmbS09P98t0j/vINsRFBzLp7oF++Xyl/2pNTxPx1+5i/dh87jxQS5LBxZecERvdK5PKOCQQ59Iq2Uqr+EZHVxpi0quXag+N1pKCU7w/m88jVHf0dilJ+0To2jAeGpfKLoe3ZkJXHR2v38a/1+/l040Giw5xc270FN/ROpE/rpths+moIpVT95lOCIyIjgJcAO/C6MeZPVdYHAzOBvkAOcLMxJlNEYoF5QD/gLWPM/TUZfE1aviMH0NczKCUi9GwVTc9W0fz22s4syTjCR2v28cGaLGat2ENCZDBXdW3G1V2bM7BtLE679uwopeqfcyY4ImIHXgaGA1nAKhFZYIzZUqnaXcBRY0x7ERkHPAvcDJQAvwO6ead6a9mOI0QGO+ieGOXvUJSqN5x2G1d0TOCKjgkUlLr4csshPtt8kA9W7+Pd7/bQJMTBsM7NuLprM4Z0iCcsSDuFlVL1gy9no/5AhjFmJ4CIzAFGAZUTnFHAk975ecA0ERFjTCGwRETa11zItWNpRg4D2sbi0L9GlapWRLCD0b0TGd07kZJyN9/+cITPNh/ky62H+GjtPoIdNoZ0iOfqrs25snMC0WFB/g5ZKdWI+ZLgJAJ7Ky1nAQPOVMcY4xKRPCAWOOJLECJyD3APQOvWrX3ZpEbtzS1iT24Rdw5KrvPvVqohCnHaGd6lGcO7NMPl9rAyM5fPNh3k8y2H+GLLIetp4CkxXN21OUM7JdAqJszfISulGpl60Z9sjJkOTAfrLqq6/v4T428GtdcH/Cl1vhx2G5e0i+OSdnE8ObIrG7Ly+GzzQT7bfJDfL9jM7xdsJiUunCGpcQxOjedH7WL1vVhKqVrny1lmH9Cq0nKSt6y6Olki4gCisAYbNwhLdxwhLiKY1IQIf4eiVINWeYDyoyM6sTO7gK+3Z/PN9mz+mZ7F28t347QLfds0ZUiHeIakxtOlRRO9K0spVeN8SXBWAakikoKVyIwDbq1SZwEwAVgOjAUWmfr2gJ0zMMawbEcOl7SLRURPskrVpLbxEbSNj+DOQSmUutyszjzK1z9k8832IxUPFoyLCOLS9lbvzuAOcSREhvg7bKVUADhnguMdU3M/8BnWbeJvGGM2i8hUIN0YswCYAbwjIhlALlYSBICIZAJNgCARGQ1cVeUOLL/64XAB2fmlDGqvt4crVZuCHfaKJyj/5ho4nF/Ckh+O8M32bL794Qjz1+0HoGOzSPqlNKV/Siz9k2NoHqUJj1Lq/DX6Jxm/tXQXT/5rC98+eoUOhFTKTzwew5YDx/l6ezYrduWyOjOXwjI3AK1iQumfHEt/b9KTHBumva1KqQr6JOMzWLojh1YxoZrcKOVHNpvQLTGKbolR3HcFuNweth7IZ8WuHFZl5vLVtsN8sCYLgLiIYAakxNAv2Up4OjaPxK5jeJRSVTTqBMfl9vDdzhyu7d7C36EopSpx2G10T4qie1IUdw9uizGGHdkFrNx1lJW7cli5K5d/bzwAQGSwg26JUfRoFUXPpGh6JEWRGB2qvTxKNXKNOsHZtP84+SUuLtHbw5Wq10SE9gmRtE+I5NYB1rOyso4WsSozl9W7j7IhK483luyi3G1dco8ND6JHUpR1R5c36YmNCPZnE5RSdaxRJzhLM6znEOr7p5RqeJKahpHUNIwbeicBUOpy8/2BfDZkHWN9Vh7r9x5j8fZsTgwzTIwOpWerKHokRdO5RRM6N48kPjJYe3qUClCNOsFZviOHTs0jidO/7JRq8IId9opn8NzuLSsodbFpX15F0rMh6xifbjxYsU1MeBAdm0XSqUUknZs3oVOLSFITIgkNsvunEUqpGtNoE5yScjerMnMruruVUoEnItjBwLaxDGx7spf2WFEZ3x/M5/sDx/n+YD5bD+YzZ+Veisutu7ZEICU2nE4tIunUvAmdmkfSPiGCVjFh+uZ0pRqQRpvgrNlzlFKXh0HtdPyNUo1JdFjQaUmPx2PYk1vE9wePs/VAPt8fPM7m/cdP6e1x2ITWsWG0i4+gbXw47eK8n/ERNA3XF4sqVd802gRnWUaO9ULAtjH+DkUp5Wc2m5AcF05yXDgjup28q7Kw1MX2Q/nsyC5kZ3YBO7IL2JldyNfbsilzeyrqNQ1zWk9tjgunXUIEybHhtIoJpXVMGJEhTn80SalGr9EmOEt3HKFHUpSefJRSZxQe7KB366b0bt30lHKX28O+Y8UVCc+JBOirbdm8vzrrlLrRYU5ax4TRKiaMVk3DvPNW8tMyOlQveylVSxplgpNfUs6GrDwmXdbW36EopRogh91Gm9hw2sSGM7TTqevyisvZm1vEntyik59Hi9my/zifbz5YcSs7gE2gRVQoSU1DaRkdSouoEFpEh9KiSQgtokNoGRVKdJhT7/RS6gI0ygRn5a5c3B6j42+UUjUuKtRJlPepzFW5PYZDx0sqkp8TCVDW0WJW7srl0PESXJ5TX58T4rTRIsqb/Hg/m0eFkBAZTEKTEOIjg4mPCCbIoT1BSlXWKBOcpRk5BDts9GnT9NyVlVKqhthtQstoq7em8iDnE9wew5GCUg7klXDgWDH7vZ8Hjlufy3Yc4dDxEjzVvEIwOsxJfEQwCU2shCc+MpiESG8CFBlMbEQQMeFBNA0L0stiqlFolAnOsh1HSEtuSohTn3WhlKo/7DahWZMQmjUJoVer6GrruNwejhSUcTi/hOz8Ug7nl5LtnU6Upe8+SnZ+KaUuT7X7aBLiIDYiuCLhiQ0PIibC+mwaZs03DQsiOtRJdJiTyBCnvu9LNTiNLsE5UlDK9wfzeeTqjv4ORSmlzpvDbqO59zLV2RhjyC91cfi4lfzkFpaRW1hKTmGZd96aso4WsSHrGLmFZaddHqusSYiDqDAn0aFBRIc5aRLqJDrUSZQ3CWoSYiVCkSEO73RyPtRp13FEqs41ugRn+Y4cQF/PoJQKbCJCkxAr8WifEHHO+sYYjpe4KhKhY0Xl5BWXc6yonGPF5RwvLudYUZlVVlzOvmPF5HnrnC0xAqtn6kSyExFsJT5NQhyEBTkID3YQEWwnLMhBRLCDsGA7EcEOwoNOzp9YFxpkJ9Rp1/FGyieNLsFZtuMIkcEOulczAFAppRorEbEGSIc6SYkL93k7YwyFZW7yisspKHGRX1JOfomL497PfG9ZQenJ+eMlLvYdK6GozEVhqZvCUlfFk6R94bAJoU47oUF2woLshDitTysBshKhMKedEKeNEKed4BPzDqtuRbnDVrEc7F0X7LAR7LAR5LDKghw2vTzXQPmU4IjICOAlwA68boz5U5X1wcBMoC+QA9xsjMn0rvsNcBfgBh4wxnxWY9FfgKUZOQxoG4tDB9kppdRFExEigq0elovh9hiKylwUlbkpKHVRVOr9LHNRUGolQsXlborLrGSoqMxNcZm7Yr7E+5lbWOydd1Fc5qbU5TnjWCRfOWziTXhOTXyCHTacdhtBdqvcaRecdhtOh7fMbsPpkIo6Tu/ksEtFXYfdhtMmlcptOGyC02HDabPKHDbB4S0/sWy3VV4+dd5mw/oUGvWlwXP+RoqIHXgZGA5kAatEZIExZkulancBR40x7UVkHPAscLOIdAHGAV2BlsCXItLBGON7ql6DTtySeeegZH98vVJKqTOwLmNZ43ia1fC+PR5DmdtDSbmbknLvp6vSfLmVCJWUuylzeShzeygtr/zprrLsoczlodTlpsxtKPdue7zEKi93W3XKXaZi/sR+zdmv5tU4u02w27xJkQh2b4JkE++n7eSnXay6thOfNsEuVJQ57NanrVI9W6X1Npt3WQQRwW4Dm1hvDLihd1LdNhzfenD6AxnGmJ0AIjIHGAVUTnBGAU965+cB08RKG0cBc4wxpcAuEcnw7m95zYR/fk6Ov9Hn3yilVGNhswkhNnu9uHPW7bGSHpfH4PImPy63weW2kjCXx3Ny3m3VcXlMxXZuj6HcY3B767k8pmJfVh2Dx1j1Xd56bg9Wfe9+TkynLBuDxzt/Ynu3oaLMbQwl5dZ3GHOiPnjMyfrGYJVXWuf2GKJCndC77n/WviQ4icDeSstZwIAz1THGuEQkD4j1ln9XZdvEql8gIvcA9wC0bl17b/fu0DySnw1OoUOzcw+4U0oppWqa1aPi/0SrMagXg4yNMdOB6QBpaWm11oHXq1X0GZ8toZRSSqnA4ctI231Aq0rLSd6yauuIiAOIwhps7Mu2SimllFI1ypcEZxWQKiIpIhKENWh4QZU6C4AJ3vmxwCJjjPGWjxORYBFJAVKBlTUTulJKKaVU9c55ico7puZ+4DOs28TfMMZsFpGpQLoxZgEwA3jHO4g4FysJwlvvn1gDkl3Aff66g0oppZRSjYeYur5n7RxEJBvYXYtfEQccqcX911eNsd3a5sZB29x4NMZ2a5vPrY0xJr5qYb1LcGqbiKQbY9L8HUdda4zt1jY3DtrmxqMxtlvbfOH0cb5KKaWUCjia4CillFIq4DTGBGe6vwPwk8bYbm1z46BtbjwaY7u1zReo0Y3BUUoppVTga4w9OEoppZQKcJrgKKWUUirgBGyCIyIjRGSbiGSIyGPVrA8Wkbne9StEJNkPYdYoH9o8UUSyRWSdd7rbH3HWJBF5Q0QOi8imM6wXEfmr92eyQUT61HWMNc2HNl8uInmVjvMTdR1jTRORViLylYhsEZHNIvJgNXUC6lj72OZAPNYhIrJSRNZ72/1UNXUC6vztY5sD7vwNICJ2EVkrIp9Us+7ijrMxJuAmrCcu7wDaAkHAeqBLlTqTgVe98+OAuf6Ouw7aPBGY5u9Ya7jdQ4A+wKYzrP8x8B9AgIHACn/HXAdtvhz4xN9x1nCbWwB9vPORwPZqfr8D6lj72OZAPNYCRHjnncAKYGCVOoF2/valzQF3/va262Hgvep+jy/2OAdqD05/IMMYs9MYUwbMAUZVqTMKeNs7Pw8YJiJShzHWNF/aHHCMMd9gvR7kTEYBM43lOyBaRFrUTXS1w4c2BxxjzAFjzBrvfD6wFUisUi2gjrWPbQ443uNX4F10eqeqd8ME1PnbxzYHHBFJAq4FXj9DlYs6zoGa4CQCeystZ3H6iaGijjHGBeQBsXUSXe3wpc0AY7zd9/NEpFU16wONrz+XQPMjb3f3f0Skq7+DqUnebureWH/lVhawx/osbYYAPNbeyxbrgMPAF8aYMx7rADl/+9JmCLzz91+ARwHPGdZf1HEO1ARHVe9fQLIxpgfwBSczYxVY1mC9m6Un8Ddgvn/DqTkiEgF8AEwxxhz3dzx14RxtDshjbYxxG2N6AUlAfxHp5ueQap0PbQ6o87eIXAccNsasrq3vCNQEZx9QObtN8pZVW0dEHEAUkFMn0dWOc7bZGJNjjCn1Lr4O9K2j2PzJl9+FgGKMOX6iu9sY8yngFJE4P4d10UTEifUf/SxjzIfVVAm4Y32uNgfqsT7BGHMM+AoYUWVVoJ2/K5ypzQF4/h4EjBSRTKwhFUNF5N0qdS7qOAdqgrMKSBWRFBEJwhqctKBKnQXABO/8WGCR8Y5kaqDO2eYq4xFGYl3TD3QLgDu8d9gMBPKMMQf8HVRtEpHmJ65Ti0h/rH/nDfrk723PDGCrMeaFM1QLqGPtS5sD9FjHi0i0dz4UGA58X6VaQJ2/fWlzoJ2/jTG/McYkGWOSsf6/WmSMGV+l2kUdZ0eNRFrPGGNcInI/8BnW3UVvGGM2i8hUIN0YswDrxPGOiGRgDdgc57+IL56PbX5AREYCLqw2T/RbwDVERGZj3UkSJyJZwO+xBuhhjHkV+BTr7poMoAi40z+R1hwf2jwW+LmIuIBiYFxDPvl7DQJuBzZ6xykA/A/QGgL2WPvS5kA81i2At0XEjpWw/dMY80kgn7/xrc0Bd/6uTk0eZ31Vg1JKKaUCTqBeolJKKaVUI6YJjlJKKaUCjiY4SimllAo4muAopZRSKuBogqOUUkqpgKMJjlJKKaUCjiY4SimllAo4/w/Y5WUkEVEFNAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "\n", "params = np.array((2,))\n", "dt = 0.1\n", "y0 = 1\n", "t = np.arange(0,4,dt)\n", "tru_y = np.exp( -params[0]*t)\n", "(t,y) = forwardEuler1D( degradation_model, t, y0=y0, params=params)\n", "\n", "plt.figure( figsize = (8,4) )\n", "plt.subplot(2,1,1)\n", "plt.plot( t, tru_y)\n", "plt.plot( t, y.T );\n", "plt.legend( ['True Y','Approximate Y'])\n", "\n", "plt.subplot(2,1,2)\n", "plt.plot( t, tru_y - y[0,:].T);\n", "plt.legend( ['True - Approximate'])\n", "\n", "plt.tight_layout();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Evaluate the effect of step size on the correctness\n", "\n", "For values of $\\Delta t$ = 0.001, 0.01, 0.1, and 0.2, simulate the model and plot the difference between the forward Euler solution and the true solution." ] }, { "cell_type": "code", "execution_count": 130, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD4CAYAAADiry33AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy86wFpkAAAACXBIWXMAAAsTAAALEwEAmpwYAABAp0lEQVR4nO3dd3wU1drA8d/Zmt4TAiSU0FHKi4igWBEEpAioYAOvAnbB+nrV13a9Xr16bderYgevivdag1IURWwgTRRpitRgIL1nk+zuef/YTQghgYRsdpLd5/v5zGdnZ87OeXYgz8yeOXNGaa0RQggRuExGByCEEKJlSaIXQogAJ4leCCECnCR6IYQIcJLohRAiwFmMDqCuhIQE3aVLF6PDEEKINmX9+vU5WuvE+ta1ukTfpUsX1q1bZ3QYQgjRpiil9jS0TppuhBAiwEmiF0KIACeJXgghAlyra6MXQgSmqqoqMjIycDgcRofSpoWEhJCSkoLVam30ZyTRCyH8IiMjg8jISLp06YJSyuhw2iStNbm5uWRkZNC1a9dGf06aboQQfuFwOIiPj5ck3wxKKeLj45v8q0gSvRDCbyTJN9/x7ENJ9D5QtmEDJStXIkM+CyFaI0n0zaRdLvbPmcu+a65l38xZVOzaZXRIQohGeuCBB3jiiSd44403+OOPPxr9Oa01N998M927d6d///5s2LCh3nLr16+nX79+dO/enZtvvrnmZDAvL4+RI0fSo0cPRo4cSX5+PgDbtm1j2LBh2O12nnjiieZ/QS9J9M1Uuno1zuxsosaNo/ynn9g5YSJZTz2Nu6zM6NCEEI3U1ES/ZMkSfvvtN3777TdeeuklrrvuunrLXXfddbz88ss1ZZcuXQrAo48+yogRI/jtt98YMWIEjz76KABxcXE8++yz3H777c3/UrVIom+movR0TJGRtP/rw3RbuoTosWPJnTeP388fR9Fnn0lzjhCtzF//+ld69uzJ8OHD2b59OwDr1q3jsssuY+DAgZSXlx9zGx9//DHTp09HKcXQoUMpKCggMzPzsDKZmZkUFRUxdOhQlFJMnz6djz76qObzM2bMAGDGjBk1y5OSkjj55JOb1HWyMaR7ZTO4y8oo+nw50eePxWS3Y7Lb6fDYo8RcfBEHHvoL+2+eQ/jw4STfew82GahNiBoPLtrMlj+KfLrNvh2iuH/8CUcts379ehYuXMjGjRtxOp0MGjSIk046icGDB/PEE08wePBgAG655RZWrFhxxOenTZvGXXfdxf79+0lNTa1ZnpKSwv79+2nfvn3Nsv3795OSknJEGYCDBw/WlE1OTubgwYPH/8UbQRJ9MxR/8QW6rIyo8eMPWx520kl0ff898t9+h+xnn2Xn+AnEXX0VCddcgyk01KBohRDffPMNkyZNIiwsDIAJEybUW+6pp57yW0xKqRbvjSSJvhkK0xdh6dCeMO9ZQG3KYiFu+hVEjRlN1hNPkPviPArT00m++24iRoyQbmYiqB3rzNtoxzqj79ixI/v27atZnpGRQceOHQ8r27FjRzIyMuot065dOzIzM2nfvj2ZmZkkJSW10DfxkDb64+TMzqb0u++IHj8BZWp4N1oSE+nw2GN0fnMB5vAIMm68iX3XXEPlngZHFBVCtJAzzjiDjz76iPLycoqLi1m0aBEAkZGRFBcX15R76qmn2Lhx4xHTXXfdBXh+CSxYsACtNatXryY6OvqwZhuA9u3bExUVxerVq9Fas2DBAiZOnFjz+fnz5wMwf/78muUtRmvdqqaTTjpJtwU5r7+ut/TqrR07djT6M+7KSp37xht626CT9NYT++msZ57RrrKyFoxSiNZjy5YtRoegtdb64Ycf1j169NCnnXaavuSSS/Tjjz+u33vvPd2zZ089YMAAXdaIv0m3262vv/56nZaWpk888US9du3amnUDBgyomV+7dq0+4YQTdFpamr7hhhu02+3WWmudk5OjzznnHN29e3c9YsQInZubq7XWOjMzU3fs2FFHRkbq6Oho3bFjR11YWHhE/fXtS2CdbiCvKt3KeoUMHjxYt4UHj+yaPAWUouv77zX5s1VZWWQ9/gRFixZh7dCB5AfuJ+KMM1ogSiFaj61bt9KnTx+jwwgI9e1LpdR6rfWR7chI081xqdixA8eWLURPrP9CzrFYk5Lo+Pjf6bRgPqbwMPZdf4M05QghWowk+uNQ+HE6mM1EjR3brO2EDxlC6quvoqxWsp/9p4+iE0KIw0mibyLtdlP4ySeEn3YqloSEZm/PmpRE3PTpFH36KY4tW3wQoRBCHE4SfROVrV2HMzOT6Am+u0oef/VVmKKjyXrqaZ9tUwghqkmib6LC9I8xhYUROeIcn23THBVFwuzZlH7zDaU/rPHZdoUQAiTRN4nb4aB42WdEjhrl8ztcYy+7FEtyMllP/kPGxxFC+JQk+iYoWbECd0kJ0RPGH7twE5lCQki88QYcP/1M8fLlPt++EOJIRg1T/N///pcTTjgBk8mEP7qTS6JvgsKP07EkJRF2yiktsv3oCy7AlpZG9tPPoJ3OFqlDCHEkfw9TfOKJJ/LBBx9whp/un5FE30jOvDxKvv2WqPHjUGZzi9ShLBYS586h8vffKfz44xapQ4hg1xqGKe7Tpw+9evXy+XdrSKMGNVNKjQaeAczAK1rrR+ustwMLgJOAXGCq1nq3UsoKvAIM8ta1QGv9Nx/G7zdFi5eA00l0A6Pd+UrkyJGE9O9P9j+fI2rcOEx2e4vWJ4QhltwFBzb5dpvJ/WDMo0ct0lqGKfa3YyZ6pZQZ+BcwEsgA1iql0rXWtTt9Xw3ka627K6WmAY8BU4GLALvWup9SKgzYopR6R2u929dfpKUVpqdj79WLkBY+CiulSLr1VvZeeSX5b71N/FV/atH6hAgmrXGYYn9ozBn9EGCH1nongFJqITARqJ3oJwIPeOffA55TnnF4NRCulLIAoUAl4NunDfhBxa5dOH7+maQ77vBLfeFDTyF8+HBy580j5qILMUdG+qVeIfzmGGfeRmvpYYr9rTFt9B2BfbXeZ3iX1VtGa+0ECoF4PEm/FMgE9gJPaK3z6laglJqtlFqnlFqXnZ3d5C/R0ooWLQKliBp3vt/qTLr1FlyFheS++qrf6hQi0LWWYYr9raUvxg4BXEAHoCtwm1IqrW4hrfVLWuvBWuvBiYmJLRxS02itKUxfRPiwoVjbtfNbvSF9+xI1dix58xdQlZXlt3qFCGSDBg1i6tSpDBgwgDFjxnDyyScDcOWVV3Lttdc2+mLs2LFjSUtLo3v37syaNYvnn3++Zt3AgQNr5p9//nlmzpxJ9+7d6datG2PGjAHgww8/JCUlhVWrVnH++edz3nnn+faL1tXQ+MXVEzAMWFbr/Z+BP9cpswwY5p23ADmAwtO2f0Wtcq8BFx+tvtY2Hn3p+vV6S6/eOv/DD/1ed8Xu3XrLCSfqzAcf9HvdQvhaaxmPPhA0dTz6xpzRrwV6KKW6KqVswDQgvU6ZdGCGd/5C4EtvxXuBcwCUUuHAUGBbk49GBir8OB0VEkLkuSP9Xretc2diLrqQ/P/8V4YxFkIct2Mmeu1pc78Rz1n7VuA/WuvNSqmHlFLVl6xfBeKVUjuAW4G7vMv/BUQopTbjOWC8rrX+2ddfoqW4KyspWrqUyHPPxRwRbkgMCdddJ8MYCyGapVH96LXWi4HFdZbdV2vegacrZd3PldS3vK0oWbkSd2HhcT9gxBeqhzHOnTeP+KuvIqRvX8NiEUK0TXJn7FEUpadjTkggfNgwQ+OIn3k1ZhnGWAhxnCTRN8BVUEDJVyuJPn8sytKoHz4txhwZSbwMYyyEOE6S6BtQtHQZuqqKqPHGNdvUJsMYCyGOlyT6BhSmp2Pr1o2QE1pHm7gMYyyE77X0MMX33HMPqampRERE+Crk4yKJvh6V+/ZRvmED0RMm4BnJoXWQYYyFaBktNUzx+PHjWbPG+OZWSfT1KPTeFh3txyEPGkOGMRai+fw1TDHA0KFDjxgawQjGXmVshbTWFH2cTtjJJ2M1aACio5FhjEUgeGzNY2zL8+29k73jevO/Q/73qGX8OUxxayKJvg7Hpk1U7tlD/KyZRodSLxnGWIjjJ8MUC8A75IHNRuSoUUaH0iAZxli0dcc68zaaL4Ypbk2kjb4WXVVF0eLFRJxzDuaoKKPDOSoZxliIpvPnMMWtiST6Wkq+/RZXfn6LPy7QF2oPY+zMyTE6HCHaBH8PU3znnXeSkpJCWVkZKSkpPPDAA77+So2iWtvNN4MHD9br1q0zpO6MW26hbNVqeny9EmWzGRJDU1Ts2sXOsecTP2sWSbfeYnQ4QhzV1q1b6dOnj9FhBIT69qVSar3WenB95eWM3stVXEzJlyuIGju2TSR5AHvXrkSOGkX+O+/gKikxOhwhRCslid6r+LPP0BUVRE8Yb3QoTRI/cybu4mIK3n3X6FCEEK2UJHqvwo/TsXXuTMiAAUaH0iSh/U4k/NRh5L0xH3dlpdHhCCFaIUn0eHrblK1bR+Soka1qyIPGip85E2d2ttwtK4SolyR6wJmbC2431o4pRodyXMKGDSPkhBPIe+VVtMtldDhCiFZGEj3gzM4GwJKUaHAkx0cpRfysWVTu2UPx8i+MDkcI0cpIoqdWok9sm4keIHLkudg6dyb35ZdlvHohGqklhykuKyvj/PPPp3fv3pxwwgk1N1sZQRI94Mxq+4lemc3EXX0Vjl9+oWz1aqPDEaJNaalhim+//Xa2bdvGjz/+yHfffceSJUt8FXKTSKLHe0avFJb4eKNDaZboiRMxJyaQ+/IrRociRKvlr2GKw8LCOPvsswGw2WwMGjSIjIwM33+hRpBBzfAkenNsLMpqNTqUZjHZ7cTPmEHWE/+g/JfNhJ54gtEhCVGvA488QsVW3w5TbO/Tm+S77z5qGaOGKS4oKGDRokXMmTOnGd/w+Emix5Po23KzTW0x06aRM+8lcl95hZSnA2uoVSGay4hhip1OJ5dccgk333wzaWlpPttuU0iiJ7ASvTkigthp08h99VUqd+/G1qWL0SEJcYRjnXkbzZfDFM+ePZsePXowd+7clgr3mKSNnsBK9ABx069AWSzkvva60aEI0ar4e5jie++9l8LCQp5++mm/fL+GBH2i1243ztzcgEr0lsREoidPovDDD6nKyjI6HCFaDX8OU5yRkcFf//pXtmzZwqBBgxg4cCCvvGJMR4mgH6bYmZvLb6cNp9099xB3xeV+q7elVe7dy++jxxB/9VUk3Xab0eEIIcMU+5AMU9xEgXCzVH1snToRNfo88t9ZiKvWT1IhRPCRRN/Ghz84mviZM3GXlJC/cKHRoQghDCSJPgDuim1ISN++hJ92GnnzF+CuqDA6HCFkeA4fOJ59KIk+QJtuqsXPmoUrJ4fCDz8yOhQR5EJCQsjNzZVk3wxaa3JzcwkJCWnS54K+H70zOxtTZCSmJu64tiLslCGE9OtH7muvEXPRhSiz2eiQRJBKSUkhIyODbO/JlTg+ISEhpKQ0bUh1SfTZ2ViSkowOo8V4hjCeyf6b51D82WdEjRljdEgiSFmtVrp27Wp0GEFJmm4C7Gap+kSOGIGtSxdyZAhjIYKSJPogSPTKbCZ+5tVUbNlK6fffGx2OEMLPGpXolVKjlVLblVI7lFJHjJ6vlLIrpd71rv9BKdWl1rr+SqlVSqnNSqlNSqlW0xiutcaZlRXwiR4gasIELElJMoSxEEHomIleKWUG/gWMAfoClyil+tYpdjWQr7XuDjwFPOb9rAX4N3Ct1voE4CygymfRN5O7qAhdWRkUid5ksxE3YwZlq1dTvmmT0eEIIfyoMWf0Q4AdWuudWutKYCEwsU6ZicB87/x7wAillAJGAT9rrX8C0Frnaq1bzdOrA71rZV0xUy/GFBUlZ/VCBJnGJPqOwL5a7zO8y+oto7V2AoVAPNAT0EqpZUqpDUqpO+urQCk1Wym1Tim1zp9dr4It0ZsjIoi99BKKP/+cip27jA5HCOEnLX0x1gIMBy7zvk5SSo2oW0hr/ZLWerDWenCiH5NusCV6gLgrrkDZbOS9/prRoQgh/KQxiX4/kFrrfYp3Wb1lvO3y0UAunrP/r7XWOVrrMmAxMKi5QftKII9z0xBLfDwxUyZT8NHHVB08aHQ4Qgg/aEyiXwv0UEp1VUrZgGlAep0y6cAM7/yFwJfa02F7GdBPKRXmPQCcCWzxTejN58zKRoWGYgoPNzoUv4r705/A5SJv/gKjQxFC+MExE723zf1GPEl7K/AfrfVmpdRDSqnqBy6+CsQrpXYAtwJ3eT+bDzyJ52CxEdigtf7U59/iOFX3ofdcNw4ettRUosaMoWDhQlyFhUaHI4RoYY0aAkFrvRhPs0vtZffVmncAFzXw2X/j6WLZ6gTDzVINiZ81k6JPPyXv3/8m8YYbjA5HCNGCgvrO2GBO9CG9exMxYgR5b8yXs3ohApwk+iBN9ACJN92Iu7iYvPnzj11YCNFmBW2id5eV4S4tDepEH9K7N5HnnUfe/AU48/ONDkcI0UKCNtEHYx/6+iTeeAPusjLyXnvd6FCEEC1EEn2QJ3p7jx5EjR1L3ltv4czNNTocIUQLkEQf5IkeIOGG69EOB7mvyt2yQgQiSfRBdFdsQ+xpaUSPH0f+22/X7BchROAI7kRvtWKOiTE6lFYh4frr0VVV5Lz8stGhCCF8LHgTfVY2loSEoLsrtiG2zp2JvmAiBQvflTFwhAgwwZvos7Ol2aaOhOuuQ7vd5M57yehQhBA+FNyJXi7EHsaWkkLM5MkU/Pe/VP3xh9HhCCF8RBK9OEzCtdcAkPPiPIMjEUL4SlAmendlJa6CAkn09bB26EDMRRdR8MEHVGZkGB2OEMIHgjLRu6QP/VHFXzMbZTKR8/wLRocihPCBoEz0crPU0VnbtSNm2lQKP/6Yyt27jQ5HCNFMQZnoqyTRH1PCrFkoq5WcF+SsXoi2LigTvZzRH5slMZHYSy+lcNEnVOzcaXQ4QohmCN5EbzJhiY83OpRWLX7m1aiQEHKe+5fRoQghmiFoE705Pg5lNhsdSqtmiYsj7vLLKVqyBMevvxodjhDiOAVtopdmm8aJ+9OVmMLC5KxeiDZMEr04KktsLHEzplP82Wc4tm41OhwhxHGQRC+OKW7GDEyRkWTLWb0QbVLQJXrtcuHKzZNE3wTm6Gji/nQlJV98QfmmX4wORwjRREGX6J25ueB2S6Jvorjp0zFFR5P93D+NDkUI0UTBl+ilD/1xMUdEEH/VVZSu/JryjRuNDkcI0QRBm+itkuibLO7yyzDHxpL9z+eMDkUI0QRBm+jljL7pTOHhxM+cSel331G2fr3R4QghGiloE71ZEv1xib30EswJCWQ/K231QrQVQZnozdHRmGw2o0Npk0yhoSTMmknZDz9QuvoHo8MRQjRCUCZ6eVZs88RMnYqlXTuy/v53tMtldDhCiGMIzkSfmGR0GG2aKSSEpDvvwLFlC/kLFxodjhDiGII00csZfXNFjR1L2LChZD/1dM11DyFE6xRUiV5rjSs7R5pufEApRfJ996ErKjj498eNDkcIcRRBlehdBQXoqio5o/cRe9euxM+aSdGiRZSuXm10OEKIBgRVopc+9L4XP3s21tRUDjz4EO7KSqPDEULUo1GJXik1Wim1XSm1Qyl1Vz3r7Uqpd73rf1BKdamzvpNSqkQpdbuP4j4uzixJ9L5mCgkh+f/upXLXLvJee93ocIQQ9ThmoldKmYF/AWOAvsAlSqm+dYpdDeRrrbsDTwGP1Vn/JLCk+eE2j5zRt4yIM84gctQocl54gcqMDKPDEULUYWlEmSHADq31TgCl1EJgIrClVpmJwAPe+feA55RSSmutlVIXALuAUl8FfbxaOtFrrckuz2Zf8b7DpuyybPon9ue8LufRJ64PSqkWqd9I7f58FyXffsvBvzxMyosvBOR3FKKtakyi7wjsq/U+AziloTJaa6dSqhCIV0o5gP8FRgINNtsopWYDswE6derU6OCbypmdjSk8HFNYmO+26Xby/MbnWbFvBRnFGThcjpp1ZmUmOTyZuJA4FmxewGu/vEZKRArndTmPUV1GBVTSt7ZvT+KNN5L1979T8sUXRJ57rtEhCSG8GpPom+MB4CmtdcnREprW+iXgJYDBgwfrlgrG133oSypLuOPrO/h2/7ec2uFUTu1wKqmRqTVT+4j2WE1WAAocBXy570s+2/0Zb2x+g1d/eZXUyFRGdR7FuLRxdI/t7rO4jBJ3xeUUfvghB/76COHDhmEKDzc6JCEEjUv0+4HUWu9TvMvqK5OhlLIA0UAunjP/C5VSfwdiALdSyqG1NmScW18m+j9K/uCGL25gV+Eu7ht2Hxf1vOio5WNCYpjcYzKTe0wm35HPl3u/5LM9nqT/2i+vcUnvS7h50M2EW9tuclRWK8kPPsCeSy8j54UXSLrd0GvvQgivxvS6WQv0UEp1VUrZgGlAep0y6cAM7/yFwJfa43StdRetdRfgaeARo5I8+C7R/5z9M5d+eikHSw/ywrkvHDPJ1xUbEsuUnlOYN3IeKy5ewbTe03hn2ztM/GgiX+79stnxGSls0CCip0wm9435OH791ehwhBA0ItFrrZ3AjcAyYCvwH631ZqXUQ0qpCd5ir+Jpk98B3Aoc0QXTaFprnyT6pbuXctWyqwi1hPLvsf9mWIdhzdpebEgsd59yN2+OfZMoexRzVsxh7oq5HCw92KztGinp9tsxh4dz4KGH0LrFWuKEEI3UqH70WuvFWuueWutuWuu/epfdp7VO9847tNYXaa27a62HVPfQqbONB7TWT/g2/MZzl5aiy8uPe/gDrTXzfprHHSvvoG98X946/y3SYtJ8Ft+AxAG8O+5d5gyaw7f7v2XixxN5Z9s7uNxtb3RIS2wsibffRvm69RR+9LHR4QgR9ILmztjm3CxV5arinm/v4bmNzzEubRyvjHqFuJA4X4eI1WRlZr+ZfDjhQ/on9OeRHx5h+pLp7Cw44rjZ6sVMmULowIFkPf44roICo8MRIqgFT6JvRh/6Vza9wqKdi7hh4A08MvwRbOaWfWhJalQq80bO45Hhj7CveB+XfHoJX+z5okXr9DVlMpH8wP24CgvJeuppo8MRIqhJoj+GPEceb2x+g5GdR3LtgGv91u9dKcX4buP57/j/0i2mG3O/mss/f/wnbu32S/2+ENK7N3GXX07Bf/5D+U8/GR2OEEFLEv0xvPzzyzhcDm4ceGNLhHVM7cLb8fro15nUfRIv/fwSN315E0WVRYbEcjwSbroJS2IimQ8+iHY6jQ5HiKAUVIle2WyYoqIa/ZnMkkze3f4uE7tN9OmF16aym+08eOqD3HPKPXy//3su/fTSNtNub44Ip93df6Ziy1by337H6HCECEpBlegtiYlNanp5/qfnUSiuG3BdC0bWOEoppvWexivnvUJxZTGXLr6UL/a2jXb7yPPOI3z4cLKfeYaqrCyjwxEi6ARdom+snQU7Sf89nam9p9I+on0LRtY0J7U7iXfHvUvXqK7MXTGX5358rtW32yulSP6/e9FVVWQ9WndgUyFESwuuRJ/U+IeC//PHfxJqCWVmv5ktGNXxSQ5P5o0xb3BB9wuY9/M87vrmLipdrfuhH7bOnYmfPZuixYsp+e47o8MRIqgEV6Jv5Bn9puxNLN+7nBknzGiR/vK+YDfbeejUh5g7aC5Ldi3h2uXXtvqLtPGzZmLt3ImDD/1FnkYlhB8FRaJ3Oxy4i4oaneif+fEZ4kLimN53egtH1jxKKa7udzV/O/1v/Jj1IzOWzOBA6QGjw2qQyW4n+f/uo3LPHnJfecXocIQIGkGR6J05OUDjulau+mMVP2T+wKx+s9rMSJLj0sbxwrkvcKD0AJctvoxf81vvYGIRw08jauxYcp5/gbINPxodjhBBITgSffXwB8cY50ZrzTMbnqF9eHsu7nWxP0LzmaHth/LG6DcAmLFkBj9k/mBsQEeR/MD9WDt0YP/cuTUHYSFEywmORN/Im6WW713O5tzNXD/w+hYf5qAl9IrrxVtj3yI5PJlrl1/Lpzs/NTqkepmjokh59hlcRUXsv/U2uZFKiBYWHIne23f7aIne6Xby7IZn6RbdjfFp4/0Vms8lhyczf8x8BiYO5K5v7uL1X15vlUMFh/TuTfID91O2Zg3ZTz9tdDhCBLTgSPTZ2WCxYI6NbbBM+u/p7C7azU2DbsJsMvsxOt+LskUxb+Q8RncZzZPrn+TxdY+3yr72MRdcQMy0qeS+8ipFn39udDhCBKyWfmZsq+DMzsYSH48y1X9cc7ldvPDTC/RP6M85qef4ObqWYTPbeOyMx4gPjefNLW+SW57Lw6c9jNVsNTq0w7S7+24cW7aSedefsXfvjr1rV6NDEiLgBM0Z/dGabX7J/YUDpQe4tM+lfhud0h9MysT/nvy/zBk0h8W7FnPTlzdRVlVmdFiHMdlspDz9FMpqZf/Nc3CXta74hAgEkuiBlftWYlZmhncc7seo/EMpxcx+M3no1IdYlbmKq5ddTZ4jz+iwDmPt0IEO/3iCih07yLzv/lZ5TUGItkwSPfB1xtcMTBpItD3aj1H516Qek3j6rKf5reA3ZiyZwR8lfxgd0mEiTjuNxJtvouiTT8h/622jwxEioAR8otdVVbjy8hpM9AdKD7A9fztnppzp58j87+xOZ/PSyJfIdeRyxeIrWt2NVfHXXEPEWWdx8LHHKPtRbqYSwlcCPtE7c3OBhrtWfp3xNUBQJHqAQe0GMX/0fACuXHIl6w6sMziiQ5TJRIfHHsWanMz+ubfU/NsJIZon8BN99tHvil2ZsZKUiBS6RgdPb48esT14c+ybxIfGc83n1/DZ7s+MDqmGOTraczNVQYHcTCWEjwRPoq/njL7cWc4PmT9wZuqZAdXbpjE6RHTgzTFv0je+L7evvJ23tr5ldEg1Qvr0Ifn++yn74Qeyn3nG6HCEaPMCP9FnNZzo12SuocJVwRkpZ/g7rFYhJiSGl0e9zNmpZ/Pomkd5ct2TrebGqpjJk4i5+GJyX36F4uXLjQ5HiDYt8BN9djYohSU+/oh1KzNWEmYJY3C7wQZE1jqEWEJ48qwnmdprKq9vfr1VPcSk3T13E3Liifxx15+p3L3b6HCEaLOCItGbY2NR1sPvCNVaszJjJad2OLVNDmDmS2aTmXtOuYc5g+awZNcSrl9+PcWVxUaHhcluJ+WZp1FmMxk33Sw3UwlxnIIi0dfXbLM9fztZZVlB22xTV/WNVY8Mf4T1B9dz5dIrOVh60OiwsHbsSId//MNzM9UDD8jNVEIch6BN9Cv3rQTg9JTT/R1Sqza+23j+de6/yCjO4PIll7M9b7vRIREx/DQSbrqRovRF5L/zjtHhCNHmBEeir+eh4F9nfE2/hH4khCYYEFXrdmqHU5k/Zj5u7Wb6kuk19xoYKeHaa4k480wO/u1RyjduNDocIdqUgE702u3GmZt7xBl9TnkOm3I2Bc1NUsejd1xv3jn/HbpEd+GmL2/izS1vGtpsokwmOvz9Mazt2pExZy6VGRmGxSJEWxPQid6Vnw9O5xGJ/tv936LRnJkqif5oksKSeGP0G5yTeg5/X/t3/rL6L1S5qwyLxxwdTcq/nkM7HOy5YjqVe/YYFosQbUlAJ/qGbpb6OuNrksKS6BXby4iw2pRQSyj/OOsfzOw3k//++l+uW34dhRWFhsUT0qsXnea/UZPsK3buMiwWIdqKoEv0Va4qvtv/HWemBN/dsMfLpEzMGTSHh097mPUH13P54svZW7TXsHhCevem84L5aLebPdOnU/Hbb4bFIkRbENiJPuvIcW7WHVxHmbNM2uePw8TuE3l55MvkV+Rz6eJLWXtgrWGx2Hv0oPOC+Sil2DN9Bo5t2wyLRYjWLrATfT1n9F9nfI3dbGdI+yFGhdWmDU4ezNtj3yYuJI5Zn81iweYFhl2ktael0fnNBSi7nb0zrqR882ZD4hCitWtUoldKjVZKbVdK7VBK3VXPertS6l3v+h+UUl28y0cqpdYrpTZ5X/36QFZnVhamqChMdjvguRv2q31fcUr7Uwi1hPozlIDSKaoTb499m7NSz+LxdY9z59d3GvaIQluXLnT+95uYwsPZ+6erKP/5Z0PiEKI1O2aiV0qZgX8BY4C+wCVKqb51il0N5GutuwNPAY95l+cA47XW/YAZwJu+Crwx6t4statoFxklGdJs4wMRtgieOusp5gyaw2d7PuOyxZexp8iYXjC2lBQ6v7kAc3Q0e/90FWUb5KElQtTWmDP6IcAOrfVOrXUlsBCYWKfMRGC+d/49YIRSSmmtf9RaVz+zbjMQqpSy+yLwxqib6L/e57nxR4Y98I3qYRNePPdFcspzmPbJNFbsXWFILNaOHen85gIsCQnsmzmTsrXGXT8QorVpTKLvCOyr9T7Du6zeMlprJ1AI1B0ucgqwQWtdcXyhNl3dRL8yYyW9YnuRHJ7srxCCwrAOw3h33Lt0iurEzStu5tkNz+Jyu/wehzU5mU5vLsCSnMze2ddQunq132MQojXyy8VYpdQJeJpzrmlg/Wyl1Dql1Lps7wXU5tJaH5boCysK+THrRzmbbyEdIjqwYMwCJveYzMubXub6L64nz5Hn9zisSUl0XjAfW0oK+665lpJvv/N7DEK0No1J9PuB1FrvU7zL6i2jlLIA0UCu930K8CEwXWv9e30VaK1f0loP1loPTmzg2a5N5S4qQldW1iT67//4Hpd2yd2wLchutvPgqQ/ywLAHWHdgHVPSp/D9/u/9HoclIYFOC+ZjS0sj47rrKP7qK7/HIERr0phEvxbooZTqqpSyAdOA9Dpl0vFcbAW4EPhSa62VUjHAp8BdWmu/nlrV7Vq56o9VRNujOTH+RH+GEZSm9JzC2+e/TbQtmmuWX8Pjax/3+8NMLLGxdH79New9e5Jx080Uf/GFX+sXojU5ZqL3trnfCCwDtgL/0VpvVko9pJSa4C32KhCvlNoB3ApUd8G8EegO3KeU2uidjhxKsgXUTfTb8rbRN64vZpPZH9UHvV5xvVg4biFTe01lwZYFXLb4MnYW7vRrDOaYGDq9/hohffuQMWcuRUuX+bV+IVqLRrXRa60Xa617aq27aa3/6l12n9Y63Tvv0FpfpLXurrUeorXe6V3+sNY6XGs9sNaU1XJf55Daid7pdvJ7we/0ipOxbfwpxBLCvUPv5dmzn+VA6QGmLprKe7++59cbrMxRUXR69VVC+/dn/223UfjJp36rW4jWImDvjK1J9EmJ7C7cTaW7kp6xPQ2OKjid3els3p/wPgOTBvLgqge55atbKHAU+K1+c0QEnV5+ibBBg/jjzjvJefFFtMv/vYKEMErgJvqsbFRoKKbwcLbne56SJGf0xkkKS2LeyHncPvh2VmasZHL6ZL7Y6792c1N4OKkvzSNqzBiyn36GvX+6iqoDB/xWvxBGCtxE7+1aqZRie/52rCYrXaO7Gh1WUDMpEzNOmMHbY98mPjSeuSvmcttXt5FTnuOf+kND6fDE47T/298o/+UXdk28gOLly/1StxBGCvhED/Br3q90i+mG1WQ1OCoB0Ce+D2+f/zZzBs3hq31fMfGjiXy842O/tN0rpYiZdAFd338Pa0oKGTfeROaDD+J2OFq8biGMEhSJfnv+dmmfb2WsJqvnYSYT/ku3mG7c+929XLv8WvaX1L1Fo2XYu3alyztvE3fVVRS8s5DdF12E49df/VK3EP4W2Ik+KZHc8lxyynPkaVKtVFp0Gm+MfoO7T7mbjVkbmfTxJN7a+pZfhlBQNhvt7ryD1JdfxplfwO4LLyLvrbcMfTauEC0hIBO9u6wMd2kplsREuRDbBpiUiUt6X8JHEz9iULtBPLrmUS759BI2HNzgl/ojTh9O2scfETb0FA7+5WEybrgRZ36+X+oWwh8CMtHX7kP/a57n57ic0bd+7SPa88KIF3js9MfIc+QxY+kM7lx5J5klmS1etyU+ntQXX6Tdn++i9Jtv2DXxAkpX/9Di9QrhDwGf6LfnbycpLImYkBhjgxKNopRibNpY0i9I59oB1/Llvi+Z8NEEXtj4AuXO8pat22QibsYMury70Psgkz+R9dTT6KqqFq1XiJYWFIlezubbnjBrGDcMvIH0C9I5M/VMnv/peSZ8NIGlu5a2eBt6SN++dH3/PaKnTCZ33jx2X345lfv2HfuDQrRSAZ3odXwMuwp2Sft8G9YhogNPnPkEr5/3OjH2GO74+g6uXHol6w+ub9F6TWFhdHj4YTo+9SSVO3ex64JJFC76pEXrFKKlBGyiV1Yru8nFqZ1yRh8ABicPZuH5C7l/2P3sLd7LlUuvZPZns/k5u2WfERs1ZgxpH32IvWdP/rjjDvZddz2O7dINU7QtgZnos7IxJybwa77nD7JnnPShDwRmk5kLe17I4smLuX3w7WzL28Zliy/jhi9uYEvulhart/oxhYm33UrZunXsuuAC9t9xJ5V797ZYnUL4UmAm+uysmvb5EHMInSM7Gx2S8KFQSygzTpjB0ilLmTNoDhuzNjL1k6nM+XJOzcHd15TFQsKsWXRf/jnxM2dS/Pnn/D72fDIfeICqg34ZkFWI4xagiT67pmtl95juMgZ9gAqzhjGz30yWTVnG9QOvZ82BNUxJn8ItK25hY9bGFqnTHB1N0m230u2zZcRefBEF773P76NGcfDxx6XvvWi1AjPRZ2VjSfD2uJELsQEvwhbBdQOuY+mUpczuP5s1B9ZwxZIruHzx5Xy+5/MWucvWmpRE8n330W3JYqJGn0fea6/z+8hR5LzwAu7SUp/XJ0RzBFyid1dW4iospCImjIKKAhnjJohE26O56X9u4vMLP+fPQ/5Mbnkut351K+M+HMdbW9+irKrM53XaUlPp8NhjdPXeWZv9zLPsGDmKvAULcFdU+Lw+IY5HwCV6l7drZVaY5yYXOaMPPmHWMC7tcymfTPqEJ896kvjQeB5d8yjnvncuT69/ukUGTgvp2ZPU556jy7sLsffsycFH/sbvo8dQ8P77aKfT5/UJ0RQBl+ir+9Dvs5UAyBl9EDObzIzsPJJ/j/03b455k6Hth/L65tcZ8/4Yrvn8GpbtXkaVy7d3vYYOGEDnN16n02uvYklIIPOee9k5YSJFS5eh3W6f1iVEY1mMDsDXqryJ/jdzDh3DOhJpizQ4ItEaDEwayMCkgWSWZPLRjo/4YMcH3L7ydmLtsUzoNoHJPSeTFp3ms/rCTz2VLsOGUbx8OdnPPMP+uXOxdetGzORJRE+YUDOEthD+oFrbkKyDBw/W69atO+7P5739Ngcf+gsP392V+JTuPHvOsz6MTgQKl9vFqsxVvP/r+3y17yuc2smgpEFc0P0CRnQeQZQtymd1aZeLok8+If+dhZRv3AhmMxHDhxM9eTKRZ5+Fstl8VpcIXkqp9VrrwfWuC7REn/XMM+TOe4lL7jQza+A13DDwBh9GJwJRTnkO6b+n88FvH7CnaA9Wk5XTOp7G6C6jOTv1bMKsYT6rq2LnLgo//JDCjz/GmZWFOSaGqPHjiZl0ASF9+/qsHhF8girR/3HvvRSs+IKps0t46qynOLfzuT6MTgQyrTWbcjaxdPdSlu1eRlZZFnaznTNSzmB0l9GcnnI6oZZQ39TlclH6/fcUfPABJcu/QFdVYe/dm5jJk4gaNw5LXJxP6hHBI6gS/d5rriFn32/MmJrN4kmLSY1K9WF0Ili4tZsfs35k6a6lfLbnM/IceYRZwjit42mckXIGwzsOJyE0wSd1uQoKKFy8mMIPPsTxyy9gtRJ51plET5pMxOnDUVZ51rE4tqBK9DsnT2aPtYi7JhSz6tJVmFTAdSwSfuZ0O1l3cB3Ldi/j631fk1WehUJxYsKJnJ5yOmemnEmfuD4opZpdl+PXXyn88CMK09Nx5eZiTkggevx4oieMx967t0/qEIEpqBL9r6efzoY0xadTO/Pm2Dd9GJkQnuadbXnbWJmxkm8yvmFTziY0msTQRE5POZ0hyUMYkjyExLDm9arRVVWUfPMthR9+QPGKr8DpxJyQQPiwYZ7p1GFYk5N986VEQAiaRK9dLrb160/6aTaqrp7CvUPv9XF0QhwutzyXb/d/y8qMlaz+YzXFVcUAdI3uypDkIZycfDInJ59MXMjxt7k78/IoWfEVpatWUbpqFa7cXABsaWmEn3oq4acOI2zIEMwRET75TqJtCppEX5WVxY4zzuTl80wMvf5+Lu51sY+jE6JhLreLbXnbWHNgDWsOrGHDwQ2UOT3DLnSP6c5J7U5iQOIA+if2p1Nkp+NqhtFuNxW//Ubpd99TumoVZWvXoh0OMJsJ7d/fk/hPO5XQfv2kbT/IBE2iL9+8md1TLuTxKSbmzH2bAYkDfBydEI1X5a5iS+4W1h5Yy5rMNfyc8zOlVZ4Bz6Lt0fRL6Ef/xP4MSBjAiYknHlfffXdlJeUbfvSc7X//vedirtaYwsMJGzLE29QzFFtaGsoso7gGsqMl+oC6M7Z6+IOCcEWPmB7GBuNyAvUdRBWYzCAX1QKe1WRlQOIABiQOYGa/mbjcLnYW7uTn7J/ZlLOJn7J/4rv936G9/09SIlLoHdebXnG96BXbi95xvUkOTz7qmb/JZiN86CmEDz0FbpmLq6CA0h/WUPq954y/ZMUKAFRoKCE9e2Lv24eQ3n0I6dMbe48emEJ9011UtG4BmejD2nf06U0ulBdAwR4o2OuZSrOhNAfK8qAsBxyFUFUGVeVQ5fDM66MMjatMYAkBix0soZ5Xa6hnmTUM7JEQEgX2qHpeow9NYfEQEgPmgPpnDFhmk5kesT3oEduDKT2nAFBSWcIvub+wKXsT2/K2sT1/O1/s/aIm+UfZomoSf9fornSN7kpadBpxIXH1HgDMMTFEnTeKqPNGAVCZkUHZmrU4tm2lYstWihZ9QsE7Cz2FTSZsaV0PJf7evQnp00f68AeggMoQ1Yk+ObXP8W2gqhwyf4IDm+DgL3DgF8j9zZPIazNZPUk2PAHC4iChJ9jCDyVqa4hnvr4zMa3BVQlOh+eg4Kw1VR8kijIgqwgqijx162MMhhUS44kjLB5Cva9hcZ7psPfV6+PALO23rUGELYKh7YcytP3QmmVlVWX8mv8r2/O2sy1/G9vztvP+b+9T7iyvKRNtjyYtOo206DS6Rnelc1RnUiNT6RjRkRBLSE05W0oKtpQUYBLg6TVUtX8/jq1bqdi6DcfWrZStX0/RJ4cefG5p146Q3r2x9+lNSO8+2Lp2wdoxBXNEeMvvENEiAirROw5mUhQKPRIbmeirymH3d7DnO9jzPexfD27vaIYh0dCuH5x4IcR2gZhOENsZYjpDaKz/ml60hspSb9IvOvTqKPD+osiFcu9rWR4UZ0LWFs/7o42/bo+qdSCocxBoaLlVfub7Q5g1rGYQtmpu7eZA6QF2Fu5kV+EudhbuZGfBTr7c+yX5FYc/2SopNImUyBRSIlNIjUwlJTKF9uHtaR/ensSwxEPJf+TIms848/Op2L4dx9ZtOLZuoWLrNkq+/RZch36ZmmNjsaamYktJwZqaijWlI7bUVKwpqViT26EsAZVOAkpAXYz9ZdYV7N+2DvvbL3BW6ln1FyrPh18/g22fwI4voKoUTBbo8D/Q+VRIHQrJ/SA6pe23o1eVe5J/7QNB9Wt5Xp0DhXeqLG54e9Ywb9KPrXNgOMpBwhbR9vdjK5fvyGdf8b6aKaM4w/NakkFW2eHPs1UoEsMSSQ5PJjksmfbh7UkKSyIxLJGE0AQSQxNJDEsk3BqOu6KCih07qNq3j8p9+6jal0FVxj4qM/ZT9ccfUHucfYsFa/v22FJTPIk/JQVbagqWdslYEuKxJCRgCvNhc6o4QtBcjHUczCQ/XHFmbJ2Hjbhd8PsK+HEBbFvsOWuPSIYBU6HXWOg0DOwB2AfZGgrRHT1TYzkrPQfDww4AdQ4G1fMF+zzz5QXUf+EZTzNXaIznF1Ldawwh1e/rW++dlwPFMcWGxBIbEkv/xP5HrHM4HfxR8gcHSg+QWZrJgbIDZJZ4Xrfnb2dlxkoqXEc+CSvUEkpiqCf5x4fGE3tiLLGDE4kN6UlcSBxx1mhiilxEZJVizyrAvf+A54CQkYFj+XJceXlHbFOFhWFJSMAS70n85oR4LPEJnmUJ1cs86+UisW8FVKJXOfmUpFpJDvfeMegogg3zYfWLnnbvsHgYMsvTHNPhf8AkwyMcwWKDyHaeqbHcLk+yr/fAkOu5zuAoPHTNoWj/oWVOx9G3rUzeg0Ak2CI910LsEZ5XW4R38s7XLA8/VNYW7v2sd94a5un1FCRCLCGkxaSRFlP/WPtaa4oqi8guyya7PJuc8hyyy7PJLjs0/3vB7+Q78imoKKi5SFxXaFIoUSlRRNujibL1JEGH077ITHyJIqrETWSJi7CiCuxFDmwF5Zh/24paUwSFRfVuzxQejjkuDnNkJKaoKMxRUZiiIjFHRmGOisQUGYU5OgpTZCTmqKhD5SIjUaGhMlREHY1K9Eqp0cAzgBl4RWv9aJ31dmABcBKQC0zVWu/2rvszcDXgAm7WWi/zWfS1aK2xF5Zj/p/2qPJ8+P5ZWPsaVBRCl9PhvL96zt4tMva3z5nMEB7vmZrKWeG95lB9MCg8NO8oPLSuosh7raLYs6woEypLvFOp5wJ3Y5ltnt5O1RfNraHeXk/Vy7yv9S4LO/QZi92zrdqTxXbkssOW21vVCYZSimh7NNH2aLrHdj9qWZfbRUFFAfmOfPIr8slz5JHvyKeosoiiiiIKKwsprCikqLKIHZX72aCKKA4pptxSDjH1b9PsMhNVBnGliiSHnUSHlfgyM7FliqiycsIqygjJziRkrwt7WRXW8iosFUd/NKO2mFER4ajQEFRoKKbQMMxh3ik0HFNYmGcKDcUUFlpTxhQWiinU+95uR9ntKJsdZbPWem/zvFqtbepgcsxEr5QyA/8CRgIZwFqlVLrWekutYlcD+Vrr7kqpacBjwFSlVF9gGnAC0AFYrpTqqfXR+h4en6r8PMwuTXgo8OxAcBTh7nsBZYOvpzi+H6UVLsoOlFFWWUyF002V002Vy02ly02l002VS1PpdHlea5a5cWtwa43brWvmtfbMu6rn3d4y2nPAcdW67qE8+xBVa0H1O6Wq11evUp5ltQpXlzlUXmFS3m2qQ59pcFnNNj21mtShbZpMni0ftsw7X3t57fW16zB5Y63etlKHlplU9XdUh8XhqbLWd8CT9xRRKBWFUp1QVjDZFCr68Drrj8Oz3OSuxOIsx+QsxVxVgtlZjtlZ6pmqSjFVlWGuKkU5yzG5HJ7J6UA5y1FOz3vldKBKCzA5M8HpfV8zHerx0mzKfGTyN1sPHRRMZs91o5rJ7Flf+33d9TXz1mOsr35v9sShTJ7JVGteqcPXedeblYl4pYivWRcGIREQdpTPmsxUaTelrgqKneUUO8spcZVT7CyjxFlOmctBaZXntcxZTqmznP3Ocn51llPuLKfcVUG5y4HD6aDc5aDcqXBVacIqINzhmcIqNBEOCHNAeAWEO9yEOYoJqSrGXgX2crAXaUKqwF4JoVXKs7xSY2nG0x1dVhNuqxm31VLzqq1mDo4aSNboQZiVGYvJgkmZMCuzZ/I+3jLc6t8eTI05ox8C7NBa7wRQSi0EJgK1E/1E4AHv/HvAc8pzuJsILNRaVwC7lFI7vNtb5ZvwD1m/YiExQFzZDr52dOOxqkvYvKEjbMgGvjyubZpNCrPJk5TMSh1KfCbPvMmbfEzqUOIymbwJCU+rtdbU/NzV2jNV0/rQD+HqctXrqz9L7c/W+ozbrQ9tv3qZ1jV1VG/LrXWtbQWjMO/UXBo7VYRQ6ZlUJTac2KjChhMrTmzK+1pnmQ0nNuXCTpXnvXJiczmxVjqx13zm0OdNuLFQhZkKLDgx48aCCzOuI+bNuLB4p8OXubHSOh5KbsVzQh/jo+25AIdSOJSiwqRwWBUVNkVFtHeZd12lUlR535coRZ53WYVSVCmoUooqtwntBO1U4FQo76RdCuUC5X011bx6JrMLzE43ZpfG7HJicYLFBVYXrD64gm/Xr2ww/l4L55LkPdWtfUpYiZVXT/qI/xvn+wfQNCbRdwT21XqfAZzSUBmttVMpVQjEe5evrvPZI64MKqVmA7MBOnXq1NjYD99GbAe29zCzr/Ol5PS+jNNtFkbbzITZLYTbzITazITbLITZzNitJqxmEzaL97XWvNWsPPMmU80Zb6DQ+sjk7/YeAWoOErXmqVX20PrDDybVv2LqW+b2Hqg89dQ5+DR4QKoVx9EOZt55qut014m9Jo76tl0rpvqWeT9LPQfVo+6vWmVrr6+93x1Aea3vU72fam/7WP+Gh72vt8yheaXdKO3EhAuT9kzK7TmYKO1G4cak3YBG4UZpjQkXaO1dr1HahQmNwlVrmaesZ9mhzyrt8qyv2a7bW1d1mTrznr3g/dWqPe917WXU3EdyqKxnx6rq+Vrva9YDJu0mBAg9dDrl+V41e06jlAYrKGutbQFKH9rWYbFVz+vDlynP/yDcWpOCZnKeG7fy/l9C4wbcyo0bTVa4lVzUoe1V/99RZjrGtMxF6FZxMVZr/RLwEni6Vx7PNoaeMwnOmeTTuAJNTRMIgXUAEyJQnN5C223MVaH9QO3HNKV4l9VbRillAaLxXJRtzGeFEEK0oMYk+rVAD6VUV6WUDc/F1fQ6ZdKBGd75C4Evtec3ZjowTSllV0p1BXoAa3wTuhBCiMY4ZtONt839RmAZnu6Vr2mtNyulHgLWaa3TgVeBN70XW/PwHAzwlvsPngu3TuCGluhxI4QQomEBNQSCEEIEq6MNgdB67twQQgjRIiTRCyFEgJNEL4QQAU4SvRBCBLhWdzFWKZUN7GnGJhKAHB+F40sSV9NIXE0jcTVNIMbVWWudWN+KVpfom0spta6hK89GkriaRuJqGomraYItLmm6EUKIACeJXgghAlwgJvqXjA6gARJX00hcTSNxNU1QxRVwbfRCCCEOF4hn9EIIIWqRRC+EEAGuTSZ6pdRopdR2pdQOpdRd9ay3K6Xe9a7/QSnVpZXEdaVSKlsptdE7zfRTXK8ppbKUUr80sF4ppZ71xv2zUmpQK4nrLKVUYa39dZ+f4kpVSq1QSm1RSm1WSs2pp4zf91kj4/L7PlNKhSil1iilfvLG9WA9Zfz+N9nIuAz5m/TWbVZK/aiU+qSedb7dX9r7gOu2MuEZKvl3IA2wAT8BfeuUuR540Ts/DXi3lcR1JfCcAfvsDGAQ8EsD68cCS/A8a3so8EMriess4BMD9ld7YJB3PhL4tZ5/S7/vs0bG5fd95t0HEd55K/ADMLROGSP+JhsTlyF/k966bwXeru/fy9f7qy2e0dc8rFxrXQlUP6y8tonAfO/8e8AI78PKjY7LEFrrr/E8J6AhE4EF2mM1EKOUat8K4jKE1jpTa73BO18MbOXIZx37fZ81Mi6/8+6DEu9bq3eq28vD73+TjYzLEEqpFOB84JUGivh0f7XFRF/fw8rr/mc/7GHlQPXDyo2OC2CK96f+e0qp1HrWG6GxsRthmPen9xKl1An+rtz7k/l/8JwN1mboPjtKXGDAPvM2Q2wEsoDPtdYN7i8//k02Ji4w5m/yaeBOwN3Aep/ur7aY6NuyRUAXrXV/4HMOHbFF/TbgGb9jAPBP4CN/Vq6UigDeB+ZqrYv8WffRHCMuQ/aZ1tqltR6I57nQQ5RSJ/qj3mNpRFx+/5tUSo0DsrTW61u6rmptMdE352Hlhsaltc7VWld4374CnNTCMTVWq3yIu9a6qPqnt9Z6MWBVSiX4o26llBVPMn1La/1BPUUM2WfHisvIfeatswBYAYyus8qIv8ljxmXQ3+RpwASl1G48TbznKKX+XaeMT/dXW0z0zXlYuaFx1WnDnYCnjbU1SAeme3uSDAUKtdaZRgellEqubpdUSg3B8/+1xZODt85Xga1a6ycbKOb3fdaYuIzYZ0qpRKVUjHc+FBgJbKtTzO9/k42Jy4i/Sa31n7XWKVrrLnjyxJda68vrFPPp/jrmw8FbG92Mh5W3grhuVkpNwPOg9Dw8V/xbnFLqHTy9MRKUUhnA/XguTKG1fhFYjKcXyQ6gDPhTK4nrQuA6pZQTKAem+eGADZ4zriuATd72XYC7gU61YjNinzUmLiP2WXtgvlLKjOfA8h+t9SdG/002Mi5D/ibr05L7S4ZAEEKIANcWm26EEEI0gSR6IYQIcJLohRAiwEmiF0KIACeJXgghApwkeiGECHCS6IUQIsD9P2ds/KO1UGXMAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "dts = np.array( [0.001,0.01,0.1,0.2])\n", "\n", "for i in range(dts.size):\n", " t = np.arange(0,4,dts[i])\n", " tru_y = y0 * np.exp( -params[0]*t)\n", " (t,y) = forwardEuler1D( degradation_model, t, y0=1, params=params)\n", " plt.plot( t, tru_y-y[0,:].T )\n", "plt.legend( ['dt=0.001','dt=0.01','dt=0.1','dt=0.2'] );" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### What happens if the step size is way too big?\n", "\n", "Use the forwardEuler method with $\\Delta t$ = 1. Plot the solution along with the true solution." ] }, { "cell_type": "code", "execution_count": 132, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAD4CAYAAADhNOGaAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy86wFpkAAAACXBIWXMAAAsTAAALEwEAmpwYAABERUlEQVR4nO29e5Qdd3Xn+9n9VrfULfXj6N2WbEnu0zZ+QFuEwIXBNmASbEOAYLFmMFwyXrNuCEmYOwGGdUNChlxyuTeQzGXNXAcIJjA2xsDYBgdjjAkkwaA2fqpPy3pYtt511C21pJb6/bt//KpOV7fO6T6Petfvs1avPqdOnapfV++q/Xvsvb+ilMJgMBgM6aUu7AYYDAaDIVyMIzAYDIaUYxyBwWAwpBzjCAwGgyHlGEdgMBgMKach7AZUQ3d3t9qyZUvYzTAYDIZY8dRTT51SSvUs3h5LR7BlyxYGBwfDbobBYDDEChF5udh2MzVkMBgMKcc4AoPBYEg5xhEYDAZDyjGOwGAwGFKOcQQGg8GQcjxxBCLyVRGxROSFEp+LiPytiOwXkedE5NWuz+4UkX32z51etMdgMBgM5ePViOBrwC1LfP52YLv9cxfw3wBEpBP4NPBaYCfwaRFZ41GbDAaDwVAGnjgCpdTPgNEldrkd+LrSPAmsFpH1wNuAx5RSo0qp08BjLO1QauLp79/Nr779eb8OHz+efwAuLPVvSxHjI/p6GDQv/wKOPxt2K6KBUvDrf4Cp8bBb4htBrRFsBA673h+xt5XafgkicpeIDIrIYD6fr64VuYfI7PkKE9Oz1X0/SZx+Gb7zYRj8StgtiQaDX9XX48zh5fdNAw99BH74ybBbEQ2OPqWvR4I7CrFZLFZK3a2UGlBKDfT0XJIhXRZdW69lszrBY88VTa5LF1Zu4e+0Y+2xf5vrwfRFGDkA1pDuDaedk8m3jaAcwVFgs+v9Jntbqe2+sOnK11Avin958hd+nSI+WEP27+Qad0UUHONQuO2IAvm9gIKLp+H8ybBbEz4psI2gHMFDwAfs6KHfAMaUUseBR4G3isgae5H4rfY2X6hb2w/AxNEXeOlUcuf7ysIx7lMvwsxUuG0Jm5lJGNmvXxvHuPAaJPjhVzYp6DR5FT56L/AL4EoROSIiHxaR/yAi/8He5RHgILAf+DvgfwNQSo0CfwHstn8+Y2/zh64rUHWN9NUf4b7dr/h2mlhg5UDqYW4GRg+E3ZpwGdmvr4PUmwcf6Gsg9fbr5D78ysa5V8YtHVSQQLyKGtqllFqvlGpUSm1SSn1FKfXflVL/3f5cKaV+Xyl1hVLqVUqpQdd3v6qU2mb//L0X7SlJfSPSvYPfXHWS7zx1hKmZOV9PF1lmZ/RI4PI36fdpv9mdv//yN+nrMpfyYIL8MGT6oa3H2Mb4iHYAzr2ST+b1iM1isWdk+tguRzl1forHcymd/zz9EsxOQvY2uxecTOMuGysHdQ36esxMwOlDYbcoXKwcZLL6J+224Tz4r363/p3Q65FCR5BlxfgRtrUr/sevUjo95Ex/bLgOOi830yFWDjqvgPXX2O9TfD0mzsLYYcj0QU9Wjw7mUjpyBtdo8c3Q0pFY20ihI9ALxh/um+Kf95/i8OiFkBsUAlYOEOi+0vT6QN/cmSz09NnvU3w98sP6d6ZfX5Op89oxpBVrCJo7oH2DviYJtY0UOoIsALdkTiPA/YMpNHJrCDq3QlOrNu7Rgzp2PI1MXdBTQZl+aGqDNVsS2+srC+dvz2QLnaakPvzKwpkmE7E7TcnMrUifI1i9BRpWsOb8ft60o4f7Bw8zM5uyoa81PH+TZ7KA0oukaeSUHTNvdxB0r2841CaFijUMjW3Q0aunhyCxC6TLotS8IwBtGxNjcO5EuO3ygfQ5gro66LkSrBx37Ozl5NlJnthbZcmKOOLEzLuNG9Lb63P+brdjHNmX3twKa0g7gLo6PSfevim9tnHuBEycWdRpIpEjxvQ5AijM9d3YlyGzqpn70rRofGofqNn5+fDOy6G+KZHGXRbWENQ366ky0AukczPzCWZpw8rpa+CQ6Uu3bcD8yMi5Lgl0jCl1BFk4f4LGyTO8d2ATT+y1OD6WkjnyxT3g+gbo3pFI4y4LKwc9O6DOTqBKcK9vWcZP6Zj5jNsRZCH/os49SRuL75W2LmjLJPJeSakjmJ8Oed9AL3MK7t99JNw2BYU1pGPmu7bNb0tz5JCVm7cHgO7t6c2tKDz43I6gX+ecnH4pnDaFiZXTSXVt3fPbnAXjhJFSRzDf6+vtauUN27q5f/Aws3PJiwa4hPwwdG2Hhqb5bZmsDhGcOBteu8JgYgzOHl344Gto1k4yn8IFY3foqEMmudMhy5LPLbQN0NcmvzdxuRXpdATtG6C5vWD4u3b2cvTMRX6+LwWLxk7MvBvnxs/vDb49YWIVefBBYnt9y2INQctqWLVuflv3lYCkzxHMzS2MrnPIZGF6HMaSta6YTkdQiAnWxv2W/rV0tTVxb9IXjafG7Zj5RY6gkEiVsoef8/c6f79DJgujL+kcgzThjpl3aGpNZ27F2Cv6gV/MNiBxjjGdjgAWJIc0NdTx7tds4vGchXVuIuyW+Udh6L/IEay+DBpbE2fcy2LloGkldGxeuL2QW5GiEZJSxUeLkOiM2pIsXih2SGinKcWOoH+B8MYdN2xmZk7xwFMJXjQuNRVSV6cNPG2JQ9aQ/rvrFt0GhWCCFK0TnDuu10wW2wbYI6QDOgclLRQcwaIRQUu77jgkzDGm1xEs8uyX96zktVs7+dbuw8wlddHYGoKGFj3UX0wae3354eI94DVbdW5Bwnp9S1IsYsghk8LcCiunk+laOi79LJNNXCchvY6gSK9v185eXh65wC8OJlN8QsfMXzkfM+8mk9Wjo4QKb1zC+TyM54v3gOsbdG5Bmhyj87f2lJgacu+TBqwiEUMOmayeNkxQboVXCmW3iMheEdkvIp8o8vkXROQZ++dFETnj+mzW9dlDXrSnLFb2QGv3gl7fLVevo2NFY3IXjRdnjbpJW12ZfImhv0NPynIrrJxOlmrruvSzrm069yQtI6TZGf2gX8o2Zqd0scaEULMjEJF64EvA24F+YJeILOhmKaX+WCl1nVLqOuC/At91fXzR+UwpdVut7amIRYlULY31/M6rN/KjPScZOZ+w+dCLZ+DcsSV6OSnr9ZVaDHTIZOHsET1vngZKLRSDzjnp2pYe2xg9qB/0S9kGJMoxejEi2AnsV0odVEpNAfcBty+x/y7gXg/OWzuZ/kuEN3bt7GVqdo7v/vpoiA3zgWLJQm5WrbeFN1Jys1tDsGINrFxb/PM05VbMzc3LU5YiTdnn+SXWS0BPryYst8ILR7ARcBf1P2JvuwQRuQzYCvzEtblFRAZF5EkReWepk4jIXfZ+g/m8R4lfmb5LhDd2rF3Fq3tXc+/uV1BJqjvurjNfDJF0LRg7yULumHk3Cez1lWTsFZi+UNo2QF+r04d0LkrScQs3FaNxhS7WmKBp1KAXi+8AHlBKudXBL1NKDQDvB74oIlcU+6JS6m6l1IBSaqCnp8eb1hR6fQsjAHbt7OVgfpzdh057c54oYOWgaRV0bCq9T4KFNxawuM58MTo26xyDNDjG5abJYD63Ig0jJLdwUykSNkLywhEcBdwZOZvsbcW4g0XTQkqpo/bvg8BPges9aFN5lEgO+e1r1rOquSFZi8ZWTo+ASvWAQS+CTZxJpPDGAs4eg8mxS7NG3RR0K1IwIihkWJfoAUOiSzBfwlJBFQ6ZLIwcgOlkJKB64Qh2A9tFZKuINKEf9pdE/4hIH7AG+IVr2xoRabZfdwOvB4K781ashvaNlxh3a1MDt1+/gUeeP87YhenAmuMrSy0GOqRlOqScHjAkrtdXEiunR0At7aX36UxJbsXMpH7Al3OvqFktYpQAanYESqkZ4CPAo0AOuF8ptUdEPiMi7iigO4D71MKJ9ywwKCLPAk8An1NKBWtpJQqM7drZy+TMHN97OgGZxufzcGGkvAcfJL/y5nLrJQ6Zfp1rMH7K/zaFyXLTZKBzT3quTL5tOMJN5dgGJKaj0ODFQZRSjwCPLNr2p4ve/1mR7/0r8Cov2lA1PX3w0s9hbnZBotVVGzq4ZlMH9+0+zJ2/uQVZakol6pT74GvrtoU3Et7ryw/DynXQ2rn0fu4CY1v/F//bFQazM1qvettNy++b6YdDP/e/TWFS7mix8wqoa0yMI0hvZrGDI7wxeqnwxh039DJ84hxPHz4TfLu8pFzjhnRMh5QzTQaJ6/UVZbmYeTeZrNZvuHjG92aFRjHhpmIkLLfCOIIl5sVvu24DrU318dc0toZgRadWW1oOp45KwoQ3ChTqzJfhCFau1fX5kzxCKlWKuxhpmDq0cvoB7xZuKkWCdCuMI1giOWRlcwO3XbuBh589zrmJGC8aO3KM5UxvJVR4o8CZQzBzsTxHkIbcCidmfqmIIYc0BBOUO1oEbRtnXobJ8/62KQCMI2hqW1J4446dvVycnuXBZ44F2y6vUKp0lc1iJL0EcyXTZGCLt+eSm1thDenkqMYVy+9byK1IqG1MjesHeyW2AYnIrTCOAJacF792Uwd961Zx3+6Y9pDPHoXJs+U7goQKbxQoVNksowcM+rpNjOl6/Umkkk5CQdkvobZRSripFAVHEP8Ro3EEsKTwhojw/tf28sLRszx/JIYFyJaqM1+MlnZdhz2p0yFWDjp6oXlVefsneTpkeqK8mHk3PX3Jtg0of0SwZovW90jA9TCOAPQ/fgnhjduv20hLYx33xnFUUMlioEOSI4fKiZl3k+SM2pEyY+bdZPrhwimdm5I0rFxp4aZiOLkVCegkGEcAywpSd6xo5LdetZ6HnjnG+GTMxCisnK4sulzMvJsECm8AMDutY+YrefC1denooSQ6gkp7wJDsEZI1BN07igs3lSIhwQTGEQB0bV9WeOP9O3s5PznDD56L2VxxpT1gsHMrpuD0pbkVsWbkAMxNV/bgg+SOkKwhnRTVWbTOY3FKFGpMBE5F2krIZPX60cV4F6g0jgB0zHDnFUve7K+5bA3bMivjNT00N6sjGqoxbkher2+5OvOlKKJbkQisYejeXl7MvMPKjM5JSZptXDy9tHBTKRISZWccgcMyvT4R4Y4bNvP0K2cYPnE2wIbVwOlDOma+kvUB0MPjhAlvAPrvkTr776uAnj5dr//My/60KyysocptoxA5lDTbWEa4qRQJibIzjsChDOGNd796E031ddz3q8Ml94kU1cwBg67D3rk19sZ9CYWY+ZbKvpfEUhOT5yuLmXfjOIIk5VYU6nFV6Bg7Nmmdj5jbhnEEDmUIb6xpa+KWq9fx3V8fYWJ6tuR+kaHSmHk3CVkEW0A16yUwf/2S5BgdO6/memSyOjflbILkXK2cTpbr2Lz8vm4SMkIyjsChzF7fHTs3c3Zihkeej8GicT4Hqy+D5pWVf9cR3iiSWxFLpid0gbVqesAt7Tr3IEkLpOVWpC1GQubFF+Ak1lVTZTgByn7GETiUKbzxusu72NLVGo/pIafGUDU4whunkiG8wakXQc1V9+CDRPT6FpAfhoYV5cfMu0nIvHgBpeDknhpsox8ujmrtiphiHIFDXT307Fi21yci3LGzl18dGmW/FeFiUzNTdsx8hXOeDklLpCpMk1V7s/fp6zkb4+KDbqwhbe+VxMw7tHZqPYek2MZ4Xj/Iq+40xd8xeuIIROQWEdkrIvtF5BNFPv+giORF5Bn75/dcn90pIvvsnzu9aE/VlDkv/u5Xb6KhTvhWlENJRw/obOlqjbtr27K5FbHCiZnvqiBm3o2TWzF60Nt2hUUto0VIVs2harLv3SQgmKBmRyAi9cCXgLcD/cAuESlmYd9SSl1n/3zZ/m4n8GngtcBO4NMisqbWNlVNmcIbPauaeUv/Wh546giTMxFdNK5lDhhs4Y3tsTbuBVg5HTZa31jd95OUW3FhVCdBVWsbYOdW7E1GbkW10XUObT3Q2hVr2/BiRLAT2K+UOqiUmgLuA24v87tvAx5TSo0qpU4DjwG3eNCm6qgga3LXzl5OX5jmR3tO+tyoKrGGQer1w7xanBLMSaDaiCGH7h06ByEJC6T5KmPm3WSyOkflzCFPmhQqVk4nya3MVPf9gm5FfG3DC0ewEXCvnB6xty3m3SLynIg8ICJOjFa530VE7hKRQREZzOd9WpSpYBHsDdu62bRmRXTLU1tDehqk0ph5N2XkVsSCyXNaaKcWR9C4QucgxLjXV6DSirTFSMB0SIFKhJtKEfPciqAWix8GtiilrkH3+u+p9ABKqbuVUgNKqYGenjIkF6uhAuGNujrhfQOb+Zf9I7w8EsEHpZWrfs7TwVkEi3vYZC0x826SUoLZykFzO7QX7XOVR1JyK5SqfbQI2jamzsHYEW/aFTBeOIKjgDsLY5O9rYBSakQp5QSkfxl4TbnfDZS6OvtmL8+43zuwmfo64b7dEQslnb5Yfcy8m6T0+mpdL3HI9OtF+OmJ2tsUJk4noZYecPNKWN0bf9sYO6If4NVG1znE/F7xwhHsBraLyFYRaQLuAB5y7yAi611vbwOcq/Uo8FYRWWMvEr/V3hYeFcSLr+to4c1XZvj24BGmZyO0aJbfC6jaH3xJEd6wcjpmfvWW2o6TyepchFMvetKsUFCqMl3epUhC9nmtC8UOMQ8hrdkRKKVmgI+gH+A54H6l1B4R+YyI3Gbv9lER2SMizwIfBT5of3cU+Au0M9kNfMbeFh4VCm/s2rmZU+cneTwXoUVjLxYDwSW8EfebfUjfqHU1mnsSSjCft2qLmXeTyeqEwzjnVjjBELVOo65YA6s2xNY2Grw4iFLqEeCRRdv+1PX6k8AnS3z3q8BXvWiHJ7g9+8o3Lbv7m3b0sK69hXt/dZhbrl6/7P6BYA1BfZNe3KyVTD8c/KfajxMm1jBsu6n243RdoXMRYtrrA6ovxV2MTL/Wdxg5UPvUSlhUI9xUihjnVpjM4sVU2OtrqK/jd2/YzM/25Tly+oKPDauAQsy8B36+p0/XaY+r8MaFUTh/ovYeH+gchO6Y51Z4ETHkkIRSE15Nk4Edbr1X64DEDOMIFrNyrR7mVWDcvzuwCYD7ByMSMeBFFIRD3AuMeTUH7BDjXh+g297apZOgaqWQWxFTx+gIN1VbdmQxmSzMTOiQ65hhHMFiCskh5Rv3pjWtvHF7D/fvPsxM2IvGE2dh7LC3vRyI78PPq4ghh0wWzryi6/nHES9i5h0aW2xlv5jaxulD+sFt7hXjCIpSRXLIrp29nDg7wT+9GHIFwkLMvEc9YEd4I6aLYDpmvgPaN3hzvMLUYWndisiilK3L69GDD+zpkBjbBnh3rxSmyuJ3PYwjKEZPX8XCGzdlM3SvbObesMtTe90DjrvwRi115osR417ffMy8l46gX+esTF/07phBUYtwUzGa2rT+RwxtwziCYlQxL95YX8d7Bzbxk+GTnBgLMeHIykFjqxZS8YpMn67XHrf0+ULMvIcRLau36JyEODrGWktxFyPTF9/cCmuoeuGmUsQ0t8I4gmJU2eu744bNzCn49mCIowJHkLzWmHk3cRXeOH9SRzt5NfQHO/v8ylj2+qrW5V2KOGfUehlU4ZDJwsg+rQcSI4wjKEaVwhuXdbXx+m1d3Lf7MHNzIfWea60zX4y4Tod4PU3mENNen46Z36Cj4ryi83KdsxI325iZ0g9sP2xjbgZG9nt7XJ8xjqAUVYYJ3nFDL0fPXOTn+0/50KhlGB+Bccsf44b4LYJ5vRjokMnq3IQL4SbBV4yXMfMO9Y06jDRutlGrcFMpnOsbs/LtxhGUopAcUlk46FuvWsua1kbu+1UI5am9zBp1E1fhDSun297W7e1x41hqYm7Wli712DYgnsEEfo0Wu7drHZCYXQ/jCEpRpfBGc0M973nNJh4bOkn+3OTyX/ASL7NG3YjoBcaYGbcnpbiLEccCY17HzLvp6dN6DxNnvT+2X1i52oWbitHQrEuRxOxeMY6gFDUsgr3vhl5m5hTf+XXAmcbWELR06NopXuPEi8clcmhuzg4d9XjoD7qOf3N7vG52v3rAEM/cCiun1zdqEW4qRQyzz40jKEUNwhvbMivZuaWT+371CirIB6c17F3W6GIy2YpzK0Jl7DBMnffnwVfIrYjR1FAhdNSPEVIM58X9iBhyyPTD6EswFZHaY2VgHEEpmlfVJLyx67WbOTRygV8cHPG4YSXwss58MeIWJujXQrGD0+uLywjJGtL6Ek1t3h979WU6dyUutuGVcFMpMllAxSq3wjiCpahhXvztV6+nvaWB+4LKND53AibO+GjcMZsXLyyc+1Qe2cmtOG/5c3yvsXyaJoOKlf1CxyvhplI4CXtxcYx45AhE5BYR2Ssi+0XkE0U+/5iIDNni9Y+LyGWuz2ZF5Bn756HF3w2VGoQ3Whrr+Z1Xb+KHL5xgdDyA5BLnJvRj6A+28Mb6+Bi3ldNz+S0d/hw/TiWYnZh5v2wD4hU55PdoMYa5FTU7AhGpB74EvB3oB3aJyOIr/DQwYIvXPwD8X67PLiqlrrN/biNKuIU3quCOnZuZmp3ju0EsGvsVMeQmTotgfk6TQbymykb2+xMz7yaT1Zncccit8FK4qRj1DdAdL2U/L0YEO4H9SqmDSqkp4D7gdvcOSqknlFLOysmTaJH66FNjRm3funau713NfbsP+79onM9BW8b7mHk3mX7Ivxh94Y3ZGd1OPx3Byh5o7Y7HAqmfEUMOmRhNh+SHvRNuKkWcRkh44wg2Au6J8CP2tlJ8GPhH1/sWERkUkSdF5J2lviQid9n7DebzAdW88UB4Y9cNvey3zjP4ss8KX35GQTg4uRVRF944/RLMTvrbA4b43OxOzHy3xzHzbgojpBiMGIO6V84eiU1uRaCLxSLyb4EB4POuzZcppQaA9wNfFJErin1XKXW3UmpAKTXQ0+OBulI5NLbo4WMNxv2Oa9ezsrmBe/3MNJ6b83cx0CEuvb4gpslgvuZQ1COH8sPQtU0nO/nFqvV6PSbqtuG1cFMpCiG18Qgx9sIRHAU2u95vsrctQERuBj4F3KaUKqTcKqWO2r8PAj8FrvegTd5Ro/BGa1MDt1+3gUeeP87YxcoXncti7BWYHvdfQLzbya2I+M1u5QCZb69fZPp0rsJYyBoUy+F1Ke5ixCX73LmXA+s0xWCEhDeOYDewXUS2ikgTcAewIPpHRK4H/j+0E7Bc29eISLP9uht4PRCtK+eB8Maunb1MTM/x4DM+JWP5HQXh0LwyHsIbhZj5Vn/PE4cF46kLOrnJb9uAeORW+B1d59DRC41t0bYNFzU7AqXUDPAR4FEgB9yvlNojIp8REScK6PPASuDbi8JEs8CgiDwLPAF8TikVradMJluz8MbVGzu4emM7/+OXPmUa+5k1uphMf/SHu36U4i5GIYQ0wjf7KZ9j5t1k+nUuy/mT/p+rWqxhnfy2+rLl962Fujo9Cot6p8nGk2VzpdQjwCOLtv2p6/XNJb73r8CrvGiDb7h7feuvrfowu3b28qnvvcCzR8a4bvNqb9rmYOWgYzO0tHt73GJksrD/MR2b3tDk//kqZWZSh0tmb/X/XCtW61yFKDuCoEaLsHA6ZNU6/89XDX4IN5Uik4UXf+T/eTzAZBYvR+flUNdYs2e/7doNrGis96c8dRBREA6O8MZodbkVvjOyH9RsgNcj4rkVVg7qm2HNVv/PFYdggqBGi6DXTMYtGA9Bm6RCjCNYDo+EN1a1NHLrtet56NljnJ+c8ahx6Jj5U3uDmRaC6JeaCCpiyKGnz9atiGhuhZXzP2beoa1b6z9E1TbGT/kj3FSKODhGG+MIysGjePFdO3u5MDXLQ88c86BRNqMHYXYquF5OV8SFN6whqGvwvs58KTL9Omdh9KVgzlcpQY4WIdq5FYVOQlCdphgEE9gYR1AOmawnwhvXbV5N37pV3Lfbw+khv1TJStHYEm3hDStnx8wHtH4R5RLME2M6qSlQR9BflbJfIAQVOuqwah20rI7uCMmFcQTl4JHwhohwxw2bee7IGC8cHfOgYczHzDv6CUEQ6V6fzzWGFtNzJSDRvB5WwA8+0Nc+qrkVfgo3FUMkHlF2GEdQHh7Oi7/r+k00N9Tx778+yOf+cZg9x8ZqCym1hqBzKzSuqLltZdOTrTm3whemxuH0y/NlgIOgqQ3WRDS3wu9S3MWIcglmZ6HYD+GmUjghpFHOrcA4gvJYvQUaVnji2TtaG/m7DwywY+0q/u7nB/ntv/1nbvp//om/fuxF9p08V/kBg4yCcHCEN6ImTeh3nflSOKUmooaV00lNHb3BnTOqwQR+CzeVItOvp+jOHQ/2vBUSQChBAvA4OeSNO3p4444eRsen+OELJ3j42WP815/s428f30ffulXceu0G3nHNei7rWkZNamZSl8juv33p/bzGvQi24bpgz70UQcbMu8lkYd+Popdb4ZSWCCJm3qGlA9o3Rc8xnjuuH8ihdJrQ/4v2DcGeuwKMIyiXTD/s/7Gnh+xsa+L9r+3l/a/txTo7wSPPH+f7zx3n84/u5fOP7uWaTR2845r1/PY1G9i4usjUz6l9wcbMFxpuC29EbYHUGtIx850BxMy7cXIrRvbD2oAfNEth5WDH24I/byYbQdsIMPvejXuqbFvRvNpIYBxBuWSy8Mw3YXwE2rq8P3x7Cx98/VY++PqtHD1zkR88d4zvP3ecv3xkmL98ZJjXXLaGW69Zz29ds57Mqhb9pbB6wFEV3rByevG2rj7Y87p7fVFxBOOnYDwfvG2Avh4v/UznuASRv1AOQeeXOLR1wcq1Nech+U1E/ksxoMcVJtj2Bl9PtXH1Cu564xXc9cYrOHRqnB88f5yHnz3Gnz08xJ9/f4jf2NrFrddu4F2jz7OirgE6i1bu9pdMH7zyZPDnXYr8MGzx939TlK5t0cutCKsHDPphOzupdSH81ECoBCsA4aZSxEDP2SwWl0tIWYJbutv4/Tdv44d/9EYe++M38gc3bufk2Qn+8/ee519+8XOO1G/igWctzk74VOK6FJmsDhGMivDGxTNw9mjwPT7Qdf67tkXTEYQ1IoBoPfzCWCh2cEJIo5hbYWMcQbm0b4DmcIU3tq9dxcfesoPH/+Ob+MFH38DAihMMz23if//2swz8xY/5918f5MFnjjLuZQmLUhRyKyIy5A06WWgxUas5ZA3pZKYwir91Ryy3Ym5O20eYtjF9Ac68HM75y8BMDZWLSGQSqUSEq7obYPIYN735Q3xv62/y/eeO8/3njvHY0ElaGuu4qW8tt167nn9zZYaWRh/mzN0jpM07vT9+pQShy7sUmX4YelDX//dbB6EcwoiZd2hq1Qv2EbhXAFu46UKw+RRu3FF2QQcylIlxBJWQ6YM9/1PHJIdxg7mxe8CS6ef63jVc37uGT/1WlsGXT/Pws8d45Pnj/OD547Q11fPWq9Zx67XrecO2HpoaPBoERk14w8pB00pdjjsMnNyKU3thQ8gie0rp6/Gq94TXhijlVoQ5TQbzWf/5HPT9VjhtWAbjCCoh0w9PfQ3OnYD2gNLUS1EkCqKuTti5tZOdWzv59K39PHlwlIefPcYP95zge08fpWNFI7dctY53XLue113eRUN9DU6hrk4beFSmQ6ycXpQLy0G7R0hhO4Jzx2FyLLzREej/xd5/1Lkufmoll0NQqmSlaGnXHZSoOMYieOIIROQW4G+AeuDLSqnPLfq8Gfg68BpgBHifUuqQ/dkngQ8Ds8BHlVKPetEmX3AXGIuCI2ho0ZKMRWior+MN27t5w/Zu/uKdV/PP+/N8/1k9SvjW4GG62pq45ep1bO1uo7Wpgbbmev27qZ7WZv17RVM9bU0NtDbX01Rfhyx+yGb6dSJVFLBycOXbwzv/mq06hyEKjjHsaTLn3GpW57qsuzq8dkCwwk2liMi0cilqdgQiUg98CXgLcATYLSIPLZKc/DBwWim1TUTuAP4KeJ+I9KM1jq8CNgA/FpEdSqloFnd3z/VdcWO4bbGGyo6Zb2qo48a+tdzYt5aJ6Vl+ujfPw88d4zu/PsLEdHmRDA11QmtTPW3NDQUH8d6pFXxg3OLj//AErOjS2y9xKNqRtDbq7zrHaG3S+9XXedCDP5+HC6fCG/qDjpfv2RGNm70QOhqmI3DdK1FwBGE6RdDnP/hTmJ3WGicRw4sRwU5gv1LqIICI3AfczkIR+tuBP7NfPwD8v6K7l7cD9ymlJoGXRGS/fbxfeNAu74mS8IY1DJe/qeKvtTTWc8vV67jl6nXMzinGp2a4MDnLhakZLkzNMj5p/7a3j9vbL0zNMG7vNz41y4XJGV45twWAyWN7eHI2W9h3dq78AlvNDXW0NtXbDkGoE6gTQezfoGeh6kQQ7G32Z3UCgnDdzLP8FfCpf53lucF/RkQvqOv9548nrmPU1envOtv1sbh01IPevuB9Ud8lfPjcWnbkn+Pj9wyW3K/YV4vvV72D/ID1c66q7+Tj3z0EHKr6OLXMstWraf6Gen70xBM89Nzl1R+oRurULF88uZcnJq/ie994KrR27Dy7gg/NTvHn9zzMiaba9JI/fetVrOto8ahlGi8cwUbAXXP2CPDaUvsopWZEZAzosrc/uei7G4udRETuAu4C6O0NsIjWYnr6wu/1XTwN547V3MuprxPaWxppb6myh3J2M/z1p/jim5th500AKKWYnJnj4tS8Eyk4l0nHqSx0LBemZplTijmlv68UzCmFwv5tb59TFLZh7zOnFNefOQ7jcHH1DnrqmwvfmbMrPs4pxdwcKPQxZucU07PKdQ7APv5i9B6u98X2sbfl5jbxutkfc+b0Kcbl0jpR5VaZrbVQZdeFAxyQXl46NV71MRb/3dVwuG4jHef2cWDufM3HqpbNs0doZJpnJtdzIB9eO2R2LR8C6k4Nc6CxtsoEUzPe5yPEZrFYKXU3cDfAwMBAeDVdM/3w9Dd0bHKQxbzcOOnqYQ79wRbe6FgwQhIRWhrraWmsZ01bQAXYHv4HuLCGv/5f3xpuNNfecbj3azzwO2ugd3FfKCDm5uD/PAavvpNH3/7GcNrg8O3XsOXY0/zoDysfuXrGnv8J34Y/+Xfv4k/CLJA4vRM++zH+jxsE3hzi9SiBF0+yo4A7Zm+Tva3oPiLSAHSgF43L+W60yGRhejxc4Y0oLAbCvPBG2COkMGPm3UQho/bMy3bMfMi2Afp/cvqQ1okIizCEm4rRuEIXa4zCtHIRvHAEu4HtIrJVRJrQi78PLdrnIeBO+/V7gJ8oPU5+CLhDRJpFZCuwHfiVB23yjyjokFo5aFoFHZvCa4ODk1EblvCGEzMfhQdfx2adyxC2bUC4C+cOhSi7ELPPrSH9AA5SuKkUEY4cqtkRKKVmgI8AjwI54H6l1B4R+YyI3Gbv9hWgy14M/hjwCfu7e4D70QvLPwR+P7IRQw5REN7ID2ujCrsHDC7hjRPhnP/sUZg8Gw1HUFen15DCLMFciJkPuQcMrk5TiI7AuVeiQKYfRg/A9ETYLbkET9YIlFKPAI8s2vanrtcTwHtLfPezwGe9aEcgtHRA+8bwPLtScHIPZG8N5/yLWSC8EUJuRZR6wKCvx4s/DO/8Vk5nfYcZM++wZovOdQmr0zQ9EY5wUykyfaDmYGQfrHtV2K1ZgCk6Vw1hDvHG83BxNDq9nLA1asMst1yMTFb/j87nwzl/fji8mjqLqauH7hBzK0ZCEm4qRRSmlUtgHEE1ZLJw6kUtvBE0UVkodmjr0nXew3QEK9dBa2c451+MO/s8aGantV1GxTYg3GCCqI0WO6+AusZILhgbR1ANmf554Y2giZpxQ7glmMOsM1+MMHt9owdhdip6tnHumNaLCBprSD94wxBuKkZDkxbqMSOChBBmmKCVg9YuneEcFTL9kN8bvPDG3Kw+b5QefCvXwoo14dzsURstQri6FdawFgxqCCifpRyiplthYxxBNYQpvBGVmHk3hdyKV4I97+lDMHMxWg++MHMrrBxInZ6XjwqhdpoiNloEvaZ25hWYDC/LuRjGEVRDU6uOiAjauJ2Y+agsjDqEJOM5r0oWtZvdLkMSdG6FldNVUKMQM+/QsUnnvARtG5PndXJdlEaL4FpD2htuOxZhHEG1ZPqDj48eOwJT56L54IPgHWOUYubdZLJaD+DssWDPG5XEOjciOoop8E6C/aCN2vWIQvZ5EYwjqJZMFkb2a+GNoIjiQjGEJ7xh5WB1LzSvCva8yxHGgvH0hE5WipptQDjz4lFcLwE7t2JF5BaMjSOoFrfwRlA4IYlRiRN3k8kGP0Jy1kuiRhghpKde1MlKUXvwgf4fXRgJNrciP7ykcFNo1NVHS9nPxjiCagljXtzKwaoNOiolamSyWq83qNyK2WnthKP44Gvt1LkNQdsGRNsxBvnwq0C4KXAy2XDrLxXBOIJq6doOdQ3BG3cURwOgoyFmp3QsexCMHIC56fBLcZci0xesbeRzOma+KyIx827CyD6P6mgR7NyK43BhNOyWFDCOoFoamnSMclCePYox826C7vVFdQ7YIejcCiunk5UiKIPIygys6AzONi6e1g/aKNsGRGpUYBxBLQS5CHb6EMxMRNe4e+zciqCMO4ox824yWa0LcOblYM4XxZh5Bye3IjDbcMKKTaepXIwjqIUghTcKc8ARvdmDFt6whnTpgEZvtVs9I8jIoclzOkkpqrYB84Uag8itKIQVR3QatX0jNLdHKnLIOIJacAwtiJ5O1KpsFiPIqqxRjJl34+Q2BOEYCzHzEe0Bg51bcVbrR/hNlISbiiFiJx2aqaFkEGSvzxqC1ZdB06Wi6JEhk9WLuH4Lb0xf1AX/ouwImldpXYBAbCMmnQQI7npERbipFGEr+y2iJkcgIp0i8piI7LN/XxLXKCLXicgvRGSPiDwnIu9zffY1EXlJRJ6xf66rpT2B07kV6psDNO4I9/hgPrdixOfciijHzLsJaoRk5XSSUtRi5t0ElX2uVLTXSxwy/VpX5LwVdkuA2kcEnwAeV0ptBx633y/mAvABpdRVwC3AF0Vktevz/6SUus7+eabG9gRLITnE55t9Zko/XONg3OD/kDfKMfNuMln9f5ud9vc8UY6Zd2jthFXr/beNgnBTDGwDIrNgXKsjuB24x359D/DOxTsopV5USu2zXx8DLCBCNZRrJIhKk6MHYG4m+sYdlPCGNQT1TXpxOspk+oPJrYjDaBGCibIrhBVHeJoMIqdWVqsjWKuUOm6/PgGsXWpnEdkJNAEHXJs/a08ZfUFEmpf47l0iMigig/l8SDKAxcj02cIbp/07R9Rj5h2CEt6wcjpsNIox826C6PVdGIXzJ6JvG+DKrZj17xxxGS2u7NG6ImEo2RVhWUcgIj8WkReK/CxQhFZKKaDkyoeIrAf+AfiQUsrJsvkk0AfcAHQCHy/1faXU3UqpAaXUQE9PhAYUQUyHWDmQev2QjTo9AWTUWsPRXhh16N6hcx38dIxRLcVdjJ4+rR9x+pB/57CGoifcVIowZTwXsawjUErdrJS6usjPg8BJ+wHvPOiLrnyISDvwA+BTSqknXcc+rjSTwN8DO734owIliAJjVk6XDmgoOWCKDpl+nUTll/DGxFktgBOHB19ji/+5FXEZLUIw0yFRFG4qRZC5FctQ69TQQ8Cd9us7gQcX7yAiTcD3gK8rpR5Y9JnjRAS9vvBCje0Jno7N0LQyAOOOwY0O8+085ZPwRhxi5t34XZXVyunkpPaN/p3DK5zcCr86TUrpax2ne2XqPIwdDrslNTuCzwFvEZF9wM32e0RkQES+bO/zu8AbgQ8WCRP9pog8DzwPdAP/pcb2BI+Iv2GC0xf1YmOcHnzg3/WIUw8Y9P9t1MfcijjEzDs0r9S5MH7ZhiPcFIdpQwinGF8JGmr5slJqBLipyPZB4Pfs198AvlHi+zfWcv7I0NMHwz/QPRKvb8j8XkDFx7jXbNF14H1zBDlobNUPlDjQ06dzHk69COuv8fbYTsx89jZvj+snfnaa4rJQ7OBENlk52PG2UJtiMou9wEkOGfchmiluxu238EYhZj4mpuvnvPj5kzpaLS62AbZuxYs6N8Zr4hI66rBijdYXicCIICZ3U8TxczokLjHzbvyMhsgPx+vB1+VjbkXUCxEWI9Ovc2JGDyy/b6VEWbipFGHIeBbBOAIv8LPXlx+G7iuhvqZZvGBxhDe8zq0YH9G94Dg9+OobdRipH4UJ4zZaBH87TfkYBVU4ZLL+51aUgXEEXuCn8EacIoYc/MqtyMewBwz+9fqsIWjt1slJcaFru86J8doRFISbYmgbs5Mw+lKozTCOwAv8ihyaOKtDy+Iy5+ngV4GxQpXNuN3sfVovYPKct8eNYyfBr9yKqAs3lSKIPKQyMI7AKxxBai+TQwpZozEa+oOuA9+0ynvHaA1Bcwe0b/D2uH5TkCb0MLdCKXu9JGYPPvCn0xS3sGKHHlfkUIgYR+AVfghvxNW4/RohOclCcYiZd+NHzaGxwzoZKW62AXZuxUGdI+MVcdBkKEZTmw65DnnB2DgCr/BjwdgahsY2LXASN7wW3ohLnflirN6i9QK8XDOJ40KxQyYLKB1G6hVWTj9QoyzcVIoI1BwyjsAr/JgXt4b0/HJcYubdeJ1bce4ETJyJ54Ovrk7/H722DYhfDxh86jTl4rd25NDTByP7/cmtKJMYPmEiSmsnrFxnjNsh47FjjFuy0GJ6PJ4qK8TMr/bumEHRebnOjfHKNuIi3FQKJ7diZH9oTTCOwEu8DBMcPwXjVryNG7x7+MV5KgT0//H8Ca0f4AVxnSYDnRPTvcM72xjZHw/hplJEQK3MOAIvyfRD/kVvkkPimDXqps0W3vDKuPM5fcy2bm+OFzReOsa5WW1ncbUN8DaYIK5BFQ7dPuVWVIBxBF6SyXonvBHX0FEHEXsRzKMF0jjGzLvxMl589CWdhBRX2wB9PcYO61yZWskPx0e4qRgNzdC1zTiCxOBlr88agpbVsGpd7ccKC6+EN+bm7NDRGD/42jfoHAivbANi7hg9zK2Ik3BTKbwOJqgQ4wi8xBHe8ORmj1Gd+VL09On68GNHajvO2CswPR7PCBkHEftm98g2YN7e4oiXUXZxXi9xyPTrmYSpC6GcviZHICKdIvKYiOyzfxct+ycisy5Rmodc27eKyC9FZL+IfMtWM4svzSthdW/txh3nmHk3Xo2Q4r5Q7OBVboU1FN+YeYfVl2ldiVptY+qCnipLgm2g/FP2W4ZaRwSfAB5XSm0HHrffF+OiUuo6+8etovFXwBeUUtuA08CHa2xP+GT6a680ee44TIwlwLg96vUVHEGMRwRg51ac1hVUayFupbiLUVenRwW12sYpW7jJdJpqolZHcDtwj/36HrTucFnYOsU3Ao6OcUXfjyxeCG/EPWLIwRHeqNUxWjlo3wQtHd60Kyy8KME8M6nDJeNuG+BNp8mKeVCFw5qtUN8c2jpBrY5grVLquP36BLC2xH4tIjIoIk+KyDvtbV3AGaXUjP3+CFBSgVtE7rKPMZjP+6AE5hVeCG/EtcpmMbzIrYh7xJCDF72+uMfMu8lk9ehofKT6YzjCTWu2eteuMCjkVvigW1EGyzoCEfmxiLxQ5Od2935KKQWUmvy8TCk1ALwf+KKIXFFpQ5VSdyulBpRSAz09Ea6/7sUimJWDtgy0dXnTpjCpVXhjdkYP/+M+LQQ6B6K1u3bbgHgvnDs4/9NaQmqtXPyEm0rhp57zMizrCJRSNyulri7y8yBwUkTWA9i/rRLHOGr/Pgj8FLgeGAFWi4jzH9wEeFi6MyS6d4DU1fYPTcJCsUMmq+vEV5tbMXoQZqeS0QOG2m92ayjeMfNuvBghJWW0CPrvOHtErw8GTK1TQw8Bd9qv7wQeXLyDiKwRkWb7dTfwemDIHkE8Abxnqe/HjsYW6LyieuOem0vGYqBDrenzcVUlK4UzL15t5JA1rJOP4hwz77BqvV73qdY2Jsb0gzNJtgGhTA/V6gg+B7xFRPYBN9vvEZEBEfmyvU8WGBSRZ9EP/s8ppZz//MeBj4nIfvSawVdqbE80qKXXN/YKTF9IjnEXpsqqNG4rB4ge/ieBTFbrCIwdru77SRot1pp97iSjmU5TzdQ0saaUGgFuKrJ9EPg9+/W/Aq8q8f2DwM5a2hBJMv2Qe1gLbzSuqOy7SYmZd6hVeMMags6t0NTqabNCwz0dsrpCnYmpcT3Fdu0uz5sVGpksvPBdPUKqNHky7hVpF9OxWeuP1BpJVQUms9gPMn2Aqi59vlBnPiE9YKitBHOcS3EXo5B9XoVjzDsx8wl58IH+306c0XoTlWLl4ivcVAw/dCvKPXXgZ0wDtSyCWTndM2hp97ZNYZLJ6nrxleZWzEzCyIHkTIWA1g9o31i9bUByRotQ23RInIWbShFS5FCCrmCEcIQ3qgmLS1IUhEO1whun9oGaTeD1qDK3Ip/TSUdxj5l3U0uSXVLvlfE8nA82V8o4Aj+ob6xOeGN2RmclJ864qyzBnMQeMNi5FVXoVlg56NmRjJh5h7ZunTNTqW2Mn9IPzCTaBnhTrrwCjCPwi2qGeEmLmXeoVnjDGoK6Bh0umSQy/VpPYPSlyr5n5ZJnG1DdvZKkxDo3znpYwCGkxhH4RU9f5cIbcRYkX4qGZl0vvpqbvWsbNMS7KO0lVJN9fvEMnD2aPNsA2xEM6xyacknqaHHVOq1DEvCCsXEEflEQ3qjAszsx80mKGHKoZl48STHzbnquBKQyxxh3xbqlyGS13sTYK+V/JwnCTcUo5FaYqaFkUM0imDWkF5orzT2IA5l+PRVSrvDG1DiceTmZD75qciuSUpG2GNVE2TnTZHEWbiqFV8p+FWAcgV9UI7yRH07mjQ4u4Y0Xy9u/0ANO6vWosASzlYOmlTq0OGkUpsrKvFeU0oupibWNLEyOwdljgZ3SOAK/qFR4Y3rCjplPYA8YKu/1JXUO2CGT1eG0M5Pl7W8NaXtKUsy8Q0u7dnDl2kZBuCnBjgACjRxKoFVFiEqiIUacmPkELgaCLbzRVL5jtHLQ0KKnUJJIJltZboWVS65tgN1pKreT4JSWSKgj6Kkht6JKjCPwk0wWxq3yhDeS3gOub9CF4yq52bt3QF29v+0Ki0rWkM7n4cKp5NoG2Mp+e3UuzXIkSbipGG1dsHKtcQSJoZIhnpWDukZdwjqpZLLlz4tbCSrFXYyu7TpHopybPWmluIuR6dc5NKfLyK2whvWDMgnCTaXwQtmvAowj8JNK5sWtnE68SlrMvJtMtrzcioun4dyxZD/4Gpp0jkS5tgHJdoyV1BxKalixmx5H2a+C3IoaMI7AT1ath+YyhTfSYNzl5lYkRZB8Ocrt9VlDsGKN7gUnlXJzKxzhpqROCzlkslqX5MzLgZzOOAI/ESlvwXjyvP6HJ964y8yoTVqd+VL0ZLW+wNT40vs5pbiTGDPv0LhC604sZxtnXk6WcFMpvJDxrICaHIGIdIrIYyKyz/69psg+bxaRZ1w/EyLyTvuzr4nIS67PrqulPZHE6fUtlRxSUFpKuHF39Or68csZd5Jj5t04uRVL6VYolcwqm8UoJ6M2DdNkUJtuRRXUOiL4BPC4Umo78Lj9fgFKqSeUUtcppa4DbgQuAD9y7fKfnM+VUs/U2J7okenXMc9LCW+kYTEQXMIby9zsTmJdknvAUN5U2dljMHk2+bYBdm7FgaVzK5x7JYllWNy0tOuOUxxGBMDtwD3263uAdy6z/3uAf1RKlVlnIAGUswhm5aBhRXJj5t0sN1WmFJzck44HX+dWrS+wnG1A8nvAoP/nalbrUJTCyukHZJKEm0pRTqfJI2p1BGuVUsft1yeA5Vaz7gDuXbTtsyLynIh8QUSaS31RRO4SkUERGczngxVtqIly4sWtId3DSWrMvJtMv51bcar45+N5uDiajgdfXb3+vy9nG5AOx1jOvHjSE+vcOMp+s9O+n2pZRyAiPxaRF4r83O7eTymlgJIT4SKyHi1i/6hr8yeBPuAGoBP4eKnvK6XuVkoNKKUGenp6lmt2dGjrhraeMow7BTc6LF9XJqmluEux3AjJyuloodbO4NoUFp1X2LkVJUZIs9PJFG4qhZNbMXrQ91Mt6wiUUjcrpa4u8vMgcNJ+wDsPemuJQ/0u8D2lVMG9KaWOK80k8PfAztr+nIiyVJjgxdO6dkqajBuWcAQpmgoB/X8/e1TrDRQjDWHFDg1NOtGulG0kVbipFLXoOVdIrVNDDwF32q/vBB5cYt9dLJoWcjkRQa8vvFBje6JJpr90ckhaYuYdHOGNUtnWVg5WdMLKTKDNCo3CgnGRyKG5Ob09LbYBdvb5cp2ElDjG7h0gdYGsE9TqCD4HvEVE9gE32+8RkQER+bKzk4hsATYD/7To+98UkeeB54Fu4L/U2J5ospTwRprmgGF54Y0k15kvxlK9vjOHYOZiemwD9P++VG6FldMPxu4dgTcrFBpX6GKNAYwIalLBVkqNADcV2T4I/J7r/SFgY5H9bqzl/LHBXU1wcWSQlYPmdmi/5PIkl0wfvPAdHSHkfuA7MfPXvi+8tgVNx2adM1HMMSa9uFoxnIXg/DBsfM3Cz6wh/WBMonBTKRwZT58xmcVBsFRGrZXTC6Np6QGDK7fi+MLtY0dg6ly6esAipXUrCgvnCY+Zd7PUGlKagiocMv0wekDrlfiIcQRB0NIB7ZsuNW6l0rUY6FBqOiRtC8UOpSKH0hQz77Bmi9ahWHw9pif0AzGNtqHmylf2qxLjCIKi2BAvTTHzbgpTZYuuRyFrNCWhow6Zfq03cH5RfoyVYOnSUpTKrRjZpx+IabseAdUcMo4gKIoJb6RtodihlPCGldMVW9MQM++mmG5F2mLm3RQLJkjraLHrCq1T4vOCsXEEQZHJXpockrZwODfF5sUdXd60USz7fOQAzE2n1zbOHdM5Ng7WkH4gdiVYuKkY9Y1ap6RcQacqMY4gKIrNi1tD0NqlM4/TRqZfG7eTWzE3m76YeYeVa7XewGLbgHQ6gsJ0iOvh5wg31TeG06YwCUCtzDiCoOguIryRtph5N4uFN04fgpmJdD74iuVWpC1m3k2pTlMabQP0333mFZg859spjCMIiqZWXW3SmQdWKp2LgQ6LSzCndQ7YwQkmcHQr8jnovDxdMfMOHZugadW8bUye1w/C1N8rS+hW1IhxBEHi7vWlMWbezWLhjULyVIpi5t1ksjA5pvUHIJ0x8w6Llf0Kwk0p7ST0LJGH5BHGEQRJT59eBJyeSGfWqJuWdp1V61wHawhW90LzynDbFRbu7PPpizqoIK22AToJ8+Se+VwbSGcgAdi5FSt8zTA2jiBIHOGNkX3p0eVdCnevz1kvSSvuefFTL6YzZt5Npl/n2Izn0yXcVIxCboUZESQDd3KIlYNVG3S0SFrJZPVDb/qido5pfvC1dsLKdfO2AcYxgn74pUm4qRTl6DnXgHEEQdK1zRbeyOnFwDQ/+GBeeGPfj2BuJt0PPpgvwWzl0hkz78YdQpofNraRycL5E3Bh1JfDG0cQJI7wxskX7Jj5tDsC++9/4TsL36eVTL9+8J3co8NG0xgz79DWo3NsXv6XdAk3laIcydsaMI4gaDJ98NLP0hsz78YR3njxUf27a3vYLQqXTJ/WHzj083SvHYFdlTWrR4tg7pViZUg8xDiCoMn0aycAxrgd4Y2ZCa1X29gSdovCxZn+MJ0ETSZr7hWH9o1atySKIwIRea+I7BGROREZWGK/W0Rkr4jsF5FPuLZvFZFf2tu/JSJNtbQnFrgNOq3hcG6c65H2Gx0W5lCkfU4c5m0ibcJNxVicW+ExtY4IXgB+B/hZqR1EpB74EvB2oB/YJSKOlf8V8AWl1DbgNPDhGtsTfZwbfM0WaGoLtSmRwLke5sEHzat0LgUYxwgu28imswzLYpxCjU72uYfU5AiUUjml1HJ5zzuB/Uqpg0qpKeA+4HZbsP5G4AF7v3vQAvbJxhHeMA8+jRkRLCTTr2PmV28JuyXh46yTGNvQZPp1RdbzJz0/dE2axWWyETjsen8EeC3QBZxRSs24tpcc/4nIXcBdAL29vf60NAjq6uFtf5neUgqL2f4WeN1HYNsl0tfp5Df/ALa/FerM8h0r1sDNfw5XvDnslkSDTQNw1bvm1008ZFlHICI/BtYV+ehTSqkHPW9RCZRSdwN3AwwMDHg/NgqSG5I/A1Y2zavgbZ8NuxXRYcsb9I9B84Y/CrsF0WHTALz3a74cellHoJS6ucZzHAU2u95vsreNAKtFpMEeFTjbDQaDwRAgQYw/dwPb7QihJuAO4CGllAKeAN5j73cnENgIw2AwGAyaWsNH3yUiR4DXAT8QkUft7RtE5BEAu7f/EeBRIAfcr5TaYx/i48DHRGQ/es3gK7W0x2AwGAyVI8qHUCS/GRgYUIODg2E3w2AwGGKFiDyllLok58uEJhgMBkPKMY7AYDAYUo5xBAaDwZByjCMwGAyGlBPLxWIRyQMvV/n1buCUh82JO+Z6zGOuxULM9VhIEq7HZUqpnsUbY+kIakFEBoutmqcVcz3mMddiIeZ6LCTJ18NMDRkMBkPKMY7AYDAYUk4aHcHdYTcgYpjrMY+5Fgsx12Mhib0eqVsjMBgMBsNC0jgiMBgMBoML4wgMBoMh5aTKEYjILSKyV0T2i8gnwm5PWIjIZhF5QkSGRGSPiPxh2G2KAiJSLyJPi8j3w25L2IjIahF5QESGRSQnIq8Lu01hISJ/bN8nL4jIvSLSEnabvCY1jkBE6oEvAW8H+oFdIpJW4eAZ4D8qpfqB3wB+P8XXws0fokulG+BvgB8qpfqAa0npdRGRjcBHgQGl1NVAPVpTJVGkxhEAO4H9SqmDSqkp4D7g9pDbFApKqeNKqV/br8+hb/KSetFpQEQ2Ab8NfDnstoSNiHQAb8TWB1FKTSmlzoTaqHBpAFaISAPQChwLuT2ekyZHsBE47Hp/hJQ//ABEZAtwPfDLkJsSNl8E/gSYC7kdUWArkAf+3p4q+7KItIXdqDBQSh0F/m/gFeA4MKaU+lG4rfKeNDkCwyJEZCXwHeCPlFJnw25PWIjIOwBLKfVU2G2JCA3Aq4H/ppS6HhgHUrmmJiJr0DMHW4ENQJuI/NtwW+U9aXIER4HNrveb7G2pREQa0U7gm0qp74bdnpB5PXCbiBxCTxneKCLfCLdJoXIEOKKUckaJD6AdQxq5GXhJKZVXSk0D3wV+M+Q2eU6aHMFuYLuIbBWRJvSCz0MhtykURETQ8785pdRfh92esFFKfVIptUkptQVtFz9RSiWu11cuSqkTwGERudLedBMwFGKTwuQV4DdEpNW+b24igQvnDWE3ICiUUjMi8hHgUfTK/1eVUntCblZYvB74d8DzIvKMve0/K6UeCa9JhojxB8A37U7TQeBDIbcnFJRSvxSRB4Bfo6PtniaBpSZMiQmDwWBIOWmaGjIYDAZDEYwjMBgMhpRjHIHBYDCkHOMIDAaDIeUYR2AwGAwpxzgCg8FgSDnGERgMBkPK+f8BwaSPXwFUoZ4AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "dt = 1\n", "y0 = 1\n", "t = np.arange(0,10,dt)\n", "tru_y = np.exp( -params[0]*t)\n", "(t,y) = forwardEuler1D( degradation_model, t, y0=y0, params=params)\n", "\n", "plt.plot(t, tru_y)\n", "plt.plot(t, y.T);" ] }, { "cell_type": "code", "execution_count": 85, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEkCAYAAAA/7cqkAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy86wFpkAAAACXBIWXMAAAsTAAALEwEAmpwYAABjS0lEQVR4nO3dd3xO1x/A8c/JkAQRhBgJSWxZImLvVXtTq2ZRtVpdtFptaVFVq7Rq7xn92ZuEqJkSKwgiJBFEEiTIPr8/nngaZHsi67xfr+fV59577rnf56l8c3PvPd8jpJQoiqIouZ9edgegKIqi6IZK6IqiKHmESuiKoih5hEroiqIoeYRK6IqiKHmESuiKoih5RL5J6EKIZkKIXRncZ7AQomwG2lcTQpwUQkQLIb5IpZ2tEOK0EOKmEGKTEKJA4nqjxOWbidttkuzzdeL660KINhn5HIqi5A/5JqFn0mAg3QkdCAPGAbPSaPcLMEdKWQkIBz5MXP8hEJ64fk5iO4QQdkAfwB5oC/whhNDPQFyKouQDeTKhCyGmCCE+TbL8M1ADKCyEcBNCXBNCrBNCiMTtk4UQZ4UQl4UQi4VGT8AVWCeE8BZCmKR1XCnlQynlWSA2ldgE0AJwS1y1Cuia+L5L4jKJ21smtu8CbJRSRkspbwM3gTrp/DoURckn8mRCB5YDAwGEEHpozm4DgZrAp4AdUAFomNh+gZSytpTSATABOkop3QAvoL+U0llK+UIIMScxub/+mpiB2MyBx1LKuMTlQMAy8b0lEACQuP1JYnvt+mT2URRFAcAguwPIClJKfyFEqBCiJlAKOA+EAmeklIEAQghvwAY4DjQXQnwFFASKA1eAncn0O/6dfABFUZRMyJMJPdFSNNfAS6M5YweITrI9HjAQQhgDfwCuUsoAIcQPgHFyHQoh5gDNk9m0UUo5I51xhQJFhRAGiWfhVkBQ4rYgoBwQKIQwAMwS279c/1LSfRRFUYC8e8kF4H9obiDWBvan0u5l8n4khCgM9EyyLQIwfbkgpRyfePnl9Vd6kzlSUw3NPclxBgHbE9/vSFwmcfuRxPY7gD6JT8HYApWBM+k9pqIo+UOePUOXUsYIIdzRXK+OT7z/mVy7x0KIJcBl4D5wNsnmlcAiIcQLoL6U8kVqxxRClEZz3b0IkJB4Y9ZOSvlUCLEHGCalvAdMADYKIX5CczloWWIXy4A1QoibaJ6Y6ZMY4xUhxGbAB4gDRksp4zP2jSiKkteJvFo+N/Fm6Dmgl5TyRnbHoyiKktXy5CWXxOe2bwKHVTJXFCW/yLNn6IqiKPlNnjxDVxRFyY/ybEIXQowTQlxNHBHaVQgxOY32HkII1wz0byOE6Jdk2VEIsfIt4u0vhLgohLgkhDghhKiRQjtVB0ZRlGTl2YQOjAJaSyn7A1+hedZcl2wAbUKXUl4CrIQQ5TPZ322gqZTSEZgKLE6hnaoDoyhKsvJkQhdCLEIztH+vEGICEC2lfCSEMBVC3BZCGCa2K5J0GeglhDgjhPAVQjRObGMjhPAUQpxLfDVIbDsDaJw49P/lCNKdJD5qmFFSyhNSyvDExVNoBg+9/rlUHRhFUVKUJxO6lHIkcA/NqM6HaB5fREoZAXgAHRKb9gH+llK+LKZlIKWsg6bey/eJ6x6iOdN3AXoD8xPXTwQ8EwcWzUlc5wW8/EXQPIW6LyfS8RE+BPYms17VgVEUJUV5dmBREmWAkCTLS9FcgtkGDAGGJ9n2d+J//0VzSQXAEFgghHBGUy6gSirHekhiuV0ppTvgnNFghRDN0ST0RhndV1GU/C0/JPQXaGqiACCl/CfxMkozQF9KeTlJ25e1XuL577sZDzxAU35XD4hK5VjGicd7mZjnJNPmuZSygRBiNP/9MmkvpbwnhHBC8wunnZQyNJl9VR0YRVFSlCcvubzmKlDptXWrgfXAinTsbwYESykTgAHAyxuKr9R5SVQFTQkBpJTuKdR9aZC4fWGSdfcSb6b+DQyQUvomF4iqA6MoSmryQ0I/BtR8OZlFonVAMWBDOvb/AxgkhLgAVAOeJa6/CMQLIS4kuSnaHNidyTgno7nm/UfitXavlxuEEHvEf1PhTQA+S6z3Ys6rdWDME9d/huYaP1LKK8DLOjD7UHVgFCXPyhcjRYUQ84CdUspDics9gS5SygE6PIYRcBRolOSmpaIoyjuTH66hA0wD6gIIIX4H2gHtdXyM8sBElcwVRcku+eIMXVEUJT/ID9fQFUVR8gWV0NNJCFFRCFEyu+NQFEVJiUro6SCEaAic5M3HHxVFUXIMldDTIIRojWZU6QAp5clsDkdRFCVF+eUpl0wRQnQBlgDdpJTHszseRVGU1Kgz9BQk1jr/C80wfJXMFUXJ8VRCT4YQYgQwE2gppfw3u+NRFEVJD3XJ5TVCiM+BMUAzKeXN7I5HURQlvVRCT5RY6+V7oC/QREoZkMYuiqIoOYpK6GiT+SygFZpk/iCbQ1IURcmwfJ/QE+fX/BNwAppLKcOyOSRFUZRMydcJPXEu0ZVoZhlqnThFnaIoSq6UbxO6EMIY2Ihmirn2UsoX2RySoijKW8mXjy0KIQoBO4EYNIOGVDJXFCXXy3cJXQhRFDgABAB9pZQx2RuRoiiKbuSrhJ5YLfEI8C8wTE3FpihKXpJvEroQwhLNFHF7gU8SJ31WFEXJM/JFQk+c7f4YsFpKOUmqaZoURcmD8nxCF0JUR5PMf5NSzsjueBRFUbJKnn5sUQhRE9gDTJBSrs7ueBRFUbJStk0SXaJECWljY5Nl/UdGRnLr1i3Kly9PsWLFsuw4iqIo79K///77SEqZ7HSY2XaGbmNjg5eXV5b0ffjwYfr27cvOnTtp27ZtlhxDURQlOwgh7qS0Lc9dQ9+5cyd9+/bFzc1NJXNFUfKVXJ/QL1y4gLu7OwAbN25k+PDh7N69myZNmmRzZIqiKO9Wrr8pOnXqVNq0acOtW7f4/vvvOXjwII6OjtkdlqIoyjuXqxN6eHg4Bw8exMXFhcWLF+Ph4UHlypWzOyxFeUVsbCyBgYFERUVldyhKLmJsbIyVlRWGhobp3ifNhC6EWA50BB5KKR2S2S6AeUB74DkwWEp5Lt0RvIXNmzdjZWXF8uXL+frrr5k1axaff/45VapUeReHV5R0CQwMxNTUFBsbGzQ/LoqSOikloaGhBAYGYmtrm+790nMNfSWQ2t3FdkDlxNcINJNFvBM///wzN27c4P79+2zduhUHBwfKly//rg6vKOkSFRWFubm5SuZKugkhMDc3z/BfdWmeoUspjwkhbFJp0gXNkHoJnBJCFBVClJFSBmcokkxwdXWlU6dO9OjRgyJFimT14RQl01QyVzIqM/9mdPGUiyWaUrQvBSauyxKBNy9zcvFY4uPi+PvvvxkyZIhK5oqSitDQUJydnXF2dqZ06dJYWlpql2NidFM9eteuXdSsWZMaNWpgZ2fHX3/9lWr7H374gVmzZqXaxtvbmz179miXd+zYwYwZuqne8eLFC5o2bUp8fMoFVy9dusTgwYN1crx35Z3eFBVCjEBzWSbTl0YCT26m/r3VnJ99k6qjN1GwkErmipIac3NzvL29AU0iLVy4MF988YV2e1xcHAYGmU8FsbGxjBgxgjNnzmBlZUV0dDT+/v5vGbUmoXt5edG+fXsAOnfuTOfOnd+6X4Dly5fTvXt39PX1U2zj6OhIYGAgd+/ezTWXcnVxhh4ElEuybJW47g1SysVSSlcppWvJksmOXE1TvQFTOF31K2o8O0nAnFY8epDsoRRFScXgwYMZOXIkdevW5auvvnrjjNnBwUGblNeuXUudOnVwdnbmo48+euOsNiIigri4OMzNzQEwMjKiatWqAPj7+9OiRQucnJxo2bIld+/efSOWZs2aaUeNP3r0CBsbG2JiYpg8eTKbNm3C2dmZTZs2sXLlSsaMGZNqv4MHD2bcuHE0aNCAChUq4ObmluznX7duHV26dAFg4MCBbNu2Tbutf//+bN++HYBOnTqxcePGDH232UkXCX0HMFBo1AOeZPX187p9J3Gp4XysY/14sagld25czMrDKUqeFBgYyIkTJ5g9e3aKba5evcqmTZv4559/8Pb2Rl9fn3Xr1r3Spnjx4nTu3Blra2v69u3LunXrSEjQTDcwduxYBg0axMWLF+nfvz/jxo1LV2wFChRgypQp9O7dG29vb3r37v3K9tT6DQ4O5vjx4+zatYuJEye+0XdMTAx+fn68rCX14YcfsnLlSgCePHnCiRMn6NChA6C5T+fp6ZmumHOC9Dy2uAFoBpQQQgQC36OZWBkp5SI01QzbAzfRPLY4JKuCTarGewO5UbwsJXcNJmFde660W4l93Vbv4tCKkmk/7ryCz72nOu3TrmwRvu9kn+H9evXqleolB9DURfr333+pXbs2oLn2bGFh8Ua7pUuXcunSJQ4dOsSsWbM4ePAgK1eu5OTJk/z9998ADBgwgK+++irDcSYntX67du2Knp4ednZ2PHjw4I19Hz16RNGiRbXLTZs2ZdSoUYSEhLB161Z69OihvQRlYWHBvXv3dBLzu5Cep1z6prFdAqN1FlEGVHZtRXCxvci1Pamwpy+nw36jbruB2RGKouQ6hQoV0r43MDDQnlUD2sflpJQMGjSI6dOnp9mfo6Mjjo6ODBgwAFtbW+1Zb1qSHlsXg6+MjIy075OrJmtiYvLGcQYOHMjatWvZuHEjK1as0K6PiorCxMTkrWN6V3L1SFGAMhUdeTrKnaDFXXA9NY6joYE06f+1ekxMyZEycyb9LtjY2LBr1y4Azp07x+3btwFo2bIlXbp0Yfz48VhYWBAWFkZERATW1tbafSMjI/Hy8qJZs2aA5mbmy+0NGjRg48aNDBgwgHXr1tG4ceNkj/3vv/9Sp06dV655m5qaEhERkWy86ek3JcWKFSM+Pp6oqCiMjY0BzbX3OnXqULp0aezs7LRtfX19cXB4YzxljpXri3MBFClZlnLjD+NjWp+mN3/h6MJRxMbFZXdYipJr9OjRg7CwMOzt7VmwYIF2tLWdnR0//fQT7733Hk5OTrRu3Zrg4FdvkUkpmTlzJlWrVsXZ2Znvv/9ee3b++++/s2LFCpycnFizZg3z5s1749hffPEFf/75JzVr1uTRo0fa9c2bN8fHx0d7UzSp9PSbmvfee4/jx49rl0uVKkX16tUZMuTVK8bu7u7a6+m5QbZNcOHq6ip1XQ9dxsdyYfFHOD/YykmTZtiNWoeZaWGdHkNRMurq1atUr149u8NQkjh37hxz5sxhzZo1ADx//hxHR0fOnTuHmZkZANHR0TRt2pTjx4+/1WOdbyO5fztCiH+llK7Jtc8TZ+gvCX1DnEcu41K1T6n/wgP/uW0IzEU3NBRFeTdcXFxo3rw58fHxHDp0iOrVqzN27FhtMge4e/cuM2bMyLZknhl56gw9qRuHlmNz/AvuUobnvTbiaK9K6irZQ52hK5mVr8/Qk6rcaigPu2zAgjAsNnfi6NHD2R2SoihKlsqzCR3AsmYbEgbvQ19fn1pH+rFj6+pkH2NSFEXJC/J0Qgcws6mB6RgPnhhb0v7iJ2xeMoOYuIS0d1QURcll8nxCBzAqXo6y490JKuZK73sz2DZ3LGGR0dkdlqIoik7li4QOIIzNsB6zm7vluvB+5FpOzunLjXth2R2Worwz27ZtQwjBtWvXsi2Ge/fu0bNnT530tW3bNnx8fDK0z4MHD+jYsaO2zO/LSo4p8ff3T9fAomnTpr2y3KBBgwzFpSv5JqEDYFCA8kNXEew8jg7xh3nwV1c8L/tld1SK8k5s2LCBRo0asWHDBp31GZfBAXxly5ZNsQJiRmUmoU+ePJnWrVtz4cIFfHx8dFZf/fWEfuLECZ30m1H5K6EDCEGZrlN53Oo36otLFNvcjY1HzqibpUqeFhkZyfHjx1m2bNkr5WA9PDxo0qQJHTp0oGrVqowcOVJbV6Vw4cKMHz8ee3t7WrZsSUhICKApd/vpp5/i6urKvHnzOHz4MDVr1sTR0ZGhQ4cSHR3N2bNncXJyIioqimfPnmFvb8/ly5dfOeNduXIlXbt2pXXr1tjY2LBgwQJmz55NzZo1qVevHmFhmr+glyxZQu3atalRowY9evTg+fPnnDhxgh07dvDll1/i7OzMrVu3uHXrFm3btqVWrVo0btw42b9EgoODsbKy0i47OTkBmtGuX375JQ4ODjg6Or4xMvVlvC/L9wJ07NgRDw8PJk6cyIsXL3B2dqZ///7a7y61fj08PGjWrBk9e/akWrVq9O/fXzc5SEqZLa9atWrJ7Pbiyl754gcLGTC5gpyzfruMiYvP7pCUPMjHxye7Q5Br166VQ4cOlVJKWb9+fenl5SWllNLd3V0aGRnJW7duybi4ONmqVSu5ZcsWKaWUgFy7dq2UUsoff/xRjh49WkopZdOmTeXHH38spZTyxYsX0srKSl6/fl1KKeWAAQPknDlzpJRSTpo0SX7++edy1KhRctq0aVJKKW/fvi3t7e2llFKuWLFCVqxYUT59+lQ+fPhQFilSRP75559SSik//fRTbT+PHj3Sfo5JkybJ+fPnSymlHDRokDZWKaVs0aKF9PX1lVJKeerUKdm8efM3vod9+/ZJMzMz2axZM/nTTz/JoKAgKaWUbm5uslWrVjIuLk7ev39flitXTt67d++NeF9+B1JK2aFDB+nu7i6llLJQoUKvHOflckr9uru7yyJFisiAgAAZHx8v69WrJz09Pd+IN7l/O4CXTCGv5p4hUFnA2K4tCR/uo+iq7gy59hHT//yRcUMHU7RggewOTcmr9k6E+5d022dpR2iX+qWDDRs28MknnwDQp08fNmzYQK1atQCoU6cOFSpUAKBv374cP36cnj17oqenp61D/sEHH9C9e3dtfy/XX79+HVtbW23tl0GDBrFw4UI+/fRTJk+eTO3atTE2Nmb+/PnJxtW8eXNMTU0xNTXFzMyMTp06AZrKjRcvauY5uHz5Mt9++y2PHz8mMjKSNm3avNFPZGQkJ06coFevXtp10dFvPvjQpk0b/Pz82LdvH3v37qVmzZpcvnyZ48eP07dvX/T19SlVqhRNmzbV/pXxNlLqt0iRItSpU0f714KzszP+/v40atTorY6XrxM6gJ5VTQqP8uDpsq5MDPmaGfPu0+/Dz6hkoWrAKHlDWFgYR44c4dKlSwghiI+PRwjBr7/+Crw5GXFKlUqTrk9aejcloaGhREZGEhsbS1RUVLL7JC11q6enp13W09PTXp8fPHgw27Zto0aNGqxcuRIPD483+klISKBo0aLaqfZSU7x4cfr160e/fv3o2LEjx44dS3MfSLnEcGYl/ez6+voZvh+RnHyf0AEoZk2RUYeJWPU+kx/8xm8LH+DSZzLNq5fK7siUvCaNM+ms4ObmxoABA16ZuLlp06bamXjOnDnD7du3sba2ZtOmTYwYMQLQJEk3Nzf69OnD+vXrkz17rFq1Kv7+/ty8eZNKlSqxZs0amjZtCsBHH33E1KlTuX37NhMmTGDBggWZij8iIoIyZcoQGxvLunXrsLTUzEGftLxukSJFsLW1ZcuWLfTq1QspJRcvXqRGjRqv9HXkyBHq1atHwYIFiYiI4NatW5QvX57GjRvz119/MWjQIMLCwjh27Bi//vrrK0nbxsaGP/74g4SEBIKCgjhz5ox2m6GhIbGxsRgaGr5yvJT6zaonjfLfTdGUFCyO6bBdPK/cmc/FWgLWj2axh6+6Warkehs2bKBbt26vrOvRo4f2aZfatWszZswYqlevjq2trbZtoUKFOHPmDA4ODhw5coTJkye/0bexsTErVqygV69eODo6oqenx8iRI1m9ejWGhob069ePiRMncvbsWY4cOZKp+KdOnUrdunVp2LAh1apV067v06cPv/76KzVr1uTWrVusW7eOZcuWUaNGDezt7bXzgib177//4urqipOTE/Xr12fYsGHUrl2bbt264eTkRI0aNWjRogUzZ86kdOnSr+zbsGFDbG1tsbOzY9y4cbi4uGi3jRgxAicnJ+1N0ZfS068u5dniXJmWkEDs/m8xPL2QA/G1OGQ3jSk9a2NsmPpUXYqSkpxcnMvDw4NZs2ZpJ7dIqnDhwkRGRmZDVMpLqjjX29LTw7DdNGTbGbTWP0efq6MZvmg/D56+/dRYiqIoWUkl9BSIeh8j3l9NDYMAfno0no/nb8E74HF2h6UoOtWsWbNkz84BdXaeC6mEnhq7zugP3omlSTTL4r7hp79W8/e5wOyOSlEUJVkqoaelfF0Mhh2iiFlx1hn+xF63ZUzd5UNcvKrYqChKzqISenqUqIT+8EMYlrHnrwJziTn5FwOXnyHsWUx2R6YoiqKlEnp6FS6J3uBd6FVpw1TDlTQLWECX34/hc+9pdkemKIoCqISeMQUKQe+14PohI/R28l3MHHr/6cHOC2oiaiXnU+VzVflcAIQQbYUQ14UQN4UQE5PZbi2EOCyEuCiE8BBCWCXXT56gbwAdfoOW3/NewnE2Gs9k0gZPZuy9RnyCGoSk5FyqfK4qn4sQQh9YCLQD7IC+Qgi715rNAlZLKZ2AKcB0XQeaowgBjT+D7kuxi7/GQbNp7Dx6iiErz/L4ubquruQ8qnyuRr4vnwvUB/YnWf4a+Pq1NleAconvBfA0rX5zQvlcnfA7KuW0cvLZtIqy0zcLZaNfDsvLQY+zOyolB1Hlc1X53JxUPtcSCEiyHAjUfa3NBaA7MA/oBpgKIcyllKGZ/D2Te9g2gaH7KLiuJ/+L/4nxMZ/T489opnd3pFvNvHvlScmcX878wrUw3V7Drla8GhPqTEi1jSqfq6HK56bPF8ACIcRg4BgQBMS/3kgIMQIYAVC+fHkdHToHKGUHww6hv64X80Km8WfxcYzflMCFgCdM6lAdQ31171nJPqp87qvye/ncIKBckmWrxHVaUsp7aM7QEUIUBnpIKR+/3pGUcjGwGDTFuTIXcg5VpCwM2YvYPIBRfrOxrzCMQSeac+XeExb2c8GiiHF2R6jkAGmdSWcFVT73P6p8LpwFKgshbIUQBYA+wI6kDYQQJYQQL/v6Gliu2zBzCeMi0G8L1OhL03tLOVp1K9eCwuj4+3H+vROW3dEp+ZQqn/sfVT4XEEK0B+YC+sByKeXPQogpaC7O7xBC9ETzZItEc8lltJTyzQtYSeTY8rm6ICW4T4NjM4ks15SeoSO5+Ri+62jHwPrWKf5Jq+RNqnyuklkZLZ+brmvoUso9wJ7X1k1O8t4N0M3DpXmBENBiEphZUXjXeHZbhPFF8W/5fscVzt0NZ3p3RwoWUJNFKYqiW+puXVaqNQj6bUI/zI/ZEZ/zc0MDdly4R5cF/3DzoTrzUbKfKp+bt6iEntUqt4YhuxFxMfS/MpztHSD0WQxdFhxn98Xg7I5OUZQ8RCX0d6FsTRh2CAqXwsl9MIfee0SV0qaMXn+Oqbt8iFWleBVF0QGV0N+VYtYwdD9YulJ870i2OHoxuL41y47fpu/iU2qKO0VR3ppK6O9SweIw4H9g1xWDw5P5wXAV8/s44RP8lA7zPTlx81F2R6goSi6mEvq7ZmgMPVdA/TFwZjGdr01kxwgXzEwM+WDZaeYfvkGCqtqoZAFVPleVz1Wygp4etPkZ2v4C13ZTaV8/dgytRucaZZl90JdBK84QGpnqY/yKkmGqfK4qn6tkpXoj4f3VcP8Shda0Y07rIkzv7sjp22F0mH+cs/5qdKmiG6p8rka+L5+bVa88Uz5XF+6cknKGjZS/2EoZcFZeDnosm8w8Iit8vVsu8rgpExISsjtC5S2o8rmqfG5OKp+rZLXydeHDg7CuB6zsiH3PZewc+x4T3C4yfe81zvqHMatXDYoWLJDdkSpv6f60aURf1e01bKPq1Sj9zTeptlHlczVU+Vzl3ShRCT48BOvfh00fUKTdTP7oP4xVJ/z5ec9VOsw/zu/9auJSvlh2R6rkMqp87qvye/lc5V0pXBIG7wK3D2HPF4gngQxu+T3O5YsxZv053l90kq/aVmVYowro6akCX7lRWmfSWUGVz/2PKp+rvFsFCkHvteD6IfwzF/4ejnMZE3aPa0xru1JM23OND1edJeyZmrtUSR9VPvc/qnxuFsnT5XN1QUo4PgcO/wg2jaH3WqSxGWtP3WHqrqsUK2TI/D41qVvBPLsjVdKgyucqmZXR8rnqDD2nEgIafwbdl8LdU7C8DeJJIAPq2/C/0Q0oWMCAvktO8fvhG8SrgUiKoqASes7n1AsG/A1P78HSVhB8EfuyZuwc24jONcry20FfBi4/zcMIVQtGyThVPjdvUQk9N7BtoinspacPK9rBzcMUNjJgTm9nZvZ04t874bSf54nH9YfZHamiKNlIJfTcopSdpgRvMRvNo43n1yGE4H3Xcuwc04gShY0YvOIsP+/2ISZOleNVlPxIJfTcpEhZGLIXbBrB9lHgMQOkpHIpU7aNbsjA+tYs8bxNjz9PcPvRs+yOVlGUd0wl9NzGuAj02wI1+oLHdNgxBuJjMTbUZ0oXB/4aUIu7Yc/pMN+Trf8GZne0iqK8Qyqh50YGBaDrn9DkKzi/Ftb3hmjNAIs29qXZ+0ljHCzN+HzLBcZv8iYiKjabA1ZyAlU+V5XPVXIqIaDFJOg0H/w8YEV7iLgPQNmiJmwYXo/PWldhu3cQHX8/zvm74dkbr5LtVPlcVT5XyelqDYJ+myD0FixtDSHXAdDXE4xrWZlNH9UnLl7Sc9FJFhxRz6znV6p8roYqn6vK5+YOQeeknFlJyunlpLx9/JVNj5/HyNHr/pXWE3bJXn+ekIHhz7MpyPxJlc9V5XNzVPlcIURbYB6gDyyVUs54bXt5YBVQNLHNRCnlnrf/daOkW9mamsca1/WENV2h2yJw6AGAmYkhv/etSfOqFkzefpm2c48xrZsjnWqUzd6Y8yHPzb48CtDtgJ0S5QrT+P0qqbZR5XM18n35XCGEPrAQaA0EAmeFEDuklEkvXn0LbJZS/imEsAP2ADZvFZmSccWsNQOQNvYHt6HwJAgajAUhEELQo5YVrjbF+HSTN2M3nMf9+kOmdHGgsJEqupmXqfK5r8rv5XPrADellH4AQoiNQBcgaUKXQJHE92bAvbeOTMmcgsVhwP/gfx/Bwe/gSQC0naEZZQpYmxdi80f1+f3wDRa438TLP5y5fZxVnfV3JK0z6aygyuf+R5XPBUsgIMlyYOK6pH4APhBCBKI5Ox+rk+iUzDE0hp4roP4Y7v+7jKnrWhAR+V9ZAEN9PT57ryobR9QnPkHSa9FJ5hz0JTZejTDNi1T53P/k+/K5QoieQFsp5bDE5QFAXSnlmCRtPkvs6zchRH1gGeAgpUx4ra8RwAiA8uXL17pz545OP4zypp37x/Nt8EFKF7RgetNZuJRyeWX706hYvt9+hf+dD6JGuaLM7e2MbYm0/5xW0k+Vz1UyKyvK5wYB5ZIsWyWuS+pDYDOAlPIkYAyUeL0jKeViKaWrlNK1ZMmS6Ti0kl4XLlxgx44dLFq0iAkTJlC3bl2uXbtGpzZzWNXyD/T0jRiyfwjzz80nNv6/gUZFjA2Z09uZBf1q4v/oGe3nebLu9B3dPEKlKMo7lZ4zdAPAF2iJJpGfBfpJKa8kabMX2CSlXCmEqA4cBixlKp2rCS5059GjR3Tp0oXixYsTGxvL0aNHEUJw584dXv7ifBb7jJlnZ/L3jb+xM7djeuPpVDCr8Eo/959E8aXbBTxvPKJlNQtm9HCipKlRcodUMiAnn6ErOZvOz9CllHHAGGA/cBXN0yxXhBBThBCdE5t9DgwXQlwANgCDU0vmim6VKFECNzc3jI2NuXnzJlZWVixbtoykfwUVMizEjw1+ZE6zOQRFBtF7Z282X9/8ypl4aTNjVg2pw/ed7PC8+Yi2c49x0OdBdnwkRVEyIV0jRaWUe6SUVaSUFaWUPyeumyyl3JH43kdK2VBKWUNK6SylPJCVQSv/SUhI4K+//sLJyYmqVasyatQoypUrR58+fZJt38q6FX93/huXUi5MPTWVMUfG8OjFI+12PT3BkIa27BrbiFJFjBm+2osJbhdVPZi3pM5vlIzKzL8ZNadoLubj48OIESOIj49nyZIlODg4YGFhgYeHB3Z2dqnumyAT2HBtA7O9ZlO4QGF+bPAjzco1e6VNTFwCcw/5sujoLcqYmTCrVw3qV1RzmGbU7du3MTU1xdzcPMVnvBUlKSkloaGhREREYGtr+8q21C65qISeC0VFRTF9+nT++OMPfvjhBz7++GP09DR/bIWFhVG8ePF093Uz/CYTPSdyPfw6var04gvXLyhoWPCVNv/eCefzzd74hz5naENbvmpbFWNDfZ1+prwsNjaWwMDAtx6IouQvxsbGWFlZvfFsu0roecjRo0f56KOPsLe3Z/78+dpBFm8jJj6GBecXsPLKSqyLWDOj8QzsS9i/0uZ5TBwz9l5j9ck7VCxZiDm9nXGyKvrWx1YUJWPe9rFFJQcICwtj2LBhfPDBB8yYMYOtW7fqJJkDFNAvwGeun7H0vaW8iHvBB3s+YPHFxcQnxGvbFCxgwJQuDqz5sA7PouPp9scJNRhJUXIYldBzOCklGzduxN7eHhMTE65cuULXrl2z5Fh1ytRha+ettLJuxe/nf2fI/iEERrw661HjyiXZ/2kTOtcoy7zDN+j+xwmu34/IkngURckYdcklB3v+/Dk9e/YkICCAJUuWUK9evXdyXCklu2/v5udTPyORfFP3GzpV6PTGDb29l4KZtO0ykVFxfNq6MiMaV8BAX50jKEpWUpdccqmEhAR69uzJuXPn3lkyB01VvY4VOrK181aqFqvKpOOT+OLoFzyJfvJKu3aOZTgwvgktq1swc991eiw6yc2H6mxdUbKLOkNXUhWfEM+KKytYeH4hxU2K83Ojn6lX5tVfLlJKdl0M5rvtl3keE8+X71VlaCNb9PXUI3qKomvqDF3JNH09fYY5DmNdh3UUMizE8APD+fXsr0TH/zd5gBCCTjXKcmB8E5pVKcnPe67y/l8n8QtRhZ0U5V3KUwndw8ODjh07ZmiflStXcu9e+su3SykZN24clSpVwsnJiXPnziXbbsOGDTg6OuLk5ETbtm159EgzGjMsLIzWrVtTuXJlWrduTXh4eIb6zS525nZs6riJ3lV7s9pnNX1398U33PeVNhamxvw1oBbz+jhz82Ek7eZ5stTTT81jqijvSJ5K6JmR0YS+d+9ebty4wY0bN1i8eDEff/zxG23i4uL45JNPcHd35+LFizg5OWmL+8+YMYOWLVty48YNWrZsqZ11PD39ZjcTAxO+rfctC1suJPRFKH139WWNzxoSklRJFkLQxdmSg+Ob0LhyCX7afZVei05w86E6W1eUrJZrE/rkyZOZO3eudnnSpElcuHCByMjIZGfSnjJlCrVr18bBwYERI0YgpcTNzQ0vLy/69++Ps7MzL168SPO427dvZ+DAgQghqFevHo8fPyY4OPiVNi8nbH327BlSSp4+fUrZsmW1+w8aNAjQzL+4bdu2dPebUzSxasLfnf+mQdkGzDw7k48OfsSDZ68W8bIoYsySga7M7e2M36NntJ/vyZ8et4hTz60rStZJafborH7VqlXrjdmsM+L27duyZs2aUkop4+PjZYUKFaSbm1uKM2mHhoZq9/3ggw/kjh07pJSaGczPnj2r3fbpp5/KGjVqvPGaPn26lFIz03fS2blbtGjxyv4vbdmyRZqamsrSpUvLxo0by7i4OCmllGZmZto2CQkJ2uX09puTJCQkyM3XN8vaa2vLBusbyP239yfb7uHTKDlyjZe0nrBLdvrdU14NfvKOI1WUvAPwkink1Vx7hm5jY4O5uTnnz5/nwIED1KxZE3Nzc+1M2np6etqZtAHc3d2pW7cujo6OHDlyhCtXriTb75w5c/D29n7jNXHixHTHFhsby59//sn58+e5d+8eTk5OTJ8+/Y12InHy5txKCEGvKr3Y3HEz5U3L8/nRz/n2+LdExrx6eaWkqRF/flCLhf1cCAp/QaffjzPv0A1i4tTZuqLoUq5N6ADDhg1j5cqVrFixgqFDhwLJz6QdFRXFqFGjcHNz49KlSwwfPjzFQknjx4/H2dn5jdfLa92WlpYEBPw3xWpgYOAbQ/BfzjxesWJFhBC8//77nDhxAoBSpUppL6UEBwdjYWGR7n5zKhszG1a3X81HTh+x028nPXf25PzD82+06+BUhoOfNaWdQxnmHPKl84LjXAp8kkyPiqJkRq5O6N26dWPfvn2cPXuWNm3apNjuZfIuUaIEkZGRuLm5abclnTkc0j5D79y5M6tXr0ZKyalTpzAzM6NMmTKvHM/S0hIfHx9CQkIAOHjwoHbWkc6dO7Nq1SoAVq1aRZcuXdLdb05mqGfImJpjWNVW89kG7xvM7+d/Jzbh1TrqxQsVYH7fmiweUIuwZzF0WXic6Xuu8iImPrluFUXJAIPsDuBtFChQgObNm1O0aFH09VMu51q0aFGGDx+Og4MDpUuXpnbt2tptgwcPZuTIkZiYmHDy5ElMTExSPWb79u3Zs2cPlSpVomDBgqxYsUK7zdnZGW9vb8qWLcv3339PkyZNMDQ0xNrampUrVwIwceJE3n//fZYtW4a1tTWbN29Os9/cxNnCGbdObsw4M4PFFxdzIugE0xtPx8bM5pV279mXpm4Fc6bvucpfx/zYd+U+07s70qDiG1PRKoqSTrl6pGhCQgIuLi5s2bKFypUr6ygyRVcO+B/gx5M/EpsQy5e1v6Rn5Z7J3jM4cesRX/99iTuhz+lbpxwT21XHzMQwmR4VRcmTI0V9fHyoVKkSLVu2VMk8h3rP5j3+7vw3NUrWYMrJKYw7Mo7QF6FvtGtQsQT7PmnCiCYV2HQ2gPfmHOXAlfvZELGi5G65+gxdyR0SZALrrq5j7r9zKVygMFMbTqWJVZNk214MfMxXbhe5dj+CDo5l+L6THRZFjN9xxIqSc+XJM3Ql99ATegywG8CGjhswNzFn9OHR/HTqJ17EvTmQy8mqKDvHNuLLNlU5ePUBLWcfZf3puySo8gGKkqZcndDnz59P9erV6d+/P9u2bWPKlCmptm/WrBkZ+avA39+f9evXa5cvXbrE4MGDMxsu165do379+hgZGTFr1qwU20kpmTRpElWqVKF69erMnz9fuz6lei+rVq2icuXKVK5cWfsUTU5TpVgVNnTYwCC7QWy6von3d77PldA3xwMY6usxunkl9n3SGIeyZnzzv0v0XqxK8ypKmlIacZTVr7cdKSqllFWrVpUBAQFSSinr168vQ0JCUm3/+qjQtLi7u8sOHTq8sq5ly5byzp07GQ9WSvngwQN55swZ+c0338hff/01xXbLly+XAwYMkPHx8dr9pJRy9+7dsm3btjIhIUGePHlS1qlTR0qpGQVra2srQ0NDZVhYmLS1tZVhYWGZivFdOXnvpGyxuYV0XuUsl1xcIuPi45Jtl5CQIDefvStr/LhfVvpmt5x94LqMik2+raLkB+TFkaIjR47Ez8+Pdu3a8csvv2BkZESJEiWIiIjA1taW2FjN889Pnz59ZXnLli3UqVOHKlWq4OnpCWjOxBs3boyLiwsuLi7aQUATJ07E09MTZ2dn5syZA0CnTp3YuHFjpmK2sLCgdu3ab8zi/bo///yTyZMno6enp90PUq73sn//flq3bk3x4sUpVqwYrVu3Zt++fZmK8V2pV6Yef3f+mxblWzDv3DyG7h9KUGTQG+2EEPRyLcehz5rSwbEM8w7foN08T077vXlzVVHyu1yb0BctWkTZsmVxd3fHwsICFxcXQDNQqFmzZuzevRuAjRs30r17d20SjYuL48yZM8ydO5cff/wR0CTMgwcPcu7cOTZt2sS4ceMATWXExo0b4+3tzfjx4wFwdXXV/iJwd3dPdlRpgwYN3uqz3bp1i02bNuHq6kq7du24ceMGAEFBQZQrV07bzsrKiqCgoBTX53RmRmbMajqLnxv9zPXw6/Tc0ZOdt3ZqC6olVaKwEXP71GTV0DrExCXQe/EpvnK7QPizmGyIXFFypnQNLBJCtAXmAfrAUinljNe2zwGaJy4WBCyklEV1GGeqgoODKVmypHZ52LBhzJw5k65du7JixQqWLFmi3da9e3cAatWqpa3zEhsby5gxY/D29kZfXx9f31frfCdlYWGhLbfbvHlz7TB/XYqOjsbY2BgvLy/+/vtvhg4dqv0lktcIIehcsTMuFi58c/wbvjn+DZ6BnkyqNwkzI7M32jetUpID45sw79ANlh6/zaGrD/mmfXV6uFjm6ro4iqILaZ6hCyH0gYVAO8AO6CuEsEvaRko5XkrpLKV0Bn4H/s6CWFNkYmLySm2Whg0b4u/vj4eHB/Hx8Tg4OGi3vaz18rLOC2iG+5cqVYoLFy7g5eVFTEzKZ31RUVHa0aRpnaEvXLhQuy4jNdetrKy0v3i6devGxYsXgZTrveTmOjAvWZlasaLNCsbVHMfBOwfpsaMHZ4LPJNu2YAEDvm5fnV1jG2FjXpAvtlygz+JTqua6ku+l55JLHeCmlNJPShkDbAS6pNK+L7BBF8GlV/Xq1bl58+Yr6wYOHEi/fv0YMmRImvs/efKEMmXKoKenx5o1a4iP19QVeb3OC4Cvr6/2F8TLM/TXXy+vwY8ePVq77mU99PTo2rUr7u7uABw9epQqVaoAKdd7adOmDQcOHCA8PJzw8HAOHDiQam2bnEpfT5/hTsNZ234tJgYmDDswjN+8fiMmPvlfsNXLFMFtZAOmdXPkavBT2s07xm8HrhMVq+rCKPlUSndLX76Anmgus7xcHgAsSKGtNRAM6KfVry6ecrG2tpYhISHy2bNn0s7OTiYkJGi3BQcHS2NjYxkeHq5dl/Qpl5CQEGltbS2llNLX11c6OjpKJycn+dVXX8lChQpJKaWMiYmRzZs3l05OTnL27NlSSilHjx6traWeUcHBwdLS0lKamppKMzMzaWlpKZ880dQGb9eunQwKCpJSShkeHi7bt28vHRwcZL169aS3t7eUUvPEx6hRo2SFChWkg4PDK0/sLFu2TFasWFFWrFhRLl++PFPx5STPYp7JKSemSIeVDrLH9h7yRtiNVNs/fBolP9lwTlpP2CWbzDwij/k+fEeRKsq7RSpPuaQ5UlQI0RNoK6Uclrg8AKgrpRyTTNsJgJWUcmwKfY0ARgCUL1++1p07dzLzOyhZn3zyCZ06daJVq1YAuLm5sX37dtasWaOzY0RHR9O0aVOOHz+OgUGurmuWaxwNOMrkE5OJjInkM9fP6FutL3oi5T8sj994xHfbL3P70TM6OpXh2w52lDZTI02VvCO1kaLpSej1gR+klG0Sl78GkFK+MWODEOI8MFpKeSKtoHQ99P/BgwecPn2azp07M3bsWPbu3cuePXu0lyt04caNGwQFBdGsWTOd9amk7dGLR3x/4nuOBR6jYdmGTGk4BYuCFim2j4qN56+jfiz0uImhnmB86yoMbmCDgX6ufahLUbTeNqEbAL5ASyAIOAv0k1Jeea1dNWAfYCvT6hRVy0XJGCklW3y38OvZXzE2MOb7+t/TyrpVqvvcCX3GDzuu4H49hGqlTfmpqwOuNsXfUcSKkjXeqpaLlDIOGAPsB64Cm6WUV4QQU4QQnZM07QNsTE8yV5SMEkLwftX32dxpM2ULl2W8x3gm/zOZZ7HPUtzH2rwQywfXZtEHtXj6Ipaei07y5ZYLhEZGv8PIFeXdUdUWlVwnNj6WPy/8ybLLyyhbqCzTG0/H2cI51X2ex8Qx//BNlnr6UcjIgC/bVKVvnfLo66ln15XcRVVbVPIUQ31DxrmMY0WbFUgkg/YNYqH3wjemu0uqYAEDJrarxt5PGlO9jCnfbrtMl4XH+fdO+DuMXFGylkroSq7lUsqFLZ220LFCRxZdWMTgvYO5+/RuqvtULmXKhuH1mN+3JiER0fT48wRfbrnAI3UZRskD1CUXJU/Y57+PqSenEpsQy4TaE+heuXuapQCeRccx/8gNlh+/jbGhPp+1rsKAetbqaRglR3urp1yyikroiq7df3afb49/y+n7p2lRrgU/NPiBYsbF0tzv5sNIftx5Bc8bj6hW2pQfO9tTt4L5O4hYUTJOXUNX8oXShUqz+L3FfOH6BZ5BnnTf0Z3jQcfT3K+SRWFWD63Dog9ciIiKo/fiU4zbcJ7gJ2/OqKQoOZk6Q1fypOth15noOZGbj2/St1pfPqv1GcYGaY8YfRETz58eN1l0zA99IRjdvCLDGlfA2FD/HUStKGlTl1yUfCk6Ppq5/85l7dW1VDCrwIzGM6huXj1d+waEPefn3VfZd+U+5Yqb8G0HO96zK6VK9CrZTl1yUfIlI30jJtSZwF+t/yIyJpJ+e/qx/PJy4hPSrsZYrnhBFg2oxbphdTEx1OejNf8yYNkZbjxQ85oqOZc6Q1fyhcdRj5lyagoH7xzEtZQr0xpNo0zhMunaNzY+gbWn7jDnoC/PYuIZWN+aT1tWwaxg6lMJKkpWUJdcFAVNPZjtt7Yz/fR09IU+k+pNokOFDunePzQymlkHfNl49i5mJoZ81roK/eqUV485Ku+USuiKkkRARADfeH6Dd4g37W3bM6neJIoUKJLu/a/ce8LUXT6c8gujSqnCfNfRjsaVS6a9o6LogEroivKauIQ4ll5ayqILiyhZsCTTGk2jduna6d5fSsn+Kw+Ytucqd8Oe06q6Bd+0r06FkoWzMGpFUQldUVJ0KeQSEz0nEhARwBCHIYxxHoOhfvqvjUfHxbPiH38WHLlJdFw8g+rbMLZlZcxM1PV1JWuohK4oqXge+5xfvX7FzdeNasWrMd3pGypZ18xQHw8joph9wJdNXgEUNTHk01ZV6Fe3PIbq+rqiY+qxRUVJRUHDgnxf/3vmN59PYZ9AIjv048K0iSQ8S7nW+ussTI2Z0cOJXWMbUb1MEb7fcYU2c49xyOcBaooA5V1RCV1RgNjYWJqXb86MXkt5WMuaAqu3c7NtWx5v3YqMT/u59Zfsy5qxblhdlg7UnEANW+1F/6WnuXLvSVaFriha6pKLku/ExsYyduxY/Pz8uHfvHnfu3CEyMpLVq1czYMAAAJ6fP8/DGb/w4sIFjKpVo9TECRSqVy9jx4lPYP3pu8w95MvjF7H0dLHiizZVKVVETVqtZJ66hq4oSSQkJPD3338TExPD8uXLuXDhAhYWFly6dAk9vf/+aJVSErF3Lw9n/UbsvXsUbt4ciy+/xKiCbYaO9+RFLAvdb7Lin9sY6OkxvLEtI5pWpLCRga4/mpIPqGvoivKakJAQPvnkE2rWrImenh5r1659JZmDZh7TIu3bU2HvHkp+9hnPz5zBr3Nn7v88jbjw9M90ZGZiyDftq3Pos6a0qG7B/CM3afarB2tP3SEuPkHXH03Jx9QZupKvXLlyhREjRgCwePFifv/9dwwNDfn999/T3Dfu0SNCfl/A4y1b0CtcmBKjPqZ4v36IAgUyFMP5u+FM23OVs/7hVCxZiIntqtOquoUq/KWkizpDV/K9qKgovvvuO5o1a8YHH3yAp6cn9vb2GBgYMHXq1HT1YVCiBGV+/AHbbf/DxNGRhzN+4VanTkQcOpShJ1lqli/G5o/q89eAWkgJw1d70XvxKbwDHmfy0ymKhjpDV/I8d3d3PvroIxwdHZk/fz6WlpZv3aeUkmeenjz4ZSYxt25RsHZtLCZOwMTePkP9xMYnsPFsAPMO+fIoMoYOTmX48r2q2JQo9NYxKnmTuimq5EuhoaF8+eWXHDx4kAULFtClSxedH0PGxfF4yxZC5v9O/OPHmHXpQsnxn2JYqlSG+omMjuOvo7dY6nmb2PgE+tUtz9gWlSlpaqTzmJXcTV1yUfIVKSXr16/HwcGBwoULc+XKlSxJ5gDCwIBifftS8cB+zD8cytPdu7nVth0hvy8g4fnzdPdT2MiAz9+rytEvm9G7djnWnb5L01/dmXPQl8jouCyJXcl70nWGLoRoC8wD9IGlUsoZybR5H/gBkMAFKWW/1PpUZ+hKVrh9+zYff/wx9+7dY8mSJdStW/edHj8mMJCHv/1GxN59GFhYUHL8eMy6dEboZezcyS8kklkHrrPn0n3MCxVgXMvK9K1TngIG6hwsv3urSy5CCH3AF2gNBAJngb5SSp8kbSoDm4EWUspwIYSFlPJhav2qhK7oWnR0NDY2NnzyySd8/vnnGBpmX4Gs5+fO8WD6DKIuXcLYzg6LiRMoVKdOhvs5fzecGXuvcfp2GOWLF+Tz96rQyaksenrqiZj86m0Ten3gByllm8TlrwGklNOTtJkJ+Eopl6Y3KJXQlawQHx+Pvn7OmNBZJiTwdPceHs6eTVxwMKatW2Hx+ecUsLHJWD9S4uEbwi97r3HtfgTVSpvyVduqNK+qHnXMj972GrolEJBkOTBxXVJVgCpCiH+EEKcSL9EoyjuXU5I5gNDTw6xTRyru3UPJTz8h8p8T3OrUmQfTZxD/JP21XYQQNK9qwZ5xjZnXx5kXsfEMXelFr0UnOXM7LAs/gZLb6OqCnAFQGWgG9AWWCCGKvt5ICDFCCOElhPAKCQnR0aEVJWfTMzamxMiRVNq/j6JduxC2ejW33mtD2Oo1yNjY9PejJ+jibMmhz5ryU1cH7oY95/2/TjJkxRlV/EsB0pfQg4BySZatEtclFQjskFLGSilvo7nmXvn1jqSUi6WUrlJK15Il1ZRd2c3Dw4OOHTtmaJ+VK1dy7969dLeXUjJu3DgqVaqEk5MT586dS7ZdTEwMI0aMoEqVKlSrVo2tW7e+sn3r1q0IIUh6mW769OlUqlSJqlWrsn///gx9juxgULIkZaZOxfZ/f2NkV50H06bh16kzEUfcMzQwyVBfjw/qWXP0y+ZMbFeNc3cf02H+ccZuOI9fSGQWfgIlp0tPQj8LVBZC2AohCgB9gB2vtdmG5uwcIUQJNJdg/HQXppJTZDSh7927lxs3bnDjxg0WL17Mxx9/nGy7n3/+GQsLC3x9ffHx8aFp06babREREcybN++VJ1Z8fHzYuHEjV65cYd++fYwaNYr4DJS5zU7G1apRfvlyrBb9CUIQOGoUd4cMJerq1Qz1Y1JAn5FNK3Lsq+aMbl6RQz4PaD3nGF+5XSAwPP2PTCp5R5oJXUoZB4wB9gNXgc1SyitCiClCiM6JzfYDoUIIH8Ad+FJKGZpVQSsZN3nyZObOnatdnjRpEhcuXCAyMpKePXtSrVo1+vfvrz1TnDJlCrVr18bBwYERI0YgpcTNzQ0vLy/69++Ps7MzL168SPO427dvZ+DAgQghqFevHo8fPyY4OPiNdsuXL+frr78GQE9PjxIlSmi3fffdd0yYMAFjY+NX+u3Tpw9GRkbY2tpSqVIlzpw5k9mv550TQmDarBkVdmyn1LffEn3tGre79+DepEnEPkz1AbE3mJkY8mWbanhOaM6g+jZs875H81keTN5+mQdPo7LoEyg5UbquoUsp90gpq0gpK0opf05cN1lKuSPxvZRSfialtJNSOkopN2Zl0ErGDR06lNWrVwOa8rEbN27EysqK8+fPM3fuXHx8fPDz8+Off/4BYMyYMZw9e5bLly/z4sULdu3aRc+ePXF1dWXdunV4e3tjYmLC+PHjcXZ2fuM1Y4ZmqEJQUBDlyv13xc7KyoqgoFev2D1+/BjQJG4XFxd69erFgwcPADh37hwBAQF06NDhlX3S029uIAwNKf5Bfyru30fxwYN5smMnt9q249Gff5KQjl+YSZUobMTkTnYc/bIZvVzLsf70XZrMdGfanquEPYvJok+g5CRqlEI+YWNjg7m5OefPn+fAgQPUrFkTc3Nz6tSpg5WVFXp6ejg7O+Pv7w9o6p/UrVsXR0dHjhw5wpUrV5Ltd86cOXh7e7/xmjhxYrpji4uLIzAwkAYNGnDu3Dnq16/PF198QUJCAp999hm//fabLr6CHE3fzIxSE76i4u5dFG7UiJB587nVrj1PduxAJmSsxG4ZMxOmdXPk8OdN6eBYhiWefjT+5Qi/HbjOk+fpvwmr5D4qoecjw4YNY+XKlaxYsYKhQ4cCYGT0X60QfX194uLiiIqKYtSoUbi5uXHp0iWGDx9OVFTyf7qndYZuaWlJQMB/T70GBga+URzL3NycggUL0r17dwB69erFuXPniIiI4PLlyzRr1gwbGxtOnTpF586d8fLySle/uVGB8uWxmj8P6zWrMTA3595XE/B/vzfPMzFmw9q8ELN7O3Pg0yY0rVqS34/cpNEvR5h7yJenUSqx50lSymx51apVSyrvVnR0tKxSpYq0tbWVcXFx0t3dXXbo0EG7ffTo0XLFihUyPDxcWlhYyOfPn8uIiAhpb28vv//+eymllB07dpRHjhxJ9zF37dol27ZtKxMSEuTJkydl7dq1k23Xu3dvefjwYSmllCtWrJA9e/Z8o03Tpk3l2bNnpZRSXr58WTo5OcmoqCjp5+en/Ux5SUJ8vHy8bZv0bdJU+lStJgPGjpPRd+9muj+fe0/kiNVnpfWEXdLx+33y98O+MiIqVocRK+8C4CVTyKtqDqx8pECBAjRv3pyiRYumOgCnaNGiDB8+HAcHB0qXLk3t2rW12wYPHszIkSMxMTHh5MmTmJiYpHrM9u3bs2fPHipVqkTBggVZsWKFdpuzszPe3t4A/PLLLwwYMIBPP/2UkiVLvtIuOfb29rz//vvY2dlhYGDAwoULc9SgIl0QenqYdemC6XvvEbpiBaFLlhLp7k6xAQMoMfIj9IsUyVB/1csU4a8BrlwOesLcQ77MOuDLsuO3GdGkIgPrW1NITYmX66nyuflIQkICLi4ubNmyhcqV3xgmoORwsQ8eEjJvHk/+9z/0zcwoMXYMxXr3RhhkLhFfCHjMnEO+eFwPwbxQAYY3qcCAeiqx53SqfK6Cj48PlSpVomXLliqZ51KGpSwoO+1nbLe6YVSlCg+m/oRfl65EHj2aoYFJL9UoV5SVQ+rw96gG2JUtwoy912g8050/PW7xTJXszZXUGbqi5EJSSiLd3Xn4y0xi7tyhUIMGWEyYgHHVKpnu89zdcOYdusFR3xCKFTRkeJMKDKxvQ2F1xp6jqDN0BYD58+dTvXp1+vfvz7Zt25gyZUqq7Zs1a0ZGfun6+/uzfv167fKlS5cYPHhwZsPl2rVr1K9fHyMjI2bNmpViu/79+1O1alUcHBwYOnQosYn1UZ48eUKnTp2oUaMG9vb2r1yXX7VqFZUrV6Zy5cqsWrUq0zFmFyEEpi1aUGHnDkp98w0vrlzhdrduBH83mbhM1klyKV+MVUPr8L9RDahRrigz912n8S9HWOh+U02ykVukdLc0q1/qKZd3r2rVqjIgIEBKKWX9+vVlSEhIqu2TPlWSHq8/NSOllC1btpR37tzJeLBSygcPHsgzZ87Ib775Rv76668pttu9e7dMSEiQCQkJsk+fPvKPP/6QUkr5888/y6+++kpKKeXDhw9lsWLFZHR0tAwNDZW2trYyNDRUhoWFSVtbWxkWFpapGHOKuPBweX/aNOlj7yCv1XSRIYv+kvEvXrxVn+fvhsvBy09L6wm7ZI0f98t5h3zl4+cxOopYySxSecpFnaHnEyNHjsTPz4927drxyy+/YGRkRIkSJYiIiMDW1lZ7Vvv06dNXlrds2UKdOnWoUqUKnp6egOZMvHHjxri4uODi4sKJEycAmDhxIp6enjg7OzNnzhwAOnXqxMaNmRs4bGFhQe3atdOcqKJ9+/YIIRBCUKdOHQIDAwHNWWxERITm8kRkJMWLF8fAwID9+/fTunVrihcvTrFixWjdujX79u3LVIw5hX7RopT6+msq7NxBwQb1CZkzh1vt2/Nk1+5MXV8HcC5XlBVD6rB9dENcrYsz+6AvjWZoBiiFq5GnOZJK6PnEokWLKFu2LO7u7lhYWODi4gKAqakpzZo1Y/fu3QBs3LiR7t27a5NoXFwcZ86cYe7cufz444+AJtEePHiQc+fOsWnTJsaNGwfAjBkzaNy4Md7e3owfPx4AV1dX7S8Cd3f3ZAchNWjQQCefMTY2ljVr1tC2raYc/5gxY7h69Sply5bF0dGRefPmoaenl2fKBiTHyNaWcgsWUH7VKvSLFuXeF1/g36cPz8+fz3SfNcoVZekgV3aPa0TjKiW0A5Sm771KSES0DqNX3pa625EPBQcHk7R88bBhw5g5cyZdu3ZlxYoVLFmyRLvt5ejNWrVqacsCxMbGMmbMGLy9vdHX18fX1zfFY1lYWGirMzZv3lz73HlWGDVqFE2aNKFx48YA7N+/H2dnZ44cOcKtW7do3bq1dlteV6huHWy3bOHJ9h2EzJnDnb79KNK+HSU/+5wCVpkbUWtf1ow/+tfC90EEC91vsuSYH6tO+NO3Tnk+alKR0mbGaXeiZCl1hp4PmZiYvDKUv2HDhvj7++Ph4UF8fDwODg7abS9LA7wsCwCa+i2lSpXiwoULeHl5EROT8p/fUVFR2sFHaZ2hL1y4ULsuIyV6AX788UdCQkKYPXu2dt2KFSvo3r07QggqVaqEra0t165dy7NlA14n9PUp2r0bFffvo8To0UQcccevfXse/vYb8ZGZr5tepZQp8/rU5NBnTenoVJbVJ+/QZKY7X/99iTuhz3T4CZSMUgk9H6pevTo3b958Zd3AgQPp168fQ4YMSXP/J0+eUKZMGfT09FizZo22DrmpqSkRERGvtPX19dX+gnh5hv766+U1+NGjR2vXlS1bNt2fZ+nSpezfv58NGzagp/ffP+ny5ctz+PBhAB48eMD169epUKECbdq04cCBA4SHhxMeHs6BAwdo06ZNuo+X2+gVLEjJsWOouG8vRdq3J3TJUm6914bwjZuQcZl/eqVCycLM6lUDjy+a0bt2ObaeC6T5LA8+2Xie6/cj0u5A0b2U7pZm9Us95fLuWVtby5CQEPns2TNpZ2cnExIStNuCg4OlsbGxDA8P165L+pRLSEiItLa2llJK6evrKx0dHaWTk5P86quvZKFChaSUUsbExMjmzZtLJycnOXv2bCmlpj7Mjh07MhVvcHCwtLS0lKamptLMzExaWlrKJ0+eSCmlbNeunQwKCpJSSqmvry8rVKgga9SoIWvUqCF//PFHKaWUQUFBsnXr1tLBwUHa29vLNWvWaPtetmyZrFixoqxYsaJcvnx5puLLrZ5fuiz9+38gfapWk7c6dpQRxzx10u+DJy/ktN0+0u67vdJ6wi45bNVZef5uuE76Vv5DKk+5qIFF+dQnn3xCp06daNWqFQBubm5s376dNWvW6OwY0dHRNG3alOPHj2OQyeHpStaQUhJx6BAPf51F7N27FGrcmFJffYmRDkYRP34ew8oT/qz4x58nL2JpUNGcUc0q0bCSOUIIHUSfv6U2sEgl9HzqwYMHnD59ms6dOzN27Fj27t3Lnj17qFIl8yMNX3fjxg2CgoJo1qyZzvpUdEvGxBC2fj2P/viThMhIir7fi5Jjx2Jgbv7WfUdGx7H+9B2Wet7mYUQ0jpZmfNysIm3sS6OvpxJ7ZqmErihKquLCw3m08A/CN2xAz8SEEiM/otiAAeglqZefWdFx8fzvXBB/HfPj9qNnVChRiI+aVqBrTUuMDPJWhcx3QSV0RVHSJdrPj4e/ziLS3R1DS0ssvvgc07ZtdXKpJD5Bsv/Kff7wuMnloKeUKmLEh41s6VfXWtWLyQCV0BVFyZBnJ0/yYMYvRF+/jknNmpSaOAGTGjV00reUkn9uhvLn0Zv8czMUU2MDPqhnzZAGNlgUUc+yp0UldEVRMkzGx/Pkf//j4bx5xIc8okjHjlh8Nh7DDDxSmpYLAY9ZfMyPvZeDMdDTo1tNS4Y3qUAli8I6O0ZeoxK6oiiZFh/5jNBlSwlbrqlWWXzwYMyHD0e/cCGdHcP/0TOWHvdji1cg0XEJtKpeipFNK+BqU1xnx8grVEJXFOWtxQYH83DOHJ7u2Il+iRKU/GQcRbt3R+hw6r/QyGhWnbzD6pP+PH4ei0v5ooxoUoHWdurJmJdUQlcURWdeXLzIgxm/8OLcOYyqVKHAR19i1aGRTo/xPCaOLV6BLPH0IzD8BdbmBfmwkS09a1lRsED+voGqErqiKDolpSRi/wFu/L6WU+WGUKJoPK69a2BboyR6OjyTjotP4IDPAxYf88M74DFFCxryQV1rBjawxsI0f95AfeuELoRoC8wD9IGlUsoZr20fDPwKvKxBukBKuTS1PlVCV5TcR0rJkydP8PX1pVq1ahgJQ87/sRffJ6WICIumSEkTarQoR/UGZTA00t2lGCkl/94JZ4mnHwd8HmCop0cX57IMa1yBqqVNdXac3OCtEroQQh/wBVoDgcBZoK+U0idJm8GAq5RyTHqDUgldUXKPbt26cenSJe7du0dCQgLR0dFMmzaNr7/+GoCEBInf+RC8D93lwe2nGBU0wKGJJY7NrShk9vaDk5Lyf/SM5f/cZrNXAFGxCTSuXIKhjWxpWlm3fx3kVG+b0OsDP0gp2yQufw0gpZyepM1gVEJXlDzrypUrBAQE8NNPP3Hv3j0KFSrEhQsXXqlu+VLwrSd4H7yL34UQ9PQFVeqUxrllOcwtdfsoYvizGNafucvqk/48eBpNxZKFGNrIlu41rTApkHdHoL7tJNGWQECS5cDEda/rIYS4KIRwE0KUS2a7oii5UHR0NFu3buWDDz6gY8eOPHv2jJUrVyabzAHKVDSj3UhH+v9YD/uGZbl59gEbp55h53xvAq6GZXpKvNcVK1SA0c0r4flVC+b2dqZgAQMm/e8yDWYc5tf913jwNCrtTvKY9Jyh9wTaSimHJS4PAOomPRsXQpgDkVLKaCHER0BvKWWLZPoaAYwAKF++fK07d+7o7pMoiqJznp6ejBgxgipVqrBgwQKmT5+OEIKFCxemu4+oyFguHwviokcgL57GYG5ZGOfW5ajsWgp9A91NySCl5Kx/OMuOa66zG+gJ2juWYUhDW5zLFdXZcbJbll9yea29PhAmpTRLrV91yUVRcq7Hjx8zYcIEdu/ezfz58+nWrRtCCP744w/69u1LsWLFMtxnfGwC18/cx/tQAOHBzyhkVgCnFuWwb1wWo4KpTwSeUXdDn7PixG22eAUSGR1HzfJFGdLQlnYOpTHUz93z+rxtQjdAc1O0JZqnWM4C/aSUV5K0KSOlDE583w2YIKWsl1q/KqErSs4jpcTNzY1PP/2ULl26MH36dMzMUj03y9Qx7vqE4X3wLoHXwjEw0seuYRlqtChHkRImOj1WRFQsW/8NZOUJf/xDn1OqiBED6lnTt055zAvr9mbtu6KLxxbbA3PRPLa4XEr5sxBiCpqZM3YIIaYDnYE4IAz4WEp5LbU+VUJXlJzl7t27jB49Gj8/PxYvXkzDhg2z/JghARF4H7rLzbMPkVJS0cUC51blKWVbRKfHSUiQePg+ZMU//njeeEQBAz261CjLoAY2OFjq9hdWVlMDixRFSVF8fDy///47P/30E+PHj+fLL7+kQIEC7zSGyPAoLroHcsXzHjEv4ihTyQznVuWxcSqh80cRbzyIYOUJf/4+F8SL2HhcrYsxqIENbXPJ5RiV0BVFSdb58+cZMWIEhQsX5q+//tLpjFWZERMVx9V/grlwOICIsCjMLExwblmOqvXLYKjjRxGfvIhli1cAa07d4U7ocyxMjehf15q+dcvl6FGoKqErivIGKSWNGjXiww8/ZMiQITlqvs+E+ARunQ/B++BdHt6JwLiQIQ5NLXFsZkXBIrr96yEhQXLUN4SVJ/w56huCob7m6ZiB9W1wKV80R30voBK6oii5lJSS4JuPOX8wAP9Lj9DX16Nq3VLUaFWe4mV0V773Jb+QSNacuoObVyAR0XHYlSnCwPrWdHYum2OKgqmErihKrvf4wXO8Dwdw7WQw8bEJWDuY49y6PJZVdH8W/Sw6jm3eQaw5eYdr9yMwNTagZy0rBtSzpkLJ7J18QyV0RVHyjBcRMVw+FsQlj0BeRMRSolxhnFuVp5KrBfo6vqkppcTrTjhrTt5h7+VgYuMljSuX4IN61rSsZoFBNtxEfduh/4qiAB4eHnTs2DFD+6xcuZJ79+6lu72UknHjxlGpUiWcnJw4d+5csu2aNWtG1apVcXZ2xtnZmYcPHwJw7NgxXFxcMDAwwM3NTdve29ub+vXrY29vj5OTE5s2bcrQ58hJTEwLULuDLQN/bkCz/lWJj03g0Aof1n57kvMH7hL9Ik5nxxJCUNumOPP71uSfiS344r0q3HwYyUdr/qXRL+7MPeTL/Sc5p8RAzrgopCh51MqVK3FwcKBsOufh3Lt3Lzdu3ODGjRucPn2ajz/+mNOnTyfbdt26dbi6vnqiVr58eVauXMmsWbNeWV+wYEFWr15N5cqVuXfvHrVq1aJNmzYULVo0U58rJzAooI99Y0vsGpblzuVQvA/d5cTfNzm75zZ2jcpSo0U5TIvr7mkVC1NjxrSozMimFTl87SHrTt9l7qEb/H7kJi2rWfBBPWsaVdL9Y5YZoRK6oiRj8uTJFC9enE8//RSASZMmYWFhQWRkJD179uTy5cvUqlWLtWvXIoRgypQp7Ny5kxcvXtCgQQP++usvtm7dipeXF/3798fExISTJ09iYpL6SMjt27czcOBAhBDUq1ePx48fExwcTJkyZdIVt42NDcAbhbOSPo5YtmxZLCwsCAkJydUJ/SWhJ7BxKoGNUwke3nmK96EALh4J5OKRQCq5lMS5dXksrHU3UMlAX4829qVpY1+au6HPWX/mLlu8Ajjg84DyxQvSr255etWyypaRqOqSi6IkY+jQoaxevRqAhIQENm7ciJWVFefPn2fu3Ln4+Pjg5+fHP//8A8CYMWM4e/Ysly9f5sWLF+zatYuePXvi6urKunXr8Pb2xsTEhPHjx2svkyR9zZihmTMmKCiIcuX+K1ZqZWVFUFDQmwECQ4YMwdnZmalTp2aoguGZM2eIiYmhYsWKmf16ciwL6yK896E9A36qT40WVvhfDmXLdC/+99s5bl98hEzQ7T3D8uYFmdiuGie+bsG8Ps6UNjNmxt5r1Jt+mDHrz3Hi1iOdVZdMD3WGrijJsLGxwdzcnPPnz/PgwQNq1qyJubk5derUwcrKCgBnZ2f8/f1p1KgR7u7uzJw5k+fPnxMWFoa9vT2dOnV6o985c+boJL5169ZhaWlJREQEPXr0YM2aNQwcODDN/YKDgxkwYACrVq1KsfxtXmBa3JiGPSvj2sEWn+P3uHgkgD1/XKRoqYI4typH1bqlMdDhQCUjA326OFvSxdmSGw8iWH/mLlv/DWTXxWBsSxSiT+1y9KhlRYksPmtXCV1RUjBs2DBWrlzJ/fv3GTp0KABGRv/9QOrr6xMXF0dUVBSjRo3Cy8uLcuXK8cMPPxAVlfyNsvHjx+Pu7v7G+j59+jBx4kQsLS0JCPhv+oHAwEAsLd+cfuDlOlNTU/r168eZM2fSTOhPnz6lQ4cO/Pzzz9Srl2rtvDzDyMSAmq3L49TCilv/PsT7UAAe665zeocfDk2tcGxqiYmpbgcqVS5lyved7JnQthp7LgWz4cxdpu+9xqwD12ljX5p+dcpTr4J5llxrVwldUVLQrVs3Jk+eTGxsLOvXr8fT0zPZdi+Td4kSJYiMjMTNzY2ePXsCmoQbERGhbZvWGXrnzp1ZsGABffr04fTp05iZmb1x/TwuLo7Hjx9TokQJYmNj2bVrF61atUq135iYGLp168bAgQO1seUn+vp6VKlTmsq1S3HP9zHnD93l7K7bnNt/h6r1NDMqFSut24FKxob6dHexoruLFTceRLDhTABbz2nO2r9uV42Pmur+kpdK6IqSggIFCtC8eXOKFi2Kvn7Kf54XLVqU4cOH4+DgQOnSpaldu7Z22+DBgxk5cmS6b4q2b9+ePXv2UKlSJQoWLMiKFSu025ydnfH29iY6Opo2bdoQGxtLfHw8rVq1Yvjw4QCcPXuWbt26ER4ezs6dO/n++++5cuUKmzdv5tixY4SGhrJy5UpA8wSOs7Nz5r+gXEgIgWXVYlhWLUZY8DMuHA7g+sn7+Hjew8apBM6tylG2su4HKlUuZcrkTnZ81bYqey8HU9fWXKf9v6QGFilKChISEnBxcWHLli1Urlw5u8NRssjzpzFcOhrIZY8gop7FYmFtinOr8lR0KYleDqy+qAYWKUoG+fj4UKlSJVq2bKmSeR5XsEgB6naqwMDpDWjaryrRL+I4sOwKa787hfehu8RE6W6gUlZTZ+iKoihJyATJ7YuP8D50l+CbTyhgrI9dY0ucmlvpdKBSZqkzdEXJhPnz51O9enX69+/Ptm3bmDJlSqrtmzVrRkZOUvz9/Vm/fr12+dKlSwwePDiz4XLt2jXq16+PkZHRGyNFk7p9+zZ169alUqVK9O7dm5iYGABmz56NnZ0dTk5OtGzZktcncX/69ClWVlaMGTMmuW7zDKEnqOBcku5f1KLnBFfK25tz4dBd1n57koPLrxByNyLtTrKJSuiKkoI//viDgwcPsm7dOmbOnMmoUaN02v/rCd3R0ZHAwEDu3r2bqf6KFy/O/Pnz+eKLL1JtN2HCBMaPH8/NmzcpVqwYy5YtA6BmzZp4eXlx8eJFevbsyVdfffXKft999x1NmjTJVGy5VSnbIrQZ7sAHU+vj2MyK2xcesXnaWbbNOY//Jd0PVHpbKqErSjJGjhyJn58f7dq145dffsHIyIgSJUoQERGBra0tsbGxgOasNenyli1bqFOnDlWqVNE+5ujv70/jxo1xcXHBxcWFEydOADBx4kQ8PT1xdnbWPs7YqVMnNm7cmKmYLSwsqF27NoaGhim2kVJy5MgR7aOLgwYNYtu2bQA0b96cggULAlCvXj0CAwO1+/377788ePCA9957L1Ox5XZFSpjQ6P3KDJregPrdKvL4wXN2L7zIhqln8PnnHnGx8dkdIqASuqIka9GiRZQtWxZ3d3csLCxwcXEBNM+VN2vWjN27dwOwceNGunfvrk2icXFxnDlzhrlz5/Ljjz8CmkR78OBBzp07x6ZNmxg3bhwAM2bMoHHjxnh7ezN+/HgAXF1dtb8I3N3dky0T0KBBg0x/rtDQUIoWLYqBgeaJ5ZRKCyxbtox27doBmqd9Pv/881Qv4+QXRgUNcWljzYCf6tNqcHX09AXua66x+psTeO25TVRkbLbGp55DV5Q0BAcHU7JkSe3ysGHDmDlzJl27dmXFihUsWbJEu6179+4A1KpVC39/fwBiY2MZM2YM3t7e6Ovr4+vrm+KxLCwstOV2mzdvjre3t+4/UBrWrl2Ll5cXR48eBTSXntq3b68teaCAvoEeVeuVoUrd0gReD8f74F1O77jNv3vvUK1+GWq0LEfRUgXfeVwqoStKGkxMTHjy5Il2uWHDhvj7++Ph4UF8fDwODg7abS9LA7wsCwCa0aGlSpXiwoULJCQkYGyc8pMSUVFR2sFH7u7u2jP3pAoWLMiJEydYuHCh9pfJnj170lWi19zcnMePHxMXF4eBgcEbpQUOHTrEzz//zNGjR7Wf5eTJk3h6evLHH38QGRlJTEwMhQsX1hYUy8+EEJSrVpxy1YoTGhSJ9+EAfE7c47JnELZOJXBuXZ4yFc3e2bykKqErShqqV6/O2rVrX1k3cOBA+vXrx3fffZfm/k+ePMHKygo9PT1WrVpFfLzmeuvrZQEAfH19tb8g0jpDHz16NKNHj87QZxFC0Lx5c9zc3OjTpw+rVq2iS5cuAJw/f56PPvqIffv2YWFhod1n3bp12vcrV67Ey8tLJfNkmFsWpuXA6tTrUoFLHoFcPhrE7QuPKGVbBOdW5angXCLLByqpa+iKkoYmTZpw/vz5V8qg9u/fn/DwcPr27Zvm/qNGjWLVqlXUqFGDa9euUaiQpmaIk5MT+vr61KhRQ3tT1N3dnQ4dOmQqzvv372NlZcXs2bP56aefsLKy4unTp4CmpMDLSzm//PILs2fPplKlSoSGhvLhhx8C8OWXXxIZGUmvXr1wdnamc+fOmYojvytkZkS9LhUZNL0hTfpU4UVkLPuXXGbt5FNcOBKQpQOV0jWwSAjRFpgH6ANLpZTJ/noWQvQA3IDaUspUH8hVA4uU3OSTTz6hU6dO2iJYbm5ubN++nTVr1ujsGNHR0TRt2pTjx49rb1oquV9CguT2hRC8DwZw3+8JRgUNaNK3ClVql85Uf6kNLErzX40QQh9YCLQGAoGzQogdUkqf19qZAp8Ayc+XpSi52DfffKOdCm7s2LHs3buXPXv26PQYd+/eZcaMGSqZ5zF6eoKKNS2oWNOC+35P8D54lyIlUi/SlllpnqELIeoDP0gp2yQufw0gpZz+Wru5wEHgS+ALdYauKIqie2879N8SCEiyHJi4LukBXIByUsrdmY5SURRFeStvfVNUCKEHzAY+T0fbEUIILyGEV0hIyNseWlEURUkiPQk9CCiXZNkqcd1LpoAD4CGE8AfqATuEEG/8SSClXCyldJVSuiYdqKEoiqK8vfQk9LNAZSGErRCiANAH2PFyo5TyiZSyhJTSRkppA5wCOqd1DV1RFEXRrTQTupQyDhgD7AeuApullFeEEFOEEOpBVUVRlBwiXc9HSSn3AHteWzc5hbbN3j4sRVEUJaPUSFFFUZQ8QiV0RVGUPCLb5hQVQoQAd9JsmLwSwCMdhqMrKq6MUXFlXE6NTcWVMW8Tl7WUMtnHBLMtob8NIYRXSiOlspOKK2NUXBmXU2NTcWVMVsWlLrkoiqLkESqhK4qi5BG5NaEvzu4AUqDiyhgVV8bl1NhUXBmTJXHlymvoiqIoypty6xm6oiiK8pocl9CFEG2FENeFEDeFEBOT2W4khNiUuP20EMImybavE9dfF0K0yQlxCSFshBAvhBDeia9F7ziuJkKIc0KIOCFEz9e2DRJC3Eh8DcpBccUn+b52vL5vFsf1mRDCRwhxUQhxWAhhnWRbdn5fqcWVnd/XSCHEpcRjHxdC2CXZlp0/j8nGld0/j0na9RBCyKRFDHXyfUkpc8wLzRR3t4AKQAHgAmD3WptRwKLE932ATYnv7RLbGwG2if3o54C4bIDL2fh92QBOwGqgZ5L1xQG/xP8WS3xfLLvjStwWmY3fV3OgYOL7j5P8f8zu7yvZuHLA91UkyfvOwL7E99n985hSXNn685jYzhQ4hqaQoasuv6+cdoZeB7gppfSTUsYAG4Eur7XpAqxKfO8GtBRCiMT1G6WU0VLK28DNxP6yO66slGZcUkp/KeVFIOG1fdsAB6WUYVLKcDSzTbXNAXFlpfTE5S6lfJ64eApNuWjI/u8rpbiyUnrieppksRDw8qZctv48phJXVkpPngCYCvwCRCVZp5PvK6cl9DRnR0raRmoqQT4BzNO5b3bEBWArhDgvhDgqhGiso5jSG1dW7JvVfRsLzUQop4QQXXUUU2bi+hDYm8l931VckM3flxBitBDiFjATGJeRfbMhLsjGn0eR8uxuOvm+1Gy0WS8YKC+lDBVC1AK2CSHsXzuDUF5lLaUMEkJUAI4IIS5JKW+9ywCEEB8ArkDTd3nctKQQV7Z+X1LKhcBCIUQ/4FtAp/cXMiuFuLLt51H8N7vb4Kw6Rk47Q09rdqRX2gghDAAzIDSd+77zuBL/hAoFkFL+i+baWJV3GFdW7JulfUspgxL/6wd4ADXfZVxCiFbAJDQTtURnZN9siCvbv68kNgJdM7nvO4krm38eU5vdTTffV1bcHHiLmwoGaG422fLfTQX719qM5tWbj5sT39vz6k0FP3R3E+Zt4ir5Mg40N0uCgOLvKq4kbVfy5k3R22hu8BVLfJ8T4ioGGCW+LwHcIJkbS1n4/7Emmh/yyq+tz9bvK5W4svv7qpzkfSfAK/F9dv88phRXjvh5TGzvwX83RXXyfb31h9D1C2gP+Cb+452UuG4KmrMSAGNgC5qbBmeACkn2nZS433WgXU6IC+gBXAG8gXNAp3ccV2001+OeoflL5kqSfYcmxnsTGJIT4gIaAJcS/3FfAj58x3EdAh4k/v/yBnbkkO8r2bhywPc1L8m/b3eSJLBs/nlMNq7s/nl8ra0HiQldV9+XGimqKIqSR+S0a+iKoihKJqmEriiKkkeohK4oipJHqISuKIqSR6iEriiKkkeohK4oipJHqISuKIqSR6iEriiKkkf8Hzi4fU2UiZHnAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Bonus cell. This is the cell Stephanie used to make a figure\n", "# in the PDF lecture notes. It illustrates the first 5 steps\n", "# taken to solve the test equation with a rate constant of 2\n", "# and step size of 0.1.\n", "params = np.array((2,))\n", "\n", "y0 = 1\n", "tru_ts = np.arange(0,0.41,0.01)\n", "tru_y1 = y0*np.exp( -params[0]*tru_ts)\n", "\n", "plt.plot( tru_ts, tru_y1 )\n", "\n", "# Illustrate forward Euler for a few steps\n", "# of size 0.1\n", "ts = np.arange(0,0.5,0.1)\n", "y = np.zeros((1,ts.size))\n", "y[0,0] = y0\n", "for i in range(1,ts.size):\n", " h = t[i]-t[i-1];\n", " yprime = degradation_model(t[i-1],y[0,i-1],params);\n", " y[0,i] = y[0,i-1] + h*yprime;\n", " plt.plot( (ts[i-1],ts[i]), (y[0,i-1], y[0,i]) )\n", " if i == 1:\n", " offsets = (0.01,0.1)\n", " else:\n", " offsets = (-0.1,-0.1)\n", " plt.annotate( \"yhat=%0.3f\\nf(yhat)=%0.3f\"% (y[0,i-1],yprime), (ts[i-1],y[0,i-1]), (ts[i-1]+offsets[0],y[0,i-1]+offsets[1]), \n", " xycoords='data', arrowprops={'arrowstyle':'->'})\n", "plt.legend(['True Solution (y)','Approximate Solution','Approximate Solution','Approximate Solution','Approximate Solution']);\n", " " ] } ], "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.8.5" } }, "nbformat": 4, "nbformat_minor": 4 }