Real-time collaboration for Jupyter Notebooks, Linux Terminals, LaTeX, VS Code, R IDE, and more,
all in one place. Commercial Alternative to JupyterHub.
Real-time collaboration for Jupyter Notebooks, Linux Terminals, LaTeX, VS Code, R IDE, and more,
all in one place. Commercial Alternative to JupyterHub.
| Download
"Guiding Future STEM Leaders through Innovative Research Training" ~ thinkingbeyond.education
Project: stephanie's main branch
Path: ThinkingBeyond Activities / BeyondAI-2024-Mentee-Projects / khurshid / ODE_methods (1).ipynb~1
Views: 1166Image: ubuntu2204
{
"nbformat": 4,
"nbformat_minor": 0,
"metadata": {
"colab": {
"provenance": [],
"toc_visible": true,
"gpuType": "T4"
},
"kernelspec": {
"name": "python3",
"display_name": "Python 3"
},
"language_info": {
"name": "python"
},
"accelerator": "GPU"
},
"cells": [
{
"cell_type": "code",
"source": [
"!pip install torchdiffeq\n",
"!pip install git+https://github.com/rtqichen/torchdiffeq\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n"
],
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "IFFc8CmQe8CX",
"outputId": "540b0c58-c082-4ce2-c4e4-7c70b163ad54"
},
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"name": "stdout",
"text": [
"Collecting torchdiffeq\n",
" Downloading torchdiffeq-0.2.5-py3-none-any.whl.metadata (440 bytes)\n",
"Requirement already satisfied: torch>=1.5.0 in /usr/local/lib/python3.10/dist-packages (from torchdiffeq) (2.5.1+cu121)\n",
"Requirement already satisfied: scipy>=1.4.0 in /usr/local/lib/python3.10/dist-packages (from torchdiffeq) (1.13.1)\n",
"Requirement already satisfied: numpy<2.3,>=1.22.4 in /usr/local/lib/python3.10/dist-packages (from scipy>=1.4.0->torchdiffeq) (1.26.4)\n",
"Requirement already satisfied: filelock in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq) (3.16.1)\n",
"Requirement already satisfied: typing-extensions>=4.8.0 in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq) (4.12.2)\n",
"Requirement already satisfied: networkx in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq) (3.4.2)\n",
"Requirement already satisfied: jinja2 in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq) (3.1.4)\n",
"Requirement already satisfied: fsspec in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq) (2024.10.0)\n",
"Requirement already satisfied: sympy==1.13.1 in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq) (1.13.1)\n",
"Requirement already satisfied: mpmath<1.4,>=1.1.0 in /usr/local/lib/python3.10/dist-packages (from sympy==1.13.1->torch>=1.5.0->torchdiffeq) (1.3.0)\n",
"Requirement already satisfied: MarkupSafe>=2.0 in /usr/local/lib/python3.10/dist-packages (from jinja2->torch>=1.5.0->torchdiffeq) (3.0.2)\n",
"Downloading torchdiffeq-0.2.5-py3-none-any.whl (32 kB)\n",
"Installing collected packages: torchdiffeq\n",
"Successfully installed torchdiffeq-0.2.5\n",
"Collecting git+https://github.com/rtqichen/torchdiffeq\n",
" Cloning https://github.com/rtqichen/torchdiffeq to /tmp/pip-req-build-haygst08\n",
" Running command git clone --filter=blob:none --quiet https://github.com/rtqichen/torchdiffeq /tmp/pip-req-build-haygst08\n",
" Resolved https://github.com/rtqichen/torchdiffeq to commit a88aac53cae738addee44251288ce5be9a018af3\n",
" Preparing metadata (setup.py) ... \u001b[?25l\u001b[?25hdone\n",
"Requirement already satisfied: torch>=1.5.0 in /usr/local/lib/python3.10/dist-packages (from torchdiffeq==0.2.5) (2.5.1+cu121)\n",
"Requirement already satisfied: scipy>=1.4.0 in /usr/local/lib/python3.10/dist-packages (from torchdiffeq==0.2.5) (1.13.1)\n",
"Requirement already satisfied: numpy<2.3,>=1.22.4 in /usr/local/lib/python3.10/dist-packages (from scipy>=1.4.0->torchdiffeq==0.2.5) (1.26.4)\n",
"Requirement already satisfied: filelock in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq==0.2.5) (3.16.1)\n",
"Requirement already satisfied: typing-extensions>=4.8.0 in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq==0.2.5) (4.12.2)\n",
"Requirement already satisfied: networkx in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq==0.2.5) (3.4.2)\n",
"Requirement already satisfied: jinja2 in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq==0.2.5) (3.1.4)\n",
"Requirement already satisfied: fsspec in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq==0.2.5) (2024.10.0)\n",
"Requirement already satisfied: sympy==1.13.1 in /usr/local/lib/python3.10/dist-packages (from torch>=1.5.0->torchdiffeq==0.2.5) (1.13.1)\n",
"Requirement already satisfied: mpmath<1.4,>=1.1.0 in /usr/local/lib/python3.10/dist-packages (from sympy==1.13.1->torch>=1.5.0->torchdiffeq==0.2.5) (1.3.0)\n",
"Requirement already satisfied: MarkupSafe>=2.0 in /usr/local/lib/python3.10/dist-packages (from jinja2->torch>=1.5.0->torchdiffeq==0.2.5) (3.0.2)\n"
]
}
]
},
{
"cell_type": "markdown",
"source": [],
"metadata": {
"id": "bVoY4cfZ2Rip"
}
},
{
"cell_type": "markdown",
"source": [
"#### **Implementations of Euler's and RK4 methods**\n",
"Euler's method is one of the oldest and simplest algorithms.\n",
"The core idea is approximating the solution function step by step using tangents lines.\n",
"\n",
"$\\frac{dy}{dt} = f(t, y)$ and $y(t_0) = y_0$, then $y_t = y_0 + f(t_0, y_0)(t-t_0)$\n",
"\n",
"\n",
"----------------------------------------------------------\n",
"\n",
"#### **Runge Kutta 4 method**\n",
"The Runge-Kutta 4th Order method is a numerical technique for solving ordinary differential equations (ODEs). It provides a balance between accuracy and computational efficiency by approximating the solution at each time step using four intermediate evaluations of the derivative.\n",
"The 4th order Runge-Kutta method approximates the curve of the function f between two points by a 4th degree polynomial. It requires 4 slope estimates at each step.\n",
"Geometrically, we can visualize the curve between two points as part of a 4th degree polynomial. The 4th order Runge-Kutta method finds 4 tangents to this curve at evenly spaced points, and uses the average of their slopes as an approximation of the entire curve between the two points.\n",
"\n",
"\n",
"\n",
"\n"
],
"metadata": {
"id": "kuxB2aig2jfF"
}
},
{
"cell_type": "markdown",
"source": [
"Here we are going to write function for each method and then utilize it initial value problem solver function. Next, set a table and graph for showing how far the hidden point of t from the exact point:"
],
"metadata": {
"id": "Sh8Hl3oY4Mqj"
}
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/",
"height": 599
},
"id": "DDj6GfzJdOBG",
"outputId": "4049e414-67a7-4954-c5c7-7e0b55e0544e"
},
"outputs": [
{
"output_type": "stream",
"name": "stdout",
"text": [
"| t | y | RK4 | Exact |\n",
"|------|-----------|-----------|-----------|\n",
"| 0.00 | 1.000000 | 1.000000 | 1.0 |\n",
"| 0.20 | 1.000000 | 1.020201 | 1.0202 |\n",
"| 0.40 | 1.040000 | 1.083287 | 1.08329 |\n",
"| 0.60 | 1.123200 | 1.197217 | 1.19722 |\n",
"| 0.80 | 1.257984 | 1.377126 | 1.37713 |\n",
"| 1.00 | 1.459261 | 1.648717 | 1.64872 |\n"
]
},
{
"output_type": "display_data",
"data": {
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
],
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAjsAAAG0CAYAAADU2ObLAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABxHUlEQVR4nO3deZxN9R/H8ded1ZgYO4OxRGRPhEKWlCXKljVrdrKVvbKVrY1ElrKLiiFJIrKElGXIj6wj2yDbjFnMcuf8/jjNMNYZZubMvfN+Ph734d5zz733c89P7vv3XW2GYRiIiIiIOCkXqwsQERERSUkKOyIiIuLUFHZERETEqSnsiIiIiFNT2BERERGnprAjIiIiTk1hR0RERJyam9UFWC02NpZz586RKVMmbDab1eWIiIhIIhiGwfXr18mbNy8uLvdvu0n3YefcuXP4+flZXYaIiIg8hNOnT5M/f/77npPuw06mTJkA82JlzpzZ4mpEREQkMUJCQvDz84v/Hb+fdB924rquMmfOrLAjIiLiYBIzBEUDlEVERMSpKeyIiIiIU1PYEREREaeW7sfsJJbdbic6OtrqMsTBubu74+rqanUZIiLpisLOAxiGwfnz57l27ZrVpYiTyJIlC3ny5NG6TiIiqURh5wHigk6uXLnImDGjfqDkoRmGQXh4OBcvXgTA19fX4opERNIHhZ37sNvt8UEne/bsVpcjTsDLywuAixcvkitXLnVpiYikAg1Qvo+4MToZM2a0uBJxJnF/nzQGTEQkdSjsJIK6riQ56e+TiEjqUjeWiIiIpAy7HbZuhaAg8PWF6tXBgu57hR0RERFJfv7+0K8fnDlz81j+/DBlCjRtmqqlqBtLnMK8efPIkiXLI7/Ppk2bsNlsWmpARORR+PtD8+YJgw7A2bPmcX//VC1HYcdJdezYEZvNdsetXr16qVbDqFGjeOqppx54Xnh4OMOGDaNIkSJkyJCBnDlzUqNGDb7//vsUra9mzZr0798/wbHnnnuOoKAgfHx8UvSzRUSclt1utugYxp3PxR3r3988L5WoG8uJ1atXj7lz5yY45unpaVE199ajRw927tzJ1KlTKVmyJJcvX2b79u1cvnw51Wvx8PAgT548qf65IiJOY+vWO1t0bmUYcPq0eV7NmqlSklp2ksgwDMLCwiy5GXdLyffh6elJnjx5EtyyZs0KmN01Hh4ebN26Nf78SZMmkStXLi5cuADA2rVrqVatGlmyZCF79uw0bNiQ48ePJ/iMM2fO0Lp1a7Jly4a3tzcVK1Zk586dzJs3j9GjR7Nv3774VqV58+bdtc5Vq1YxfPhwGjRoQKFChahQoQJvvvkmnTt3jj/n6tWrtG/fnqxZs5IxY0bq16/P0aNH7/ndO3bsSOPGjRMc69+/PzX/+w+rY8eObN68mSlTpsTXd/Lkybt2Yy1fvpxSpUrh6elJoUKF+PjjjxO8b6FChRg3bhydO3cmU6ZMFChQgFmzZt2zNhERpxYUlLznJQO17CRReHg4jz32mCWfHRoaire3d7K8V1wXTrt27di3bx8nTpzg3Xff5bvvviN37twAhIWFMXDgQMqWLUtoaCjvvfceTZo0ISAgABcXF0JDQ6lRowb58uVj1apV5MmThz179hAbG0vLli05cOAAa9eu5ZdffgG4Z9dQnjx5WLNmDU2bNiVTpkx3Padjx44cPXqUVatWkTlzZoYMGUKDBg04ePAg7u7uSf7+U6ZM4ciRI5QuXZoxY8YAkDNnTk6ePJngvN27d9OiRQtGjRpFy5Yt2b59O7169SJ79ux07Ngx/ryPP/6YsWPHMnz4cJYtW0bPnj2pUaMGxYsXT3JtIiIOLbGrw6fmKvJGOhccHGwARnBw8B3PRUREGAcPHjQiIiLij4WGhhqAJbfQ0NBEf68OHToYrq6uhre3d4LbBx98EH9OZGSk8dRTTxktWrQwSpYsaXTt2vW+7/nvv/8agPHXX38ZhmEYM2fONDJlymRcvnz5ruePHDnSKFeu3ANr3bx5s5E/f37D3d3dqFixotG/f3/jt99+i3/+yJEjBmBs27Yt/tilS5cMLy8v49tvvzUMwzDmzp1r+Pj4JPj+r776aoLP6devn1GjRo34xzVq1DD69euX4Jxff/3VAIyrV68ahmEYbdq0MV588cUE5wwaNMgoWbJk/OOCBQsar7/+evzj2NhYI1euXMYXX3xx1+97t79XIiJOIybGMHx9jVizw+rOm81mGH5+5nmP4H6/37dTy04SZcyYkdDQUMs+Oylq1arFF198keBYtmzZ4u97eHiwePFiypYtS8GCBfn0008TnHv06FHee+89du7cyaVLl4iNjQXg1KlTlC5dmoCAAMqXL5/gPR/G888/z4kTJ/j999/Zvn07GzZsYMqUKYwePZp3332XQ4cO4ebmRuXKleNfkz17dooXL86hQ4ce6bMf5NChQ7z66qsJjlWtWpXJkydjt9vjt3soW7Zs/PM2m408efLE74ElIpKuuLoSnDUrPkFBGECCZVTjFlWdPDlV19tR2Ekim82WbF1JKc3b25uiRYve95zt27cDcOXKFa5cuZLguzVq1IiCBQsye/Zs8ubNS2xsLKVLlyYqKgq4uc9TcnB3d6d69epUr16dIUOG8P777zNmzBiGDBnyUO/n4uJyxxinlNye4fauNJvNFh8ORUTSE2P5cnwOHiQGCPfyInNExM0n8+c3g47W2ZHUcvz4cQYMGMDs2bOpXLkyHTp0iP+Bvnz5MocPH+add97hhRdeoESJEly9ejXB68uWLUtAQABXrly56/t7eHhgf8iphSVLliQmJoYbN25QokQJYmJi2LlzZ/zzcfWVLFnyrq/PmTMnQbcNfgsICEhyfSVKlGDbtm0Jjm3bto1ixYppE08Rkdtdu8aNrl0B+MTNjfCjR+HXX+Hrr80/AwNTPeiAwo5Ti4yM5Pz58wluly5dAswd3V9//XXq1q1Lp06dmDt3Lvv374+faZQ1a1ayZ8/OrFmzOHbsGBs3bmTgwIEJ3r9169bkyZOHxo0bs23bNk6cOMHy5cvZsWMHYM5SCgwMJCAggEuXLhEZGXnXOmvWrMnMmTPZvXs3J0+eZM2aNQwfPpxatWqROXNmnnjiCV599VW6du3Kb7/9xr59+3j99dfJly/fHV1McWrXrs2uXbtYsGABR48eZeTIkRw4cCDBOYUKFWLnzp2cPHkyQTfdrd566y02bNjA2LFjOXLkCPPnz+fzzz/n7bffTtr/GCIi6UDsoEF4Xb3KYSDi7bfJky+fOb28dWvzT6v+T+IjjQ5yAkkdoOwoOnTocNdBzsWLFzcMwzBGjx5t+Pr6GpcuXYp/zfLlyw0PDw8jICDAMAzDWL9+vVGiRAnD09PTKFu2rLFp0yYDMFasWBH/mpMnTxrNmjUzMmfObGTMmNGoWLGisXPnTsMwDOPGjRtGs2bNjCxZshiAMXfu3LvWOm7cOOPZZ581smXLZmTIkMF4/PHHjb59+yao7cqVK0a7du0MHx8fw8vLy6hbt65x5MiR+OdvH6BsGIbx3nvvGblz5zZ8fHyMAQMGGH369EkwQPnw4cNGlSpVDC8vLwMwAgMD7xigbBiGsWzZMqNkyZKGu7u7UaBAAePDDz9M8DkFCxY0Pv300wTHypUrZ4wcOfKu39eR/16JiNzTpk3xg5AbZs6cqIHDjyIpA5RthpHExVucTEhICD4+PgQHB5M5c+YEz924cYPAwEAKFy5MhgwZLKpQnI3+XomI07lxg9gyZXA5dowZQMQnnzBgwIAU/cj7/X7fTgOURURE5NGMHYvLsWOcAz7Pn59dPXtaXVECGrMjIiIiD2/fPoxJkwDoBQx6//0012qtsCMiIiIPJyYGunTBFhPDMuB46dK8/vrrVld1B3VjiYiIyMP57DPYtYtrwJvA7PHj0+SyHGrZERERkaQ7cQLeeQeAt4Fizz/Pyy+/bG1N96CWHREREUkaw4AePSAigl+Br4A/PvoIm832oFdaQi07IiIikjQLF8L69US6uNANaNWqFc8884zVVd2TWnZEREQk8S5ehP/W0HkvNpZTHh6sGzfO4qLuTy07clebNm3CZrNx7do1q0sREZG0pH9/uHKFQxky8Anw5ptvUrhwYaurui+FnVRit8OmTbBkifnnQ+6PmWgdO3bEZrPdcatXr17KfrCIiDivH3+EJUuItdl4/cYNMmXNyogRI6yu6oHUjZUK/P2hXz84c+bmsfz5YcqUlN38tV69esydOzfBMU9Pz5T7wNtERUXh4eGRap8nIiIp6Pp1c1AyMDNjRvaEhfHJu++SNWtWiwt7MLXspDB/f2jePGHQATh71jzu759yn+3p6UmePHkS3LJmzcrJkyex2WwEBATEn3vt2jVsNhubNm265/v99ttvVK9eHS8vL/z8/Ojbty9hYWHxzxcqVIixY8fSvn17MmfOTLdu3VLuy4mISOoaPhzOnOFK1qy8FRbG448/Tq9evayuKlEUdpLIMCAsLHG3kBDo29d8zd3eB8wWn5CQxL2flVu2Hj9+nHr16tGsWTP279/PN998w2+//UafPn0SnPfRRx9Rrlw59u7dy7vvvmtRtSIikqy2b4dp0wBoHxFBBDB+/PhU7S14FGkq7GzZsoVGjRqRN29ebDYbK1eufOBrIiMjGTFiBAULFsTT05NChQoxZ86cFKsxPBweeyxxNx8fswXnXgzDbPHx8Unc+4WHJ63W1atX89hjjyW4jXvIEfPjx4+nbdu29O/fnyeeeILnnnuOzz77jAULFnDjxo3482rXrs1bb71FkSJFKFKkyEN9loiIpCGRkdClCxgG24sV48cbN6hcuTKvvfaa1ZUlWpoasxMWFka5cuXo3LkzTRM5mKVFixZcuHCBr776iqJFixIUFERsbGwKV+oYatWqxRdffJHgWLZs2QgJCUnye+3bt4/9+/ezePHi+GOGYRAbG0tgYCAlSpQAoGLFio9WtIiIpC0TJsChQ8Rky8YrR48C8PHHH6fZBQTvJk2Fnfr161O/fv1En7927Vo2b97MiRMnyJYtG2COG7mfyMhIIiMj4x8n9Yc/Y0YIDU3cuVu2QIMGDz5vzRp4/vnEfXZSeHt7U7Ro0TuOh/73BYxb+sWio6Pv+16hoaF0796dvn373vFcgQIFEnymiIg4iYMH4YMPAPiwYEEuX7lCkyZNqFq1qsWFJU2aCjtJtWrVKipWrMikSZNYuHAh3t7evPLKK4wdOxYvL6+7vmb8+PGMHj36oT/TZoPE/p6/9JI56+rs2buPt7HZzOdfeglSc9+0nDlzAhAUFET58uUBEgxWvpunn36agwcP3jU8iYiIE4qNNbuvoqO5WLkyw3fuxM3NjQkTJlhdWZI5dNg5ceIEv/32GxkyZGDFihVcunSJXr16cfny5TumXMcZNmwYAwcOjH8cEhKCn59fitTn6mpOL2/e3Aw2twaeuNa/yZNTLuhERkZy/vz5BMfc3NzIkSMHVapUYcKECRQuXJiLFy/yzn+bud3LkCFDqFKlCn369KFLly54e3tz8OBB1q9fz+eff54yX0BERKzzxRewYwdGpky0unIFMBcQLFasmMWFJV2aGqCcVLGxsdhsNhYvXkylSpVo0KABn3zyCfPnzyciIuKur/H09CRz5swJbimpaVNYtgzy5Ut4PH9+83hKrrOzdu1afH19E9yqVasGwJw5c4iJiaFChQr079+f999//77vVbZsWTZv3syRI0eoXr065cuX57333iNv3rwp9wVERMQap07B0KEAbGnQgF+PHiV79uwOO8vWoVt2fH19yZcvHz4+PvHHSpQogWEYnDlzhieeeMLC6m5q2hRefRW2boWgIPD1herVU7brat68ecybN++ez5coUYLt27cnOHbrGJ6aNWsmeAzwzDPPsG7dunu+58mTJx+qVhERSUMMA3r2hNBQYipVovn69QCMGTPGIRYQvBuHDjtVq1blu+++IzQ0lMceewyAI0eO4OLiQv78+S2uLiFXV6hZ0+oqREREHuCbb8yZMx4efFKiBJf++IOSJUs69EKxaaobKzQ0lICAgPjBsoGBgQQEBHDq1CnAHG/Tvn37+PPbtGlD9uzZ6dSpEwcPHmTLli0MGjSIzp0733OAsoiIiNzD5cvmarjApe7dGfHfciOffPIJbm6O2z6SpsLOrl27KF++fPwMoYEDB8aPDQFz9lBc8AF47LHHWL9+PdeuXaNixYq0bduWRo0a8dlnn1lSv4iIiEN76y34918oVYrugYHExMTQoEED6tata3Vlj8Rm3D4wI50JCQnBx8eH4ODgOwYr37hxg8DAQAoXLkyGDBksqlCcjf5eiUiatH69uRaKzcYfkydTuV8/XF1dOXDgAE8++aTV1d3hfr/ft0tTLTsiIiJigbAw6N4dgNjevXlj9mwAevfunSaDTlIp7IiIiKR3I0dCYCAUKMC8okU5cOAAWbNmZeTIkVZXliwcd7SRiIiIPLo//4RPPwUg7OOPGdqrFwCjRo2K34rJ0allR0REJL2Kjja3hIiNhTZtGLVzJ//++y/FixenZ8+eVleXbNSyIyIikl599BHs3w/Zs3O0d28m16gBmFPN3d3dLS4u+ahlRxzeqFGjeOqpp5L9fTdt2oTNZuPatWvJ/t4iIpY7cgT+2xjb+OQT3hwzhpiYGBo2bEiDBg0sLi55KeykFrsdNm2CJUvMP+32FP24jh07YrPZsNlsuLu7U7hwYQYPHsyNGzcSnGez2Vi5cmX84+joaFq3bk2+fPk4cOBAgnMjIyN56qmnsNlsD9wlPaXcXq+IiDyE2Fjo1g0iI+Gll1idJQs///wzHh4efPrf+B1nom6s1ODvD/36wZkzN4/lz29uiZ6CO4HWq1ePuXPnEh0dze7du+nQoQM2m42JEyfe9fzw8HCaNWvG0aNH+e233yhcuHCC5wcPHkzevHnZt29fitUsIiKp4KuvYPNmyJiRyM8+o/9/LTkDBw6kaNGiFheX/NSyk9L8/aF584RBB+DsWfO4v3+KfbSnpyd58uTBz8+Pxo0bU6dOHdb/t6Hb7a5du8aLL77IuXPn7hp0fvrpJ9atW8dHH32UqM+22WzMnDmThg0bkjFjRkqUKMGOHTs4duwYNWvWxNvbm+eee47jx48neN3333/P008/TYYMGXj88ccZPXo0MTExABQqVAiAJk2aYLPZ4h/HWbhwIYUKFcLHx4dWrVpx/fr1+OciIyPp27cvuXLlIkOGDFSrVo0///wzwevXrFlDsWLF8PLyolatWtrYVESc07lzMGiQef+DD/h4+XJOnDhB3rx5GTFihLW1pRCFnaQyDHPxpcTcQkLMPUbutkh13LF+/czzEvN+j7DY9YEDB9i+fTseHh53PHf+/Hlq/DcobfPmzeTJkyfB8xcuXKBr164sXLiQjBkzJvozx44dS/v27QkICODJJ5+kTZs2dO/enWHDhrFr1y4Mw6BPnz7x52/dupX27dvTr18/Dh48yMyZM5k3bx4ffPABQHw4mTt3LkFBQQnCyvHjx1m5ciWrV69m9erVbN68mQkTJsQ/P3jwYJYvX878+fPZs2cPRYsWpW7duly5cgWA06dP07RpUxo1akRAQABdunRh6NChif6uIiIO4803ITgYnnmGM02axP8bO2nSpPhNtZ2Okc4FBwcbgBEcHHzHcxEREcbBgweNiIiImwdDQw3DjB2pfwsNTfT36tChg+Hq6mp4e3sbnp6eBmC4uLgYy5YtS3AeYHh4eBhPPvmkERYWdsf7xMbGGvXq1TPGjh1rGIZhBAYGGoCxd+/e+34+YLzzzjvxj3fs2GEAxldffRV/bMmSJUaGDBniH7/wwgvGuHHjErzPwoULDV9f3wTvu2LFigTnjBw50siYMaMREhISf2zQoEFG5cqVDcMwjNDQUMPd3d1YvHhx/PNRUVFG3rx5jUmTJhmGYRjDhg0zSpYsmeB9hwwZYgDG1atX7/tdk+quf69ERFLD8uXm74mbm2Hs22e0bt3aAIyqVasasbGxVleXJPf7/b6dWnacWK1atQgICGDnzp106NCBTp060axZszvOa9iwIUeOHGHmzJl3PDd16lSuX7/OsGHDkvz5ZcuWjb+fO3duAMqUKZPg2I0bNwgJCQFg3759jBkzhsceeyz+1rVrV4KCgggPD7/vZxUqVIhMmTLFP/b19eXixYuA2eoTHR1N1apV4593d3enUqVKHDp0CIBDhw5RuXLlBO/57LPPJvk7i4ikWdeuQe/e5v0hQ9gaHMySJUuw2WxMnToVm81maXkpSQOUkypjRggNTdy5W7ZAYqbvrVkDzz+fuM9OAm9v7/iBZnPmzKFcuXJ89dVXvPHGGwnOa9euHa+88gqdO3fGMAwGDhwY/9zGjRvZsWMHnp6eCV4Tt8v8/Pnz7/n5t67REPcf0d2OxcbGAhAaGsro0aNpepdB2w/aMPP29SBsNlv8+4qICDB4MJw/D8WKYR82jDf/+z+A3bp1o3z58hYXl7IUdpLKZgNv78Sd+9JL5qyrs2fvPt7GZjOff+klcHVN3jpv4+LiwvDhwxk4cCBt2rTBy8srwfMdOnTAxcWFTp06ERsby9tvvw3AZ599xvvvvx9/3rlz56hbty7ffPPNHS0hj+rpp5/m8OHD950J4O7ujj2J0/aLFCmCh4cH27Zto2DBgoA5xf7PP/+kf//+AJQoUYJVq1YleN3vv/+etC8gIpJWbdoE/23uyezZzFqwgH379pElS5YE/8Y7K4WdlOTqak4vb97cDDa3Bp645sLJk1M86MR57bXXGDRoENOmTYsPM7dq164dLi4udOjQAcMwGDRoEAUKFEhwTtzgtSJFipA/f/5kre+9996jYcOGFChQgObNm+Pi4sK+ffs4cOBA/H+MhQoVYsOGDVStWhVPT0+yZs36wPf19vamZ8+eDBo0iGzZslGgQAEmTZpEeHh4fCtXjx49+Pjjjxk0aBBdunRh9+7dzJs3L1m/n4iIJSIizDV1ALp353KpUrzTpAlgTiTJkSOHhcWlDo3ZSWlNm8KyZZAvX8Lj+fObx1NwnZ3bubm50adPHyZNmkRYWNhdz2nbti0LFy5k2LBh91yPJ6XUrVuX1atXs27dOp555hmqVKnCp59+Gt8aA/Dxxx+zfv16/Pz8ktTsOmHCBJo1a0a7du14+umnOXbsGD///HN8WCpQoADLly9n5cqVlCtXjhkzZjBu3Lhk/44iIqlu7Fg4ehTy5oWJExkxYgRXrlyhTJky9OjRw+rqUoXNMB5hPrMTCAkJwcfHh+DgYDJnzpzguRs3bhAYGEjhwoUfOGbkgex22LoVgoLA1xeqV0+1Fh1JW5L175WIyP3s2wcVKpi/QStX8mfevFSuXBnDMNi8eTPPJ2a8aBp1v9/v26kbK7W4ukLNmlZXISIi6UVMjLmjud0OzZtjb9iQnv8FnXbt2jl00EkqdWOJiIg4oylTYNcuyJIFpk5l1qxZ7N69m8yZM/Phhx9aXV2qUtgRERFxNidOwLvvmvc/+oiL/83IBXj//ffj1z5LL9SNJSIi4kwMA7p3N2dh1aoFnTszpHNnrl27xlNPPUXPnj2trjDVqWUnEdL5GG5JZvr7JCIpasEC+OUXyJABZs5k2/bt8UtpTJ8+HTe39NfOkf6+cRLErcobHh5+xyJ8Ig8rbuuL21d9FhF5ZBcuwIAB5v1Ro4gpXJhezZsD8MYbb6TbbXAUdu7D1dWVLFmyxO+xlDFjRqfeO0RSlmEYhIeHc/HiRbJkyYKrlh4QkeTWvz9cvQpPPQUDBzJt2jT2799PtmzZmDBhgtXVWUZh5wHy5MkDEB94RB5VlixZ4v9eiYgkm9WrYelSc6mTr74i6NIl3v1vkPL48ePTxUrJ96Kw8wA2mw1fX19y5cpFdHS01eWIg3N3d1eLjogkv+vXIW7g8cCB8PTTvN22LdevX6dSpUp06dLF2vosprCTSK6urvqREhGRtGn4cDhzBh5/HEaNYuPGjXz99dfYbDamTZuGi0v6no+Uvr+9iIiIo9u+HaZNM+/PmsUNF5f4Pa969uxJxYoVLSwubVDYERERcVSRkeaWEIYBnTrBCy8wfvx4jh49iq+vrzY0/o/CjoiIiKMaPx4OHYJcueCjj/j777/jZ11NmTIFHx8fiwtMGzRmR0RExBH9738Q13IzdSpG1qz0aNqUqKgo6tevT/P/1tcRteyIiIg4Hrvd7L6KjoZGjeC111iwYAGbN2/Gy8uLadOmaV24WyjsiIiIOJovvoDff4dMmWD6dC5dvsxbb70FwKhRoyhcuLDFBaYt6sYSERFxJKdOwbBh5v2JEyF/fgZ37szly5cpU6YMA+K2i5B4atkRERFxFIZhLh4YGgpVq0L37mzevJm5c+cCMHPmTO27dxcKOyIiIo5i6VJYswY8PGD2bCKjo+PX1OnevXu63ejzQRR2REREHMGlS9C3r3n/nXegRAkmTZrE33//Te7cuRk/fry19aVhCjsiIiKO4K23zMBTqhQMGcLhw4f54IMPAPj000/JmjWrxQWmXQo7IiIiad26dbBgAdhs8OWXxLq50bVrVyIjI6lbty6tWrWyusI0TWFHREQkLQsLg+7dzftvvglVqjBr1iy2bt2Kt7c3M2fO1Jo6D5Cmws6WLVto1KgRefPmxWazsXLlyvuev2nTJmw22x238+fPp07BIiIiKe299+DkSShQAD74gLNnzzJ48GAAxo0bR8GCBa2tzwGkqbATFhZGuXLlmBa3e2siHT58mKCgoPhbrly5UqhCERGRVPTnnzB5snl/xgwMb2969erF9evXqVy5Mr1797a0PEeRphYVrF+/PvXr10/y63LlykWWLFmSvyARERGrREebW0LExkLbtlC/Psu++45Vq1bh7u7OV199haurq9VVOoQ01bLzsJ566il8fX158cUX2bZt233PjYyMJCQkJMFNREQkzfnwQ9i/H7Jnh08/5cqVK/Tp0weA4cOHU6pUKYsLdBwOHXZ8fX2ZMWMGy5cvZ/ny5fj5+VGzZk327Nlzz9eMHz8eHx+f+Jufn18qViwiIpIIhw/DmDHm/cmTIWdO3n77bS5evEiJEiUYFrddhCSKzTAMw+oi7sZms7FixQoaN26cpNfVqFGDAgUKsHDhwrs+HxkZSWRkZPzjkJAQ/Pz8CA4OJnPmzI9SsoiIyKOLjYVatWDLFqhbF376iV82bODFF1/EZrOxbds2rZSM+fvt4+OTqN/vNDVmJzlUqlSJ33777Z7Pe3p64unpmYoViYiIJMGXX5pBJ2NGmDGD8IgIunXrBkDv3r0VdB6C04WdgIAAfH19rS5DREQk6c6dg0GDzPsffACFCvHe228TGBiIn58f48aNs7Y+B5Wmwk5oaCjHjh2LfxwYGEhAQADZsmWjQIECDBs2jLNnz7JgwQIAJk+eTOHChSlVqhQ3btzgyy+/ZOPGjaxbt86qryAiIvLw+vSBkBCoVAnefJPff/+dTz/9FIAZM2aQKVMmiwt0TGkq7OzatYtatWrFPx44cCAAHTp0YN68eQQFBXHq1Kn456Oionjrrbc4e/YsGTNmpGzZsvzyyy8J3kNERMQh+PvDihXg5gZffsmN6Gg6depEbGws7du3p0GDBlZX6LDS7ADl1JKUAU4iIiIp4to1KFECzp83dzQfO5ahQ4cyceJEfH19+d///qeNPm+TlN9vh556LiIi4hQGDzaDTvHiMGIEf/75Jx9++CFgdl8p6DwahR0RERErbdoEs2eb92fPJtJmo2PHjsTGxtK2bVteeeUVS8tzBmlqzI6IiEi6EhEBXbua93v0gOrVGTNiBAcPHiR37txMmTLF2vqchFp2RERErDJmDBw7BnnzwoQJ7N69m4kTJwLwxRdfkD17dosLdA4KOyIiIlYICDD3vwKYPp0oLy86deqE3W6nZcuWNGnSxNLynInCjoiISGqLiTF3NLfboXlzePVVPvjgA/766y9y5szJ1KlTra7QqSjsiIiIpLYpU2D3bsiSBaZOJSAgIH515GnTppEzZ05r63MyGqAsIiKSmk6cgHffNe9//DGRWbPS/qWXiImJoVmzZrz22mvW1ueEFHZERERSi2FA9+7mLKzataFTJ0YNHx7ffTV9+nSrK3RKCjsiIiKpZf58+OUXyJABZs5kx++/M2nSJABmzpxJrly5LC4wedntsHUrBAWBry9Urw6urqlfh8KOiIhIarhwAf7b85HRownz9aVDgwbExsbSrl07p5t95e8P/frBmTM3j+XPbw5Xato0dWvRAGUREZHU0K8fXL0K5cvDwIEMHTqUo0ePkj9/fj777DOrq0tW/v7mJLNbgw7A2bPmcX//1K1HG4FqI1AREUlpP/wAr7xi9uH88Qcbrl6lTp06AKxbt44XX3zR4gKTj90OhQrdGXTi2GxmC09g4KN1aWkjUBERkbQiJAR69TLvDxxIcJEidOrUCYBevXo5VdABc4zOvYIOmGO0T582z0stGrMjIiKSkoYPN3/9ixSBUaPo37s3p0+fpkiRIvGDk51JUFDynpcc1LIjIiKSUrZtg7jp5LNmseqXX5g3bx42m4358+fj7e1tbX0pwNc3ec9LDmrZERERSQmRkeaO5oYBnTvzb5kydC1dGoBBgwZRtWpViwtMGdWqQcaMEB5+9+fjxuxUr556NSnsiIiIpIRx4+DQIcidG2PSJLp16cLFixcpVaoUY8aMsbq6FDNmzP2DDsDkyam73o66sURERJLbgQMwfrx5f+pU5qxcycqVK/Hw8GDx4sV4enpaW18KmTcPxo417/fsabbg3Cp/fli2LPXX2VHLjoiISHKy280dzaOj4ZVXOFquHP2efhqADz74gHLlyllcYMrYsMHstQMYNsxs2Jo6VSsoi4iIOJ/p02HnTsiUiejJk3m9VSvCwsKoVasWA+NWUHYy//sfNGsGMTHQqhW8/7553NUVata0tDRA3VgiIiLJ59Qps1kDYOJEPpg/nz/++AMfHx/mz5+Pi4vz/eyePw8NGkBwMFStCnPnQlr7mmrZERERSQ6GYQ5UCQuDatXYUbYs79eoAcCMGTPw8/OzuMDkFxYGjRqZGe+JJ2DlSnOP07RGYUdERCQ5LF0Ka9aAhwdhkyfTrmVL7HY7bdu2pVWrVlZXl+zsdmjbFnbtguzZza+eI4fVVd1dGmtoEhERcUCXLkHfvub9d9+l7/TpHD9+nAIFCvD5559bW1sKeest+P578PQ0/yxa1OqK7k1hR0RE5FENHGgGntKlWfHEE8yZMwebzcbChQvJkiWL1dUlu88+gylTzPsLFphjddIydWOJiIg8ip9/hoULwWbj3/Hj6dqxIwBDhgzh+eeft7a2FLBqFfTvb94fPx5atLC0nERRy46IiMjDCg2F7t0BiH3zTVp+8gmXL1+mfPnyjB492uLikt+uXdC6tTkWu2tXGDLE6ooSR2FHRETkYb33HvzzDxQowCdZs/Lrr7/i7e3N0qVL8fDwsLq6ZPXPP9CwobkVxEsvwbRpN7d/SOvUjSUiIvIw/vgjfuDKof79GTpoEACff/45xYoVs7KyZHftGrz8Mly4AGXKwHffgbu71VUlnlp2REREkio62twSIjaWqBYtaPDZZ9jtdlq3bk2HDh2sri5ZRUVB8+bmKsl588KPP0LmzFZXlTQKOyIiIkn14Yfw118YOXLQJzqakydPUqhQIb744gtsjtK3kwiGAT16mPteeXvD6tXgiGsjKuyIiIgkxeHDMGYMAFubNmX2ihW4urqyZMkSfHx8LC4ueY0bd3P7h2+/hfLlra7o4SjsiIiIJFZsLHTrBpGRhFavToNFiwAYO3YsVapUsbi45PX11/DOO+b9qVPN/a8clcKOiIhIYs2eDVu2YHh70/LKFcLCw6lduzaDBw+2urJktWULdOpk3n/rLejVy9p6HpXCjoiISGKcPQv/hZqVFSuy5n//I3v27CxcuBBXV1eLi0s+hw9D48bmwORmzWDSJKsrenQKOyIiIg9iGNC7N4SEcPWJJ2i+eTMA8+bNI2/evBYXl3z+/dfsrrp6FSpXNheGdnGCpOAEX0FERCSF+fvD999juLnR8MIFYoH+/fvTsGFDqytLNhER8MorcOIEFC5sbgvh5WV1VclDYUdEROR+rl6FPn0AmJs7N9tDQnjmmWeYOHGixYUln9hYaN8efv8dsmSBNWsgVy6rq0o+CjsiIiL3M3gwnD/PhWzZ6Hn2LFmyZOGbb75xqu0ghg6FZcvMVZFXroQnn7S6ouSl7SJERETu5ddf4csvAWh25QpRwDdz51K4cGFr60pGM2eaayQCzJkDNWpYW09KSFMtO1u2bKFRo0bkzZsXm83GypUrE/3abdu24ebmxlNPPZVi9YmISDoSEWGuqQPM8fRkG9CvXz8aN25saVnJ6aefzHHXAKNHw+uvW1tPSklTYScsLIxy5coxbdq0JL3u2rVrtG/fnhdeeCGFKhMRkXRn9Gg4doyLHh4MiIzkmWeeYZIzzMP+z7590KIF2O3QoQO8+67VFaWcNNWNVb9+ferXr5/k1/Xo0YM2bdrg6uqapNYgERGRu9q7Fz76CIAuUVHYfHycapzOmTPmLuahoVC7NsyaBU60pdcd0lTLzsOYO3cuJ06cYOTIkYk6PzIykpCQkAQ3ERGReDEx5o7mdjvfAj8Ac+bMcZpxOtevQ8OG5hqJJUvC8uXgJBnunhw67Bw9epShQ4eyaNEi3NwS10g1fvx4fHx84m9+jrh9q4iIpJzJk2HPHq7abPQF+vbtS9OmTa2uKlnExJhdV/v2Qe7c8OOP5lRzZ+ewYcdut9OmTRtGjx5NsWLFEv26YcOGERwcHH87ffp0ClYpIiIO5fhxjPfeA2CgYVDAicbpGIa5XNDateZigT/8AIUKWV1V6khTY3aS4vr16+zatYu9e/fS57/FnmJjYzEMAzc3N9atW0ft2rXveJ2npyeenp6pXa6IiKR1hgHdu2OLiGAD8EO2bOxZtsxpfjM++sicZm6zmTuaP/OM1RWlHocNO5kzZ+avv/5KcGz69Ols3LiRZcuWOU3fqoiIpJL582HDBiKA7sDir7+mQIECVleVLL77Ln4PUz75xNzoMz1JU2EnNDSUY8eOxT8ODAwkICCAbNmyUaBAAYYNG8bZs2dZsGABLi4ulC5dOsHrc+XKRYYMGe44LiIicl8XLhDTrx9uwHtAu1GjqFu3rtVVJYsdO6BdO/P+m29Cv37W1mOFNBV2du3aRa1ateIfDxw4EIAOHTowb948goKCOHXqlFXliYiIk4rq2ROPkBD2AP978UVWO8miM8ePm5t7RkZCo0bw6afOPcX8XmyGYRhWF2GlkJAQfHx8CA4OJnPmzFaXIyIiqcz4/ntsjRsTA7ySOzcLDhwgR44cVpf1yC5fhueegyNHoEIF2LwZvL2trir5JOX3O0217IiIiKQKux22boUTJ4jo25eMwGQXF0atWuUUQScyEpo0MYNOgQLmzCtnCjpJpbAjIiLpi7+/OXDlzBkAMgLRQNkOHahUqZKlpSUHw4DOnc0slzmzuZaOr6/VVVnLYdfZERERSTJ/f2jePD7oxHEDXpw3z3zewb33njm13M3NXB1Zc3YUdkREJL2w280WnbsMVbX9d6N/f/M8BzVnDrz/vnl/5kyoU8faetIKhR0REUkftm69o0UnAcOA06fN8xzQL79A9+7m/REjzK4sMSnsiIhI+hAUlLznpSEHDkCzZubeV61bw9ixVleUtijsiIhI+pA7d+LOc7DRvEFB8PLLEBIC1avD3Lnpcy2d+1HYERER52cYhH71lXn3XufYbODnZyYGBxEWZi4WeOoUPPEErFgBTrKVV7JS2BEREedmGMS8+SaPff01sXGHbm/6iHs8eTK4uqZmdQ/Nbje7rHbvhhw5YM0ayJ7d6qrSJoUdERFxXoaBMXQobtOmAdDP25uL06djy5cv4Xn588OyZdC0qQVFPpyBA83FAj094fvvoWhRqytKu7SooIiIOK8xY7BNmgRAb5uNJitXkrtOHejWzZx1FRRkjtGpXt1hWnQApkyBzz4z7y9caG4LIfemsCMiIs5p4kQYNQqAAUDhSZOoE7fwjKsr1KxpVWWP5PvvYcAA8/7EifDaa9bW4wgUdkRExPlMmQJDhwIwDLjQujWfvPWWtTUlgz//NMfpGIa5ps6gQVZX5BgUdkRExLnMmGGuhAyMBn4qV47tX36JzcHnY588ac68ioiAevXg8881xTyxFHZERMR5zJsHPXsCMBH4PHt2/lixgowZM1pa1qO6dg0aNIALF6BcOfj2W3PvK0kcXSoREXEOS5bAG28AMAV4182NjStXUrhwYWvrekRRUebqyIcOQd68sHo1ZMpkdVWORWFHREQcn78/tGsHsbHMBPoDc2bNolq1ahYX9mjixuZs3AiPPQY//mjOkpek0To7IiLi2FavhlatwG5noYsLPYG3336bTp06WV3ZI3v/fbNnztXV7Lp66imrK3JMCjsiIuK41q83+3iio1mRIQMdY2N5uWFDJkyYYHVlj2zRInjvPfP+559D/frW1uPIFHZERMQxbd4Mr74KUVFs9PGhxY0blCxdmsWLF+PqQAsE3s3mzdC5s3l/0CDo0cPaehydwo6IiDie7dvNrb4jItidJw/1g4PJkiMHP/zwA5kzZ7a6ukfy99/QpAlER0Pz5uAEjVSWU9gRERHHsmuX2acTFsbxwoWpev48hrs7K1asoFChQlZX90guXjSnmF+9ClWqwIIF4KJf6kemSygiIo5j3z546SUICeHik09SNjCQSGCWE8y8ioiAV16BwEB4/HFYtQq8vKyuyjko7IiIiGM4eBBefBGuXiW0dGlKnjhBOObMq44dO1pd3SOJjTVnzu/cCVmzwpo1kDOn1VU5D4UdERFJ+44ehRdegH//JbJ0acqfP8/lqCgaOsnMqyFDYPly8PCAlSuheHGrK3IuCjsiIpK2BQZC7dpw/jwxJUtSMzKSY5cuUa5cOb7++muHn3n1xRfw0Ufm/Tlz4Pnnra3HGSnsiIhI2nX6tBl0zpwhtnhxmvv48PvRo+TLl48ff/yRTA6+b8KaNdCnj3l/7Fho29baepyVwo6IiKRNQUFm19XJkxhFi9K3ZEm+37GDTJkysWbNGvLly2d1hY9k715o0cIcr9OxI4wYYXVFzkthR0RE0p6LF82gc/QoFCzIx/XrM23FClxdXVm2bBlly5a1usJHcuYMNGwIYWHm15w5E2w2q6tyXkkKOzVq1OD3339PqVpERETgyhVz1tWhQ5AvH9/17MmgqVMBmDlzJi+99JLFBT6akBBzPcRz56BkSVi2zByYLCknSWHn+vXrVK1alSZNmnDo0KGUqklERNKr4GBzHZ39+yF3braNHUubd94BYMSIEbzxxhsWF/hooqPNrqv9+yFPHnPMTpYsVlfl/JIUdvbs2cOiRYs4cOAAZcuWpXPnzpw+fTqlahMRkfTk+nVzZeTduyFHDo588QX1+/UjJiaGNm3aMHbsWKsrfCSGAb17w88/Q8aM8MMPULCg1VWlD0kes9O6dWsOHTrE5MmT+emnnyhWrBhvvfUWly9fTon6REQkPQgPh0aNYMcOyJKFi4sX80Lfvly/fp0aNWowZ84cbA4+qGXSJJg92xybs2QJVKxodUXpx0MNUHZzc6N3794cP36c4cOH89VXX1GkSBHef/99wsPDk7tGERFxZjduQOPG5lbfmTIRtmIFdYcM4cyZMzz55JOsWLECT09Pq6t8JN9+C0OHmvcnTza3hZDU80izsTJmzMi7777LX3/9RZkyZRg5ciSPP/4406ZNIyYmJrlqFBERZxUVZW7tvX49eHsT88MPNJswgYCAAHLlysWaNWvImjWr1VU+ku3boX17836/ftC3r7X1pEc2wzCMpLwgODiYvXv3smfPHvbu3cvevXs5fPgwdrsdgEKFCnH69GkKFy7M9OnTqVOnTooUnlxCQkLw8fEhODiYzJkzW12OiEj6ERMDLVuCvz9kyEDsjz/SYe5cFi1aRMaMGfn111+pVKmS1VU+kmPHzN3LL1+GV181t4Rw8AWf04yk/H67JeWNCxcuzKlTpwAwDIN8+fLxzDPP0LZtWypWrEjFihXJmjUrJ06cYNiwYdSrV49Zs2bRuXPnh/82IiLifOx2s7nD39+cd/399wxes4ZFixbFr6Xj6EHn8mVo0MD8s2JFWLxYQccqSWrZadCgAZUqVaJixYo888wz5M6d+77nd+/enbVr1/LPP/88cqEpRS07IiKpLDYW3ngD5s0DNzfw9+ejw4cZNGgQAPPnz6d9XL+Pg7pxw1wq6LffzBlXv/9uTjWX5JNiLTtr1qxJUiE1a9Zk9uzZSXqNiIg4sbj51/Pmmc0cS5ey4OrV+KDz4YcfOnzQiY2FTp3MoOPjAz/+qKBjtSSFnaRq0KABy5YtS8mPEBERR2EYMGAAzJhhzr9esIAfM2Sgc8uWALz11lu8/fbbFhf56N59F5YuNRutli+HUqWsrkhSNOz4+PjQtGnTlPwIERFxBIYBw4bBlCnm46++4vfHH+e12rWx2+28/vrrTJo0ydoak8FXX8G4ceb92bPNfa/EemlqI9AtW7bQqFEj8ubNi81mY+XKlfc9/7fffqNq1apkz54dLy8vnnzyST799NPUKVZERBJvzBiYONG8P306h6pU4eWXXyYiIoJ69eoxZ84cXFzS1E9Skq1fD927m/fffdfcyVzShhRt2UmqsLAwypUrR+fOnRPVIuTt7U2fPn0oW7Ys3t7e/Pbbb3Tv3h1vb2+6deuWChWLiMgDTZwIo0aZ9z/5hDONGlH3uee4cuUKlSpV4rvvvsPd3d3SEh/VgQPmckF2O7RtC6NHW12R3CrJ6+ykFpvNxooVK2jcuHGSXte0aVO8vb1ZuHBhos7XbCwRkRQ0ZQr072/eHz+eK926Ub16dQ4ePEjx4sX57bffyJEjh6UlPqpz58y1dE6fhuefh3XrwMEXfHYISfn9duw2w9vs3buX7du3U6NGjXueExkZSUhISIKbiIikgBkzbgad994jtE8fGjZsyMGDB8mbNy8///yzwwed0FBzS6/Tp6F4cVixQkEnLXKKsJM/f348PT2pWLEivXv3pkuXLvc8d/z48fj4+MTf/Pz8UrFSEZF0Yu5c6NnTvD94MDeGDqVx48bs2LGDLFmysHbtWgo6+Jbfdju0bg179kCOHOYU82zZrK5K7sYpws7WrVvZtWsXM2bMYPLkySxZsuSe5w4bNozg4OD42+nTp1OxUhGRdGDJEnPRQIC+fYkeO5aWrVqxYcMGvL29+emnnyhTpoy1NT4iwzAbrVavhgwZYNUqKFLE6qrkXtLUAOWHVbhwYQDKlCnDhQsXGDVqFK1bt77ruZ6eng6/e66ISJrl7w/t2plpoHt37B9/TMcOHVi1ahWenp788MMPVKlSxeoqH9mUKfD55+b9hQvh2WetrUfuzynCzq1iY2OJjIy0ugwRkfRn9Wpo1crs3+nQAWPaNHr16sXXX3+Nm5sby5Yto1atWlZX+chWrICBA837H35ozsKStC1NhZ3Q0FCOHTsW/zgwMJCAgACyZctGgQIFGDZsGGfPnmXBggUATJs2jQIFCvDkk08C5jo9H330EX379rWkfhGRdGvdOmjWDKKjoVUrjC+/ZNCQIcyaNQubzcaiRYto2LCh1VU+sj/+MKeWGwb06AFvvWV1RZIYaSrs7Nq1K0HqH/hfdO7QoQPz5s0jKCgoftd1MFtxhg0bRmBgIG5ubhQpUoSJEyfSPW5VJxERSXmbNkHjxhAVBU2awIIFvD9+PB9//DEAs2fPpuV/W0I4ssBAc+ZVRATUrw9Tp5q7Xkjal2bX2UktWmdHROQRbN8OL70EYWHw8svg78/k6dMZMGAAAJ9++in946afO7CrV6FqVTh0CJ56CrZsgUyZrK4qfUu36+yIiEgq2rXLbOIIC4M6dWDZMr5auDA+6IwZM8Ypgk5UlNlDd+gQ5MtnDk1S0HEsCjsiIpJ0+/aZLTohIeaywd9/zzfff0/Xrl0BePvtt3nnnXcsLvLRGQZ07Qq//gqPPWaupZMvn9VVSVIp7IiISNIcPGi25Fy9au6TsHo1y3/6ibZt22IYBt27d2fSpEnYnGBAy9ixsGABuLrCd99BuXJWVyQPQ2FHREQS7+hReOEFuHQJKlSAn35i5YYNtGrVCrvdTrt27Zg+fbpTBJ2FC2HkSPP+9OlQr5619cjDU9gREZHECQyE2rXh/HkoWxZ+/pkftm6lRYsWxMTE0KZNG+bOnYuLi+P/tGzadHMR6MGDoVs3S8uRR+T4fyNFRCTlnT5tBp0zZ6BECVi/nh9//51mzZoRHR1Nq1atmD9/Pq6urlZX+sgOHTJn0EdHw2uvwfjxVlckj0phR0RE7i8oyAw6J09C0aLwyy+s3bOHpk2bEh0dzWuvvcbChQtxc0tTS7c9lAsXoEEDuHbN3AJi/nxwgoaqdE//E4qIyL1dvGiO0Tl2DAoVgo0bWf+//9G4cWOioqJo2rQpixcvdoqgEx4Or7xiZroiReD778HLy+qqJDko7IiIyN1duQIvvmj26+TPDxs3suHIEV555RUiIyN59dVXWbJkCe7u7lZX+shiY839S//4A7JlgzVrIGdOq6uS5OL4UVxERJJfcLC5js7+/ZAnD2zYwK8nT9KoUSNu3LhBo0aN+Pbbb/Hw8LC60odit8PWrWYPna8vrFplbtju4QErV0KxYlZXKMlJYUdERBK6ft1cGXn3bsiRA375hS3nz9OwYUMiIiJo0KAB3333ncMGHX9/6NfPHGt9u3nzoHr1VC9JUpjCjoiI3BQeDg0bwo4dkDUr/PILGy9coFGjRoSHh1OvXj2WL1+Op6en1ZU+FH9/aN7cXBn5bhz0a8kDaMyOiIiYbtwwdy/fsgUyZ4aff+bn8+d5+eWXCQ8Pp27duvj7+5MhQwarK30odrvZonOvoGOzQf/+5nniXBR2RETE3O2yeXNYvx68veGnn1h94QKvvPJK/BidlStX4uXA05O2br1711UcwzCXE9q6NfVqktShsCMikt5FR0OrVuYulxkymHtdBQXRpEkToqKiaNasGcuWLXPYFp04QUHJe544DoUdEZH0zG6H9u1hxQpzKtL337MkKIiWLVvGbwGxdOlShx2MfKvETiX39U3ZOiT1KeyIiKRXsbHQpQssXQpubrB8OfODgmjbti12u52OHTuyYMECp1gw8OzZm5t63ovNBn5+mo3ljBR2RETSI8OAXr3MudaurrB0KbPOnaNTp04YhkG3bt346quvnGKvq19/haefhu3bIWNG89jtm7LHPZ482bwc4lwUdkRE0hvDgAEDYOZM81d+wQKmnjtH9+7dMQyDvn37MmPGDIffvdwwYOJEqFPH3PWibFnYtw+WL4d8+RKemz8/LFsGTZtaU6ukLMdvmxQRkcQzDBg2DKZMMR9/9RUfnTvHoEGDABg0aBATJ07EdnvTh4MJDoaOHc3VkAE6dIDp082WnaJF4dVXE66gXL26WnScmcKOiEh6MmaM2dwBGNOn8+7x43zwwQcAvPvuu4wePdrhg87+/dCsmbl3qYcHfP65OTTp1q/l6go1a1pWoqQyhR0RkfRiwgQYNQqA2E8+oc9ff/HFF18AMH78eIYOHWphccljwQLo0QMiIqBgQbNrqmJFq6sSqynsiIikB5Mnm91XQMz779Pujz9YunQpNpuNGTNm0K1bN2vre0SRkebqxzNmmI/r1YNFiyB7dkvLkjRCYUdExNnNmGEOSAaihw/n1W3b+Omnn3B3d2fhwoW0bNnS4gIfzT//mIs/79pldlWNGgXvvAMOPr5akpHCjoiIM5s7F3r2BOBGv37U2bSJbdu34+Xlhb+/P/Xq1bO4wEezdi20bQtXrkC2bPD111C3rtVVSVqj3Csi4qy+/hreeAOA8C5dqPLrr2zbvp0sWbKwfv16hw46sbEwejQ0aGAGnWeegT17FHTk7tSyIyLijJYvN7eBMAyut2nD07/+yrHjx8mdOzfr1q2jbNmyVlf40C5fhtdfN1t1wGy4+vRT8PS0ti5JuxR2RESczerV5saedjvXGjemzKZNnDl3jkKFCrF+/XqKFi1qdYUP7c8/zfE5p06Bl5e5LmK7dlZXJWmdwo6IiDNZt85cZCYmhn/r1KHk5s1cunqVUqVK8fPPP5Pv9qWDHYRhwOzZ8OabEBVlLgy4fLm5KrLIg2jMjoiIs9i0CRo3hqgozlWuzONbt3Lp6lUqV67M5s2bHTbohIdDp07QvbsZdBo3NmdeKehIYinsiIg4g+3boWFDiIjgZOnSPP7HH4RGRtKwYUM2bNhAdgddcObYMXj2WZg/35xKPnEi+PuDj4/VlYkjUdgREXF0u3ZB/foQFsbRQoV48sABIg2Drl27smLFCry9va2u8KF8/z1UqGBu/5ArF2zYAIMH37ljuciDKOyIiDiyffvgpZcgJIS/c+em3MmTRAJjxoxh5syZuLk53tDMmBgYOtTsrgoJgapVYe9e7WUlD8/x/isQERHTwYNQpw5cvcr/fHyocuECUa6ufDVrFp07d7a6uody4QK0bg2//mo+HjDA7Lpyd7e2LnFsCjsiIo7oyBF44QW4dImDXl48FxxMbMaM/LBsGfXr17e6uoeybRu0aAHnzsFjj8GcOfDaa1ZXJc5A3VgiIo4mMNAMOufPc8jdneoREWTIlYvNmzc7ZNAxDJgyxeymOncOSpSAP/5Q0JHko5YdERFHcvo01K4NZ85w2MWFGtHRZCtalLVr11KkSBGrq0uy0FDo0gW++cZ83KqVuZ7OY49ZW5c4F7XsiIg4iqAgM+icPMlRoGZsLIUrVWL79u0OGXQOHYJKlcyg4+YGn31mbueloCPJTWFHRMQRXLyI8cILcOwYgUBtoPprr7Fp0yZy5sxpdXVJ9s035uadhw5BvnywebO5OrKmlUtKUNgREUnrrlwhtk4dbIcOcRoz6HQYMYKlS5fi5eVldXVJEhUF/fub3VVhYWZD1Z498NxzVlcmzkxjdkRE0rLgYKJr1cL9r78IAuq5uTHqyy/p0KGD1ZUl2dmz5myr7dvNx8OGwdix4OpqbV3i/NJUy86WLVto1KgRefPmxWazsXLlyvue7+/vz4svvkjOnDnJnDkzzz77LD///HPqFCsiktzsdnN/qyVLzD+vXSOsRg3c9+/nX6CZjw9fbNjgkEFn40Z4+mkz6Pj4mKsjjxunoCOpI02FnbCwMMqVK8e0adMSdf6WLVt48cUXWbNmDbt376ZWrVo0atSIvXv3pnClIiLJzN8fChWCWrWgTRuoVQt7zpx479vHFeCNggWZ/+efPP/881ZXmiSGYS4K+OKLcPEilCsHu3fDK69YXZmkJzbDMAyri7gbm83GihUraNy4cZJeV6pUKVq2bMl777131+cjIyOJjIyMfxwSEoKfnx/BwcFkzpz5UUoWEXk4/v7QvLmZDG5jADMLF6bFrl1ky5Yt9Wt7BNeuQceOZisOmPenTwcHG2YkaVRISAg+Pj6J+v1OUy07jyo2Npbr16/f9x+E8ePH4+PjE3/z8/NLxQpFRG5jt0O/fncNOnG6R0eTzcG2+d63DypWNIOOhwfMmmWuiKygI1ZwqrDz0UcfERoaSosWLe55zrBhwwgODo6/nT59OhUrFBG5zdatcObMPZ+2AbYzZ8zzHMT8+VClChw/DgULmttAdO2qaeViHaeZjfX1118zevRovv/+e3LlynXP8zw9PfH09EzFykRE7iMoKHnPs9CNG2Yj1axZ5uP69WHRInCw3jdxQk7RsrN06VK6dOnCt99+S506dawuR0Qk8a5eTdx5vr4pW8cjOnkSqlc3g47NBmPGwOrVCjqSNjh8y86SJUvo3LkzS5cu5eWXX7a6HBGRxAkPJ3bECJg8GRfMgch37eWx2SB/fjNJpFFr10LbtnDlihluvv4a6ta1uiqRm9JU2AkNDeXYsWPxjwMDAwkICCBbtmwUKFCAYcOGcfbsWRYsWACYXVcdOnRgypQpVK5cmfPnzwPg5eWFj4MN5hORdGTTJuydO+MaGAjAr0BNwLDZsN06UDlukMvkyWlyQZrYWLMFZ8wYc3z1M8/Ad9+Z43RE0pI01Y21a9cuypcvT/ny5QEYOHAg5cuXj59GHhQUxKlTp+LPnzVrFjExMfTu3RtfX9/4W79+/SypX0TkvkJCoGdPqFUL18BATgNNPDy4uHQptuXLseXLl/D8/Plh2TJo2tSScu/n8mV4+WUYPdoMOj17mmOoFXQkLUqz6+yklqTM0xcReWhr10K3bvDfDNAZwNT8+Vm0alX8/8HDbjcTQ1CQOUanevU02aLz55/mskCnTplTyWfOhHbtrK5K0puk/H6nqW4sERGnc+UKDBxozscGjgNdgNjnn2fTsmUJdyx3dYWaNa2oMlEMwxyA3LevuaFn0aLmeohlylhdmcj9paluLBERp7JiBZQsCfPnEwt8CpQFKrz1Fr/88kvCoJPGhYebKyD36GEGnSZNYNcuBR1xDGrZERFJbhcvQp8+5mhd4KirK+3tdg489hjz586lefPmFheYNEePmt1W+/eDiwtMmABvv61FAsVxKOyIiCQXwzB3LO/bFy5fJtbFhYmGwWi7nSIlS7LL35/ixYtbXWWSrFwJHTqYY6tz54alS9N0T5vIXakbS0QkOZw9a27l3bYtXL7MySxZqBgby3DDoEmrVuzcudOhgk5MDAwZYnZXhYRAtWqwZ4+CjjgmhR0RkUdhGPDll+bYnNWriXV3Z0rOnDxx7Rp/ubnx2Wef8fXXX/PYY49ZXWmiXbgAL74IkyaZjwcMgI0bIW9ea+sSeVjqxhIReViBgeYOlxs2AHCpSBHqnzvHrn//JW/evHz33Xc899xzFheZNNu2wWuvmbPfH3vM3Kn8tdesrkrk0ahlR0QkqWJjYepUKF0aNmzAyJCBpc88Q+7jx9kVEUGtWrXYu3evQwUdwzAXaq5Z0ww6JUua6+ko6IgzUNgREUmKw4fh+efNQcjh4YRVrEj9fPlo/eef4OLCyJEjWb9+Pbly5bK60kS7fh1atTK7q2JioHVr2LkTnnzS6spEkoe6sUREEiMmBj7+GEaOhMhIjMceY9srr/DS8uVEREbi6+vL119/TU0HG8F78CA0awZ//w1ubvDpp9C7t6aVi3NR2BEReZD9+6FzZ9i9G4CYF16gn5cX07/+GoD69eszf/58h1okEMxp5F26QFgY5MtnLgv07LNWVyWS/NSNJSJyL1FRMGoUVKhgBp0sWQgcOZLiJ04wffVq3NzcmDRpEqtXr3aooBMVBf36md1VYWFQu7Y5rVxBR5yVWnZERO7mzz/N1pwDBwAwGjdm9lNP0eeDD4iOjqZgwYIsXbqUKlWqWFxo0pw5Ay1awI4d5uPhw2HMmDS536hIslHLjojIrSIiYPBgqFLFDDo5c3Jt5kwaRETQfdQooqOjadq0KXv37nW4oLNxIzz9tBl0fHxg1Sr44AMFHXF+atkREYmzdSu88Ya5GRRAmzb8VK8e7QcO5NKlS2TIkIGPPvqIXr16YXOgEbyxseYCgSNGmPefegqWLYMiRayuTCR1qGVHRCQ0FN5805xSfvQo5M3LjW+/pZu3Nw3at+fSpUuUK1eOXbt20bt3b4cKOteumVs+DBtmBp1OnWD7dgUdSV/UsiMi6dv69eYqyP/8Yz5+4w12t25N6549OXr0KDabjbfffpuxY8fi6elpba1JFBBg7lZ+/Dh4esLnn5sNVw6U1USShcKOiKRP167BW2+Z+yEAFCyIfcYMxu/axai6dbHb7eTPn58FCxZQq1YtS0t9GPPmQc+ecOMGFCpkdltVqGB1VSLWUNgRkfRn1SozCZw7ZzZz9OnDyW7deL1HD7Zt2wZAy5Yt+eKLL8iaNavFxSbNjRvm4s6zZ5uPGzSAhQshWzZr6xKxksbsiEj68e+/0KYNvPqqGXSKFcPYvJmvypWj7HPPsW3bNjJlysSCBQtYsmSJwwWdkyehWjUz6Nhs5pTyH35Q0BFRy46IOD/DgG+/hT594NIlcHGBt9/mTJcudO3bl7Vr1wJQrVo1Fi5cSKFChayt9yH89BO0bQtXr0L27PD11/DSS1ZXJZI2qGVHRJxbUJA5HalVKzPolC6NsWMH80uWpPQzz7B27Vo8PT358MMP2bRpk8MFHbvd3K7r5ZfNoFOpkrkasoKOyE1q2RER52QY5ijdgQPNwcju7jBiBOc6dKD7m2+yevVqACpVqsS8efMoUaKEpeU+jEuXzNacdevMxz17mht5OtikMZEUp5YdEXE+//wD9eqZ2z1cuwYVK2Ls2sWiIkUo/fTTrF69Gg8PD8aPH8+2bdscMuj8+ac5u2rdOvDyggULYPp0BR2Ru1HLjog4j9hYmDEDhgwxFwr09IQxYzjfpg09+vTh+++/B6BChQrMnz+fUqVKWVxw0hkGzJxpbuQZFQVPPAHLl0OZMlZXJpJ2qWVHRJzD0aNQqxb07m0GnapVMQICWJQ3L6XKleP777/H3d2d999/nx07djhk0AkPhw4dzO6qqChzKNKffyroiDyIWnZExLHZ7TB5MrzzjrnIjLc3TJhAYP369Ozdm59//hmA8uXLM2/ePMqWLWttvQ/p6FFo1gz++svcuHPCBHNNRK2GLPJgatkREcf1v//Bc8/B22+bQadOHez79vFpdDSly5bl559/xtPTkw8++ICdO3em+aBjt8OmTbBkifmn3W4eX7kSKlY0g07u3LBhg/mVFXREEkctOyLieKKjzaaNsWPN+z4+8PHH7KtQgS6tWrFr1y4AatSowaxZsyhWrJjFBT+Yv785DufMmZvH8uc3Q87KlebjatXM5YJ8fS0pUcRhqWVHRBzLnj3wzDPw3ntm0GnUiIhduxh27BgVKlZk165d+Pj4MHv2bDZu3OgwQad584RBB8zHcUFn4EDYuFFBR+RhqGVHRBzDjRvm/geTJpn9O9mzw9Sp/Jo7N93q1+fYsWMANG/enM8++wxfB0kFdrvZomMY9z4ne3bza7u6pl5dIs5ELTsikvZt3w7ly8P48WY6aNmSfzdvpuPPP1P7hRc4duwY+fLlY+XKlXz33XcOE3QAtm69s0Xndpcvm+eJyMNRy46IpF1hYTBiBHz2mdn0kScP9mnTmHXhAsOrVePatWvYbDZ69OjBhAkTyJw5s9UVJ1lQUPKeJyJ3UtgRkbRp40bo0gUCA83HHTuy5/XX6T50aPwA5PLly/PFF19QuXJlCwt9eAEBMGtW4s51oMYqkTRH3VgikrYEB0P37vDCC2bQ8fPj+nff0TNDBiq++GL8AOTPP/+cP//80yGDzvbt5sad5cubU8zvx2YDPz+oXj1VShNxSmrZEZG048cfzaBz9iwARs+eLCpdmoE9e3Lp0iUA2rVrx4cffkju3LmtrDTJDAN++QU++AA2bzaPubhAy5bm9PK33755Xpy4dXQmT9bgZJFHoZYdEbHe5cvQrh00bGgGnaJFOTp7NtX376d9795cunSJUqVKsXnzZhYsWOBQQSc2FlasgEqV4KWXzKDj7m720P39N3z9tTmtfNkyyJcv4Wvz5zePN21qTe0izkItOyJirWXLzP2sLl4EFxfCu3dn8I0bTO/WDcMw8Pb2ZtSoUfTr1w93d3erq020mBhYutScQHbwoHnMy8tsuHrrLTPI3KppU3j1VXPWVVCQOUanenW16IgkB4UdEbHG+fPQp4+5ZTcQW6IES+rUodf8+YSEhADQsmVLPvroI/LfngzSsBs3YN48c12cuLHVPj7mV+3XD3LmvPdrXV2hZs3UqFIkfVHYEZHUZRiwaJH5y3/1KoabG0ebNqXJ7t0cnDoVgAoVKjB58mSqVatmcbGJFxoKM2fCxx/fnCaeMycMGAC9epmBR0SsobAjIqnn9GmzH+ennwCIKFGCAT4+zPz2WwBy587N+PHj6dChAy4ujjGk8MoV+PxzmDLFvA9mF9XgwfDGG5Axo7X1iUgaG6C8ZcsWGjVqRN68ebHZbKyM2xTmHoKCgmjTpg3FihXDxcWF/v37p0qdIpJEsbFms0epUvDTTxgeHqyqUoWshw8z8/ff8fDwYMiQIRw5coROnTo5RNA5f94MNAULwsiRZtB54gn46is4fhzefFNBRyStSFP/ooSFhVGuXDmmTZuWqPMjIyPJmTMn77zzDuXKlUvh6kTkoRw/DnXqQI8ecP06ZwsUoLKHB6/+/juRsbE0adKEgwcPOswKyCdPmuOpCxWCDz80u6/KljUHIx86BJ07g4eH1VWKyK3SVDdW/fr1qV+/fqLPL1SoEFOmTAFgzpw5iXpNZGQkkZGR8Y/jBkKKSDKz22HqVBg+HCIiiPbw4AMvL8aeOkUsULZsWT799FNq165tdaWJ8vffMGECLF5szrQCqFLF3M3i5ZdvrokjImlPmmrZSQ3jx4/Hx8cn/ubn52d1SSLO59Ahc970gAEQEcEf3t4Uj4pidHAw+fz8mD9/Pnv27HGIoLNnDzRvDiVLwvz5ZtB58UX49VdzJeSGDRV0RNK6dBd2hg0bRnBwcPzt9OnTVpck4jyio2HcOHjqKdixgzBXV7oBlcPCuJolC5MmTeLIkSO0b98e1zS+gMzWrVC/PlSoYM6ONwxo3Bh27oR168wp4go5Io4hTXVjpQZPT088PT2tLkPE+QQEmANW9u4F4Eegh93ORQ8P3nrzTYYPH062bNksLfFBDAN+/tnc0uG338xjLi7QujUMG2aOrxYRx5Puwo6IJLPISHj/fYwJE7DFxHAZ6A8sAl5//XXGjh1LoUKFLC3xQeK2dBg3zuy2AnOQcadOMGgQFClibX0i8mgUdkTkwez2u+9jsHMn0e3b437kCDbgO6AP8NRLL7FnwgTKly9vceH3Fx1t7k01YYI5ABnM6eI9epj7Vd2+V5WIOKY0FXZCQ0M5duxY/OPAwEACAgLIli0bBQoUYNiwYZw9e5YFCxbEnxMQEBD/2n///ZeAgAA8PDwoWbJkapcv4pz8/c3Vjs+cuXksXz4iS5fGfd063A2DC0Av4NLzz7Ps/fepXr26VdUmSkQEzJ1rbunwzz/msSxZzLVx+vaFHDksLU9EkpnNMAzD6iLibNq0iVq1at1xvEOHDsybN4+OHTty8uRJNm3aFP+c7S4jBAsWLMjJkycT9ZkhISH4+PgQHBzsEGt8iKQqf39zKtJt/0wYQNx/eQuARU8/zeCJE3nhhRfu+t9kWnH9OsyYYW7pcOGCeSxXLrMVp2dP0D8BIo4jKb/faSrsWEFhR+Qe7HZz5bxbW3RuYQDXXF3Z7u9Pg0aN0nTIuXwZPvvMXPbn6lXzWIEC5grInTubu5GLiGNJyu93murGEpE0ZOvWewYdMFt2strtvJw5c5qdg33uHHzyidmaExZmHiteHIYOhbZtwd3d2vpEJHUo7IhIQoYBO3cSOWoUiVqkIW6L7zQkMNAcjzNnDkRFmceeespc7bhJE3NstYikHwo7ImI6fhwWLSJm3jzcTp5MXNABc3ZWGnHwIIwfD0uWmL1wAFWrmiGnXr002wAlIilMYUckPbt8Gb79FhYtMvc+wPxHIQxYAbzs5kaWmBjumhFsNsif35yGbrFdu8w1clasuHmsbl1zW67nn7euLhFJGxR2RNKbyEj48UdYuND8MzoaADvwC+ZigP9Wrcrbo0aRJTgY22uvma+7dS5DXBPJ5MmW9QkZBmzZYoacdetuHm/a1Aw5FSpYUpaIpEEKOyLpgWHAtm1mwPn2W7h2Lf6pvZgBZwlQ5qWXGDFiBM/f2hyybNmd6+zkz28GnaZNU6f+WxgG/PSTuaXDf41RuLpCmzbmwGMtsSUit1PYEXFmR46YAWfRIrhl7akLbm7Mi4lhIfC3qyutWrXip0GDKFeu3J3v0bQpvPrq3VdQTkV2u7kh57hxsG+feczT05w6PmgQFC6cquWIiANR2BFxNv/+C998Y4acP/6IPxzl6ckKV1dmhoezOSYGL29vunXrxpr+/SlQoMD939PV1dzm2wJRUbB4sbmlw5Ej5jFvb3MRwIED09T4aBFJoxR2RJxBRAT88IMZcNauhZgYAAxXVw76+fHx+fMsvXGDCCBPnjy837cvPXr0IGvWrNbWfR8REfDll/Dhh3D6tHksa1ZzO4e+fSGNb6AuImmIwo6Io4qNNUfoLlxojqsJCYl/KviJJ1jq5sZ7hw5x8b/uqyeffJK3336b119/HU/PRE8sT3UhITB9Onz6KVy8aB7Lkwfeegu6d4dMmaytT0Qcj8KOiKM5eNAcg7N4MZw6FX841s+PPSVK8N6RI/x09Gj88Zdffpm+fftSp04dXFxcrKg4US5dgilTzC0dgoPNY4UKmVs6dOoEGTJYWp6IODCFHRFHcOGCuVLewoWwZ8/N4z4+BNety8LYWN75+WeC/5uDnSlTJjp37kzv3r154oknLCo6cc6eNTfmnDkTwsPNY08+CcOGQevW2tJBRB6dwo5IWhUeDitXmgFn/fqbSwK7uRFbty47ixdn7N69/PTtt/EveeKJJ+jbty8dOnQgUxrv7zl2zNzSYd68+KV+ePppc7Xjxo0hDTdCiYiDUdgRSUvsdvj1VzPg+PtDaOjN5ypX5lL9+sy4epXPly7lwo8/AmCz2ahbty59+/albt26abqrCuDAAXNLh6VLzWFHYK5yPHw4vPSStnQQkeSnsCOSFuzff3MczrlzN48XLkxs27ZsyJOHj3/4gXWjR2P8t5Jx7ty56dy5M127dqWwAywy88cf5ho5339/81j9+mbIqVbNurpExPkp7IhY5dw5+PprsxVn//6bx7NmhRYtOFu7NjP272fOnDmcuyUA1alTh+7du/Pqq6/insYHtBiG2VA1bhxs2GAes9mgWTMz5JQvb219IpI+KOyIpKbQULN7auFC89c/br8pd3do2JCI5s355vp15ixezNaZM+NfljNnTjp16kTXrl0pWrSoRcUnnmHA6tVmyPn9d/OYmxu8/joMGWIOQBYRSS0KOyIpLSbGDDYLF5rbcsdNOQKoWpXYtm35zdeXL/39Wd61K+H/Pe/i4sKLL75Ip06daNKkCR4eHhZ9gcSz2+G778yQ89df5jFPT+jSxdzSoWBBa+sTkfRJYUckJRgGBASYAWfJEjh//uZzTzwB7drxT7VqzNm0ifkTJ/LPP//EP128eHE6duxIu3btyJcvX+rXfhd2+/23xoqKMr/qhAnmLCuAxx6DXr1gwABzUUAREaso7Igkp9OnzUHGCxeai//FyZ4dWrXiUv36LD56lCVLl7Lzvffin/bx8aFVq1Z07NiRypUrY0tDU5L8/e++6fmUKVC3rrmlw0cf3Xw+Wzbo3x/69DGHH4mIWE1hR+RRhYSY2zUsWgSbNt0ch+PpCa+8QmiTJiy7fp1F337Lr198Qex/861dXFyoU6cOHTt2pHHjxnh5eVn3He7B3x+aN7/5leKcPWsOMs6UCa5fN4/5+sLbb0O3bmarjohIWqGwI/IwoqNh3TqzBef77+HGjZvP1ahBZMuWrMmQgQWrVrGmY0eioqLin65SpQqtW7emRYsW5EnD/Tt2u9mic3vQgZvHrl+HwoVh6FDo0MHMdyIiaY3CjkhiGQbs2mUGnKVL4d9/bz735JNEtmjBupw5WbhlC2sGDSIsLCz+6VKlStGmTRtatWrF448/bkHxSbd1a8Kuq3uZPRteeCHl6xEReVgKOyIPcvKk2UW1aBEcPnzzeK5c3GjShPW+vszetYt1EycSGRkZ/3TBggVp06YNrVu3pkyZMqlf90O4dAm2b4dt2xIu/nc/cTuTi4ikVQo7Indz9ao5DmfhQrOJI46XFxH16vFrvnxM/ftvfvnqK2JiYuKfLlq0KM2aNaNZs2ZUrFgxTQ00vp1hmNlt27abAefWLJdYvr7JX5uISHJS2BGJExUFP/1kBpwffjAfA4bNRmilSmzKl4/Jp07x68qV8Vs2AJQuXTo+4JQuXTrNBpwbN+DPP28Gm+3b4fLlO88rUQKqVoUqVeCdd8wN1+82bsdmM2dlVa+e8rWLiDwKhR1J3wzDXOJ34UL45hu4ciX+qesFC/KLry/j//mHP3fuTPCyChUqxAecYsWKpXbViXLhws1gs20b7N59c3fxOBkyQKVK8NxzZsB59llzlnycrFnN2Vg2W8LAE5fnJk9OuN6OiEhapLAj6dPx4zfH4cStggeE+fiwLkcOJp49y85//oH/Fvvz8vKiTp06NGzYkJdffjnNLPYXJzYWDh26GWy2bTO/4u1y5zZDTdytfHm438LMTZuavXl3W2dn8mTzeRGRtM5mGHdroE4/QkJC8PHxITg4mMyZM1tdjiTVg5b2vdXly/Dtt2Yrzo4d8Ycj3d35ycuL6SEhbABi/zueP39+GjZsSMOGDalVqxYZM2ZM8a+TWOHh5i7iccFmxw64di3hOTYblCqVMNwULnyzVSYpknKZRURSQ1J+v9WyI47rfkv7xjU5REaaO1IuXIixZg22//px7MAGm40FhsHK6GjCoqNxdXWlcqVKNGjQgIYNG1KuXLk0M/7m3LmbwWb7dti719xy61YZM0LlyjeDTZUqkCVL8ny+qyvUrJk87yUiktoUdsQx3W9p3+bNYcwYYk+dInbpUtz+W+LXBuwFFgJLgPOGQZEiRWj/0ku8+OKL1KpViyzJlQ4egd0OBw4knCV18uSd5+XLdzPYPPcclCtnbp4uIiIJqRtL3ViOx26HQoXuueKdgRls4pwGFgOLgDM+PtSuXZuX/gs4RYoUSfFyH+T6ddi582aw+f13cweKW7m4QNmyNwcSV60KBQo8XJeUiIgzUDeWOLe1a++7tG/c7/8aYLqXF9SowfO1avFVjRpUqFABNzdr/9qfPp1wIPG+feYA41s99pjZDRUXbCpXBmVxEZGHo7AjaUtUlNkVdfo0nDoVfzNOnybq2DFsp0/jERGRqLcqPnYsK4cOtTTcxMTA/v0Jw83dclqBAgkHEpcuDRZnMhERp6F/TiX1GIa5H0FciLkt0HD6tDnd5y49qzYgqXtMFqlWLdUTQ3Cw2Q0VF2x27oRbtsgCzMG+Tz2VcLxN/vypWqaISLqisCPJJyzMDCy3hpjb79+6O/g93MAcZ3Pqv1vc/XOurmQqVYrHq1Xj3W++wevKFWwWLu1rGObA4VtnSf31151ZzcfHXKwvLthUqmR2U4mISOpQ2JHEsdvNVpf7tcrcbe+BuwjPkoVLXl78YxgcDA3lYGhoglDzL2Cz2ShWrBgVKlSgQoUKdK5ShaeffpoMGTKYb/LCC6m+tG90tDnl+9ZZUkFBd573+OM3g03VquZaNy4uyVqKiIgkgcKOmGHh2rV7h5hTp8xxNHb7g98rUyYMPz/CcuTgXy8vTsbE8L/r1/nj/Hl2nDnDqdhYoq5dS7ACns1mo3jx4lSoUIG2/4Wb8uXLkylTpnt/Tios7Xv1qrlYX1zLzR9/wO3Dhdzc4OmnE3ZJaWNMEZG0RVPP08PU88hIMxDcr1UmNPTB7+PmZoaJAgXMQJMtG+fc3DgRHc2BkBD+vHCBgMBATpw4kWAn8FtlyZKFMmXKxN9Kly5NuXLl7h9s7ieZlvY1DHN7hVsHEh88eOd5WbMmnP79zDPg5fVwpYuIyMPT1PM0wB5l56/pWwk/HkTGIr6U6VUdV48UWF8/NhYuXrx/q8yFC4l7rxw5zGlB/90ic+fmoqcn/xgGR2/c4H+XL3Pin384efIkx3/4gZDbF4O5hbe3N8WKFaN06dIJwk3evHmTdVViO65spSZBgC9QHUjMVY6MhD17Eo63uXjxzvOeeCLhLKnixdUlJSLiaNJU2NmyZQsffvghu3fvJigoiBUrVtC4ceP7vmbTpk0MHDiQ//3vf/j5+fHOO+/QsWPHVKn3Xn4f7E+BT/rxlP1m98q5t/NzauAUqkxKYvfK9ev3HvB76pTZYhMV9eD38fK6GWT8/IjMk4erjz3GeXd3/jEMjkdG8s+//3L27FlOnz7Nye3buXi3X/9buLi4UKhQIYoVK0bx4sXjb8WKFSNfvnwpvtVCYnaLiHPpUsIdwHftMgPPrTw8oGLFm8Hm2WchV64U/QoiIpIK0lTYCQsLo1y5cnTu3JmmiRhzERgYyMsvv0yPHj1YvHgxGzZsoEuXLvj6+lK3bt1UqPhOvw/2p9KHzTHX8b0pj/0seT5szu8suxl4oqPNTY/u1ypz++6Od+PigpEnD9G+voTnyEGwjw+XM2YkyN2dU8CJ6GhOXr/OpcuXuXDmDGd37rxvq8ytfHx8KFy4cPytUKFC8feLFi2Kp2dSJ4QnjwftFjF5Mnh732y1OXz4zvfIkSPhQOIKFSBu/LOIiDiPNDtmx2azPbBlZ8iQIfz4448cOHAg/lirVq24du0aa9euTdTnJOeYHXuUnQsZC5HHfoa79XQYQARe3Cj+BF6XzpPh6iVsty+dexfhnl5cypiZfz0zcs4tA2dc3DkZ68KJaIMjN2I4GnaDiJjoJNebMaM3efLkIXfu3OTOnZs8efLEP/bz88PPzw8fH58kv29Ks9uhWrW7z4S6nxIlEo63eeIJbbcgIuKo0s2YnR07dlCnTp0Ex+rWrUv//v3v+ZrIyEgib+m/SGwLR2L8NX1rgq6r29mAjESQ8fD+m/XgwRnyc4oCnKIAp/GLvx/3ODQyE0Te820fWng4nDhh3pxRmTLw8ss3u6SyZ7e6IhERsYJDh53z58+TO3fuBMdy585NSEgIEREReN1lmsz48eMZPXp0itQTfjxxTQ2f0ocltOA0flwkF4Ytrh3IdktLw837rkTFj3950J8Jt8B0Tna72QP4IMOGQevWKV+PiIikbQ4ddh7GsGHDGDhwYPzjkJAQ/Pz8kuW9MxZJ3AIrtT5txoD+Kbu6rzPbtAlq1XrweVrvRkREgLsOLXEYefLk4cJt06ovXLhA5syZ79qqA+Dp6UnmzJkT3JJLmV7VOeean9h7tK7EYuOsqx9leinoPIrq1c1ZV/cab2OzgZ9fiu8WISIiDsKhw86zzz7Lhg0bEhxbv349zz77rCX1uHq4cmrgFIA7Ak/c49MDJ6fMejvpiKurOb0c7gw8KbhbhIiIOKg0FXZCQ0MJCAggICAAMKeWBwQEcOrUKcDsgmrfvn38+T169ODEiRMMHjyYv//+m+nTp/Ptt98yYMAAK8oHoMqkpvwxaBnnXfMlOB7kmp8/Bi1L+jo7cldxu0XkS3iZyZ/fPJ4Mu0WIiIiTSFNTzzdt2kStuwzG6NChA/PmzaNjx46cPHmSTZs2JXjNgAEDOHjwIPnz5+fdd99N0qKCKbVdRKqtoJzOJdNuESIi4mCS8vudpsKOFdLF3lgiIiJOJim/32mqG0tEREQkuSnsiIiIiFNT2BERERGnprAjIiIiTk1hR0RERJyawo6IiIg4NYUdERERcWoKOyIiIuLUFHZERETEqblZXYDV4haQDgkJsbgSERERSay43+3EbASR7sPO9evXAfDz87O4EhEREUmq69ev4+Pjc99z0v3eWLGxsZw7d45MmTJhs9mS9b1DQkLw8/Pj9OnT2ncrBek6pw5d59Sh65x6dK1TR0pdZ8MwuH79Onnz5sXF5f6jctJ9y46Liwv58+dP0c/InDmz/kNKBbrOqUPXOXXoOqceXevUkRLX+UEtOnE0QFlEREScmsKOiIiIODWFnRTk6enJyJEj8fT0tLoUp6brnDp0nVOHrnPq0bVOHWnhOqf7AcoiIiLi3NSyIyIiIk5NYUdEREScmsKOiIiIODWFHREREXFqCjuPaNq0aRQqVIgMGTJQuXJl/vjjj/ue/9133/Hkk0+SIUMGypQpw5o1a1KpUseWlOs8e/ZsqlevTtasWcmaNSt16tR54P8uYkrq3+c4S5cuxWaz0bhx45Qt0Ekk9Tpfu3aN3r174+vri6enJ8WKFdO/HYmU1Gs9efJkihcvjpeXF35+fgwYMIAbN26kUrWOZ8uWLTRq1Ii8efNis9lYuXLlA1+zadMmnn76aTw9PSlatCjz5s1L8Tox5KEtXbrU8PDwMObMmWP873//M7p27WpkyZLFuHDhwl3P37Ztm+Hq6mpMmjTJOHjwoPHOO+8Y7u7uxl9//ZXKlTuWpF7nNm3aGNOmTTP27t1rHDp0yOjYsaPh4+NjnDlzJpUrdyxJvc5xAgMDjXz58hnVq1c3Xn311dQp1oEl9TpHRkYaFStWNBo0aGD89ttvRmBgoLFp0yYjICAglSt3PEm91osXLzY8PT2NxYsXG4GBgcbPP/9s+Pr6GgMGDEjlyh3HmjVrjBEjRhj+/v4GYKxYseK+5584ccLImDGjMXDgQOPgwYPG1KlTDVdXV2Pt2rUpWqfCziOoVKmS0bt37/jHdrvdyJs3rzF+/Pi7nt+iRQvj5ZdfTnCscuXKRvfu3VO0TkeX1Ot8u5iYGCNTpkzG/PnzU6pEp/Aw1zkmJsZ47rnnjC+//NLo0KGDwk4iJPU6f/HFF8bjjz9uREVFpVaJTiOp17p3795G7dq1ExwbOHCgUbVq1RSt01kkJuwMHjzYKFWqVIJjLVu2NOrWrZuClRmGurEeUlRUFLt376ZOnTrxx1xcXKhTpw47duy462t27NiR4HyAunXr3vN8ebjrfLvw8HCio6PJli1bSpXp8B72Oo8ZM4ZcuXLxxhtvpEaZDu9hrvOqVat49tln6d27N7lz56Z06dKMGzcOu92eWmU7pIe51s899xy7d++O7+o6ceIEa9asoUGDBqlSc3pg1e9gut8I9GFdunQJu91O7ty5ExzPnTs3f//9911fc/78+buef/78+RSr09E9zHW+3ZAhQ8ibN+8d/4HJTQ9znX/77Te++uorAgICUqFC5/Aw1/nEiRNs3LiRtm3bsmbNGo4dO0avXr2Ijo5m5MiRqVG2Q3qYa92mTRsuXbpEtWrVMAyDmJgYevTowfDhw1Oj5HThXr+DISEhRERE4OXllSKfq5YdcWoTJkxg6dKlrFixggwZMlhdjtO4fv067dq1Y/bs2eTIkcPqcpxabGwsuXLlYtasWVSoUIGWLVsyYsQIZsyYYXVpTmfTpk2MGzeO6dOns2fPHvz9/fnxxx8ZO3as1aXJI1LLzkPKkSMHrq6uXLhwIcHxCxcukCdPnru+Jk+ePEk6Xx7uOsf56KOPmDBhAr/88gtly5ZNyTIdXlKv8/Hjxzl58iSNGjWKPxYbGwuAm5sbhw8fpkiRIilbtAN6mL/Pvr6+uLu74+rqGn+sRIkSnD9/nqioKDw8PFK0Zkf1MNf63XffpV27dnTp0gWAMmXKEBYWRrdu3RgxYgQuLmofeFT3+h3MnDlzirXqgFp2HpqHhwcVKlRgw4YN8cdiY2PZsGEDzz777F1f8+yzzyY4H2D9+vX3PF8e7joDTJo0ibFjx7J27VoqVqyYGqU6tKRe5yeffJK//vqLgICA+Nsrr7xCrVq1CAgIwM/PLzXLdxgP8/e5atWqHDt2LD5MAhw5cgRfX18Fnft4mGsdHh5+R6CJC5mGtpFMFpb9Dqbo8Gcnt3TpUsPT09OYN2+ecfDgQaNbt25GlixZjPPnzxuGYRjt2rUzhg4dGn/+tm3bDDc3N+Ojjz4yDh06ZIwcOVJTzxMhqdd5woQJhoeHh7Fs2TIjKCgo/nb9+nWrvoJDSOp1vp1mYyVOUq/zqVOnjEyZMhl9+vQxDh8+bKxevdrIlSuX8f7771v1FRxGUq/1yJEjjUyZMhlLliwxTpw4Yaxbt84oUqSI0aJFC6u+Qpp3/fp1Y+/evcbevXsNwPjkk0+MvXv3Gv/8849hGIYxdOhQo127dvHnx009HzRokHHo0CFj2rRpmnruCKZOnWoUKFDA8PDwMCpVqmT8/vvv8c/VqFHD6NChQ4Lzv/32W6NYsWKGh4eHUapUKePHH39M5YodU1Kuc8GCBQ3gjtvIkSNTv3AHk9S/z7dS2Em8pF7n7du3G5UrVzY8PT2Nxx9/3Pjggw+MmJiYVK7aMSXlWkdHRxujRo0yihQpYmTIkMHw8/MzevXqZVy9ejX1C3cQv/76613/vY27rh06dDBq1Khxx2ueeuopw8PDw3j88ceNuXPnpnidNsNQ25yIiIg4L43ZEREREaemsCMiIiJOTWFHREREnJrCjoiIiDg1hR0RERFxago7IiIi4tQUdkRERMSpKeyIiIiIU1PYEREREaemsCMiTsswDMaMGcPWrVutLkVELKSwIyJO68iRI4wcOZKgoCCrSxERCynsiIjT2r17NwBPP/20xZWIiJW0EaiIOKVKlSrx559/Jjjm4+PDtWvXrClIRCzjZnUBIiIpYciQIYwaNYrIyEjee+89ALJkyWJtUSJiCbXsiIjTKliwILVr12bu3LlWlyIiFtKYHRFxSsHBwZw6dYqyZctaXYqIWExhR0Sc0v79+wEUdkREYUdEnFNc2ClXrpzFlYiI1RR2RMQp7d+/H19fX3LkyGF1KSJiMYUdEXFKp06dIn/+/FaXISJpgKaei4hTKly4MBs3bmTSpEnkzZuXEiVKUKFCBavLEhELaOq5iDilc+fO0blzZ7Zt20ZoaCifffYZb775ptVliYgFFHZERETEqWnMjoiIiDg1hR0RERFxago7IiIi4tQUdkRERMSpKeyIiIiIU1PYEREREaemsCMiIiJOTWFHREREnJrCjoiIiDg1hR0RERFxago7IiIi4tQUdkRERMSp/R8ZITJVpzRiWQAAAABJRU5ErkJggg==\n"
},
"metadata": {}
}
],
"source": [
"def solveIVP(f, tspan, y0, h, solver):\n",
"\n",
" # Initialise t and y arrays\n",
" t = np.arange(tspan[0], tspan[1] + h, h)\n",
" y = np.zeros((len(t),len(y0)))\n",
" y[0,:] = y0\n",
"\n",
" # Loop through steps and calculate single step solver solution\n",
" for n in range(len(t) - 1):\n",
" y[n+1,:] = solver(f, t[n], y[n,:], h)\n",
"\n",
" return t, y\n",
"\n",
"\n",
"def euler(f, t, y, h):\n",
" return y + h * f(t, y)\n",
"\n",
"def rk4(f, t, y, h):\n",
" k1 = f(t, y)\n",
" k2 = f(t + 0.5 * h, y + 0.5 * h * k1)\n",
" k3 = f(t + 0.5 * h, y + 0.5 * h * k2)\n",
" k4 = f(t + h, y + h * k3)\n",
" return y + h / 6 * (k1 + 2 * k2 + 2 * k3 + k4)\n",
"\n",
"# Define the ODE function\n",
"def f(t, y):\n",
" return t * y\n",
"\n",
"def exact(t):\n",
" return np.exp(t**2/2)\n",
"\n",
"\n",
"# Define IVP parameters\n",
"tspan = [0, 1] # boundaries of the t domain\n",
"y0 = [1] # initial values\n",
"h = 0.2 # step length\n",
"\n",
"# Calculate the solution to the IVP\n",
"t, yEuler = solveIVP(f, tspan, y0, h, euler)\n",
"t, yRK4 = solveIVP(f, tspan, y0, h, rk4)\n",
"\n",
"# Print table of solution values\n",
"print(\"| t | yEuler | RK4 | Exact |\")\n",
"print(\"|------|-----------|-----------|-----------|\")\n",
"for n in range(len(t)):\n",
" print(f\"| {t[n]:4.2f} | {yEuler[n,0]:9.6f} | {yRK4[n,0]:9.6f} | {exact(t[n]):9.6} |\" )\n",
"\n",
"#calculate exact solution\n",
"tExact=np.linspace(tspan[0], tspan[1], 200)\n",
"yExact=exact(tExact)\n",
"# Plot solution\n",
"fig, ax = plt.subplots()\n",
"plt.plot(tExact, yExact, \"k-\", label=\"Exact Solution\")\n",
"plt.plot(t, yEuler, \"bo-\", label=\"Euler\")\n",
"plt.plot(t, yRK4, \"r-o\", label=\"RK4 method\")\n",
"plt.xlabel(\"$t$\", fontsize=12)\n",
"plt.ylabel(\"$y$\", fontsize=12)\n",
"plt.legend()\n",
"plt.show()"
]
},
{
"cell_type": "code",
"source": [
"# Define ODE function and the exact solution\n",
"def f(t,y):\n",
" return t * y\n",
"\n",
"\n",
"def exact(t):\n",
" return np.exp(t ** 2 / 2)\n",
"\n",
"\n",
"# Define IVP parameters\n",
"tspan = [0, 1]\n",
"y0 = [1]\n",
"hvals = [ 0.2, 0.1, 0.05, 0.025 ] # step length values\n",
"\n",
"# Loop through h values and calculate errors\n",
"errors = []\n",
"errorsRk4 = []\n",
"print(f\"Exact solution: y(1) = {exact(1):0.6f}\\n\")\n",
"print(\"| h | Euler | Error | RK4 | Error |\")\n",
"print(\"|:-----:|:--------:|:--------:|:-------:|:--------:|\")\n",
"for h in hvals:\n",
" t, y = solveIVP(f, tspan, y0, h, euler)\n",
" t, yrk4 = solveIVP(f, tspan, y0, h, rk4)\n",
" errors.append(abs(y[-1,0] - exact(t[-1])))\n",
" errorsRk4.append(abs(yrk4[-1,0] - exact(t[-1])))\n",
" print(f\"| {h:0.3f} | {y[-1,0]:0.6f} | {errors[-1]:0.2e} | {yrk4[-1,0]:0.6f} | {errorsRk4[-1]:0.2e} |\")\n",
"\n",
"# Plot errors\n",
"fig, ax = plt.subplots()\n",
"plt.plot(hvals, errors, \"bo-\")\n",
"plt.plot(hvals, errorsRk4, \"ro-\")\n",
"plt.xlabel(\"$h$\", fontsize=12)\n",
"plt.ylabel(\"$E(h)$\", fontsize=12)\n",
"plt.show()"
],
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/",
"height": 599
},
"id": "ix0mHSwbsKmB",
"outputId": "e1769dc3-87e5-42e3-f2d7-3af78063b7c6"
},
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"name": "stdout",
"text": [
"Exact solution: y(1) = 1.648721\n",
"\n",
"| h | Euler | Error | RK4 | Error |\n",
"|:-----:|:--------:|:--------:|:-------:|:--------:|\n",
"| 0.200 | 1.459261 | 1.89e-01 | 1.648717 | 4.59e-06 |\n",
"| 0.100 | 1.547110 | 1.02e-01 | 1.648721 | 2.64e-07 |\n",
"| 0.050 | 1.595942 | 5.28e-02 | 1.648721 | 1.55e-08 |\n",
"| 0.025 | 1.621801 | 2.69e-02 | 1.648721 | 9.33e-10 |\n"
]
},
{
"output_type": "display_data",
"data": {
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
],
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAk4AAAG0CAYAAADEuKgXAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABR60lEQVR4nO3deXxU1f3/8dckkEQKSWTLAoGAIIusssRYEMR8BVwpsUWwEpCCVhYhLoClYLVtKFCkCkJtBWwVQRSxooVKAAsYFhMg7AUKApKExV8SZEkgOb8/jkwdCDAJk8wkeT8fj3mYe+fMnc+dDJm355w512GMMYiIiIjIdfl5uwARERGR8kLBSURERMRNCk4iIiIiblJwEhEREXGTgpOIiIiImxScRERERNyk4CQiIiLipireLqCiKSws5NixY9SoUQOHw+HtckRERMQNxhhOnz5NZGQkfn5X71dScPKwY8eOERUV5e0yREREpASOHDlC/fr1r3q/gpOH1ahRA7AvfHBwsJerEREREXfk5uYSFRXl/By/GgUnD7s0PBccHKzgJCIiUs5cb5qNJoeLiIiIuEnBSURERMRNCk4iIiIiblJwEhEREXGTgpOIiIiImxScRERERNyk4CQiIiLiJgUnERERETcpOImIiIi4SSuHi4iIiM8rKIC1ayEjAyIioGtX8Pcv+zoUnERERMSnLVkCzzwDR4/+b1/9+vCnP0HfvmVbi4bqRERExGctWQKPPOIamgC++cbuX7KkbOtRcBIRERGfVFBge5qMufK+S/tGj7btyoqCk4iIiPiktWuv7Gn6IWPgyBHbrqwoOImIiIhPysjwbDtPUHASERERn3PxIqxf717biIjSreWH9K06ERER8Snr18PTT0N6+rXbORz223Vdu5ZNXaAeJxEREfERx4/DoEHQpYsNTTVrwlNP2YDkcLi2vbQ9Y0bZruek4CQiIiJeVVAAb7wBzZrB22/bfUOHwt69MHs2fPAB1Kvn+pj69e3+sl7HSUN1IiIi4jWbNtlhudRUu92+vQ1LMTH/a9O3Lzz8sFYOFxERkUrq1CkYPx7++le7rEBICPzud3ZorqhA5O8P3buXeZlXUHASERGRMlNYCHPnwrhxNjwBJCTAH/4AYWHerc0dCk4iIiJSJrZsscNyGzbY7datYdassv1W3I3S5HAREREpVdnZMHIkdOxoQ1ONGjB9up3XVJ5CE6jHSUREREqJMfD3v8Pzz9ulBgD694dp0yAy0ru1lZSCk4iIiHjc9u12WG7dOrvdvLkdluvRw7t13SgN1YmIiIjH5OZCYqJdVmDdOqhWDSZPhm3byn9oAvU4iYiIiAcYA4sW2dB06aK78fF2LlODBt6tzZMUnEREROSG7N4NI0bAqlV2u0kTeP116NXLu3WVBg3ViYiISImcOWPXY2rb1oamoCB4+WU7v6kihiZQj5OIiIgUkzHw0UcwejQcOWL3PfAAvPYaNGrk1dJKnYKTiIiIuG3fPrsm04oVdjs62gamBx/0alllRkN1IiIicl3nzsHEidCqlQ1NAQEwYQLs3Fl5QhOox0lERESuY9kyGDUKDh602/feCzNnQtOm3q3LG9TjJCIiIkU6eBAeesj2KB08CPXrwwcfwPLllTM0gYKTiIiIXCYvD377W2jZEj75BKpUgRdesMsOxMeDw+HtCr3HZ4PTrFmziI6OJigoiJiYGDZt2nTVtjt37iQ+Pp7o6GgcDgczZsy4os2l+y6/DR8+3Nmme/fuV9z/1FNPlcbpiYiI+KR//Qtat4Zf/xrOn4e774b0dPjDH6B6dW9X530+GZwWLVpEYmIikyZNIi0tjbZt29KzZ0+OX7pC4GXOnj1L48aNmTx5MuHh4UW22bx5MxkZGc7b559/DsBPf/pTl3ZDhw51aTdlyhTPnpyIiIgPOnIEHnkEeva035yLiIAFCyA5GVq08HZ1vsMng9P06dMZOnQogwcPpmXLlsyZM4dq1aoxd+7cItt36tSJqVOn8uijjxIYGFhkmzp16hAeHu68LVu2jFtuuYVu3bq5tKtWrZpLu+DgYI+fn4iIiK/Iz4cpU2w4+vBD8Pe36zPt2QP9+1fuYbmi+Fxwys/PJzU1lbi4OOc+Pz8/4uLiSElJ8dhzvPPOOzzxxBM4LntHvPvuu9SuXZtWrVoxfvx4zp49e81j5eXlkZub63ITEREpD1avhnbtYOxYuwp4ly6QlgavvgrqNyiazy1HcPLkSQoKCggLC3PZHxYWxp49ezzyHEuXLiU7O5tBgwa57B8wYAANGzYkMjKS9PR0xo4dy969e1myZMlVj5WUlMRvfvMbj9QlIiJSFjIy4Nln4b337HadOjB1KgwcqB6m6/G54FQW3nrrLXr37k1kZKTL/mHDhjl/bt26NREREdxzzz0cOHCAW265pchjjR8/nsTEROd2bm4uUVFRpVO4iIjIDbh40a6/NHEinD5tQ9Ivf2m/QXfzzd6urnzwueBUu3Zt/P39ycrKctmflZV11YnfxfH111+zcuXKa/YiXRITEwPA/v37rxqcAgMDrzqvSkRExFesXw9PP22/IQfQuTO88QZ06ODdusobn5vjFBAQQIcOHUhOTnbuKywsJDk5mdjY2Bs+/rx586hbty7333//ddtu3boVgIiIiBt+XhEREW84fhwGDbLzl9LToWZNePNNSElRaCoJn+txAkhMTCQhIYGOHTvSuXNnZsyYwZkzZxg8eDAAAwcOpF69eiQlJQF2sveuXbucP3/zzTds3bqV6tWr06RJE+dxCwsLmTdvHgkJCVSp4nrqBw4cYMGCBdx3333UqlWL9PR0xowZw1133UWbNm3K6MxFREQ8o6AA/vxn+NWvIDvb7hs6FH7/e6hd26ullWs+GZz69evHiRMnmDhxIpmZmbRr147ly5c7J4wfPnwYP7//dZYdO3aM9u3bO7enTZvGtGnT6NatG2vWrHHuX7lyJYcPH+aJJ5644jkDAgJYuXKlM6RFRUURHx/PhAkTSu9ERURESsGmTXZYLjXVbrdvD7Nnw/czUOQGOIwxxttFVCS5ubmEhISQk5OjNaBERKRMnToF48fDX/8KxkBICPzud/DUU3Z9Jrk6dz+/fbLHSURERNxXWAhz58K4cTY8ASQk2MukXLa6j9wgBScREZFyLC3NDstt3Gi3W7eGWbOga1fv1lVR+dy36kREROT6srNhxAjo1MmGpho1YPp0O69Joan0qMdJRESkHDEG/v53eP55u9QA2GvKTZsGl63rLKVAwUlERKSc2L7dDsutW2e3mze3w3I9eni3rspEQ3UiIiI+LjcXEhPtsgLr1kG1ajB5MmzbptBU1tTjJCIi4qOMgUWLbGjKyLD74uPtXKYGDbxbW2Wl4CQiIuKDdu+2k79XrbLbTZrA669Dr17erauy01CdiIiIDzlzxq7H1LatDU1BQfDyy3Z+k0KT96nHSURExAcYAx99BKNHw5Ejdt8DD8Brr0GjRl4tTX5AwUlERMTL9u2DkSNhxQq7HR1tA9ODD3q1LCmChupERES85Nw5mDgRWrWyoSkgACZMgJ07FZp8lXqcREREvGDZMhg1Cg4etNv33gszZ0LTpt6tS65NPU4iIiJl6OBBeOgh26N08CDUrw8ffADLlys0lQcKTiIiImUgLw9++1to2RI++QSqVIEXXrDLDsTHg8Ph7QrFHRqqExERKWUrVtg1mfbvt9vdu9tLpbRs6dWypATU4yQiIlJKjhyBRx6x6y/t3w8REbBggV2fSaGpfFJwEhER8bD8fJgyBVq0gA8/BH9/uz7Tnj3Qv7+G5cozDdWJiIh40OrVMHy4nbsE0KWLHZZr08a7dYlnqMdJRETEAzIyYMAA6NHDhqY6dWD+fPj3vxWaKhIFJxERkRtw8SLMmAHNmsF779lhuKefhr17ISFBw3IVjYbqRERESmj9ehuS0tPtdufO8MYb0KGDd+uS0qMeJxERkWI6fhwGDbLzl9LToWZNePNNSElRaKroFJxERETcVFBge5SaNYO337b7hg61w3JDh4KfPlUrPA3ViYiIuGHTJjssl5pqt9u3h9mzISbGu3VJ2VI2FhERuYZTp2DYMLjjDhuaQkLsxXg3b1ZoqozU4yQiIlKEwkKYOxfGjbPhCey35P7wBwgL825t4j0KTiIiIpdJS7PDchs32u3Wre0ill27ercu8T4N1YmIiHwvO9tejLdTJxuaqleH6dPtEJ1Ck4B6nERERDAG/v53eP55u9QA2GvKTZsGkZHerU18i4KTiIhUatu322G5devsdvPmdliuRw/v1iW+SUN1IiJSKeXmQmKiXVZg3TqoVg0mT4Zt2xSa5OrU4yQiIpWKMbBokQ1NGRl2X3y8ncvUoIF3axPf57M9TrNmzSI6OpqgoCBiYmLYtGnTVdvu3LmT+Ph4oqOjcTgczJgx44o2L730Eg6Hw+XWvHlzlzbnz59n+PDh1KpVi+rVqxMfH09WVpanT01ERLxk926Ii7PzlzIyoEkT+Oc/4YMPFJrEPT4ZnBYtWkRiYiKTJk0iLS2Ntm3b0rNnT45fmrF3mbNnz9K4cWMmT55MeHj4VY972223kZGR4bytuzSg/b0xY8bwySefsHjxYr744guOHTtG3759PXpuIiJS9s6csesxtW0Lq1ZBUBC8/LKd39Srl7erk/LEJ4fqpk+fztChQxk8eDAAc+bM4dNPP2Xu3LmMGzfuivadOnWiU6dOAEXef0mVKlWuGqxycnJ46623WLBgAT2+H9yeN28eLVq0YMOGDdxxxx03eloiIlLGjIGPPoLRo+HIEbvvgQfgtdegUSOvlibllM/1OOXn55OamkpcXJxzn5+fH3FxcaSkpNzQsfft20dkZCSNGzfmscce4/Dhw877UlNTuXDhgsvzNm/enAYNGlzzefPy8sjNzXW5iYiI9+3bB7172/lLR45AdDT84x/wyScKTVJyPhecTp48SUFBAWGXrWcfFhZGZmZmiY8bExPD/PnzWb58ObNnz+bgwYN07dqV06dPA5CZmUlAQAChoaHFet6kpCRCQkKct6ioqBLXKCIiN+7cOZg4EVq1ghUrICAAJkyAnTvhwQe9XZ2Udz45VFcaevfu7fy5TZs2xMTE0LBhQ95//32GDBlS4uOOHz+exMRE53Zubq7Ck4iIl3zyCYwaBYcO2e1777UX5G3a1KtlSQXic8Gpdu3a+Pv7X/FttqysrGtO/C6u0NBQbr31Vvbv3w9AeHg4+fn5ZGdnu/Q6Xe95AwMDCQwM9FhdIiJSfAcPwjPP2OAEUL8+zJgBffuCw+HV0qSC8bmhuoCAADp06EBycrJzX2FhIcnJycTGxnrseb777jsOHDhAREQEAB06dKBq1aouz7t3714OHz7s0ecVERHPycuD3/4WWra0oalKFXjhBbvsQHy8QpN4ns/1OAEkJiaSkJBAx44d6dy5MzNmzODMmTPOb9kNHDiQevXqkZSUBNgJ5bt27XL+/M0337B161aqV69OkyZNAHjuued48MEHadiwIceOHWPSpEn4+/vTv39/AEJCQhgyZAiJiYnUrFmT4OBgRo4cSWxsrL5RJyLig1assBfk/X7ggO7d7aVSWrb0allSwflkcOrXrx8nTpxg4sSJZGZm0q5dO5YvX+6cMH748GH8/P7XWXbs2DHat2/v3J42bRrTpk2jW7durFmzBoCjR4/Sv39/Tp06RZ06dejSpQsbNmygTp06zse9+uqr+Pn5ER8fT15eHj179uSNN94om5MWERG3HDkCY8bAhx/a7YgI+OMf4dFH1cMkpc9hjDHeLqIiyc3NJSQkhJycHIKDg71djohIhZGfb+ctvfyyXdDS3x9GjoTf/Ab051ZulLuf3z7Z4yQiIvJDq1fD8OF27hJAly52WK5NG+/WJZWPz00OFxERuSQjAwYMgB49bGiqUwfmz4d//1uhSbxDwUlERHzOxYt2WK5ZM3jvPTt36emnYe9eSEjQXCbxHg3ViYiIT1m/3oak9HS73bkzvPEGdOjg3bpEQD1OIiLiI44fh0GD7Pyl9HSoWRPefBNSUhSaxHcoOImIiFcVFNgepWbN4O237b6hQ+2w3NCh4KdPKvEhGqoTERGv2bjRDsulpdnt9u1h9myIifFuXSJXoxwvIiJl7tQpGDYMYmNtaAoJsRfj3bxZoUl8m3qcRESkzBQWwty5MG6cDU9gvyX3hz/A9xeHEPFpCk4iIlIm0tLssNzGjXa7dWu7iGXXrt6tS6Q4NFQnIiKlKjvbXoy3UycbmqpXh+nTITVVoUnKH/U4iYhIqTAG/v53eP55u9QAQP/+MG0aREZ6tzaRklJwEhERj9u+3Q7LrVtnt5s3t8NyPXp4ty6RG6WhOhER8ZjcXEhMtMsKrFsH1arB5MmwbZtCk1QM6nESEZEbZgwsWmRDU0aG3Rcfb+cyNWjg3dpEPEnBSUREbsju3Xby96pVdrtJE3j9dejVy7t1iZQGDdWJiEiJnDlj12Nq29aGpqAgePllO79JoUkqKvU4iYhIsRgDS5bAmDFw5Ijd98AD8Npr0KiRd2sTKW0KTiIi4rZ9+2DkSFixwm5HR9vA9OCDXi1LpMxoqE5ERK7r3DmYOBFatbKhKSAAJkyAnTsVmqRyUY+TiIhc0yefwKhRcOiQ3b73XntB3qZNvVqWiFeox0lERIp08CA89JC9HToE9evDBx/A8uUKTVJ5KTiJiIiLvDz47W+hZUvb21SlCrzwgl12ID4eHA5vVyjiPRqqExERpxUr7JpM+/fb7e7d7aVSWrb0alkiPkM9TiIiwpEj8Mgjdv2l/fshIgIWLLDrMyk0ifyPgpOISCWWnw9TpkCLFvDhh+DvD6NHw5490L+/huVELqehOhGRSmr1ahg+3M5dAujSxQ7LtWnj3bpEfJl6nEREKpmMDBgwAHr0sKGpTh2YPx/+/W+FJpHrUXASEakkLl6EGTOgWTN47z07DPf007B3LyQkaFhOxB0aqhMRqQTWrbMhaft2u925M7zxBnTo4N26RMob9TiJiFRgx4/DoEHQtasNTTVrwptvQkqKQpNISSg4iYhUQAUFtkepWTN4+2277xe/sMNyQ4eCn/76i5SIhupERCqYjRvtsFxamt1u396GqDvu8G5dIhWB/p9DRKSCOHUKhg2D2FgbmkJC7MV4N29WaBLxFJ8NTrNmzSI6OpqgoCBiYmLYtGnTVdvu3LmT+Ph4oqOjcTgczJgx44o2SUlJdOrUiRo1alC3bl369OnD3r17Xdp0794dh8Phcnvqqac8fWoiIh5VWAh//asdlvvLX8AY+y25vXvtOk3+/t6uUKTi8MngtGjRIhITE5k0aRJpaWm0bduWnj17cvz48SLbnz17lsaNGzN58mTCw8OLbPPFF18wfPhwNmzYwOeff86FCxe49957OXPmjEu7oUOHkpGR4bxNmTLF4+cnIuIpaWlw55123tKpU9C6tV2Paf58CAvzdnUiFY/DGGO8XcTlYmJi6NSpEzNnzgSgsLCQqKgoRo4cybhx46752OjoaEaPHs3o0aOv2e7EiRPUrVuXL774grvuuguwPU7t2rUrssfqavLy8sjLy3Nu5+bmEhUVRU5ODsHBwW4fR0SkOLKzYcIEmD3b9jhVrw4vv2wv0Fu1qrerEyl/cnNzCQkJue7nt8/1OOXn55OamkpcXJxzn5+fH3FxcaSkpHjseXJycgCoWbOmy/53332X2rVr06pVK8aPH8/Zs2eveZykpCRCQkKct6ioKI/VKCJyOWPgb3+zw3KzZtnQ1L+/HZYbM0ahSaS0+dy36k6ePElBQQFhl/Uxh4WFsWfPHo88R2FhIaNHj+bHP/4xrVq1cu4fMGAADRs2JDIykvT0dMaOHcvevXtZsmTJVY81fvx4EhMTnduXepxERDxt+3b7bbl16+x28+Y2PPXo4d26RCoTnwtOZWH48OHs2LGDdZf++nxv2LBhzp9bt25NREQE99xzDwcOHOCWW24p8liBgYEEBgaWar0iUrnl5sJLL8Frr9n1mapVg4kTbQ9TQIC3qxOpXHxuqK527dr4+/uTlZXlsj8rK+uqE7+LY8SIESxbtozVq1dTv379a7aNiYkBYP/+/Tf8vCIixWUMLFxoe5ZefdWGpvh4e2HesWMVmkS8weeCU0BAAB06dCA5Odm5r7CwkOTkZGJjY0t8XGMMI0aM4KOPPmLVqlU0atTouo/ZunUrABERESV+XhGRkti9G+Li7PyljAxo0gT++U/44ANo0MDb1YlUXj45VJeYmEhCQgIdO3akc+fOzJgxgzNnzjB48GAABg4cSL169UhKSgLshPJdu3Y5f/7mm2/YunUr1atXp0mTJoAdnluwYAEff/wxNWrUIDMzE4CQkBBuuukmDhw4wIIFC7jvvvuoVasW6enpjBkzhrvuuos2bdp44VUQkcrozBl45RWYPh0uXICgIHjxRXj+efuziHiZ8VGvv/66adCggQkICDCdO3c2GzZscN7XrVs3k5CQ4Nw+ePCgAa64devWzdmmqPsBM2/ePGOMMYcPHzZ33XWXqVmzpgkMDDRNmjQxzz//vMnJySlW3Tk5OQYo9uNEpHIrLDTmgw+MiYoyxg7SGfPAA8b897/erkykcnD389sn13Eqz9xdB0JE5JJ9+2DkSFixwm5HR9uJ4A8+6NWyRCqVcruOk4hIZXHunP12XKtWNjQFBNhFLXfuVGgS8VU+OcdJRKSi++QTGDUKDh2y2/feay/I27SpV8sSketQj5OISBk6eBAeesjeDh2C+vXtN+WWL1doEikPFJxERMpAXh789rfQsqXtbapSBV54wS47EB8PDoe3KxQRd2ioTkSklK1YYS++e2kt3e7d7aVSWrb0alkiUgLqcRIRKSVHjsAjj0CvXjY0RUTAggWwapVCk0h5peAkIuJh+fkwZQq0aAEffgj+/jB6NOzZY1cC17CcSPmloToREQ9avRqGD7dzlwC6dLHDcroAgUjFoB4nEREPyMiAAQOgRw8bmurUgfnz4d//VmgSqUgUnEREbsDFizBjBjRrBu+9Z4fhnn4a9u6FhAQNy4lUNBqqExEpoXXrbEjavt1ud+4Mb7wBHTp4ty4RKT3qcRIRKabjx2HQIOja1YammjXhzTchJUWhSaSiU3ASEXFTQYHtUWrWDN5+2+77xS/ssNzQoeCnv6giFZ6G6kRE3LBxox2WS0uz2+3b2xB1xx3erUtEypb+/0hE5BpOnYJhwyA21oamkBB7Md7NmxWaRCoj9TiJiBShsBDmzoVx42x4AvstuT/8AcLCvFubiHiPgpOIyGXS0uyw3MaNdrt1a7uIZdeu3q1LRLxPQ3UiIt/LzrYX4+3UyYam6tVh+nRITVVoEhFLPU4iUukZA3//Ozz/vF1qAOw15aZNg8hI79YmIr5FwUlEKrXt2+2w3Lp1drt5czss16OHd+sSEd+koToRqZRycyEx0S4rsG4dVKsGkyfDtm0KTSJydepxEpFKxRhYuBCefdZemBcgPt7OZWrQwLu1iYjvU3ASkUpj9247+XvVKrvdpAm8/jr06uXdukSk/NBQnYhUeGfO2PWY2ra1oSkoCF5+2c5vUmgSkeJQj5OIVFjGwJIlMGYMHDli9z3wALz2GjRq5N3aRKR8UnASkQpp3z4YORJWrLDb0dE2MD34oFfLEpFyTkN1IlKhnDsHEydCq1Y2NAUEwIQJsHOnQpOI3Dj1OIlIhfHJJzBqFBw6ZLfvvddekLdpU6+WJSIViHqcRKTcO3gQHnrI3g4dgvr14YMPYPlyhSYR8Sz1OIlIuVBQAGvX2rWXIiLsteMuXoSpU+F3v4Pz56FKFbuo5a9/ba8zJyLiaQpOIuLzliyBZ56Bo0f/t692bRuUMjPtdvfu9lIpLVt6pUQRqSQ8FpzOnz+Pw+EgMDDQU4cUEWHJEnjkEbu0wA+dPGn/GxoKb7wBjz4KDkeZlycilUyJg9OaNWv4+OOPWb9+Pbt27eLcuXMAVKtWjRYtWnDnnXfSp08funfv7qlaRaSSKSiwPU2Xh6Yf+tGP4Gc/U2gSkbJRrMnhFy5cYObMmTRu3JgePXrw97//ndDQUH7+85/zwgsv8PzzzzNgwABCQ0N555136NGjB40aNWLmzJlcuHChWIXNmjWL6OhogoKCiImJYdOmTVdtu3PnTuLj44mOjsbhcDBjxowSHfP8+fMMHz6cWrVqUb16deLj48nKyipW3SLiOWvXug7PFeWbb2w7EZGyUKwepyZNmpCfn09CQgI/+9nPuP3226/ZPjU1lcWLF/P73/+eadOmcejSd4SvY9GiRSQmJjJnzhxiYmKYMWMGPXv2ZO/evdStW/eK9mfPnqVx48b89Kc/ZcyYMSU+5pgxY/j0009ZvHgxISEhjBgxgr59+7J+/Xq36hYRz7peaLrk0sV6RURKnSmGOXPmmPPnzxfnIcYYY/Ly8sycOXPcbt+5c2czfPhw53ZBQYGJjIw0SUlJ131sw4YNzauvvlrsY2ZnZ5uqVauaxYsXO9vs3r3bACYlJcXt2nNycgxgcnJy3H6MiFxp505jbr3VGDtQd+3b6tXerlZEyjt3P7+LNVT35JNPlmjyd0BAAE8++aRbbfPz80lNTSUuLs65z8/Pj7i4OFJSUor93O4eMzU1lQsXLri0ad68OQ0aNLjm8+bl5ZGbm+tyE5GSKyiAKVOgfXv4z3+uPXfJ4YCoKLs0gYhIWfC5BTBPnjxJQUEBYWFhLvvDwsLIvPS941I4ZmZmJgEBAYSGhhbreZOSkggJCXHeoqKiSlSjiMDevdClC4wdC/n5cN998Oc/24B0eYC6tD1jBvj7l3mpIlJJ3XBw2rlzJ++88w7vv/8+27ZtK/Yk8PJu/Pjx5OTkOG9HLl2CXUTcVlAA06dDu3awYQMEB8PcubBsGQwdalcBr1fP9TGXVgfv29crJYtIJXVD6zi99tprJCYmUlhYCIDD4aBKlSrceuuttGnTxnnr3bu328esXbs2/v7+V3ybLSsri/Dw8BLV6c4xw8PDyc/PJzs726XX6XrPGxgYqLWrRG7Avn0weDBc+g7GvffCX/9qh+Au6dsXHn74ypXD1dMkImXthnqcpkyZwu23386uXbvYt28f//jHP5g0aRK33XYbaWlpTJgwgQceeKBYxwwICKBDhw4kJyc79xUWFpKcnExsbGyJ6nTnmB06dKBq1aoubfbu3cvhw4dL/LwicnWFhfD669C2rQ1N1avDm2/a68sVNeLt729XB+/f3/5XoUlEvOGGepxyc3OZNGkSzZs3B+CWW27h/vvvd95//vx5duzYUezjJiYmkpCQQMeOHencuTMzZszgzJkzDB48GICBAwdSr149kpKSADv5e9euXc6fv/nmG7Zu3Ur16tVp0qSJW8cMCQlhyJAhJCYmUrNmTYKDgxk5ciSxsbHccccdJX+RROQK//0vPPEEfPGF3e7Rww7NNWzo3bpERK7nhoJTly5dOHz48FXvDwoKomPHjsU+br9+/Thx4gQTJ04kMzOTdu3asXz5cufk7sOHD+Pn97/OsmPHjtG+fXvn9rRp05g2bRrdunVjzZo1bh0T4NVXX8XPz4/4+Hjy8vLo2bMnb7zxRrHrF5GiFRbayd7PPw9nzkC1avYivU89BX4+91UVEZErOYy51sUMXM2YMYPbbruNli1bUq9ePbZs2cKDDz7Ipk2biIyMLM06y43c3FxCQkLIyckhODjY2+WI+Iyvv4YhQ+DSaPhdd8G8edC4sXfrEhEB9z+/ixWc/H8wqSA4OJgWLVqQm5vLiRMneOWVV+jTp0+RK3tXJgpOIq6MsZO9ExPhu+/gpptg8mQYMUK9TCLiO0olOJ05c4YdO3awfft2l9upU6fswRwO6tWrR5s2bWjbtq3zW3UtWrS48TMqJxScRP7nyBH4xS/gX/+y2z/+se1latrUu3WJiFyuVILT1WRmZpKenu4Spnbv3s358+dxOBwUFBTc6FOUGwpOIraXaf58GD0acnMhKAh+9zt45hl9G05EfJO7n983NDn8kvDwcMLDw7n33nud+woLC9m3bx/bt2/3xFOISDlx7JhdtPKzz+x2TIwNUd9/+VZEpFzzSHAqip+fH82aNaNZs2al9RQi4kOMgXfegVGjIDsbAgLglVfg2WfVyyQiFUexpma2bNmSv/3tb+Tn57v9mLy8PObNm0fLli2LXZyIlA+ZmdCnDwwcaENTx46QlgYvvKDQJCIVS7F6nAYNGkRiYiLPPPMMDz30EHFxcdx+++00atSIatWqAXYC+cGDB/nqq69YuXIln3zyCQEBATz//POlcgIi4j3GwMKF9hty334LVavCSy/ZwFSl1PqzRUS8p9iTw0+fPs1bb73F/PnzSU9Px/H9JcqrfP9X8uLFiwAYY2jVqhVPPPEETzzxRKWZKK3J4VJZHD8Ov/wlLFlit9u3t3OZ2rTxalkiIiVSJt+qO3ToEF9++SV79uxxLklQq1YtmjdvTmxsLI0aNSrpocstBSepDD74wIamkydtz9Kvfw3jx9seJxGR8qjUvlX3zjvv0Lt3b2rVqkV0dDTR0dE3UqeIlCMnT9phuUWL7Hbr1vD227a3SUSkMij2ur0JCQmsWLHCuX3x4kX279/v0aJExPcsXQq33WZDk78/TJgAX32l0CQilUuxg9PlI3s5OTk0a9aMVatWXdF2x44dvPvuuyWvTkS87ttv4ec/h5/8xM5ratkSNmywSw0EBHi7OhGRsuWRK0VdbZrUtm3bGDhwoCeeQkS8YNky28v07rv2unLjxkFqql1uQESkMtIXhkXkCtnZ9nIpb79tt5s1sz/HxHizKhER79O1yUXExfLl0KqVDUoOBzz3HGzZotAkIgIl7HHKzs72cBki4m25ufbyKH/9q91u2hTmzYMf/9i7dYmI+JIS9TiNHDmSmjVrcs899zBx4kQcDgfHjh1zLn4pIuXLypW2l+lSaHrmGdi6VaFJRORyxe5xWr58Odu2bSM9PZ1t27axbt06jDEkJCQwZMgQGjVqRIsWLWjevDknTpwojZpFxENOn7aXR5kzx243bgxz50K3bt6tS0TEV93QyuEAFy5cYNeuXc4glZ6eTnp6OsePH7dP4HBQUFDgkWLLA60cLuXF6tXwxBNw6JDdHj4cJk+G6tW9WpaIiFeU2srhl6tatSpt27albdu2PP744879WVlZziAlIr7jzBm7rMDMmXa7YUPby9Sjh3frEhEpD264x0lcqcdJfNnatTBoEPz3v3b7ySdh6lSoUcOrZYmIeJ27n99ajkCkEjh7FhIT7dyl//4XoqJgxQo7t0mhSUTEfVoAU6SCS0mxvUz/+Y/dfuIJmD4dQkK8WpaISLmkHieRCur8efuNuS5dbGiKjIRPP4W33lJoEhEpKfU4iVRAmzZBQgLs2WO3Bw6EGTPg5pu9WpaISLmnHieRCiQvD158EWJjbWgKD4ePP7aXT1FoEhG5cepxEqkgUlNtL9POnXZ7wAB47TWoVcu7dYmIVCTqcRIp5/LzYeJEexHenTuhTh348EN4912FJhERT1OPk0g5tm2b7WXats1u//SnMGuWDU8iIuJ56nESKYcuXIBXXoGOHW1oqlULFi2C999XaBIRKU3qcRIpZ3bssL1MaWl2+yc/gdmzISzMu3WJiFQG6nESKScuXoSkJLj9dhuabr7ZzmP68EOFJhGRsqIeJ5FyYNcuu/r35s12+8EH4c9/hogIr5YlIlLp+HSP06xZs4iOjiYoKIiYmBg2bdp0zfaLFy+mefPmBAUF0bp1az777DOX+x0OR5G3qVOnOttER0dfcf/kyZNL5fxErqegwF6E9/bbbWgKCbFrMn38sUKTiIg3+GxwWrRoEYmJiUyaNIm0tDTatm1Lz549OX78eJHtv/zyS/r378+QIUPYsmULffr0oU+fPuzYscPZJiMjw+U2d+5cHA4H8fHxLsd6+eWXXdqNHDmyVM9VpCh790LXrvayKXl50Lu3XW5g4EBwOLxdnYhI5eQwxhhvF1GUmJgYOnXqxMyZMwEoLCwkKiqKkSNHMm7cuCva9+vXjzNnzrBs2TLnvjvuuIN27doxZ86cIp+jT58+nD59muTkZOe+6OhoRo8ezejRo0tUd25uLiEhIeTk5BAcHFyiY0jlVlBgF6588UV7vbkaNezlUgYPVmASESkt7n5++2SPU35+PqmpqcTFxTn3+fn5ERcXR0pKSpGPSUlJcWkP0LNnz6u2z8rK4tNPP2XIkCFX3Dd58mRq1apF+/btmTp1KhcvXrxqrXl5eeTm5rrcREpq/37o3h0SE21o+r//s9+ie+IJhSYREV/gk5PDT548SUFBAWGXfVUoLCyMPZeuWnqZzMzMIttnZmYW2f7tt9+mRo0a9O3b12X/qFGjuP3226lZsyZffvkl48ePJyMjg+nTpxd5nKSkJH7zm9+4e2oiRSostAtXjh0L585B9erwxz/C0KEKTCIivsQng1NZmDt3Lo899hhBQUEu+xMTE50/t2nThoCAAJ588kmSkpIIDAy84jjjx493eUxubi5RUVGlV7hUOP/9r+1R+uILu3333TB3LkRHe7UsEREpgk8Gp9q1a+Pv709WVpbL/qysLMLDw4t8THh4uNvt165dy969e1m0aNF1a4mJieHixYscOnSIZs2aXXF/YGBgkYFK5HqMsUsKPPccnDkD1arZb9A99RT4+eQguoiI+OSf54CAADp06OAyabuwsJDk5GRiY2OLfExsbKxLe4DPP/+8yPZvvfUWHTp0oG3bttetZevWrfj5+VG3bt1inoXI1X39Ndx7L/zylzY0de0K6enw9NMKTSIivswne5zADpklJCTQsWNHOnfuzIwZMzhz5gyDBw8GYODAgdSrV4+kpCQAnnnmGbp168Yf//hH7r//fhYuXMhXX33Fm2++6XLc3NxcFi9ezB//+McrnjMlJYWNGzdy9913U6NGDVJSUhgzZgw///nPufnmm0v/pKXCMwbeestO/j59Gm66ya4GPnKkApOISHngs8GpX79+nDhxgokTJ5KZmUm7du1Yvny5cwL44cOH8fvBJ82dd97JggULmDBhAi+++CJNmzZl6dKltGrVyuW4CxcuxBhD//79r3jOwMBAFi5cyEsvvUReXh6NGjVizJgxLnOYRErq6FH4xS9gxQq7feedMG8e3Hqrd+sSERH3+ew6TuWV1nGSyxljV/sePRpyciAwEH73O7vt7+/t6kREBNz//PbZHieRiuDYMRg2DD791G537mxDVPPm3q1LRERKRrMqREqBMfDOO3DbbTY0BQTA5Mmwfr1Ck4hIeaYeJxEPy8qCJ5+0F+IF6NDB9jLddpt36xIRkRunHicRDzEGFi2yAenjj6FqVXjlFUhJUWgSEako1OMk4gEnTtg1mD74wG63awfz54MbS4WJiEg5oh4nkRv04Ye2R+mDD6BKFZg0CTZuVGgSEamI1OMkUkKnTsGIEbBwod1u3dr2Mt1+u1fLEhGRUqQeJ5ES+Phj28u0cKFdi+lXv4LNmxWaREQqOvU4iRTD//t/MGqUXWoAoGVL28vUqZNXyxIRkTKiHicRN336qe1leucde125sWMhNVWhSUSkMlGPk8h1ZGfbi/LOm2e3b73V9jLFxnqzKhER8Qb1OIlcw4oVdtL3vHngcNgAtXWrQpOISGWlHieRIuTmwnPPwV/+YrebNLHhqUsX79YlIiLepR4nkcusXGl7mS6FplGjbC+TQpOIiKjHSeR7330HL7wAs2fb7UaNbC9Tt27erUtERHyHepxEgDVroE2b/4Wmp5+G9HSFJhERcaXgJJXamTN2KO7uu+HgQWjQwA7VzZoF1at7uzoREfE1GqqTSmvdOhg0CA4csNvDhsHUqRAc7NWyRETEh6nHSSqdc+fssgJ33WVDU/36sHw5/PnPCk0iInJt6nGSSiUlxfYy/ec/dvuJJ2D6dAgJ8WpZIiJSTqjHSSqF8+ftN+a6dLGhKTLSXkLlrbcUmkRExH3qcZIKb/NmSEiA3bvt9uOPw5/+BDff7N26RESk/FGPk1RYeXnwq1/Zy6Ps3g1hYbB0KfztbwpNIiJSMupxkgopLc32Mu3YYbf794fXX4datbxbl4iIlG/qcZIKJT8fJk2Czp1taKpTBz74ABYsUGgSEZEbpx4nqTC2bbO9TNu22e1HHoE33rDhSURExBPU4yTl3oUL8Mor0LGjDU21asGiRbB4sUKTiIh4lnqcpFwoKIC1ayEjAyIioGtX8Pe3w3GDBkFqqm3Xpw/MmWMngouIiHiagpP4vCVL4Jln4OjR/+2rX99egHfxYjuv6eab7eTvAQPA4fBerSIiUrEpOIlPW7LEzlUyxnX/0aPw7rv25wcesJdLiYws+/pERKRyUXASn1VQYHuaLg9NP1SzJnz0EVTRO1lERMqAJoeLz1q71nV4rijffgvr1pVNPSIiIgpO4rMyMjzbTkRE5EYpOInPql7dvXYREaVbh4iIyCU+HZxmzZpFdHQ0QUFBxMTEsGnTpmu2X7x4Mc2bNycoKIjWrVvz2Wefudw/aNAgHA6Hy61Xr14ubb799lsee+wxgoODCQ0NZciQIXz33XcePze5ti++gKefvnYbhwOiouzSBCIiImXBZ4PTokWLSExMZNKkSaSlpdG2bVt69uzJ8ePHi2z/5Zdf0r9/f4YMGcKWLVvo06cPffr0Yceli5V9r1evXmRkZDhv7733nsv9jz32GDt37uTzzz9n2bJl/Pvf/2bYsGGldp7iKj8fXnwR7r7bzm+6tB7T5UsMXNqeMcOu5yQiIlIWHMZc6ztL3hMTE0OnTp2YOXMmAIWFhURFRTFy5EjGjRt3Rft+/fpx5swZli1b5tx3xx130K5dO+bMmQPYHqfs7GyWLl1a5HPu3r2bli1bsnnzZjp27AjA8uXLue+++zh69CiRRXzfPS8vj7y8POd2bm4uUVFR5OTkEBwcXOLzr4z+8x947DH46iu7/cQT8Kc/wb/+deU6TlFRNjT17euVUkVEpILJzc0lJCTkup/fPtnjlJ+fT2pqKnFxcc59fn5+xMXFkZKSUuRjUlJSXNoD9OzZ84r2a9asoW7dujRr1oxf/vKXnDp1yuUYoaGhztAEEBcXh5+fHxs3bizyeZOSkggJCXHeoqKiin2+lZ0x8NZb0L69DU0332wXtnzrLTvPqW9fOHQIVq+2F+tdvRoOHlRoEhGRsueTq9+cPHmSgoICwi67bkZYWBh79uwp8jGZmZlFts/MzHRu9+rVi759+9KoUSMOHDjAiy++SO/evUlJScHf35/MzEzq1q3rcowqVapQs2ZNl+P80Pjx40lMTHRuX+pxEvecOgXDhtmFLsEO0f3tb3Zl8B/y94fu3cu8PBERERc+GZxKy6OPPur8uXXr1rRp04ZbbrmFNWvWcM8995TomIGBgQQGBnqqxEolORkGDoRjx+wClr/7HTz7rOYsiYiI7/LJobratWvj7+9PVlaWy/6srCzCw8OLfEx4eHix2gM0btyY2rVrs3//fucxLp98fvHiRb799ttrHkeKJz8fXngB/u//bGi69VbYsMHuU2gSERFf5pPBKSAggA4dOpCcnOzcV1hYSHJyMrGxsUU+JjY21qU9wOeff37V9gBHjx7l1KlTRHy/EFBsbCzZ2dmkpqY626xatYrCwkJiYmJu5JTke3v2wB13wNSpdm7Tk09CWhp06ODtykRERNxgfNTChQtNYGCgmT9/vtm1a5cZNmyYCQ0NNZmZmcYYYx5//HEzbtw4Z/v169ebKlWqmGnTppndu3ebSZMmmapVq5rt27cbY4w5ffq0ee6550xKSoo5ePCgWblypbn99ttN06ZNzfnz553H6dWrl2nfvr3ZuHGjWbdunWnatKnp37+/23Xn5OQYwOTk5HjolagYCguNmT3bmJtuMgaMqVXLmI8+8nZVIiIilruf3z47x6lfv36cOHGCiRMnkpmZSbt27Vi+fLlzAvjhw4fx8/tfh9mdd97JggULmDBhAi+++CJNmzZl6dKltGrVCgB/f3/S09N5++23yc7OJjIyknvvvZdXXnnFZY7Su+++y4gRI7jnnnvw8/MjPj6e1157rWxPvoI5cQJ+8Qv4xz/s9v/9H8yfD0Ws7iAiIuLTfHYdp/LK3XUgKot//QsSEiAzEwICICkJRo8GP58cJBYRkcrK3c9vn+1xkvLt/HkYP94uUgnQooVdg6ldO29WJSIicmMUnMTjdu6EAQMgPd1uP/20nQxerZp36xIREblRGjARjzEGZs6Ejh1taKpTBz75BGbNUmgSEZGKQT1O4hFZWfbacp99Zrd79YJ580DLX4mISEWiHie5YZ99Bm3a2P8GBtoL8372mUKTiIhUPOpxkhI7d86u9j1zpt1u1cpOAG/d2rt1iYiIlBb1OEmJpKdDp07/C02jRsHmzQpNIiJSsSk4SbEUFtolBjp1st+eCwuzw3J/+hMEBXm7OhERkdKloTpxW0YGDB4MK1bY7QcegLfegrp1vVuXiIhIWVGPk7jlH/+wE8BXrLA9S2+8YfcpNImISGWiHie5prNn4dlnYc4cu922rZ0A3rKld+sSERHxBvU4yVVt2QIdOvwvND37LGzcqNAkIiKVl4KTXKGwEKZNg5gY2LMHIiLsxXqnTbPrNImIiFRWGqoTF998AwkJkJxstx9+GP76V6hd27t1iYiI+AL1OInTRx/ZCeDJyfbacm++afcpNImIiFjqcRK++w7GjLE9S2DnNb37LjRr5t26REREfI16nCq5r76C22+3ocnhgLFj4csvFZpERESKoh6nSqqgAKZOhV//Gi5ehHr14O9/h7vv9nZlIiIivkvBqRI6cgQefxy++MJuP/II/PnPULOmd+sSERHxdRqqq2QWL7YTwL/4An70I5g7F95/X6FJRETEHepxqiROn4ZRo2D+fLvdqZNdAbxJE6+WJSIiUq6ox6kS2LgR2rWzocnhgF/9CtavV2gSEREpLvU4VWAFBZCUBC+9ZH9u0MBOAL/rLm9XJiIiUj4pOFVQhw7ZCeDr1tntRx+F2bMhNNSbVYmIiJRvGqqrgN57D9q2taGpRg3429/sfCaFJhERkRujHqcKJCcHRoyAd96x27Gx9ufGjb1bl4iISEWhHqcKYv16OwH8nXfAzw8mTYJ//1uhSURExJPU41TOXbwIr7wCv/0tFBZCdLS9ztydd3q7MhERkYpHwakcKCiAtWshIwMiIqBrV/D3h//+Fx57DDZssO0efxxefx1CQrxbr4iISEWl4OTjliyBZ56Bo0f/t69+ffjJT+y6TKdP26A0ezb07++1MkVERCoFBScftmSJvY6cMa77jx61PUsAXbrYeU0NG5Z9fSIiIpWNgpOPKiiwPU2Xh6YfCgmB5GQICCi7ukRERCozfavOR61d6zo8V5ScHPjyy7KpR0RERHw8OM2aNYvo6GiCgoKIiYlh06ZN12y/ePFimjdvTlBQEK1bt+azzz5z3nfhwgXGjh1L69at+dGPfkRkZCQDBw7k2LFjLseIjo7G4XC43CZPnlwq53ctGRmebSciIiI3zmeD06JFi0hMTGTSpEmkpaXRtm1bevbsyfHjx4ts/+WXX9K/f3+GDBnCli1b6NOnD3369GHHjh0AnD17lrS0NH7961+TlpbGkiVL2Lt3Lw899NAVx3r55ZfJyMhw3kaOHFmq51qUiAjPthMREZEb5zDmWrNovCcmJoZOnToxc+ZMAAoLC4mKimLkyJGMGzfuivb9+vXjzJkzLFu2zLnvjjvuoF27dsyZM6fI59i8eTOdO3fm66+/pkGDBoDtcRo9ejSjR48uUd25ubmEhISQk5NDcHBwiY4Bdo5TdDR8803R85wcDvvtuoMH7dIEIiIiUnLufn77ZI9Tfn4+qampxMXFOff5+fkRFxdHSkpKkY9JSUlxaQ/Qs2fPq7YHyMnJweFwEHrZRdwmT55MrVq1aN++PVOnTuXixYtXPUZeXh65ubkuN0/w94c//cn+7HC43ndpe8YMhSYREZGy5JPB6eTJkxQUFBAWFuayPywsjMzMzCIfk5mZWaz258+fZ+zYsfTv398lWY4aNYqFCxeyevVqnnzySX7/+9/zwgsvXLXWpKQkQkJCnLeoqCh3T/O6+vaFDz6AevVc99evb/f37euxpxIRERE3VMrlCC5cuMDPfvYzjDHMnj3b5b7ExETnz23atCEgIIAnn3ySpKQkAgMDrzjW+PHjXR6Tm5vr8fD08MNFrxwuIiIiZcsng1Pt2rXx9/cnKyvLZX9WVhbh4eFFPiY8PNyt9pdC09dff82qVauuOw8pJiaGixcvcujQIZo1a3bF/YGBgUUGKk/y94fu3Uv1KURERMQNPjlUFxAQQIcOHUhOTnbuKywsJDk5mdjY2CIfExsb69Ie4PPPP3dpfyk07du3j5UrV1KrVq3r1rJ161b8/PyoW7duCc9GREREKgqf7HECO2SWkJBAx44d6dy5MzNmzODMmTMMHjwYgIEDB1KvXj2SkpIAeOaZZ+jWrRt//OMfuf/++1m4cCFfffUVb775JmBD0yOPPEJaWhrLli2joKDAOf+pZs2aBAQEkJKSwsaNG7n77rupUaMGKSkpjBkzhp///OfcfPPN3nkhRERExGf4bHDq168fJ06cYOLEiWRmZtKuXTuWL1/unAB++PBh/Pz+12F25513smDBAiZMmMCLL75I06ZNWbp0Ka1atQLgm2++4R//+AcA7dq1c3mu1atX0717dwIDA1m4cCEvvfQSeXl5NGrUiDFjxrjMYRIREZHKy2fXcSqvPLWOk4iIiJSdcr2Ok4iIiIgvUnASERERcZOCk4iIiIibFJxERERE3KTgJCIiIuImBScRERERNyk4iYiIiLhJwUlERETETQpOIiIiIm5ScBIRERFxk4KTiIiIiJsUnERERETcpOAkIiIi4iYFJxERERE3KTiJiIiIuEnBSURERMRNCk4iIiIiblJwEhEREXGTgpOIiIiImxScRERERNyk4CQiIiLiJgUnERERETcpOImIiIi4ScFJRERExE0KTiIiIiJuUnASERERcZOCk4iIiIibFJxERERE3KTgJCIiIuImBScRERERNyk4iYiIiLhJwUlERETETVW8XYC4oaAA1q6FjAyIiICuXcHf39tVlS29BqL3gEjl5iN/A3y6x2nWrFlER0cTFBRETEwMmzZtumb7xYsX07x5c4KCgmjdujWfffaZy/3GGCZOnEhERAQ33XQTcXFx7Nu3z6XNt99+y2OPPUZwcDChoaEMGTKE7777zuPn5rYlSyA6Gu6+GwYMsP+Njrb7Kwu9BqL3gEjl5kt/A4yPWrhwoQkICDBz5841O3fuNEOHDjWhoaEmKyuryPbr1683/v7+ZsqUKWbXrl1mwoQJpmrVqmb79u3ONpMnTzYhISFm6dKlZtu2beahhx4yjRo1MufOnXO26dWrl2nbtq3ZsGGDWbt2rWnSpInp37+/23Xn5OQYwOTk5JT85C/58ENjHA5jwPXmcNjbhx/e+HP4Or0GoveASOVWRn8D3P38dhhjTNnHteuLiYmhU6dOzJw5E4DCwkKioqIYOXIk48aNu6J9v379OHPmDMuWLXPuu+OOO2jXrh1z5szBGENkZCTPPvsszz33HAA5OTmEhYUxf/58Hn30UXbv3k3Lli3ZvHkzHTt2BGD58uXcd999HD16lMjIyOvWnZubS0hICDk5OQQHB5f8BSgosGn66NGi73c4IDwc1q2ruMMVBQXQpYvtli1KZXgNKjt33wNr1+o9UFn45keWlJaCAjskl5lZ9P0OB9SvDwcP3vDfAHc/v31yjlN+fj6pqamMHz/euc/Pz4+4uDhSUlKKfExKSgqJiYku+3r27MnSpUsBOHjwIJmZmcTFxTnvDwkJISYmhpSUFB599FFSUlIIDQ11hiaAuLg4/Pz82LhxIz/5yU+ueN68vDzy8vKc27m5uSU65yusXXv10AT2j0dGBtxyi2eerzzSayCX3gNNmni7EhHxBmPgyBH7mdm9e5k8pU8Gp5MnT1JQUEBYWJjL/rCwMPbs2VPkYzIzM4tsn/l9Sr303+u1qVu3rsv9VapUoWbNms42l0tKSuI3v/mNm2dWDFf7P+zLVa1acf9Pu6AALly4fruK/BpUdsV5D1TxyT9nUhocDm9XIGXl4kXIz79+O3c/Mz1Af2lu0Pjx4116unJzc4mKirrxA0dEuNfuX/8qs5Rd5tassRMAr6civwaVnd4DIpWbu38D3P3M9ACf/FZd7dq18ff3Jysry2V/VlYW4eHhRT4mPDz8mu0v/fd6bY4fP+5y/8WLF/n222+v+ryBgYEEBwe73Dyia1c7bnu1/7NyOCAqyrarqPQaiN4DIpWbD/4N8MngFBAQQIcOHUhOTnbuKywsJDk5mdjY2CIfExsb69Ie4PPPP3e2b9SoEeHh4S5tcnNz2bhxo7NNbGws2dnZpKamOtusWrWKwsJCYmJiPHZ+bvH3hz/9yf58+Rvm0vaMGRV7iEqvgeg9IFK5+eLfAI98h68ULFy40AQGBpr58+ebXbt2mWHDhpnQ0FCTmZlpjDHm8ccfN+PGjXO2X79+valSpYqZNm2a2b17t5k0aVKRyxGEhoaajz/+2KSnp5uHH364yOUI2rdvbzZu3GjWrVtnmjZt6r3lCIyxX7OsX9/1K5hRUZXrK9h6DUTvAZHKrQz+BpT75QgAZs6cydSpU8nMzKRdu3a89tprzp6f7t27Ex0dzfz5853tFy9ezIQJEzh06BBNmzZlypQp3Hfffc77jTFMmjSJN998k+zsbLp06cIbb7zBrbfe6mzz7bffMmLECD755BP8/PyIj4/ntddeo3r16m7V7LHlCH7IR1ZL9Sq9BqL3gEjlVsp/A9z9/Pbp4FQelUpwEhERkVLl7ue3T85xEhEREfFFCk4iIiIiblJwEhEREXGTgpOIiIiImxScRERERNyk4CQiIiLiJgUnERERETcpOImIiIi4ScFJRERExE1VvF1ARXNpIfbc3FwvVyIiIiLuuvS5fb0Lqig4edjp06cBiIqK8nIlIiIiUlynT58mJCTkqvfrWnUeVlhYyLFjx6hRowYOh8Njx83NzSUqKoojR45U2mvgVfbXoLKfP+g10PlX7vMHvQalef7GGE6fPk1kZCR+flefyaQeJw/z8/Ojfv36pXb84ODgSvmP5Ycq+2tQ2c8f9Bro/Cv3+YNeg9I6/2v1NF2iyeEiIiIiblJwEhEREXGTglM5ERgYyKRJkwgMDPR2KV5T2V+Dyn7+oNdA51+5zx/0GvjC+WtyuIiIiIib1OMkIiIi4iYFJxERERE3KTiJiIiIuEnBSURERMRNCk5laNasWURHRxMUFERMTAybNm26ZvvFixfTvHlzgoKCaN26NZ999pnzvgsXLjB27Fhat27Nj370IyIjIxk4cCDHjh1zOUZ0dDQOh8PlNnny5FI5v+vx5PkDDBo06Ipz69Wrl0ubb7/9lscee4zg4GBCQ0MZMmQI3333ncfPzR2ePv/Lz/3SberUqc42vvT7h+K9Bjt37iQ+Pt55DjNmzCjRMc+fP8/w4cOpVasW1atXJz4+nqysLE+elts8ff5JSUl06tSJGjVqULduXfr06cPevXtd2nTv3v2K98BTTz3l6VNzm6dfg5deeumK82vevLlLm4r8Hijq37jD4WD48OHONr70HijO+f/lL3+ha9eu3Hzzzdx8883ExcVd0d4Yw8SJE4mIiOCmm24iLi6Offv2ubTx+OeAkTKxcOFCExAQYObOnWt27txphg4dakJDQ01WVlaR7devX2/8/f3NlClTzK5du8yECRNM1apVzfbt240xxmRnZ5u4uDizaNEis2fPHpOSkmI6d+5sOnTo4HKchg0bmpdfftlkZGQ4b999912pn+/lPH3+xhiTkJBgevXq5XJu3377rctxevXqZdq2bWs2bNhg1q5da5o0aWL69+9fqudalNI4/x+ed0ZGhpk7d65xOBzmwIEDzja+8vs3pvivwaZNm8xzzz1n3nvvPRMeHm5effXVEh3zqaeeMlFRUSY5Odl89dVX5o477jB33nlnaZ3mVZXG+ffs2dPMmzfP7Nixw2zdutXcd999pkGDBi6/427dupmhQ4e6vAdycnJK6zSvqTReg0mTJpnbbrvN5fxOnDjh0qYivweOHz/ucu6ff/65Aczq1audbXzlPVDc8x8wYICZNWuW2bJli9m9e7cZNGiQCQkJMUePHnW2mTx5sgkJCTFLly4127ZtMw899JBp1KiROXfunLONpz8HFJzKSOfOnc3w4cOd2wUFBSYyMtIkJSUV2f5nP/uZuf/++132xcTEmCeffPKqz7Fp0yYDmK+//tq5r2HDhkX+YytrpXH+CQkJ5uGHH77qc+7atcsAZvPmzc59//znP43D4TDffPNNCc+kZMri9//www+bHj16uOzzld+/McV/DX7oaudxvWNmZ2ebqlWrmsWLFzvb7N692wAmJSXlBs6m+Erj/C93/PhxA5gvvvjCua9bt27mmWeeKUnJHlcar8GkSZNM27Ztr/q4yvYeeOaZZ8wtt9xiCgsLnft85T1wI+dvjDEXL140NWrUMG+//bYxxpjCwkITHh5upk6d6myTnZ1tAgMDzXvvvWeMKZ3PAQ3VlYH8/HxSU1OJi4tz7vPz8yMuLo6UlJQiH5OSkuLSHqBnz55XbQ+Qk5ODw+EgNDTUZf/kyZOpVasW7du3Z+rUqVy8eLHkJ1MCpXn+a9asoW7dujRr1oxf/vKXnDp1yuUYoaGhdOzY0bkvLi4OPz8/Nm7c6IlTc0tZ/P6zsrL49NNPGTJkyBX3efv3DyV7DTxxzNTUVC5cuODSpnnz5jRo0KDEz1tatXpCTk4OADVr1nTZ/+6771K7dm1atWrF+PHjOXv2rMee012l+Rrs27ePyMhIGjduzGOPPcbhw4ed91Wm90B+fj7vvPMOTzzxxBUXmff2e8AT53/27FkuXLjgfH8fPHiQzMxMl2OGhIQQExPjPGZpfA7oIr9l4OTJkxQUFBAWFuayPywsjD179hT5mMzMzCLbZ2ZmFtn+/PnzjB07lv79+7tc+HDUqFHcfvvt1KxZky+//JLx48eTkZHB9OnTb/Cs3Fda59+rVy/69u1Lo0aNOHDgAC+++CK9e/cmJSUFf39/MjMzqVu3rssxqlSpQs2aNa/6OpaGsvj9v/3229SoUYO+ffu67PeF3z+U7DXwxDEzMzMJCAi44n8mrvValobSOP/LFRYWMnr0aH784x/TqlUr5/4BAwbQsGFDIiMjSU9PZ+zYsezdu5clS5Z45HndVVqvQUxMDPPnz6dZs2ZkZGTwm9/8hq5du7Jjxw5q1KhRqd4DS5cuJTs7m0GDBrns94X3gCfOf+zYsURGRjqD0qXf37X+VpbG54CCUwVw4cIFfvazn2GMYfbs2S73JSYmOn9u06YNAQEBPPnkkyQlJZX7JfsfffRR58+tW7emTZs23HLLLaxZs4Z77rnHi5WVvblz5/LYY48RFBTksr8i//7F1fDhw9mxYwfr1q1z2T9s2DDnz61btyYiIoJ77rmHAwcOcMstt5R1mR7Xu3dv589t2rQhJiaGhg0b8v777xfZA1uRvfXWW/Tu3ZvIyEiX/RXhPTB58mQWLlzImjVrrvg7V9Y0VFcGateujb+//xXf4sjKyiI8PLzIx4SHh7vV/lJo+vrrr/n8889depuKEhMTw8WLFzl06FDxT6SESvP8f6hx48bUrl2b/fv3O49x/PhxlzYXL17k22+/veZxPK20z3/t2rXs3buXX/ziF9etxRu/fyjZa+CJY4aHh5Ofn092drbHnrckSuP8f2jEiBEsW7aM1atXU79+/Wu2jYmJAXD+Oykrpf0aXBIaGsqtt97q8negMrwHvv76a1auXOn23wEo2/fAjZz/tGnTmDx5Mv/6179o06aNc/+lx13vb4CnPwcUnMpAQEAAHTp0IDk52bmvsLCQ5ORkYmNji3xMbGysS3uAzz//3KX9pdC0b98+Vq5cSa1ata5by9atW/Hz87ui67I0ldb5X+7o0aOcOnWKiIgI5zGys7NJTU11tlm1ahWFhYXOPxxlobTP/6233qJDhw60bdv2urV44/cPJXsNPHHMDh06ULVqVZc2e/fu5fDhwyV+3tKqtSSMMYwYMYKPPvqIVatW0ahRo+s+ZuvWrQDOfydlpbReg8t99913HDhwwHl+Ff09cMm8efOoW7cu999//3XbeuM9UNLznzJlCq+88grLly93macE0KhRI8LDw12OmZuby8aNG53HLJXPgRJNKZdiW7hwoQkMDDTz5883u3btMsOGDTOhoaEmMzPTGGPM448/bsaNG+dsv379elOlShUzbdo0s3v3bjNp0iSXr6Pn5+ebhx56yNSvX99s3brV5WumeXl5xhhjvvzyS/Pqq6+arVu3mgMHDph33nnH1KlTxwwcOLDcn//p06fNc889Z1JSUszBgwfNypUrze23326aNm1qzp8/7zxOr169TPv27c3GjRvNunXrTNOmTb22HIEnz/+SnJwcU61aNTN79uwrntOXfv/GFP81yMvLM1u2bDFbtmwxERER5rnnnjNbtmwx+/btc/uYxtivojdo0MCsWrXKfPXVVyY2NtbExsaW3Ym7WWtJzv+Xv/ylCQkJMWvWrHH5G3D27FljjDH79+83L7/8svnqq6/MwYMHzccff2waN25s7rrrrrI9+e+Vxmvw7LPPmjVr1piDBw+a9evXm7i4OFO7dm1z/PhxZ5uK/B4wxn47rUGDBmbs2LFXPKcvvQeKe/6TJ082AQEB5oMPPnB5f58+fdqlTWhoqPn4449Nenq6efjhh4tcjsCTnwMKTmXo9ddfNw0aNDABAQGmc+fOZsOGDc77unXrZhISElzav//+++bWW281AQEB5rbbbjOffvqp876DBw8aoMjbpfU7UlNTTUxMjAkJCTFBQUGmRYsW5ve//71LsChLnjz/s2fPmnvvvdfUqVPHVK1a1TRs2NAMHTrU5QPTGGNOnTpl+vfvb6pXr26Cg4PN4MGDXf7RlSVPnv8lf/7zn81NN91ksrOzr7jP137/xhTvNbjae7xbt25uH9MYY86dO2eefvppc/PNN5tq1aqZn/zkJyYjI6M0T/OqPH3+V/sbMG/ePGOMMYcPHzZ33XWXqVmzpgkMDDRNmjQxzz//vNfWcTLG869Bv379TEREhAkICDD16tUz/fr1M/v373d5zor8HjDGmBUrVhjA7N2794rn87X3QHHOv2HDhkWe/6RJk5xtCgsLza9//WsTFhZmAgMDzT333HPF6+DpzwGHMcaUrK9KREREpHLRHCcRERERNyk4iYiIiLhJwUlERETETQpOIiIiIm5ScBIRERFxk4KTiIiIiJsUnERERETcpOAkIiIi4iYFJxGR6/juu+/w8/Nj+vTp3i5FRLxMwUlE5Dp27NiBMYbbbrvN26WIiJcpOImIXMf27dsBaNmypZcrERFvU3ASEbmO7du3ExwcTFRUlLdLEREvU3ASEbmO7du306JFC9LS0ujduzc1atSgXr16/OlPf/J2aSJSxhzGGOPtIkREfFnt2rUJCwvj//2//8fgwYNp0KABf/nLX0hLS2Pbtm20bt3a2yWKSBmp4u0CRER8WUZGBqdOncLhcJCWluYcrrvrrrto2bIlW7ZsUXASqUQ0VCcicg3p6ekAvPzyyy5znKpWrQpAQECAV+oSEe9QcBIRuYZL36j7yU9+4rJ/z549ADRr1qzMaxIR71FwEhG5hu3bt1OvXj3Cw8Nd9m/bto0qVapoiQKRSkbBSUTkGrZv306bNm2u2J+ens6tt95KYGCgF6oSEW9RcBIRuYqCggJ2795N27Ztr7hv27ZtRQYqEanYFJxERK5i3759nD9//oqAdO7cOfbv36/gJFIJKTiJiFzFpYnhlwekHTt2UFBQoOAkUglpAUwRERERN6nHSURERMRNCk4iIiIiblJwEhEREXGTgpOIiIiImxScRERERNyk4CQiIiLiJgUnERERETcpOImIiIi4ScFJRERExE0KTiIiIiJuUnASERERcZOCk4iIiIib/j+2do/hrQ3tKAAAAABJRU5ErkJggg==\n"
},
"metadata": {}
}
]
}
]
}