diff --git a/examples/whm_gr_test/init_cond.py b/examples/whm_gr_test/init_cond.py index c2580cddb..de15165c1 100755 --- a/examples/whm_gr_test/init_cond.py +++ b/examples/whm_gr_test/init_cond.py @@ -8,7 +8,7 @@ sim.param['DU2M'] = swiftest.AU2M sim.param['T0'] = 0.0 sim.param['DT'] = 0.25 * swiftest.JD2S / swiftest.YR2S -sim.param['TSTOP'] = 1000.0 +sim.param['TSTOP'] = 100.0 sim.param['ISTEP_OUT'] = 1461 sim.param['ISTEP_DUMP'] = 1461 sim.param['CHK_QMIN_COORD'] = "HELIO" diff --git a/examples/whm_gr_test/param.swifter.in b/examples/whm_gr_test/param.swifter.in index 789250f41..6addd694c 100644 --- a/examples/whm_gr_test/param.swifter.in +++ b/examples/whm_gr_test/param.swifter.in @@ -1,6 +1,6 @@ ! VERSION Swifter parameter file converted from Swiftest T0 0.0 -TSTOP 1000.0 +TSTOP 100.0 DT 0.0006844626967830253 ISTEP_OUT 1461 ISTEP_DUMP 1461 diff --git a/examples/whm_gr_test/param.swiftest.in b/examples/whm_gr_test/param.swiftest.in index ace6f3cad..c9b7462f0 100644 --- a/examples/whm_gr_test/param.swiftest.in +++ b/examples/whm_gr_test/param.swiftest.in @@ -1,6 +1,6 @@ ! VERSION Swiftest parameter input T0 0.0 -TSTOP 1000.0 +TSTOP 100.0 DT 0.0006844626967830253 ISTEP_OUT 1461 ISTEP_DUMP 1461 diff --git a/examples/whm_gr_test/swiftest_relativity.ipynb b/examples/whm_gr_test/swiftest_relativity.ipynb index fe0e42cc9..cd203cdd4 100644 --- a/examples/whm_gr_test/swiftest_relativity.ipynb +++ b/examples/whm_gr_test/swiftest_relativity.ipynb @@ -23,9 +23,9 @@ "output_type": "stream", "text": [ "Reading Swifter file param.swifter.in\n", - "Reading in time 1.000e+03\n", + "Reading in time 1.000e+02\n", "Creating Dataset\n", - "Successfully converted 1001 output frames.\n", + "Successfully converted 101 output frames.\n", "Swifter simulation data stored as xarray DataSet .ds\n" ] } @@ -46,9 +46,9 @@ "output_type": "stream", "text": [ "Reading Swiftest file param.swiftest.in\n", - "Reading in time 1.000e+03\n", + "Reading in time 1.000e+02\n", "Creating Dataset\n", - "Successfully converted 1001 output frames.\n", + "Successfully converted 101 output frames.\n", "Swiftest simulation data stored as xarray DataSet .ds\n" ] } @@ -116,16 +116,26 @@ "text": [ "Mean precession rate for Mercury long. peri. (arcsec/100 y)\n", "JPL Horizons : 571.3210506300043\n", - "Swifter GR : 571.6748624042897\n", - "Swiftest GR : 571.6748624128945\n", - "Obs - Swifter : -0.3538117742853119\n", - "Obs - Swiftest : -0.35381178289026777\n", - "Swiftest - Swifter: 8.604956747149117e-09\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,) " ] }, { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEGCAYAAAB2EqL0AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAA40klEQVR4nO3deZxN5R/A8c/3zj72rUWytSjLGIwt2WXPloqUpYRSaUeL5PdTKhVtZFfJElkiS0mhnyyDQgihhiIiDGOW+/39cS+NMcO9zMyZ5ft+ve7LPc89zznf52K+c57nnOcRVcUYY4y5GJfTARhjjMkeLGEYY4zxiSUMY4wxPrGEYYwxxieWMIwxxvgk0OkAMlLRokW1dOnSTodhjDHZRnR09CFVLZbaZzk6YZQuXZp169Y5HYYxxmQbIrI3rc+sS8oYY4xPLGEYY4zxiSUMY4wxPsnRYxipSUhIICYmhri4OKdDyXVCQ0MpUaIEQUFBTodijLkEuS5hxMTEkC9fPkqXLo2IOB1OrqGqHD58mJiYGMqUKeN0OMaYS5DruqTi4uIoUqSIJYtMJiIUKVLEruyMycZyXcIALFk4xL53Y7K3XJkwjDEmp/rs6/d45ZP7SUpKSPdjW8JwQN68edmzZw9hYWFERkZSvnx5+vTpg9vtZs+ePVSsWPGC9QcPHszw4cPPKStdujSHDh3yK46WLVty9OhRf8M3xmRBcaeO0HtCLYbs+5Bv49Zw7MTRdD9Hrhv0zkquu+46Nm7cSGJiIo0aNWLOnDlUrVo1w8+rqqgqX375ZYafyxiT8WZ89TZz93zMT8EJ1DhemG51BlCoQKqze1wWu8LIAgIDA7nlllvYuXNnuhzvrbfeomLFilSsWJERI0YAsGfPHm6++WYefvhhqlatyu+//372qmT06NFERkYSGRlJmTJlaNiwIQBTp06lUqVKVKxYkf79+589ft68eXn++eepXLkytWrV4sCBAwB89tlnVKxYkcqVK1OvXr10aYsxJm0790bTd1xd/rN/Aj8FJ9A0viQf9PmGelVaZMj5cvUVxstfbOHn/cfS9Zjli+fnpdsr+FXn5MmTLF26lCFDhvhc5+233+aTTz45u71//34AoqOjmThxIqtXr0ZVqVmzJvXr16dQoUJs376diRMn8sEHH5xzrD59+tCnTx8SEhJo1KgRTz75JPv376d///5ER0dTqFAhmjZtypw5c2jXrh2xsbHUqlWLoUOH8uyzzzJ27FheeOEFhgwZwuLFi7nmmmusq8uYDOROSuT7TUt5Zd1T7AuEqKNX8uhtz1DlpqaIK+OuA+wKw0G7du0iMjKSOnXq0KpVK1q08P23gieeeIKNGzeefRUvXhyAlStX0r59e/LkyUPevHnp0KEDK1asAKBUqVLUqlUrzWP269ePRo0acfvtt7N27VoaNGhAsWLFCAwMpEuXLixfvhyA4OBgWrduDUC1atXYs2cPAHXq1KF79+6MHTuWpKSkS/lKjDEXsf6X9XSbUIuHf3yamCChV772vNJ9PlXLN8/QZAG5/ArD3yuB9HZmDCM9qWqan+XJkyfNzyZNmsTevXt57733LnqcoKCgs7fIBgQEkJiYCMDo0aNZvXo1CxYsIDIyko0bN1KkSJFLaYYxJoWkxHg+/3Y0Q/aNhWCoGVuIxmXr07n5fzItBrvCyGHq1avHnDlzOHnyJLGxscyePZu6detesE50dDTDhw/nk08+weX9DaVmzZp89913HDp0iKSkJKZOnUr9+vUveJxdu3ZRs2ZNhgwZQtGiRfn999/TrV3G5Faqyu69K2n4URWG7BtLmNvNXa47+E+XxZmaLCATrzBEpBwwPVlRWWAQsAwYDeQF9gBdVPW8gQURaQ6MBAKAcao6LKNjzgiJiYmEhIRccJ/t27dTokSJs9tvv/02d955p0/Hr1q1Kt27d6dGjRoA9OzZkypVqpztNkrNe++9x99//312sDsqKopx48bx6quv0rBhQ1SVli1b0rZt2wue+5lnnmHHjh2oKo0bN6Zy5co+xWyMSduTY2/n65C9BLqESrEhdCx3Lx0aP+5ILHKhrocMO6lIALAPqAnMBJ5W1e9E5H6gjKq+mMr+vwC3ATHAWqCzqv58ofNERUVpygWUtm7dys0335xubfHXjz/+yIMPPsiaNWsci8FJTn//xmQX7826nzlH13Ag0NP9e7erBs/fOy7DZ0wQkWhVjUrtM6fGMBoDu1R1r/fKY7m3/CtgMfBiiv1rADtV9VcAEZkGtAUumDCymtGjR/POO++cvdXVGGNS2n9wN8PnPMpXIXshUKh2Ih+v3TOXKwul/3MV/nIqYXQCpnrfbwbaAHOBO4FrU9n/GiB5h3gMnquT84hIL6AXQMmSJdMp3PRx5vZVY4xJKSExgfFfPMe8wwv5PURwqfLmTS9Qp+qdhAUHOB0e4MCgt4gE40kQn3mL7gf6ikg0kA+IT61aKmWp9qWp6hhVjVLVqGLFnM/IxhhzMT/+PJv+E9vw/rFF/B4kNDlVjqm3vEGTWp2yTLIAZ64wWgDrVfUAgKpuA5oCiMiNQKtU6sRw7pVHCWB/BsdpjDEZKjHJzY9b59Jj3YtosFAqHu4u2ZfWdbpTKG+o0+Gdx4mE0Zl/u6MQkStU9aCIuIAX8NwxldJa4AYRKYNnsLwTcE9mBGuMMRkh7tQRBn7cheVBv+ESqB97E0+0GUzpEs4+H3YhmZowRCQcz51OvZMVdxaRvt73nwMTvfsWx3P7bEtVTRSRR/AMiAcAE1R1SyaGbowx6WbCFw/z9t8rIARujg+kXpFOPNy1Py5X1l4zJlPHMFT1pKoWUdV/kpWNVNUbva8B6r3PV1X3q2rLZPt96d3nOlUdmplxp7ehQ4dSoUIFIiIiiIyMZPXq1T7VGzRoEF9//TUAK1asoEKFCkRGRrJq1ap0mXn2wIED3HPPPZQtW5Zq1apRu3ZtZs+eDcC3335LgQIFqFKlCjfddBNPP/30ZZ/PmNwm9uRxnh3bkw8OeW4Mbe9uyEfd1vJIhwFZPllALp8axAmrVq1i/vz5rF+/npCQEA4dOkR8fGrj/OdLPjnhlClTePrpp+nRoweTJk1i3bp1tGzZ8gK1z5WYmEhg4L9//apKu3bt6NatG59++ikAe/fuZd68eWf3qVu3LvPnz+fUqVNUqVKF9u3bU6dOHZ/PaUxu5XYrz33Ukm/1N2KDXRRPgGcrvUbjKN//z2YFljAy2R9//EHRokXPPu1dtGhRANasWcOwYcP4/PPPmTt3Lp06deKff/7B7XZTvnx5fv31V7p3707r1q05evQoM2bMYPHixSxZsoTvv/+eU6dOsXLlSgYOHEjr1q159NFH2bRpE4mJiQwePJi2bdsyadIkFixYQFxcHLGxsXzzzTdn4/rmm28IDg4+57bfUqVK8eijj57XhjMLP+3bty+Dvy1jsr9V66cya/1UFgfF4ELoHNCCgfcNy/CJAjNC7k4YCwfAn5vS95hXVYIWac9a0rRpU4YMGcKNN95IkyZNuPvuu6lfvz5Vq1Zlw4YNgKe7qWLFiqxdu5bExERq1jz3kZOePXuycuVKWrduTceOHc9eYZyZOPC5556jUaNGTJgwgaNHj1KjRg2aNGkCeK5wfvrpJwoXLnzOMbds2eLz4k1Hjhxhx44dtuaFMRew7+/DfLroZT5KWAZBcNPpAMZ1WUmBPHmdDu2SZb8Ul83lzZuX6OhoxowZQ7Fixbj77ruZNGkSgYGBXH/99WzdupU1a9bw5JNPsnz5clasWHHRyQNTWrJkCcOGDSMyMpIGDRoQFxfHb7/9BsBtt912XrJITd++falcuTLVq1c/W7ZixQoiIiK46qqraN26NVdddZV/jTcmF3AnJbJ951f0n9GMjxKWEeZ283ihdkzplr2TBeT2K4wLXAlkpICAABo0aECDBg2oVKkSkydPpnv37tStW5eFCxcSFBREkyZN6N69O0lJSeet330xqsqsWbMoV67cOeWrV69Oc4rzChUqMGvWrLPb77//PocOHSIq6t8pZc6MYfzyyy/ceuuttG/fnsjISL9iMyYnU7eb/p80YxEHIQwax5WiZ/P3qFiqtNOhpQu7wshk27dvZ8eOHWe3N27cSKlSpQDP1OQjRoygdu3aFCtWjMOHD7Nt2zYqVLjwfdn58uXj+PHjZ7ebNWvGu+++e3ZNizNdXRfSqFEj4uLiGDVq1NmykydPprrvjTfeyMCBA3nttdcuelxjcouV0bPo8mF1FnGQ4gnKYwVaMaL3/ByTLMASRqY7ceIE3bp1o3z58kRERPDzzz8zePBgwLMGxYEDB86ODURERBAREXHR2SkbNmzIzz//TGRkJNOnT+fFF18kISGBiIgIKlasyIsvppzL8Xwiwpw5c/juu+8oU6YMNWrUoFu3bmkmhT59+rB8+XJ2797t3xdgTA5z9PghBkxqwZM/DWJ3SBytk67nsy7reLBdtlyB4YIcmd48s2TF6c1zO/v+TU7y9aoxvL/pPXaGKBVjg+lZYyiNo5o7HdZlyYrTmxtjTLZ1Ii6Odz7vwtSEXyAE7g2qQa9uoyiUJ9jp0DKUJQxjjPGRut0s/t9ohm1/n8OBLirFhtCtcm+a1XnQ6dAyhSUMY4zxQULCSR6a3JDVQScJd0HngIa0bv0iESVzzzIKljCMMeYipix6k89+m8yuEKVO7JW0ibiXlrd2dzqsTGcJwxhj0rD91xUMWvooPwcnERyk3Ha6LIO6z6RgeM4eq0iLJQxjjEkhKSmB92Y9x9zjX3IqAKr9U5THmg2jarnaTofmKHsOwwFZaXrzo0eP8sEHH6T5uU15bnKbiV+8QKfx1Rl3ahF/BbroUqABEx75NtcnC7CEkemST2/+008/8fXXX3PttddevCKe6c3PTCJ4ZnrzjRs3sn379gxJGGemPK9Xrx6//vor0dHRTJs2jZiYmLP71K1blw0bNrBhwwbmz5/P999/f0lxGOO0TVvnMWLGI7z191y2hSRR63gxplQfTN/272aLtSoyQ6YlDBEpJyIbk72OicjjIhIpIj94y9aJSI006j8hIltEZLOITBWRrLfgrQ9Sm968ePHirFmzhg4dOgAwd+5cwsLCiI+PJy4ujrJlywLQvXt3Zs6cybhx45gxYwZDhgyhc+fODBo0iOnTp5990js2Npb777+f6tWrU6VKFebOnQt4ZqStUaMGkZGRREREsGPHDgYMGMCuXbuIjIzkmWeeOSdWm/Lc5AYJCYnMXzGBB38YyPhT3wHwXKnX+G/XL4kof0e2nIY8o2TaGIaqbgciAUQkAM/a3LOBscDLqrpQRFoCrwMNktcVkWuAx4DyqnpKRGbgWdd70uXE9Nqa19j297bLOcR5bip8E/1r9E/zcyenNx89ejT9+vWjS5cuxMfHk5SUxLBhw9i8eTMbN248L1ab8tzkdIcObePteS8zL2AzuFzUii1MmwoduL1+9lrYKLM4NejdGNilqntFRIH83vICwP406gQCYSKSAIRfYL8s7cz05itWrGDZsmXcfffdDBs2jO7du6c6vXlSUtIlTW8+b968s7PcnpnevHbt2gwdOpSYmBg6dOjADTfc4Ndx+/bty8qVKwkODmbt2rXAv1Oeb9++nQEDBtiU5ybbGP/FQ4z4eyUEeLa7BtelX5f3CA60K4q0OJUwOgFTve8fBxaLyHA8XWS3pNxZVfd5P/8NOAUsUdUlqR1YRHoBvQBKlix5wSAudCWQkZya3vzmm2+mZs2aLFiwgGbNmjFu3Liz3V2psSnPTU508tRRnp58ByvCDhLmdlMlrhiPNXuXCmUrOR1alpfpqVREgoE2wGfeooeAJ1T1WuAJYHwqdQoBbYEyQHEgj4jcm9rxVXWMqkapalSxYlnvCUwnpzf/9ddfKVu2LI899hht2rThp59+Oq9ucjbluclpPvi8Fx0+vZUVYQfJl+Tmxeue4MOHvrVk4SMnrr1aAOtV9YB3uxvwuff9Z0Bqg95NgN2q+peqJnj3P+9KJDtwcnrz6dOnU7FiRSIjI9m2bRtdu3alSJEi1KlTh4oVK5436G1TnpucYtuv0Tw5pi2jjq9iX6BQ59gNLO+2idvr9XQ6tGwl06c3F5FpwGJVnejd3go8pKrfikhj4HVVrZaiTk1gAlAdT5fUJGCdqr57oXPZ9OZZj33/JjMlJLn5/NtxfLhnJH8Furg6QRnWcC43X1uasOAAp8PLkrLM9OYiEg7cBvROVvwgMFJEAoE4vOMPIlIcGKeqLVV1tYjMBNYDicAGYExmxm6MyT7U7WbZ6neZu+lLvgnZT1CA0CmwA/06PknePAWcDi/bytSEoaongSIpylYC1VLZdz/QMtn2S8BLGR2jMSZ7S3Ir0xcP4NWDCyEEbo4LoP0Ng7ijXju7A+oy5cq5pFT1ouMCJv3l5NUdjfMSE+L4J/YQz0+5n1V59pPXrdwR3pLeHZ8nX56CToeXI+S6hBEaGsrhw4cpUqSIJY1MpKocPnyY0NBs+YC+yQYe/bgRKwOOQ16IjAulR7XBNIpq7XRYOUquSxglSpQgJiaGv/76y+lQcp3Q0FBKlCjhdBgmh1m8cgSzNs9iVZjn9vD7Q1vzRLdXHY4qZ8p1CSMoKIgyZco4HYYx5jIdPHKAEXN78YX8CmFQ7XQo73b5yrqfMlCuSxjGmOwtKTGepavfZeKWKWwOSaD0aXi40n9pUbut06HleJYwjDHZRlx8PM9+0oRlAUcgBNq7q/FI59Fckc/GxjKDJQxjTLawY+86nl/Uk62hSVQ8FUybMnfQudlzToeVq1jCMMZkaYcO7eK/8x5gpRwiIFi5y1WDR+8bTcE8IU6HlutYwjDGZFlfrhjFh9s+4PcgpeLx/LSp+BAdG97ndFi5liUMY0yWs/ePnYxd8iRz2Q3BcH/YrTzebZQ9O+UwSxjGmCwjMSGOyQuf472/l5AoQsTxMJ5pNobIGyKdDs1gCcMYk0X8fSKWZ6Y3Yk3gSRDhwTxtaNLkGcpfU9Dp0IyXJQxjjKNOnjjIjG8/5MOD0zgR6OKW2GJ0qd6NetW7OR2aScEShjHGMbv3rOSeZb054XKRH+igNejT+V2uLhTudGgmFZYwjDGZTt1uXph0L/MCNhGMEHWsMHdX60XzW7o4HZq5AEsYxphMNWp2L2b8/T2HvGtT3BNWnad6THQ4KuOLTEsYIlIOmJ6sqCwwCPgWGA2E4llN72FVXZNK/YLAOKAioMD9qroqY6M2xqSXVRsmM3vDNBYGxECgixrHrmB4j/kUyhvmdGjGR5mWMFR1OxAJICIBwD5gNjAWeFlVF4pIS+B1oEEqhxgJLFLVjiISDFgnpzHZgDspiXHzX2LUkTkkBgiiyuvlXuDWqh3JG2KdHNmJU39bjYFdqrpXRBTI7y0vAOxPubOI5AfqAd0BVDUeiM+cUI0xl2rtxo9YsnkR05I2gQitEqpzZ83mVKt0l9OhmUvgd8IQkTxAnKomXcZ5OwFTve8fBxaLyHDABdySyv5lgb+AiSJSGYgG+qlqbCrx9QJ6AZQsWfIyQjTGXKq4hCQWr3iTF37/GIDCiW4eKdOfNvW6EBIY4HB05lJddEV0EXGJyD0iskBEDgLbgD9EZIuIvCEiN/hzQm93UhvgM2/RQ8ATqnot8AQwPpVqgUBVYJSqVgFigQGpHV9Vx6hqlKpGFStWzJ/QjDHppP/Elrzw+8cEqtL4VDleq/Eudzbqaskim/PlCmMZ8DUwENisqm4AESkMNASGichsVf3Ex3O2ANar6gHvdjegn/f9Z3gGtlOKAWJUdbV3eyZpJAxjjHPGzXuETw4t43CIi9Kn4YGbnqddvc5Oh2XSiS8Jo4mqJqQsVNW/gVnALBEJ8uOcnfm3Owo8Yxb18dwt1QjYkcq5/hSR30WknHfwvDHwsx/nNMZkoF/2rGf4gqdZH36QAm6oF1eSwZ0/plihwk6HZtLRRRNGasniUvYBEJFw4Dagd7LiB4GRIhIIxOEdfxCR4sA4VW3p3e9RYIq3S+tXoIcv5zTGZBx3UhKvT+vFwtM/cCyPcF18AIMafkjEDbWcDs1kAJ8HvUXkyVSK/wGiVXWjL8dQ1ZNAkRRlK4Fqqey7H2iZbHsjEOVrvMaYjDV+Xl+i/9zCipDDhAm0dd3B4F4vOx2WyUD+3CUV5X194d1uBawF+ojIZ6r6enoHZ4zJek7GJfDRopd5/5/lEAI3nQ7gv62XcePVRS5e2WRr/iSMIkBVVT0BICIv4Rl8rofnNldLGMbkYAmnY1m/ZRofrh3L2tBYCiW6eeDq7nRq1JuQ0PwXP4DJ9vxJGCU592G5BKCUqp4SkdPpG5YxJqvpP7UVX8lhCIX6p4rTs/FrtrBRLuNPwvgU+EFE5nq3bwemeh/kszuWjMmhFq54i0+3TGVjWBxlTkOrK++gd7fBTodlHOBzwlDV/4jIl8CtgAB9VHWd92Obk9iYHOZUXBwvT2vHAtkHYdAssTjPdZ5O4XwFnQ7NOMSfu6QEuBkooKpDRKSkiNRIbWZZY0z25U5K5OtV7/HRlk/4MfQ0lU+E0bPuBzSIsJsUczt/uqQ+ANx4Hq4bAhzH8+Be9QyIyxjjgCOxsbw0rRXLAj1jFR1d5RjQe7pN6WEA/xJGTVWtKiIbAFT1iPchOmNMDrBj91qeW/Ig20KTqHoiL3dHdKNFnV6I66JTzplcwp+EkeBdx0IBRKQYnisOY0w2djz2CE9ObcEPQbEQCp0DbuWBe9/kygK25Iw5lz8J4x08Cx5dKSJDgY7ACxkSlTEmU8xc+gHjf/2AmGCh+vHCNL6+PV2aP+50WCaL8ucuqSkiEo1n4j+Adqq6NWPCMsZkpF92r2Hk0qdYHnSU8EDlvqDaPN7nQ4IDrfvJpO2iCSONOaQAWohIC1V9K51jMsZkkPjTJ/jvtPuYzU4CApUqx/LTq8Fb3FrJJgs0F+fLFUY+75/l8NwRNc+7fTuwPCOCMsakvy+/n8KnP7/Fj8GeCRseLticu+56hYLhdu+K8Y0v05u/DCAiS/DMJXXcuz2Yf1fNM8ZkUVu2f8nnaz5iVtJmkoKFW/65jhfveJES15w3SbQxF3Q5c0nFA6XTNRpjTLpxu5Wf9/xMz/89Q6zLBSL0LdiFtu36cXWBMKfDM9mQPwnjY2CNiMzGc2tte2ByhkRljLks8aeP8/zHXVgUtBtcLm45fjVtKnekVd1eTodmsjF/7pIaKiILgbreoh6qusHX+iJSDpierKgsMAjP0qyjgVAgEXg4relGvM+BrAP2qWprX89tTG7y+dcv8dK+zyEIiia6aZuvLo/dNwqXS5wOzWRzvtwlJaqqAKq6Hlh/oX3S4l2LO9K7fwCwD89zHWOBl1V1oYi0xLOuRoM0DtMP2ArY5PvGpPBbzFreWvgCS4P341Ll1tgy9L9rMiWL2braJn34ctP1MhF5VERKJi8UkWARaSQik4Fufp63MbBLVffi6d46kwAKAPtTqyAiJfCs8jfOz3MZk+ON/KwXrZbez9Lg/YgqT1/Riff7fmHJwqQrX7qkmgP341n7ogxwFAjDk2yWAG/7uqZ3Mp2Aqd73jwOLRWS495i3pFFnBPAs/97mmyoR6QX0AihZsuSFdjUm29u4ZS5z1s1gFj8BUO9EY17r2p+8+a52ODKTE8lFepLO3VkkCCgKnFLVo5d0Qs+EhfuBCqp6QETeAb5T1VkichfQS1WbpKjTGmipqg+LSAPgaV/GMKKionTdunUX282YbOdY7ElGzxvItNNLSRAhb5Kb0fWmUq5UBUKDbGZZc+lEJFpVU53L3p+7pFDVBOCPy4ynBbBeVQ94t7vhGZsAz3MdqXU51QHaeMc4QoH8IvKJqt57mbEYk+38sT+aN758nq+C9hGmSlupT7tbOlH5+ginQzM5nF8JI5105t/uKPBcbdTHc7dUI2BHygqqOhAYCJDsCsOShclVDh6LY+mqN3j1j+lokFD+dCCP1XqHOhF1L17ZmHSQqQlDRMKB24DeyYofBEaKSCAQh3f8QUSKA+NUtWVmxmhMVnT02J88OaU1P4afJr9bqasR9G33Gtdeca3ToZlcxJ8lWh9Q1fGXczJVPQkUSVG2EjhvjgJV3Q+clyxU9Vs8VyPG5AoDJrZiges3XGFK1bhwelZ/jrpV2zkdlsmF/LnCeFNEuuB5uG4NMFVVt2RMWMaYCV88wvJ9a4kOO0nRRDeNg1vxfK/XELEH8Iwz/EkYh4H/AsF4HsCbISLvqOqHGRGYMbnV4X+OMWb+U3ya+AOEQaXTQYy5dyl5wws5HZrJ5fxJGP+o6jfe94tEZCSwGrCEYUw6SEqMZ/ayQXy2ezE/hyRybbzySOVRtKxhg9oma/B70FtE+uN5FqMAcDzdIzImF0pyK8NmdGJawg4IgY4SQZ87x3JlQVtX22Qdl3KX1Cw8U3u0BV5J33CMyX2Wr53IqPXvsjk0gRviXHS/6THa1H/A6bCMOY8/CaOQiFyrqjuBnSIyFtgALMiY0IzJ2ZKSEnjhk9uZzz5CgpXb9ToevWsyVxcq6HRoxqTKn4SRH/hWRA4BPwMFgaSMCMqYnEzdbhasGMmMbZ+yITSOqOMF6XLrSJpEVnU6NGMuyJ+E0RDYDNTEs763YlcXxvgl5uA+XplzDytC/oZQaJ50Da/2WUBgoM3/ZLI+fxZQ+sn7dpX3ZYzxw+qNn/LSulfYFyJE/VOM3g0HULNCE8TlyyoDxjjPibmkjMlVdh84zJC5LVgXcgqChAfzNOXejq9SOE+w06EZ4xdLGMZkoE8XvcGbf0wmPkSoebIAjcu2pnOzAU6HZcwl8WcuqUeAKap6JAPjMSZH2P1bNI8t6c6eICjsVu7I15z7Ow0jb4j9jmayL3/+9V4FrBWR9cAEYPHF1vE2JrdxJyXyzKSWLAn8g+AApcbREtxZsw/Na7VzOjRjLps/g94viMiLQFOgB/CeiMwAxqvqrowK0JjsYvy8/nx5YBG/BLsBuC9vTR7vcVkTPBuTpfi74p6KyJ/An3hmrS0EzBSRr1T12YwI0Jis7pddSxnz3TAWB/wJwVD96FW89cB8CuYNcTo0Y9KVP2MYj+FZTvUQnmVUn1HVBBFx4Vkl74IJQ0TKAdOTFZUFBuFZ22I0nqVXE4GHVXVNirrXAh/h6RZzA2NUdaSvsRuTEZISk5i5bDRjfxvFgUAhQJW3K75K1fLNKRAe5HR4xqQ7nxKGeCbgrwx0UNW9yT9TVbeItL7YMVR1O55p0RGRAGAfMBsYC7ysqgu9a3a/DjRIUT0ReEpV14tIPiDae1Xzsy/xG5Petu5YyNT/TWA22yBQaBhbgR4NulKlvC0QaXIunxKGtyuqSspkkezzrX6etzGwS1X3iojimXYEPDPg7k/l+H8Af3jfHxeRrcA1eKYoMSbTqCrrN83g/vX/wS1CuNtNryKtueOOIdYFZXI8f8YwVolIdVVdmw7n7QRM9b5/HFgsIsMBF3DLhSqKSGmgCp61OIzJNPGnj/PUpA6sCv4Dt0uo808lHm//Gjdda+tqm9zB37mk+ojIHiAWEDwXHxH+nFBEgoE2wEBv0UPAE6o6S0TuAsYDTdKomxfP9OqPq+qxNPbpBfQCKFmypD+hGZOmkZ/1ZHzsD2ioUCIB2hZsTe+ur9hyqSZXEV8fpRCRUqmVp9VNdYHjtAX6qmpT7/Y/QEFvt5fgWdkvfyr1goD5eJ7/eMuXc0VFRem6dev8Cc+Yc+zcvYJx34zga9d2glWpHleHV3t+QHiwPYBnciYRiVbVqNQ+8+dffbc0yof4GU9n/u2OAs+YRX08d0s1wnPH1Tm8iWQ8sNXXZGHM5UhISOSZia1ZGrIPAiHUDa9FDKZu1J1Oh2aMY/xJGLHJ3ocCrQG/BrtFJBy4DeidrPhBYKSIBAJxeLuTRKQ4ME5VWwJ1gPuATSKy0VvvOVX90p/zG+OLb1aNZO6WuXwT8hcAd7o68vSdfQkPL+pwZMY4y+cuqfMqioQA81S1WfqGlH6sS8r4Y9+R40xb/BKTTn8FwHXxwuiOy7nKVsAzuUh6dUmlFI7n4TtjsjV1u9m99zuGLBpIdHgseZPcdM7fmnvbPU7hAgWdDs+YLMOfJ7034VllDyAAKIb/4xfGZDnPftSMRfInhEON+Dw8dMswoio0cDosY7Icf64wkj/NnQgcUNXEdI7HmEyzbPX7TNowmfVhp7jxlIsGRVvxSNehdqusMWnwZ7Zav26fNSarOnU6gTemd2Vu0ibiw4QG8UX4b5c5FMhX0OnQjMnSfF5MWEQmi0jBZNuFRGRChkRlTAZZtGI490+qwWe6mVKnA3i2xKu8++C3liyM8YE/XVIRqnr0zIaqHhGRKukfkjHp72RcAsNm3MVs3QmhcHfgTTzZ/VPCQ21WWWN85U/CcIlIoTNLtIpIYT/rG5Pp1O1mwfJXGf/LdHaGKJEnwulxyxs0qlLP6dCMyXb8+YH/JvA/EZmJ526pu4ChGRKVMekgKSmB/h+1YLHrAOFBbjoQwQOdPqRkkfNmnjHG+MCf9TCWAevwTN8heNbGsOnFTZY0+5uRzPrlY34MO03N2EJ0jHqR5jVuczosY7I1f9bDmKOq1bA1KEwWtnffVl6e34O1obEQBk2TrmBoz0WEBttYhTGXy58uqR/ScT0MY9JVYkIcExY8y/gjS0kIEar9XYKnOoymUqnSTodmTI6R6ethGJPe5i7/iFlb32VDaBy4XPTOW5eHu76Py2UP4BmTnvxJGC0yLApjLsGBAz/x1fqpvHngCxJDhVrHr+Demt2oX72r06EZkyP5kzB+A7oAZVV1iIiUBK4C7Alwk+l+27+ZuxZ3JtblQoCBVz1Aq7v6UiDMxiqMySj+JIwPADeeu6SGAMfxLJdaPQPiMiZNL05uzxx2gstF/RPlaF65Da1vtasKYzKaPwmjpqpWFZENcPZJ7+AMisuY88xY/AJv7Z9NrMvF1QlK8/xRPNltktNhGZNr+JMwEkQkAO8U5yJSDM8Vh09EpBwwPVlRWWAQnqVZR+NZxS8ReFhV16RSvzkwEs/U6uNUdZgfsZts7K+/tvPK3Ef5OugPcLmocawoT7WfQfkSxZwOzZhcxZ+E8Q4wG7hSRIYCHYEXfa2sqtuBSABv4tnnPd5Y4GVVXSgiLYHXgQbJ63r3fx/P8q4xwFoRmWcPDuZ8oz5/idlHZvFHkCCqvFjiAe7s9oTTYRmTK/kzvfkUEYkGGnuL2qrqtks8b2Ngl6ruFREFzszVUADYn8r+NYCdqvorgIhMA9piDxHmWGt/msGUNR+yNOggBAq1j5bnjR4jKZD/KqdDMybXumjCEJF5KYu8fzYTEVS1zSWctxMw1fv+cWCxiAzHM936Lansfw3we7LtGKBmGvH2AnoBlCxZ8hJCM06KT3Tz7cblDNr0MrFBLq5McPN8lbeocGN9CuQLdTo8Y3I1X64wauP5YT0VWM2/CeOSeAfK2wADvUUPAU+o6iwRuQsYDzRJWS2VQ2kqZajqGGAMQFRUVKr7mKxH3W7W/jiBmRvmsDBgL7hc3J7YnD6tH6HkVaWcDs8Yg28J4yo8YwedgXuABcBUVd1yiedsAaxX1QPe7W5AP+/7z4BxqdSJAa5Ntl2C1LuuTDb16eKXGHZwDgTAdaeF26+6i+6tnyfAntY2Jsu4aMJQ1SRgEbBIRELwJI5vRWSIqr57CefszL/dUeD5wV8fz91SjYAdqdRZC9wgImXwDJZ3wpO8TDb318FtvPR5X1aFHiAQaJHUiMfuHMxVhQo7HZoxJgVfpzcPAVrh+WFfGs8dU5/7ezIRCcdztdI7WfGDwEgRCQTi8I4/iEhxPLfPtlTVRBF5BFiM57baCZdxhWOyAFWlz9hb+V/IMYJClWsShN439uH2+n2dDs0YkwZRvXA3v4hMBioCC4Fpqro5MwJLD1FRUbpu3TqnwzApTFv8NF/tXcGakJPkTXLTLKAtL3UdimfZFWOMk0QkWlWjUvvMlyuM+/DMTnsj8Fiy/9RnZqu15cuMTw4fO8a7sx9lFushBErFw5RO31EgX1GnQzPG+MCXMQxXZgRici53UiLL133A5B8/YV3IKQolunn0uiF0uPV2AgJtdhljsgt/nvQ2xm9xCUm88ukdzOZXCIEm8cV5psPHFC9yhdOhGWP8ZAnDZJhdu1fw8uKn2BB2ipvjAmh1bTfua/64LWxkTDZlCcOku4TEBIZM6cR83U5gCLTR63i2y0cUyFvA6dCMMZfBEoZJV9+vn8H7a19jU2g8EbFhNL3xWbrddqfTYRlj0oElDJMuTp1OYPj0Psx2ryYhVOgg1zP4oVmIy+6ZMCansIRhLos7KZEJC/oy+dBKjga4uPlUEE/Un0Tt8pWdDs0Yk84sYZhLFp+QyPOftGCR608IcHGXqxJd7nifslcUcjo0Y0wGsIRh/JaUGM+S/73PqG0T2B0CNWMLcF+Np6gf1d7p0IwxGcgShvHLb39sofeCu4kJEsKD3HRwXU//B2YQHhLidGjGmAxmCcP4xJ2UyMtTOrI4aQfxgULNv8tzX+P+1I+o6nRoxphMYgnDXNSsb95l1s6JbApJAJeLrsFVeeaJyU6HZYzJZJYwTJq2/fIlM9aMYXbiThJDhBrHi/DafZ9TtICtVWFMbmQJw6Rq2bp5vPDTQI4FuHABr18/gAbVOxEWHOB0aMYYh1jCMOc4fmwfM759hxFHvoQAF/WOV+Lu2h2oV62j06EZYxyWaQlDRMoB05MVlQUGAbWBct6ygsBRVY1Mpf4TQE9AgU1AD1WNy8CQc53layfTb8sbJIrgUqVfkYb06PqOLWxkjAEyMWGo6nYgEkBEAvCszT1bVUec2UdE3gT+SVlXRK4BHgPKq+opEZmBZ13vSRkeeC6QlBjPUxNvZ2nwfhCh7omy3FXvvzSoVMnp0IwxWYhTXVKNgV2quvdMgXh+jb0LaJRGnUAgTEQSgHBgf4ZHmQuMnfcKn/31KX8EC1cnKPdd3Z77uv3H6bCMMVmQUwmjEzA1RVld4ICq7ki5s6ruE5HhwG/AKWCJqi5J7cAi0gvoBVCyZMl0DTonid40k1H/e5PVoScgUKhxpASv9JjFlQXCnQ7NGJNFZXrCEJFgoA0wMMVHnTk/iZypUwhoC5QBjgKfici9qvpJyn1VdQwwBiAqKkrTL/KcQVV57/MhfHpsBidCXYS73Qyv/B+qV2pLaJDdAWWMSZsTVxgtgPWqeuBMgYgEAh2AamnUaQLsVtW/vPt/DtwCnJcwTNpmff0Cy/d8zzdBh8Dlok1iU57rNIA8eYo5HZoxJhtwImGkdiXRBNimqjFp1PkNqCUi4Xi6pBoD6zIuxJwl9nQiX62Zy+B9cyEISsVDv4iXaFz9Dlsu1Rjjs0xNGN4f+LcBvVN8dN6YhogUB8apaktVXS0iM4H1QCKwAW+3k0lbUmI8O35dzBvfvMKasBOEu93cFX47XVr146qiVzsdnjEmmxHVnNvNHxUVpevW5d4LkUEf385s9x4Aap3Kzx0VH6D5rfc7G5QxJksTkWhVjUrtM3vSOwf65n9vMGnLdH4MigMROgc05tkH3yIwwJZLNcZcOksYOUhikptXpjzI5+7VJAULlU+H8J+WH1OmRHmnQzPG5ACWMHKIH9aPZ3L0RFYG/0OJBGhX4j/0btHO6bCMMTmIJYxs7mjsaYbPuJe5rm0QDHVPF2ZE9yUEB9sKeMaY9GUJIxv7+ZclDF02kJ9C44mIDeHOSgNpW6c94rKxCmNM+rOEkQ2dOp3Aq9PuZDa7CAhR7nZF0rfrOArlDXU6NGNMDmYJI5tZ+sNEPvxxJFtDk6h8PJwm1z9C9+b3OR2WMSYXsISRTfx19C9en9mDrwL3kBQqdAy4jhf6zCIg0OZ/MsZkDksYWVxiQhyvz7yHqfE7IAgijofTu/6H1KtUxenQjDG5jCWMLGz3n78z8suuLA04BMADYTW5p/17XJHfxiqMMZnPEkYWdOL4H8z87g3GH1zM0QAXtU4U49Emg4go18Dp0IwxuZgljCzmyD/76DizKQcDXQS5hJ7hNeh59wfkCbXnKowxzrKEkUW4kxJ5/pPbmU8MBLpoGFuBFtW60aJmC6dDM8YYwBJGlvD1qgmM2PI2e4Mg3O2mReANDH54mtNhGWPMOSxhOOj3mNWM/OoFFrv+hCCodrQo/e+czc3FCzodmjHGnMcShkOmLH6fj34fxf4gz4p3/Qq2oMd9rxFgK+AZY7KoTEsYIlIOmJ6sqCwwCKgNlPOWFQSOqmpkKvULAuOAioAC96vqqoyLOGPE7FvNnB/G8uGJ1RAk1DsewdNt+lGmZA2nQzPGmAvKtIShqtuBSAARCQD2AbNVdcSZfUTkTeCfNA4xElikqh1FJBgIz9CA01lcQhKzl89i5G8vE+tyEepWht7Un1pVOpE/NMjp8Iwx5qKc6pJqDOxS1b1nCkREgLuARil3FpH8QD2gO4CqxgPxmRJpOog7dYRXZvRlNpvA5aLJyeo0rnw7TWu3dzo0Y4zxmVMJoxMwNUVZXeCAqu5IZf+ywF/ARBGpDEQD/VQ1NuWOItIL6AVQsmTJdA3aX6rKjKVjeCXmXdwiFEpy07dER+5u+rKjcRljzKXI9IUTvN1JbYDPUnzUmfOTyBmBQFVglKpWAWKBAantqKpjVDVKVaOKFSuWTlH7b9fuFXQaU5X/7nsPtwh1j1/H2KZfWbIwxmRbTlxhtADWq+qBMwUiEgh0AKqlUScGiFHV1d7tmaSRMJymqrw0qRtfazTHQ12UjoeHyvWhZb2+TodmjDGXxYmEkdqVRBNgm6rGpFZBVf8Ukd9FpJx38Lwx8HMGx+m3CV/04qs/17I5OBFRoVlcdd7oNR7P8IwxxmRvmZowRCQcuA3oneKj88Y0RKQ4ME5VW3qLHgWmeLu0fgV6ZHC4PktISGT03MGMiV0FwVA8QZlyxyKKFirhdGjGGJNuMjVhqOpJoEgq5d1TKdsPtEy2vRGIysDw/KZuN6s3jGXUhrGsDzpNqFt5sFg/ujW5m5DQ/E6HZ4wx6cqe9L5EqsqoOQ8x6vj/IAganCpOn6ZDqVA2S+U0Y4xJN5YwLkH0T9MYuepNNoWc4spEpX3hbjzY5WmCAzP9pjNjjMk0ljD89PzE9iyQHSSFCvXii/BY0/cpV6qi02EZY0yGs4Tho1XrJzAheiw/BJ8g4ng4rSLe4J569Z0OyxhjMo0ljIs4ciKOoZ/ex1fBW3EHC7VO5+Odnt8QFmrrahtjchdLGBewesNkRqwdyeaQBCJiw+hUZRitazZEXDZWYYzJfSxhpOLw8ThGzLqXObIdQqA91/FMjxnkCwt2OjRjjHGMJYwUjsTGc8eMahwOdFHxVBC3XXsf97d6wumwjDHGcZYwUsgX7KaqXEWJsGt4rMt4AgNtrQpjjAFLGOcJDArlrfuXOh2GMcZkOTZ6a4wxxieWMIwxxvjEEoYxxhifWMIwxhjjE0sYxhhjfGIJwxhjjE8sYRhjjPGJJQxjjDE+EVV1OoYMIyJ/AXsvsXpR4FA6hpMdWJtzvtzWXrA2+6uUqhZL7YMcnTAuh4isU9Vctd6qtTnny23tBWtzerIuKWOMMT6xhGGMMcYnljDSNsbpABxgbc75clt7wdqcbmwMwxhjjE/sCsMYY4xPLGEYY4zxiSWMFESkuYhsF5GdIjLA6XjSi4hcKyLLRGSriGwRkX7e8sIi8pWI7PD+WShZnYHe72G7iDRzLvpLJyIBIrJBROZ7t3N0ewFEpKCIzBSRbd6/79o5ud0i8oT33/RmEZkqIqE5sb0iMkFEDorI5mRlfrdTRKqJyCbvZ++IiPgchKray/sCAoBdQFkgGPgRKO90XOnUtquBqt73+YBfgPLA68AAb/kA4DXv+/Le9ocAZbzfS4DT7biEdj8JfArM927n6PZ62zIZ6Ol9HwwUzKntBq4BdgNh3u0ZQPec2F6gHlAV2JyszO92AmuA2oAAC4EWvsZgVxjnqgHsVNVfVTUemAa0dTimdKGqf6jqeu/748BWPP/Z2uL5AYP3z3be922Baap6WlV3AzvxfD/ZhoiUAFoB45IV59j2AohIfjw/WMYDqGq8qh4lZ7c7EAgTkUAgHNhPDmyvqi4H/k5R7Fc7ReRqIL+qrlJP9vgoWZ2LsoRxrmuA35Ntx3jLchQRKQ1UAVYDV6rqH+BJKsAV3t1ywncxAngWcCcry8ntBc/V8V/ARG9X3DgRyUMObbeq7gOGA78BfwD/qOoScmh7U+FvO6/xvk9Z7hNLGOdKrS8vR913LCJ5gVnA46p67EK7plKWbb4LEWkNHFTVaF+rpFKWbdqbTCCebotRqloFiMXTVZGWbN1ub599WzzdLsWBPCJy74WqpFKWbdrrh7TaeVntt4Rxrhjg2mTbJfBc3uYIIhKEJ1lMUdXPvcUHvJepeP886C3P7t9FHaCNiOzB07XYSEQ+Iee294wYIEZVV3u3Z+JJIDm13U2A3ar6l6omAJ8Dt5Bz25uSv+2M8b5PWe4TSxjnWgvcICJlRCQY6ATMczimdOG9E2I8sFVV30r20Tygm/d9N2BusvJOIhIiImWAG/AMlmULqjpQVUuoamk8f4/fqOq95ND2nqGqfwK/i0g5b1Fj4Gdybrt/A2qJSLj333hjPONzObW9KfnVTm+31XERqeX9vromq3NxTo/8Z7UX0BLPHUS7gOedjicd23UrnkvPn4CN3ldLoAiwFNjh/bNwsjrPe7+H7fhxJ0VWewEN+PcuqdzQ3khgnffveg5QKCe3G3gZ2AZsBj7Gc2dQjmsvMBXPOE0CniuFBy6lnUCU97vaBbyHd8YPX142NYgxxhifWJeUMcYYn1jCMMYY4xNLGMYYY3xiCcMYY4xPLGEYY4zxiSUMYy5CRIqIyEbv608R2ed9f0JEPsigcz4uIl0vss80EbkhI85vTGrstlpj/CAig4ETqjo8A88RCKzHM7tw4gX2qw/cq6oPZlQsxiRnVxjGXCIRaZBsnY3BIjJZRJaIyB4R6SAir3vXHVjknZblzFoE34lItIgsPjOtQwqNgPWqmigi14nI+mTnvEFEzsyPtQJo4k0wxmQ4SxjGpJ/r8Eyn3hb4BFimqpWAU0Arb9J4F+ioqtWACcDQVI5TB4gGUNVdwD8iEun9rAcwyfuZG8+01ZUzqD3GnMN+MzEm/SxU1QQR2YRnMa5F3vJNQGmgHFAR+Mq7yFkAnqkeUroaz3xIZ4wDeojIk8DdnLt+w0E8s7T6OiuvMZfMEoYx6ec0eH7zF5EE/XeA0I3n/5oAW1S19kWOcwoITbY9C3gJ+AaIVtXDyT4L9e5vTIazLiljMs92oJiI1AbPdPMiUiGV/bYC15/ZUNU4YDEwCpiYYt8bgS0ZE64x57KEYUwmUc+yvx2B10TkRzwzBt+Syq4L8SyzmtwUPLMNLzlTICJXAqfUu+KaMRnNbqs1JgsSkdnAs6q6w7v9NFBAVV9Mts8TwDFVHe9QmCaXsTEMY7KmAXgGv3d4k8d1eG63Te4onvUfjMkUdoVhjDHGJzaGYYwxxieWMIwxxvjEEoYxxhifWMIwxhjjE0sYxhhjfPJ/ZMu6uD4iMSIAAAAASUVORK5CYII=\n", + "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", "text/plain": [ "
" ] diff --git a/src/gr/gr.f90 b/src/gr/gr.f90 index fd441c03b..8831a93d5 100644 --- a/src/gr/gr.f90 +++ b/src/gr/gr.f90 @@ -101,6 +101,29 @@ module pure subroutine gr_pseudovel2vel(param, mu, xh, pv, vh) return end subroutine gr_pseudovel2vel + module pure subroutine gr_pv2vh_body(self, param) + !! author: David A. Minton + !! + !! Wrapper function that converts from pseudovelocity to heliocentric velocity for swiftest bodies + implicit none + ! Arguments + class(swiftest_body), intent(inout) :: self !! Swiftest particle object + class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters + ! Internals + integer(I4B) :: i + real(DP), dimension(:,:), allocatable :: vh !! Temporary holder of pseudovelocity for in-place conversion + + associate(n => self%nbody) + if (n == 0) return + allocate(vh, mold = self%vh) + do i = 1, n + call gr_pseudovel2vel(param, self%mu(i), self%xh(:, i), self%vh(:, i), vh(:, i)) + end do + call move_alloc(vh, self%vh) + end associate + return + end subroutine gr_pv2vh_body + module pure subroutine gr_vel2pseudovel(param, mu, xh, vh, pv) !! author: David A. Minton !! @@ -177,4 +200,28 @@ module pure subroutine gr_vel2pseudovel(param, mu, xh, vh, pv) return end subroutine gr_vel2pseudovel + module pure subroutine gr_vh2pv_body(self, param) + !! author: David A. Minton + !! + !! Wrapper function that converts from heliocentric velocity to pseudovelocity for Swiftest bodies + implicit none + ! Arguments + class(swiftest_body), intent(inout) :: self !! Swiftest particle object + class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters + ! Internals + integer(I4B) :: i + real(DP), dimension(:,:), allocatable :: pv !! Temporary holder of pseudovelocity for in-place conversion + + associate(n => self%nbody) + if (n == 0) return + allocate(pv, mold = self%vh) + do i = 1, n + call gr_vel2pseudovel(param, self%mu(i), self%xh(:, i), self%vh(:, i), pv(:, i)) + end do + call move_alloc(pv, self%vh) + end associate + return + end subroutine gr_vh2pv_body + + end submodule s_gr \ No newline at end of file diff --git a/src/modules/swiftest_classes.f90 b/src/modules/swiftest_classes.f90 index d6d147634..0b5fdf4c8 100644 --- a/src/modules/swiftest_classes.f90 +++ b/src/modules/swiftest_classes.f90 @@ -159,6 +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 :: 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 @@ -407,6 +409,12 @@ module pure subroutine gr_pseudovel2vel(param, mu, xh, pv, vh) real(DP), dimension(:), intent(out) :: vh !! Heliocentric velocity vector end subroutine gr_pseudovel2vel + module pure subroutine gr_pv2vh_body(self, param) + implicit none + class(swiftest_body), intent(inout) :: self !! Swiftest particle object + class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters + end subroutine gr_pv2vh_body + module pure subroutine gr_vel2pseudovel(param, mu, xh, vh, pv) implicit none class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters @@ -416,6 +424,12 @@ module pure subroutine gr_vel2pseudovel(param, mu, xh, vh, pv) real(DP), dimension(:), intent(out) :: pv !! Pseudovelocity vector - see Saha & Tremain (1994), eq. (32) end subroutine gr_vel2pseudovel + module pure subroutine gr_vh2pv_body(self, param) + implicit none + class(swiftest_body), intent(inout) :: self !! Swiftest particle object + class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters + end subroutine gr_vh2pv_body + module subroutine io_dump_param(self, param_file_name) implicit none class(swiftest_parameters),intent(in) :: self !! Output collection of parameters diff --git a/src/modules/whm_classes.f90 b/src/modules/whm_classes.f90 index b0f176a0a..00d73a471 100644 --- a/src/modules/whm_classes.f90 +++ b/src/modules/whm_classes.f90 @@ -30,21 +30,19 @@ module whm_classes !! Note to developers: If you add componenets to this class, be sure to update methods and subroutines that traverse the !! component list, such as whm_setup_pl and whm_spill_pl contains - procedure, public :: h2j => whm_coord_h2j_pl !! Convert position and velcoity vectors from heliocentric to Jacobi coordinates - procedure, public :: j2h => whm_coord_j2h_pl !! Convert position and velcoity vectors from Jacobi to helliocentric coordinates - procedure, public :: vh2vj => whm_coord_vh2vj_pl !! Convert velocity vectors from heliocentric to Jacobi coordinates - procedure, public :: drift => whm_drift_pl !! Loop through massive bodies and call Danby drift routine - procedure, public :: fill => whm_fill_pl !! "Fills" bodies from one object into another depending on the results of a mask (uses the MERGE intrinsic) - procedure, public :: accel => whm_getacch_pl !! Compute heliocentric accelerations of massive bodies - procedure, public :: accel_gr => whm_gr_getacch_pl !! Acceleration term arising from the post-Newtonian correction - procedure, public :: p4 => whm_gr_p4_pl !! Position kick due to p**4 term in the post-Newtonian correction - procedure, public :: vh2pv => whm_gr_vh2pv_pl !! Converts from heliocentric velocity to psudeovelocity for GR calculations - procedure, public :: pv2vh => whm_gr_pv2vh_pl !! Converts from psudeovelocity to heliocentric velocity for GR calculations - procedure, public :: setup => whm_setup_pl !! Constructor method - Allocates space for number of particles - procedure, public :: set_mu => whm_util_set_mu_eta_pl !! Sets the Jacobi mass value for all massive bodies. - procedure, public :: set_ir3 => whm_setup_set_ir3j !! Sets both the heliocentric and jacobi inverse radius terms (1/rj**3 and 1/rh**3) - procedure, public :: step => whm_step_pl !! Steps the body forward one stepsize - procedure, public :: spill => whm_spill_pl !!"Spills" bodies from one object to another depending on the results of a mask (uses the PACK intrinsic) + procedure, public :: h2j => whm_coord_h2j_pl !! Convert position and velcoity vectors from heliocentric to Jacobi coordinates + procedure, public :: j2h => whm_coord_j2h_pl !! Convert position and velcoity vectors from Jacobi to helliocentric coordinates + procedure, public :: vh2vj => whm_coord_vh2vj_pl !! Convert velocity vectors from heliocentric to Jacobi coordinates + procedure, public :: drift => whm_drift_pl !! Loop through massive bodies and call Danby drift routine + procedure, public :: fill => whm_fill_pl !! "Fills" bodies from one object into another depending on the results of a mask (uses the MERGE intrinsic) + procedure, public :: accel => whm_getacch_pl !! Compute heliocentric accelerations of massive bodies + procedure, public :: accel_gr => whm_gr_getacch_pl !! Acceleration term arising from the post-Newtonian correction + procedure, public :: gr_pos_kick => whm_gr_p4_pl !! Position kick due to p**4 term in the post-Newtonian correction + procedure, public :: setup => whm_setup_pl !! Constructor method - Allocates space for number of particles + procedure, public :: set_mu => whm_util_set_mu_eta_pl !! Sets the Jacobi mass value for all massive bodies. + procedure, public :: set_ir3 => whm_setup_set_ir3j !! Sets both the heliocentric and jacobi inverse radius terms (1/rj**3 and 1/rh**3) + procedure, public :: step => whm_step_pl !! Steps the body forward one stepsize + procedure, public :: spill => whm_spill_pl !!"Spills" bodies from one object to another depending on the results of a mask (uses the PACK intrinsic) end type whm_pl !******************************************************************************************************************************** @@ -57,14 +55,12 @@ module whm_classes !! component list, such as whm_setup_tp and whm_spill_tp contains private - procedure, public :: drift => whm_drift_tp !! Loop through test particles and call Danby drift routine - procedure, public :: accel => whm_getacch_tp !! Compute heliocentric accelerations of test particles - procedure, public :: accel_gr => whm_gr_getacch_tp !! Acceleration term arising from the post-Newtonian correction - procedure, public :: p4 => whm_gr_p4_tp !! Position kick due to p**4 term in the post-Newtonian correction - procedure, public :: vh2pv => whm_gr_vh2pv_tp !! Converts from heliocentric velocity to psudeovelocity for GR calculations - procedure, public :: pv2vh => whm_gr_pv2vh_tp !! Converts from psudeovelocity to heliocentric velocity for GR calculations - procedure, public :: setup => whm_setup_tp !! Allocates new components of the whm class and recursively calls parent allocations - procedure, public :: step => whm_step_tp !! Steps the particle forward one stepsize + procedure, public :: drift => whm_drift_tp !! Loop through test particles and call Danby drift routine + procedure, public :: accel => whm_getacch_tp !! Compute heliocentric accelerations of test particles + procedure, public :: accel_gr => whm_gr_getacch_tp !! Acceleration term arising from the post-Newtonian correction + procedure, public :: gr_pos_kick => whm_gr_p4_tp !! Position kick due to p**4 term in the post-Newtonian correction + procedure, public :: setup => whm_setup_tp !! Allocates new components of the whm class and recursively calls parent allocations + procedure, public :: step => whm_step_tp !! Steps the particle forward one stepsize end type whm_tp !******************************************************************************************************************************** @@ -179,34 +175,6 @@ module pure subroutine whm_gr_p4_tp(self, param, dt) real(DP), intent(in) :: dt !! Step size end subroutine whm_gr_p4_tp - module pure subroutine whm_gr_pv2vh_pl(self, param) - use swiftest_classes, only : swiftest_parameters - implicit none - class(whm_pl), intent(inout) :: self !! Swiftest particle object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters - end subroutine whm_gr_pv2vh_pl - - module pure subroutine whm_gr_pv2vh_tp(self, param) - use swiftest_classes, only : swiftest_parameters - implicit none - class(whm_tp), intent(inout) :: self !! WHM test particle object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters - end subroutine whm_gr_pv2vh_tp - - module pure subroutine whm_gr_vh2pv_pl(self, param) - use swiftest_classes, only : swiftest_parameters - implicit none - class(whm_pl), intent(inout) :: self !! Swiftest particle object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters - end subroutine whm_gr_vh2pv_pl - - module pure subroutine whm_gr_vh2pv_tp(self, param) - use swiftest_classes, only : swiftest_parameters - implicit none - class(whm_tp), intent(inout) :: self !! WHM test particle object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters - end subroutine whm_gr_vh2pv_tp - !> Reads WHM massive body object in from file module subroutine whm_setup_pl(self,n) implicit none diff --git a/src/whm/whm_gr.f90 b/src/whm/whm_gr.f90 index 88018057e..3cf159504 100644 --- a/src/whm/whm_gr.f90 +++ b/src/whm/whm_gr.f90 @@ -17,20 +17,19 @@ module subroutine whm_gr_getacch_pl(self, param) !! author: David A. Minton real(DP), dimension(:, :), allocatable :: aj real(DP) :: beta, rjmag4 - associate(n => self%nbody, mu => self%muj, c2 => param%inv_c2, & - ah => self%ah, xj => self%xj, GMpl => self%Gmass, eta => self%eta) - if (n == 0) return - allocate(aj, mold = ah) - do i = 1, n - rjmag4 = (dot_product(xj(:, i), xj(:, i)))**2 - beta = - mu(i)**2 * c2 - aj(:, i) = 2 * beta * xj(:, i) / rjmag4 + associate(pl => self, npl => self%nbody, inv_c2 => param%inv_c2) + if (npl == 0) return + allocate(aj, mold = pl%ah) + do i = 1, npl + rjmag4 = (dot_product(pl%xj(:, i), pl%xj(:, i)))**2 + beta = -pl%muj(i)**2 * inv_c2 + aj(:, i) = 2 * beta * pl%xj(:, i) / rjmag4 end do suma(:) = 0.0_DP - ah(:, 1) = ah(:, 1) + aj(:, 1) - do i = 2, n - suma(:) = suma(:) + GMpl(i) * aj(:, i) / eta(i) - ah(:, i) = ah(:, i) + aj(:, i) + suma(:) + pl%ah(:, 1) = pl%ah(:, 1) + aj(:, 1) + do i = 2, npl + suma(:) = suma(:) + pl%Gmass(i) * aj(:, i) / pl%eta(i) + pl%ah(:, i) = pl%ah(:, i) + aj(:, i) + suma(:) end do end associate return @@ -51,13 +50,12 @@ module subroutine whm_gr_getacch_tp(self, param) integer(I4B) :: i real(DP) :: rjmag4, beta - associate(n => self%nbody, mu => self%mu,& - c2 => param%inv_c2, ah => self%ah, xh => self%xh, status => self%status) - if (n == 0) return - do i = 1, n - rjmag4 = (dot_product(xh(:, i), xh(:, i)))**2 - beta = - mu(i)**2 * c2 - ah(:, i) = ah(:, i) + beta * xh(:, i) / rjmag4 + associate(tp => self, ntp => self%nbody, inv_c2 => param%inv_c2) + if (ntp == 0) return + do i = 1, ntp + rjmag4 = (dot_product(tp%xh(:, i), tp%xh(:, i)))**2 + beta = - tp%mu(i)**2 * inv_c2 + tp%ah(:, i) = tp%ah(:, i) + beta * tp%xh(:, i) / rjmag4 end do end associate return @@ -113,100 +111,4 @@ module pure subroutine whm_gr_p4_tp(self, param, dt) return end subroutine whm_gr_p4_tp - module pure subroutine whm_gr_pv2vh_pl(self, param) - !! author: David A. Minton - !! - !! Wrapper function that converts from pseudovelocity to heliocentric velocity for massive bodies - !! in a WHM object - implicit none - ! Arguments - class(whm_pl), intent(inout) :: self !! Swiftest particle object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters - ! Internals - integer(I4B) :: i - real(DP), dimension(:,:), allocatable :: vh !! Temporary holder of pseudovelocity for in-place conversion - - associate(pl => self, npl => self%nbody) - if (npl == 0) return - allocate(vh, mold = pl%vh) - do i = 1, npl - call gr_pseudovel2vel(param, pl%mu(i), pl%xh(:, i), pl%vh(:, i), vh(:, i)) - end do - call move_alloc(vh, pl%vh) - end associate - return - end subroutine whm_gr_pv2vh_pl - - module pure subroutine whm_gr_pv2vh_tp(self, param) - !! author: David A. Minton - !! - !! Wrapper function that converts from pseudovelocity to heliocentric velocity for test particles bodies - !! in a WHM object - implicit none - ! Arguments - class(whm_tp), intent(inout) :: self !! Swiftest particle object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters - ! Internals - integer(I4B) :: i - real(DP), dimension(:,:), allocatable :: vh !! Temporary holder of pseudovelocity for in-place conversion - - associate(tp => self, ntp => self%nbody) - if (ntp == 0) return - allocate(vh, mold = tp%vh) - do i = 1, ntp - call gr_pseudovel2vel(param, tp%mu(i), tp%xh(:, i), tp%vh(:, i), vh(:, i)) - end do - call move_alloc(vh, tp%vh) - end associate - return - end subroutine whm_gr_pv2vh_tp - - module pure subroutine whm_gr_vh2pv_pl(self, param) - !! author: David A. Minton - !! - !! Wrapper function that converts from heliocentric velocity to pseudovelocity for massive bodies - !! in a WHM object - implicit none - ! Arguments - class(whm_pl), intent(inout) :: self !! Swiftest particle object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters - ! Internals - integer(I4B) :: i - real(DP), dimension(:,:), allocatable :: pv !! Temporary holder of pseudovelocity for in-place conversion - - associate(pl => self, npl => self%nbody) - if (npl == 0) return - allocate(pv, mold = pl%vh) - do i = 1, npl - call gr_vel2pseudovel(param, pl%mu(i), pl%xh(:, i), pl%vh(:, i), pv(:, i)) - end do - call move_alloc(pv, pl%vh) - end associate - return - end subroutine whm_gr_vh2pv_pl - - module pure subroutine whm_gr_vh2pv_tp(self, param) - !! author: David A. Minton - !! - !! Wrapper function that converts from heliocentric velocity to pseudovelocity for teset particles - !! in a WHM object - implicit none - ! Arguments - class(whm_tp), intent(inout) :: self !! Swiftest particle object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters of on parameters - ! Internals - integer(I4B) :: i - real(DP), dimension(:,:), allocatable :: pv !! Temporary holder of pseudovelocity for in-place conversion - - associate(tp => self, ntp => self%nbody) - if (ntp == 0) return - allocate(pv, mold = tp%vh) - do i = 1, ntp - call gr_vel2pseudovel(param, tp%mu(i), tp%xh(:, i), tp%vh(:, i), pv(:, i)) - end do - call move_alloc(pv, tp%vh) - end associate - return - end subroutine whm_gr_vh2pv_tp - end submodule s_whm_gr \ No newline at end of file diff --git a/src/whm/whm_step.f90 b/src/whm/whm_step.f90 index fb84fe49e..aac704f67 100644 --- a/src/whm/whm_step.f90 +++ b/src/whm/whm_step.f90 @@ -54,9 +54,9 @@ module subroutine whm_step_pl(self, system, param, t, dt) call pl%set_beg_end(xbeg = pl%xh) call pl%kick(dth) call pl%vh2vj(cb) - if (param%lgr) call pl%p4(param, dth) + if (param%lgr) call pl%gr_pos_kick(param, dth) call pl%drift(system, param, dt) - if (param%lgr) call pl%p4(param, dth) + if (param%lgr) call pl%gr_pos_kick(param, dth) call pl%j2h(cb) call pl%accel(system, param, t + dt) call pl%kick(dth) @@ -93,9 +93,9 @@ module subroutine whm_step_tp(self, system, param, t, dt) tp%lfirst = .false. end if call tp%kick(dth) - if (param%lgr) call tp%p4(param, dth) + if (param%lgr) call tp%gr_pos_kick(param, dth) call tp%drift(system, param, dt) - if (param%lgr) call tp%p4(param, dth) + if (param%lgr) call tp%gr_pos_kick(param, dth) call tp%accel(system, param, t + dt, lbeg=.false.) call tp%kick(dth) end associate