{ "cells": [ { "cell_type": "markdown", "id": "c6776669-1b39-4076-a156-9fe383220bc3", "metadata": {}, "source": [ "# Mercury's relativistic precession\n", "\n", "One of the early successes of [general relativity](https://en.wikipedia.org/wiki/General_relativity) (GR) was the explanation of an [excess in the precession of Mercury's perihelion](https://en.wikipedia.org/wiki/Mercury_(planet)#Advance_of_perihelion) that could not be accounted for via Newtonian gravity. In this tutorial\n", "we will show how we can reproduce this classical result with heyoka.py.\n", "\n", "In the weak gravitational field regime, the formulation of GR can be greatly simplified via the formalism of [post-Newtonian approximations](https://en.wikipedia.org/wiki/Post-Newtonian_expansion) (PN). In the specific case of Mercury's relativistic perihelion precession, a first-order post-Newtonian expansion of the full GR equations is sufficient. The 1PN Hamiltonian of a test particle orbiting a spherical mass reads, in Cartesian coordinates:\n", "\n", "$$\n", "\\mathcal{H}_\\mathrm{1PN} \\left(v_x, v_y, v_z, x, y, z \\right) = \\frac{1}{2}v^2-\\frac{\\mu}{r} + \\varepsilon \\left(\\frac{\\mu^2}{2r^2} -\\frac{1}{8}v^4-\\frac{3}{2}\\frac{\\mu v^2}{r} \\right),\n", "$$\n", "\n", "where $\\varepsilon = 1/c^2$, $v=\\sqrt{v_x^2+v_y^2+v_z^2}$, $r=\\sqrt{x^2+y^2+z^2}$ and $\\mu=GM$ is the gravitational parameter of the system. In other words, $\\mathcal{H}_\\mathrm{1PN}$ is a Keplerian one-body Hamiltonian augmented by a velocity-dependent perturbation of magnitude $\\varepsilon$.\n", "\n", "Let us write down the Hamiltonian using heyoka.py's expression system:" ] }, { "cell_type": "code", "execution_count": 1, "id": "36bead79-c512-4fcf-935a-83bdfef1d98f", "metadata": {}, "outputs": [], "source": [ "import heyoka as hy\n", "\n", "# Create the symbolic variables.\n", "vx, vy, vz, x, y, z = hy.make_vars(\"vx\", \"vy\", \"vz\", \"x\", \"y\", \"z\")\n", "\n", "# mu and eps constants.\n", "mu = 0.01720209895 * 0.01720209895 * 365 * 365\n", "eps = 2.5037803127808595e-10\n", "\n", "# Auxiliary quantities.\n", "v2 = vx*vx + vy*vy + vz*vz\n", "r2 = x*x + y*y + z*z\n", "r = hy.sqrt(r2)\n", "\n", "# Define the Hamiltonian.\n", "Ham = 1./2*v2 - mu/r + eps * (mu**2/(2*r2) - 1/8.*v2*v2 - 3./2.*mu*v2/r)" ] }, { "cell_type": "markdown", "id": "f671aa05-e1a5-4f60-8ec5-ba76edd8e78d", "metadata": {}, "source": [ "As usual, we are using Solar masses, astronomical units and years as units of measure.\n", "\n", "Let us now define the initial conditions for Mercury's orbit. Semi-major axis, eccentricity and inclination are those of Mercury's current orbit. The initial argument of perihelion $\\omega_0$ is set to a value of 1:" ] }, { "cell_type": "code", "execution_count": 2, "id": "f27796a6-84fb-4309-b00b-5006ed5c153f", "metadata": {}, "outputs": [], "source": [ "import pykep as pk\n", "import numpy as np\n", "\n", "# Mercury's initial position and velocity.\n", "r0, v0 = pk.par2ic([0.387098, 0.205630, 7*2*np.pi/360, 0., 1., 0.], mu)" ] }, { "cell_type": "markdown", "id": "e38dd5e8-117a-4b4b-adfa-71199252a148", "metadata": {}, "source": [ "We can now proceed to the creation of the integrator object. We will be using the {func}`~heyoka.hamiltonian()` function to automatically generate the equations of motion from the Hamiltonian:" ] }, { "cell_type": "code", "execution_count": 3, "id": "f265e089-4fee-47a5-a6c8-20f378c00c0e", "metadata": {}, "outputs": [], "source": [ "ta = hy.taylor_adaptive(\n", " # Hamilton's equations.\n", " hy.hamiltonian(Ham, [x, y, z], [vx, vy, vz]),\n", " # Initial conditions.\n", " r0 + v0\n", ")" ] }, { "cell_type": "markdown", "id": "1ea4d7a4-ea43-481a-b201-70c5b17871b0", "metadata": {}, "source": [ "Let us now numerically integrate for 100 years:" ] }, { "cell_type": "code", "execution_count": 4, "id": "cf09022e-764b-410b-9bd2-d4cd86e591f1", "metadata": {}, "outputs": [], "source": [ "tgrid = np.linspace(0, 100, 1000)\n", "\n", "_, _, _, _, _, out = ta.propagate_grid(tgrid)" ] }, { "cell_type": "markdown", "id": "869777a8-17bb-4d37-9181-8b0fc7772097", "metadata": {}, "source": [ "We can now convert the result of the numerical integration to Keplerian orbital elements:" ] }, { "cell_type": "code", "execution_count": 5, "id": "238cf32d-c3c3-482e-bd67-2d6f143f2a7f", "metadata": {}, "outputs": [], "source": [ "kep_out = np.array([pk.ic2par(_[0:3], _[3:], mu = mu) for _ in out])" ] }, { "cell_type": "markdown", "id": "cd44da3d-32cd-49ea-b03f-59d6f470efa9", "metadata": {}, "source": [ "We are now ready to visualise the results. GR predicts for Mercury a relativistic precession of $42.98^{\\prime\\prime}$ per century, or roughly $2.08373\\times 10^{-4}\\,\\mathrm{rad}$ per century. Let us then plot the time evolution of the argument of perihelion and compare it to the theoretical prediction:" ] }, { "cell_type": "code", "execution_count": 6, "id": "710fb755-116b-4ead-a445-5a0f75caa783", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABKUAAAJOCAYAAABm7rQwAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAADFcElEQVR4nOzdZ3gV1f728e/OTu+dBBJ67xIUASkiCqJHFFRURLCACNIhsddzDiSAINJUBBtSjqgUUUFAaqT3XgIhEAgEkpBe9jwv+LMftykETEiA+3NdXDLzm7VmrY1ocmfNGpNhGAYiIiIiIiIiIiI3kF1ZD0BERERERERERG4/CqVEREREREREROSGUyglIiIiIiIiIiI3nEIpERERERERERG54RRKiYiIiIiIiIjIDadQSkREREREREREbjiFUiIiIiIiIiIicsMplBIRERERERERkRtOoZSIiIiIiIiIiNxwCqVEREREREREROSGUyh1E/rss89o3749np6emEwmkpKSynpIIiIiIiIiIiLXRKFUOdW+fXu+/PLLAmvp6el07tyZN95448YOSkRERERERESkhNiX9QDk2g0dOhSAP/74o0zHISIiIiIiIiJyvbRSSkREREREREREbjiFUiIiIiIiIiIicsMplCon/vvf/+Lu7m79tXbtWvr375/vnIiIiIiIiIjIrUB7SpUT/fv358knn7Qe9+zZk+7du9OtWzfruUqVKpXF0ERERERERERESpxCqXLC19cXX19f67GLiwuBgYHUrFmzDEclIiIiIiIiIlI6FErdhM6cOcOZM2c4cuQIALt378bDw4PKlSvbBFsiIiIiIiIiIuWV9pS6CU2fPp077riDvn37AtC2bVvuuOMOFi1aVMYjExEREREREREpHpNhGEZZD0JERERERERERG4vWiklIiIiIiIiIiI3nEIpERERERERERG54bTReRmzWCycPn0aDw8PTCZTWQ9HREREREREROQfMQyDS5cuUbFiRezsCl8PpVCqjJ0+fZrQ0NCyHoaIiIiIiIiISIk6efIkISEhhdYVSpUxDw8P4PIflKenZxmPRkRERERERETkn0lJSSE0NNSaeRRGoVQZu/LInqenp0IpEREREREREbllXG2bIm10LiIiIiIiIiIiN5xCKRERERERERERueEUSomIiIiIiIiIyA2nPaVuEnl5eeTk5JT1METKFQcHB8xmc1kPQ0RERERERK6DQqlyzjAMzpw5Q1JSUlkPRaRc8vb2Jigo6Kob6ImIiIiIiEj5olCqnLsSSAUGBuLq6qpvvEX+j2EYpKenk5CQAEBwcHAZj0hERERERESuhUKpciwvL88aSPn5+ZX1cETKHRcXFwASEhIIDAzUo3wiIiIiIiI3EW10Xo5d2UPK1dW1jEciUn5d+fuhPddERERERERuLgqlbgJ6ZE+kcPr7ISIiIiIicnNSKCUiIiIiIiIiIjecQimRv6hatSoTJ04ssf7at2/P0KFDC62/9957NG3atMTuV54cP34ck8nEjh07ynooIiIiIiIiUg4plJIS16dPH0wmE2PGjLE5/9NPP5X7R602b95Mv379btj9Ro4cyYoVK66pTUkHZyWhT58+PProozbnQkNDiY+Pp2HDhmUzKBERERERESnXFEpJqXB2diYyMpKLFy+W9VCKJTs7G4CAgIAburG8u7t7uX6z4j/ZPNxsNhMUFIS9vV7yKSIiIiIiIvkplJJS0bFjR4KCghg9enSh1xT06NrEiROpWrWq9fjKCpz//ve/VKhQAW9vb95//31yc3MZNWoUvr6+hISEMHPmTJt+Tp06RY8ePfDx8cHPz4+uXbty/PjxfP2OHj2aihUrUrt2bSD/KqSkpCT69etHhQoVcHZ2pmHDhixZsgSAxMREnn76aUJCQnB1daVRo0bMmTPnmj6nv38GV8Y1btw4goOD8fPzY+DAgdZwqH379pw4cYJhw4ZhMplsVp5t2LCBtm3b4uLiQmhoKIMHDyYtLc1aj4+P56GHHsLFxYVq1arx3Xff5ZuvyWRi+vTpdO3aFTc3N/7973+Tl5fHiy++SLVq1XBxcaFOnTp8/PHHNnP46quvWLhwoXVMf/zxR4GP761evZq77roLJycngoODee2118jNzbXW27dvz+DBgwkPD8fX15egoCDee++9a/pMRURERERE5OagJQw3EcMwyMjJK5N7uziYr+nRO7PZzH//+1+eeeYZBg8eTEhIyHXfe+XKlYSEhLBmzRrWr1/Piy++SHR0NG3btmXjxo3MmzeP/v37c//99xMaGkp6ejr33nsvbdq0Yc2aNdjb2/Pvf/+bzp07s2vXLhwdHQFYsWIFnp6eLF++HMMw8t3XYrHw4IMPcunSJb799ltq1KjBvn37MJvNAGRmZhIWFkZERASenp78/PPP9OrVi+rVq9OiRYvrnu+qVasIDg5m1apVHDlyhB49etC0aVP69u3LDz/8QJMmTejXrx99+/a1ttm9ezedOnXiww8/5IsvvuDcuXO8+uqrvPrqq8yaNQuA5557jvPnz/PHH3/g4ODA8OHDSUhIyHf/d999l9GjRzNhwgTMZjMWi4WQkBDmz5+Pv78/GzZsoF+/fgQHB/Pkk08ycuRI9u/fT0pKivVevr6+nD592qbfU6dO0aVLF/r06cPXX3/NgQMH6Nu3L87OzjbB01dffcXw4cPZuHEj0dHR9OnTh9atW3P//fdf92cqIiIiIiIi5Y9CqZtIRk4e9d/5rUzuve+DTrg6Xtu/Lo899hhNmzbl3Xff5Ysvvrjue/v6+jJp0iTs7OyoU6cOUVFRpKen88YbbwDw+uuvM2bMGNavX89TTz3F3LlzsbOzY8aMGdYgbdasWXh7e/PHH3/wwAMPAODm5saMGTOsIdXf/f7772zatIn9+/dbV1JVr17dWq9UqRIjR460Hg8aNIhff/2V//3vf/8olPLx8WHy5MmYzWbq1q3LQw89xIoVK+jbty++vr6YzWY8PDwICgqythk7dizPPPOMdVP1WrVqMWnSJNq1a8e0adM4fvw4v//+O5s3b6Z58+YAzJgxg1q1auW7/zPPPMMLL7xgc+7999+3/r5atWps2LCB+fPn8+STT+Lu7o6LiwtZWVk2Y/q7qVOnEhoayuTJkzGZTNStW5fTp08TERHBO++8g53d5YWbjRs35t1337XOY/LkyaxYsUKhlIiIiIiIyC1GoZSUqsjISDp06MCIESOuu48GDRpYAwuAChUq2GyebTab8fPzs6762bp1K0eOHMHDw8Omn8zMTI4ePWo9btSoUaGBFMCOHTsICQmxBlJ/l5eXx5gxY5g3bx6nTp0iKyuLrKws3NzcrmueVzRo0MC6GgsgODiY3bt3F9nmypxnz55tPWcYBhaLhZiYGA4dOoS9vT3NmjWz1mvWrImPj0++vq6EVn81ffp0ZsyYwYkTJ8jIyCA7O/ua3xq4f/9+WrZsabPirnXr1qSmphIXF0flypWBy6HUXwUHBxe4oktERERERERubgqlbiIuDmb2fdCpzO59Pdq2bUunTp1444036NOnj03Nzs4u32NzBW2s7eDgYHNsMpkKPGexWIDLj92FhYXZBDRXBAQEWH9/tfDIxcWlyPr48eOZMGECEydOpFGjRri5uTF06FDrpunXq6i5FcZisfDyyy8zePDgfLXKlStz8ODBAtsV9Nji3z+X+fPnM2zYMMaPH0/Lli3x8PBg7NixbNy48WpTyXevvz8CeuX+fz1/PfMXERERERGRm49CqZuIyWS65kfoyoMxY8bQtGnTfCuOAgICOHPmjE1Y8ddNsa9Xs2bNmDdvHoGBgXh6el53P40bNyYuLo5Dhw4VuFpq7dq1dO3alWeffRa4HAwdPnyYevXqXfc9i8PR0ZG8PNu9xZo1a8bevXupWbNmgW3q1q1Lbm4u27dvJywsDIAjR46QlJR01futXbuWVq1aMWDAAOu5v644K2xMf1e/fn0WLFhg8+e9YcMGPDw8qFSp0lXHISIiIiIiIrcWvX1PSl2jRo3o2bMnn3zyic359u3bc+7cOaKiojh69ChTpkzhl19++cf369mzJ/7+/nTt2pW1a9cSExPD6tWrGTJkCHFxccXup127drRt25bu3buzfPlyYmJi+OWXX/j111+By4+/LV++nA0bNrB//35efvllzpw584/HfzVVq1ZlzZo1nDp1ivPnzwMQERFBdHQ0AwcOZMeOHRw+fJhFixYxaNAg4HIo1bFjR/r168emTZvYvn07/fr1w8XF5aob2NesWZMtW7bw22+/cejQId5++202b96cb0y7du3i4MGDnD9/vsAVbwMGDODkyZMMGjSIAwcOsHDhQt59912GDx9u83imiIiIiIjI7erAgQPX/FTKzUzfCcoN8eGHH+Z7VKxevXpMnTqVKVOm0KRJEzZt2mSzcfj1cnV1Zc2aNVSuXJlu3bpRr149XnjhBTIyMq555dSCBQu48847efrpp6lfvz7h4eHWFUFvv/02zZo1o1OnTrRv356goCAeffTRfzz+q/nggw84fvw4NWrUsD6O2LhxY1avXs3hw4dp06YNd9xxB2+//TbBwcHWdl9//TUVKlSgbdu2PPbYY/Tt2xcPDw+cnZ2LvF///v3p1q0bPXr0oEWLFiQmJtqsmgLo27cvderUoXnz5gQEBLB+/fp8/VSqVImlS5eyadMmmjRpQv/+/XnxxRd56623SuBTERERERERuXlt3ryZbt26Ub9+fV555ZUCt1q5FZmM22Wm5VRKSgpeXl4kJyfnC0wyMzOJiYmhWrVqVw0ORK5VXFwcoaGh/P7779x3331lPZzrpr8nIiIiIiJysxs1ahTjxo0D4JFHHuGbb775R9vRlLWiso6/uvk2KBKR67Jy5UpSU1Np1KgR8fHxhIeHU7VqVdq2bVvWQxMREREREblt5OTkMG/ePGrVqkWLFi0AGDp0KImJiYwcOZL69euX8QhvHIVSIreJnJwc3njjDY4dO4aHhwetWrVi9uzZ+d52JyIiIiIiIiUvLS2NGTNm8NFHHxEbG0vnzp2t+ypXqlSJmTNnlvEIbzyFUiK3iU6dOtGpU6eyHoaIiIiIiMht5fz583zyySdMnjyZCxcuABAYGEibNm1s3k5+O1IoJSIiIiIiIiJSCiIjI3n//ffJyMgAoEaNGowcOZLevXvj4uJSxqMrewqlRERERERERERKyF9XP/n6+pKRkUFYWBgRERF069YNs9lcxiMsP+zKegAiIiIiIiIiIjczwzBYtWoVnTt3ZsaMGdbzzz33HCtWrGDz5s088cQThQZSmTl5fLhkH1P/OEJqVu6NGnaZ00opEREREREREZHrkJeXx48//khUVBSbN28G4MSJE7z00kuYTCacnJzo0KFDoe3TsnJ59bttbD+ZRFJ6DgCfrznGy+1q8HzrqjjZ39qrqrRSSkRERERERETkGmRmZvLZZ59Rr149nnjiCTZv3oyzszMDBgzg559/Ltbm5WeSM5mzKZZVB89ZAykvFwcupucwf8tJzLfBBuhaKSUiIiIiIiIicg1efPFFvvvuOwB8fHwYOHAggwYNIjAwsNA2eRaDDxbvJcDDieSMHL5YF4PF+P91LxcH1kXcy297z+Ln5oi9+dZfR3Trz1DKnePHj2MymdixY0dZD+WaVK1alYkTJ5ZYf+3bt2fo0KEl1t8/1adPHx599FHrcUmMr7zNUURERERE5HqcOnWKc+fOWY/79etHSEgIH330EbGxsXz44YdFBlInL6SzYFscX0WfYNyyQ3y+1jaQ6tumGp/2CsPD2YHHw0K4t27hfd1KtFJKStTVlij27t2b995778YM5jp9+eWXDB06lKSkJJvzmzdvxs3NrWwGVQZ++OEHHBwcinXtH3/8wb333svFixfx9va+rj5ERERERETKmwMHDjB27Fi++eYbBg8ezLhx4wBo27Ytx44dK/L7ndjEdMxmE7l5FjpNXENmjsVaszPBiAfqsHjnae6s6subD9Uv9bmURwqlpETFx8dbfz9v3jzeeecdDh48aD3n4uLCxYsXy2Jo5OXlYTKZsLO7vgWCAQEBJTyikpednY2jo2OJ9OXr61su+hAREREREbnRoqOjiYyMZOHChdZz+/fvxzAMTCYTJpOp0EDq5IV0Ei5l0XvmJjJz8sj965IooHfLKjzUuCJ3VfNl4L01S3Ue5Z0e35MSFRQUZP3l5eWFyWTKd+6KY8eOce+99+Lq6kqTJk2Ijo626WvDhg20bdsWFxcXQkNDGTx4MGlpadb6xYsXee655/Dx8cHV1ZUHH3yQw4cPW+tffvkl3t7eLFmyhPr16+Pk5MSJEyfIzs4mPDycSpUq4ebmRosWLfjjjz+Ayyt+nn/+eZKTk63/obmysuvvj+8lJSXRr18/KlSogLOzMw0bNmTJkiUAJCYm8vTTTxMSEoKrqyuNGjVizpw51/RZvvfeezRt2pRPP/2U0NBQXF1deeKJJ2xWcF155G706NFUrFiR2rVrA5eXlvbo0QMfHx/8/Pzo2rUrx48ft7bLy8tj+PDheHt74+fnR3h4OIZh+x/Kvz96l5WVRXh4OKGhoTg5OVGrVi2++OILjh8/zr333gtcfpbaZDLRp0+fAvso7p/Zb7/9Rr169XB3d6dz5842YaeIiIiIiEhpWbZsGW3btqVVq1bWQKpr166sX7++WBuYJ6fn0GXSWrpP20BqVm6+QKpRJS/ee6QBd1XTD/BBodRNKS0trdBfmZmZxb42IyOjWNeWljfffJORI0eyY8cOateuzdNPP01ubi4Au3fvplOnTnTr1o1du3Yxb9481q1bx6uvvmpt36dPH7Zs2cKiRYuIjo7GMAy6dOlCTk6O9Zr09HRGjx7NjBkz2Lt3L4GBgTz//POsX7+euXPnsmvXLp544gk6d+7M4cOHadWqFRMnTsTT05P4+Hji4+MZOXJkvrFbLBYefPBBNmzYwLfffsu+ffsYM2YMZvPl13VmZmYSFhbGkiVL2LNnD/369aNXr15s3Ljxmj6jI0eOMH/+fBYvXsyvv/7Kjh07GDhwoM01K1asYP/+/SxfvpwlS5aQnp7Ovffei7u7O2vWrGHdunXWcCc7OxuA8ePHM3PmTL744gvWrVvHhQsX+PHHH4scy3PPPcfcuXOZNGkS+/fvZ/r06bi7uxMaGsqCBQsAOHjwIPHx8Xz88ccF9lHcP7Nx48bxzTffsGbNGmJjYwv8MxARERERESlpP//8M2vXrsXBwYHnn3+effv28dNPP9GqVasCr7+Qls2Z5Exizqcx+pf9dJm0lkuZudZ6pwYV6NasEm89VI99H3TihwGtivVmvtuGIWUqOTnZAIzk5OR8tYyMDGPfvn1GRkaGzXmg0F9dunSxudbV1bXQa9u1a2dzrb+/f4HXXa9Zs2YZXl5e+c7HxMQYgDFjxgzrub179xqAsX//fsMwDKNXr15Gv379bNqtXbvWsLOzMzIyMoxDhw4ZgLF+/Xpr/fz584aLi4sxf/586/0BY8eOHdZrjhw5YphMJuPUqVM2fd93333G66+/XuS4q1SpYkyYMMEwDMP47bffDDs7O+PgwYPF/jy6dOlijBgxwnrcrl07Y8iQIYVe/+677xpms9k4efKk9dwvv/xi2NnZGfHx8YZhGEbv3r2NChUqGFlZWdZrvvjiC6NOnTqGxWKxnsvKyjJcXFyM3377zTAMwwgODjbGjBljrefk5BghISFG165dCxzfwYMHDcBYvnx5gWNdtWqVARgXL160Of/XPq7lz+zIkSPWa6ZMmWJUqFCh0M+psL8nIiIiIiIiRUlNTTUmTpxobNq0yXru+PHjxogRI2y+DytMdm6e0XrMCqNKxJICf7WLWmlk5uSW5hTKraKyjr/SnlJSZho3bmz9fXBwMAAJCQnUrVuXrVu3cuTIEWbPnm29xjAMLBYLMTExHD58GHt7e1q0aGGt+/n5UadOHfbv32895+joaHOfbdu2YRiG9TG3K7KysvDz8yv22Hfs2EFISEi+fq7Iy8tjzJgxzJs3j1OnTpGVlUVWVtY1b5ReuXJlQkJCrMctW7bEYrFw8OBBgoKCAGjUqJHNPlJXPjsPDw+bvjIzMzl69CjJycnEx8fTsmVLa83e3p7mzZvne4Tvr/M1m820a9fumsb/V/v37y/Wn5mrqys1atSwHgcHB5OQkHDd9xUREREREfmrc+fOMXnyZCZPnsyFCxd45JFHrI/qValSxbqZeVFmrY9h9sZY4i5m5Ks5mu34fXg7Kng54WRvLvHx30oUSt2EUlNTC61deXzsiqK+mf/7ht9/3XPoRvjrpnBXli9aLBbrP19++WUGDx6cr13lypU5dOhQgX0a/7fp3BUuLi42xxaLBbPZzNatW/N9Vu7u7sUeu4uLS5H18ePHM2HCBCZOnEijRo1wc3Nj6NCh1sfnrteVufx1Tn8PuiwWC2FhYTaB3hXXu1n71eZbHIUFXn//M/v7ZoEmk6nQtiIiIiIiIsUVExNj3crkynY2NWvW5OGHH873fUlBPlyyj+2xF+nWLIT3F++zqVX1c2Vuv5acSkoHoLKfa+lM4hajUOomdC2rbUrr2tLWrFkz9u7dS82aBb+JoH79+uTm5rJx40brs72JiYkcOnSIevXqFdrvHXfcQV5eHgkJCbRp06bAaxwdHcnLyytyfI0bNyYuLo5Dhw4VuFpq7dq1dO3alWeffRa4HBQdPny4yLEVJDY2ltOnT1OxYkXg8hsg7OzsCl2hBZc/u3nz5hEYGIinp2eB1wQHB/Pnn3/Stm1bAHJzc9m6dSvNmjUr8PpGjRphsVhYvXo1HTt2zFe/slKrqM/tev/MRERERERE/qnw8HDGjx9vXQjRvHlzIiIieOyxx/ItWPi77zbGsmjnKf48dgGAbbFJNvX/9W/JHaHe2JvtCPJyLpXx36q00bmUSxEREURHRzNw4EB27NjB4cOHWbRoEYMGDQKgVq1adO3alb59+7Ju3Tp27tzJs88+S6VKlejatWuh/dauXZuePXvy3HPP8cMPPxATE8PmzZuJjIxk6dKlwOW37KWmprJixQrOnz9Penp6vn7atWtH27Zt6d69O8uXLycmJoZffvmFX3/9Fbicti9fvpwNGzawf/9+Xn75Zc6cOXPNn4OzszO9e/dm586drF27lsGDB/Pkk09aH90rSM+ePfH396dr166sXbuWmJgYVq9ezZAhQ4iLiwNgyJAhjBkzhh9//JEDBw4wYMAAm7f6/V3VqlXp3bs3L7zwAj/99BMxMTH88ccfzJ8/H7i8xNVkMrFkyRLOnTtX4Gq+6/0zExERERERuVZXtn+5okqVKlgsFh544AFWrFjBpk2bePzxx4sMpFIyc0hIyeSNH3dbA6krKno5M6pTHSI616V5FR/szYpXroc+NSmXGjduzOrVqzl8+DBt2rThjjvu4O2337buPQUwa9YswsLCePjhh2nZsiWGYbB06dJ8j3/93axZs3juuecYMWIEderU4ZFHHmHjxo2EhoYC0KpVK/r370+PHj0ICAggKiqqwH4WLFjAnXfeydNPP039+vUJDw+3rhR6++23adasGZ06daJ9+/YEBQXx6KOPXvPnULNmTbp160aXLl144IEHaNiwIVOnTi2yjaurK2vWrKFy5cp069aNevXq8cILL5CRkWFdOTVixAiee+45+vTpQ8uWLfHw8OCxxx4rst9p06bx+OOPM2DAAOrWrUvfvn2tb2esVKkS77//Pq+99hoVKlSweUviX13vn5mIiIiIiEhx5OXl8f3339OiRQu+/vpr6/nnn3+ebdu28dtvv9GhQ4dCH9X7dU88ry3YxQtfbqbxe8u4678rbOqTn7mDAx92Zl1EBwbeW5NX2tfQ2/T+AZNRxpu1TJ06lbFjxxIfH0+DBg2YOHFioY9VAaxevZrhw4ezd+9eKlasSHh4OP3797e5ZsGCBbz99tscPXqUGjVq8J///CffN9xF3TcnJ4e33nqLpUuXcuzYMby8vOjYsSNjxoyxPkYFlzfHHjlyJHPmzCEjI4P77ruPqVOn2mxMfTUpKSl4eXmRnJyc71GrzMxMYmJiqFatGs7OWgJ4u3nvvff46aef2LFjR1kPpVzT3xMREREREcnMzOSrr75i3LhxHDlyBICwsDC2bNlSrPY5eRYupmXTasxKci35Y5LGIV60runPiPtra1VUMRSVdfxVmX6S8+bNY+jQobz55pts376dNm3a8OCDDxIbG1vg9TExMXTp0oU2bdqwfft23njjDQYPHsyCBQus10RHR9OjRw969erFzp076dWrF08++SQbN24s9n3T09PZtm0bb7/9Ntu2beOHH37g0KFDPPLIIzbjGTp0KD/++CNz585l3bp1pKam8vDDD191PyIRERERERER+eeSkpIYPXo0VatWpX///hw5cgQfHx/efvttfvnllyLbnryQzv+2nCQzJ4++X2/hrv+usAmk+rapxvOtq3J3dV++fakFEZ3rKpAqYWW6UqpFixY0a9aMadOmWc/Vq1ePRx99lNGjR+e7PiIigkWLFtm8Pr5///7s3LmT6OhoAHr06EFKSorNv3ydO3fGx8eHOXPmXNd9ATZv3sxdd93FiRMnqFy5MsnJyQQEBPDNN9/Qo0cPAE6fPk1oaChLly6lU6dOxfoMtFJKCqOVUsWjvyciIiIiIrevRx55hMWLFwMQGhrKiBEjePHFF4t8u3p2roVzqVkMmbOdLScu5qtX8naha9OKDL6vFs4ORW+CLgUr9yulsrOz2bp1Kw888IDN+QceeIANGzYU2CY6Ojrf9Z06dWLLli3k5OQUec2VPq/nvgDJycmYTCa8vb0B2Lp1Kzk5OTb9VKxYkYYNGxbZT1ZWFikpKTa/RAry3nvvKZASERERERH5iwMHDnDu3Dnr8aBBg2jYsCFff/01R48eZciQIUUGUgDvLtpD6zErCwyk7O1MfNorjPDOdRVI3QD2ZXXj8+fPk5eXR4UKFWzOV6hQodC3lJ05c6bA63Nzczl//jzBwcGFXnOlz+u5b2ZmJq+99hrPPPOMNeE7c+YMjo6O+Pj4FLsfgNGjR/P+++8XWhcRERERERERW9HR0URGRrJw4UJef/11/vvf/wLQsWNHdu3aVehm43kWg7MpmZxKymDW+hi8XByYs+mkte7maObZllVoXcOfVjX8uJCeTaCHnsC4UcoslLri7//iGIZR5M71BV3/9/PF6bO4983JyeGpp57CYrFc9a1nxRn/66+/zvDhw63HKSkp1re+FdWniBRMfz9ERERERG5NV97WHRkZydq1a63nT58+bf391d58N/H3Q3yy8kih9cnPNOPeuoHWYwVSN1aZhVL+/v6YzeZ8q4oSEhLyrWK6IigoqMDr7e3t8fPzK/KaK31ey31zcnJ48skniYmJYeXKlTbPQQYFBZGdnc3FixdtVkslJCTQqlWrQuft5OSEk5NTofW/cnBwAC5vvO7i4lKsNiK3m/T0dOD//30REREREZGb37x58/j3v//Nnj17gMtf7/fq1YtRo0ZRt27dQttdWSgScz6NH7bF5QukAj2cSLiUxUdPNqFpqDfVA4p+1E9KV5mFUo6OjoSFhbF8+XIee+wx6/nly5fTtWvXAtu0bNnSuoHZFcuWLaN58+bWb0hbtmzJ8uXLGTZsmM01V4Ki4t73SiB1+PBhVq1aZQ29rggLC8PBwYHly5fz5JNPAhAfH8+ePXuIioq6no8kH7PZjLe3NwkJCQC4urpeNQUWuV0YhkF6ejoJCQl4e3tjNut5bxERERGRW8Xq1avZs2cPHh4evPzyywwdOpRKlSoV2Wbp7nhe/2E3T90VysLtpzmTkmlTH/t4Yx5uXJHYC+nUCfIozeFLMZXp43vDhw+nV69eNG/enJYtW/LZZ58RGxtL//79gcuPup06dYqvv/4auPymvcmTJzN8+HD69u1LdHQ0X3zxhfWtegBDhgyhbdu2REZG0rVrVxYuXMjvv//OunXrin3f3NxcHn/8cbZt28aSJUvIy8uzrqzy9fXF0dERLy8vXnzxRUaMGIGfnx++vr6MHDmSRo0a0bFjxxL7jIKCggCswZSI2PL29rb+PRERERERkZvPuXPn+OSTT+jatSthYWEAjBw5ksqVK9O/f3/rC8cKE3cxnU9XH+PbjScwDPh09TGb+tN3hTKgfU1CfV0BFEiVI2UaSvXo0YPExEQ++OAD4uPjadiwIUuXLqVKlSrA5ZVHsbGx1uurVavG0qVLGTZsGFOmTKFixYpMmjSJ7t27W69p1aoVc+fO5a233uLtt9+mRo0azJs3jxYtWhT7vnFxcSxatAiApk2b2ox51apVtG/fHoAJEyZgb2/Pk08+SUZGBvfddx9ffvllia7YMJlMBAcHExgYaH3DoIhc5uDgoBVSIiIiIiI3qWPHjjF+/HhmzpxJZmYm+/bt4/vvvwegevXqvPbaa4W2zbMY/HEwgQYVvXh/8T6W7ztrU29Z3Y/nWlZhxYEEhtxXmyAv7RVVHpkM7RJcplJSUvDy8iI5OdlmzyoRERERERGRW9GOHTuIjIxk/vz5WCwWAJo3b87rr79Ot27dimx7PjWLCcsPcTYlk9/353+i6O7qvkzrGYa3q4O2vylDxc06yvzteyIiIiIiIiJye3jxxReZOXOm9bhTp05ERETQvn37q4ZIFovBJysOM3tjbL6ah5M9d1bzZcQDtfFxcyzxcUvpUCglIiIiIiIiIqUiLy8PwLrtRqNGjTCbzfTo0YNRo0bl2zLnr7JzLXwdfZwOdQNZfegcUb8eJCMnz1qv5O3C2Mcbs3hXPP3aVqeav1upzkVKnh7fK2N6fE9ERERERERuNZmZmXz11VeMGzeO9957j549ewKQlpZGQkIC1apVK7p9Th7f/nmCf/+8v8B6x3qBvNK+BmFVfEt87PLP6fE9EREREREREbmhkpKSmDZtGh9//DFnz17efHzGjBnWUMrNze2qgdTWExd4+rONZOdZrOdMJmhfO4BVB8/Rv10NXnuwbulNQm4YhVIiIiIiIiIi8o+cOnWKCRMm8Omnn5KamgpAaGgoI0aM4MUXXyyybXp2LqeTMuk9cxPVA9zYHptkE0g92DCIJ5uHcm/dQOKTMwj00Jv0bhUKpURERERERETkH+nduzcrVqwAoGHDhoSHh/PUU0/h4OBQZLvDZy/x8CfryMq9HEKdSsqwqT/YMIhpz4ZZj4O9XEp45FKWFEqJiIiIiIiIyDWJjo6mVq1a+Pv7AzB8+HBycnKIiIjgwQcfLPRNele2tc7KtbBwxyl+2HbKGkgBhPq6cDY5i+daVqFL42BqBbqX/mSkzGij8zKmjc5FRERERETkZmCxWFi6dCmRkZGsW7eOd999l/feew+4HDYVFkRdce5SFl0nr8PNyZ7UrFzikzNt6ndW9WH+yy0xmUzF6k/KL210LiIiIiIiIiL/WE5ODnPmzCEqKoq9e/cC4ODgQFpamvWaogIkwzCYv+Ukv+09y+m/BVEAPq4OzOxzJ9X83az9KJC6PSiUEhEREREREZECTZ06lTFjxnDy5EkAPDw86N+/P0OHDqVixYqFtjMMg9kbY3F1NOPsYCZiwW6beqivC189fxff/hlL65p+3FHZp1TnIeWTQikRERERERERKdD27ds5efIkFSpUYOjQofTv3x9vb+8i2/xvy0n2nErmq+gTBdYXv3oPtSq44+xg5p1/1S+FUcvNQqGUiIiIiIiIiHDs2DHGjx9P3759adq0KQCjRo3izjvv5LnnnsPZ2fmqfew5lcyo73flO29ngmaVfWhTK4BGIV4lPXS5SSmUEhEREREREbmNbd++ncjISP73v/9hsVi4ePEi3333HQC1a9emdu3ahbZdeeAsuXkGJxLTmb/lJIcTUm3qs56/k9TMXOoFe1JTb9KTv1EoJSIiIiIiInKbMQyDlStXEhkZyfLly63nO3XqRL9+/YrVx5nkTF76agsWI3/N0d6Oni0qc2+dwJIastyCFEqJiIiIiIiI3Ga6du3K4sWLATCbzfTo0YPw8HCaNGlSaJtLmTlk51rwcHag98xNRB9LtKnfWdWHA2cuUcHTmV+HtMHebFeqc5Cbn0IpERERERERkVtcRkYGjo6OmM1mAO655x5+//13XnzxRUaMGEHVqlWLbJ+bZ+GRyeuJT86gopcLx86nWWtmOxO9W1blzYfqkZWbh9nOpEBKisVkGEYBC+3kRklJScHLy4vk5GQ8PT3LejgiIiIiIiJyC7l48SJTp05l0qRJTJ48mSeeeAKA1NRUMjIyCAgIKLK9YRi89NUWVhxIKLDuZG/H+tc64O/uVOJjl5tXcbMOrZQSERERERERucXExcUxYcIEPvvsM1JTL28+Pnv2bGso5e7ujrt74RuPx11MZ1PMBZLSc/IFUndX96VJiDePNK2IYaBASq6bQikRERERERGRW8S+ffsYO3Yss2fPJicnB4BGjRoRERHBk08+edX2r/+wm2V7z5CYlp2vZmeCZcPaUjPQo8THLbcnhVIiIiIiIiIit4i+ffuyYcMGANq1a0dERASdO3fGZDIV2e50UgYbjiYyZ1NsgfV3/1Wf+sGeCqSkRCmUEhEREREREbkJWSwWfv75Z+655x58fHwAiIiI4MsvvyQiIoIWLVoU2X7RztP8vu8s/dpWp/fMTflWR335/J1k5lg4m5LJcy2rXDXYErlW2ui8jGmjcxEREREREbkW2dnZzJkzh6ioKPbt28e///1v3nzzzWK3P34+jT+PJfLuor1k5Vry1V+6pxqdGgZxZ1Xfkhy23Ea00bmIiIiIiIjILeTSpUt8/vnnTJgwgbi4OAA8PT2xty/et/ZX1qT0/3YrB85csqn5uzvSoKIXqVm5DOpQCy9Xh5IdvEgBFEqJiIiIiIiIlHMffPABEyZMICkpCYCgoCCGDRvGyy+/jJeXV6Ht4pMzWHf4PEcSUpm1/jjZebYro7o3C+HxsBDuqOyNs4O5NKcgko9CKREREREREZFy7ujRoyQlJVG7dm1GjRpFr169cHJyumq78O93sfbw+QJrVf1cGXZ/LUJ8XEt6uCLFoj2lypj2lBIREREREZG/2rZtG1FRUbz11ls0bNgQgIMHD7J37166du2K2VzwiqbMnDwOnLlEkxAvvt8ax5RVRziemG6tV/Z15Y7K3uw8mcTMPndSPcD9hsxHbj/aU0pERERERETkJmEYBitWrCAyMpLff/8dACcnJ7766isA6tSpQ506dYrsY8wvB/hyw/ECa9UD3Jj+bBi1K3iU6LhF/gmFUiIiIiIiIiJlJDc3lwULFhAVFcW2bdsAMJvNPPXUU4wYMaJYfSzbe4b3F+/jVFJGgfUPH21Ir7urlNiYRUqKQikRERERERGRMmAYBu3atWPDhg0AuLi48NJLLzF8+HCqVq1aZNvj59OY9sdRWtbw4/UfdpORk2et1Q/2pOfdlencIIjoY4k82DC4NKchct0USomIiIiIiIjcIBcvXsTLyws7OztMJhNdunThwIEDDBo0iFdffRV/f/8i2284ep7Xf9jNif/bK2relpM29f7tavDag3Wtxw83rljykxApIdrovIxpo3MREREREZFbX1xcHBMmTOCzzz7jm2++4dFHHwUgNTUVk8mEm5tbke1z8yzsjEvmzR93c+DMpXz1J8JCqBfsSY87Q3Fz0voTKVva6FxERERERESkjO3bt4+xY8cye/ZscnJyAPjpp5+soZS7e+FvwIs5n8ZbP+2me7MQvt8ax4ajiTb1x+6oxDsP12f7yYu0rRWAvdmu1OYhUhoUSomIiIiIiIiUsPXr1xMZGcnixYut59q3b09ERASdOnUqsq3FYvDHoQTmbT7J+iOJrD9iG0YFejjxaoeadG1aCS8XBzrUrVAqcxApbQqlREREREREREqQYRgMGTKErVu3YjKZeOyxxwgPD6dFixZFtttzKplQX1f+t+Uk//55v00t2MuZl9pUZ+6mWN7oUo976waW5hREbgiFUiIiIiIiIiL/QHZ2NnPmzOHRRx/Fy8sLk8nEm2++ydKlSxk5ciR16tQptK1hGKw8kMDJC+m8t3gfbo5m0rLzbK755Ok7aFXDDz93J168p1ppT0fkhtFG52VMG52LiIiIiIjcnC5dusTnn3/OhAkTiIuLIzIykvDw8GvqY+GOUwyZu6PAmoPZxGsP1lMQJTcdbXQuIiIiIiIiUgrOnj3LJ598wpQpU0hKSgIgKCgIb2/vq7bdHZdMZV9X9sYn8030CX7Zc8amHtG5LjtOXuS5llVpWd0POztTKcxApHxQKCUiIiIiIiJSDIZhMHjwYGbMmEFmZiYAtWvXZtSoUfTq1QsnJ6ci2285foHHp0cXWh/UoSavtK9RomMWKc8USomIiIiIiIgUg8lkIjExkczMTFq0aEFERASPPPIIZrO5yHbJGTmM+t9Olu07W2D97uq+zOl7NyaTVkXJ7UV7SpUx7SklIiIiIiJS/hiGwe+//87YsWP55JNPrJuVHzx4kDNnztC2bdsiQ6SUzBwGfLuNIC9nLBaDH7afsqk/17IK7/6rAWsOn6NBRU8CPZxLdT4iN1Jxsw6FUmVMoZSIiIiIiEj5kZuby4IFC4iKimLbtm0AvPTSS3z++efFap+Zk8fgOdtZdTCBnLz8325X8HRiw2v3YdZeUXIL00bnIiIiIiIiIsWUkZHBrFmzGD9+PMeOHQPA1dWVl156iWHDhl21vWEYxF3MYOWBhAIf02tTy597avoTVsVHgZTI/1EoJSIiIiIiIrc1wzAICwtj//79APj5+TFo0CAGDhyIv79/oe0sFoO3F+4hz2Kw53Qye06l2NRdHMxsevPyqigXB7P2jBL5G4VSIiIiIiIicts5deoUFStWxGQyYTKZeOKJJ/jqq68YMWIEL7zwAm5ubkW2P5OcycaYRGZvjC2w/nzrqrSvE4iHs0NpDF/klqA9pcqY9pQSERERERG5cfbu3UtUVBTfffcdCxcupEuXLgCkp6fj4OCAg0PhIdKaQ+fIzrVQJ8iDBz9eS2pWrrVmZ4Lpz4Yxb/NJ6gZ7MKpT3VKfi0h5pT2lRERERERERP7PunXriIyMZMmSJdZzK1assIZSrq6uhbY9nZTBqaQMXvhyM7mW/Os6ut1RiYcaB3NfvQo80CCo5AcvcotSKCUiIiIiIiK3JMMwWLx4MZGRkWzYsAEAk8nEY489RkREBHfddddV+0jLyuVfn6wjMS27wHrdIA/GPtFEm5eLXAeFUiIiIiIiInLLevvtt9m1axeOjo707t2bESNGUKdOnUKvP5+aRXpWHuuPnmfdkfPsO51iE0h1bVqRBhU9ubdOIKG+rjiY7RRIiVwnhVIiIiIiIiJyS0hJSeGLL76gb9++uLu7YzKZeOedd9i8eTNDhgwhODi4yPZ5FoMnpkcTcz6twHqwlzP/fawRbk76VlqkJOhvkoiIiIiIiNzUzp49y8cff8zUqVNJTk7GZDIxdOhQALp370737t0LbWuxGNjZmfhp+ym++fOETSDl+H+roEwm+HVIWwI9nXB2MJf2dERuGwqlRERERERE5KZ05MgRxo0bx5dffklWVhYAderUoVKlSsVqP3TudtYePk+XRsF88+cJm5q/uyNLBrXBwCA710Jlv8I3QheR66NQSkRERERERG4qeXl5PPvss8yfPx+LxQJAixYtiIiIoGvXrtjZ2RXZ/sftcfxvSxwbjiYC5Aukvnz+TlrW8MPJXquiREqTQikRERERERG5qZjNZnJycrBYLHTp0oWIiAjatGmDyVT4huOJqVnsPZ1CjUB3hs3baVNzcTDzWLNKNK/iQ8KlLNrVDiiyLxEpGSbDMIyyHsTtLCUlBS8vL5KTk/H09Czr4YiIiIiIiJQrubm5fP/994wfP565c+dSo0YNAA4ePEh2djaNGjUqsv3yfWf5eddpNsVc4HRyZr762Mcb071ZCHZ6g55IiSlu1qGVUiIiIiIiIlLuZGRkMGvWLMaNG0dMTAwAEydO5JNPPgEu7x1VFMMwSM3KZdCcbWTmWPLVm1fxoXGIN12bVlIgJVJGFEqJiIiIiIhIuXHhwgWmTJnCJ598wrlz5wDw9/dn0KBBDBw4sMi2h89e4pc9Z3j6rspELNjFygMJNvVBHWpSyduFY+fTGH5/bb1JT6SMKZQSERERERGRciEvL48mTZoQFxcHQNWqVRkxYgQvvPACrq6Fv/0uz2KQcCmT13/YzZYTF/lo+SGbegVPJx5qVJFX2tfA1VHfBouUF/rbKCIiIiIiImXm4MGD1K5dG5PJhNlsplevXixdupSIiAieeOIJ7O2v/m3rf37ez8z1MTbnXBzM2NuZSM/J49NezWka6l1KMxCR66WNzsuYNjoXEREREZHbjWEYrFu3jsjISH7++WeWLVvG/fffD0BWVhaOjo6Fvv3OYjE4lZTBppgLjP5lP/fXr8CcTSetdWcHO15oXY2ed1ehopczl7Jy8XR2uCHzEpHLtNG5iIiIiIiIlCsWi4XFixcTGRlJdHQ0ACaTiU2bNllDKScnpyL7mPrHEcYt+/+P5/01kAKY8GRTHmwUbD1WICVSfimUEhERERERkVKVl5fHV199xdixYzlw4ABwOXzq3bs3I0eOpFatWoW2zcrNw8nezNFzqSzacZqPVxy21lwczDSs5MmOk0l81qs5ob4u1Az0KPX5iEjJUCglIiIiIiIipcpkMjF+/HgOHDiAl5cXr7zyCkOGDCEoKKjIdot3nmbw3O00quTFkYRU0rPzbOpvPFSPXndXIc9iYLYr+HE/ESm/FEqJiIiIiIhIiTpz5gzTp08nPDwcV1dX7Ozs+OCDDzh27Bgvv/zyVffTTUjJZOofR/nmzxMYBuyKS7apPxEWwivta1DN3w1AgZTITUqhlIiIiIiIiJSIw4cPM27cOL766iuysrIICAhg4MCBAHTv3r3ItimZOcxYc4x7agXwVfRxft4Vb1N/IiyEfm2rM3fzSfq1rU4FT+dSm4eI3BgKpUREREREROQf2bJlC5GRkSxYsIArL3i/++67qV279lXbJqVnM37ZIXafSmbHySQmrTxiU28c4sXXL9yFt6sjAG8/XL/kJyAiZUKhlIiIiIiIiFyX7OxsunTpwooVK6znHnroISIiIrjnnnswmYp+rM4wDKb932N6fxfq60LtQA9e7VDTGkiJyK1FoZSIiIiIiIgUm2EY1rDJ0dERNzc37O3tefrppwkPD6dhw4aFts3KzWP6H8e4s5oPS3bF88O2ODJzLNa6v7sj819uyZYTF+nUIAgvF4dSn4+IlB2TcWVtpZSJlJQUvLy8SE5OvupmfyIiIiIiImUlPT2dWbNmMWnSJJYtW0aVKlWAy/tIOTk5Ubly5SLb5+RZ+G5jLO8u2ltgvV3tAPq2qc49tfxLfOwicmMVN+vQSikREREREREp1IULF5gyZQqTJk3i/PnzAEydOpXIyEgAatWqVWjbuIvp+Ls7cfhsKk9+Gk1GTp61ZjJBl0bBLNt7hn5tqzOqU93SnYiIlDsKpURERERERCSf2NhYPvroI2bMmEFaWhoA1apVY+TIkfTp06fItpk5eRw4c4nu0zZQ0duZkxcybOod6wXy7N1VaF8nkIzsPFwczaU1DREpxxRKiYiIiIiIiI2srCzuuOMOLly4AEDTpk2JiIjg8ccfx96+6G8jTySm0eXjtaRlX14V9fdA6t46Aczofaf1WIGUyO1LoZSIiIiIiMhtzjAMtm7dSlhYGCaTCScnJ1544QW2bdtGREQE999/f6Fv0rNYLm9TfPRcKn8eS+SXPWesgRRA/WBPvF0d6N4shNoVPKjs53pD5iQi5Z9CKRERERERkduUxWJh8eLFREZGEh0dzapVq2jfvj0AY8aMwWwuehXThbRsHvx4DWdTsgqsN6joyZJB9xQaaInI7U2hlIiIiIiIyG0mKyuL2bNnM3bsWA4cOACAk5MT+/bts4ZSRQVShmGwcMdpftkTX2Ag5evmyOfPNaeKn6sCKREplEIpERERERGR20RWVhaffPIJEyZM4PTp0wB4e3szYMAABg0aRFBQUKFtLRaDdxbtwcneTOMQL4bO22FTr+TtwpJB9/Dj9lM0CvEirIpPaU5FRG4BCqVERERERERuE/b29nz22WecPn2aSpUqMWzYMPr164eHh0eR7RbtPE300UTmbIotsP7dSy2oHeSBj5sjL9xTrTSGLiK3IIVSIiIiIiIit6jDhw8zbdo0/vvf/+Ls7IzZbGbMmDEkJyfTs2dPHB0dr97H2UsMnrM93/m6QR5U8nahQUVPWtX0L43hi8gtTqGUiIiIiIjILWbz5s1ERUWxYMECDMOgXr169O3bF4Bu3boV2XbVgQTOp2bx/dY4dsUlk5GTZ1Of/VILKng6U83fDbOd9osSkeunUEpEREREROQWYBgGy5YtIzIyklWrVlnPP/TQQzRp0qRYfSSmZtH36y3kWox8NSd7Ox67oxKtavhp83IRKREKpURERERERG5y6enp3HPPPWzffvkxO3t7e5555hlGjRpFw4YNC213KimDpPRsaga6M3TuDn7Zc8am3r1ZCAfOpODmZM+cvndrZZSIlCiFUiIiIiIiIjehvLw8zGYzAK6urgQFBeHm5kbfvn0ZNmwYlStXLrJ9bp6FJ6dHcyopo8B6t2aViOzeCHuzXYmPXUQEFEqJiIiIiIjcVBITE5kyZQqffvopGzduJCQkBIDJkyfj7e2Nr6/vVfsYPGc7i3aeLrDmaG/H6lHtCfZyKdFxi4j8nUIpERERERGRm8CJEyf46KOPmDFjBunp6QDMnDmTd955B4Dq1asX2X75vrOsP3KeCp7O+QKpBxsG8VzLqgR5OZOZk6dASkRuCIVSIiIiIiIi5dju3buJiopizpw55OVdfhNe06ZNiYiI4PHHH79q+3cX7mHO5pNk51oKrC8ZdA8NK3mV6JhFRIpDoZSIiIiIiEg5lZqaSsuWLUlLSwOgQ4cOREREcP/991/1DXjxyRlsirnAV9EnbM6b7UxYDIMx3RoR4uOqQEpEyoxCKRERERERkXLCYrGwevVq7r33XgDc3d3p168fJ0+eJDw8nDvvvLPI9nM3xfLDtlP0vLsy7y3ay8X0HJv6zD7NqernxsX0bMKqXH3vKRGR0mQyDMMo60HczlJSUvDy8iI5ORlPT8+yHo6IiIiIiJSBrKwsvv32W8aOHcvBgwdZt24drVu3BsAwjKuuijqVlMH6I+cJ/35XgfXnWlbhgfpB3FPLv8THLiLyd8XNOrRSSkREREREpIykpKTw6aefMnHiRE6fvrz5uJeXF8ePH7eGUoUFUpk5eRw4c4mGFT0Z8O1WdsYl29Tb1wmgbpAnO08mMfz+2ni7OpbuZERErpFCKRERERERkRssLS2NDz/8kGnTppGSkgJApUqVGDZsGH379i1yZUFCSiYrDySw6mACv+09SyVvF04lZVjrjzSpSM8Wlbmrmu9VV1iJiJQlhVIiIiIiIiI3mJOTE/PnzyclJYW6desSHh5Oz549cXS8+mqmN3/aw/J9Z63Hfw2kQnxcGHZ/bar5u5XKuEVESpJCKRERERERkVK2adMmPv/8c6ZMmYKjoyP29vZ89NFHmEwm/vWvf2FnZ1dgu4zsPKKPnad1TX++iT7B3M0nOZKQaq2HVfGhfe0ADpy5xLuP1CfQw/lGTUlE5B9TKCUiIiIiIlIKDMPgt99+IzIykj/++AOA1q1b06dPHwAeffTRq/Yx9reDzFwfU2AtxMeF9/7VgEYhXiU0YhGRG0uhlIiIiIiISAnKzc1l/vz5REVFsXPnTgDs7e155plnaNGiRbH6WHv4HO8s3EvM+bQC6x90bcBzLauW1JBFRMqEQikREREREZEScvHiRZo1a8bx48cBcHNzo1+/fgwdOpTKlSsX2XbriYu8+eNuOjcM4tPVx8jIybPW6gV7MuL+2rSq6ceW4xe5p6Z/aU5DROSGUCglIiIiIiLyD2RmZuLsfHkvJx8fH2rUqEFaWhqDBw9mwIAB+Pr6Ftl+64mLDJ+/gxOJ6QAcOHPJpt6nVVXee6SB9bht7YASnoGISNlQKCUiIiIiInIdTpw4wfjx45k9ezZ79+4lKCgIgJkzZxIQEICLi0uR7bNy89gdl8y/f95vDaQAzHYm8iwGQzvWwsvFge5hIaU6DxGRsqJQSkRERERE5Brs2rWLqKgo5s6dS17e5Ufs5s6dy9ChQwGKfEzv+Pk0hszbQcOKnvy29yznU7Ns6g81CmbyM3eQkpGLl6tDqc1BRKQ8UCglIiIiIiJyFYZhsGbNGiIjI/nll1+s5++77z4iIiLo2LHjVdtvOJrIdxtj2XkyiZ0nk2zqFb2c6de2Og83qYjJZFIgJSK3BYVSIiIiIiIiV3HhwgU6d+5MZmYmdnZ2PP7444SHhxMWFlZom9w8C7M3xnJXNV82H7/AOwv32tRrBrrz7r/qM331UV5pV5N7amnzchG5vdiV9QCmTp1KtWrVcHZ2JiwsjLVr1xZ5/erVqwkLC8PZ2Znq1aszffr0fNcsWLCA+vXr4+TkRP369fnxxx+v+b4//PADnTp1wt/fH5PJxI4dO/L10b59e0wmk82vp5566to+ABERERERKXeysrJYtGiR9djPz48BAwbQv39/Dh48yLx584oMpFYfOseUVUd5d9FeHvx4bb5A6uOnmjL/5Za0qRXA7JfuViAlIrelMg2l5s2bx9ChQ3nzzTfZvn07bdq04cEHHyQ2NrbA62NiYujSpQtt2rRh+/btvPHGGwwePJgFCxZYr4mOjqZHjx706tWLnTt30qtXL5588kk2btx4TfdNS0ujdevWjBkzpsg59O3bl/j4eOuvTz/99B9+KiIiIiIiUlaSk5OJjIykatWqdO3alc2bN1tr48ePZ9q0adSsWbPIPpbvO0vvmZuY8PuhfDU3RzOjOtWha9NK+Lo5lvj4RURuJibDMIyyunmLFi1o1qwZ06ZNs56rV68ejz76KKNHj853fUREBIsWLWL//v3Wc/3792fnzp1ER0cD0KNHD1JSUmye8+7cuTM+Pj7MmTPnmu97/PhxqlWrxvbt22natKlNrX379jRt2pSJEyde92eQkpKCl5cXycnJeHp6Xnc/IiIiIiJy/eLj45k4cSLTp08nJSUFgEqVKjF9+nQefvjhIttui71IBU9nftwWx664ZJbtO2tT/89jDcnKsfBAgwpU8nbBZDKV2jxERMqD4mYdZbZSKjs7m61bt/LAAw/YnH/ggQfYsGFDgW2io6PzXd+pUye2bNlCTk5Okddc6fN67luU2bNn4+/vT4MGDRg5ciSXLl0q8vqsrCxSUlJsfomIiIiISNm4ePEiffv2pWrVqkRFRZGSkkK9evWYNWsWx44du2ogteNkEt2mbqD1mJWMW3YoXyD1fOuq9GxRhRfuqUaIj6sCKRGRvyizjc7Pnz9PXl4eFSpUsDlfoUIFzpw5U2CbM2fOFHh9bm4u58+fJzg4uNBrrvR5PfctTM+ePalWrRpBQUHs2bOH119/nZ07d7J8+fJC24wePZr333//mu4jIiIiIiKlw93dnd9++43s7GxatWpFREQEDz/8MHZ2Bf/8PjfPgp3JRHpOHq//sJvFO0/b1Ct5u3AqKYPWNf345oUW2NkphBIRKUyZv33v7z8pMAyjyJ8eFHT9388Xp89rvW9B+vbta/19w4YNqVWrFs2bN2fbtm00a9aswDavv/46w4cPtx6npKQQGhp6TfcVEREREZFrZxgGv/32G19//TVfffUVDg4OODg4MGXKFHx9fWndunWR7S+kZdN54hp8XB3xcnFg0/ELNvVud1Ri/JNNOJKQSpCXswIpEZGrKLNQyt/fH7PZnG91UkJCQr5VTFcEBQUVeL29vT1+fn5FXnOlz+u5b3E1a9YMBwcHDh8+XGgo5eTkhJOT0z+6j4iIiIiIFF9OTg7z588nKiqKXbt2AfDQQw/Rs2dPAP71r38V2T4718LI/+1k8a7TGAYkXMqyqfu6ObLxjftwMF9eXVWrgkcpzEJE5NZTZntKOTo6EhYWlu9Rt+XLl9OqVasC27Rs2TLf9cuWLaN58+Y4ODgUec2VPq/nvsW1d+9ecnJyCA4O/kf9iIiIiIjIP5eWlsakSZOoVasWzz77LLt27cLd3Z3hw4fTrl27q7aPOZ/Gn8cS+Wn7KRbtvBxIXVHVz5V3/1Wf9x9pwPRnw6yBlIiIFF+ZPr43fPhwevXqRfPmzWnZsiWfffYZsbGx9O/fH7j8qNupU6f4+uuvgctv2ps8eTLDhw+nb9++REdH88UXX1jfqgcwZMgQ2rZtS2RkJF27dmXhwoX8/vvvrFu3rtj3Bbhw4QKxsbGcPn35GfGDBw8Cl1diBQUFcfToUWbPnk2XLl3w9/dn3759jBgxgjvuuOOqy35FRERERKR0nTlzhoYNG5KYmAhAQEAAQ4YMYcCAAfj4+BTazjAM3vppDycvZvDn0USy8yw2dbOdiW1v3Y+Xq0Opjl9E5HZQpqFUjx49SExM5IMPPiA+Pp6GDRuydOlSqlSpAlx+LWtsbKz1+mrVqrF06VKGDRvGlClTqFixIpMmTaJ79+7Wa1q1asXcuXN56623ePvtt6lRowbz5s2jRYsWxb4vwKJFi3j++eetx0899RQA7777Lu+99x6Ojo6sWLGCjz/+mNTUVEJDQ3nooYd49913MZvNpfaZiYiIiIhIwVJSUqyvHg8KCqJBgwbExcUxcuRI+vTpg4uLS9HtM3PYevwiszfGFlh/8Z5qtKjmq0BKRKSEmAzjr4tQ5UZLSUnBy8uL5ORk6/9ARURERESk+Hbt2kVUVBSLFy/m6NGj+Pv7A5dXSwUEBBT5Q+OFO05xPjWbdrX9eXx6NEnpOdaayQTf92/Jb3vPEurjQq+WVUt7KiIit4TiZh1l/vY9ERERERGRa2UYBqtXryYyMpJff/3Ven7JkiX06dMHuLxaqjAJlzI5kZjOsHk7sBjw4d/qDzUO5qFGwYRV8SWsim8pzEBERBRKiYiIiIjITSMvL4+FCxcSGRnJpk2bALCzs+OJJ54gPDy80Ldg/1VmTh5dJ68nPjnT5ryHkz0mEwR7uTCxR1NtXi4iUsoUSomIiIiIyE3j3LlzPP3002RnZ+Ps7Mzzzz/PiBEjqFGjRqFtElOzSMrIYcqqI2w5fpGs3DzOpmRZ6x3rVaBzwyAebhyMvZ0Js50Jk8l0I6YjInJbUyglIiIiIiLlVnJyMr/88ov1xUNBQUEMHjwYJycnBg8eTGBgYJHtDcPgmc83cvDspQLr/u5OjH+iiTYvFxEpAwqlRERERESk3ImPj2fixIlMnz6dlJQU6tatS9OmTQEYO3Zsoe3yLAZp2bl4OjuwcMcp5myKtQmk/NwccXOyx2IY/K9/SzydHXBz0rdFIiJlQf/1FRERERGRcuPQoUOMHTuWr7/+muzsbADq169PcnJysdq/+t02VuxPoHaQO3tOpdjUvFwcWPhqa0J8XMnNs2CvPaNERMqUQikRERERESlzZ8+eZcCAAfz4448YhgFA69atiYiI4KGHHsLOrugA6dc98Xy36SRrDp0DyBdIfdYrjDa1AnBxNAMokBIRKQcUSomIiIiISJnz9vYmOjoawzB45JFHCA8Pp3Xr1kW22XMqmd/2nuGB+kH0/3abTc3DyZ4PHm2Av7sTh8+mcn/9Ctq8XESknDEZV34MIWUiJSUFLy8vkpOT8fT0LOvhiIiIiIiUupycHObPn8///vc/FixYgNl8efXSL7/8QpUqVahfv36R7f84mMD3W+NYtu8s2bmWfPUPH23I03eGajWUiEgZKW7WoZVSIiIiIiJyQ6SlpfHFF1/w0UcfceLECQB++OEHnnjiCQAefPDBItsbhkF6dh5D5u4gOSMHADsTWP7vx+z31Q2korcLT4SFKJASEbkJKJQSEREREZFSdf78eSZPnszkyZNJTEwEIDAwkCFDhtCxY8ci2+6PT2H+lpN0bhDE2N8OsuXERZv6y+1q0LlBELEX0nm4cbAe0RMRuYkolBIRERERkVITExNDgwYNyMjIAKBGjRqMHDmS3r174+LiUmg7i8XgfGoW7y7ay6aYC8xaf9ym7u/uRKcGFXilfQ08nR1oEupdirMQEZHSoFBKRERERERKVEJCAoGBgQBUrVqVJk2akJOTQ0REBN26dbPuIfV3hmFwOjmTil7OjF9+kCmrjtrUawW6E+Ljwoajiczo3ZymCqJERG5qCqVEREREROQfMwyDP/74g8jISNavX09sbCw+Pj6YTCZ+/vln6+8Laxt3MYMlu+KJ/PUAtQLdOZyQaq2bTPBKuxr0aV2VQA9ncvMs2jNKROQWoFBKRERERESuW15eHj/99BORkZFs3rwZADs7O1atWkW3bt0A8PX1LbKPmeuP8+GSfdbjvwZSAOOfaEK3ZiHWYwVSIiK3BoVSIiIiIiJyzbKysvj6668ZN24chw4dAsDZ2Znnn3+eESNGUKNGjULbJqVn4+XiwPojiWyKSWTSyiM29W7NKhF3MYP3/tUAR3s7agS4lepcRESkbCiUEhERERGRa5aQkMCAAQPIzc3Fx8eHgQMHMmjQIOteUoVZujueAbO3FVp/o0td+rUtPNASEZFbh0IpERERERG5qtOnT/Prr7/ywgsvABAaGsqIESOoUKECffv2xd3dvcj2F9Kymb76KF9HHy+w3r1ZCP3bVadmYNH9iIjIrcNkGIZR1oO4naWkpODl5UVycjKenp5lPRwRERERERsHDx5k7NixfPPNN2RnZ7N7924aNmxYrLZxF9N5/Yfd/KtxRdYcPseSXfE29T6tqjLg3hr8uO0U3cNC8Hd3Ko0piIjIDVbcrEMrpUREREREJJ8///yTyMhIFi5cyJWfY99zzz1kZ2dftW1qVi7jfjvIxpgL7I9PYe3h8zb1OhU8+K5vC/z+L4R6uZ0e1xMRuR0plBIREREREasTJ07w3HPPsWbNGuu5Rx55hIiICFq1anXV9oZh8PmaY3y54Xi+WrvaATg72PFSm+rWQEpERG5fCqVERERERMQqMDCQAwcO4ODgwLPPPsuoUaOoV69eoddn51r4ZOVhgr1cmL76KGeSM8nOs1jrro5m1oTfS0Z2HiE+LphMphsxDRERuQkolBIRERERuU2lpaUxY8YMli5dyi+//IKdnR0uLi5899131KlTh5CQkCLb51kMFmyL45OVRwqst6nlT++WVbVXlIiIFEihlIiIiIjIbebcuXNMnjyZyZMnc+HCBQCWLFnCI488AsB9991XaNtNMReo5ONCUno2T332J5cyc23qI+6vzXebYnmuZVVeaa+9okREpHAKpUREREREbhMxMTGMHz+emTNnkpGRAUCNGjUYNWoUDzzwQJFts3Mt7IpL4slPo7EzgeVv7/BuXsWHPq2r8nDjigy6r1ZpTUFERG4hCqVERERERG4De/fupUmTJuTl5QEQFhZGREQE3bp1w2w2F9k2PjmDBz9eS1J6DpA/kGpd04/ZL91dKuMWEZFbl0IpEREREZFbkGEYnDhxgqpVqwJQv3597rjjDnx8fIiIiKBDhw6Fbjqem2fBAGatj+FCWg57TydbAymAOhU8aFPLn2fvrkJWroUgL+cbMCMREbnVKJQSEREREbmF5OXl8eOPPxIZGcmhQ4eIjY3Fy8sLk8nEH3/8gZubW5HtkzNy6DxxDfHJmQXWawS48cuQNtjZ6S16IiLyzyiUEhERERG5BWRmZvL1118zbtw4Dh8+DICzszObNm3i/vvvB7hqILXywFl+3nWmwEAq2MuZKT2bUcnbRYGUiIiUCIVSIiIiIiI3seTkZKZOncrHH3/M2bNnAfDx8WHgwIEMGjSIwMDAQttaLAYDv9tGrsUgrIoPY345YFMP9HBi5cj2HIhPoYKnM6G+rqU6FxERub2YDMMwrn6ZlJaUlBS8vLxITk7G09OzrIcjIiIiIjeZmJgYatWqRV5eHqGhoQwfPpyXXnoJd3f3Itst23uGPw6d47uNsQXWv3z+TuoGeWq/KBERuWbFzTq0UkpERERE5Cayf/9+Vq5cycCBAwGoVq0ar7/+OrVq1eLpp5/GwcGh0LYJlzJJzczFbGfi5W+38tcfTzuYTfRuWZWkjBwqervQvk7hK6xERERKglZKlTGtlBIRERGR4oiOjiYyMpKFCxdiMpnYv38/derUKVbb1YfOcexcKtP+OErCpax89cnP3MFdVX0J9NSqKBER+ee0UkpERERE5CZnGAZLly4lMjKStWvXWs8/8sgjFPdny8kZOfT9egvZuZZ8NQ9ne+6rG0iXhsHavFxERG44hVIiIiIiIuXQgQMHeOKJJ9izZw8ADg4O9OrVi5EjR1KvXr1C2x1JSOXkxXSahfrw2g+7+GXPGZv68PtrcyYlExcHM292qacwSkREyoxCKRERERGRcsIwDEymyyFRaGgo8fHxeHh48PLLLzN06FAqVapUZPs8i8FzX2zkdHImfm6OJKZl29S7NAqif7saONrbldocREREikuhlIiIiIhIGTt37hyffPIJa9euZeXKlZhMJtzc3Pjxxx9p1KgR3t7eV+3j9R92MWfTSetxYlo2QZ7OZOXmkZadx4rh7Qj1dS3FWYiIiFwbhVIiIiIiImXk2LFjfPTRR8ycOZOMjAwAli9fzgMPPABAmzZtimw/7Y+j/LonnlY1/W0CKYCHGgfzzsP1MQGpWbkKpEREpNxRKCUiIiIicoNt376dqKgo5s+fj8VyeQPy5s2bExERwX333XfV9qOX7ufTNcesxzvjkm3q3/dvSfOqvtbjwBIat4iISElSKCUiIiIicgNt2rSJFi1aWI87depEREQE7du3t+4nVZi4i+lsj02yCaQAvFwcyM61MKFHE5wdzDaBlIiISHmlUEpEREREpBTl5eVx8OBB6tevD8Cdd95JWFgYtWvXJjw8nKZNmxbZfvbGE3y5/jhNQ71ZuOM02XkWm/q0ns3oUC+QzBwLXi4OpTUNERGREqdQSkRERESkFGRmZvLVV18xbtw4EhMTiY2Nxd3dHZPJxIYNG3B0dCyy/blLWaw9fI43f9wDwOGEVJv6My0q06FOIPfVC8RkMuFkby61uYiIiJQGhVIiIiIiIiUoKSmJadOm8fHHH3P27FkAfH192b17Ny1btgQoNJBKTM1i6Z4zPNQomFe+3cqWExdt6s+3rkr9YE9+33+W8E518HYtOtgSEREpz0yGYRhlPYjbWUpKCl5eXiQnJ+Pp6VnWwxERERGR65SQkMDYsWP59NNPuXTpEgCVK1dm+PDhvPjii7i7uxfaNjE1i+X7zvLj9lNsjLmQr96pQQVeaF2NO6v6YmdX9L5TIiIiZa24WYdWSomIiIiIlICkpCTGjx+PYRg0bNiQiIgIevTogYPD1fd5en/xPhbtPJ3vvMkEFb1cGH5/HeoEeZTGsEVERMqMQikRERERkesQHR1NdHQ0w4cPB6B27dq89957NG/enAcffLDQN+llZOexfP9ZmoR48dmaY2yKuWCzX1TdIA9eaV+DxNRsnmlRGWcH7RUlIiK3Jj2+V8b0+J6IiIjIzcNisbB06VIiIyNZt24ddnZ2HDp0iBo1ahS7j9FL9/PpmmMF1ip4OjHlmWY0r+pbUkMWERG54fT4noiIiIhICcnJyWHOnDlERUWxd+9eABwcHOjVq9dVH8/LzMnD2cHMppgLvP3THg6evWRTD/RwIuFSFh8+2pBed1cptTmIiIiUNwqlRERERESKsHXrVh577DFOnjwJgIeHBy+//DJDhw6lUqVKRbZdtvcML3+7lQ51AllxIMGmVruCO693qUe7WgEcOZdKrcDCN0IXERG5FSmUEhERERH5G4vFgp2dHQC1atUiJSWFChUqMHToUPr374+3t3eR7feeTmbQd9s5dj4NIF8g9fRdoYzu1th6XLuCNjEXEZHbj0IpEREREZH/c+zYMcaPH8/OnTtZu3YtJpMJT09Pli9fTqNGjXB2di6y/Z5TyaRm5fLJysPWQAogwMOJAHcnXu1Qk8TULP7VpGJpT0VERKTcUyglIiIiIre97du3ExUVxfz587FYLABs2LCB1q1bA3DnnXcW2vbkhXT6f7sVgL2nU/LV76sbyBd9Cm8vIiJyu1IoJSIiIiK3JcMwWLlyJZGRkSxfvtx6vlOnTkRERNCqVaur9rH3dDLfRJ8oMIyqHuDGsy2q8GCjoBIdt4iIyK1CoZSIiIiI3JZWrlxJx44dAbCzs6NHjx6Eh4fTtGnTQttk5uQx7reDtKrpR2xiOu8t3mdTr13BnWnPhjFj7TGebB7KHZV9SnMKIiIiNzWTYRhGWQ/idpaSkoKXlxfJycl4enqW9XBEREREblkZGRns37+fZs2aAZc3M2/RogV33303w4cPp1q1akW2jz6aSPTR80xaeaTA+oddG9C+TiChvq4lPnYREZGbSXGzDq2UEhEREZFb2sWLF5k2bRoff/wxFouFEydO4Orqip2dHRs3brS+Za8o6w6f59kvNuY7H+jhhMWA51pWoVfLqqUwehERkVuXQikRERERuSWdOnWKCRMm8Omnn5KamgpA5cqVOXLkCI0bNwYoNJDaFnsRezsTb/20h5w8g/3xtntGfdC1AXWDPGlexQc7O1PpTkREROQWpVBKRERERG4pJ06c4P333+fbb78lJycHgEaNGhEeHk6PHj1wcHAosv2+0yl0n7aBwja56NE8lOe0KkpEROQfUyglIiIiIreU9PR0Zs2aBUDbtm2JiIjgwQcfxGQqeEVTUno29mY7TMB/lu5nzqZYm0CqXe0Atp64SMd6gYx/silmrYwSEREpEQqlREREROSmZbFYWLp0KXv27OG1114DoF69eowZM4Z27dpx9913F9k+MTWLDuNXk56dS4iPKzHn02zqnRsEMe3ZZlgMMIEe1RMRESlBevteGdPb90RERESuXXZ2NnPmzGHs2LHs3bsXs9nM0aNHqVKlSrHa51kM3vhhN//behJLAV8Ne7k4sPGN+3B2MJfwyEVERG59evueiIiIiNxyUlNT+fzzz/noo4+Ii4sDwMPDg/79++Pq6nrV9gt3nCIxNRt3Z3vmbTlpUwur4sN7/2rAobOXqOrvqkBKRESklCmUEhEREZGbwpo1a3j00Ue5ePEiAEFBQQwdOpT+/fvj5eVVaDvDMHh/8T42H7/A3tMpBV6z+c2OBHg4AdAopPC+REREpOQolBIRERGRcisnJ8f6trzGjRuTm5tLrVq1GDVqFL169cLZ2bnI9mlZuew4mcSXG47nqzmYTbzUpjoNKnpaAykRERG5cRRKiYiIiEi5s23bNqKiojh58iTr1q3DZDLh7e3Nhg0bqFevHmZz4Y/WfbcxlhOJabSpFcDA77aRnJFjU1/86j2cSkrHz92JO6v6lvZUREREpBAKpURERESkXDAMgxUrVhAZGcnvv/9uPb9r1y6aNGkCQMOGDQttfyEtmxOJabzx424APl1zzKbeqUEFOjUIolGIlx7RExERKQcUSomIiIhImcrLy2PBggVERUWxdetWAMxmMz169CA8PNwaSBUkPTuX9Ow8PJ0d6DZ1PccT023qDzUK5nRyBrl5BpOevgMne21eLiIiUl4olBIRERGRMrV06VJ69OgBgIuLCy+99BLDhw+natWqhba5mJbN+dQsIhbsYvvJJOztTOTkGdb6vXUCeLhxRR67oxImE5hMptKehoiIiFwjhVIiIiIickNdvHiR/fv306pVKwC6dOlCy5YteeCBB3j11Vfx9/cvsr1hGPSZtYmdccnWc38NpLxdHYh6vIk2LxcRESnnFEqJiIiIyA0RFxfHhAkT+Oyzz3Bzc+P48eM4OztjNptZv359oauZsnLzuJiWg4+bA//bEsdve8/YBFK1At2pHuCGYUBk98Y42Nvh7qQvc0VERMo7/d9aRERERErVvn37GDt2LLNnzyYn5/Kb8KpVq0ZcXBw1a9YEin68bti8HSzdfabAmpujmc+ea041f7eSH7iIiIiUKoVSIiIiIlIq9u3bx2uvvcbixYut59q3b094eDidO3e+6j5Pqw4mMPvPWH7ff7bA+vRnm9G6pj8ezg4lOm4RERG5MRRKiYiIiEipyM3NZfHixZhMJh577DHCw8Np0aJFkW1WHUjg2z9P0KtlFfp/s5WsXIu15uFsz7SeYfi4ObA//hKdGgRpA3MREZGbmMkwDOPql0lpSUlJwcvLi+TkZDw9Pct6OCIiIiLXJTs7m++++47Y2Fjeeecd6/mPP/6Yzp07U6dOnSLbRx9N5NuNJ/h5V3yB9bceqkevllVwsjeX6LhFRESk5BU361AoVcYUSomIiMjN7NKlS3z++edMmDCBuLg4HBwciImJoVKlSsVqbxgG6dl5dPxoNfHJmTY1kwl6t6xKTp6FNx+qh6ujFvmLiIjcDIqbdej/7CIiIiJyzc6ePcukSZOYOnUqSUlJAAQHBzN06NCr/qBtf3wKs9bHUCvQg0/XHOV8arZN/cV7qjGgfQ0ycy1U8nYprSmIiIhIGVMoJSIiIiLXZOnSpXTr1o2srCwAateuzahRo+jVqxdOTk6FtjMMg8S0bD5YvI/oY4n56v7uTtxbJ4CB99bE182x1MYvIiIi5YNCKRERERG5qvT0dFxdXQFo0aIFZrOZFi1aEBERQdeuXbGzsyuwXWZOHn8cTKBNrQC+WBfDR8sP2dSbhHrTtpY/P2w7xfRnw2gU4lXqcxEREZHyQaGUiIiIiBTIMAxWrFhBZGQkmZmZrF27FgA/Pz92795NtWrVCn37nWEYnE7O5Ovo43y6+hhO9nY2b9IDePbuygy8tybBXi6MeKDojdBFRETk1qNQSkRERERs5ObmsmDBAqKioti2bRsAZrOZI0eOULNmTQCqV69eZB9zN5/k9R92W4//HkiN7taIp++qXMIjFxERkZuJQikRERERASAjI4Mvv/yScePGcezYMQBcXV156aWXGD58OFWqVCm0bcKlTBzNdny+9nK7KauO2tT7ta2Or5sjz95dhfOXsqji51p6ExEREZGbgkIpEREREQFg8eLFDBgwALj8iN6gQYN49dVX8fPzK7Ld8n1n6fv1lkLrwzrWZkjHWtZjdyd9CSoiIiIKpURERERuWydPnuTYsWO0a9cOgG7dunHvvffy2GOP8cILL+Dm5lZk+5TMHL798wQz18XYnHc025GdZ+GpO0Pp07oqtQI9Sm0OIiIicvNSKCUiIiJym9m3bx9RUVHMnj2bChUqcOzYMRwdHbG3t2flypVFtj12LpUXvtxMpwZB7ItPYe3h8zb1ni0q8/bD9Vl96Byta/prVZSIiIgUSl8liIiIiNwm1q9fT2RkJIsXL7aeq1WrFgkJCYSEhBTZNjMnjwnLD7HqYALHE9P5dM0xm3r1ADdm9bmTEB9XzHYmOjUIKpU5iIiIyK3jmkOp0aNHs2PHDs6ePYubmxv16tXjscceo3Xr1qUxPhERERH5hzZt2sSwYcPYsGEDACaTiW7duhEeHs5dd91VZNvj59NwdTLzw7ZT+YIof3cnXm5bnb2nk3nqrspU8Sv6cT8RERGRvzIZhmFcS4Nq1apRr149fH19uXTpErt27eLEiRN07NiR+fPn4+3tXUpDvTWlpKTg5eVFcnIynp6eZT0cERERuQVt3bqV5s2b4+joSO/evRk5ciS1a9cu9PrcPAsfrzhMdp6Fz9ccw/K3rxZNJtj4xn34ujpib7Yr5dGLiIjIzaa4Wcc1r5SKiYnJd27Tpk3079+fgQMHMnv27GvtUkRERERKyKVLl/j8889JTk7m/fffByAsLIxPP/2Uf/3rXwQHBxfZ3mIxWLIrnk9WHimw3qaWPz3uDCXQw7nExy4iIiK3l2teKVWYHTt20KZNGy5dulQS3d02tFJKRERESsLZs2eZNGkSU6dOJSkpCScnJ44fP05Q0NX3dvp1TzwVvV2wGPDKt1uJT860qc/s05z5m+N4uEkwDzeuWFpTEBERkVtEqa2U+qtZs2bh7u6Oo6MjP/30EwEBAf+kOxERERG5RkePHmXcuHHMmjWLrKwsAGrXrs2oUaPw8fEpsm1unoUdJ5Po/+22AuuNKnnxXMsqdKhbgQ51K5T42EVEROT29o9CqY0bN/K///2PpKQkunTpwqJFi0pqXCIiIiJyFXPnzqVnz55YLBYAWrRoQUREBF27dsXOrui9ns6nZvHgx2s5dymrwPpdVX2Z379liY9ZRERE5Ip/FEpNnz6dadOm8euvvzJy5Eg2b95Mw4YNS2psIiIiIvIXhmGQlJRkXQHVoUMHnJycuPfee4mIiKBNmzaYTKYC2+bmWcjKtTBs3g5yLQaezvY2gVTNQHeeCAvh+dbVOJuSiY+b4w2Zk4iIiNy+rnlPqbZt2xIVFcXdd99tc37Hjh106dKF06dPl+gAb3XaU0pERESuJjc3l++//56oqCi8vb1ZuXKltXbmzJmr7huVnp3LAxPWEHcxo8B6iI8Lq0fdi9mu4EBLRERE5FqU2p5SjRs35p577uGuu+6ie/fuNGrUCHd3d+bMmUNGRsFf6IiIiIjItcvIyGDWrFmMHz+eY8eOAeDq6kpcXBwhISEARQZShmHwx6Fz/Lr7jE0g1STEiz2nU6gV6E7U443xc3dSICUiIiI33DWHUpMnT2bAgAGMHTuWDz74wPq2PZPJxH//+98SH6CIiIjI7ebChQtMnTqVSZMmce7cOQD8/PwYNGgQr776Kn5+foW2tVgMes/aRGJqNhW9nfl9f4JN3cfVgW9faoHFAAezCVfHf7Sbg4iIiMh1u+bH9/4qLy+Po0ePkpSURJUqVahQQW9luVZ6fE9ERET+7uuvv6Z3794AVKlShREjRvDCCy/g5uZWZLvVh87x+76zfPPniQLrnz/XnLpBHoT6upb4mEVERESuKLXH9/7KbDZTu3btf9KFiIiIyG1v7969xMfH07FjRwCefvpp5s2bR8+ePXniiSdwcHAotO3WExeJu5hO01BvXvxyM7mW///zRid7O6b2bMbm4xfxcLbn/vr6AaKIiIiUH/9opZT8c1opJSIicnsyDIN169YRFRXFkiVLqFKlCocPHy4ygPqrDUfOszMumUkrDpORk5evPqFHE1rX9CfQw7mkhy4iIiJSpOJmHXY3cEwFmjp1KtWqVcPZ2ZmwsDDWrl1b5PWrV68mLCwMZ2dnqlevzvTp0/Nds2DBAurXr4+TkxP169fnxx9/vOb7/vDDD3Tq1Al/f39MJhM7duzI10dWVhaDBg3C398fNzc3HnnkEeLi4q7tAxAREZHbisViYeHChbRu3Zq2bduyZMkSTCYTzZs3JykpqVh9pGfn0vfrLUT+eiBfIBXq68JDjYL5V+OKCqRERESkXCvTUGrevHkMHTqUN998k+3bt9OmTRsefPBBYmNjC7w+JiaGLl260KZNG7Zv384bb7zB4MGDWbBggfWa6OhoevToQa9evdi5cye9evXiySefZOPGjdd037S0NFq3bs2YMWMKHf/QoUP58ccfmTt3LuvWrSM1NZWHH36YvLz8P60UERERWblyJQ0bNuTRRx8lOjoaR0dH+vbty4EDB/j+++8JCAgosN3BM5dYuOMU22Mv8tzMTdR/5zfSsv//1xvhnevwn8caMunpO1gb3oEpPZthby7znz2KiIiIFKlMH99r0aIFzZo1Y9q0adZz9erV49FHH2X06NH5ro+IiGDRokXs37/feq5///7s3LmT6OhoAHr06EFKSgq//PKL9ZrOnTvj4+PDnDlzrvm+x48fp1q1amzfvp2mTZtazycnJxMQEMA333xDjx49ADh9+jShoaEsXbqUTp06Fesz0ON7IiIit48///yTli1b4unpySuvvMKQIUMIDg4uso3FYtB27CriLmYUWO9QN5CpPZvh7GAujSGLiIiIXLNy//hednY2W7du5YEHHrA5/8ADD7Bhw4YC20RHR+e7vlOnTmzZsoWcnJwir7nS5/XctyBbt24lJyfHpp+KFSvSsGHDa+pHREREbk1nz57ljTfe4K233rKeu/vuu/nmm2+IjY1lzJgxBQZShmGQlpULwIdL9lH9jaU2gVQVP1cczXY4O9ixNvxeZva5U4GUiIiI3JT+0dv3/onz58+Tl5dHhQq2b4GpUKECZ86cKbDNmTNnCrw+NzeX8+fPExwcXOg1V/q8nvsWNhZHR0d8fHyuqZ+srCyysrKsxykpKcW+p4iIiJR/R44cYdy4cXz55ZdkZWXh4uLC0KFD8ff3B+DZZ58tsn3497v4accp2tUO4Pf9CTa1BxsG8e9HG+LiaCY1M5dAT+0ZJSIiIjevMgulrjCZTDbHhmHkO3e16/9+vjh9Xut9i+tq/YwePZr333//H99HREREypctW7YQGRnJggULrF+f3H333URERODr63vV9hN/P8TE3w9bj/8eSH37YgvuqeVvPXZ1LPMv40RERET+kTL7asbf3x+z2ZxvVVFCQkK+VUxXBAUFFXi9vb09fn5+RV5zpc/ruW9hY8nOzubixYs2q6USEhJo1apVoe1ef/11hg8fbj1OSUkhNDS02PcVERGR8mf69Om88sor1uMuXboQERFBmzZtivxhlWEYbDlxkTPJmTaBFECDip442dvx5kP1Sc7ItgmkRERERG4FZbanlKOjI2FhYSxfvtzm/PLlywsNdVq2bJnv+mXLltG8eXMcHByKvOZKn9dz34KEhYXh4OBg0098fDx79uwpsh8nJyc8PT1tfomIiMjNJTc3l7Nnz1qPH3nkEVxdXXn22WfZtWsXP//8M23bti00kPpuYyytx6yk5eiVPDE9mkFzttvUJ/Rows+D2/DDgNaEVfGhQ93i/+BMRERE5GZRpuu+hw8fTq9evWjevDktW7bks88+IzY2lv79+wOXVxWdOnWKr7/+Grj8pr3JkyczfPhw+vbtS3R0NF988YX1rXoAQ4YMoW3btkRGRtK1a1cWLlzI77//zrp164p9X4ALFy4QGxvL6dOnATh48CBweYVUUFAQXl5evPjii4wYMQI/Pz98fX0ZOXIkjRo1omPHjqX+2YmIiMiNl56ezqxZsxg/fjy1atXit99+Ay6/7OT06dN4eXkV2T4pPZvVh87xxo+7C6w/e3dlWlb3p0ujoBIfu4iIiEh5U6ahVI8ePUhMTOSDDz4gPj6ehg0bsnTpUqpUqQJcXnkUGxtrvb5atWosXbqUYcOGMWXKFCpWrMikSZPo3r279ZpWrVoxd+5c3nrrLd5++21q1KjBvHnzaNGiRbHvC7Bo0SKef/556/FTTz0FwLvvvst7770HwIQJE7C3t+fJJ58kIyOD++67jy+//BKzWW/AERERuZUkJiYydepUJk2axPnz5wG4dOkS58+ft25gXlggdeBMCrPWHefFNtV468c9bDp+waY+uENN7qzmy+/7zjKqc13cnbRXlIiIiNweTMaVnTilTKSkpODl5UVycrIe5RMRESlnYmNj+eijj5gxYwZpaWkAVK1alZEjR/L888/j6upaaNuUzBx+2R3PvM0n2RablK/eqoYfvVtV5f56FbCz++cvWxEREREpL4qbdehHcSIiIiKFWLZsGR9//DEATZo0ISIigieeeAJ7+6t/CfWfJfuZt+WkzTmTCeoHe3IpM5c3utSjYaWiH/cTERERuZUplBIRERHh8pvw1q5dS3p6Op07dwagV69eLFu2jBdffJEHHnig0I3LM3PyWLIrHg9neyJ/PQDAsXNp1nr1ADf+3bUhQV7OVA9wL/3JiIiIiNwEFEqJiIjIbc1isbBo0SIiIyP5888/qVmzJgcOHMBsNuPk5MT8+fOv2seUVUf4ZOWRAmv+7o68/0gDWtX0L+mhi4iIiNzUFEqJiIjIbSkrK4vZs2czduxYDhy4vLrJycmJ++67j7S0tCL3PziSkEo1fzd2xiUR+csBNsbYbl7+UKNglu6JZ+zjTXg8LKRU5yEiIiJys1IoJSIiIredn376iYEDB3L69Gng8pvzBgwYwODBgwkKCiqy7ZJdp3n1u+14OtuTmWshO9dirfm7O/Jh14Y82CiYzJw8nB30Rl4RERGRwiiUEhERkduOn58fp0+fpmLFigwbNox+/fpd9S24RxIuMXD2dg6evQRASmauTb1bs0p89GRT67ECKREREZGiKZQSERGRW9rhw4cZP348fn5+/Oc//wHgnnvu4YcffqBLly44OTkV2jYrN48v1x+ngqczP24/ZQ2kALxdHXiwYRD929VgU8wFOtarUOpzEREREbmVmAzDMMp6ELezlJQUvLy8SE5OvupPaEVERKT4tmzZQmRkJAsWLMAwDNzd3YmLi8PLy+uqbeOTM3jxyy3EnE8jIycvX71VDT++63t3aQxbRERE5KZX3KxDK6VERETklmEYBsuWLSMyMpJVq1ZZzz/00ENEREQU6wdA++NTmLspln3xKflqdYM8eLJ5qFZFiYiIiJQAhVIiIiJyyxg7diwREREA2Nvb88wzzzBq1CgaNmxYaJuM7Dze/HE3YVV9OJucyaSVR2zq1QPc+PGV1vyyJ5576wZSwdO5VOcgIiIicrvQ43tlTI/viYiIXL/09HQuXrxIpUqVADh58iSNGzemT58+DBs2jMqVKxfZfuuJi6w+dI5JKw4XWH/n4fq0rxNA9QD3Eh+7iIiIyK1Kj++JiIjILSsxMZEpU6bwySefcPfdd7N48WIAQkNDOX36NC4uLoW2vZCWjZO9HfviU3hierRNzcFs4vGwEDYeu8Cjd1TihXuqleo8RERERG5nCqVERETkphEbG8tHH33E559/Tnp6OgD79u3j0qVLeHh4ABQaSO04mURKRg4DZm8DIMDD9q177/6rPvfU9KdWBY9SnIGIiIiIXKFQSkRERMq9ffv2MXr0aObMmUNe3uW34TVt2pSIiAgef/xx7O2L/pLm6LlUuk/bQJ7l/+9akJqVa/39v5pU5PnWWhUlIiIiciMplBIREZFyb82aNXz77bcA3HfffURERNCxY0dMJlOB159OysDebMLOZGLi74f49s9Ym/qgDjXZGHOBltX9GNShJvZmu1Kfg4iIiIjYUiglIiIi5YrFYmHRokU4ODjw0EMPAdC7d282b97MK6+8QvPmzYtsn5SeTaeJa7iUmVtgvW3tAIbfX7vQQEtEREREbgy9fa+M6e17IiIil2VlZTF79mzGjh3LgQMHqFOnDvv27cPOrnirmAzD4P3F+5i98QQ5efm/vPF0tufPN+7D1VE/kxMREREpTXr7noiIiNwUUlJS+PTTT5k4cSKnT58GwMvLi27dupGVlVXkm/QMwyDy14MkXMqkaag3X244blO/u7ov03qGsSMuiQB3JwVSIiIiIuWIvjITERGRMvPVV18xZMgQkpOTAahUqRLDhg2jb9++V11BPPqX/fy25wzHEy+/he+Hbads6utf60Al78uB1r11Akth9CIiIiLyTyiUEhERkRvKMAzrfk4hISEkJydTt25dwsPD6dmzJ46OjkW2T8nM4dCZS3y6+li+mq+bIy/eU41K3i7WQEpEREREyieFUiIiInJDbN68mcjISOrUqcN//vMfADp06MDy5cvp0KFDkXtHfbXhOPtOp+DqZGb2xliycy029Z8Gtsbf3REXBzN+7k6lOg8RERERKRna6LyMaaNzERG5lRmGwbJly4iMjGTVqlXA5f2i4uPji9wr6oqUzBxiE9N5+JN1BdY71qvAvXUD6NmiSomOW0RERESunzY6FxERkTKTm5vL/PnziYqKYufOnQDY29vzzDPPMGrUqCIDqdjEdA4nXOKeWv48Pm0Dh86m2tRfvbcmx86nkpKRy+Rn7sDZwVyqcxERERGR0qFQSkRERErce++9Z31Ez83NjX79+jF06FAqV65caJvkjBzOpmTS/9utHDuXlq/ernYADzUO5omwEOueVCIiIiJy81IoJSIiIv9YYmIiqampVKly+TG6F198kZkzZzJgwAAGDBiAr6/vVft4+Zst/HnsQoE1D2d7xnRvRLCXNi8XERERuVUolBIREZHrduLECT766CNmzJhB586dWbBgAQDVqlUjNjYWe/uCv9TIzMkjPjkTs8nEtxtPEHcx3SaQqh/sSfs6Afi7O/FMi8rkWQzcnPRli4iIiMitRF/diYiIyDXbvXs3UVFRzJkzh7y8POByQJWVlYWT0+W33xUWSAGEf7+LRTtPF1hztLfjox5NqBukF4CIiIiI3MoUSomIiEix/fnnn3z44YcsXbrUeu6+++4jIiKCjh07XnWvp+ijiXzz53GW7j5jc97FwUxGTh6fP9ec5lV88HFzLJXxi4iIiEj5oVBKREREim3Tpk0sXboUOzs7unfvTnh4OM2bNy+yzeyNJ/hkxRF6tazCtD+OkpqVa605mu348vk7qRnozsGzl2hTK6C0pyAiIiIi5YRCKRERESlQVlYW3377LRUqVODhhx8GLm9gfuzYMV599VVq1qxZZPttsReZuS6GJbviARj720Gb+sgHatOndTXc/2+vqEBP51KYhYiIiIiUVybDMIyyHsTtLCUlBS8vL5KTk/H01N4ZIiJS9pKTk/n000+ZOHEi8fHxNGzYkF27dl310TwAi8Xg2Pk0Atyd6DZtPUfPpdnU6wd78uzdVThwJoXwznWtgZSIiIiI3DqKm3XoK0EREREBID4+no8//php06aRkpICQKVKlejTpw+5ubk4ODgU2vbw2UtMWXWE1Kw8ft9/Nl/96bsq894j9XGyN5fa+EVERETk5qJQSkRERPj4448JDw8nOzsbgLp16xIeHk7Pnj1xdCx60/H07Fw+/Hk/aw6dy1er5O1Ci2q+DOtYS4GUiIiIiNhQKCUiInKbysvLw2y+HBTVqVOH7OxsWrVqRUREBA8//DB2dnYFtruUmcO8zSd5qHEwC7bGMW7ZIZt6y+p+PHVXKN9vjePdf9WnZqBHqc9FRERERG4+2lOqjGlPKRERuZEMw+DXX38lKiqKe+65hw8//NB6fvPmzdx1111Ftj+TnMkX647x+dqYAus9mofyaoeahPq6lvjYRUREROTmUNysQ6FUGVMoJSIiN0Jubi7z5s0jKiqKXbt2ARAYGEhcXFyRe0X91cIdpxgyd0eh9Q+7NqBXy6olMFoRERERuZlpo3MREREhLS2NmTNnMn78eE6cOAGAm5sb/fr1Y9iwYUUGUmdTMrmUmcOo73dRzd+NH7adsqn3aVWVe+sG0raWP2dTsqjg6VSqcxERERGRW4tCKRERkVvYG2+8waRJkwAICAhg8ODBDBgwAF9f3yLbrTl0judmbrIeb49Nsqm/3K46rz9Yz3oc5OVccoMWERERkduCQikREZFbyPHjx7FYLFSvXh2AgQMHsnTpUoYPH06fPn1wcXEpsJ3FYmBnZyI5I4fvNsbydfRxm3p1fzeOnU9jQPsaPNK0IjUC3Et7KiIiIiJyi1MoJSIicgvYtWsXkZGRzJs3j+7duzNv3jwAateuzaFDhzCZTIW2PXoulcemrKdt7QD2x6dw9FyaTf3xsBDGPt6YhEtZBHo4FdmXiIiIiEhxKZQSERG5SRmGwerVq4mMjOTXX3+1nk9OTiY3Nxd7+8v/my8sRMrOtTB55WGW7TtLSmYuS3bF29RrBLgxo/edVPJ2wWQyUcFTj+iJiIiISMlRKCUiInIT+u2333jnnXfYtOnyvk92dnY8/vjjhIeHExYWVmg7wzBYtPM0vm6OHDqbyqT/1959h0dV5u8fv2fSOySBhFATmqFDKNIElaIgRXFFUQRRhAWBEEpsrK5+BRMEFZAqq6Ii4LLWRQUEKRoF6b0GQguQEDLpZeb8/uDHrCNJiEoSQt6v68qlcz7PZ85zcuWAuX3OM+uOOtSbVvfTtAeaavnWU+rbIkShgV4leh0AAACouAilAAAoh/bs2aMtW7bIzc1NTzzxhCZMmKB69eoVOt5qMzR73RGdSMrQ5zvPSpJ+v4Bqy/N3K9DbTWazSU2q+5Xk9AEAAABCKQAAbnapqamaP3++GjVqpD59+kiSnn76aaWlpWnUqFEKCgoqst8wDK3Zn6i31h753fEr/+zasIr6NAtRVR7PAwAAQCkyGcbV/yRFWbBYLPLz81Nqaqp8fX3LejoAgJvI2bNn9dZbb2n+/PlKS0tTy5YttW3btmJtNL5y22kFeLsqNStPMd8c1NnUbIf6f8d20qYjSWofFqDmNSuV0BUAAACgIipu1sFKKQAAbjKHDh3S9OnT9eGHHyo3N1eS1KhRI40bN042m01OTk6F9tpshnafSdWET3cVWG8Y5KPHbq+lxiF+ahzCI3oAAAAoO4RSAADcRF566SW9+uqrurqQuWPHjoqOjlbv3r1lNpsL7MnJt8rVyazUrDzdN3uzTqdkOdSbVvfTgXMWtapVWStGti/xawAAAACKg1AKAIAyZBiG8vLy5OrqKkmKiIiQYRjq27evJk+erI4dOxbaa7UZSs/O1/1zf1RmrlUuziaHQKpOgKceb19HT3SsI0tWvtxdCw61AAAAgLJAKAUAQBnIy8vTihUrFBsbqwcffFBTpkyRJN133306cOCAbrvttiL7c/KtuvftTTp+MaPAelUfN60e30WuzleCKD9Plxt7AQAAAMBfRCgFAEApysjI0OLFizVz5kydPHlSkpSWlqYXXnhBZrNZZrO50EDKMAxl59m0/tAFbT6a5BBIuTiZVCfAS41CfDWyS115uznbAykAAADgZkQoBQBAKUhKStKcOXM0Z84cJScnS5KqVq2qcePG6e9//3uh+0VJV8Koxxb/ov1nLcrJtykz1+pQ93Zz1uboO1XJ07VErwEAAAC4kQilAAAoBc8995zeffddSVLdunU1ceJEDRkyRB4eHkX2/Xw8WWv3n9ePR5MLrM9/rJUaBPkQSAEAAKDcIZQCAKAE7Nq1S97e3qpbt64kKSoqSjt37tSkSZM0YMAAOTk5Fdr76a+ntO+sRQPb1NTj/9qi3HybvebqbNa/R7bXzlOXZTaZdE+TaiV+LQAAAEBJMBlXP3MaZcJiscjPz0+pqany9fUt6+kAAP4CwzC0YcMGxcTE6Ntvv9Vjjz2mDz/8sNj9v564pJ+PJ+uttUeUb7v2r+fYB5upc/1AVfMrenUVAAAAUJaKm3WwUgoAgL/IarXqiy++UExMjLZs2SJJ9j2iDMOQyWS67nvk5Fs1fMmvSsnMczju6mRWi1qV5O/pqv4tqrN5OQAAAG4ZhFIAAPwFy5Yt00svvaTDhw9Lktzd3TVs2DBNmDBBYWFhhfYdTLRo+8nLupyVqxVbT+lEcqZDPfqe29SqViWFVPJQTX/PEr0GAAAAoCwQSgEA8BfEx8fr8OHDqly5skaPHq0xY8aoatWqRfYYhqFRH23X8aSMAusd6gZocPva8nbjr2kAAADcuvivXQAAiuns2bN6++23dccdd6h3796SpL///e/y8PDQU089JW9v7wL7snKtSriUqYbBPpr+3UG9s/6YQ/22YB+5uzgpPilDq8Z1VvVK7BkFAACAWx8bnZcxNjoHgJvfoUOHNH36dH344YfKzc1Vu3btFBcXV6y9oiRp9NLt+u/ucwr2dVeiJduh1rFegN4a2FKB3q7KsxrsGQUAAIByj43OAQD4i3755RfFxMTo888/19X/h9OpUydNnjy5WP0LNx7TtG8O6ur//vl9IPXe0Da687b/Pern6ly8kAsAAAC4FRBKAQBQgDFjxmjOnDn213379lV0dLQ6dOhQZN/plEx9ueusqlfy0NRVBx1q7cMC1Kl+oPq1CNHRC+nq2rDovacAAACAWxmhFAAAkvLy8pSfny8Pjyv7OXXt2lULFizQY489pkmTJik8PLzI/k9/PaWYbw8pKT2nwPrU+5tqULta9tc1KvOJegAAAKjYCKUAABVaRkaGFi9erBkzZujvf/+7nn32WUlS//79FR8fr+rVqxfZn5mbr1/iL+m5/+xRvu3abRqHdQxVsxp+6ts8pETmDwAAAJRXbHRextjoHADKRlJSkubMmaM5c+YoOTlZktS8eXPt2LHjuhuY7zx1WW98d0jD7wjTm2sOa+epyw71MXfV0wOtamjbyRTd37K6nMzsFQUAAICKg43OAQAowIkTJzRjxgwtXrxYWVlZkqS6detq4sSJGjJkSJGBVGZuvv67+5yWxJ3UnjOp2nw0yaHepk5lDelQR/c0Dpazk1mhgV4lei0AAABAeUYoBQCoUKZMmaKPPvpIkhQREaHo6Gg98MADcnJyKrTngiVbAd5uivnmoD6IO+lQq+LjpnubBOvn48n6x32N1bSGX4nOHwAAALhVEEoBAG5ZhmHohx9+UK1atVS3bl1J0qRJk3T+/HlFR0frrrvuKnRlVG6+TZ/vPKPLmbmauuqgfNydlZadb68H+7rr9QFN1aaOv7zc+OsUAAAA+KP4r2gAwC3HarXq888/V0xMjLZu3aphw4Zp8eLFkqRmzZpp9erV132PRZuOa/p3h+yvfxtIBXi56pV+jdW1YdUbP3kAAACggiCUAgDcMrKzs/Xhhx9q+vTpOnLkiCTJ3d1dfn5+Mgyj0FVRhmFo64kUNQ7x1Z4zqZr1/RH9dCzZYcyE7g303f5EPX9vuDrUCyzxawEAAABudYRSAIBbwrx58/TKK68oMTFRklS5cmWNHj1aY8aMUdWqRa9o+nr3OY35ZIdqB3jq7OUs5Vn/98G0Pu7O+mffxnqgVQ2Nubt+iV4DAAAAUJEQSgEAbgkXL15UYmKiatasqaioKD311FPy9vYusufUpUyN+ni79pxJlSSdTM6UdGW/qERLtvq1CNHbD7cs8bkDAAAAFRGhFACg3Dl48KCmT5+uAQMGqFevXpKk0aNHq06dOnrkkUfk4uJSaG9KRq6e/2yP6lbx1v5zFnsgJV1ZFfVUpzD9vWtdrT90QbeHBpT4tQAAAAAVlckwDOP6w1BSLBaL/Pz8lJqaKl9f37KeDgDc1H7++WfFxMToiy++kGEY6tSpkzZt2lSs3otpORr2/laHEOq3WtaqpM9GdbyR0wUAAAAqpOJmHayUAgDc1AzD0KpVqxQbG6uNGzfaj/fr10/R0dHFeo+DiRat3Ha6wECqU71AdQuvqjsaVLlhcwYAAABwfYRSAICb2qOPPqpPPvlEkuTi4qLHHntMkyZNUnh4eKE9WblWjV++UzUqeyg+KUPfH7zgUK/p76H1E7rq2MUM1fL3lIerU4leAwAAAIBrEUoBAG4qGRkZMplM8vT0lCT16dNHX331lUaOHKnIyEhVr169yP69Z1K1/uAFfbsvscD6C73C1blBoJydzGoY7HPD5w8AAACgeNhTqoyxpxQAXHHx4kXNnj1b77zzjl544QVFRUVJkvLz85WWlqbKlSsX2rs9IUWerk7KyrVqwLyfZPvN32yerk6KGdBMS+JO6O7wII3sUrekLwUAAACo0IqbdRBKlTFCKQAVXXx8vGbMmKF//etfysrKkiR17dpV69evv27v3jOpOp2SqdFLd8gwDIcwSrqyKqprwyqqH8SKKAAAAKC0sNE5AOCmtnPnTsXGxmrFihWyWq2SpNatWys6Olr333//dftPp2Tqgbk/KddqK7DevVGQht8RdkPnDAAAAODGIZQCAJSJadOmacWKFZKkHj16KDo6WnfeeadMJlOB4xOSM2XI0JmULK3cfkb/2XFav13r+1SnUDmZTYqoXVl3hwfJyVzw+wAAAAC4ORBKAQBKnNVq1WeffaZWrVopLOzK6qXJkyfLbDZr8uTJatmyZZH9adl5um/2Jlmy8wusR9SurOd7hctMEAUAAACUG+wpVcbYUwrArSw7O1sffPCB3njjDR09elQjRozQ/Pnzi91vGIZivj2kD+NOKCPXaj/u5mxWTr5NlTxdtHHynfJxcy50hRUAAACA0sWeUgCAMnP58mXNmzdPb7/9ts6fPy9Jqly5smrVqnXdXsMw9NQHv+pEcoba1w3QRz8nONQbh/hq+Yj2OnUpUx4uTvJ1dymRawAAAABQsgilAAA31P/93/8pNjZWaWlpkqSaNWtqwoQJevLJJ+Xt7V1k71trD+vf207rdMqVT+E7djHDob5uQheFVbnyHuHVWF0KAAAAlGeEUgCAGyo7O1tpaWlq0qSJJk+erIcfflguLoWvZsrJt+rUpUxl59n01tojDjUvVyc1q1FJfVuEyNXJbA+kAAAAAJR/hFIAgD8tLi5OMTExGjVqlHr06CFJGjt2rNq3b69evXoVuc/TkrgT+vl4svaftehEcuY19WVP367bwwJKbO4AAAAAyhahFADgD7HZbFq1apViYmK0efNmSVJaWpo9lKpatap69+5daH+e1aYzKVn6xxf7Cqz3aBSkdmEBahfqf+MnDwAAAOCmQSgFACiWvLw8ffLJJ4qNjdW+fVcCJRcXFw0ePFiTJk0qsnfP6VT9Ep+sv7Wuqcfe/UV7zqQ61F/p11gXLDk6l5qtaQ80lauzucSuAwAAAMDNgVAKAFAsffv21bfffitJ8vHx0YgRIxQZGanq1asX2pORk6+zl7M0eul2JVzK1P/994BDvWWtSurdtJoebVdbTubCH/UDAAAAcOshlAIAFOjixYvy8vKSp6enJOmRRx7Rjh07FBkZqZEjR6pSpUrXfY8xn+zQuoMXCqx5uzlr1sMtVdPf80ZOGwAAAEA5YTIMwyjrSVRkFotFfn5+Sk1Nla8vH28OoOwdP35cM2fO1L/+9S/FxMRozJgxkq48vme1WuXu7l5gX3aeVQmXMnXgnEUf/HRC3u4u2nj4or1ev6q3HoyoobvDgxQa6KU8q03uLk6lck0AAAAASk9xsw5WSgEAJEk7duxQbGysVqxYIZvNJknasGGDPZRycXGRi4tLof1TPt+rT7edLrDmbDZpxkPN1axGJfsxJzOBFAAAAFCREUoBQAW3bt06vf7661qzZo39WM+ePRUdHa2uXbsW2peTb5Wbs5N2nrqs5VtPXRNIhVXx0qlLmfpgWFs1CPJRoLdbSV0CAAAAgHKIUAoAKrjZs2drzZo1cnJy0sCBAzVp0iS1aNGiyJ7Fm+M1bdUB9WwSrB8OXlBGrtVeM5uk959oq471ApWSmUsYBQAAAKBAhFIAUIFkZ2frgw8+UM+ePVWnTh1J0rPPPqsaNWooKipKoaGhRfbvP2vRgo3H9MXOs5Kk/+4+51Afe1c9DesUqkqerpJEIAUAAACgUIRSAFABXL58WfPmzdPbb7+t8+fP65lnntHs2bMlSe3atVO7du0K7U3PydfnO86ofd0ATfr3Lu07a3GoP9ymph5oVUNrD5zX8DvC5ONe+L5TAAAAAHAVoRQA3MLOnDmjN998UwsWLFB6erokqVatWmrSpMl1e+OTMjRzzWEdOZ+mg4lp19Tvb1ldMQOaydXZLElqG+p/YycPAAAA4JZGKAUAt6jIyEjNnTtXeXl5kqSmTZtq8uTJGjhwYJGfoidJ2XlWTVt1QKv3n7+mVq+qtxqH+GpCjwb2QAoAAAAA/ihCKQC4Rbm5uSkvL0933HGHoqOjde+998pkMhU4NjUrT/M3HFP/FtX1n+2ntWDjcYd62zr+er53uH44dEFDO9Sx7xkFAAAAAH9Wmf8v7rlz5yo0NFTu7u6KiIjQpk2bihy/YcMGRUREyN3dXWFhYZo/f/41Y1auXKlGjRrJzc1NjRo10mefffaHz2sYhl5++WWFhITIw8NDXbt21b59+xzGdO3aVSaTyeHr4Ycf/hPfBQD482w2m77++mt17txZ69atsx+PiopSXFycNmzYoF69ehUaSCWl52ju+qOa98Mx9Xxr4zWB1IBWNTT1gaZqUbOSIrs1IJACAAAAcEOUaSi1fPlyRUZG6oUXXtCOHTvUuXNn3XvvvUpISChwfHx8vHr16qXOnTtrx44dev755zV27FitXLnSPiYuLk4DBw7U4MGDtWvXLg0ePFgPPfSQfvnllz903tjYWM2cOVNz5szR1q1bFRwcrO7duystzXFfleHDh+vcuXP2rwULFtzg7xIAFCwvL09LlixRs2bN1KdPH23evFkzZsyw14OCgnT77bcX2JuTb5UkfbcvUa3/b61DEOXqZFbf5iGSpH/2bawZDzVXvareJXglAAAAACoik2EYRlmdvF27dmrVqpXmzZtnPxYeHq7+/ftr2rRp14yPjo7Wl19+qQMHDtiPjRw5Urt27VJcXJwkaeDAgbJYLPrmm2/sY+655x5VrlxZn3zySbHOaxiGQkJCFBkZqejoaElSTk6OgoKCFBMToxEjRki6slKqRYsWeuutt/7098BiscjPz0+pqany9fX90+8DoOJIT0/XokWL9Oabb+rUqVOSJB8fH40cOVKRkZEKCQkptPe8JVvxSRl6/F9b1DjEVzsSLjvUh7Svrftb1VCLmpWUkZMvLzee8gYAAADwxxQ36yizlVK5ubnatm2bevTo4XC8R48e+umnnwrsiYuLu2Z8z5499euvv9o38i1szNX3LM554+PjlZiY6DDGzc1NXbp0uWZuH3/8sQIDA9W4cWNNnDjxmpVUv5eTkyOLxeLwBQB/RI8ePRQVFaVTp04pKChI06ZNU0JCgmJjY4sMpLbEX9Lt077Xwwt/Vm6+rcBA6p/9mqhFzUqSRCAFAAAAoESV2W8cSUlJslqtCgoKcjgeFBSkxMTEAnsSExMLHJ+fn6+kpCRVq1at0DFX37M45736z4LGnDx50v760UcfVWhoqIKDg7V3714999xz2rVrl9asWVPodU+bNk3//Oc/C60DwO8dP35cISEhcnd3lyQ9+eSTSkpK0qRJkzR48GD78d/LzrPK1cmspPQcrdx+Rp9sSdBv18Z2qheo+KQMjbqzrtrU8VftAM/SuBwAAAAAkHQTfPre7zfeNQyj0M14Cxv/++PFec8bMWb48OH2f2/SpInq16+v1q1ba/v27WrVqlWB83/uuecUFRVlf22xWFSzZs0CxwKo2Hbs2KGYmBh9+umnmjt3rv3R4SFDhmjo0KFycnIqtPf4xXT1nrVZTmaTcq025ebbHOq9m1XTO4MK/nMKAAAAAEpDmYVSgYGBcnJyumZV1IULF65ZoXRVcHBwgeOdnZ0VEBBQ5Jir71mc8wYHB0u6smKqWrVqxZqbJLVq1UouLi46cuRIoaGUm5ub3NzcCn0PABWbYRhat26dYmJiHFZd7tixw/7vzs6F/9Gdb7VpwcbjWrXnnLLyrNfU61bx0oLBrVWjsseNnTgAAAAA/EFltqeUq6urIiIirnnUbc2aNerQoUOBPe3bt79m/OrVq9W6dWu5uLgUOebqexbnvFcfyfvtmNzcXG3YsKHQuUnSvn37lJeX5xBkAUBxffrpp2rTpo26deumNWvWyMnJSYMGDdLOnTs1f/78QvvyrDbFfHtQS+JOaEncSU3/7pD2nf3ffnUtalbSry920/huDTTzoRaqV9Vb7i6Fr7ICAAAAgNJQpo/vRUVFafDgwWrdurXat2+vhQsXKiEhQSNHjpR05VG3M2fOaMmSJZKufNLenDlzFBUVpeHDhysuLk6LFy+2f6qeJI0bN0533HGHYmJi1K9fP33xxRdau3atNm/eXOzzmkwmRUZGaurUqapfv77q16+vqVOnytPTU4MGDZIkHTt2TB9//LF69eqlwMBA7d+/XxMmTFDLli3VsWPH0voWAriFvP/++9q2bZs8PDz05JNPKioqSqGhoYWONwxD8zYc0/aTl7X2wPkCx/z47F0K9nWXk9mkcd3ql9TUAQAAAOAPK9NQauDAgUpOTtYrr7yic+fOqUmTJlq1apVq164tSTp37pwSEhLs40NDQ7Vq1SqNHz9e77zzjkJCQjRr1iwNGDDAPqZDhw5atmyZXnzxRU2ZMkV169bV8uXL1a5du2KfV5ImT56srKwsjRo1SikpKWrXrp1Wr14tHx8fSVdWXH3//fd6++23lZ6erpo1a6p379566aWXitznBQAkKSUlRfPmzdPjjz+uGjVqSJJeeOEFtW7dWs8884yqVKly3ff44fBFxX576Jrjgd6uiqhdWXfdVlXVK/GYHgAAAICbk8kwfvtZTChtFotFfn5+Sk1Nla+vb1lPB0AJO336tN58800tXLhQ6enpioqK0owZM4rVu2LrKbk4m7QlPkXf7UvUpYxch/o34zrLkpWn+kE+8vdyLYnpAwAAAMB1FTfrKPNP3wOAimD//v2aPn26Pv74Y+Xl5UmSmjZtqttvv/26vYZh6MC5NE1eubvAelgVLw1sXVPh1Qi2AQAAAJQfhFIAUMIeffRRLV261P66S5cuio6O1j333COTyVRgz+mUTBmG5Ofpovvf+VHHLmY41Ae0qqFNRy7qtmq++uCJNoW+DwAAAADcrAilAOAGs9lsMpv/9+GmVatWlclkUv/+/RUdHe2wx921vYZSs/LUe9ZmZeVa5eHqpNSsPHs9wMtVj7evo3Hd6stqM+RkJowCAAAAUD6xp1QZY08p4NaRm5urpUuXavr06VqwYIE6deokSUpMTFRqaqoaNmxYZH++1aY+c37UgXOWAuv+Xq76MfouebjyYQoAAAAAbl7sKQUApSQtLU2LFi3Sm2++qdOnT0uSZs2aZQ+lgoODFRwcfE2fYRiy2gztOn1ZW0+k6Lwl2yGQ8vNwUef6gbqnSbDCq/nK1clMIAUAAADglkEoBQB/0oULFzRr1iy98847unz5sqQrAVRkZKRGjhxZZK9hGHr8X1u06UhSgXV3F7PWRN2hqj7uN3raAAAAAHBTIJQCgD/BMAzdfffd2rt3rySpfv36mjRpkgYPHix396KDpAPnLFq973yhgdSCwREKDfQikAIAAABwSyOUAoBi2r59uxo3biw3NzeZTCaNGjVK77//vqKjo9WvXz85ORX+aN2CDcf007FkDe1YR3//aJuy82z2mrPZpC+e6ajk9Fzl5tvUrVFQaVwOAAAAAJQpNjovY2x0DtzcDMPQ999/r5iYGK1du1aLFy/WsGHDJElWq1Vms1kmU+GfgLf79GVtOpKk6d8dKrD+Qq9wdW1YRfWDfEpk/gAAAABQ2tjoHAD+AqvVqpUrVyomJkbbt2+XJDk5Oen48eP2MYWtjMq32mTJzpePu7OGL/lV5y05DvVmNfx0W7CP0nPyNbh9bbm7sHk5AAAAgIqHUAoAfsMwDC1YsEDTp0+3B1AeHh566qmnFBUVpTp16hTaeygxTXHHkrT5aLLWHjh/Tf2ZO+upT/MQNQjyLnJ1FQAAAABUBIRSAPAbJpNJ//nPf3T8+HEFBARozJgxGj16tAIDA6/bG7l8pw6csxRYi6hdWU91DlUlT9cbPWUAAAAAKJcIpQBUaKdOndJbb72liRMnqlq1apKkl156SX369NGwYcPk5eVVYF9qZp52nEpR5/pVNHvdEc1Zd1T5tv9t0dctPEhBvm46fjFDcwa1VIC3W6lcDwAAAACUF4RSACqk/fv3KzY2Vh9//LHy8/Pl4uKi119/XZLUsWNHdezYscj+5z7brVV7EgusRdSurFf7N1Y1P48bPm8AAAAAuFUQSgGoUDZv3qzY2Fh99dVX9mNdu3ZV9+7di9X/0c8n9crX+5WbbyuwvmBwhHo2Dr4hcwUAAACAWxmhFIAKwTAM9ejRQ2vXrpV0Ze+o+++/X9HR0Wrbtm2Rvb+euKQ564+qR6Ngvfj5Xodap3qBGtElTM1qVNL+sxa1rxtQYtcAAAAAALcSQikAt6z8/Hw5O1/5Y85kMqlhw4bauHGjhgwZogkTJqhhw4ZF9n+566xe+mKvUjLzJEk/HLroUJ9yXyM92SnU/ppACgAAAACKz2QYhnH9YSgpFotFfn5+Sk1Nla+vb1lPB7glpKWladGiRXrzzTf173//W+3atZMkJSYmyjAM+4bmhcnIydeWE5c0ZukOpefk248HeLkqOSNXE7o3UEglD/VrESJnJ3OJXgsAAAAAlDfFzTpYKQXglnH+/HnNmjVLc+fO1eXLlyVJCxYssIdSwcGF7/W07WSKpny+V/c2CdbK7ad1IjnToT6iS5jG3lVfJ5Iz1DjEr8SuAQAAAAAqCkIpAOXesWPH9MYbb+i9995TTk6OJKlhw4aaNGmSHnvssSJ7c/Kt+mZPohZtOq795yzaf87iUG9Zq5Ieb19bvZpWk5uzE4EUAAAAANwghFIAyjWbzaYePXro+PHjkqR27dopOjpa/fr1k9lc8KN1Npuh7w9eUPMaflq06bgWbYp3qDev4acnOobq419Oasp9jdSsRqWSvgwAAAAAqHAIpQCUK4ZhaP369ercubNcXFxkNps1btw4ffvtt4qOjtYdd9whk8lUYG++1abPd57VgXMWLd4cf03d191Zswe1Uts6/vJwdVL/ltVL+nIAAAAAoMJio/MyxkbnQPHk5+dr5cqVio2N1fbt2/Xhhx/aH80zDKPQIOq3Fm+O16tf7y+wFujtppf6NFKf5iE3dN4AAAAAUNGw0TmAW0JWVpbee+89zZgxw/6Inqenpy5cuGAfU1AgZRiGfjh8UU2r+2nt/vP6fOcZ/Xz8ksOYN/7WXHtOX9bg9nVUr6p3yV4IAAAAAMABoRSAm5LNZtPUqVM1a9YsXbx4UZIUEBCgsWPHavTo0QoICCiy/5u9iRr18fYCa+4uZr3Qu5EejKihByNq3PC5AwAAAACuj1AKwE3JbDbr+++/18WLF1WnTh1NmDBBw4YNk6enZ4Hjrz7Cdy41S2OW7tCvJ1Mc6iaTZBhSvxYhevvhlqVxCQAAAACAIhBKAbgp7N+/XzNmzNDrr7+uKlWqSJJeffVVJSQk6KGHHpKzc+F/XJ1OydT9c39S8xp+Ons5W/vPWew1N2ezxndvoBF3hGnHqctqGORT4tcCAAAAALg+QikAZWrz5s2KjY3VV199JUmqXr26XnnlFUlSp06diuy9nJmrJz/4Vdv+/6qotQcuONTDq/nqm3Gd7a9b1ap8I6cOAAAAAPgLCKUAlDqbzaavv/5aMTEx+umnnyRd2az8gQceUN++fa/bn5mbryPn0/XdvkR7IHWVj5uzHmpTU3UCvdQ+zL9E5g8AAAAA+OsIpQCUKqvVqrZt22r79iubkLu6umrIkCGaOHGiGjRoUGhfdp5VzyzdLpsh7TmTqotpOQ71IF83/Rh9l5ydzCU6fwAAAADAjUEoBaDEZWVlycPDQ5Lk5OSk1q1b6+jRoxo1apTGjh2ratWqFdl/5Hya1h28cM3jeVc9e+9t6lA3gEAKAAAAAMoRk2EYRllPoiKzWCzy8/NTamqqfH19y3o6wA11/vx5zZo1S3PnztXatWsVEREhSbpw4YLc3d2L/JlfsfWU3FzMqh3gpQfn/aR82//+qPJ0ddKyp2/Xez+eUNtQfz3StlaJXwsAAAAAoHiKm3WwUgrADXf06FG98cYbev/995WTc+UxuyVLlthDqapVqxbaeygxTYfOp2nyyt0F1sfcVU933lZVzWpU0psDW9zwuQMAAAAASgehFIAbZtu2bYqJidHKlStls9kkSe3atVN0dLT69et33f4Llmz1e2ezsvNs19Q8XZ3ULtRfE3o0vOHzBgAAAACUPkIpADdEfn6++vXrpzNnzkiSevXqpejoaHXu3Fkmk6nAnpPJGcrJt2nF1lPafSZV206myPqbx/TG3lVPoVW81KVBVVXycJHZXPD7AAAAAADKH0IpAH9Kfn6+vvzyS/Xt21fOzs5ydnbW5MmTtXXrVk2ePFlNmzYtsj8r16r+7/yolMy8Auvh1Xw1rlsDORFEAQAAAMAtiVAKwB+SlZWl9957T2+88Ybi4+O1bNkyDRw4UJI0duzYQvvyrTblWm3ycHHSgo3H9cFPJxwCqWp+7krPzpenm5PWRHWRj5tzoSusAAAAAADlH6EUgGK5dOmS3nnnHc2ePVsXL16UJAUEBCgjI+O6vYZhaODCn7X/rEXt6wZo3cELDvWwKl768plOsloNyST5uruUyDUAAAAAAG4ehFIAipSbm6vo6GgtWrTIHkDVqVNHEyZM0LBhw+Tp6Vlk/8KNx/TxLwk6mZwpSdcEUt+M66zbgn1YFQUAAAAAFQyhFIAiubq6asuWLcrIyFDz5s01efJkPfTQQ3J2LvyPj6MX0rTzVKrqVfXW1FUHHWrVK3noiY51VNPfU5m5+Qqv5lvSlwAAAAAAuAkRSgGwMwxDmzdv1ltvvaVFixbJ399fkhQbG6u0tDT17NmzyBVNS39J0Hf7EvXj0STl/+ZT9K56/4k26tqwaonNHwAAAABQfhBKAZDNZtNXX32lmJgYxcXFSZJatGihKVOmSJI6duxYZL/VZuhcapae/2xPgfU+zUPUqJqvujSocmMnDgAAAAAotwilgAosNzdXH330kaZPn66DB688Zufq6qqhQ4faP1GvML8cT9aa/ef1SLtaembpDh04Z3Go/1//Jqrs6aqk9Bw93r42e0YBAAAAABwQSgEVVE5OjsLDwxUfHy9J8vX11ahRozRu3DgFBwcX2pedZ9XplCxNXrlbJ5Mz9e7meId6sxp+6tk4WA+3qSlnJ3OJXgMAAAAAoPwilAIqkMuXL6tSpUqSJDc3N3Xp0kXZ2dkaP368RowYIV/fgjcdNwxD+TZDLk5mTVixS//dc86hXr2Sh7zdnHU2NUuzH2mp2gFeJX0pAAAAAIByzmQYxrW7EaPUWCwW+fn5KTU1tdBAAPirjhw5ojfeeENLlizRli1b1LRpU0lScnKyvL295ebmVmBfTr5Vxy5k6Nt9iZq97ojuvi1Iaw+ct9fDAr30SNtaerhtTXm7OctmSE5mHtMDAAAAgIqsuFkHK6WAW9ivv/6qmJgYrVy5Ulfz55UrV9pDqYCAgCL7X/vvAS2JO2l//dtAymySXh/QTG1D/e3HnMijAAAAAADFRCgF3GIMw9Dq1asVExOj9evX24/37t1bkydPVufOnQvtS87IVaC3m7bEX9JXu87qw5//F0jV8vdUeDUf7T9n0fzHIlTFx01VfdxL/HoAAAAAALcmQingFpOTk6MhQ4bo/PnzcnZ21iOPPKJJkybZV0cVZvHmeP3ffw/I09VJmbnWa+ov9g5Xj8aFb4AOAAAAAMAfQSgFlHOZmZlavny5Hn/8cTk5Ocnd3V3PPfec4uPjFRUVpVq1ahXZf+xiuuauP6aV209feb/fBVKjutbVsE6hCvQueN8pAAAAAAD+DEIpoJxKTk7WO++8o9mzZyspKUl+fn564IEHJEnjxo0rsvd0SqZmf39U9zYN1qzvj2h7wmWH+sQeDdSneYi+3n1OQzrUkbcbf1QAAAAAAG4sftMEypmEhATNnDlTixYtUmZmpiSpTp06xeo9nZKpmG8PaduJSzqbmq3lv55yqN/TOFhvPdxC7i5OkqTRd9a7oXMHAAAAAOAqQimgnMjMzNSIESP0ySefyGq98ohdixYtFB0drQcffFDOzkXfztl5Vr3x3SF9tevsNbX2YQEK8nXTuG4N7IEUAAAAAAAliVAKKCc8PDx06NAhWa1W3XXXXYqOjlb37t1lMpkKHJ+alaeZqw+pbWiAPtmSoM1HkxzqzWr46d3HW+vw+XS1C/OXi5O5NC4DAAAAAABJhFLATclms+mLL77Q3Llz9emnn6pSpUoymUx6++235ezsrDZt2hTZfzkzV4s3x+uDuJP6IO7kNfX+LUI0/I4wVfV1V1Vf95K6DAAAAAAACkUoBdxEcnJy9NFHH2n69Ok6dOiQJGn+/Pl69tlnJUnt27cvtHd7QorqVvHW9oQUPfn+VtmM/9Uqe7roodY1tfSXBD3XK1yD2hX9iXwAAAAAAJQ0QingJmCxWLRgwQK99dZbOnv2yp5Pfn5+GjVqlIYOHVpkb1J6jradTNGID7cpwMtVyRm5DvVH29XSw21qqWkNPz17722FPu4HAAAAAEBpIpQCylh6errCwsKUnJwsSQoJCdH48eP19NNPy9fXt8jeXacua8C8n5T//5dF/T6QGti6pl67v6n9NYEUAAAAAOBmQSgFlIHExEQFBwdLkry9vdWzZ09t375dkydP1qBBg+Tm5lZgX3pOvtydzfrxWLL2nL6sr3adswdSktSrabDcnZ00sE1NBfq4qXolj1K5HgAAAAAA/ihCKaAUbd26VTExMfr888+1d+9e3XbbbZKkuXPnysfHR2Zz4Z+Al5CcqXve3qjMXGuB9btuq6q5j0aUyLwBAAAAALjRCKWAEmYYhlavXq2YmBitX7/efvy7776zh1J+fn6F9lttht77MV5f7T5XYCBVr6q35j/WStX8WBUFAAAAACg/CKWAEpKfn68VK1YoNjZWu3btkiQ5Oztr0KBBmjRpkpo0aVJob3aeVX//aJuC/dxVzc9DM9ccdqg3qe6rf4/soP/uPqcWtSqpbhXvEr0WAAAAAABuNJNhGMb1h6GkWCwW+fn5KTU19bqbWqN8yczMVO3atZWUlCQvLy8NHz5c48ePV61atQrtMQxD//rxhDYfuaj1hy4WOGbN+DsUUslDXm5kygAAAACAm09xsw5+qwVukOTkZC1dulSjR4+W2WyWp6enpkyZIovFolGjRsnf37/Q3tx8m/KsNu0+napXv95/Tb1X02C5uzipdW1/1Q/yKcnLAAAAAACgVBBKAX9RQkKCZs6cqUWLFikzM1N16tRRnz59JEljx44tsnfF1lNKzcrTv36M17nU7GvqXz3TSbX8PeXn6VIicwcAAAAAoKwQSgF/0p49ezR9+nR98sknys/PlyS1bNlSnp6e1+01DEPHkzI0eeXuAuthVbzUt3mImtYofAN0AAAAAADKM0Ip4A9KTU3VoEGDtGrVKvuxu+++W9HR0erWrZtMJlOBffvOpiojx6oGQd564v2t2pFw2aH+Yu9w/XQsWTUre+if/QrfBB0AAAAAgFsBoRTwB/n6+urs2bMym80aMGCAoqOjFRERUeh4wzBkycrXwAU/Kz0n/5q6j5uzHr29tp7sFKqnOoeV5NQBAAAAALhpEEoBRcjJydFHH32kd999V6tXr5aPj49MJpMWLFggf39/1atXr8h+m83QQwvi9OvJlALrlTxdtDn6LnnzSXoAAAAAgAqG34SBAlgsFi1YsEBvvfWWzp49K0lavHixIiMjJUlt27YtsM8wDOXbDC2JO6ntJ1Pk5+niEEi5Opk1sE1NPdU5VC5OZkkikAIAAAAAVEj8Ngz8xrlz5zRr1izNmzdPqampkqTq1atr/PjxevLJJ6/bP3zJr1p74EKBNVcns9ZN7KIala+/EToAAAAAALc6Qing/7t06ZLq1aunzMxMSVJ4eLgmT56sQYMGydXVtcjeA+cs+uHQxWsCKVdns8wmad5jEari7UYgBQAAAADA/0cohQrt2LFjqlu3riTJ399fvXv31pkzZxQdHa377rtPZrO50N6Zaw5r9b5E9WkeorfWHlae1bDXzCbp89EdFezrruw8m2oFEEYBAAAAAPBbhFKocAzD0HfffaeYmBht3LhRhw4dsm9Y/v7778vTs+gA6WDilVVRs74/8v9fH3KoT+rZUF0aVFGT6n4lcwEAAAAAANwCCKVQYeTn52vFihWKjY3Vrl27JEnOzs768ccf7aFUYYFUUnqOtp1MUcd6gXrqg191OiXLoT76zrrydXfRwcQ0PdkpVO4uTiV7MQAAAAAAlHOEUrjlZWdn691339WMGTN04sQJSZKXl5eefvppRUZGqlatWoX2HruYrvUHL+g/289o/zmLvN2clZ6Tb6+P6BKmB1vVUP0gn5K+DAAAAAAAbimEUrjl5eXlacqUKbp8+bKqVKmisWPHatSoUfL3979u78RPd2lHwmX766uBlLPZpMYhvnqqU5iq+LiV1NQBAAAAALhlEUrhlnPy5EktXbpUzz77rEwmk3x8fPTKK6/IyclJTzzxhDw8PArsS07P0aYjSepcP1Bz1h/Vsi2nlJVntdf7NA/R7WH+SkjOVGS3BvJw5RE9AAAAAAD+LEIp3DJ2796t2NhYLVu2TFarVa1atVLPnj0lSWPGjLlu/z++3Kf/7j5XYK1JdV9N6tGQT9EDAAAAAOAGIZRCuWYYhjZu3KiYmBh988039uN33323KleuXKz3+M/203rpi31K+81eUU5mk6w2Q5I0/7EI3dMk+MZOHAAAAACACo5QCuXWhQsX1LdvX/3yyy+SJLPZrAcffFCTJ09WREREkb3f7DmnFz/fq97NqmlJ3EmHWqd6gZrUs6Fuq+ajoxfS1TjEr8SuAQAAAACAiopQCuWKYRgymUySpMDAQKWnp8vNzU1PPPGEJkyYoHr16hXZv2b/eU3+9y6lZOZJ0jWB1KSeDTX6zv+9B4EUAAAAAAAlg1AK5UJqaqoWLFigjz76SHFxcfLy8pLZbNaSJUtUvXp1BQUFFdmflJ6j3acv6/nP9tgDKUmq6e+h5PRcvXZ/E9lsVzYzBwAAAAAAJY9QCje1c+fO6a233tL8+fNlsVgkSR999JFGjBghSWrVqlWhvTsSUhS5fKdqB3gp7liS8qyGQ31ohzp6uW/jkps8AAAAAAAoFKEUbkqHDx/W9OnTtWTJEuXm5kqSwsPDNXnyZA0aNKjI3nyrTd/tO693Nx/XyeRMnUzOdKi3DfXX3yJqqHezaiU2fwAAAAAAUDRCKdx0zp49q0aNGslqtUqSOnTooOjoaN13330ym80F9qTn5Gvu+qPq0ThY3+5N1PwNxxzq7cMCFH3vbXp303GNvrOewqv5lvh1AAAAAACAwpkMwzCuPwwlxWKxyM/PT6mpqfL1rZhBiWEY2r17t5o3b24/9sADDyg/P1/R0dHq2LFjob02m6Evd53VL/HJ+mTLqWvqrs5mLRwcodvDAuTu4lQi8wcAAAAAAP9T3KyDlVIoM3l5eVqxYoViY2O1b98+HT16VHXq1JEkLV++XC4uLtd9j2VbT+n5z/Zcc9zZbFIVHzdN7NFQXRtWvdFTBwAAAAAAfxGhFEpdRkaGFi9erJkzZ+rkyZOSJG9vb+3YscMeShUUSBmGobUHLiikkrvmrj+mC2nZ2noixWHMu4+3Vp7Vps4NqsjbjR9vAAAAAABuVvzWjlJjsVg0c+ZMzZkzR8nJyZKkKlWqaNy4cRo1apQqV65cZP/3By5o+JJfC6y5OJk0qWdDdWsUdMPnDQAAAAAAbjxCKZQawzD05ptvymKxKCwsTBMnTtTQoUPl4eFR4PjLmbnycnNWSkauxq/YqR+PJjvUW9WqpO0Jl/W3iBqa/rfmBb4HAAAAAAC4ORFKocTs2rVLn376qV599VWZTCb5+fnp9ddfl7+/vwYMGCBn58J//E4kZeietzeqsqerXJ3NOpmc6VB/5s56mtizoS6kZcvf07WkLwUAAAAAANxgfPpeGbvVPn3PMAxt2LBBMTEx+vbbbyVJ33//ve66665i9afn5OvpJb/qp2PJBdbDqnjp+6guMplMN2zOAAAAAADgxuHT91CqrFarPv/8c8XExGjr1q2SJLPZrAcffFDBwcHX7f/paJIuZ+Xp8Pm0awKptqH+GtklTJcy8tS8hh+BFAAAAAAAtwBCKfxlCQkJ6t69uw4fPixJcnd31xNPPKEJEyaobt26hfbl5ts06uPtik9K17GLGdfUK3u66Ofn75abs1OJzR0AAAAAAJQNQin8ZTVq1JDJZFKlSpU0evRojR07VlWrVi2y5+zlLG04fFFrD5y/pubmbNb47g3UunZlAikAAAAAAG5RhFL4y8xms/7973+rdu3a8vHxKXTcvB+OKd9qU0Ttyhr63lblWm32mpuzWd9F3qF1By+oYbCPOtYLLI2pAwAAAACAMkIohRuiSZMmhdaOX0zXvrMWxXx7sMD68M6hujs8SHUCvTSsU2hJTREAAAAAANxECKVQImy2Kx/qmJqVp75zflR6Tr5DvVE1X2XlWVU7wFPP9wpn83IAAAAAACoYc1lPYO7cuQoNDZW7u7siIiK0adOmIsdv2LBBERERcnd3V1hYmObPn3/NmJUrV6pRo0Zyc3NTo0aN9Nlnn/3h8xqGoZdfflkhISHy8PBQ165dtW/fPocxOTk5GjNmjAIDA+Xl5aW+ffvq9OnTf+K7cOtISM7UntOpun/eT2r2z9Vq+eoah0Dq8fa19cGwtvrv2E5aP7Gr3n+iLYEUAAAAAAAVUJmGUsuXL1dkZKReeOEF7dixQ507d9a9996rhISEAsfHx8erV69e6ty5s3bs2KHnn39eY8eO1cqVK+1j4uLiNHDgQA0ePFi7du3S4MGD9dBDD+mXX375Q+eNjY3VzJkzNWfOHG3dulXBwcHq3r270tLS7GMiIyP12WefadmyZdq8ebPS09N13333yWq1lsB36+aXm2/TgPk/qc+czdp16vI1q6PCqnjphd7h6tKgCkEUAAAAAAAVnMkwDKOsTt6uXTu1atVK8+bNsx8LDw9X//79NW3atGvGR0dH68svv9SBAwfsx0aOHKldu3YpLi5OkjRw4EBZLBZ988039jH33HOPKleurE8++aRY5zUMQyEhIYqMjFR0dLSkK6uigoKCFBMToxEjRig1NVVVqlTRhx9+qIEDB0qSzp49q5o1a2rVqlXq2bNnsb4HFotFfn5+Sk1Nla+vb3G/dTcVq83Qok3H9dHPJ3U6Jct+vHXtysq12uTu4qT3hraRh4uTzGbCKAAAAAAAbmXFzTrKbKVUbm6utm3bph49ejgc79Gjh3766acCe+Li4q4Z37NnT/3666/Ky8srcszV9yzOeePj45WYmOgwxs3NTV26dLGP2bZtm/Ly8hzGhISEqEmTJoXOX7oSblksFoev8iwzN19//2ibXv/moEMgVb2ShxYPaaMvn+mkFSPay8vNmUAKAAAAAADYlVkolZSUJKvVqqCgIIfjQUFBSkxMLLAnMTGxwPH5+flKSkoqcszV9yzOea/+83pjXF1dVbly5WLPX5KmTZsmPz8/+1fNmjULHVse5OTZdPh8msOxL5/pqE2T75Sfp0sZzQoAAAAAANzsyvzT936/t5BhGEXuN1TQ+N8fL8573qgxv3e9Mc8995yioqLsry0WS7kOpip7uWrx0DZKyciV1WboUkaumtWoVNbTAgAAAAAAN7kyC6UCAwPl5OR0zaqiCxcuXLNC6arg4OACxzs7OysgIKDIMVffszjnDQ4OlnRlNVS1atUKHZObm6uUlBSH1VIXLlxQhw4dCr1uNzc3ubm5FVovj+pW8ZaqlPUsAAAAAABAeVJmj++5uroqIiJCa9ascTi+Zs2aQkOd9u3bXzN+9erVat26tVxcXIocc/U9i3Pe0NBQBQcHO4zJzc3Vhg0b7GMiIiLk4uLiMObcuXPau3dvkaEUAAAAAAAAyvjxvaioKA0ePFitW7dW+/bttXDhQiUkJGjkyJGSrjzqdubMGS1ZskTSlU/amzNnjqKiojR8+HDFxcVp8eLF9k/Vk6Rx48bpjjvuUExMjPr166cvvvhCa9eu1ebNm4t9XpPJpMjISE2dOlX169dX/fr1NXXqVHl6emrQoEGSJD8/Pz355JOaMGGCAgIC5O/vr4kTJ6pp06bq1q1baX0LAQAAAAAAyqUyDaUGDhyo5ORkvfLKKzp37pyaNGmiVatWqXbt2pKurDxKSEiwjw8NDdWqVas0fvx4vfPOOwoJCdGsWbM0YMAA+5gOHTpo2bJlevHFFzVlyhTVrVtXy5cvV7t27Yp9XkmaPHmysrKyNGrUKKWkpKhdu3ZavXq1fHx87GPefPNNOTs766GHHlJWVpbuvvtuvf/++3JycirJbxsAAAAAAEC5ZzKu7hSOMmGxWOTn56fU1FT5+vqW9XQAAAAAAAD+kuJmHWW2pxQAAAAAAAAqLkIpAAAAAAAAlDpCKQAAAAAAAJQ6QikAAAAAAACUOkIpAAAAAAAAlDpCKQAAAAAAAJQ6QikAAAAAAACUOkIpAAAAAAAAlDpCKQAAAAAAAJQ6QikAAAAAAACUOkIpAAAAAAAAlDpCKQAAAAAAAJQ6QikAAAAAAACUOkIpAAAAAAAAlDpCKQAAAAAAAJQ6QikAAAAAAACUOueynkBFZxiGJMlisZTxTAAAAAAAAP66qxnH1cyjMIRSZSwtLU2SVLNmzTKeCQAAAAAAwI2TlpYmPz+/Qusm43qxFUqUzWbT2bNn5ePjI5PJVNbT+VMsFotq1qypU6dOydfXt6ynA5Q57gnAEfcE4Ih7AnDEPQE4uhXuCcMwlJaWppCQEJnNhe8cxUqpMmY2m1WjRo2ynsYN4evrW25vGKAkcE8AjrgnAEfcE4Aj7gnAUXm/J4paIXUVG50DAAAAAACg1BFKAQAAAAAAoNQRSuEvc3Nz00svvSQ3N7eyngpwU+CeABxxTwCOuCcAR9wTgKOKdE+w0TkAAAAAAABKHSulAAAAAAAAUOoIpQAAAAAAAFDqCKUAAAAAAABQ6gil8JfNnTtXoaGhcnd3V0REhDZt2lTWUwJK3LRp09SmTRv5+PioatWq6t+/vw4dOuQwxjAMvfzyywoJCZGHh4e6du2qffv2ldGMgdI1bdo0mUwmRUZG2o9xT6CiOXPmjB577DEFBATI09NTLVq00LZt2+x17glUJPn5+XrxxRcVGhoqDw8PhYWF6ZVXXpHNZrOP4Z7ArWzjxo3q06ePQkJCZDKZ9PnnnzvUi/Pzn5OTozFjxigwMFBeXl7q27evTp8+XYpXceMRSuEvWb58uSIjI/XCCy9ox44d6ty5s+69914lJCSU9dSAErVhwwaNHj1aP//8s9asWaP8/Hz16NFDGRkZ9jGxsbGaOXOm5syZo61btyo4OFjdu3dXWlpaGc4cKHlbt27VwoUL1axZM4fj3BOoSFJSUtSxY0e5uLjom2++0f79+zVjxgxVqlTJPoZ7AhVJTEyM5s+frzlz5ujAgQOKjY3V9OnTNXv2bPsY7gncyjIyMtS8eXPNmTOnwHpxfv4jIyP12WefadmyZdq8ebPS09N13333yWq1ltZl3HgG8Be0bdvWGDlypMOx2267zXj22WfLaEZA2bhw4YIhydiwYYNhGIZhs9mM4OBg4/XXX7ePyc7ONvz8/Iz58+eX1TSBEpeWlmbUr1/fWLNmjdGlSxdj3LhxhmFwT6DiiY6ONjp16lRonXsCFU3v3r2NYcOGORx74IEHjMcee8wwDO4JVCySjM8++8z+ujg//5cvXzZcXFyMZcuW2cecOXPGMJvNxrfffltqc7/RWCmFPy03N1fbtm1Tjx49HI736NFDP/30UxnNCigbqampkiR/f39JUnx8vBITEx3uDzc3N3Xp0oX7A7e00aNHq3fv3urWrZvDce4JVDRffvmlWrdurb/97W+qWrWqWrZsqUWLFtnr3BOoaDp16qTvv/9ehw8fliTt2rVLmzdvVq9evSRxT6BiK87P/7Zt25SXl+cwJiQkRE2aNCnX94hzWU8A5VdSUpKsVquCgoIcjgcFBSkxMbGMZgWUPsMwFBUVpU6dOqlJkyaSZL8HCro/Tp48WepzBErDsmXLtH37dm3duvWaGvcEKprjx49r3rx5ioqK0vPPP68tW7Zo7NixcnNz0+OPP849gQonOjpaqampuu222+Tk5CSr1arXXntNjzzyiCT+nkDFVpyf/8TERLm6uqpy5crXjCnPv38TSuEvM5lMDq8Nw7jmGHAre+aZZ7R7925t3rz5mhr3ByqKU6dOady4cVq9erXc3d0LHcc9gYrCZrOpdevWmjp1qiSpZcuW2rdvn+bNm6fHH3/cPo57AhXF8uXL9dFHH2np0qVq3Lixdu7cqcjISIWEhGjIkCH2cdwTqMj+zM9/eb9HeHwPf1pgYKCcnJyuSWUvXLhwTcIL3KrGjBmjL7/8UuvXr1eNGjXsx4ODgyWJ+wMVxrZt23ThwgVFRETI2dlZzs7O2rBhg2bNmiVnZ2f7zz33BCqKatWqqVGjRg7HwsPD7R8Gw98TqGgmTZqkZ599Vg8//LCaNm2qwYMHa/z48Zo2bZok7glUbMX5+Q8ODlZubq5SUlIKHVMeEUrhT3N1dVVERITWrFnjcHzNmjXq0KFDGc0KKB2GYeiZZ57Rf/7zH61bt06hoaEO9dDQUAUHBzvcH7m5udqwYQP3B25Jd999t/bs2aOdO3fav1q3bq1HH31UO3fuVFhYGPcEKpSOHTvq0KFDDscOHz6s2rVrS+LvCVQ8mZmZMpsdf/10cnKSzWaTxD2Biq04P/8RERFycXFxGHPu3Dnt3bu3XN8jPL6HvyQqKkqDBw9W69at1b59ey1cuFAJCQkaOXJkWU8NKFGjR4/W0qVL9cUXX8jHx8f+fzX8/Pzk4eEhk8mkyMhITZ06VfXr11f9+vU1depUeXp6atCgQWU8e+DG8/Hxse+pdpWXl5cCAgLsx7knUJGMHz9eHTp00NSpU/XQQw9py5YtWrhwoRYuXChJ/D2BCqdPnz567bXXVKtWLTVu3Fg7duzQzJkzNWzYMEncE7j1paen6+jRo/bX8fHx2rlzp/z9/VWrVq3r/vz7+fnpySef1IQJExQQECB/f39NnDhRTZs2veYDZsqVMvvcP9wy3nnnHaN27dqGq6ur0apVK2PDhg1lPSWgxEkq8Ou9996zj7HZbMZLL71kBAcHG25ubsYdd9xh7Nmzp+wmDZSyLl26GOPGjbO/5p5ARfPVV18ZTZo0Mdzc3IzbbrvNWLhwoUOdewIVicViMcaNG2fUqlXLcHd3N8LCwowXXnjByMnJsY/hnsCtbP369QX+/jBkyBDDMIr385+VlWU888wzhr+/v+Hh4WHcd999RkJCQhlczY1jMgzDKKM8DAAAAAAAABUUe0oBAAAAAACg1BFKAQAAAAAAoNQRSgEAAAAAAKDUEUoBAAAAAACg1BFKAQAAAAAAoNQRSgEAAAAAAKDUEUoBAAAAAACg1BFKAQAAAAAAoNQRSgEAAJSwl19+WS1atCiz80+ZMkVPP/10mZ2/uObMmaO+ffuW9TQAAEApMRmGYZT1JAAAAMork8lUZH3IkCGaM2eOcnJyFBAQUEqz+p/z58+rfv362r17t+rUqVPq5/8jcnJyVKdOHX366afq1KlTWU8HAACUMOeyngAAAEB5du7cOfu/L1++XP/4xz906NAh+zEPDw95e3vL29u7LKanxYsXq3379mUeSFmtVplMJpnNhS/Ud3Nz06BBgzR79mxCKQAAKgAe3wMAAPgLgoOD7V9+fn4ymUzXHPv943tDhw5V//79NXXqVAUFBalSpUr65z//qfz8fE2aNEn+/v6qUaOG/vWvfzmc68yZMxo4cKAqV66sgIAA9evXTydOnChyfsuWLXN4JG7JkiUKCAhQTk6Ow7gBAwbo8ccft7/+6quvFBERIXd3d4WFhdnnd9XMmTPVtGlTeXl5qWbNmho1apTS09Pt9ffff1+VKlXS119/rUaNGsnNzU0nT57UDz/8oLZt28rLy0uVKlVSx44ddfLkSXtf37599fnnnysrK6tY338AAFB+EUoBAACUgXXr1uns2bPauHGjZs6cqZdffln33XefKleurF9++UUjR47UyJEjderUKUlSZmam7rzzTnl7e2vjxo3avHmzvL29dc899yg3N7fAc6SkpGjv3r1q3bq1/djf/vY3Wa1Wffnll/ZjSUlJ+vrrr/XEE09Ikr777js99thjGjt2rPbv368FCxbo/fff12uvvWbvMZvNmjVrlvbu3asPPvhA69at0+TJkx3On5mZqWnTpundd9/Vvn375O/vr/79+6tLly7avXu34uLi9PTTTzs8Atm6dWvl5eVpy5Ytf/2bDAAAbmqEUgAAAGXA399fs2bNUsOGDTVs2DA1bNhQmZmZev7551W/fn0999xzcnV11Y8//ijpyoons9msd999V02bNlV4eLjee+89JSQk6IcffijwHCdPnpRhGAoJCbEf8/Dw0KBBg/Tee+/Zj3388ceqUaOGunbtKkl67bXX9Oyzz2rIkCEKCwtT9+7d9eqrr2rBggX2nsjISN15550KDQ3VXXfdpVdffVUrVqxwOH9eXp7mzp2rDh06qGHDhrJarUpNTdV9992nunXrKjw8XEOGDFGtWrXsPVdXUF1vBRgAACj/2FMKAACgDDRu3Nhhf6WgoCA1adLE/trJyUkBAQG6cOGCJGnbtm06evSofHx8HN4nOztbx44dK/AcVx+Bc3d3dzg+fPhwtWnTRmfOnFH16tX13nvvaejQofYVS9u2bdPWrVsdVkZZrVZlZ2crMzNTnp6eWr9+vaZOnar9+/fLYrEoPz9f2dnZysjIkJeXlyTJ1dVVzZo1s7+Hv7+/hg4dqp49e6p79+7q1q2bHnroIVWrVs1hfh4eHsrMzCzeNxIAAJRbrJQCAAAoAy4uLg6vTSZTgcdsNpskyWazKSIiQjt37nT4Onz4sAYNGlTgOQIDAyVdeYzvt1q2bKnmzZtryZIl2r59u/bs2aOhQ4fa6zabTf/85z8dzrNnzx4dOXJE7u7uOnnypHr16qUmTZpo5cqV2rZtm9555x1JV1ZHXeXh4XHNpxO+9957iouLU4cOHbR8+XI1aNBAP//8s8OYS5cuqUqVKtf7FgIAgHKOlVIAAADlQKtWrbR8+XJVrVpVvr6+xeqpW7eufH19tX//fjVo0MCh9tRTT+nNN9/UmTNn1K1bN9WsWdPhXIcOHVK9evUKfN9ff/1V+fn5mjFjhn211+8f3StKy5Yt1bJlSz333HNq3769li5dqttvv12SdOzYMWVnZ6tly5bFfj8AAFA+sVIKAACgHHj00UcVGBiofv36adOmTYqPj9eGDRs0btw4nT59usAes9msbt26afPmzQW+35kzZ7Ro0SINGzbMofaPf/xDS5Ys0csvv6x9+/bpwIEDWr58uV588UVJV8Ku/Px8zZ49W8ePH9eHH36o+fPnX/ca4uPj9dxzzykuLk4nT57U6tWrdfjwYYWHh9vHbNq0SWFhYapbt+4f+fYAAIByiFAKAACgHPD09NTGjRtVq1YtPfDAAwoPD9ewYcOUlZVV5Mqpp59+WsuWLbM/BniVr6+vBgwYIG9vb/Xv39+h1rNnT3399ddas2aN2rRpo9tvv10zZ85U7dq1JUktWrTQzJkzFRMToyZNmujjjz/WtGnTinUNBw8e1IABA9SgQQM9/fTTeuaZZzRixAj7mE8++UTDhw//A98ZAABQXpkMwzDKehIAAAAoGYZh6Pbbb1dkZKQeeeQRh1r37t0VHh6uWbNmldHsHO3du1d33323Dh8+LD8/v7KeDgAAKGGslAIAALiFmUwmLVy4UPn5+fZjly5d0rJly7Ru3TqNHj26DGfn6OzZs1qyZAmBFAAAFQQrpQAAACqYOnXqKCUlRVOmTNHEiRPLejoAAKCCIpQCAAAAAABAqePxPQAAAAAAAJQ6QikAAAAAAACUOkIpAAAAAAAAlDpCKQAAAAAAAJQ6QikAAAAAAACUOkIpAAAAAAAAlDpCKQAAAAAAAJQ6QikAAAAAAACUOkIpAAAAAAAAlLr/Bwa8xMOkC3O7AAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "%matplotlib inline\n", "from matplotlib.pylab import plt\n", "\n", "fig = plt.figure(figsize=(12, 6))\n", "\n", "plt.plot(tgrid, kep_out[:, 4], label=\"Numerical integration\")\n", "# NOTE: offset slightly the line for the theoretical\n", "# prediction in order to improve the readability.\n", "plt.plot(tgrid, 1.0000035 + 2.08373e-4 * tgrid/100., 'k--', label=\"Theoretical prediction\")\n", "\n", "plt.xlabel(\"Time (years)\")\n", "plt.ylabel(r\"$\\omega$\")\n", "\n", "plt.legend()\n", "plt.tight_layout()" ] }, { "cell_type": "markdown", "id": "62bbd751-e951-4b24-90a7-129f5b9c3f46", "metadata": {}, "source": [ "Note that we added a slight offset to the theoretical prediction, because otherwise the two lines would have been overlapping, resulting in an unreadable graph.\n", "\n", "We can see how the numerically-computed precession matches very well the theoretical prediction. Note also how, in addition to the long-term linear behaviour, the numerical integration also highlights small short-term periodic oscillations which are also due to relativistic effects." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.10.13" } }, "nbformat": 4, "nbformat_minor": 5 }