From af307f768d01ccefe3f84a05b5b7542a07ce0896 Mon Sep 17 00:00:00 2001 From: David A Minton Date: Mon, 12 Jul 2021 16:55:19 -0400 Subject: [PATCH] Streamlined pseudovelocity conversion now that the methods are attached to the swiftest_body class --- .../whm_gr_test/swiftest_relativity.ipynb | 25 +++++------------ src/io/io.f90 | 28 ++++--------------- src/modules/swiftest_classes.f90 | 4 +-- src/modules/whm_classes.f90 | 4 +-- src/whm/whm_setup.f90 | 22 ++++++--------- src/whm/whm_step.f90 | 2 +- 6 files changed, 26 insertions(+), 59 deletions(-) diff --git a/examples/whm_gr_test/swiftest_relativity.ipynb b/examples/whm_gr_test/swiftest_relativity.ipynb index cd203cdd4..4026960a7 100644 --- a/examples/whm_gr_test/swiftest_relativity.ipynb +++ b/examples/whm_gr_test/swiftest_relativity.ipynb @@ -8,7 +8,6 @@ "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", - "import pandas as pd\n", "import swiftest\n", "from astroquery.jplhorizons import Horizons" ] @@ -76,7 +75,7 @@ "outputs": [], "source": [ "obj = Horizons(id='1', id_type='majorbody',location='@sun',\n", - " epochs={'start':'2021-01-28', 'stop':'3021-02-05',\n", + " epochs={'start':'2021-01-28', 'stop':'2121-02-05',\n", " 'step':'1y'})\n", "el = obj.elements()\n", "t = (el['datetime_jd']-el['datetime_jd'][0]) / 365.25\n", @@ -115,27 +114,17 @@ "output_type": "stream", "text": [ "Mean precession rate for Mercury long. peri. (arcsec/100 y)\n", - "JPL Horizons : 571.3210506300043\n", + "JPL Horizons : 573.8351991142854\n", "Swifter GR : 579.6804682138315\n", - "Swiftest GR : 579.6804681921402\n" - ] - }, - { - "ename": "ValueError", - "evalue": "operands could not be broadcast together with shapes (1000,) (100,) ", - "output_type": "error", - "traceback": [ - "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", - "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)", - "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 11\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34mf'Swifter GR : {np.mean(dvarpi_swifter)}'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 12\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34mf'Swiftest GR : {np.mean(dvarpi_swiftest)}'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 13\u001b[0;31m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34mf'Obs - Swifter : {np.mean(dvarpi_obs - dvarpi_swifter)}'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 14\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34mf'Obs - Swiftest : {np.mean(dvarpi_obs - dvarpi_swiftest)}'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 15\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34mf'Swiftest - Swifter: {np.mean(dvarpi_swiftest - dvarpi_swifter)}'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", - "\u001b[0;32m~/.local/lib/python3.7/site-packages/numpy/ma/core.py\u001b[0m in \u001b[0;36m__sub__\u001b[0;34m(self, other)\u001b[0m\n\u001b[1;32m 4151\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_delegate_binop\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mother\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 4152\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mNotImplemented\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 4153\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0msubtract\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mother\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 4154\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 4155\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0m__rsub__\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mother\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", - "\u001b[0;32m~/.local/lib/python3.7/site-packages/numpy/ma/core.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, a, b, *args, **kwargs)\u001b[0m\n\u001b[1;32m 1013\u001b[0m \u001b[0;32mwith\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0merrstate\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 1014\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mseterr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdivide\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m'ignore'\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0minvalid\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m'ignore'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1015\u001b[0;31m \u001b[0mresult\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mf\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mda\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdb\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 1016\u001b[0m \u001b[0;31m# Get the mask for the result\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 1017\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mma\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmb\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mgetmask\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mgetmask\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mb\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", - "\u001b[0;31mValueError\u001b[0m: operands could not be broadcast together with shapes (1000,) (100,) " + "Swiftest GR : 579.4895342016788\n", + "Obs - Swifter : -5.845269099546055\n", + "Obs - Swiftest : -5.654335087393189\n", + "Swiftest - Swifter: -0.19093401215286576\n" ] }, { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEGCAYAAAB2EqL0AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAA37UlEQVR4nO3dd3iUVfbA8e9JJ4QaAkgJXVqEAKFLEbCACFYEUUFXsa0F14ZdVnfVxYKVH90VpVdhFewUMUAgSI1AaKG3BAgJaef3xwwYQgIzkMmknM/z5GHeO++d97wBcnLLe6+oKsYYY8zF+Hg7AGOMMUWDJQxjjDEusYRhjDHGJZYwjDHGuMQShjHGGJf4eTsAT6pUqZLWrl3b22EYY0yRERMTc1hVw3J7r1gnjNq1a7Nq1Spvh2GMMUWGiOzM6z3rkjLGGOMSSxjGGGNcYgnDGGOMS4r1GEZu0tPTSUhIIDU11duhlDhBQUHUqFEDf39/b4dijLkEJS5hJCQkUKZMGWrXro2IeDucEkNVOXLkCAkJCdSpU8fb4RhjLkGJ65JKTU0lNDTUkkUBExFCQ0OtZWdMEVbiEgZgycJL7PtuTNFWIhOGMcYUV9HxRxi3dDue2LrCEoYXhISEsGPHDkqVKkVkZCRNmjTh4YcfJisrix07dhAREXHB+q+//jojRow4p6x27docPnzYrTh69epFYmKiu+EbYwqpubF7GPJlDJN+30lKema+f36JG/QuTOrVq0dsbCwZGRl069aNOXPm0LJlS49fV1VRVf73v/95/FrGGM9bsf0oU1fuZubqBGqFBvPZwJYEB+T/j3drYRQCfn5+dOjQga1bt+bL573//vtEREQQERHBhx9+CMCOHTto3Lgxjz76KC1btmT37t1nWyWjRo0iMjKSyMhI6tSpwzXXXAPA5MmTueqqq4iIiOD5558/+/khISG89NJLNG/enHbt2nHgwAEApk+fTkREBM2bN6dz5875ci/GmLwlnUrn05+30u//ljP/j70MbBvOoqGdaXxFWY9cr0S3MN74ZgMb9x7P189sUq0sr93U1K06p06d4scff2T48OEu1/nggw+YNGnS2eO9e/cCEBMTw4QJE4iOjkZVadu2LV26dKFChQrExcUxYcIEPvvss3M+6+GHH+bhhx8mPT2dbt268fTTT7N3716ef/55YmJiqFChAtdddx1z5szh5ptvJjk5mXbt2vHWW2/x3HPPMWbMGF5++WWGDx/OwoULqV69unV1GeNhm/cf56kpsWzef4IyQX78+HQXKpcN8ug1rYXhRdu2bSMyMpKOHTty44030rNnT5frDh06lNjY2LNf1apVA2Dp0qXccsstlC5dmpCQEG699VaWLFkCQK1atWjXrl2en/nkk0/SrVs3brrpJlauXEnXrl0JCwvDz8+PgQMHsnjxYgACAgLo3bs3AK1atWLHjh0AdOzYkcGDBzNmzBgyM/O//9QYA1sPnuC1ueu54cMlxB9O5t+3XsWioZ09niyghLcw3G0J5LczYxj56UIzI0qXLp3nexMnTmTnzp188sknF/0cf3//s1NkfX19ycjIAGDUqFFER0ezYMECIiMjiY2NJTQ09FJuwxiTi5U7jnLn/y0nS6FR1TKM7N+ChlXLFNj1rYVRzHTu3Jk5c+Zw6tQpkpOTmT17Np06dbpgnZiYGEaMGMGkSZPw8XH8k2jbti2//vorhw8fJjMzk8mTJ9OlS5cLfs62bdto27Ytw4cPp1KlSuzevTvf7suYkkpVycpShn+zkTtGLadi6QD+dnUdxg9uXaDJAgqwhSEiDYGp2YrqAq8CPwOjgBBgBzBQVc8bWBCRG4CRgC8wVlXf9nTMnpCRkUFgYOAFz4mLi6NGjRpnjz/44APuuOMOlz6/ZcuWDB48mDZt2gDwwAMP0KJFi7PdRrn55JNPOHr06NnB7qioKMaOHcu///1vrrnmGlSVXr160bdv3wte+9lnn2XLli2oKt27d6d58+YuxWyMyduzM/5gRkwCABVLBzB2UGsia5b3SiziiYc7LnpREV9gD9AWmAE8o6q/isj9QB1VfSWX8/8ErgUSgJXAAFXdeKHrREVFac4NlDZt2kTjxo3z7V7ctXbtWh588EFWrFjhtRi8ydvff2OKivV7khi7JJ45sXupXzmEW1pU59Gu9Ty+YoKIxKhqVG7veWsMozuwTVV3Olsei53l3wMLgVdynN8G2Kqq8QAiMgXoC1wwYRQ2o0aN4qOPPjo71dUYY3JKPp3Ba/M2nG1VVCkbyPSH2lOhdICXI/NewugPTHa+Xg/0AeYCdwA1czm/OpC9QzwBR+vkPCIyBBgCEB4enk/h5o8z01eNMSanzCzlp80HmRu7h/l/7KOUvy8zH+lAnUqlKRXg6+3wAC8kDBEJwJEghjmL7gc+EpFXgXlAWm7VcinLtS9NVUcDo8HRJXXZARtjjIedzsjk3//bzMTfdgDQo3FlPrmrJUH+hSNRnOGNFkZPYLWqHgBQ1c3AdQAiciVwYy51Eji35VED2OvhOI0xxqMyMrNITsvkxVnrWLBuHy3Dy/O3q+vSsX5ooUsW4J2EMYC/uqMQkcqqelBEfICXccyYymkl0EBE6uAYLO8P3FUQwRpjjKc8N+MPZq3ZA0CHeqF8elfLQjFWkZcCTRgiEoxjptND2YoHiMhjztezgAnOc6vhmD7bS1UzROTvOAbEfYHxqrqhAEM3xph8sycxhb9/vZo1uxK5tkkV2tcNZXCH2vj4FO49Ywr0wT1VPaWqoaqalK1spKpe6fx6QZ3zfFV1r6r2ynbe/5zn1FPVtwoy7vz21ltv0bRpU5o1a0ZkZCTR0dEu1Xv11Vf54YcfAFiyZAlNmzYlMjKS5cuX58vKswcOHOCuu+6ibt26tGrVivbt2zN79mwAfvnlF8qVK0eLFi1o1KgRzzzzzGVfz5iSJi0ji8cnr6H7e78QuzuRgW3D+XhAC+6/uk6hTxZQwpcG8Ybly5czf/58Vq9eTWBgIIcPHyYtLbdx/vNlX5zwq6++4plnnuG+++5j4sSJrFq1il69el2g9rkyMjLw8/vrr19Vufnmmxk0aBBff/01ADt37mTevHlnz+nUqRPz588nJSWFFi1acMstt9CxY0eXr2lMSaWqTPp9J8u2HuG7DfupXzmEEXc099oDeJfKEkYB27dvH5UqVTr7tHelSpUAWLFiBW+//TazZs1i7ty59O/fn6SkJLKysmjSpAnx8fEMHjyY3r17k5iYyLRp01i4cCGLFi1i2bJlpKSksHTpUoYNG0bv3r15/PHHWbduHRkZGbz++uv07duXiRMnsmDBAlJTU0lOTuann346G9dPP/1EQEDAOdN+a9WqxeOPP37ePZzZ+GnPnj0e/m4ZU/QdTU7j/e/jmPT7LgAe7VqP525o5OWoLk3JThjfvgD71+XvZ1a9CnrmvWrJddddx/Dhw7nyyivp0aMHd955J126dKFly5asWbMGcHQ3RUREsHLlSjIyMmjb9txHTh544AGWLl1K7969uf3228+2MM4sHPjiiy/SrVs3xo8fT2JiIm3atKFHjx6Ao4Xzxx9/ULFixXM+c8OGDS5v3nTs2DG2bNlie14YcwEHT6QyL3YvU1fuZsvBk3RtGMand7WkdGDR/bFbdCMvokJCQoiJiWHJkiX8/PPP3Hnnnbz99tsMHjyY+vXrs2nTJlasWMHTTz/N4sWLyczMvOjigTktWrSIefPmnd3GNTU1lV27HL/dXHvttecli9w89thjLF26lICAAFauXAk4ElmzZs2Ii4vjhRdeoGrVqm7evTElQ3pmFveOW8Hm/ScI8PVh/OAormlY2ePLenhayU4YF2gJeJKvry9du3ala9euXHXVVXzxxRcMHjyYTp068e233+Lv70+PHj0YPHgwmZmZ5+3ffTGqysyZM2nYsOE55dHR0Xkucd60aVNmzpx59vjTTz/l8OHDREX9taTMmTGMP//8k6uvvppbbrmFyMhIt2Izprj7Je4gw7/ZSPzhZG5oWpVnrr+S+pULdlVZT7HlzQtYXFwcW7ZsOXscGxtLrVq1AMfS5B9++CHt27cnLCyMI0eOsHnzZpo2vfC+HWXKlOHEiRNnj6+//no+/vjjs3tanOnqupBu3bqRmprK559/frbs1KlTuZ575ZVXMmzYMN55552Lfq4xJcXuo6e4+p2fGDxhJQiMvTeKUfe0KjbJAixhFLiTJ08yaNAgmjRpQrNmzdi4cSOvv/464NiD4sCBA2fHBpo1a0azZs0u2oy95ppr2LhxI5GRkUydOpVXXnmF9PR0mjVrRkREBK+8knMtx/OJCHPmzOHXX3+lTp06tGnThkGDBuWZFB5++GEWL17M9u3b3fsGGFPMpKZn8sVvOxg4NppjyWk80a0+cx/rSI8mVbwdWr7zyvLmBaUwLm9e0tn33xQnexJTeHRSDGsTkigT6MfE+9vQqlYFb4d1WQrj8ubGGFNknUrLYPXORB6fvJr0TOWtWyLoGXEFFQvxsh75wRKGMca4YeeRZG4ftZxDJ05TsXQA8/7egbphId4Oq0BYwjDGGBdNWLadd77bTKCfL0/1aECPxlVKTLIASxjGGHNR2afKtgwvz+t9mtKsRnlvh1XgLGEYY0wesm+XWrlMILe2qM4rvZsU6iXIPckShjHG5KCqTFm5mzFL4ok/lExwgC9ThrQrUd1PubHnMLygMC1vnpiYyGeffZbn+7bkuSlplm45zL3jVzBs1joOnTjNuEFRrH/9+hKfLMASRoHLvrz5H3/8wQ8//EDNmjUvXhHH8uZnFhE8s7x5bGwscXFxHkkYZ5Y879y5M/Hx8cTExDBlyhQSEhLOntOpUyfWrFnDmjVrmD9/PsuWLbukOIzxttT0TOas2cPd46JZsuUwtUODiX31Oro3rlIk9qooCAWWMESkoYjEZvs6LiJPiUikiPzuLFslIm3yqD9URDaIyHoRmSwiQQUVe37KbXnzatWqsWLFCm699VYA5s6dS6lSpUhLSyM1NZW6desCMHjwYGbMmMHYsWOZNm0aw4cPZ8CAAbz66qtMnTr17JPeycnJ3H///bRu3ZoWLVowd+5cwLEibZs2bYiMjKRZs2Zs2bKFF154gW3bthEZGcmzzz57Tqy25LkpCTIys4jdncg/pq3lqamxlCvlz/v9mjNlSHt8LVGco8DGMFQ1DogEEBFfHHtzzwbGAG+o6rci0gt4F+iava6IVAeeAJqoaoqITMOxr/fEy4npnRXvsPno5sv5iPM0qtiI59s8n+f73lzefNSoUTz55JMMHDiQtLQ0MjMzefvtt1m/fj2xsbHnxWpLnpuS4N2FcYxeHA9ARPWyfHhnZLFa/yk/eatLqjuwTVV3AgqUdZaXA/bmUccPKCUifkDwBc4r1M4sbz569GjCwsK48847mThxIn5+frkub75kyZJLWt787bffJjIykq5du55d3rx9+/b861//4p133mHnzp2UKlXKrc997LHHaN68Oa1btz5bdmbJ86pVq9K7d29b8twUGXsSUxg49ndGL46ndmgwb/RpyqxHOlqyuABvzZLqD0x2vn4KWCgiI3AksA45T1bVPc73dwEpwCJVXZTbB4vIEGAIQHh4+AWDuFBLwJO8tbx548aNadu2LQsWLOD6669n7NixZ7u7cmNLnpviKCtLeWhSDN9vPABAh3qhvH1rM8JDg70cWeFX4C0MEQkA+gDTnUWPAENVtSYwFBiXS50KQF+gDlANKC0id+f2+ao6WlWjVDUqLCzME7dwWby5vHl8fDx169bliSeeoE+fPvzxxx/n1c3Oljw3xc0fCYkMm7WO7zceoG6l0sx+tANfP9jOkoWLvNEl1RNYraoHnMeDgFnO19OB3Aa9ewDbVfWQqqY7zz+vJVIUeHN586lTpxIREUFkZCSbN2/m3nvvJTQ0lI4dOxIREXHeoLcteW6Ki4PHUxk6NZabP13G1FW7aVA5hEVDO9MivGivLFvQCnx5cxGZAixU1QnO403AI6r6i4h0B95V1VY56rQFxgOtcXRJTQRWqerHF7qWLW9e+Nj33xSk9Mwslm09zJgl8SzbeoQGlUP44M5I6oWFUCrA19vhFUqFZnlzEQkGrgUeylb8IDDSOZidinP8QUSqAWNVtZeqRovIDGA1kAGsAUYXZOzGmKLlVFoGL89Zz6zVjqneD3Wuyz+ua0iAnz1+dqkKNGGo6ikgNEfZUqBVLufuBXplO34NeM3TMRpjirbMLGXPsRT6/d9y9h9PpfEVZXm4S116RlxhyeIylci1pFT1ouMCJv8V590dTeGQlpHFnaOXs2ZXIkH+Prx8Y2PubleLIH/rfsoPJS5hBAUFceTIEUJDQy1pFCBV5ciRIwQFFckH9E0RMH3Vbr6K3kXs7kS6NarMk90b0LxmeW+HVayUuIRRo0YNEhISOHTokLdDKXGCgoKoUaOGt8MwxczB46m8OHs9P2w6QFiZQF7p3YS/XV3H22EVSyUuYfj7+1Onjv1jMqaoO5Gazpe/72Tayt3sPpbC7a1q8ObNEdb95EElLmEYY4q+k6czeGLyGn6OO0SZQD++fqAtbeuGXryiuSyWMIwxRUZ6ZhZz1uzhvUV/cvBEKg9cXYchXepSuYyNjRUESxjGmCIh8VQaA8ZEs2nfcepUKs20h9oTVbuit8MqUSxhGGMKtYzMLCav2MUXy3ey68gpnr2+IXe3q0W5Uv7eDq3EsYRhjCm0/jxwgqemxLJx33GCA3z5/O6WdG9cxdthlViWMIwxhc6Rk6f5dv1+3lywkZBAP0b2j6RP82r27JSXWcIwxhQqMTuP8cTkNexJTMHXR5gypD31K4d4OyyDJQxjTCFxPDWdGasSGD5/I2WC/Hjrlgha1KxgyaIQsYRhjPG6nzcfZMiXq0jPVK6qXo4x90ZRtZxNlS1sLGEYY7wmIzOLYbPWMT0mgfCKwfSMqMrgjrUtWRRSljCMMV7xxjcbmBmTwPHUDCqFBPDpXS25qkY5b4dlLsAShjGmQG3ce5zh8zfwe/xRKoUE8tnAZvS66gpvh2VcUGAJQ0QaAlOzFdUFXgV+AUYBQTh203tUVVfkUr88MBaIABS4X1WXezZqY0x+ST6dwdgl2xm/bDtJKelUKRvIoqe6UC7YHsArKgosYahqHBAJICK+wB5gNjAGeENVvxWRXsC7QNdcPmIk8J2q3i4iAUBwQcRtjLk8qsr3Gw8w8bcd/LbtCCGBfsx5rCP1wkpTJsiSRVHirS6p7sA2Vd0pIgqUdZaXA/bmPFlEygKdgcEAqpoGpBVMqMaYS5WWkcW8tXt5ZvpaAG5pUZ1/3hxBSKD1hhdFbv+tiUhpIFVVMy/juv2Byc7XTwELRWQE4AN0yOX8usAhYIKINAdigCdVNTmX+IYAQwDCw8MvI0RjzKVKTc9kb2IKT09bS+zuRK4oF8SbN0dwdYNKBPrZfhVF1UV3RBcRHxG5S0QWiMhBYDOwT0Q2iMh/RKSBOxd0dif1AaY7ix4BhqpqTWAoMC6Xan5AS+BzVW0BJAMv5Pb5qjpaVaNUNSosLMyd0Iwx+eSl2evp9t6vxO5OpGvDML78W1u6N65iyaKIc6WF8TPwAzAMWK+qWQAiUhG4BnhbRGar6iQXr9kTWK2qB5zHg4Anna+n4xjYzikBSFDVaOfxDPJIGMYY79m8/zj/92s8s9fsIbR0AB/cGUnnK+0Xt+LClYTRQ1XTcxaq6lFgJjBTRNwZuRrAX91R4Biz6IJjtlQ3YEsu19ovIrtFpKFz8Lw7sNGNaxpjPOhochp3j41m477jBPr5EFWrAmMHRVE+OMDboZl8dNGEkVuyuJRzAEQkGLgWeChb8YPASBHxA1Jxjj+ISDVgrKr2cp73OPCVs0srHrjPlWsaYzxHVZn0+05+2nyQjfuO06lBJd65rRnVypfydmjGA1we9BaRp3MpTgJiVDXWlc9Q1VNAaI6ypUCrXM7dC/TKdhwLRLkarzHGs7YfTuajH7cwe80eAO5uF86bN1/l5aiMJ7kzSyrK+fWN8/hGYCXwsIhMV9V38zs4Y0zhk5KWyQ+bDvDeojh2HDlF2zoVeeuWCOqF2aqyxZ07CSMUaKmqJwFE5DUcg8+dcUxztYRhTDGXmaX8/evV/Lj5IIF+Poy+pxXXNqliGxuVEO4kjHDOfVguHailqikicjp/wzLGFDY/bDzAiEVxbN5/gmsahvHSjU1sr4oSxp2E8TXwu4jMdR7fBEx2PshnM5aMKaYOHE9lyH9XsTYhibqVSjOyfyR9I6t7OyzjBS4nDFX9p4j8D7gaEOBhVV3lfHugJ4IzxnhPWkYWM1cnMHpxPAeOp/JEt/o82Lmurf9UgrkzS0qAxkA5VR0uIuEi0ia3lWWNMUXbkZOneXzyGn7bdoSyQX5MGNyatnVDL17RFGvudEl9BmTheLhuOHACx4N7rT0QlzHGC46npjNx2Q5GL47ndEYm79x2FTe3qG5LehjAvYTRVlVbisgaAFU95nyIzhhTDCSdSmfAmN/ZuO84VcoGMvORDjSsWsbbYZlCxJ2Eke7cx0IBRCQMR4vDGFOEpWVkMWZJPGOWxHPqdCbPXt+QO1rVoHJZ21fbnMudhPERjg2PqojIW8DtwMseicoYUyCi44/wxJQ1HDh+murlS/FBv0iuaVTZ22GZQsqdWVJfiUgMjoX/AG5W1U2eCcsY40lJp9IZtzSeMUu2U7lsIG/eHEG/qJoE+F10xwNTgl00YeSxhhRATxHpqarv53NMxhgPUVW+/H0n45ZuZ+eRUwT6+TCyfwsia5b3dmimCHClhXFm1KshjhlR85zHNwGLPRGUMSb/xe5OZOrKXUxesZtKIQFMvK81kTXL2xLkxmWuLG/+BoCILMKxltQJ5/Hr/LVrnjGmkEpNz2Tm6gRenrMeVahTqTTfD+2Mn691Pxn3XM5aUmlA7XyNxhiTb7KylN3HTvGfhXHM/2MfZQL9eKV3EzpdWcmShbkk7iSML4EVIjIbx9TaW4AvPBKVMeay/XPBRiYs2wFAo6pleL9fJE2qlfVuUKZIc2eW1Fsi8i3QyVl0n6qucbW+iDQEpmYrqgu8imNr1lFAEJABPJrXciPO50BWAXtUtber1zamJNmflMo/pseybOsRImuWZ0jnutzQtCo+PrYEubk8rsySElVVAFVdDay+0Dl5ce7FHek83xfYg+O5jjHAG6r6rYj0wrGvRtc8PuZJYBNgvyYZk0N6ZhavzdvA19G7AIiqVYGRA1pQ3bZLNfnElRbGzyIyE5irqrvOFDqXBbkaGAT8DEx047rdgW2qulNElL8SQDlgb24VRKQGjl3+3gLymuprTIn0c9xB3l/0J+v2JHFFuSBG3NGcjvUreTssU8y4kjBuAO7HsfdFHSARKAX4AIuAD1zd0zub/sBk5+ungIUiMsL5mR3yqPMh8Bx/TfPNlYgMAYYAhIeHuxmWMUXLydMZTFi6nfe+/xOAVrUqMOPh9rYDnvEIV6bVpuJYqfYzEfEHKgEpqpp4KRd0tkz6AMOcRY8AQ1V1poj0A8YBPXLU6Q0cVNUYEel6kXhHA6MBoqKiLthNZkxRlXw6g8krdvF19C7iDydTKSSAqQ+1p3r5UpYsjMe4M0sKVU0H9l3mNXsCq1X1gPN4EI6xCXA81zE2lzodgT7OMY4goKyITFLVuy8zFmOKpPcW/cn4ZdsBeLBTHe7rWIdqNlZhPMythJFPBvBXdxQ4xiy64Jgt1Q3YkrOCqg7D2SJxtjCesWRhSpqDJ1LZm5jK+9//yeI/D3F1/Uo8f0MjIqqXtVaFKRAFmjBEJBi4FngoW/GDwEgR8QNScY4/iEg1YKyq9irIGI0pjNIysrj98+XsOnqK4ABfBrSpydPXNiSsTKC3QzMliDtbtP5NVcddzsVU9RQQmqNsKdAql3P3AuclC1X9BUdrxJgSYfTibfy0+SC7jp6iY/1Q3ugTQf3KId4Oy5RA7rQw3hORgTgerlsBTFbVDZ4Jyxiz43AyL85ex2/bjgAwpHNdhvVsZN1PxmvcSRhHgDeBABwP4E0TkY9U9f88EZgxJVXSqXSmx+xmRkwCm/efoE/zarzXrzn+tv6T8TJ3EkaSqv7kfP2diIwEogFLGMbkk+Op6Tw0aRW/xx/F31cYPziKbo2qeDssY4BLGPQWkedxPItRDjiR7xEZUwJlZim/xx9h2Kx17ElMsX21TaF0KbOkZuJY2qMv8K/8DceYkufwydMMnrCC9XuOU6VsIFOHtCOqdkVvh2XMedxJGBVEpKaqbgW2isgYYA2wwDOhGVP8TVu5m89/3ca+pBSeu6Eht7e0VoUpvNxJGGWBX0TkMLARKA9keiIoY4q7HYeTeW3eBn798xCVywQy5t4oOjUI83ZYxlyQOwnjGmA90BbH/t6KtS6MccvR5DT+/b9NTI9JINDPh+dvaMSQznXxtb0qTBHgzgZKfzhfLnd+GWPccDw1nXvGRbNh73EqBPsz85EO1A2zB/BM0eGNtaSMKVF2Hz3Fhz9sYebqBPx8hPfuaM41jSpTsXSAt0Mzxi2WMIzxoN+2HWbw+JWkZWYRWbM8w3o2om3d0ItXNKYQcmctqb8DX6nqMQ/GY0yxcPJ0Bv/8ZiNTV+0mvGIwj3erzw0RVSkT5O/t0Iy5ZO60MKoCK0VkNTAeWHixfbyNKYlG/rCF8cu2k5SSTtWyQbzXrzmt7bkKUwy4M+j9soi8AlwH3Ad8IiLTgHGqus1TARpTVCzdcpjJK3axYN0+6lQqzcj+kXRtWNnbYRmTb9zdcU9FZD+wH8eqtRWAGSLyvao+54kAjSns0jOz+HL5Tv65YCOqUCkkkDmPdqRcsHU/meLFnTGMJ3Bsp3oYxzaqz6pquoj44Ngl74IJQ0QaAlOzFdUFXsWxt8UoHFuvZgCPquqKHHVrAv/F0S2WBYxW1ZGuxm6MJ2RlKUu2Hub7jfuZ9PsuSgf4Mv3hDlQvX8qShSmWXEoY4liAvzlwq6ruzP6eqmaJSO+LfYaqxuFYFh0R8QX2ALOBMcAbqvqtc8/ud4GuOapnAP9Q1dUiUgaIcbZqNroSvzH57XRGJuOX7uCd7zYD0K5uRT4b2MqmyppizaWE4eyKapEzWWR7f5Ob1+0ObFPVnSKiOJYdAccKuHtz+fx9wD7n6xMisgmojmOJEmMKjKqSkp7JoPErWLnjGPXCSvNGnwiuql7OWhWm2HNnDGO5iLRW1ZX5cN3+wGTn66eAhSIyAvABOlyooojUBlrg2IvDmAL1xJRYvlnr+J2mQeUQRt3Tinr2tLYpIdxdS+phEdkBJAOCo/HRzJ0LikgA0AcY5ix6BBiqqjNFpB8wDuiRR90QHMurP6Wqx/M4ZwgwBCA8PNyd0IzJ0y9xB/nPwjg27D1O3bDSDO1xJb2bXWHbpZoSRVx9lEJEauVWnlc31QU+py/wmKpe5zxOAso7u70Ex85+ZXOp5w/Mx/H8x/uuXCsqKkpXrVrlTnjGnCMlLZN3vtvMxN92EBLoR0T1sowb1JrSgbZIgimeRCRGVaNye8+df/WD8igf7mY8A/irOwocYxZdcMyW6oZjxtU5nIlkHLDJ1WRhzOXIyMzi+Znr2LjvOJv2HSe8YjBf3N+GOpVKezs0Y7zGnYSRnO11ENAbcGuwW0SCgWuBh7IVPwiMFBE/IBVnd5KIVAPGqmovoCNwD7BORGKd9V5U1f+5c31jXHEqLePsYoEAgzvU5vU+Tb0clTHe53KX1HkVRQKBeap6ff6GlH+sS8q4Y19SCt+t38/kFbv488BJompVYNQ9ragUEujt0IwpMPnVJZVTMI6H74wp8lSVB/+7ivV7jhMc4Ms7t11F38jqBPn7ejs0YwoNd570Xodjlz0AXyAM98cvjCl0pq3czcTfdrBx33EGtAnnkS71CA8N9nZYxhQ67rQwsj/NnQEcUNWMfI7HmAJzPDWdv01cycodxyjl78u/b72K/q1r2lRZY/Lgzmq1bk2fNaawSk3P5KvoXSxcv5+VO45xX8faPHd9I0oFWPeTMRfiTpfUF8CTqproPK4AvKeq93soNmPy3dHkNB74YiWrdyXi5yO8e1sz+rWu6e2wjCkS3OmSanYmWQCo6jERaZH/IRmT/1LSMpm1JoFPf9rK4eQ0PrmrBd0bVbFWhTFucCdh+IhIhTNbtIpIRTfrG+MVx5LTuGd8NOv3HKdskB9fP9CWKNsBzxi3ufMD/z3gNxGZgWO2VD/gLY9EZUw+UFW+it7Fhz/8yfHUDJ67oSF9mlejRgWbAWXMpXBnP4yfgVU4lu8QHHtj2PLiplBal5DEK3PXE7s7kYZVyvDRgBZ0qFfJ22EZU6S5sx/GHFVthe1BYQqxpFPpPD/zD77bsJ/KZQJ5/aYmDGgbTqCfjVUYc7nc6ZL6PR/3wzAm38XsPMaLs9YRd+AEQf4+THqgLVdWKePtsIwpNgp8Pwxj8lt0/BE++2Ubv/55CD8f4Yv729CpfiV8fOwBPGPykzsJo6fHojDmEm3ad5x7xq0gLTOLOpVK89UDbalWvpS3wzKmWHInYewCBgJ1VXW4iIQDVQF7AtwUuBOp6bw5fxNTV+2mQrA/o+9sRYvwCpQrZftqG+Mp7iSMz4AsHLOkhgMncGyX2toDcRmTpwnLtvPGN465F1G1KjCsVyNa1bLnKozxNHcSRltVbSkia+Dsk94BHorLmPPsOnKKd77bzIJ1+2hQOYSXbmxM14aVvR2WMSWGOwkjXUR8cS5xLiJhOFocLhGRhsDUbEV1gVdxbM06CscufhnAo6q6Ipf6NwAjcSytPlZV33YjdlOEZWUp//52ExOW7SAjS6lRoRSf3NWShlVtBpQxBcmdhPERMBuoIiJvAbcDr7haWVXjgEgAZ+LZ4/y8McAbqvqtiPQC3gW6Zq/rPP9THNu7JgArRWSePThY/M1ancC8tXv5Je4Qvj7CzEc60KpWBW+HZUyJ5M7y5l+JSAzQ3VnUV1U3X+J1uwPbVHWniChQ1lleDtiby/ltgK2qGg8gIlOAvthDhMXW8dR0Pvt5G6N+3QZArdBgfni6C/6+Pl6OzJiS66IJQ0Tm5Sxy/nm9iKCqfS7huv2Byc7XTwELRWQE4AN0yOX86sDubMcJQNs84h0CDAEIDw+/hNCMN6VlZLFq51E+/GELK7YfpVq5IEbd04qq5YIsWRjjZa60MNrj+GE9GYjmr4RxSZwD5X2AYc6iR4ChqjpTRPoB44AeOavl8lGaSxmqOhoYDRAVFZXrOabw+ujHLXzy81YAbo6sxsu9m1ApJNDLURljwLWEURXH2MEA4C5gATBZVTdc4jV7AqtV9YDzeBDwpPP1dGBsLnUSgOy73NQg964rU0St35PEQ1/GsCcxhXphpfnHdQ25vmlVfO1pbWMKjYsmDFXNBL4DvhORQByJ4xcRGa6qH1/CNQfwV3cUOH7wd8ExW6obsCWXOiuBBiJSB8dgeX8cycsUcVlZyv1frOSXuEOUCfLj1hbVeenGxoRaq8KYQsfV5c0DgRtx/LCvjWPG1Cx3LyYiwThaKw9lK34QGCkifkAqzvEHEamGY/psL1XNEJG/AwtxTKsdfxktHFMIqCovzl7P2t2JbNx3nAaVQxjZvwVNqpW9eGVjjFeI6oW7+Z17eUcA3wJTVHV9QQSWH6KionTVqlXeDsPkcCw5jXcXbmbyCsc8hrvbhfPPvhE4tl0xxniTiMSoalRu77nSwrgHx+q0VwJPZPtPfWa1WvuV0LgkKSWdMYvj+W7DfrYePEmb2hX58oE2tleFMUWEK2MYNpfRXLa0jCxembOeeWv34ucjjLq7JTdEXOHtsIwxbnDnSW9j3Jaansnc2D2M+jWe7YeTuaddLZ7q0cAGtY0pgixhGI9JTc/kjlHLWbcniVqhwXw+sCXXN61qGxsZU0RZwjD5LjNLGbskntlr9hB34AQv39iYu9qGExxg/9yMKcrsf7DJV/uTUnns69XE7DxGkL8P79zajH6ta168ojGm0LOEYfJFanomYxbH88EPf1LK35e3b72K/m1sLS9jihNLGOay7U1M4aEvY1i3JwlfH+GrB9sRWbO8t8MyxuQzSxjmkqVnZvF19C7eXLARXx/htZua0K1RZWqFlvZ2aMYYD7CEYS7J3sQUBo6NZvvhZCKql+WDfpE0qGI74BlTnFnCMG45eTqDR79azeI/D1Em0I8P7mxOz4grCPK3p7WNKe4sYRiXTV+1m7FLthN34ATVy5figzsjaVOnorfDMsYUEEsY5qJidh7j45+28EvcIcqV8ufTu1pyYzNb1sOYksYShslTZpby/cYDPDtjLSdSM6hWLojvhnambJC/t0MzxniBJQyTqw17k5iyYjdf/r4Tf19hzmMdaVilDKUCbKzCmJLKEoY5z9Ith7l7XDQAzWuW54N+zakbFuLlqIwx3lZgCUNEGgJTsxXVBV4F2gMNnWXlgURVjcyl/lDgAUCBdcB9qprqwZBLnIPHU/n75DWs2H6U8IrBjBsURf3KIbaxkTEGKMCEoapxQCSAiPji2Jt7tqp+eOYcEXkPSMpZV0SqA08ATVQ1RUSm4djXe6LHAy8hXpy9jq+jdwGOVsXzNzS05yqMMefwVpdUd2Cbqu48UyCOX2P7Ad3yqOMHlBKRdCAY2OvxKEuAeWv38s3avXy/8QA1K5big36RRNW2qbLGmPN5K2H0BybnKOsEHFDVLTlPVtU9IjIC2AWkAItUdVFuHywiQ4AhAOHhtvhdXk6ezuCZaWv5bsN+AEJLBzDzkQ5ULhPk5ciMMYVVgScMEQkA+gDDcrw1gPOTyJk6FYC+QB0gEZguIner6qSc56rqaGA0QFRUlOZf5MWDqjIpeheLNuxnyZbDhAT68c3jV3NFuSB7WtsYc0HeaGH0BFar6oEzBSLiB9wKtMqjTg9gu6oecp4/C+gAnJcwTN6OnDzN+GXb+fTnbQDc1rIGI+5oZoPaxhiXeCNh5NaS6AFsVtWEPOrsAtqJSDCOLqnuwCrPhVi8JJ/OYOWOo7w0ez17ElNofEVZ/nN7M5pcUdaShTHGZQWaMJw/8K8FHsrx1nljGiJSDRirqr1UNVpEZgCrgQxgDc5uJ3NhqspjX6/ml7hDBPn78PKNjbmjVU3KBdvT2sYY94hq8e3mj4qK0lWrSm5DZP4fe/nkp61s3n+Cbo0q8+z1DWl8RVlvh2WMKcREJEZVo3J7z570LoayspT3v/+TT37eSoVgf57oVp8nujfAz9fH26EZY4owSxjFSEZmFh//tJVf4g6yNiGJbo0q894dzalQOsDboRljigFLGMVEWkYWIxbFMXpxPACv9m7C/VfX8XJUxpjixBJGEZeUks64pduZsWo3e5NSubVldf5ze3N8fWz2kzEmf1nCKMJOpWVw7/gVrN2dSNWyQUy4rzVdrwyzqbLGGI+whFEEpaZnMj0mgS9+20H8oZO8flMTbm5RnfLBNlZhjPEcSxhFzM4jyQwav4IdR04RHODL27c1o19UTW+HZYwpASxhFBHHU9MZ+cMWxi3dToVgf0bf04oejavgY2MVxpgCYgmjCFi96xjDZq4j7sAJAnx9GH1vFK1tCXJjTAGzhFGI7UtKYcEf+3hzwSZK+fsy9t4omtUoR+WytgS5MabgWcIopDbtO87Nny7jdEYWDauUYdIDbQkrE+jtsIwxJZgljEImMWkf//r+INNWJVAh2J8J97WmRc0KlAqwvSqMMd5lCaOQyMrM4MVJvVmdnsC++OfpUK8Bj3drQPt6od4OzRhjAEsYhcIPy8fz7sb32ecntE6tyBO3t6d3y/reDssYY85hCcOLdu1azkc/vsIv7CfNF1omVeDZ2+fRuFp5b4dmjDHnsYThJV8t/JTxCZ9x0M8HEJ4s35P77n3H1oAyxhRaBZYwRKQhMDVbUV3gVaA90NBZVh5IVNXIXOqXB8YCEYAC96vqcs9F7Bnx23/i+7Uz+SRpMfj50PlEc4be+BD163TydmjGGHNBBZYwVDUOiAQQEV9gDzBbVT88c46IvAck5fERI4HvVPV2EQkAgj0acD5LTc9kztJ5fLjjZZJ9fAjMUv555T/o2OpuygbZdqnGmMLPW11S3YFtqrrzTIE4lljtB3TLebKIlAU6A4MBVDUNSCuQSPPB6dQk3pz8N+b6xIGPD91PNaFHs3707Hibt0MzxhiXeSth9Acm5yjrBBxQ1S25nF8XOARMEJHmQAzwpKom5zxRRIYAQwDCw8PzNWh3qSr//d9bvHdoCuojhKfBoFq30e+6170alzHGXIoC3+TZ2Z3UB5ie460BnJ9EzvADWgKfq2oLIBl4IbcTVXW0qkapalRYWFg+Re2+bduXcNuYSEYcnoqK0PlEfd7vudCShTGmyPJGC6MnsFpVD5wpEBE/4FagVR51EoAEVY12Hs8gj4ThbarKqxPv5UddzYlAH+qkwSON/07Pqx/ydmjGGHNZvJEwcmtJ9AA2q2pCbhVUdb+I7BaRhs7B8+7ARg/H6bbx3zzET/tWsTYwDT8Vrk+N4j9DxtsOeMaYYqFAE4aIBAPXAjl/3T5vTENEqgFjVbWXs+hx4Ctnl1Y8cJ+Hw3VZenoGo7/5F6NO/AaBUC1dmXzH91Qsd4W3QzPGmHxToAlDVU8B5y2OpKqDcynbC/TKdhwLRHkwPLdpVhbRa8Ywes04VvqnEJSlPBD2BIN79CcwqKy3wzPGmHxlT3pfIlVl9LzHHQ/g+UOnlDAe7v4mzRp08HZoxhjjEZYwLsHaDbP4YOnbrA84RViWcluFQTw48BkC/Ap80pkxxhQYSxhuennCbXwjcWQFCR3TyvNkj09oXCfS22EZY4zHWcJw0W8xY/hy9X9ZGpBI05OB3NT4nwzs1uviFY0xppiwhHERR0+e4q2vB/J9wBY0QIg6XYqPB39PSOly3g7NGGMKlCWMC4he8wX/XvUe2wKV5ieD6R/1Nje27or42FiFMabksYSRiyMnUhkx4w7m++yAAOjvG8kT94+jTKkAb4dmjDFeYwkjh2PJadwyrRXH/HxonhJIt5oDuf/God4OyxhjvM4SRg5lArKIkqqEB9fk8YHj8PXz9XZIxhhTKFjCyMHPP4j37//R22EYY0yhY6O3xhhjXGIJwxhjjEssYRhjjHGJJQxjjDEusYRhjDHGJZYwjDHGuMQShjHGGJdYwjDGGOMSUVVvx+AxInII2HmJ1SsBh/MxnKLA7rn4K2n3C3bP7qqlqmG5vVGsE8blEJFVqlqo9hD3NLvn4q+k3S/YPecn65IyxhjjEksYxhhjXGIJI2+jvR2AF9g9F38l7X7B7jnf2BiGMcYYl1gLwxhjjEssYRhjjHGJJYwcROQGEYkTka0i8oK348kvIlJTRH4WkU0iskFEnnSWVxSR70Vki/PPCtnqDHN+H+JE5HrvRX/pRMRXRNaIyHzncbG+XwARKS8iM0Rks/Pvu31xvm8RGer8N71eRCaLSFBxvF8RGS8iB0VkfbYyt+9TRFqJyDrnex+JiLgchKral/ML8AW2AXWBAGAt0MTbceXTvV0BtHS+LgP8CTQB3gVecJa/ALzjfN3Eef+BQB3n98XX2/dxCff9NPA1MN95XKzv13kvXwAPOF8HAOWL630D1YHtQCnn8TRgcHG8X6Az0BJYn63M7fsEVgDtAQG+BXq6GoO1MM7VBtiqqvGqmgZMAfp6OaZ8oar7VHW18/UJYBOO/2x9cfyAwfnnzc7XfYEpqnpaVbcDW3F8f4oMEakB3AiMzVZcbO8XQETK4vjBMg5AVdNUNZHifd9+QCkR8QOCgb0Uw/tV1cXA0RzFbt2niFwBlFXV5erIHv/NVueiLGGcqzqwO9txgrOsWBGR2kALIBqooqr7wJFUgMrO04rD9+JD4DkgK1tZcb5fcLSODwETnF1xY0WkNMX0vlV1DzAC2AXsA5JUdRHF9H5z4e59Vne+zlnuEksY58qtL69YzTsWkRBgJvCUqh6/0Km5lBWZ74WI9AYOqmqMq1VyKSsy95uNH45ui89VtQWQjKOrIi9F+r6dffZ9cXS7VANKi8jdF6qSS1mRuV835HWfl3X/ljDOlQDUzHZcA0fztlgQEX8cyeIrVZ3lLD7gbKbi/POgs7yofy86An1EZAeOrsVuIjKJ4nu/ZyQACaoa7TyegSOBFNf77gFsV9VDqpoOzAI6UHzvNyd37zPB+TpnuUssYZxrJdBAROqISADQH5jn5ZjyhXMmxDhgk6q+n+2tecAg5+tBwNxs5f1FJFBE6gANcAyWFQmqOkxVa6hqbRx/jz+p6t0U0/s9Q1X3A7tFpKGzqDuwkeJ737uAdiIS7Pw33h3H+Fxxvd+c3LpPZ7fVCRFp5/x+3ZutzsV5e+S/sH0BvXDMINoGvOTtePLxvq7G0fT8A4h1fvUCQoEfgS3OPytmq/OS8/sQhxszKQrbF9CVv2ZJlYT7jQRWOf+u5wAVivN9A28Am4H1wJc4ZgYVu/sFJuMYp0nH0VL426XcJxDl/F5tAz7BueKHK1+2NIgxxhiXWJeUMcYYl1jCMMYY4xJLGMYYY1xiCcMYY4xLLGEYY4xxiSUMYy5CREJFJNb5tV9E9jhfnxSRzzx0zadE5N6LnDNFRBp44vrG5Mam1RrjBhF5HTipqiM8eA0/YDWO1YUzLnBeF+BuVX3QU7EYk521MIy5RCLSNds+G6+LyBciskhEdojIrSLyrnPfge+cy7Kc2YvgVxGJEZGFZ5Z1yKEbsFpVM0SknoisznbNBiJyZn2sJUAPZ4IxxuMsYRiTf+rhWE69LzAJ+FlVrwJSgBudSeNj4HZVbQWMB97K5XM6AjEAqroNSBKRSOd79wETne9l4Vi2urmH7seYc9hvJsbkn29VNV1E1uHYjOs7Z/k6oDbQEIgAvnducuaLY6mHnK7AsR7SGWOB+0TkaeBOzt2/4SCOVVpdXZXXmEtmCcOY/HMaHL/5i0i6/jVAmIXj/5oAG1S1/UU+JwUIynY8E3gN+AmIUdUj2d4Lcp5vjMdZl5QxBScOCBOR9uBYbl5EmuZy3iag/pkDVU0FFgKfAxNynHslsMEz4RpzLksYxhQQdWz7ezvwjoisxbFicIdcTv0Wxzar2X2FY7XhRWcKRKQKkKLOHdeM8TSbVmtMISQis4HnVHWL8/gZoJyqvpLtnKHAcVUd56UwTQljYxjGFE4v4Bj83uJMHvVwTLfNLhHH/g/GFAhrYRhjjHGJjWEYY4xxiSUMY4wxLrGEYYwxxiWWMIwxxrjEEoYxxhiX/D/hGs5n3/VMywAAAABJRU5ErkJggg==\n", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEGCAYAAABPdROvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAABCrElEQVR4nO3dd3iURdfA4d/ZFFIooYROCKFDgACRItKLUkRBRUFfxYbYy2fBhr5YXlRUbKCIAooiIlJERRClCdJ7M3RCINTQQiDJnu+PXTCGgFlSNuXc18WV3XnamZDsyTMzz4yoKsYYY8zlcng7AGOMMfmbJRJjjDFZYonEGGNMllgiMcYYkyWWSIwxxmSJr7cD8IYyZcpoeHi4t8Mwxph8ZcWKFYdUNTR9eaFMJOHh4SxfvtzbYRhjTL4iIrsyKremLWOMMVliicQYY0yWWCIxxhiTJYWyjyQjycnJxMbGkpSU5O1QCqWAgAAqV66Mn5+ft0MxxnjIEolbbGwsxYoVIzw8HBHxdjiFiqpy+PBhYmNjqVatmrfDMcZ4yJq23JKSkihdurQlES8QEUqXLm13g8bkU5ZI0rAk4j32vTcm/7JEYowxhcDxY3t47ZsenDi+N9vPbYkkDylatCg7d+4kMDCQqKgo6tWrx8CBA3E6nezcuZPIyMhLHv/yyy8zbNiwf5SFh4dz6NAhj+Lo1q0bCQkJnoZvjMmj1m6YyE2Tu/Jd0k6mLPgs289vne15UPXq1Vm9ejUpKSl06NCBqVOn0qRJkxy/rqqiqvz00085fi1jTM5zpqYw9qcBfHB4KaWdSuSxPlx706Bsv47dkeRhvr6+XHnllWzdujVbzvfOO+8QGRlJZGQkw4cPB2Dnzp3UrVuXBx54gCZNmrBnz57zdzEff/wxUVFRREVFUa1aNdq3bw/AhAkTaNCgAZGRkTzzzDPnz1+0aFGef/55GjVqRIsWLYiPjwdg0qRJREZG0qhRI9q0aZMtdTHGXNqePYu5+8sWvHtkGfUTAzm7/0Vev+MpSgb7Z/u17I4kA//9YQMb445n6znrVSzOS9fW9+iYxMRE5syZw5AhQzJ9zLvvvsv48ePPv4+LiwNgxYoVjBkzhiVLlqCqNG/enLZt21KyZEm2bNnCmDFjGDFixD/ONXDgQAYOHEhycjIdOnTgiSeeIC4ujmeeeYYVK1ZQsmRJunTpwtSpU7n++us5deoULVq04LXXXuPpp5/m008/5YUXXmDIkCH88ssvVKpUyZrMjMlh6nQyafbjDIubgw9w1fEmzDt0C9/e14oqpYJy5Jp2R5IHbdu2jaioKFq1akX37t3p2rVrpo99/PHHWb169fl/FStWBGDhwoX06tWL4OBgihYtSu/evVmwYAEAVatWpUWLFhc956OPPkqHDh249tprWbZsGe3atSM0NBRfX19uvfVW5s+fD4C/vz89evQAoGnTpuzcuROAVq1a0b9/fz799FNSU1Mv51tijMmkSb8+wSv7fyPKEUT/CkP5Oe5mPuzXlMhKJXLsmnZHkgFP7xyy27k+kuykqhfdFhwcfNFtY8eOZdeuXXz44Yf/eh4/P7/zw3h9fHxISUkB4OOPP2bJkiX8+OOPREVFsXr1akqXLn051TDGXELC0R28v/dXmkkAH/9nMZ2GL6BxmB8d6pTL0evaHUkh0aZNG6ZOnUpiYiKnTp1iypQptG7d+pLHrFixgmHDhjF+/HgcDtePSvPmzZk3bx6HDh0iNTWVCRMm0LZt20ueZ9u2bTRv3pwhQ4ZQpkwZ9uzZk231Msb87f1f7uekwLOtX+fPHUfZfvAUtzWvmuPXtTuSPCIlJYUiRYpccp8tW7ZQuXLl8+/fffddbrrppkydv0mTJvTv359mzZoBcM8999C4cePzzU8Z+fDDDzly5Mj5Tvbo6GhGjx7N//73P9q3b4+q0q1bN6677rpLXvupp54iJiYGVaVjx440atQoUzEbYzJvw+bv+S4pltuCq1Ojehfe+WoFIUF+dG9YIcevLZdqqiiooqOjNf3CVps2baJu3bpeigjWrFnDvffey9KlS70Wg7d5+//AmPzKmZrCf764gr2azA83/sJpSnPl0N+4s1U4z3evl23XEZEVqhqdvtyatvKAjz/+mL59+/Lqq696OxRjTD6jTifDp9zEWkcKT4T3pFjxSnyzbA8pTqVfLjRrgTVt5QnnhtkaY0x6MVtnMnXVSIr7F6dkQCjFgiLo1OxO/IoEk5pyllcm9WDy2X30KVKJa9u9SkqqkwlLd9O6Zhmqlbn4QJrslGuJRERqAxPTFEUAg4GWQG13WQiQoKpRGRwfAowGIgEF7lLVxWm2Pwm8BYSqqmdzghhjTB60a9cC7lnwJMcFUhL/nti0WMxIrnSU46zjDL/rce4tVpeHr/8GcTgYPW8b+44l5ero01xLJKq6BYgCEBEfYC8wRVWHn9tHRN4Gjl3kFO8BM1X1RhHxB4LSHFcF6AzszpHgjTEmlx06uIn75jyAE2h08gEWxZUluvwxGpfdzvqjc1gQGE+iOrgzoCmP9B4LwBeLdzL05810a1CeLvVydshvWt5q2uoIbFPVXecKxPUAQh+gQ/qdRaQ40AboD6CqZ4GzaXZ5F3gamJZzIRtjTO44eWIf98+4hSOihMb2YqOjFq/fWIfejSvhcAipzleZuWYLn/22mPcPlmWXrqJBpRK8+uMmOtUtx3u3NMbhyL2lGbyVSG4BJqQraw3Eq2pMBvtHAAeBMSLSCFgBPKqqp0SkJ7BXVddcak0LERkADAAICwvLhioYY0zOGP7TPWyVVKrGtqNYxR58cXMUpYv+/XiAj0Po3rgOnRrUZMTv2xgxdyvTVsfRtlYoH93aGD+f3B1HleujttzNUj2BSek29eXC5HKOL9AEGKmqjYFTwCARCQKex9XXckmqOkpVo1U1OjQ09LLjz0mvvfYa9evXp2HDhkRFRbFkyZJMHTd48GB+/fVXABYsWED9+vWJiopi8eLF2TKTb3x8PP369SMiIoKmTZvSsmVLpkyZAsDcuXMpUaIEjRs3pk6dOjz55JNZvp4xhdnJE/v4IXEn9U8UpU3zRxh7Z7N/JJG0ivj68HjnWvz4SGuevqY2n/ynKUV8fXI5Yu8M/+0KrFTV+HMFIuIL9OafnfFpxQKxqnruk/U7XImlOlANWCMiO4HKwEoRKZ9DseeYxYsXM2PGDFauXMnatWv59ddfqVKlSqaOHTJkCJ06dQLgq6++4sknn2T16tVs2bLF40RyblqTc1SV66+/njZt2rB9+3ZWrFjBN998Q2xs7Pl9WrduzapVq1i1ahUzZszgjz/+8Oiaxpi/ffLT8yQ6HDStcDtPdK6FTyaaqGqVK8YD7WoQ4Jf7SQS8k0gyuvPoBGxW1dgM9kdV9wN73CO/wNXHslFV16lqWVUNV9VwXAmniXv/fGXfvn2UKVPm/NPtZcqUoWLFiixdupTevXsDMG3aNAIDAzl79ixJSUlEREQA0L9/f7777jtGjx7Nt99+y5AhQ+jbty+DBw9m4sSJREVFMXHiRE6dOsVdd93FFVdcQePGjZk2zdWlNHbsWG666SauvfZaunTp8o+4fvvtN/z9/f8xPLlq1ao8/PDDF9Th3IJce/dm/wpsxhQG8cdOMfv4UmqeER7ulX8eCcjVPhJ3U1Rn4L50my7oMxGRisBoVe3mLnoY+MrdNLYduDPHAv15EOxfl73nLN8Aug696OYuXbowZMgQatWqRadOnbj55ptp27YtTZo0YdWqVYCr2SoyMpJly5aRkpJC8+bN/3GOe+65h4ULF9KjRw9uvPFGxo4dy/Lly89PuPjcc8/RoUMHPv/8cxISEmjWrNn5O5nFixezdu1aSpUq9Y9zbtiwIdOLah09epSYmBhbc8SYy6CqvDXxdfYWEZ4JvQbfXO7nyIpcTSSqmghcMO2rqvbPoCwO6Jbm/Wrggkfz0x0TntUYvaVo0aKsWLGCBQsW8Pvvv3PzzTczdOhQ+vfvT40aNdi0aRNLly7liSeeYP78+aSmpv7rpIvpzZo1i+nTp59fjjcpKYndu10jpjt37nxBEsnIgw8+yMKFC/H392fZsmWAK8E1bNiQLVu2MGjQIMqXz3cti8Z43dTVe4lNnklJXyd9Or7o7XA8Yk+2Z+QSdw45ycfHh3bt2tGuXTsaNGjAuHHj6N+/P61bt+bnn3/Gz8+PTp060b9/f1JTUy9Yn/3fqCqTJ0+mdu3a/yhfsmTJRaeSr1+/PpMnTz7//qOPPuLQoUNER/+d01u3bs2MGTP466+/uOqqq+jVqxdRUVEexWZMYXbo5BlG/jSN/ZXOcE/x+vgXKebtkDySf+6dCrgtW7YQE/P3yOfVq1dTtaprnpw2bdowfPhwWrZsSWhoKIcPH2bz5s3Ur3/pJ1eLFSvGiRMnzr+/+uqr+eCDD86vKXKuyexSOnToQFJSEiNHjjxflpiYmOG+tWrV4tlnn+WNN9741/MaY/728vQNlC/2NQ6gT8tnvR2Ox+yOJI84efIkDz/8MAkJCfj6+lKjRg1GjRoFuNYAiY+PP9/30LBhQ8qWLculnpsBaN++PUOHDiUqKopnn32WF198kccee4yGDRuiqoSHhzNjxoxLnkNEmDp1Ko8//jhvvvkmoaGhBAcHXzRZDBw4kGHDhrFjxw6qVat2Gd8JY/KvQ4c2s2DN56w7uJb1p/fjAJoWrUpoYFNKl7+ero0jL+j7mL0xnrgdo4ipfJz/BFWnfIXG3gk+C2waeTebwtz77P/A5GeJJw/Q69uOxPlAMadS3xHI6ZRUNjnOctYhlEh1UvtYE65s+jLXNq5E0SK+JKcofd77glPl3qG8+PJ1v4V5ulnrYtPI2x2JMcZkgxE/DyDOB0bUvZdWTR+g/7iVLIw5SIPyPvSquZ5f93/F0lKrObqxN1PmdmZXUiNOOosTXfUDjojwZrt383QSuRRLJMYYk0Wbtkxj/Kmt3BhQidbNHmHpjiPM/+sg/9e5Fg91qIHINfRPfYypvw/inT0/c6zqL8AvVEx1stnHwUsVuxBR7YJpBvMNSyTGGJMFqSlnGbLoZUIUHuv6KQAf/BZDmaL+3NM64nxfpsPHl96dhtEx4SFWbPyWHYc3sePEHsoGlOKGjp6NwMxrLJEYY0wWTPz1cdY7Ungz/EZKlAhj1e6jLIg5xKCudQj0v3DKkhIh4XS48mkvRJpzbPivMcZcplMn9/PRvnm0JJBrWrvmjv3gt62EBPlxW4vcWeY2L7BEYowxl+nbuc9z3CE8csVTiMPB+r3H+G3zAe5uVY2iRQpPg48lkjwkL00jn5CQwIgRIy663aaWN4Vd0umjjDuwhJYEElnvJgDemxNDsQBf7mgV7t3gcpklkjwir0wjf86lEolNLW8MTJ03mMM+wr0NXXPQLog5yOyN8dzXJoLiAX5eji53WSLJI7w5jfyGDRto1qwZUVFRNGzYkJiYGAYNGsS2bduIioriqaee+kesNrW8KQxSncrM9ftYtPUQyanOf2xLTk5kTNxcotSP6EZ3cjbFyUvTN1C1dBD3tI7wUsTeU3ga8TzwxtI32Hxkc7aes06pOjzT7JmLbvfmNPIff/wxjz76KLfeeitnz54lNTWVoUOHsn79elavXn1BrDa1vCnoFm07xKvTV+Fz6icOJVXjlG9t2tcpS6vqZWgcFsL69a8S5wPP1/kP4nDw2fxtbD94ijF3XuG1xaW8yRJJHuHNaeRbtmzJa6+9RmxsLL1796ZmzZoendemljcFhdOpPPflJ8Qd/ZrjJQ5zrJQDH1Vang1h247r2LwlgbIh84kpfoSIVB/OBN/C7sOJvD8nhs71ytG+dllvV8ErLJFk4FJ3DjnJW9PI161bl+bNm/Pjjz9y9dVXM3r06PPNZhmxqeVNQTX5j9n86vwQZwi09SlJt2rdWLnvTybrdk5X/AKAw06lcXJx9h26kXvHr8Qh4OfjYHCPet4N3otyLZG4l8lNuyZ7BDAYaAmc+2QLARJUNSqD40OA0UAkoMBdqrpYRN4CrgXOAtuAO1U1IUcqkYO2bNmCw+E4fzeQfhr522+/ndtvv/38NPL79++/7GnkP/jgA0SEVatW0bhxY7Zv305ERASPPPII27dvZ+3atTRq1Ogfx6bVoUMHnnvuOUaOHMn9998PZG5q+QkT0q+wbEzekepUZq59ldRgmNrhU6qGtQRcS7oOOLKNKX+8RkhgKTo3e5xixSuRnOpk/l8H+WFNHFfWKEOVUkHerYAX5Vpnu6puUdUod5JoCiQCU1T15jTlk4HvL3KK94CZqloHaARscpfPBiJVtSHwF5D/JvPHNY38HXfcQb169WjYsCEbN27k5ZdfBjKeRr5hw4aZmkZ+48aN5zvbX3zxRZKTk2nYsCGRkZG8+KJrFbaJEycSGRlJVFQUmzdv5vbbb6d06dK0atWKyMjICzrbz00tP2/ePKpVq0azZs244447Ljm1/Pz589mxY0cWv0vG5JzvFs5mZfAROkv580nknJKlqnPXtZ/Tu9MwihWvBLjuQjrWLcfwWxrTJzpzIywLKq9MIy8iXYCXVLVVmjIBdgMdVDUm3f7FgTVAhF4iYBHpBdyoqrde6vo2jXzeZP8HxltSncqAkW1ZGXyEGZ3HUqnSJVf1LrQuNo28t4b/3gKkb+doDcSnTyJuEcBBYIyIrBKR0SKS0dqwdwE/Z3RBERkgIstFZPnBgwezErsxpoD5bsGs83cjlkQ8l+uJRET8gZ7ApHSb+nJhcjnHF2gCjFTVxsApYFC68z4PpABfZXQCVR2lqtGqGh0aGpqFGhhjCpJUpzJr/WsAPNppqJejyZ+8cUfSFVipqvHnCkTEF+jNPzvj04oFYlX13Jwh3+FKLOeOvwPoAdx6qaavf1MYV4vMK+x7b7xl4txprAw+Qhe7G7ls3kgkGd15dAI2q2psBvujqvuBPe6RXwAdgY0AInIN8AzQU1UzHjqUCQEBARw+fNg+0LxAVTl8+DABAQHeDsUUMqlO5actr+KvyhPXDPd2OPlWrj5HIiJBuEbT3Zdu0wV9JiJSERitqt3cRQ8DX7mbxrYDd7rLPwSKALPdo5j+VNWBeKhy5crExsZi/SfeERAQQOXKlb0dhilkxvz4BmuCznCrXyPKlYv0djj5lldGbXlbRqO2jDH5V0pyEn+sGEGVclGEh7XB4fPvfyMnnTlF7y+akyIw/dZFBAQWz4VI87eLjdqyJ9uNMfneu1Nv5ovE7bAFijuVho6itC7XlI6N76NcuYao08nuPX/w154FlCpWifBKLfl89uvs8RceDultSSSLLJEYY/K1mfP/yxeJ27nBvwKNQhux6sAalp/ex8L4+fxv5nxqpjqIl1SOO9I8wLva9aV+ki93dx/slbgLEkskxph8a+u2WQzeNoko8ef53lOITxRGjFlGzIGTVC2yjrDi8zkdFEdDRymurNKEJtVasyVuO7+sX8JpPUjv6Kfx8S18s/VmN0skxph86eSJfTw+7/8IAt7u9gVbDqVw59hlJCWn8vU9zalSqj1bD97F4m2HGbtoJ3NjhR6nKjBtTRAlgxoz/ObGtKxe2tvVKBAskRhj8qUv5zzJTh/4vNETHEipQp9PFlMyyJ+v7mlOrXLFAKhSKoj2tcvynxZVef2nTUxaEUv72qEMu6kRpYsW8XINCg5LJMaYfCcx8RBfH11DW59i1K17Oz3eX0CJQD++f+BKyhW/8HmkKqWCGHlbU+ISTlO+eAAOx6UnPDWesURijMl3pswbTIJDuDvqQV6atoHdRxKZcG+LDJNIWhVDAnMpwsLF1mw3xuQrycmJjNs3n8bqz25neyavjOWhDjVpHmH9Hd5iicQYk2dk5gHpmQtfZ5+P0LV8b56fsp7oqiV5pEONXIjOXIw1bRlj8oSRc7fx7q9/0SQshCurl6FtrVAaVQn5xz7O1BQ+2zGdcCc8NzeSiiF+DL8lCl8f+5vYmyyRGGO8bvfhRN799S9qlSvKiaQU3v31L96Z/RdXhJfkoQ41aVOzDBv3HWfKnEFs81HqH2jEXVfV5InOtQguYh9j3mb/A8YYr1JVXv5hAzUCVtKm3ArqlKlD1c5tWXkonFELY7nj86UUD/ClsmMO+8MWUPesL0/f8hZNwit4O3TjZonEGONVszfGs33HbFKqfs2Xpx2wJwb2/EBJp3J3/SsILPUoa7atZlHqj1R0Ovi0z1RKhFgSyUsskRhjvOb02VQ+nDEFCfsKf4Qf248kOfk0G3b9zk975vDe0eWUOnwbASqICCM6jqRESLi3wzbpWCIxxnjN6NmzSSr9EWcdMLb1MMLCrgKgevXO9ARWr/uKD1cOZ63zNJ80HXR+u8lbPE4kIhIMJKlqqofH1eafS+lGAIOBlsC5lQ9DgARVjcrg+BBgNBAJKHCXqi4WkVLu84YDO4E+qnrUk9iMMbnvbHIKC/c8y9Ei8FmzwdSscc0F+0Q1uJXRDW7l7JkT+Bcp5oUoTWb865g5EXGISD8R+VFEDgCbgX0iskFE3hKRmpm5kKpuUdUod5JoCiQCU1T15jTlk4HvL3KK94CZqloHaARscpcPAuaoak1gjvu9MSaPGz39GTYGptAvuDkN6998yX0tieRtmRl8/TtQHXgWKK+qVVS1LNAa+BMYKiK3eXjdjsA2Vd11rkBc6+T24cL13BGR4kAb4DMAVT2rqgnuzdcB49yvxwHXexiLMSaXHUvYyYRjM6lxRnjw+hHeDsdkUWaatjqpanL6QlU9gusOYrKI+Hl43QvWaMeVmOJVNSaD/SOAg8AYEWkErAAeVdVTQDlV3eeOaZ+IlM3ogiIyABgAEBYW5mG4xpjs9L8f7uKkQ7ip7GP4+9ssvPndv96RZJRELmefc0TEH+gJTEq3qS8Z3I24+QJNgJGq2hg4hYdNWKo6SlWjVTU6NDTUk0ONMdlo/pJ3+dF5kIYJ5bi1063eDsdkg0x3tovIExkUHwNWqOpqD67ZFVipqvFpzu0L9MbVd5KRWCBWVZe433/H34kkXkQquO9GKgAHPIjFGJNLjh3bzfCf7+W7M3FUPquUKTvY1gQpIDyZoCYaGAhUcv8bALQDPhWRpz04T0Z3Hp2Azaoam9EBqrof2OMe+QWuPpaN7tfTgTvcr+8ApnkQizEmi3bsnMuWv2Zccp9ZC17lusnd+D5pL921Mru2P0e/Vg1zKUKT0zwZ/lsaaKKqJwFE5CVcdwZtcPVZvPlvJxCRIKAzcF+6TRf0mYhIRWC0qnZzFz0MfOVuGtsO3OkuHwp8KyJ3A7uBmzyokzEmC5asHMXDa97ntEO4avGr3Nf0EaIi+53ffuzYbl6fcTs/pRymnvjyRtSL3De9GPWrFqVJWIj3AjfZSjIzbTOAiGwCGqnqWff7IsBqVa0rIqvcfRf5QnR0tC5fvtzbYRiTr/226E2e+usLwpwOri7dkK+PrOaoQ4hIFSr7BlPOvwTzEmM54oABIY24t/tnPPrtRmZvjOenR6+iRlkb0pvfiMgKVY1OX+7JHcnXwJ8icq7p6FpggvsBxY0XP8wYU9B8+dNzvH1gOrVSfXmkxVh8AqpQIj6WTX+9wT7dxoGURNamnqCc+PB+88HUr3sDM9fv48d1+3jq6tqWRAqYTN+RAIhIU+AqQICFqpov/6y3OxJjLt+OvX/Rd1YvKpzx4a/dz3DKGXJ+m0PAqeDv4+DqyPJcU788FUICKB7gyy2jllCueBGmPtgKP1s/JF/K8h2J+4HBukAJVR0iImEi0kxVl2ZnoMaYvO2DHx8n0V+4u8ELhPXowuFTZ/FzOKhSKpCKIYFsP3iKCUt38/3KWH5YE3f+OF+HMO6uKyyJFECeNG2NAJxAB2AIcALXA4lX5EBcxpg8aMH6FSz020mL5BL0aNUnw31qly/Gyz3rM6hrHbYdPEn88ST2HztDxZAA6lcskcsRm9zgSSJprqpNRGQVgKoedY+gMsYUAqlOZfy8ZzgTLDzS9uV/3T/Az4f6FUtY8igEPLnHTBYRH1wz7yIiobjuUIwxhcAXv//GiqD9tNNSRNbq4u1wTB7iSSJ5H5gClBOR14CFwOs5EpUxJk+JP57EvA2vkCLweIdXvR2OyWMy3bSlql+JyApcT5UDXK+qmy51jDEmf0g8m0Kgnw+uMTUXGjrpK9YUP0RX33KEV22Ty9GZvO5fE8lF5tgC6CoiXVX1nWyOyRiTS46dTubjedv4fOEOGoeF8GG/JpRJN//Vz+v2sS9pBP5ByhNd3vdSpCYvy0zTVjH3v2jgfv6ea2sgUC/nQjPG5KQv5v1Jt7e+Y+TcGK6qUYZVuxPo8f5CVu7+e4HRY4nJfDlzGBuKnuWe0tGElq3vxYhNXuXJFCmzgBtU9YT7fTFgkqpeuD5mHmcPJJrCasnKUcyJmcaSU3vY7uP63fdVpaQT2hapyoK9d7PrWACtapQhyN+HuKPHOVPkYdRHmNZvIUUCbARWYZYdU6SEAWfTvD+La510Y0w+sG7DJO5Z9wGBTqVWShGaJpTnqhoVOJmcwJ7T8UxO3kVomcHcVLoTW0504HTKHkJ9f+RPf2F4zdssiZiL8iSRfAksFZEpuIYA9+LvJW6NMXmYOp0MX/YGpVT5rMsPXPPRJu64Mpx7evzdOr1uwySGLH2VH3zmQIk558tbEEiHFk95I2yTT3gyaus1EfkZ15K4AHeq6qqcCcsYk50Wr/yYpXKGQeVa8+2aJBS4s1X4P/ZpUP8mJtS6lh/mvUhC0hFKB5WjVNEKRDe4DXHYtCbm4jIzakvU3ZGiqiuBlZfaxxiTtzhTUxi+bhSVFK6+8jXavbOM7g0qULlk0AX7+voF0KvTW16I0uRnmfkz43cReVhEwtIWioi/iHQQkXH8vUKhMSaPmfXH62xypPJgtZ58v+YoJ8+kcG/rCG+HZQqQzDRtXQPchWvtkWpAAhCIKwnNAt7NzJrt7mVyJ6YpigAGAy2Bc0vohgAJqhqVwfE7cU0UmQqknBs5ICJRwMdAAJACPGAzEhvjkpycyIdbJ1FDHHRs+QId311Ey4jSNKhsHecm+/xrIlHVJFwz/44QET+gDHBaVRM8uZCqbgGiANxzdu0Fpqjq8HP7iMjbwLFLnKa9qh5KV/Ym8F9V/VlEurnft/MkNmMKqm9mP84uH/ioTn8mLN/PvmNJvHmjrZVuspdHPWiqmqyq+zxNIhnoCGxT1V3nCtzrnfQh3drtmQkLKO5+XQKIu8S+xhQahw/9xYj4P2hFEA0jH+CD32JoXbMMrWuGejs0U8B4Mvw3O93ChQmjNRCvqjEXOUaBWSKiwCeqOspd/hjwi4gMw5UYr8zoYBEZAAwACAsLy2gXYwqUD2Y/RJLAM23fYOT87Rw7ncygrnW8HZYpgHJ9TJ97DZOewKR0m/py6buRVqraBOgKPCgi52aOux94XFWrAI8Dn2V0sKqOUtVoVY0ODbW/yEz+dOx0Mp/M28aBE0mX3G/D5u/5/kwctwbXwD+kOWP+2EmvxpVsbRCTIzxZavduVc3wQ9pDXYGVqhqf5ty+QG+g6cUOUtU499cD7ocimwHzcY0Ye9S92yRgdDbEaEyek+pUnvhqNsknhvHbxgQSipzliED3gIoMbPcG5Ss0BiA+fi2vL36VUgoDrxnJ4J+3APB/XWpf6vTGXDZPmrbeFpFbcY2MWgpMUNUNl3HNjO48OgGbVTU2owNEJBhwqOoJ9+suuJb7BVefSFtgLq5lgC/WNGZMvvb2rC2cOf0K60qeIDxZKJMYRLi/P9Mkjh9++Q/di5Rnx5kjrJZkcMDTod157PtYft10gPvaRlApJNDbVTAFlCeJ5DDwKuCPa/TVtyLyvqp+ktkTiEgQ0Bm4L92mC/pMRKQiMFpVuwHlgCnutRJ8ga9VdaZ713uB99x3NUm4+0GMKUh+WrePRcs/ZFuVk9xbrC4De07kg99iGDF3G+V8Y4gI/Ybp7KdqstDRWZ1UR1cGL6hMsSJHeOrq2vbciMlRnsz+u9LdR3HufSCwRFXz3VhCm/3X5Cc7Dp3ipg+mUSLsNYqKD9/2W3B+AsWjp86yOjaBNXsS2BB7iKOn4URSCqeTU+nRsAL3talOiSA/L9fAFBTZMfvvuRM9g+tZkhK4HhA0xuSgj+duo17pEaz2Ed6Jfu4fs/CWDPanfe2ytK9dFqjlvSBNoXY5o7YmA1uBytia7cbkqEMnz7B583hWlTjGf4Kr0yjyFm+HZMwFPLkjKSkiVVR1K7BVRD4FVgE/5kxoxpjxi7ZD6M+UTVUe7J4dgyaNyX6eJJLiwFwROQRsxDUvVmpOBGWMgaTkVFatfYvtofBK5e4EBZXxdkjGZMiTRNIeWA80xzXJomJ3I8ZkmjqdHq3rMXV5DHtDllAjxYdr276Sg5EZkzWZ/qlW1bWq6lTVxao6VlXHZTCBojEmA8MmXU/nsQ3ZvGV6pvZXVeYtf554PwdPNrwPH1//HI7QmMtny54Zk8OmznmGcYnbOOqAuxY9x6q14//1mEmLl7AmeDPRKQG0uuLBXIjSmMtnicSYHLRh02Re2fMjzQlgWqfRlFZhwIqhTJ87jKTkjLsYv5k/hy/W38tph/Bkq5dyOWJjPOfJA4kPAV+p6tGcDSnn2QOJJjccObKVm6dejwATr59GyVLVOXRoC3dPvYntfkqDEwHomb6UrtCO6KolaVatNIvW/MiEva9ywgfeqHMf7Vs+4u1qGHNedjyQWB5YJiIrgc+BX2yddmMubvgvAznigC9bvErJUtUB+H6jL9u2P03XiG+YH7yD00U/p96x8fyy2sEvq5VdgUn4OeDTZq/QqH5vL9fAmMzxpLP9BaAmrmna+wMxIvK6iFTPodiMybeOH9vDz0n76RlQiXp1rgdg0bZDDP15M23q1+XNu39g5nXT6V+sFqeLpJAYeJbTAcnU0kDGdhxtScTkKx5NkaKqKiL7gf24ZgEuCXwnIrNV9emcCNCY/OjHxW+Q5BBujHLNT7rv2Gke/noVEaFFefPGRogIJUtV54kbvucJL8dqTFZ5sh7JI7jW/jiEa82Pp1Q1WUQcuKZut0RiDK7nRSbtW0A9fKhfpzdOp/LExDUkJafy8W1NKVrEWwuTGpMzMvUT7V5PvRHQO+066wCq6hSRHjkRnDH50bpNk4hxOBlcoRMAXy/dzeLth/lf7wbUKFvUy9EZk/0y1Ufi7lRvnD6JpNm+KVujMiYfm7TmUwKdSrcrn2HPkUT+99MmWtcswy1XVPF2aMbkCE/usReLyBWquuxyLiQitYGJaYoigMFAS1xTroBr/q4EVY3K4PiduKatTwVS0g5BE5GHgYdw9dv8aP01xltOHN/LzDP76R5QkaDgcgz4ZgkiwtAbGuJemM2YAsfTubYGuj/QTwGC62YlUwtbqeoWXCsrIiI+wF5giqoOP7ePiLwNHLtUDOmnZRGR9sB1QENVPSMiZTNbIWOy2w+LXifJIfRudB9vz/qLP7Ye5vVeDWyZW1OgeZJIumbjdTsC29I2lbn7YfrgWnfdE/cDQ1X1DICqHsi2KI3xwM6d8/hg3zwaqT/PzynHur1buT6qIn2bWZOWKdg8SSR3XKR8yGVc94I12oHWQLyqxlzkGAVmiYgCn6jqKHd5LaC1iLyGa832JzNqfhORAbjXcw8LC7uMkE1hdiIpmTdnbuF4UjJ3X1WNhpVD/rE9MfEQj//2MD4oO3bexRn/03zUrwndGpS3Ji1T4HmSSE6leR0A9AA87mQXEX+gJ/Bsuk19uTC5pNVKVePcTVezRWSzqs7HVYeSQAvgCuBbEYlI/9S9O/GMAtcUKZ7GbQqvtbEJPDxhFXuOnKSYP0xbHUerGqW5rXlVWkSUJiTQlxe/u4FtDifVYjsRXq8VL19bn5LBNmOvKRwynUhU9e2070VkGJC5ObH/qSuwUlXj05zLF+gNNL3E9ePcXw+IyBSgGTAfiAW+dyeOpSLixLWm/MHLiM2Yf5iwdDeDp62nXvFtlKr5KTsdThqnOPA9E8yEXwP52ucMqf5JrAs8S/TRqtzQ8xl6NKzo7bCNyVVZeTIqCNfIK09ldOfRCdisqrEZHSAiwYBDVU+4X3fh7ya1qbj6VeaKSC3AH9dDk8Zkye7Dibw0bQPdKi9gecAPCNA3uDrbE/ezxXGCE0EnCFYlyCm0SSrL8/0nUrGkPSdiCh9Pnmxfh6ufAsAHCMXD/hERCQI6A/el23RBn4mIVARGq2o3oBwwxd3W7At8raoz3bt+DnwuIuuBs8AdNpmkyQ5DZ26iZalRzA2Mobr68H7nT6hcuYW3wzImz/HkjiTt0+spuDrGUzy5mKomAqUzKO+fQVkc0M39ejuuJ+szOudZ4DZP4jDm3yzZtp+TB59gVZkjdHCU4H+9phBU1EaWG5MRT/pIMnyq3ZiCJvHUET6a3ZNVJU9za2A1nr7hexw+Nj+WMReT6WnkRWSciISkeV9SRD7PkaiM8ZINm7/nP9+0Z3VAIrf5NWZQn+mWRIz5F578hjRU1YRzb1T1qIg0zv6QjMl9x4/t4cOZA5l4ehclUFqd6MBTD77v7bCMyRc8SSQOESl5bqldESnl4fHG5EnHj+2h7+RuxDqU1mfLsmD3nTx4f08cDnuQ0JjM8CQRvA0sEpHvcI3e6gO8liNRGZNL1Olk8PRbiHMoT5a/mxfn1uKRjjWJrFTC26EZk294sh7J78ByXM9sCK61STbmYGzGXDZ1Otn411R+WDuGeSd30jo4jCev/wb/IsX+sd/XvzzIHOdxHinVnPdXNKBuhSI81L6Gl6I2Jn/KVCJxL7E7VVWbApY8TJ6WmHiIeyd2Zq0jBT9VGjoCmJC0m3Vft+Htqz+lYkXXCgTrNkxiWPwC2kgx/jwygKOnDjP2zivw9830GBRjDJ41bf2ZlfVIjMktX/36BGsdKTwV2pLrWr1AiRJhzPljKC/8NZ4+v/SnniOIuNTT7HUoZVJhbdwD7E08xPPd61G/ojVpGeOpXFuPxJjccCxhJ2MOr6SdT3Fu7zbqfHmrZk/x0JlazNj6CgdJolRKIGVTirH10DVUrFSDT66tT72Kxb0YuTH5l7fWIzEmR3w25wlOCjxy5WAAklOdfPDbVsb/uYsjp4oQXvod6lcqgb+PAz8f4ab2Zbkm0qZ6NyYrPEkku4FbgQhVHSIiYUB5wJ54N3lCfPxavj7+Fz38y1KzxjUcT0rmwa9WsiDmEJ3rleOOluFcWb20Des1Jpt5kkhGAE5co7aG4Fo/fTKuNUCM8bqPf3+KVIEH2r7OniOJ3DV2GTsOneLNGxvSJ9pWKTQmp3iSSJqrahMRWQXnn2y3lXtMnrAl5kemJO3lpoCq/LgtlI/nLSTVqXxxdzOurF7G2+EZU6B5kkiSRcQH91TyIhKK6w7FGK9KTTnLywtfoLgqv8Tcwu7Vm2ldswwv96xP9VBbH8SYnOZJInkfmAKUc6+PfiPwYo5EZYwHxv/yEOsdKdSNu4LASrUZ3qkmTcJKejssYwoNT6aR/0pEVgAd3UXXqermnAnLmMzZs2cxHx5YRL3T/lQIf4R3b2lsI7CMyWX/mkhEJP267Od+S68WEVS1Z2YuJCK1gYlpiiKAwUBLoLa7LARIUNWoDI7fiauDPxVIUdXodNufBN4CQlXVltotBFJTzvLyr4/gQPFNeZShNzayJGKMF2TmjqQlsAfXUrhL+DuReERVtwBRAO6+lr3AFFUdfm4fEXkbOHaJ07TPKEmISBVcS/juvpzYTN52JiWVZTuOEhEaTMWQQAD2H9zKc9P7ssw3iSsOR/LqXX0I8PPxcqTGFE6ZSSTlcX1I9wX6AT8CE1R1Qxau2xHYlnbVRffEkH1wDS/21LvA08C0LMRk8pgzKal8uzyWEb9vZd+xJADCSgURXWola/ico75CsyMNeeCGD88nGGNM7vvXRKKqqcBMYKaIFMGVUOaKyBBV/eAyr3sLrjuctFrjWgc+5mKhALNERIFPVHUUgIj0BPaq6ppLNWuIyABgAEBYWNhlhm1yy8a449z7xXL2JpymadWSvNijHju2TWV1/Hjm+h6jZKrwQtX/47rb++NjDxga41WZnUa+CNAdVxIJxzWC6/vLuaD72ZOewLPpNvXlwuSSVitVjRORssBsEdmMa1r754Eu/3Zdd+IZBRAdHa2XE7vJHfuOneausa65QcfeEcWZQ58zbvn3rHOkEOSv3BAYxv2d36dUKZvu3Zi8IDOd7eOASOBn4L+quj6L1+wKrFTV+DTX8AV6A00vdpCqxrm/HhCRKUAz4ChQDTh3N1IZWCkizVR1fxbjNLng7JkTLFkzhvo1ulGqVA1OJCVz55hlOJJjuafxCt5e/Bg7fSBM4dlybejZ6gWKFqvg7bCNMWlk5o7kP7hm+60FPJKm+ejc7L+eTpma0Z1HJ2CzqsZmdICIBAMOVT3hft0FGKKq64CyafbbCUTbqK38QZ1OXpjUg59TjyAbR1FXfSlxphi+QQmcKAnvHoW64sPbETfSseXT+PjaRArG5EWZ6SPJtlV+RCQIV8f9fek2XdBnIiIVgdGq2g0oB0xxJzFf4GtVnZldcRnv+OLn+/g59Qi3B0VQ3L84vx9cz9oiR6jlDKBnyXo0j+hOw3p9EIctNGVMXubJk+1ZpqqJQOkMyvtnUBYHdHO/3g40ysT5w7McpMkVi5aP4J2Di+nsE8KTN0zh978O8r/Fy7m1eRiv9Wrg7fCMMR6wP/VMrtu1awFPrxtBhNPBq72+Z/+JM/zft2uoW6E4L/ao5+3wjDEeskRictXGzVO5/bf7EeC9TiPwDyzDIxNWcTbFyUf9GttDhcbkQ5ZITK75c8Un3Ln4BYoojGv3PiFlm3PflytYtvMor/duQITN1GtMvpSrfSSmcDl0aDN/rB3LjqNb2Zm4j3mpxwjHh4+7f0miT3V6ffQHuw4n8sp19bkuqpK3wzXGXCZLJCbbxWydyZfLhjHjzH6SRfBVpbLTQbeAsjzdfSxrDwTy4FcL8fVx8OXdzWlZ/YLxF8aYfMQSick2Tqfy7Nhr+clnFwFOpXdQGDc1eZCI8Pb4+QUBMGHpbl6YuoyaZYvy6e3RVCkV5OWojTFZZYnEZItTZ1J4cUxfZgfuovHxYmw9ej8tbuxC7ZrlAEhJdTJs1l98PG8bbWqF8lG/xhQL8PNy1MaY7GCJxGRZ7NFEXv3yIf4osZmrUooxqO9MHpywhnu+WE63BuXZfyyJjfuOk5TspF/zMIb0rI+vj43zMKagsERisiTVqbwy/r8sLb6UJqkBvHfbL/gXKcZ3A69kyIyN/LRuHzXLFqVfs6q0iChF53rlbPEpYwoYSyQmS8b8NpdNQTMIS/Xhoz4/4F+kGAABfj683qsBr9tT6sYUeJZIzGWLPXKKOX89yalAYXTrN21WXmMKKWuoNpdFVXn32/tZH3yWASVbUKvGNd4OyRjjJZZIzGX5du4U5hdZSeNkf+699mNvh2OM8SJr2jIem7HwGz7Z/gp+AkO7fYbDx36MjCnM7BPAeGTCL6/wbtw3BAPvNv0vFStGeTskY4yX5VoiEZHawMQ0RRHAYKAlUNtdFgIkqGpUBsfvBE4AqUCKqka7y98CrgXOAtuAO1U1ISfqUJidSEpm1LR7GHd6BRHJPgztPI461Rp7OyxjTB6Qa4lEVbcAUQAi4gPsBaao6vBz+4jI28CxS5ymfQbL6M4GnlXVFBF5A3gWeCYbQy+0kpJT+ej3rcxYu5cqPq+yqtRBohKDeLnXVKpXqOjt8IwxeYS3mrY6AttUdde5AnE9pdYH6ODJiVR1Vpq3fwI3ZkuEhdyirYd4bso6dh8+Tufwt1kUeJTuUp5Bd0wlJDjY2+EZY/IQbyWSC9ZoB1oD8aoac5FjFJglIgp8oqqjMtjnLv7ZfHaeiAwABgCEhYVdVtAFzZxFb7B4zzxikw6zN/U04b5Fuafpk0z8qyZfL9nNlaGrqFprKot8znBvsTo8fP1EWz/dGHMBUdXcvaCIPxAH1FfV+DTlI4Gtqvr2RY6rqKpxIlIWV3PWw6o6P83254FooLf+S6Wio6N1+fLl2VCb/Ck15SzDp97M2FNbKeZUKuNLBd9glqcc47hDiDxZhGB/H5b4J1LCqTxSqRN9ugz3dtjGGC8TkRXn+qfT8sYdSVdgZbok4gv0Bppe7CBVjXN/PSAiU4BmwHz38XcAPYCO/5ZECrtTJ/fzzPe9mKcn6RdYlad6fYfDpwjfr9rLnz8s4IqQsWwosRcBHijRkP90fNueWDfGXJI3EklfLmzW6gRsVtXYjA4QkWDAoaon3K+7AEPc267B1bneVlUTcy7s/C/h6A7umXo9WyWVFyp0oE+X95iz6QDDZm1h8/4TNK1alVf7TiXE9wQAQUXLejliY0x+kKuJRESCgM7Afek2XdBnIiIVgdGq2g0oB0xxzxrrC3ytqjPdu34IFAFmu7f/qaoDc6wS+dSxY7sZMLUXOySVD+vdR/GK/+HGjxezYtdRwksH8d4tUVzbsCIOhwCB3g7XGJOP5HofSV5Q2PpIThzfy4DJPdgiybxe/S5+i+/IpBWxlA4uwv91qcWNTSvjZ+uDGGP+RV7qIzG5aO2Giby69H/ESAovhd3Ki/MbkpC4lwGtI3ioQw1bpdAYk2WWSAqoXbsW8N7855jtTKAUyusRt/H28pacSU5i2oNXUa9icW+HaIwpICyRFACnz6byw9o4pq3eyy3RlTm1/w3e2v87Pgr3hzTglvZvM3DiLnYfTuDLu5tZEjHGZCtLJPmY06m8PXsL4//czbHTyZQLOMKEMw+xrugZrpKiDOn2GcVK1uGJb1ezdMcR3u/bmOYRpb0dtjGmgLEe1nxszKKdfPT7NlpElGJ41wOEhL3BpuAkog9H8PjVMzntU40bRi7ip3X7eaF7XXo2svmxjDHZz+5I8qmY+BO8MXMzHeuUoWO58QzZ8StlgOH1BvF/s8rRf8wKTp5JAeCzO6LpWLecdwM2xhRYlkjyoeRUJ49/u5pSRZII9X2WV/YdpJUEM7TXN4SUrMboCgncMmoxNcoWZeStTalSKsjbIRtjCjBLJPnQB79t5cjBlVSs+jk/pjoZWDySgdeOw8fXH4CoKiH88UwHigf62fMhxpgcZ4kkn/lj6yHmL/0ICf+V/SJ8WPsu2rZ44oL9Shct4oXojDGFkSWSfGT7gRN89kN/dlTaSYTTh/c6jiAs7Cpvh2WMKeQskeQT+w7G8tLkXqwulUQnKclrN062SRWNMXmCJZJ8IG7/Ou6d0Y+9gcrdQU159IYxtsCUMSbPsE+jPO5Iwh4GzujHQR8n95e8ncduGmdJxBiTp9gnUh52KvEoD0zqyR5f5dai/bjvuqe9HZIxxlzAEkkedfpMIg+Nv5oN/in08+/Eozc97+2QjDEmQ5ZIvESdTlasGcfT49vS+fNIvv/1yfPbjhw/zMCx7Vhe5DR9HI15qt9w7wVqjDH/Itc620WkNjAxTVEEMBhoCdR2l4UACaoalcHxO4ETQCqQcm5xFREp5T5vOLAT6KOqR3OgCtlm4+apvLD4JWIcToo5lcriy0t7f2HDN+u5rfX7PDvjZjYEpHCrX1MG9Rvr7XCNMeaSvLJCooj4AHuB5qq6K03528AxVR2SwTE7gWhVPZSu/E3giKoOFZFBQElVfeZS1/fmComr13/NA8teJ1jhgbCuXHPls/j5BfH+9H6MORlDEaeSInBf8c7c3/tdr8RojDEZyWsrJHYEtqVLIgL0ATp4eK7rgHbu1+OAucAlE4m3LF01modWDydUhdHXjKVCxabnt/XpMJ71Yx5jX7FF/KdKH/pdM9iLkRpjTOZ5K5HcAkxIV9YaiFfVmIsco8AsEVHgE1Ud5S4vp6r7AFR1n4hk+JSeiAwABgCEhYVlNX6PqNPJ1N8H8drun6isDj7tMYHQsvXPb49LOE3fT/8kIfEGvu7zFg0ql8jV+IwxJityvWlLRPyBOKC+qsanKR8JbFXVty9yXEVVjXMnitnAw6o6X0QSVDUkzX5HVbXkpWLIzaat/ftX899fBrKQU9ROciBJzxJYLIKSQX74iIDA0h1HOJaYzPh7mtOoSkiuxGWMMZ7KS01bXYGV6ZKIL9AbaHqxg1Q1zv31gIhMAZoB84F4EangvhupABzI0eg9MPfPt3lu0xiSgUaHanEy8FFKFg/kyKmz7Dx0CqcqqlAswJcP+ja2JGKMyZe8kUj6cmGzVidgs6rGZnSAiAQDDlU94X7dBTjXIT8duAMY6v46LUei9tD0355j8O7phKc4OLinP41adeOJzrXxcYi3QzPGmGyVq4lERIKAzsB96TZd0GciIhWB0araDSgHTHH1x+MLfK2qM927DgW+FZG7gd3ATTlXg8z58qf7ePPgIuom+bJt3yBevaktXRtU8HZYxhiTI3I1kahqIlA6g/L+GZTFAd3cr7cDjS5yzsO4RoHlCZ9Ov533j66iwaki7D76IuMGtCOyknWeG2MKLpv99zKpKruPJLI34TT7EpKISzhN7I4X+MlvI41OBHE4+TUmPdSKCiUCvR2qMcbkKEskl+HIqbM8MmEVC7f+/Wxkq9KfsbZsDE1OB1G5wkeM7N6IYgF+XozSGGNyhyUSD63fe4z7vlzBwZNnGNS1Dg0rlWDLpv/yzpEY2ktx3r57Nn5+Qd4O0xhjco0lEg9MWRXLoMnrKB3sz3cDW9KwcgjfzHyYd47MpY0UZdjNMy2JGGMKHZv91wO7dv9E24pz+fjmokRWLMZXPz/Aa/FzaSfFePfmWfgXKebtEI0xJtfZHYkHjiR+x+LgfSxe8BMB85Qkh9DRUZy3+szCr0iwt8MzxhivsETigae7jOCG3fPYdmAtW49tp6hfUe7t8Zk1ZxljCjVLJB4oWroGjUrXyPiBFmOMKaSsj8QYY0yWWCIxxhiTJZZIjDHGZIklEmOMMVliicQYY0yWWCIxxhiTJZZIjDHGZIklEmOMMVkiqurtGHKdiBwEdl3m4WWAQ/+6V8FidS4crM6FQ1bqXFVVQ9MXFspEkhUislxVo70dR26yOhcOVufCISfqbE1bxhhjssQSiTHGmCyxROK5Ud4OwAuszoWD1blwyPY6Wx+JMcaYLLE7EmOMMVliicQYY0yWWCLxgIhcIyJbRGSriAzydjzZTUSqiMjvIrJJRDaIyKPu8lIiMltEYtxfS3o71uwmIj4iskpEZrjfF+g6i0iIiHwnIpvd/98tC0GdH3f/XK8XkQkiElDQ6iwin4vIARFZn6bsonUUkWfdn2dbROTqy72uJZJMEhEf4COgK1AP6Csi9bwbVbZLAf5PVesCLYAH3XUcBMxR1ZrAHPf7guZRYFOa9wW9zu8BM1W1DtAIV90LbJ1FpBLwCBCtqpGAD3ALBa/OY4Fr0pVlWEf37/YtQH33MSPcn3Mes0SSec2Araq6XVXPAt8A13k5pmylqvtUdaX79QlcHy6VcNVznHu3ccD1Xgkwh4hIZaA7MDpNcYGts4gUB9oAnwGo6llVTaAA19nNFwgUEV8gCIijgNVZVecDR9IVX6yO1wHfqOoZVd0BbMX1OecxSySZVwnYk+Z9rLusQBKRcKAxsAQop6r7wJVsgLJeDC0nDAeeBpxpygpynSOAg8AYd3PeaBEJpgDXWVX3AsOA3cA+4JiqzqIA1zmNi9Ux2z7TLJFknmRQViDHTotIUWAy8JiqHvd2PDlJRHoAB1R1hbdjyUW+QBNgpKo2Bk6R/5t0LsndL3AdUA2oCASLyG3ejcrrsu0zzRJJ5sUCVdK8r4zr1rhAERE/XEnkK1X93l0cLyIV3NsrAAe8FV8OaAX0FJGduJorO4jIeAp2nWOBWFVd4n7/Ha7EUpDr3AnYoaoHVTUZ+B64koJd53MuVsds+0yzRJJ5y4CaIlJNRPxxdVJN93JM2UpEBFe7+SZVfSfNpunAHe7XdwDTcju2nKKqz6pqZVUNx/V/+puq3kbBrvN+YI+I1HYXdQQ2UoDrjKtJq4WIBLl/zjvi6gMsyHU+52J1nA7cIiJFRKQaUBNYejkXsCfbPSAi3XC1p/sAn6vqa96NKHuJyFXAAmAdf/cXPIern+RbIAzXL+RNqpq+Qy/fE5F2wJOq2kNESlOA6ywiUbgGF/gD24E7cf1hWZDr/F/gZlyjE1cB9wBFKUB1FpEJQDtcU8XHAy8BU7lIHUXkeeAuXN+Tx1T158u6riUSY4wxWWFNW8YYY7LEEokxxpgssURijDEmSyyRGGOMyRJLJMYYY7LEEokxWSAipUVktfvffhHZ6359UkRG5NA1HxOR2/9ln29EpGZOXN+Y9Gz4rzHZREReBk6q6rAcvIYvsBJooqopl9ivLXCbqt6bU7EYc47dkRiTA0SkXZq1TV4WkXEiMktEdopIbxF5U0TWichM97Q0iEhTEZknIitE5Jdz01qk0wFYqaopIlJdRFamuWZNETk3Z9gCoJM78RiToyyRGJM7quOaqv46YDzwu6o2AE4D3d3J5APgRlVtCnwOZDRzQitgBYCqbgOOuZ9SB9fT6WPd25y4pgVvlEP1MeY8+2vFmNzxs6omi8g6XFPszHSXrwPCgdpAJDDbNRUUPrimO0+vAv9cgGs0cKeIPIFr+o+060kcwDXTbWGa2dh4gSUSY3LHGXDdKYhIsv7dOenE9XsowAZVbfkv5zkNBKR5PxnXfEq/AStU9XCabQHu/Y3JUda0ZUzesAUIFZGW4JrOX0TqZ7DfJqDGuTeqmgT8AowExqTbtxawIWfCNeZvlkiMyQPcyzffCLwhImuA1bjWy0jvZ1zL5Kb1Fa4FiWadKxCRcsDpcyvjGZOTbPivMfmMiEwBnlbVGPf7J4ESqvpimn0eB46r6mdeCtMUItZHYkz+MwhXp3uMO6lUxzUsOK0E4MtcjssUUnZHYowxJkusj8QYY0yWWCIxxhiTJZZIjDHGZIklEmOMMVliicQYY0yW/D8/PAcVLa7mvgAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] diff --git a/src/io/io.f90 b/src/io/io.f90 index dc479da5f..245c6a052 100644 --- a/src/io/io.f90 +++ b/src/io/io.f90 @@ -975,12 +975,8 @@ module subroutine io_write_discard(self, param) call util_exit(FAILURE) end select lfirst = .false. - if (param%lgr) then - select type(discards) - class is (whm_tp) - call discards%pv2vh(param) - end select - end if + if (param%lgr) call discards%pv2v(param) + write(LUN, HDRFMT) t, nsp, param%lbig_discard do i = 1, nsp write(LUN, NAMEFMT) sub, dname(i), dstatus(i) @@ -993,11 +989,8 @@ module subroutine io_write_discard(self, param) if (param%lgr) then allocate(pltemp, source = pl) - select type(pltemp) - class is (whm_pl) - call pltemp%pv2vh(param) - allocate(vh, source = pltemp%vh) - end select + call pltemp%pv2v(param) + allocate(vh, source = pltemp%vh) deallocate(pltemp) else allocate(vh, source = pl%vh) @@ -1007,7 +1000,6 @@ module subroutine io_write_discard(self, param) do i = 1, npl write(LUN, PLNAMEFMT) name(i), GMpl(i), Rpl(i) write(LUN, VECFMT) xh(1, i), xh(2, i), xh(3, i) - write(LUN, VECFMT) vh(1, i), vh(2, i), vh(3, i) end do deallocate(vh) @@ -1189,16 +1181,8 @@ module subroutine io_write_frame_system(self, iu, param) call io_write_hdr(iu, param%t, pl%nbody, tp%nbody, param%out_form, param%out_type) if (param%lgr) then - associate(vh => pl%vh, vht => tp%vh) - select type(pl) - class is (whm_pl) - call pl%pv2vh(param) - end select - select type(tp) - class is (whm_tp) - call tp%pv2vh(param) - end select - end associate + call pl%pv2v(param) + call tp%pv2v(param) end if if (param%out_form == EL) then ! Do an orbital element conversion prior to writing out the frame, as we have access to the central body here diff --git a/src/modules/swiftest_classes.f90 b/src/modules/swiftest_classes.f90 index 0b5fdf4c8..b0a804610 100644 --- a/src/modules/swiftest_classes.f90 +++ b/src/modules/swiftest_classes.f90 @@ -159,8 +159,8 @@ module swiftest_classes procedure(abstract_step_body), public, deferred :: step procedure(abstract_accel), public, deferred :: accel ! These are concrete because the implementation is the same for all types of particles - procedure, public :: vh2pv => gr_vh2pv_body !! Converts from heliocentric velocity to psudeovelocity for GR calculations - procedure, public :: pv2vh => gr_pv2vh_body !! Converts from psudeovelocity to heliocentric velocity for GR calculations + procedure, public :: v2pv => gr_vh2pv_body !! Converts from velocity to psudeovelocity for GR calculations using symplectic integrators + procedure, public :: pv2v => gr_pv2vh_body !! Converts from psudeovelocity to velocity for GR calculations using symplectic integrators procedure, public :: initialize => io_read_body_in !! Read in body initial conditions from a file procedure, public :: read_frame => io_read_frame_body !! I/O routine for writing out a single frame of time-series data for the central body procedure, public :: write_frame => io_write_frame_body !! I/O routine for writing out a single frame of time-series data for the central body diff --git a/src/modules/whm_classes.f90 b/src/modules/whm_classes.f90 index 00d73a471..6107a719d 100644 --- a/src/modules/whm_classes.f90 +++ b/src/modules/whm_classes.f90 @@ -162,7 +162,7 @@ end subroutine whm_gr_getacch_tp module pure subroutine whm_gr_p4_pl(self, param, dt) use swiftest_classes, only : swiftest_parameters implicit none - class(whm_pl), intent(inout) :: self !! Swiftest particle object + class(whm_pl), intent(inout) :: self !! WHM massive body object class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters real(DP), intent(in) :: dt !! Step size end subroutine whm_gr_p4_pl @@ -178,7 +178,7 @@ end subroutine whm_gr_p4_tp !> Reads WHM massive body object in from file module subroutine whm_setup_pl(self,n) implicit none - class(whm_pl), intent(inout) :: self !! Swiftest test particle object + class(whm_pl), intent(inout) :: self !! WHM massive body objectobject integer(I4B), intent(in) :: n !! Number of test particles to allocate end subroutine whm_setup_pl diff --git a/src/whm/whm_setup.f90 b/src/whm/whm_setup.f90 index b0b25dc55..f9a28478f 100644 --- a/src/whm/whm_setup.f90 +++ b/src/whm/whm_setup.f90 @@ -58,15 +58,14 @@ module subroutine whm_util_set_mu_eta_pl(self, cb) ! Internals integer(I4B) :: i - associate(pl => self, npl => self%nbody, GMpl => self%Gmass, muj => self%muj, & - eta => self%eta, GMcb => cb%Gmass) + associate(pl => self, npl => self%nbody) if (npl == 0) return call util_set_mu_pl(pl, cb) - eta(1) = GMcb + GMpl(1) - muj(1) = eta(1) + pl%eta(1) = cb%Gmass + pl%Gmass(1) + pl%muj(1) = pl%eta(1) do i = 2, npl - eta(i) = eta(i - 1) + GMpl(i) - muj(i) = GMcb * eta(i) / eta(i - 1) + pl%eta(i) = pl%eta(i - 1) + pl%Gmass(i) + pl%muj(i) = cb%Gmass * pl%eta(i) / pl%eta(i - 1) end do end associate @@ -81,20 +80,15 @@ module subroutine whm_setup_system(self, param) ! Arguments class(whm_nbody_system), intent(inout) :: self !! Swiftest system object class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters of on parameters + call io_read_initialize_system(self, param) ! Make sure that the discard list gets allocated initially call self%tp_discards%setup(self%tp%nbody) call self%pl%set_mu(self%cb) call self%tp%set_mu(self%cb) if (param%lgr) then - select type(pl => self%pl) - class is (whm_pl) - call pl%vh2pv(param) - end select - select type(tp => self%tp) - class is (whm_tp) - call tp%vh2pv(param) - end select + call self%pl%v2pv(param) + call self%tp%v2pv(param) end if end subroutine whm_setup_system diff --git a/src/whm/whm_step.f90 b/src/whm/whm_step.f90 index aac704f67..55e7611b1 100644 --- a/src/whm/whm_step.f90 +++ b/src/whm/whm_step.f90 @@ -16,7 +16,7 @@ module subroutine whm_step_system(self, param, t, dt) real(DP), intent(in) :: t !! Current simulation time real(DP), intent(in) :: dt !! Current stepsize - associate(system => self, cb => self%cb, pl => self%pl, tp => self%tp, ntp => self%tp%nbody) + associate(system => self, cb => self%cb, pl => self%pl, tp => self%tp) call pl%set_rhill(cb) call pl%step(system, param, t, dt) call tp%step(system, param, t, dt)