From 19757e6714f59fbe79a77ed1661ceae6dabab538 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Mon, 2 Sep 2024 13:42:08 -0400 Subject: [PATCH] add an example of comparing the compressible solvers (#236) also fix __repr__ of Pyro --- docs/source/compressible-rt-compare.ipynb | 198 ++++++++++++++++++++++ docs/source/index.rst | 1 + pyro/pyro_sim.py | 5 +- 3 files changed, 203 insertions(+), 1 deletion(-) create mode 100644 docs/source/compressible-rt-compare.ipynb diff --git a/docs/source/compressible-rt-compare.ipynb b/docs/source/compressible-rt-compare.ipynb new file mode 100644 index 000000000..2137afee7 --- /dev/null +++ b/docs/source/compressible-rt-compare.ipynb @@ -0,0 +1,198 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "5132ed09-12eb-4c88-9356-eb4ae6425001", + "metadata": {}, + "source": [ + "# Comparing the Compressible Solvers" + ] + }, + { + "cell_type": "markdown", + "id": "c6fb6238-510a-451a-ae84-8ab8e132dc5d", + "metadata": {}, + "source": [ + "Here we'll run the same problem (the single-mode Rayleigh-Taylor instability)\n", + "with 3 different compressible solvers." + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "id": "2151f491-17f2-4d7f-bc75-1f4f2fa93bdd", + "metadata": {}, + "outputs": [], + "source": [ + "from pyro import Pyro" + ] + }, + { + "cell_type": "markdown", + "id": "818bcb7a-92b1-4580-b3f4-d5862ee26006", + "metadata": {}, + "source": [ + "We'll append the `Pyro` objects to the `runs` list as we do the simulations." + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "id": "2940380a-67f3-4dff-8bd8-af8079464d30", + "metadata": {}, + "outputs": [], + "source": [ + "runs = []\n", + "solvers = [\"compressible\", \"compressible_rk\", \"compressible_fv4\"]" + ] + }, + { + "cell_type": "markdown", + "id": "fb6e7f9f-a693-4776-874d-fd47c5e17986", + "metadata": {}, + "source": [ + "We'll override the default number of cells for the `rt` problem to make it smaller (and run faster)." + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "id": "51f5f2d9-2e42-4257-aa1f-b7e6e4079393", + "metadata": {}, + "outputs": [], + "source": [ + "problem = \"rt\"\n", + "params = {\"mesh.nx\": 32, \"mesh.ny\": 96, \"driver.cfl\": 0.8}" + ] + }, + { + "cell_type": "markdown", + "id": "11247f6f-7b2f-4599-8281-c36f526bf01c", + "metadata": {}, + "source": [ + "Now we loop over the solvers, setup the `Pyro` object and run the simulations." + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "id": "78825cbd-1643-471b-b6fc-1efc77d833df", + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/home/zingale/development/pyro2/pyro/compressible/problems/rt.py:72: RuntimeWarning: invalid value encountered in divide\n", + " 0.5*(xmom[:, :]**2 + ymom[:, :]**2)/dens[:, :]\n", + "/home/zingale/development/pyro2/pyro/compressible_rk/problems/rt.py:72: RuntimeWarning: invalid value encountered in divide\n", + " 0.5*(xmom[:, :]**2 + ymom[:, :]**2)/dens[:, :]\n", + "/home/zingale/development/pyro2/pyro/compressible_fv4/problems/rt.py:72: RuntimeWarning: invalid value encountered in divide\n", + " 0.5*(xmom[:, :]**2 + ymom[:, :]**2)/dens[:, :]\n" + ] + } + ], + "source": [ + "for s in solvers:\n", + " p = Pyro(s)\n", + " p.initialize_problem(problem_name=problem, inputs_dict=params)\n", + " p.run_sim()\n", + " runs.append(p)" + ] + }, + { + "cell_type": "markdown", + "id": "62349272-56b3-4a6c-9012-c00bf7e9a9db", + "metadata": {}, + "source": [ + "Finally, we'll plot them. We'll just look at the density and put them all on the same set of axes.\n", + "Since our data is stored as $\\rho(x, y)$, `imshow()` will want to put $x$ on the vertical axis, so\n", + "we transpose the data before plotting." + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "id": "9afb841c-f08f-4554-8fcb-06939d0fa131", + "metadata": {}, + "outputs": [], + "source": [ + "import matplotlib.pyplot as plt\n", + "from mpl_toolkits.axes_grid1 import ImageGrid" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "id": "e208db0e-267e-4499-85e8-1e9f915d2a51", + "metadata": { + "scrolled": true + }, + "outputs": [ + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 7, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAm8AAAIdCAYAAABr3JpCAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy80BEi2AAAACXBIWXMAABJ0AAASdAHeZh94AACBdElEQVR4nO39eXxURb4//r+6O/sGISaEAAlLlD1DAgq4sBkRcOaKl4FhZITMgOgICpdrRNCvIHjdUD7gcOM4IzJwxxkUF0DgZ4YIgUEEHRYRRSRgEgYQsgAJCYR0+vz+gHSfqu70vp3068kjj0dXn9N1qs8puqvrvKtKpyiKAiIiIiLSBH2gC0BEREREzmPjjYiIiEhD2HgjIiIi0hA23oiIiIg0hI03IiIiIg1h442IiIhIQ9h4IyIiItIQNt6IiIiINISNNyIiIiINYeONiIiISEPYeCMiIiLSEDbeKKj85S9/gU6nQ3Fxsd3n7Bk+fDi6dOniw1KSVrWW+qXT6ZCXlxfQMpB7tFQHDx06hLvvvhuJiYnQ6XRYtGiRz49JzgkLdAGIiIgouBiNRowfPx6NjY1YsmQJ2rZti6ysLI/yPHz4MAYMGACj0Yj169fjl7/8pdfKG2rYeKOg8tBDD2HSpEmIiIgIdFGoFWL9okDTSh08efIkTp48iddffx2zZs3yOD+TyYSHH34YUVFRuHz5slfKGMp425QEjY2NuHr1asCObzAYEBUVBb2eVbM1Yv1yX6DPXWsR6POolTr4008/AQDatWvnlfz+8Ic/4Ntvv0V+fr5X8gt1wV17gti1a9fw6quvon///oiJiUGbNm0wcOBArFy5UtivtLQUDz30ENq3b4/IyEh0794dCxYsQH19vbDfokWLoNPp8N1332HOnDno0KEDYmJicPfdd+PYsWMAgI8++gg5OTmIjo5Gly5d8Kc//cmqXM2xMEVFRRg8eDBiYmKQmpqK2bNnW/3aaT7mt99+i7lz56JTp06IiorC3r17AQANDQ148cUX0adPH0RFRaFt27b4xS9+gYMHDwr5mEwmLF++HFlZWYiPj0dCQgJ69OiBadOmobGx0bzfnj17MGbMGKSmpiIqKgodO3bE2LFjzceDg9gPo9GIRYsWISMjA5GRkcjKysK6deucvmbHjx/HQw89hA4dOiAiIgJdunRBfn4+6urqnM7DX1i/LEKxfjk6d7YcOHAAqamp6N27N8rLy10+pox10CLU6uDw4cMxbNgwAMBvf/tb6HQ66HQ6HD16FDqdDnPnzrX5ul//+teIiIhARUWF8PypU6fw7LPPYtGiRUhPT3epLGQbb5u64dq1a7j33ntRXFyMUaNG4Te/+Q2ioqLwzTff4KOPPjJ3MZeVleG2227DpUuX8Nhjj+Hmm29GcXExXnrpJXz++ef47LPPEBYmXoKpU6ciLi4OCxYsQEVFBV5//XXce++9WLJkCZ566in8/ve/x+9+9zusWrUKjzzyCHr37o0777xTyOPAgQP44IMP8PDDD2PKlCnYsWMH3njjDRw5cgTbtm2z+sU3efJkREdH47//+7+h0+nQoUMHNDY2YvTo0dizZw8eeughzJo1C5cuXcKf//xn3HHHHdi1axcGDhwIAPif//kfPPfcc/jFL36BRx99FAaDAT/++CM2bdqEhoYGhIeH49ixY7jnnnvMH7Lt27fHuXPnsHv3bnz99dcYPHiww/M+b9481NXV4bHHHgMArF69Gr/+9a9x9epVh8Hb+/fvx8iRI9G2bVs88sgj6NixI77++mu88cYb+Pzzz7Fz506Eh4c7WQN8i/WL9cveubOlsLAQv/zlL5GVlYVPPvnE494S1sHQroPPPPMM7rjjDrz44ouYMWMG7rrrLgBAeno6br31Vvztb3/D0qVLYTAYzK+pqanBxo0bMWbMGCQnJwv5/f73v0e3bt0wZ84c/PWvf3WqDOSAQi575ZVXFADK/PnzrbY1NTWZHz/44IMKAGXLli3CPk8++aQCQHn77bfNzy1cuFABoPz85z9XTCaT+fkVK1YoAJT4+HilvLzc/Pz58+eVyMhIZdKkSULeABQAyscffyw8/8QTTygAlL///e9Wxxw2bJjS2Ngo7L9s2TIFgPLpp58Kz1+6dEnp3LmzMmzYMPNz2dnZSq9eveyes+b3sW/fPrv7rV69WgGg7Nixw+q59PR05eLFi+bnL168qKSnpyuJiYlKfX29+flhw4YpGRkZQr5ZWVlKjx49lJqaGuH5jz76SAGgrF692m65/In1i/XL3rlTblyHqVOnKoqiKGvXrlXCw8OV+++/XyinJ1gHWQd37Nhh83UrV660ec3ffvttBYDy4YcfCs+vW7dO0el0yp49e4T3un79epfKQyLeNnXDu+++i8TERDz33HNW25p/8ZlMJmzatAnZ2dkYO3assM/8+fOh1+vx8ccfW73+iSeegE6nM6ebf/H8x3/8Bzp37mx+Pjk5GT169MDx48et8ujRowfGjRsnPPf0008DgM1jzpkzx+rX8V//+lf07NkTAwYMQGVlpfnv2rVruOeee7B7925cuXIFANCmTRucPn0au3fvbvGctWnTBgCwceNGt+NNfv/735vzac7z0UcfxYULF+wOsf/mm29w+PBhPPjgg2hoaBDez5133onY2Fj84x//cKtMvsD6xfrVzNa5U3v55ZcxdepU/O53v8OHH36I6Ohot44jYx1kHWxJ863RtWvXCs+vXbsW7dq1w89//nPzcxcuXMDs2bPx8MMPY8iQIV45Pl3Hxpsbjh8/jp49eyIqKqrFfSoqKnD58mX06dPHalu7du3QoUMHnDx50mpbt27dhHRiYiIAoGvXrlb7JiYmoqqqyur5Xr16WT3XoUMHtG3b1uYxb7nlFqvnjh49iu+//x7JyclWf++88w6amppQWVkJAHjxxRcRFRWFu+66Cx07dsTkyZPxt7/9DdeuXTPnN2nSJOTm5uLFF19Eu3btMHLkSLzyyisoKyuzOnZLbL2v3r17AzdGRrXk6NGjAICFCxdavZeUlBTU1dXh3LlzTpfD11i/WL+a2Tp3zT766CPMnz8f06dPxx//+EfhFpanWAdZB1vS3EDbuHEjampqgBtxj//85z+tRtHm5+dDURS8/PLLXjk2WTDmLci09AHc0vPX7yJ4JiYmxma+/fr1w7Jly1p8XXNcw5AhQ3DixAkUFhZix44d2LFjB/72t7/hhRdewO7du9GuXTtERkZi27Zt+PLLL1FYWIhdu3bhueeew6JFi/C3v/0NDzzwgMfvoyXN5+i///u/MXr0aJv7NH+BtHasX97ny/pl69w1u+2221BaWooPPvgAM2bMMMdnBTvWQe/z92fclClT8NFHH+H999/H9OnT8X//939QFAVTp04173PgwAG88847eP7551FVVWVuhJ8/fx64MZq1pKQEnTt3RmRkpNfKFirYeHPDLbfcgu+//x4NDQ0tVrrk5GTEx8fj22+/tdp24cIFnD17Fv379/dJ+Zp/hamdPXsWFy9etPrV25Kbb74ZFRUVGDlypFND2uPi4jB+/HiMHz8eAFBQUICZM2di1apVwtDw2267DbfddhtwYwRSdnY2nn32Wac+2I4ePYr7779feO67774DbPyal98Lbnw55ObmOjxOoLF+WWP9stapUyesWbMGI0eORG5uLj799FOnguKdwTpojXXQYuzYsbjpppuwdu1ac+OtZ8+e5vcNAOXl5VAUBc8995zN2++PP/44AOCrr77SzA+PYMLbpm6YPHkyLly4gBdeeMFqW/MvIL1ebx5y/umnnwr7vPzyyzCZTD77JXbs2DFs2LBBeO6VV14BAKs4kZZMmTIFP/30U4u/StVd8M23FtRycnIAANXV1S3u06lTJyQnJ5v3ceTNN9/EpUuXzOlLly7hj3/8I9q2bWse1m5LdnY2+vbtiz/+8Y82bz0YjUany+APrF+sX87q2LEjdu7cibS0NIwaNQqff/65V/JlHWQdtCc8PBwPPvggdu/ejb/97W84fvy40OuGG43Y9evXW/3NnDkTuNFLuH79enTv3t1r5Qol7Hlzw+zZs/HJJ5/ghRdewFdffYVRo0YhKioK3377LY4dO4aioiLgRpzEtm3bMG7cODz22GPIzMzErl278N5772Ho0KFWld1b+vXrh9/85jd4+OGHcfPNN2PHjh344IMPMGzYMPzqV79y+j1u27YN+fn52L59O0aOHImEhASUl5fjs88+Q1RUFHbs2AHciNMYPHgwBg0ahLS0NJw9exZ/+tOfEBERgUmTJgEAXnjhBfzjH//Az3/+c3Tt2hWKouCTTz7B999/j6eeesqpMt10000YNGgQfvvb3wI3htGXl5fj7bfftnt7SafT4f/+7/8wcuRIZGVl4Xe/+x369OmD+vp6lJSU4KOPPsJLL70UNGtFsn6xfrkiNTUVxcXFyM3NxejRo7F582a7X/TOYB1kHXRk6tSpeOONN/D73/8eer0ev/nNb4TtaWlpNpe/ap6Lb/DgwVweyxOBHu6qVVeuXFFeeOEFpXfv3kpkZKTSpk0bZeDAgcr//u//CvudPHlS+c1vfqMkJycr4eHhSteuXZX58+crdXV1wn7NQ9p//PFH4fkff/xRAaAsXLjQqgy2hos3TyGwbds25bbbblOioqKUlJQUZdasWVZDyFs6ZrPGxkZlxYoVysCBA5WYmBglJiZGyczMVB588EGlsLDQvN9LL72k3HXXXUpycrISERGhdOrUSfnlL3+p7N+/37zPjh07lIkTJyoZGRlKVFSUkpiYqNx2223Kn//8Z2HaAHvD6Ldt26Y899xzSufOnZWIiAilb9++yrvvvuvUeVEURSktLVUeeeQRJSMjQwkPD1fatWun5OTkKE8//bQwRUEwYP0K7frl6NyppwppVllZqfTv31+JiYlRioqKXDqeLayDoV0HW5oqRK1v374KACU3N9fpfDlViHfoFG9Eg1LQ0Ol0mDp1Kv7yl78EuijUCrF+UaCxDhIx5o2IiIhIUzyKefvDH/6AnTt34ptvvsH58+dx+fJlJCUlYeDAgXj00UeFyfqcZTQasXz5cqxduxYlJSWIjo7GoEGD8PTTT2Po0KGeFJeIKGCuXLkiBKO3JDU11S/lodDDOth6eNR4e+WVV1BRUYF+/fqhR48eiIqKwsmTJ7FlyxZs2bIF//Vf/2V3Dh1ZY2Mjxo4di6KiIiQlJeG+++5DVVUVCgsLUVhYiNWrV2PKlCmeFJmIKCDee+89cyC6PYxkIV9hHWw9PIp5+/zzz5GdnW01Cmb37t0YPXo06urqsGfPHqeXxXjxxRfxzDPPIDs7G5999pl5UsGioiKMGTMGBoMBP/zwA9LT090tMhFRQJw9e9bmnGgyLcxFSNrEOth6+GzAwvTp07Fq1So8//zzNifokxmNRqSmpqKqqgpffPGF1WSTjz76KN566y3MnTsXr7/+ui+KTERERBT0fDZgoXkRYGeXvdizZw+qqqrQpUsXm7OEN8+ls3HjRi+XlIiIiEg7fNJ4O3ToEN577z3o9XqMGTPGqdccPHgQADBgwACb25ufP3HiBGpra71YWiIiIiLt8MoKC2+++Sb27duHhoYGlJWVYe/evQgPD0dBQQGysrKcyqOsrAwA0LlzZ5vb4+PjkZCQgJqaGpSVlaFv374t5nX+/HlUVFQIz9XU1OCHH35Av379uAgu+UxDQwNOnTqFYcOGoW3btsI21ksKFNZLCkb26mVLzpw5gwsXLtjclpiYiLS0NC+XMkh5Y6bfX/3qVwoA819MTIzy1ltvKUaj0ek8Hn74YQWA8swzz7S4T1pamgJA2bNnj928mmfV5h//AvW3YcMG1kv+Bd0f6yX/gvHPVr205fTp00p8nL7FfOLj45XTp087lZfWeaXnbd26dVi3bh3q6upw/PhxrFixAo888gg++OADbNiwwe6abL7w2GOPYcKECcJz3333HSZOnIjspF8gJty5Fr5maHWqZVOgC+B99Y0XcbDqE5s9yM7XS8WHJdT5MG8HgrmeBrQu+v56u10vk1vh5yUFjfrGizhYYbte2nLhwgXUXjbho9UdkNk1XNhW8mMj/vO3Z3HhwoWQ6H3z6sL0sbGx6N+/P1avXg29Xo933nkHS5cuxcKFCx2+Ni4uDgBQV1fX4j7NC9rGx8fbzSslJQUpKSk2t8WEt0V8RJLD8miKPoBfyJ4w+fJLK7Bs3Wpyul76co4lXSAbb0FcTwNZF/14vd2rlzf5rHhEcGFgY7NuXcPQq4fYeDP59EdQ8PHZb+GpU6cCLowOzcjIAACcOnXK5vba2lrU1NQI+xIRERGFGq/2vKklJycDgFUgbEuys7MBAPv377e5vfn5bt26Oex5IyIiotbJpJjQpJisngslPut5Ky4uBgBkZmY6tf/tt9+OpKQklJaWYu/evVbb161bBwAYN26cl0tKREREpB1uN94+//xzbN26FU1NTVbbPvnkEzzzzDPAjZUW1KZMmYKePXti5cqVwvNhYWGYO3cuAGDmzJm4ePGieVtRURFWrVqFyMhIzJ49290iExERkcaZoNj8CyVu3zY9fvw4fvvb36Jdu3bIyclBSkoKLl68iGPHjuHEiRMAgDlz5mDy5MnC68rLy3Hs2DFUVlZa5Zmfn48dO3agqKgImZmZGDFiBKqrq1FcXAxFUfD2229zXVMiIiIKaW433oYNG4Znn30Wu3btwtGjR/HPf/4Ter0eaWlpeOihh/Dwww/jrrvucinP8PBwbN26FcuXL8eaNWuwefNmREVFYdSoUZg/fz6GDh3qbnGJiIioFTDd+Cc/F0rcbrx17doVS5Yscfl1zbFwLQkPD0d+fj7y8/PdLRoRERG1UiYFaJKm2GnFM0/ZFMzTZhIRERGRxGdThRARERF5m60BCqE2YIE9b0REREQawp43IiIi0gwTFDSx542IiIiItII9b0RERKQZjHljzxsRERGRprDnjYiIiDTDpCg25nkLrZ43Nt6IiIhIM0w3/uTnQglvmxIRERFpCHveiIiISDOabEwVIqdbO/a8EREREWkIe96IiIhIM64vTG/9XChhzxsRERGRhrDnjYiIiDSDo03Z80ZERESkKex5Iwo2Op2YDrHJJ8kDct3RwnFZv4NTEF9TE3Rogs7quVDCxhsRERFphkmxHqDAAQtEREREFLTY80ZERESawdumbLwRtW6BioEi5zHGkYhcxMYbERF5Rm5w2vvRwMapNgTxNW2y0fMmp1s7xrwRERFRyGpsbMS2bdswZ84c9O/fH3FxcYiMjETXrl0xbdo0HD161K18jUYjXnvtNWRlZSEmJgZJSUkYO3Ysdu3a5XGZ2XgjIiIizTBBB5Mi/XnQ87Zz506MGjUKK1asQFVVFXJzc/Hzn/8cAPDOO+8gOzsbGzdudCnPxsZGjBkzBvn5+Thz5gzuu+8+/OxnP0NhYSFGjBiBtWvXul1esPFGREREoUyv12PChAn44osvcOrUKWzYsAEffvghSkpKMG/ePDQ0NGDq1KmoqqpyOs+lS5eiqKgI2dnZOH78ONavX4/t27ejsLAQer0eM2bMQHl5uftldvuVROQfOp34R9QsWOuGorT8pzXNk4p5+09rguiaNo82Vf950vM2cuRIvP/++xg8eLDwvMFgwEsvvYQePXrg0qVL2LJli1P5GY1GLFu2DABQUFCAxMRE87bc3FxMmzYNDQ0NWLFihdtlZuONiIiINON6g00v/fnmx4tOp0NWVhYA4PTp0069Zs+ePaiqqkKXLl2sGoQAMGnSJABw+VasGhtvRERERC0oKSkBAKSmpjq1/8GDBwEAAwYMsLm9+fkTJ06gtrbWrTJxqhAiIiLSDOXGIAX5OagaWmrJyclISUlx61hFRUU4ePAgIiMjMXr0aKdeU1ZWBgDo3Lmzze3x8fFISEhATU0NysrK0LdvX5fLxcYbkdaoY5tcmYuJtMHepL28vr7nr3g0+Th6XltvGDdunNVzCxcuxKJFi1zOq7KyEtOmTQMAPPnkk+jQoYNTr7t8+TIAIDY2tsV94uLiUFNTw563kNJa/pPL70OLQbxERORX9ibp3bBhAzIzM4VtycnJLh/j6tWrGD9+PMrLyzF06FAsXLjQw1J7FxtvRERE1CpkZmaiT58+HuVhNBoxceJE7Nq1Czk5Odi0aRPCw8Odfn1cXBwAoK6ursV9mnvn4uPj3SojG29ERESkGSZFjyZFb/WcNzQ1NWHy5Mn45JNP0KtXLxQWFqJNmzYu5ZGRkQEAOHXqlM3ttbW1qKmpEfZ1FUebEhERUchTFAXTpk3D+++/j+7du6OoqAg33XSTy/lkZ2cDAPbv329ze/Pz3bp1c7vnjY03Ii0L1klayXt4fb0vGCfPDcYyBSkTdDBBL/15/v9j1qxZWLNmDdLT07F9+3akpaW5lc/tt9+OpKQklJaWYu/evVbb161bB7QwuMJZbLwRERGRZnh7hQUAeOqpp1BQUIC0tDRs374d6enpDl8zZcoU9OzZEytXrhSeDwsLw9y5cwEAM2fOxMWLF83bioqKsGrVKkRGRmL27Nlul5cxb0RERBSyNm3ahKVLlwI3bmUuWbLE5n533nknpk+fbk6Xl5fj2LFjqKystNo3Pz8fO3bsQFFRETIzMzFixAhUV1ejuLgYiqLg7bffdqqB2BI23oiIiEgzmmwMWJDTrqiurjY/3r17N3bv3t3ivurGmz3h4eHYunUrli9fjjVr1mDz5s2IiorCqFGjMH/+fAwdOtTt8oKNNyIiIgpleXl5yMvLc/l1xcXFdreHh4cjPz8f+fn5HpTONjbeiIiISDNMNmLcvDFgQUs4YIGIiIhIQ9jzRkRERJpxfbSp3uq5UMKeNyIiIiINYc8bERERaYYvl8fSCjbeiIiISDOaV1WQnwslofVuiYiIiDSOPW9ERESkGU0K0KTorJ4LJWy8EZFv6DU0+ksuKxcDb92C/Xpr6f8OBQQbb0RERKQZJuhtTBUSWlFgofVuiYiIiDSOPW9ERESkGSboraYGCbWeNzbeiIgotKljzOT4N1/Fn/nrONQqsfFGREREmsHlsdh4IyIiIg1pUnQ2pgoJrcZbaN0kJiIiItI49rwRERE181fsGWPc3KbYWB5LCbG+qNB6t0REREQax543LQiVX2jBPus5EREFXJOiR5M0VYicbu1C690SERERaRx73oiIiEgzTDamBjEFrDSB4XbPW2NjI7Zt24Y5c+agf//+iIuLQ2RkJLp27Ypp06bh6NGjLueZl5cHnU7X4t/o0aPdLS4RERFRq+B2z9vOnTsxatQoAECnTp2Qm5sLg8GAAwcO4J133sG7776L9957D/fff7/Led9xxx3IzMy0er5fv37uFpeIiIhaAZONmDd5uazWzu3Gm16vx4QJEzB37lwMHjzY/HxTUxOeeeYZvPLKK5g6dSpOnDiBpKQkl/KePn068vLy3C0aERERUavldlN15MiReP/994WGGwAYDAa89NJL6NGjBy5duoQtW7Z4o5xEREREaILe5l8o8cm71el0yMrKAgCcPn3aF4cgIiKiEKQoOpikPyXElsfy2WjTkpISAEBqaqrLr92xYwe+/vpr1NfXIzU1FcOHD8eIESN8UEoiIiIibfFJ462oqAgHDx5EZGSkWyNE165dK6QXL16MQYMG4b333kNGRobD158/fx4VFRXCc82NSaJAYb2kYMR6SVrTBJ3VbdImsOfNI5WVlZg2bRoA4Mknn0SHDh2cfm3//v1x66234u6770Z6ejqqq6vxxRdfYMGCBdi3bx9yc3Nx8OBBxMXF2c2noKAAzz//vMfvhcibWC8pGLFeEmmPVxtvV69exfjx41FeXo6hQ4di4cKFLr1+zpw5QjomJgYTJkzAvffei5ycHJSUlODNN99Efn6+3Xwee+wxTJgwQXiupKQE48aNc6k8RN7EeknBiPWStEZR9FZTgyicKsQ9RqMREydOxK5du5CTk4NNmzYhPDzcK3knJCRg9uzZeOKJJ7B161aHjbeUlBSkpKR45dhE3sJ6ScGI9ZJIe7zSeGtqasLkyZPxySefoFevXigsLESbNm28kbVZjx49AABnzpzxar5ERESkHddj3nRWz4USj/sZFUXBtGnT8P7776N79+4oKirCTTfd5J3SqVRXVwOAw3g3IiIiotbM4563WbNmYc2aNUhPT8f27duRlpbmnZJJ1q9fDwAYOHCgT/InIiKi4GeCzirmTV6ovrXzqOftqaeeQkFBAdLS0rB9+3akp6c7fM2UKVPQs2dPrFy5Unj+0KFD2LJlC5qamoTnr1y5ggULFuCjjz6CwWDAzJkzPSkyERERaZjpxm1T9V+oNd7c7nnbtGkTli5dCgDo1q0blixZYnO/O++8E9OnTzeny8vLcezYMVRWVgr7lZaW4oEHHkBSUhJycnKQnJyMiooKHDp0CBUVFYiIiMBbb71lXrmBiIiIKBS53XhrjkEDgN27d2P37t0t7qtuvLUkKysLTzzxBL766iscOXIEVVVVMBgMSE9Px/jx4/H444+jd+/e7haXiIiIWgGTjalC5HRr53bjLS8vD3l5eS6/rri42Obz3bp1w4oVK9wtDhEREVFI8NnapkRERETeZlL0aArxnrfQerdEREREGseeNyIiItIMk42pQUwBK01gsOeNiIiISEPY80ZERESawZg3Nt6Igp+iuP9aXWhNXNkq8HoT2WVSdDApOqvnQkloNVWJiIiINI49b0RERKQZTdCjSep7ktOtXWi9WyIiIiKNY88bERF5l0kVt6cPrVikkBHAa6wo1jFunoSKahEbb0TBxpufQnJeDGgPPrzeROQiNt6IiIhIM0zQwyRFfcnp1i603i0RERGRxrHnjYiIvItxbq1fAK+xCTo0yfO8IbTqHBtvREFB8U/ErS9jolrTF7b8Xkw+jEvzJfOxQiyam1o1TtLL26ZEREREmsKeNyIiItKM6z1v8tqm7HkjIiIioiDFnjciIvKMJzGBrSlWsjUJ4mtqgg5N4IAFIqLgZ29wRahNr05EIY2Nt1Dn6WhDfmkSEZEfcbQpY96IiIiINIWNNyIiItIMk6K3+eeJAwcO4NVXX8XEiRPRtWtX6HQ66HQ6HDlyxO08T5w4gRkzZqB79+6IjIxETEwMevfujfz8fFRUVHhUXt42JQplwbSQuSfHdvTaQN7eZ2gBkVcpNgYoePq/bPHixdi4caOHuVjs2bMHo0aNQl1dHbp164af//znaGhowL59+/Daa6/h3Xffxe7du9GtWze38mfPm5qiuP8XzHS6lv+COW9faq3XmoiIXDZkyBA8++yz+Pjjj3Hq1ClkZGR4lN8jjzyCuro6PPXUUzh+/Dg+/PBDbN68GaWlpRg9ejTOnj2Lp59+2u382fNGREREmtGkWK9tKqddNW/ePA9LZVFVVYUjR47AYDBg0aJF0Ost/WSxsbF47rnn8Omnn2Lv3r1uH4ONNyIico0313pV58U53wLLW9dVzifErmtkZCQAmOPmWpKUlOT2MXjblIiIiDRDgfVgBSWImjNxcXG44447YDQasXDhQphMJvO2uro6LF68GAAwbdo0t48R2j1v3oxfCqbAb7LGax18/Hne5GMxdpGoVSopKbF6Ljk5GSkpKX4tx5///GeMHj0ar776Kj744ANkZ2ejoaEBe/fuhdFoxEsvvYRZs2a5nX9oN96IiIhIU+xN0jtu3Dir/RcuXIhFixb5rXwA0KtXL+zZswe/+tWv8Pnnn+PkyZPmbbm5uRg6dKhH+bPxRkRERK3Chg0bkJmZKTyXnJzs93IUFxdj/PjxaN++PT799FMMGjQI9fX12Lx5M5566ikMHz4cH374IX7xi1+4lT8bb0RE5Jg3Byk4e4wQC3T3O39cU/k4XrimJuis5nlrTmdmZqJPnz4eH8MT1dXVGD9+PBoaGvDpp58iPT0dANC2bVvMmDEDbdu2xa9+9Ss8/vjjGDNmDMLCXG+KBU+EHxEFHue08xzPIZFPKYrl1mnzXzD9V9uyZQuqq6sxePBgc8NN7T//8z8RERGBsrIy4XaqK0Kw581PH6iBDGoPlgD6QAaJB+xYQfQJQkREfvfvf/8bAJCQkGBze1hYGGJjY3Ht2jVcuHDBrWOEYOONiIiItMrWWqaerm3qTR06dAAAHDx4EEaj0eq26A8//GButHXp0sWtYwTPuyUiouBiuhGv5K/YKKvjK+JfINhbAtCTv0AIhvMplCEwRfCWKVOmoGfPnli5cqXw/JgxYxAdHY3S0lLMmzcPRqPRvK2yshIPP/wwAGDYsGFo3769W8dmzxsRERFphr2pQty1ZcsWLFmyxJw+e/YsAGDy5MmIjo4GAOTk5KCgoMC8T3l5OY4dO4bKykohr/bt2+ONN97AI488gmXLluGDDz5ATk4Orly5gn379uHixYto3749/vSnP7ldXjbe/MWXMXCu5OXrX3z2Ys18GQMXTNGqZJurdc+TuuqoPnDSXiJSqaiowL59+6yeP3z4sPlxVFSU0/lNnz4d/fr1w/Lly/H5559jy5YtCAsLQ9euXTF9+nTk5+d7NHEwG29ERESkGYqNqUIUeNYxkZeXh7y8PJdeU1xcbHf7oEGD8Pe//92jcrWEjTciItIGfyxi7694NH/1/gYqto18io03IiIi0gxfxLxpDRtvRMHGl/FX3uxVcNTz4c8YN0d5eRoDp36v3u7J0Mr1JgoSJhuNNY0PXHUZG2+hIFATBDMInIiIyOvYeCMiIiLNUGzcNlV425SIiCiEBMPtZd61IBew8UZERESawQELod54C2SgsPrY3g7strdd2qYYPFshTdckhYnae1+OJip25Xq4eu2CPShcUfzza9uXk0W7SAkziE9EhIvbI8SPJ8XQcll1TeL70l0zijtcaxS3G5tcLK0X+bNXpflY7MkhalVCu/FGREREmmKyMUmvnG7t2HgjIiLvcmWZPHfJU7a4MmlvMMS42ePJXQlfTcrrj2tKTmPjjYiIiDSDo03ZeCMKba7GIboSTxkdKaTrMhOF9IVbxI+f+lTxWE3xUlxauJ1pOBvF2E1DrRhPF/OTWLbEH8SYuNiSC0Jad6VBzN9er4OjXhHGmxF5FQcshGLjTQuB4a6u2acXv7ia2sQK6SudLOm6VPFL7VqC9B/AQY3QS3HgETXi+4z9yfKFG/3vOmGb4ZKYhsnBnNjyeXDldkAggsL9fVwiIgpJodd4IyIi13nrh4kr+TCWyrc0ek1Nio3lsULsd7Nn80QQERERkV+x540o2Di6newJve9+r13tniyky0ZHCOnMAeVC+oH23wrpnpFnhXQ7w2UhbUDLP62bpGkCqpvihPT3DR2E9Kfn+gjpH/anC+mMT68J6aiS8y0e22Mavd5EgaLAxoAFThUSwjz5EHX0IelKDJyDIfBKTJSQrhicJKQrh4iBaf1usXxpjmwrfoG2D78kpKP04mSmsqsmcSLVc41thPT+i5YvwW9+6Cxsu+mLBCGdvLdKSOvqrogH82WMmy+vNRERkQ+x8UZERESaoSg6q6lBOFUIERERAMBPo/NbPLz37li4TeppN8VZ7nyYYiJsvED10nrL7Xf95aviRm/dLg/WUfhOlSEIyqNRbLwRERGRZnB5LDbeiIKDyeT8L3FPlqmRjyHH77nQ03FhSEcx6ymVQnrZzVuF9K2RYtB/iiFGSBt0cixhONwnTrJ7T/RJIf2r+CNC+qvOKUL6xT5jhfSVtZb3mvj5v+0f2lHvhqs9Lt643r4cFEFEfud2462xsRHFxcXYsmULiouLUVJSgsbGRqSlpWHkyJF48skn0atXL5fzNRqNWL58OdauXYuSkhJER0dj0KBBePrppzF06FB3i2ubNz/QHH0p2uPgC9PYvq2QLvm1OAnvL4d/IaQntv1SSGeGWybOjdOJs95bf2G6pkk5JaQvJ+23lDNdnBD4/YG3CekP+g0Wy/l3cRLfsLPirPcedfkHy7UmIiKPcHksDxpvO3fuxKhRowAAnTp1Qm5uLgwGAw4cOIB33nkH7777Lt577z3cf//9TufZ2NiIsWPHoqioCElJSbjvvvtQVVWFwsJCFBYWYvXq1ZgyZYq7RSYiIm/wVeyUo55E9XFdmezVhX2VSDGOrV5a1u3fwy1fm6ZkcUoZmb7C8kO7U7E4C0CMvCRbg528PFmo3lf7usLLky1zwIIHk/Tq9XpMmDABX3zxBU6dOoUNGzbgww8/RElJCebNm4eGhgZMnToVVVVVTuR23dKlS1FUVITs7GwcP34c69evx/bt21FYWAi9Xo8ZM2agvLzciZyIiIiIWie3e95GjhyJkSNHWj1vMBjw0ksvYcOGDTh27Bi2bNniVG+Z0WjEsmXLAAAFBQVITLT82snNzcW0adPw1ltvYcWKFXj99dfdLTZR8PPk16+j18q/gF28BVw1tJP5cbfHjgnblnT6REh3CZNj2sSJc/1JDg/oECaW5T5pQuA+ff4ipP+/x39hfnzS0EPYlrRTDB2w4uiWvT+vN1ErwOWxfDRgQafTISsrC8eOHcPp06edes2ePXtQVVWFLl26YPDgwVbbJ02ahLfeegsbN270rPFmLzDckw9RR1+KMvWXpPTaxs7ipLsnHhO/eP405M9CekikOLltjN7+8HVvkr8U2+iizY8HiOF16JUsxuKN+o9vhPSMZLGR371APA/h5aqAePlaeTMI3BF715qB4URE5GM+i7QuKSkBAKSmpjq1/8GDBwEAAwYMsLm9+fkTJ06gtrbWa+UkIiIi7WiOeZP/QolPet6Kiopw8OBBREZGYvTo0U69pqysDADQuXNnm9vj4+ORkJCAmpoalJWVoW/fvi3mdf78eVRUVAjPNTcmiQKF9ZKCkdP10h8TvLoyKa8n+cqbVUsOXurbTthWO7lGSB8auNr8OE4vLlUou2yyTMx7+82/FbY1viveWWhzpNr8WFcvTejrK/6atNfdgSbUIq833iorKzFt2jQAwJNPPokOHTo4fA0AXL58PeYkNja2xX3i4uJQU1PjsOetoKAAzz//vEvlJvI11ksKRqyXpDVcmN7LjberV69i/PjxKC8vx9ChQ7Fw4UJvZu+0xx57DBMmTBCeKykpwbhx467/AnDy14ZiZz+dq0O37cRJmdqLv/SO/06cnPSDOwqE9IBIOabNfzFunpBj8e6ObhLS79/xlpD+5dXHhHTPlW3Mj/XnqmGXi78ovXat7Wxztl7aK4unHH68STF7df3FnvB20y2jvV/tLA5Q6BQWuAEJnpJjN7uHi+9F/V5nTI8WttVVi5MVxx6SBjA46vXxx/V2t14SUVDyWuPNaDRi4sSJ2LVrF3JycrBp0yaEhzs/Q3pc3PUPy7q6uhb3ae6di4+Pt5tXSkoKUlJS7O5D5G+slxSMWC9Ja2z1wQTDsq3+5JXGW1NTEyZPnoxPPvkEvXr1QmFhIdq0aePEKy0yMjIAAKdO2R52X1tbi5qaGmFfIiLyoZbuVPhqVLUHy7VZUc8dYRDzbUoUOwCqf5Zgflz/H2KM25cD1wjpGAdxbmrqmLi9Uj63maQYuE03mR+3+1osg+GCFCrU5Ob59/bSbc5SX1d1GdxscXFtUy+MNlUUBdOmTcP777+P7t27o6ioCDfddJMTrxRlZ2cDAPbv329ze/Pz3bp1c9jzRkRERNRaedzzNmvWLKxZswbp6enYvn070tLS3Mrn9ttvR1JSEkpLS7F3716rud7WrVsHAIzDoNZJsRP75Mnsk3p5CRkHeSWJ6+iWPSBu3tj1I/NjLce4uUr9Xl9WnQMAuP+Bx4V0r3LprkOltEavI7643iF2S4laNy6P5WHj7amnnkJBQQHS0tKwfft2pKenO3zNlClT8OWXX2LWrFmYNWuWpSBhYZg7dy6eeeYZzJw5E5999hnatr3+RVJUVIRVq1YhMjISs2fP9qTI9r8kXcnGxTysqpVq0EH5WPELc8UwcXb3/hE+mdEl6Mjvc8Wwvwnp+SfzzI+7vCvdRpDWBfRmELhLefFLkoiIfMztVsGmTZuwdOlS4MatzCVLltjc784778T06dPN6fLychw7dgyVlZVW++bn52PHjh0oKipCZmYmRowYgerqahQXF0NRFLz99ttONRCJiMjL7MVDeWvFEkfLtdk7jhQP15Ri6QGtGJAgbKvOFke6/+LWf5kf/78O+4RtBp13RvPLo+2/HvR/Qvq/0geZH3/yVbawrd1BsfzJ+y0xcYbzl8QDuRLJ749rKh/HwRJ8zlAUG1OFsOfNOdXVlqkadu/ejd27d7e4r7rxZk94eDi2bt2K5cuXY82aNdi8eTOioqIwatQozJ8/H0OHDnW3uEREREStgtuNt7y8POTl5Tmxp6i4uNju9vDwcOTn5yM/P9/dohEREVErxalCfLQ8VlBTFM8Cgt09rBQMdfVmy7xKvcb8IGwbFl0lpA06cVLQ1kqeKFU+D+rzVP2lOF1M1NflCArufoLYq5eKB8P3TfJtJ/t5VQ0QJ4yeNniHkM4M89lyyJohn4Npg/8ppDd9PkJI37TFOkREoHPtGtklX29zniH2zUbUyoVe442IiIg0S7ER4xZqP0/YeCMiIttMJkuwuQu9d+oR2lbLy1nvbP/4Lb0sQZyq5qcRyUJ6wG8Pmx//v/bbhG2dpd5TcYF5//Quy3ca3kj7yvz4xV+IvbmnxojnYdm5e8yP96/OEral7qgQ0rqayy0XwoMeWbvX2Nlr6uakwJwqxF+1lIiIiIi8gj1vRMHA1HT9z9sU+3nq4sXei/O3i/uPiT+MllxoqhfSZUaDkK5RIoX0zWFiD0A7g7g9Uuf8WsiualAahXR1U4OQPm4Uz0OCTtyeEWY5L5E68WNTPkerhoij4pN3xgpppVbqCXFwjVzSUl6+qFtEAWKyMVWInG7tQq/xppicDwh2YU4hR3SR4hfVmbss6ac77BS2xTiYV6jJQfnl7vhg4Wq55fMwQ3WenrjrYWFb96PiF7/SIH75OuSta+1JsDkREZETQq/xRkREzjEpUG7EJanjmlxZdcSVfeXYKUWOiYq39GL++74UYdNfn1gmpLMi1HFsMU6XIRiIMXhAL+n3/J87f25+fPjpz4Rtv4mZK6Q7bVL1utbWCduszrebMXCOXifUHSHmzd3R+TZ+b4fYiIXg7KIhIiIiIpvY80ZERESawdGmodh4szU1czOXunClffX2K44pSVybLibbMgFt5zBxbbpzTeKtgu31XYT05zU3C+lLjeIkvve0+05I3xVzwvz4lnAxeNrbfmi0dMv/s767sG1bdW8h3Sb8ipAe2kaerLhMSHcOs8Sxqc8fbJxf3enz9gvqq2vt7tB7U2Amj27sKE7K+7NeZS3uCwAHr1k+Mv7f6bHCtsOn08S8r4hxiP26nRbSUzvsEdL3xojXTL515IrLpqtCurBevMW25uztQvqbkx2FdHi0OMAhq+MZ8+P/6vgPYVuUziik5XNY37GDkA77rtaJd+BlAahbROQ7odd4IyIi5+gs8UqKUdVIdTQoyt7AHTuvVaTRsrow8Svqck/LXG59xx8VtokxbqFDft/yeTn9g+XHfty/xBHiwjWFg+vq5jWFdF2Fa+pmZxl73th4IyIiIg1RbIxPCLW+ZQ5YICIiItIQ9rwRBQOTCWi6cWvBQfykZ8cRf5/WdBFvuWQnnBPSPzWJsYRv/tuy6Prxf4qxmG1+FA8VXi8e68fO3YT00hHixLgde74nliVSjDvT2/mtaYJ4S+fINTHebmnJKCFdt0OMgetQLr6+MVacm+GbrpZbTy/dJeb9+047hHQP6Rzu7NJVSLf7Rpow1x/X281liIiCkQIbt03dvQerUaHXeFN/SXqTgyzrM8QvwVtTvzE/Ptckfomtu9RLSK89MFhIR5WKE/6GiWEM2HezOFBgTLblWPkpRcK29DBx/iNHE/zKE+2WG8WDLz+fa378/zvYT9gWd1z80jNKUy9t6yoOaHgoe6+QHhFnieW4NbVc2PZ9hnis2PKz9t6GZ+xda35JEhGRj4Ve442IiJwTGwMl8sYI9QuqUfGerCTiwmuVeHF0fHVvy1fWc6m77L5W/UMzWFedcZe99zZNOi9zels6A2KPSbMNXBBnOnD7urp7TRuuAJXuHI9Bb62rRhMRERG1cux5IyIiIs24Pl2rPFVIwIoTEKHdePPlxJVSEHJNuniqkyMumx9/3yBO4vnXb24T0u2+EIOnE8quCemwejEI65q0SPunDT8zPx6QWyps+1W8mA6HwcabsWiU5mHaXp8pHmuv5VjpheL5jbgkTsprjBGPVXNWfJ9/DRPPQ4eci+bH6vMHAF9K59dqKmI/Xmv38tADhhvnQ4rJdHe9QdhYu9B8jBvqU+13vh+50llIf3ck3fw49YhYrtifGoS0oV6cQyqhRNz/VLg4aGBzWn8h3SNpn5AOt3Prq1G6bbO5ZpCQrt0tHqvzjhoxA+k8NcWI9Sm8zhJn+l1iurDtSJJ4jmTyOW4nXQO/XG+FN1mIWpPQbrwREVGLTDGRMEVfX8FFX3nB/8ePE1ePuXyL5YdrvwixAX7ZJP0wUf3QPGUUf9DKo6hjdZZ8B0aKg7A8WenDHnkVkH81WEZw1SniD9lUg/heO4dZRmLH6MT33S9C/DGgPmfy+Qz0NTUZIu3u2yIuTM/GGxEREWkHV1jggAUiIiIiTWHPG1EwMOiv/wFQrl1zuLuz5DsJugjx9tG1NuIe1dfEaMGyenHh+riTlls0UVVijJuuUYw70zVIk+zWijGPyYfEWzj/uKOnkJ7R7gshHWPnvki9dA/lH6fFvJIPiWUxVImLw5vixbLowsXftVFVltfHnRRv9RzMFmPe4sPE8yKf4+br3Mwv15sxb9SaKLrrf/JzHjhw4ACKiorwr3/9C1999RVKS6/Hg3/zzTfo27ev2/levnwZy5cvx4cffogTJ06gqakJHTp0wODBg7FgwQL07t3biVyshWDjTWdeRFdRjA73dvsoBvHUNiSK26sbLV+S5xrE+Iuo78QvkoRy8YsnvEb6IroibT8nxke0T7As5vxh3xxh2y/iTgjpSAeBAw1SYPiHP4n5tf/C8h8o9gdxAh8lRvzS00lxKAnl4n++q9J52H+zZUb/cL0Y19HQViqoFBSuNDXCV3RCzElodd1T63atXRQabjRso9X/Pxulz05PBu2oBxOFS5+bKeJnQJ/M0+bHTVKDfU9DvJD+vO4W8+N13w8QtjVWiPkqEZbPtZFZ4sLuCzp8KqS7h4uTqjvrRKM4yOrFs6OF9PbDlvnYdNfExnZ4svjDZ1LP/ebHd8T+IGzrFyHGsanPWV1KR2FbtPSZK1xXb11TiNdVfU2v1fomntAdixcvxsaNG72a5/Hjx3HPPfegrKwMHTt2RG5uLnQ6HX788Uf8/e9/x+jRo9l4IyIiotZPsTE1iKfjFYYMGYJ+/fphwIABGDhwIO68806UlZW5nV9NTQ3uuecelJeX4/XXX8ecOXOg11sa5WfOnIHR6H4HEhtvREREFNLmzZvn1fyWLFmCsrIyzJ49G3PnzrXanpaW5lH+bLwRBQNFNfbdn7NNSnEiFVfFW0JlF8X7/dEVlrLppGI2xou3wXXSrRP9pToxr3Ix7qzslHis2j7irW/5VrmwrzRNRIWUV89yaV43aQ1aU4w4NYP8XvRGy3tRnwMAOFYpziGX0VaafsFRLI4vr3cg6hSRrwX58lgNDQ1YtWoVAGDOnDk+OUboNd70Osu9fF8uIi5P+hkpB4Zb5vQ5eTFJ2BZ7WtzX0CBN4hkmBTxLx9JfE+O72hyzfEke/VFs7V/sJha7nd7+ObkobZbz66E6lk4qhylWjG+Q34f8PmNPi9Xz60rLsbq1rRK2NUVJ/3PlyUp9ea3VsSHemLCXiIg0a//+/bhw4QI6deqELl264Msvv8TGjRtRWVmJtLQ0/PznP8eAAQOcyKllodd4IyIip9SnhEFJvN4LGRNu6Y00ST/MdJ6sU6/qFdSHiz2elzuK6V+3/9r8+KQxRti2rHyUkC7fnmF+nP5PMeA/7JLY66seGf31nf2EbS9OE8u7olOR+XG0TuyxlV1RLCOJ5QEKX68Sj9Nrt2WAlxIpvm9jG3Gw1/q7hpkf7xvZRdj2/3X5REj/h+qc/bmjuG/M1+Jx1NfVW9cU0nW9nGZ5fOWCe00Qe/O8lZSUWO2fnJyMlJQUq+d95ciRIwCAjh07Ys6cOVixYoWwfdGiRZg2bRr++Mc/IizMvXPA8eNERESkHUoLfwDGjRuHvn37Cn8FBQV+LV51dTVwY/qRN954A/PmzcPJkydRVVWFd999F4mJiVi1ahWef/55t4/BnjciIiJqFTZs2IDMTHHN7eTk5Bb39wXTjTCdxsZG/O53v8PLL79s3vbggw8iJiYGDzzwAJYtW4Z58+YhLs716WdCr/FmUixz0PgyiLdJ7HM2XJMCw69YLlbFaXGSsrR68bWNseJlqksVg7PjzkiLsleLwdmGC5ZA8Ygz4pxyF01it3+MTpxgVHbRJHbfR5wRu90NF1Rzu0nn92oH8TbH5TSp3D+JMW/h0nk4ozpP8ZHSIujS+ZXPv0+vtTowX57fyOk8TFYLlPuENDTdIC6viOqr4jW6eFasLx0aLe/vWhvx+tV0FutlmzKxYz/+lPj+9HXirayISmnAglQ34/Xi/vb2jagUyyIfC0axLFdSxXjMSxnie0tQld3QKF5j+RzFR0l1UzrH8jXwqeY65cuYTyI/s3fbNDMzE3369AlQya6Lj7fMOTh9+nSr7ePGjUNycjIqKirw5ZdfYuTIkS4fI/Qab0RE5JSryTooyTe+JNWxOVJj0Ooni95ORI69hqQU/1OfKn5Bp4dbBiptuCgGfJd/liGkO29TDZ6SfszprkirWpw9b36YskNsXO+8RYxNK5uwzfw4I8z+6hhlqlHKO/8p5nPLjp/Enc9b3puugxifZTCI56HzPyxlLFfE971hvHhe7k74zvxYPp/y+VZfG69dU+k49R0sZWgIa50DvLp06WLzsbxPRUUFfvrpJ5vbHWHMGxEREWmHnZi3YJCdnW1+XFVVZXOf5ufduWUKNt6IiIiIvKdTp0649dZbAQA7duyw2n7ixAnz6g3uThnC26ZEQUBpNEIx3bgV4sP4JEVakzKqSvy5eumKtNagXtx+LdZym8OYLP72u5oi7htdKf02lGP6GsVbJmF10tq2ihhPCbQc8ybvK+dltRanVJZrsfbfi041CXBYvfQTX2//HMrnWL4G/rjeSpMf4+yIfE5nYx1p/9+CnTJlCr788kvMmjULs2bNErYtWLAADzzwAF544QUMHToUP/vZzwAAly5dwiOPPIKmpiaMGzcOHTt2bCF3+0Ku8aY0NUGBH74kpaDkSGnS9Zqrqg94k1jpGmOkL5J24msvdxbTYVfE4OxYKRhbHbwfUSMeq14agACD/QEL8v5yflDP/yQPWEgUy1krhQIo0mLyVoHeqsEAwvmzcX7l8++va634Y9ABkZ9c6WBEU9qN+h2haiDL9VyOh3Ll/4H6/2aE2AivTxf/H9eYLP/vPy3rJWxr9714zKZYS14NiWK+4ZfFhemj1fO+SYNboirFz7jjjZaRi8mGczbekHrf9i3mo5MH0cRYynQlQxzE1hgnfjZGXrB8zsrvWz4vt/Y5aX4sn0/5fAvXzVvXVDqOugzXwoLnR8WWLVuwZMkSc/rs2bMAgMmTJyM6+vq1ycnJEaYdKS8vx7Fjx1BZWWmV37hx4/DEE0/gjTfewKBBgzB48GC0adMGX3zxBSoqKtCzZ0+89dZbbpc35BpvREREpHFejnGrqKjAvn37rJ4/fPiw+XFUVJTVdntWrFiBO+64A//7v/+LgwcPoqGhAd26dcOjjz6KJ598EgkJCU7kYhsbb0RERKQdPljbNC8vD3l5eS69pri42OE+EydOxMSJEz0omW0csEBERESkIaHX82Y0AqZGJ3b0kBQfEH9KvLd/tsYyGaohQZwrqKabOFHqtXZiXjop0NskzZUjx3vpVIu066W33iQFeTqKDLumSJOfyqdSHf8lxbzJ5TTGiNsv9BXTEdVSDJzqPF2sEc9Rh1P2A9J9qlF1EkxuxnCYTAD8UGbpvCSUiRewtCpWSIfFidsv3Wz5yGiKlybdrRd/Cyp6B/VSKppeDtWU9miy88ta3lfOS54YVy6LXNYmacnKmp6qSXprxfcpn6Ma6Rx2KZP/0/mxbnKSXmqNFADSJL3BNFWIP4Re442IiJwSFtuIsBs/mkwxlngf+YeZzoPGoTov9TEAmI/d7McGy+S1dafjhW1J0oosV25SLYB+k9jgjgoXv/ij4i0NbtOP5cK2iJquQvpik+WH4zUHK7eo942oEfdtqhCD3PVd082PryaKX81X20k/sg2W9xZ+WXzf8nn5MdNyzuTzKZ9v9bXw1jWVj6Mug+mSHzpSWik23oiIiEgzFMV6xUNfroAYjBjzRkRERKQh7HkjCgKKoljdavCH6NKLQjr+O3FNxctZ0pqQHVWT75nkeavEOQDjTou3aJQGaR1IaV6/JmnKwSid87dU5H3lvORzK5dFLmt1PynoTRXXZkqQJ90VfwPHfye+Nrr0vJA2BeA6B6JuEfmMD0abak3INd7sfUmqA/vdyVcgBSXH/lgjpMN+sMy823hLvbDNmClO3qjXSV8W58U4hegqKQD6mvhFpqgmSFQ87Gs16KQFqeWYUfV5kMphVU5pEkh9kjhBsLFdy9cj7AdxwELsj9ViORwEhXv1Wju5jUhrdAYFesON//Nhqv+vJukzyZPPFXVeYdJngkH8vKltsnz26YzSjwerMlm2N8bKA1LEdHySan3JEvGzI/yyWNw61UTljQ7+u6v3lfORP6OaVGWoSxXPg0maSzeq2lJ++X3L50V9zuTzKZ9v9bXw2jWVjqMug87Az0t3hVzjjYiIiDRM0dnoOfD/8liBxMYbERERaYYOgHRDKgArmwYWBywQERERaQh73vxEd06MyUr5l2VNs58M4iLJ19qLE4jq68XA7sQj4m+M2B8vCWnlmhh8rQuzXOZrbcSfK231Unyd3Xdhvf+1ttLPH3XMhFQOuZyJR9oJ6Qt9xVg+U4wYExJx3vI+Uv4lTfgqnd/WHElhL17P5Zi7SqleHhDX2jOFi9ekvqPlmoRJk9Xe9LV47KgfxXmsTFKMj/w+GhLFeJx4vVh/DHZ+Wsv7ynnJx5LjjeSyJu/vIKQrf2aJXTLGi3nHnBb/f6YcuCqk5XPsKq9eb6LWgAMW2HgjIiLbGuvDobt8ffSsrtEyoEhRpMB3k/SzT2+npS0Hs6vzMor5XqsTR+42qlZ4MUVJ+yaIjeiwq5bteqP9AQCKwbJdLrlOaiA3SqvM2KPeV85Hpi6DXD69tHCL+r3J71s+L+oyyOcTRmkEuPpaeOuaAtA1Wn4sXbtsGWzWWC+9UXIaG29ERESkHRywwJg3IiIiIi1hzxtRkHNlTjqr2C5Ht2rqxDkGI4+eFtIdr7QX0pfTLbc8wmvFezkxJy+IeVeJaSsxYqwnksV5/qJ00sL3drKS95XzsjpWTa3dsrY9JN7OiahJND9ujBdvU8WVixN4hZWeE9Im6Rw74svrTdQqMOaNjTd/UWrFD/i4b86aH3eqF2e1v9RVjEuIrBVrZfwJccJfnBODreVF3nSRlvyutZOCrXXiF7CjaA55fzk/9bGUujq75bzpoFj9ImrFBZWvxYnb2/xoic+IOiHNWi+dXyLyXMS5METiekNWd9mFRqgcA+Ukfa34mRH570QhfSCjs/mxIV6cBLyhrTi4Ju6M5bPKcFUqT6Q0wW+T6nNMXlRdGikTrvqh4GjpdmFfecSN/DmtKoNeWlzE0CDuqzda0vUp4qe2IV4cMHOg2nLOIv8t/iiRz7fdqc3dvKaQ6k7kadVAtXNsgriLZ46IiIi0gz1vjHkjIiIi0hL2vBEREZF2sOcttBtvnixO7igvq8BhkxgdYTpXYX4ceUGcvLb9UXHRdcUoRSI0NEjbpUmAJEqsJVg7Jk2MDYvVO4raEMn7y/mpj4VqKYC9XpoQ+PgpId22PFJI68LEWA51cL1JmgDYEb9eay/n7828HJXVdFGsi/rDYmB/m2OqayTN5aQ0SvXQwbGa2sQK6cy0CiEdJa1/Yy8eU95XzqupTRshrftJKqr8f+jfZ4VkjOr/K3TiDQtF+v9oMrn2fyqY62bMWSD6RgiVIn1O+YJ8jLiyTkL6dG/LdTSEiZ+L1+LF93451fL1djVZ3Bbzk3hews5b4ojlT9OGttIi9qqJyh3dulLvK+cjU5chvF78v1GfKr5W36iaeF163/J5OX3Bcs7iysRj+uOaysdRX1ODgzFNLbMxVUiILZDF26ZEREREGuJR4+3AgQN49dVXMXHiRHTt2hU6nQ46nQ5HjhxxK7+8vDxzHrb+Ro8e7UlxiYiISON0iu2/UOLRbdPFixdj48aN3ivNDXfccQcyMzOtnu/Xr5/Xj0VERESkJR413oYMGYJ+/fphwIABGDhwIO68806UlZU58Ur7pk+fjry8PI/zIdKK5t5lfxzHJVL8lnLlSou7Qu+gI98gRq3Vp4txPaOT/yVmJ7/cTtmbpNiuYcnHhfSn6cOEdGyJFEEnLVQvx+spV6VJfz3gj+ssH8ufxyTyOQ5Y8KzxNm/ePO+VxE+C8UtSkYLv5bQVF78kG5PjzI9vTTshZiW/1EG55S/JW9PKhfS/k7ubH4f92/6klFZB4vJ5sFsS+wLxBenv4xL5WuKxOsRFXx+0Yqp3baUId8jHuOlrccDMiZstE3krncXJaOvTxB8a6kEKYdJ84UlfiwOtTBVVloRe/PysSxfz7RjufJS9el85H/k46jIkfS1OWF7fXkzXqG5MmcKlfI1Svj9afiTJ59Mf11Q+jroMl6/UtfAKcoQDFoiIiIg0JCinCtmxYwe+/vpr1NfXIzU1FcOHD8eIESMCXSwiIiKigAvKxtvatWuF9OLFizFo0CC89957yMjIcPj68+fPo6JCnOeppKTE6+UkckWo10v1urcAcOFm8eNnYMxJrx1LzuvvN98tpOM+F8siz0EYSkK9XpL22BpdytGmAdS/f3/ceuutuPvuu5Geno7q6mp88cUXWLBgAfbt24fc3FwcPHgQcXFxdvMpKCjA888/77dyB5wUZ6ULEy9rTVfLxLnD2h4TtskTnzq6jy7vL+f3Zte+5sdJ30pfkHIsn5cntw12IVcvSRPs1Ut96VnoDddjy5pUC5PrDPamTHaRatJjRRo4YvhBjKlN33az+XFV72hhmyJ9m4XVW8p700Ex1kvO16QaiGNIThK26dPEhn283hJrJ681L1PvK+cjH8dUVd1i+TojXUhXZlu+A40x4qe2rkw8L0nfWQbb6KV8m6TF5r12XeWJrFXXVf3e9E1uztKr2Jik12rS3tYtqBpvc+bMEdIxMTGYMGEC7r33XuTk5KCkpARvvvkm8vPz7ebz2GOPYcKECcJzJSUlGDdunE/KTeQM1ksKRqyXRNoTVI23liQkJGD27Nl44oknsHXrVoeNt5SUFKSkpPitfETOYL2kYMR6SZrDqUK0M9q0R48eAIAzZ84EuihEREREAaOJnjcAqK6+Hg/gKN6tVXA0j5uDucR0MWLMw4Velse9I08L2xzN6yaT95fzUx/rpu1STIoc8yYf29H7dnHBb/IxV+tp2wQhebm3WB86h4mLZLsSjynvK+clH0suC66Ic4ZZxWOq32sI1UPT1asw6a7P0aXTq66n3jfxRTopbkmehyzy0I/mx2llbYRtCBe/znRXLLFepioxtspkZy5NpX07Id2tfaWQjtJJEzrbod5Xzkc+DlQxbyapPuqPlgrp1POJlnyiI8V8GqX5M1WLwsvnUydfRz9cV3UZTMrVFl7hhBDraZNppudt/fr1AICBAwcGuihEREREAeP3xtuUKVPQs2dPrFy5Unj+0KFD2LJlC5qk0UZXrlzBggUL8NFHH8FgMGDmzJl+LjEREREFCy5M7+Ft0y1btmDJkiXm9NmzZwEAkydPRnT09VtmOTk5KCgoMO9TXl6OY8eOobJS7EIuLS3FAw88gKSkJOTk5CA5ORkVFRU4dOgQKioqEBERgbfeegtZWVmeFJmIiIhI0zxqvFVUVGDfvn1Wzx8+fNj8OCoqyqm8srKy8MQTT+Crr77CkSNHUFVVBYPBgPT0dIwfPx6PP/44evfu7UlxichVjuLa1OT5BqU5oxrSxRif7MwyIR2r914smZyXfKxL6eK8WZHnq4S0PN+YEAPH2EyiwOJoU88ab3l5ecjLy3PpNcXFxTaf79atG1asWOFJcbTLlS9IJ15vukkM5I3qfdH8OFnfAJH4hWuA/YBVk/Q/RM5PfSy5HKgRJ8p0+UsuRIPGiYKCzg9RNvIxFPH/ualGtaD8ZQeLmpta/jbXGaRJZFXHrc8QB7N0ixEHZYXDUiZHZ0S9b4eYGmHbyYyOQjr6qCU3nfRDSLkqfs42nf7JknA0yMDOefDLNZWPo3jhs5uNN+0MWCAiIiIiDU0VQkRERMS1TdnzRkRERKQp7HkLFA8Cwa02S/Ebtd3EmI170780Pw738hyMcn73pn9vfry3223CtvgyKc5EjluT36e9hevl88cYOO/wNP5SLSJcSF64JUJI/+dN38Nf7paOteqWTCGdelgsK644PxGrldZUN/V6QO/FRegdHk/+gJKOrY7fkgeVyNRxVlaT0cqxdZZ8628Sj9kuXIytM7jQxaPeV87niHScaGES5JbLB3hwHry18LyrhPOvKoPJzc8bLkzPnjciIiIiLWHPGxEREWlLiMW4ydjzRkRERKQh7HkjIt+QFn+/1FP8qdwz8oxL2XkSOSYfSy5LqqOF6ik4tBQ7ZYMwV5qDWE5dmOWr8GqyGDt1c/Q5IR3uQpePel85n39Ix1GXQY5xkye8hs6k2jW0Yr2A671uVqGHIdYTx8ZbK6CLjRXSVX3F/+h3xB/3W1nUx9rSd4iwLeFLsZzKxUt+KxcREbUSnKSXt02JiIiItIQ9b0RERKQZnKSXPW9EREREmsKeNy2Sgm+bUhOFtLHPZSHdJazKL8WSjyWXQy6nztOF6imoyEHV1zqL1zulR4WQTtLX+6Vcto4ll0Uua8R58f+MYjT6sHTkLHnBdq/lG2GZQLquo/g5lBxWa+MVrpPzkY+jLoPSIC5Eb48r50SxN/G5ljDmjT1vRERERFrCnjciIiLSDva8seeNiIiISEvY80ZE7pFibXSxMUK6ukekkB7b4QchHa7zX4yjfKy7pbJs7XGnkO7wvfheFHV8ZmuJG3KCTqczx1QFIl7Koxg3dWywnI88CW5iG/Pj8FQxPrKtQVxQ3l1yPvJx1GVQfjovvlguv/q9uRArLJ/PQF9Td6+vzsbo0lCbqpiNNw3SRUYI6Ys94oX0iG5fC+kYvfPB1k0e9j2rj5XbXfyCPNijv5BOLBPfh8JZ7YmIiBzibVMiIiIiDWHPGxEREWkHByyw542IiIhIS9jzpgXSpLxomyAkq/qJoZqj2n7rj1I5dHeb74R0Ub9sIZ34pfg+0HBNTDc1+axs1AI5+Fmue3a2KYni9bzUQ/wpPCTuuPhy6adyk/TLWS9FIItTAIvkmmKVF+yX5e897hDSqdJ7wWVVsLmjetlKJ5tWB5f7KtDdawMUHJGu0eU+KebHt3YWY3WTvTRgQc7n1s7lQvpUn1vMj2PO/CS+2GCn9svvOwgHMHh7cmUuj8WeNyIiIgpxBw4cwKuvvoqJEyeia9eu5pHWR44c8Ur+iqJg6NChXsuXPW9ERESkHT6IeVu8eDE2btzoWSZ2rFy5Ev/85z+h0+m80sPJnjciIiIKaUOGDMGzzz6Ljz/+GKdOnUJGRobX8v7xxx8xf/58jB07Funp6V7Jkz1vROQWXUS4kK7v2kZId+4jxu2kGmo8Ol6TF+Nx5LLIZa3vmiykY85Vmh8rVxiL6SiGyV7Pgq8Wl3eJFCdW2c/yVXh/ghiLZrDTpeNKdKOcT450nIP9epsfp38WmH4Ve9cmqK6pD3re5s2b51kGLVAUBdOnT4dOp8Obb76JoUOHeiVfNt4CRR1U6iDQVhcmXqaGjHZCuk2/KiHdJbxSSKs/NLwdSi3npz6WXA65nPL7iKy6KKQVe4HhrTQonIiI7NPSgIU//elP2L59O9544w2v9bqBt02JiIiIvO/UqVN46qmnMGTIEMycOdOrebPnjYiIiLSlhZ62kpISq+eSk5ORkpJic39fmjFjBq5evYq3334belemsnECG29ERORVfomBcvRlqC6DvBB9RkchHTGo2vy4d9Rp8TBemrpfzkc+jroMcvmUMnFfe+/Nk3nf7AmKWEUnjBs3zuq5hQsXYtGiRX4tx1/+8hd8+umnWLhwIXr37u3EK1zDxlswkv6T6NqIC8+fuzVSSP+uyw4hHa5z/j+rVRC4g/+grgSNy+X4dZd/Cel3bh0tpDNOiO9TqZAm7fXRBJLBzpNh5V7/wFVPFpooDlA4N1AcwDC701dCWq4PBj8GqcjHCpeiNR/s/KWQ/n+3il8AXY+p3uu1RjFzL08mHVTXmygI2Yt527BhAzIzM4VtycniACRfO3v2LP7rv/4LvXv3xoIFC3xyDDbeiIiIqFXIzMxEnz59AlqGRx99FDU1Ndi6dSsiIiJ8cgw23oiIiEg7gnxh+k2bNiE2Nhbz58+32vbTT9enJfrtb3+L2NhY5OXlIS8vz+VjsPFGRERE5EV1dXXYuXNni9v/9a/rYUTDhw93K/+Qbrz5ahFeeBh7oosSY9rqe7YX0jF3VQjpQTHi6Bp7k0r6k1wOuZzrpPdRf0h8n9G1l4W0cuWq22UJ1mttizfLKuflclnl+MuYaPPjmv7i6K2MYWVCun+UOAlpMMuKPCWkM4aK7+XSt53Mj9vsEhcYVy57tnB5UF3v1kR1LuQ5I8vHiTFQv+663fw42VArbPNWbKacj3yc8V0PmR//fdxIYVvnZeL/JZ29hepDQZD3vNn7P92lSxeUlZXhm2++Qd++fd0+Bud5IyIiIu1QLIMWzIMXAtB4mzJlCnr27ImVK1f6/dgh3fNGREREtGXLFixZssScPnv2LABg8uTJiI6+fschJycHBQUF5n3Ky8tx7NgxVFZW2sjRt9h4IyIiIu3wwW3TiooK7Nu3z+r5w4cPmx9HRUV5dhAvYuONKMj4NU5JntDTKsYtRkhf622J/frpPxuEbSvSPxPS4Trvzn/my+Xg5bI+Lr2X2Q9MMj+OruwkbIs4IsbLKfX1YuZy/Is0aWpIx6V5k3RelSZL+sro/sK27qNPCunbYk6YH7syT6Yn5OOoy/Dl6C7CtktHfiakoz+1xMfpDNL/YS/P5B8q3Bn1WVxc7PJxSktLXX6NLWy8BYrqP5guQpzcVOmcKqTLxor/GRd1FyfljdIZfVJEb5PLOUt6H4vGjhfSt5zrIKT1pWfMj5UGseHAheqJiEJIEA1QCAQ20YmIiIg0hD1vREREpBn2lscKFWy8ERGRRwKxHqs8d5usabBliaQrj14Qtj3RYY+QbmuQ4hQDQF2GqVL5XnlUXPc54qLlvRm+/E7MSDovbp9frrEb1EKv8abXm+PNdL6Mk5KDRqVJFXWRlvXOlE7i5LQ//qe44Pf04duFdPeI82LW0s1/vQvBAFYffw7+w7oSNC6XQy6n/D6mDxeDP9deFCeq7PqR5fW6f58TtikN0iL28geY3661e5EIuvBw6AzX64RLH5qO3pdUD60+VKV4S11crJCu6yXWzbKJluO9NGCjsK2t3v4XoHz9HU0mLb8ze2fW0dV1dGxpqXmr9/KC6r0uuPafwrb098Tg8tijUt2UJ/GVFra3ut72rqmLwejN11vXFA40ONydSBuCfJJefwi9xhsRERFpFm+bcsACERERkaaw542IiIi0g7dN2XgjCgZKp/ZQoq8v9K67IsXveRKvJ8VIKVFijFtDezHGrbpnhJDW51YJ6Vd6fmp+nBp2SdjmKK7MlVhMW7wZtegoHlNOp4VbAt5fvv0DYdsrSfcK6QtF4iS+SUfF6xl5Xoyn8+n1jr4RR3klAvje/WzhYQC7u/naC3zXR4uz3dfffouQvjDjsvnxf0tzSiaH1Qhp9fV2NTbTWY7yVafl8slzYr7+3/eYHyf+KUvYFrPnByGtXG052DHYrik5L+Qab7q4WOjDEwAASqMUpmyvIpukbXqpAsoz04eJp1ZJiBPSV7ommh+fukfcN2+U+B/19tjjQtrbM9cHivw+5PdpGiee07/EDDc/7rxNnPk/+kdxNJmu5rKQVozSRMZWQeIOrq+QuXStwy0NIl1jI3Cl5ZcSEZGH2PPGmDciIiIiLQm5njciIiLSLt2NP/m5UMLGGxERaYJ6HsJr/TKEbWX/IX59T+5yxPxYjs+MkGasVMebeRqb6Sx7cZdy+eTy/0L13t79jyHCtu5XugvpiG/KzI+t5hwkzQq5xtulwR1xre31hd+jLoj/QfQNYqCwet4YnRQTpUhxT02R4h3oq+3ESXlruorbwwZaYrSevGWXsO3myJ+EdLiDhef10gQ3Bg8mvPFmNJ1cDqsPRSkpv88hUgxcyn2WIN7/zRwmbGvcL04mm3AyWUhHVUsf1tekay1fX1XMmyL9pDPJ1zrRcq2vXIwCPoHLfnwgDhEdrsdiRlaL+euleHbhtMqXWiprkzj+ANcSxfcdnVErpO/r+i8hPTRejHKP0lniRA3SEAK9Tk579iXYJL3cYOentbyvq+Syhkv/E5pUlSDJIMZTLu61SUjv6tRTSG8t7S2k60sThHTEBfF6G6TrLVxj6RxY1c0Wrve1s/XAyyBqPUIsxk0Wco03IiIi0i5O0ssBC0RERESawp43IiLSBtUUP8ZoMTQlsp04R096pGWOwii9OC2Ut2/xe4O6DHopJEEuv/q9ye/bKM1/F+GjudwCilOFsPFGFAweuf0zpN18PRi70hgvbKuXApmaFEuHuUkKepK/hGKkgLmbwsUYtxRpMlB1TBtsxCGq5+aTJxmNkL5wHE1u2mQ1PkyOkZT2d+HDWZ7m1vpYIrmsJjvB5I7mWRwRf1RID+lbIqTP9xRj3iob7V9v9TW2jm8V36nV9Q67fr3PHK/DIrulJiItCbnGW/sJ5Ujodj3g+MeL7YRt9VcjhbTJpApaN0kTs+qlL8ko8UMzo604aewD7U4K6V5Rp82PY/XiDNjyl4OnX5JqVl+A0neao/voVl+KLnyhOnofJuno8nnoHnHO/HhJn43CtqPdOwrpL6q7Cemyi4lCuv6q+AVp7/rqra61eL26tq02P645WY0jbgxYICIiJ7HnjTFvRERERFoScj1vREREpF0cbeph4+3AgQMoKirCv/71L3z11VcoLS0FAHzzzTfo27evW3kajUYsX74ca9euRUlJCaKjozFo0CA8/fTTGDp0qCfFJSIiF5gajTBJcZD+ZBXloVgmmVXPrwgAPdufF9Idwy3hDPYm5bWVdpan8wu2VAbrEBOx/Or3Jr/vnxK7Cuko1cS8SqP9OUP9Qf3OTErgy6NVHjXeFi9ejI0bNzqxp3MaGxsxduxYFBUVISkpCffddx+qqqpQWFiIwsJCrF69GlOmTPHoGA8kH0B6h+ujcRpTxbevDgSHgyBnq//8UuCwHOgt76/ebh3jJuclTSbsw58YckybL7kyMSogxufJcYI/iy4X0n07nhLSjWn+udbltVfxWYs5tSwlrAYdw6/HTaaGX3QjB+fI70Ue1RZhVdek7apr5ChwvzUR3quDuFFH5zg9vEpIq7+Iva25LMawqz47BpHfMebNs5i3IUOG4Nlnn8XHH3+MU6dOISMjw4lXtWzp0qUoKipCdnY2jh8/jvXr12P79u0oLCyEXq/HjBkzUF5e7kRORERERK2TRz1v8+bN81pBjEYjli1bBgAoKChAYqJldGBubi6mTZuGt956CytWrMDrr7/uteMSERGRtoRajJssaAYs7NmzB1VVVejSpQsGDx5stX3SpEl46623sHHjRjbeiIj8wBAXA4MhDgCgi7JM/qrEx4o7honxZ2iyE4ChkxfgtdyS1l2uFzYpV8XbvUpGmvlxxRgxdOJ3Kd8I6Xi95bXy7X97MW6erA3tCvk46vn85PI1SeVXv7ex0vt+ZUyakE48Ypk2SV92RtimvqYAoMTFqAooXVN7k/0apJt4RjGkQlerirtTXVNDUyMgTjXpHN42DZ7G28GDBwEAAwYMsLm9+fkTJ06gtrYW8fHxNvcj0qIwXZM5dk6v+G4GHzlezyo+y06Mm7zd+rXe/fT0Z/ylzOq9qJLyF2mEtOs1qW0irRWPJumd+eN6h4VQfCJRKAiaxltZWRkAoHPnzja3x8fHIyEhATU1NSgrK7M7mvX8+fOoqKgQnispuT7LeYTOaJ5FXg6QtzcJuzyBrByELHPlS9LRAAVHI6H0HvzkkEdLGexPRO/R6Cq5nPIgAetBHdL1UW92ECQu5+3oWrt6fYVjq65lhK7l0VP26iVRoLBektZwqpAgarxdvnx91YPY2NgW94mLi0NNTQ1qa2tb3Ac3Yuaef/55r5eRyBOslxSMWC+JtCdoGm/e9Nhjj2HChAnCcyUlJRg3blzAykTEeknByF69LM/rjsj2qQAARRUCZYyVeqatFqK1dG1b3YEOU1rcN6xOWqZO6jAP620JkHq85w5hm3r5PEi3/D25xW9vDV5XlhN0tL6umtUatlK8mfq9ye/78WzxvPz5uTvMj43f9Ra2KVJYmzFWdRyDeEydUSy/MJOTtK984yKsrr0lH9U1bTj3E/AGXMeYt+BpvMXFXQ+Kraura3Gf5t45R/FuKSkpSElJ8XIJiTzDeknBiPWSSHuCpvHWPEfcqVOnbG6vra1FTU2NsK87DFDMv8IM0s86+xO1OohxcxSXZicQ3JWRUM5st0d+F1Y/mF3M2pOgcvl9mBzEwKl/ssm/ROUgcUhB4I6utaPra7dcTm6z5/qAhSabebjyi91ReazrmmsDGDype444ioH0JC9vsq6X9gcwyAMc9NL/On9cbw5YoFaFPW/BszB9dnY2AGD//v02tzc/361bN440JSIiopAVNI2322+/HUlJSSgtLcXevXuttq9btw4AGB9EREQUwppHm8p/ocTvt02nTJmCL7/8ErNmzcKsWbMsBQkLw9y5c/HMM89g5syZ+Oyzz9C2bVsAQFFREVatWoXIyEjMnj3b30UmIgpJ94/Yh/aZ1+901JssM9bJUyHJGk2WSPgmqY/AlTVxY/TXhHRm1E/mxykGcdYBeT1pe3MSWq8/a0nLayrbm7TXk7AR+Tj2pnwy2Zm0NwqNwrb+UeISkgv7XTA/Lrk5VdimvqaONEqjG9QhF+F6+9dUvZa0+pqeK6nFSncGLCD0bpPKPGq8bdmyBUuWLDGnz549CwCYPHkyoqOjAQA5OTkoKCgw71NeXo5jx46hsrLSKr/8/Hzs2LEDRUVFyMzMxIgRI1BdXY3i4mIoioK3334b6enpnhRZmAxVJseieMKVWCKHHywe/qRQx9BYx5mJXBk95ehYrrI3MaoVF+OM5C8QT660vfg4d2OLrsdiXs9X/mL0pF66GqvpKMbNlzFvMldi4HwZ4yZzeA7kho3VJLzSOfbD9fbndSMi3/Oo8VZRUYF9+/ZZPX/48GHz4yhp+Q17wsPDsXXrVixfvhxr1qzB5s2bERUVhVGjRmH+/PkYOnSoJ8UlIiIijdMpCnTSoDU53dp51HjLy8tDXl6eS68pLi62uz08PBz5+fnIz8/3pGhERERErVLQTBVCRETBJTu6DF1iIwHp1nSjYv+rw/60Sy33kMghLdZL5rW8/Jy9ZQUdhZ6opymyWsZPik1z9w69VT52yiCzO2mvg1CNjmGWmLeUWDFOUL5O9q6ru9cU0nVT71sa3WD3dS3iVCHBM9qUiIiIiBwLuZ439WSoJnkiV3nxcg/Igd8y9a8PbwcT2xst5SgI3NXRU/Z+jTn6pekqu+fJKuhbeq0P4yHU19rdAQt61eTRMrleevIL2Oq4HkzCa9Ub4OKxXRnVBw8HJbhaF60n4rU8NDnIy9EkvtYDGOTXe/962xvJSKQ1XJiePW9EREREmhJyPW9ERESkYYx5Y+ONiIhs0+lMllvr6tu9dgYOAEC4m8dzNO+gvX3tsbrVLSXVeckDBxwOYHCTvQEK8q1yR7fq1eytlWt1PuVb+PYGhDhdAvtlUpdB5yC8qCW8bRqCjTchtsjNiuMMTxaX9zSWyBX+nNzUEXtxRpA+wFw+J3661u7GFunVX5KO2ImZcjqPGzypp3Sdo4XlHcbAOeKF6+1qvSCi4BZyjTciIiLSuBD/XckBC0REREQawp43IiJySLj16mC6E3enN3FliiVPprmR2V//2T+hJe6GsLgyNY08kZHDGDgny+fo3Hv9tr2NmLdQ64lj440oCNib503+0HTlg9DVLzTHH8I+jL90NH+anWN7e05Be+RzIAeTuxwDJ/HF9eY8b0StS8g13gx+CsL2ZLJTb1N/sTmaCNWbx/I29Xly9IXpT+KEy0RE5FOcKoQxb0RERERaEnI9b0RE5BwDFBhsLZrnoIPd3V4Bm8cir3F4fu1cV096elo6rrt3TDjPGxtvREREpCWKcv1Pfi6EhHTjzZ9xUq7Etfl0Ul4/BnZ7yv5M4fbPUTDFxDlDHYvpaoC7a8cJ3HlwNHO9I96su/4aQQg3Bij44nozFpOodQnpxhsRERFpi87GbVLtdEt4BwcsEBEREWkIe96IiMgm9fyD4kS29gPfm+z0C7gyKMFXt/nlsAp1GEYwrPfsykL0rnA8B2HL18YX19Tt+Qc5VQgbb0TBJpji83w5B6GnMXCeHMubHE3a60gwXW8i0oaQa7zpdYrPJ8V1lacf3v78EvSEp1+g6vPk6NdxoK5xsNUtIqLWRmeyXgXM2ytwBTvGvBERERFpCBtvRETkkOFG/Jsvb/Oqj+HP28kmRWf+C5RAlMFf59vrx1Ba+PPAgQMH8Oqrr2LixIno2rUrdDoddDodjhw54lI+jY2N2LZtG+bMmYP+/fsjLi4OkZGR6Nq1K6ZNm4ajR496VtAbQu62KREREWmXL1ZYWLx4MTZu3OhZJgB27tyJUaNGAQA6deqE3NxcGAwGHDhwAO+88w7effddvPfee7j//vs9Og4bbwHiy185/pyANFAcjZoi5wRTjF5rqbeuTiBNRIE3ZMgQ9OvXDwMGDMDAgQNx5513oqyszOV89Ho9JkyYgLlz52Lw4MHm55uamvDMM8/glVdewdSpU3HixAkkJSW5XV423oiIiEhDbCyP5WGHyLx58zx6fbORI0di5MiRVs8bDAa89NJL2LBhA44dO4YtW7ZgypQpbh+HMW9ERORQE3TmP38cw9Zx9DqT+c9X1LFnvuwl9ddx7J0zR+fbW/xxDC3Q6XTIysoCAJw+fdqjvNjzRkRERJrhi5g3fykpKQEApKamepRPyDXe9G6OePHl4tHkOU+vh7eur7tz7AXj/IPkG65eZ3WvjLt1hHWLQkVz40gtOTkZKSkpASmPWlFREQ4ePIjIyEiMHj3ao7xCrvFGREREGmZneaxx48ZZ7b5w4UIsWrTIP2VrQWVlJaZNmwYAePLJJ9GhQweP8mPjjYiIiFqFDRs2IDMzU3guOTk5YOUBgKtXr2L8+PEoLy/H0KFDsXDhQo/zZOONiIg8Ym/Rckf7urKouUmxvNZXgxb8dYvZ0zVxnaU+Z67w1zV1i42Yt+aet8zMTPTp08e3x3eB0WjExIkTsWvXLuTk5GDTpk0IDw/3OF823oiIiEg7FBtThVhNHRJ4TU1NmDx5Mj755BP06tULhYWFaNOmjVfyZuPNSRyg0Lrx+nrO0TQAwXyOg3kKAw42INIeRVEwbdo0vP/+++jevTuKiopw0003eS1/Nt6IiIhIM7QwVcisWbOwZs0apKenY/v27UhLS/Nq/pykl4iIbDK1MLlqE/TCnye8lY/MlcXQm6fqCWQvpytl8NVi8r64pnJe6gl7tb4k3pQpU9CzZ0+sXLlSeP6pp55CQUEB0tLSsH37dqSnp3v92Ox5IyIiIm3xcjt7y5YtWLJkiTl99uxZAMDkyZMRHR0NAMjJyUFBQYF5n/Lychw7dgyVlZXm5zZt2oSlS5cCALp16ybkqXbnnXdi+vTpbpeXjTci8otgjisjotBWUVGBffv2WT1/+PBh8+OoqCiH+VRXV5sf7969G7t3725xXzbeiIiIKCT4IuYtLy8PeXl5Lr2muLjYK/m4g403IiLyKntzi7kyP5uvliUM9hG8vpoDzpPeb29dU/IONt6IiIhIO0zK9T/5uRDCxhsRERFph521TUMFG28hIFCB4sE8KSsREZFWsfFGREREmqGFSXp9jY03IiJyyN7Era4sfi7vqw5299UC544GKATbXQL5bonvBjD4/prKx/H5ovUhgo03IiIi0pYgXIjen9h4C5BQmLDUn+8x2H45ExER+Qobb0RERKQZjHlj442IiEJMsPfUy+ULhTs15Bo23oiIiEg7OM8bG29ERESkHTpFgU4asCCnW7uQa7yZoGuxC9qbXens5vYvR+fbX9fWxOtOREQ+FnKNNyIi8i55Xq/WtIi53kf344L9h15QX1PTjT/5uRDi/Cx8RERERBRw7HkjIiIiDbGOeQu1EQtsvKm0lji1YB0GH8jz21quLRERERtvREREpB2cKoSNNyIisk0PRdWT33JEuLzAubsB7PKi5YG4i+CrAQqOjuOPAQzW57PlBeRlrlxTe4vPq8vgr3PdGrHxRkRERNqhKNYL04fYPG8cbUpERESkIR433oxGI1577TVkZWUhJiYGSUlJGDt2LHbt2uVyXnl5edDpdC3+jR492tPitgqGG7cyWvoLVlotNxERBY/mhenlv1Di0W3TxsZGjB07FkVFRUhKSsJ9992HqqoqFBYWorCwEKtXr8aUKVNczveOO+5AZmam1fP9+vXzpLhEROQD9mKcvMlbk8HKo8/VPx4dxZ65G6flSkybN0fHq8+ZvYl24cfr6DHeNvWs8bZ06VIUFRUhOzsbn332GRITEwEARUVFGDNmDGbMmIHhw4cjPT3dpXynT5+OvLw8T4pGRERE1Cq5fdvUaDRi2bJlAICCggJzww0AcnNzMW3aNDQ0NGDFihXeKSkRERGFPJ0C6EzSX2h1vLnfeNuzZw+qqqrQpUsXDB482Gr7pEmTAAAbN270rIQUMrFhofI+iYiIPOH2bdODBw8CAAYMGGBze/PzJ06cQG1tLeLj453Oe8eOHfj6669RX1+P1NRUDB8+HCNGjHC3qERERNRaMObN/cZbWVkZAKBz5842t8fHxyMhIQE1NTUoKytD3759nc577dq1Qnrx4sUYNGgQ3nvvPWRkZDh8/fnz51FRUSE8V1JS4vTxiXyB9ZKCkb16qdeZLAHvDoLd1VwJuPdHL7tJEcuj9+Aem3rggaPBC96aeFcuv7fYO/e+uobqARTeGoASitxuvF2+fBkAEBsb2+I+cXFxqKmpQW1trVN59u/fH7feeivuvvtupKeno7q6Gl988QUWLFiAffv2ITc3FwcPHkRcXJzdfAoKCvD888+7+I6IfIv1koIR6yVpDpfHCq4VFubMmSOkY2JiMGHCBNx7773IyclBSUkJ3nzzTeTn59vN57HHHsOECROE50pKSjBu3DiflNvbGO91nXwetL64vNbrJbVOrJdE2uN2462596uurq7FfZp751yJd7MlISEBs2fPxhNPPIGtW7c6bLylpKQgJSXFo2MSeRvrJQUj1kvSmuuT8ipWz4UStxtvzbFnp06dsrm9trYWNTU1wr6e6NGjBwDgzJkzHudFRESucW2yV/e+Sb0ZA6XuqffV3Qx/LCZvi7t3IRydX/V19eY583psGwcsuD9VSHZ2NgBg//79Nrc3P9+tWzePe94AoLq6GlD1+BERERGFIrcbb7fffjuSkpJQWlqKvXv3Wm1ft24dAHgtbmL9+vUAgIEDB3olPyIiItIgBYBJ+gutjjf3G29hYWGYO3cuAGDmzJm4ePGieVtRURFWrVqFyMhIzJ49W3jdlClT0LNnT6xcuVJ4/tChQ9iyZQuampqE569cuYIFCxbgo48+gsFgwMyZM90tclAJlgXa9VA8+vOnYDlnREREgeTRaNP8/Hzs2LEDRUVFyMzMxIgRI1BdXY3i4mIoioK3337bal3T8vJyHDt2DJWVlcLzpaWleOCBB5CUlIScnBwkJyejoqIChw4dQkVFBSIiIvDWW28hKyvLkyITEZGHWtv8XP6Ij3OlDIGileuqUxQbAxZC60e8R4238PBwbN26FcuXL8eaNWuwefNmREVFYdSoUZg/fz6GDh3qdF5ZWVl44okn8NVXX+HIkSOoqqqCwWBAeno6xo8fj8cffxy9e/f2pLhEREREmufxPG/h4eHIz893OH1Hs+LiYpvPd+vWjYvYExERkX0cbRpck/S2Zv7sivdnLJorx/LlsPrWNqEvefZ/htefiFozNt6IiIhIO9jzxsYbERGFFm8uVO8LvlqIvtVonh5Efi6EuD1VCBERERH5H3veiIiISDN0sDFVSIjN9xlyjTe9nUldPQlyDuREsf6eLNddcjn9OYBB5qtrrZVrEQj+/D/i6rE4wIGItCTkGm9ERNT6yA1wdQPe3jZH+cjc/RHiyQ8Ee68NyR8eio0BCiH2u5kxb0REREQawp43IiIi0g5OFcLGm1qwLnDeWuOoHL2vQMbEkeu0fE59GSNJRORtbLwREZFDJsX/UTb+WihdPa+aoznffNWQD8Tcbpq9pux5Y+ONiIiINIST9HLAAhEREZGWsOeNiIiINMR6kt5QmyuEjbcg1FoHKLjKn5P6Bpq9yaMp8LR+bfiZQtS6sPFGREQ2mRR9QILa1cdXcyXYXT2wQG5825u011+L1tsboOBoUIS7gyYCeS1tlcHt8nDAAmPeiIiIiLSEPW9ERESkHSbl+p/8XAhh4y0IMB7FOerz1Jrj34iIiOxh442IiDTBXoyUvXg4RwvT24uP81YMnKNJeL21+HwwxLX5HBemZ+ONiIiINIQDFjhggYiIiEhL2HgjIiIiDVEsvW/mXjjPet4OHDiAV199FRMnTkTXrl2h0+mg0+lw5MgRt/M0Go147bXXkJWVhZiYGCQlJWHs2LHYtWuXR2VFKN421UNpcYCAN4PgOQjBt1w9v/66trzuFAp8tTi7J5MhuzInnL0YOEfxcd4SjDFuwXhd/WXx4sXYuHGj1/JrbGzE2LFjUVRUhKSkJNx3332oqqpCYWEhCgsLsXr1akyZMsXt/EOu8UZEREQa5oOpQoYMGYJ+/fphwIABGDhwIO68806UlZW5nd/SpUtRVFSE7OxsfPbZZ0hMTAQAFBUVYcyYMZgxYwaGDx+O9PR0t/Jn442IiIhC2rx587yWl9FoxLJlywAABQUF5oYbAOTm5mLatGl46623sGLFCrz++utuHYMxb0RERKQdisn2X5DYs2cPqqqq0KVLFwwePNhq+6RJkwDAo9u07HlTYbxS68Vrqz2eXDNO4kxEgXLw4EEAwIABA2xub37+xIkTqK2tRXx8vMvHYOONiIgc8lUwu71jeHMAg5o8mMHeJL3W+Vr2dTRhr6OJeVsqg3U+vrlJ5o9rKh/HK4MX7EzSW1JSYrV7cnIyUlJSPD+uk5pj5Tp37mxze3x8PBISElBTU4OysjL07dvX5WOw8UZERETaodgYsHCjMTdu3Dir3RcuXIhFixb5q3S4fPkyACA2NrbFfeLi4lBTU4Pa2lq3jsHGGxEREbUKGzZsQGZmpvBccnJywMrjK2y8ERERkXbYWR4rMzMTffr0CUy5boiLiwMA1NXVtbhPc++cO/Fu4GhTIiIiIu/JyMgAAJw6dcrm9traWtTU1Aj7uoo9b0RERKQdQb4wfXZ2NgBg//79Nrc3P9+tWzf2vBEREREF2u23346kpCSUlpZi7969VtvXrVsHtDC4wllsvBEREZF2yIvS2+qJ84MpU6agZ8+eWLlypfB8WFgY5s6dCwCYOXMmLl68aN5WVFSEVatWITIyErNnz3b72LxtqgHyhKOtdcLZUHmf5BxOtBt4JujcmgusyU6/gAHOz4Tvr0XiQ4Un87r54poG0//xLVu2YMmSJeb02bNnAQCTJ09GdHQ0ACAnJwcFBQXmfcrLy3Hs2DFUVlZa5Zefn48dO3agqKgImZmZGDFiBKqrq1FcXAxFUfD222+7va4p2HgjIiIiTTEpgMlk/ZwHKioqsG/fPqvnDx8+bH4cFRXldH7h4eHYunUrli9fjjVr1mDz5s2IiorCqFGjMH/+fAwdOtSj8rLxRkRERBpi6zapZ423vLw85OXlufSa4uJiu9vDw8ORn5+P/Px8j8pmC2PeiIiIiDSEPW9ByFEcQDDFCfgSY+CIiMhKkE8V4g9svBERkVe5EsDuD/JC9DJXBkI4Woy+pX0dLVKvLoM8sEAuv68Wqrcn2K5pqGPjjYiIiLTDZGNheg8HLGgNY96IiIiINIQ9b0RERKQdigmKYrJ6LpSw8RYEQmUAgqfU54mDF/zDk0k9tTyhaqi+byLSBjbeiIiISDsUGzFuIfabiY03IiIi0g5OFcIBC0RERERawp63AGGcm2c4gS8RUYgymQB57j55rdNWjj1vRERERBrCnjciIiLSDsa8seeNiIiISEvY80ZERESaoSgKFCnGTWHPGxEREREFK/a8ERERkXYw5o09b0RERERawp43IiIi0g6TYr0elrxcVivHxhsRERFph2K6/ic/F0J425SIiIhIQzxuvBmNRrz22mvIyspCTEwMkpKSMHbsWOzatSso8iMiIqLWQ1EAxaSIf6F119SzxltjYyPGjBmD/Px8nDlzBvfddx9+9rOfobCwECNGjMDatWsDmh8RERFRa+NR423p0qUoKipCdnY2jh8/jvXr12P79u0oLCyEXq/HjBkzUF5eHrD8iIiIqJVpjnmT/0KI2403o9GIZcuWAQAKCgqQmJho3pabm4tp06ahoaEBK1asCEh+RERERK2R2423PXv2oKqqCl26dMHgwYOttk+aNAkAsHHjxoDkR0RERK2PVbzbjb9Q4nbj7eDBgwCAAQMG2Nze/PyJEydQW1vr9/yIiIiIWiO353krKysDAHTu3Nnm9vj4eCQkJKCmpgZlZWXo27ev3/I7f/48KioqhOe+++47AEB5qdHBO/MPE3SBLkKropcnbAyQ5vrV0NBgtU0L9VLmST0Nlmvijtb2vt2tl/8ubfRTCZ3j7rnV6+zHQ9nL12CVl3tlMClinWqyt6+D+mdS3Ot3Cbbvneb6Zate2lOv1FrFuNWjzqtlC3ZuN94uX74MAIiNjW1xn7i4ONTU1DjVU+bN/AoKCvD888/b3PbEwxccloXIU9988w1ycnKE5+zVyzkPV/upZBTKXK2Xzzxyzk8lo1Bmq17akpiYiPj4eByu3WNze3x8vBAv35q1yhUWHnvsMUyYMEF47tChQ/jNb36D999/H7179w5Y2bSkpKQE48aNw4YNG5CZmRno4mjCd999h4kTJ+KWW26x2sZ66R2sl65jvfQ91kvX2auXtqSlpeH777/HhQu2O2ESExORlpbm5VIGJ7cbb3FxcQCAurqWuyqbe9Pi4+P9ml9KSgpSUlJsbuvduzf69OnjsDxkkZmZyXPmooSEBKvnWC+9i/XSdayXvsd66Tpb9bIlaWlpIdNAs8ftAQsZGRkAgFOnTtncXltbi5qaGmFff+ZHRERE1Bq53XjLzs4GAOzfv9/m9ubnu3Xr5lTPm7fzIyIiImqN3G683X777UhKSkJpaSn27t1rtX3dunUAgHHjxgUkPyIiIqLWyO3GW1hYGObOnQsAmDlzJi5evGjeVlRUhFWrViEyMhKzZ88WXjdlyhT07NkTK1eu9Ep+zkpOTsbChQuRnJzs1utDEc+Z61w9ZzzHruM5cx3rpe/xnLmO58x9OkVR3J6UqLGxEWPHjkVRURGSkpIwYsQIVFdXo7i4GIqiYPXq1Zg6darwmuHDh2Pnzp1YuHAhFi1a5HF+RERERKHEo4Xpw8PDsXXrVrz66qtITU3F5s2bceDAAYwaNQrFxcUuN7S8nR8RERFRa+NRzxsRERER+ZdHPW9ERERE5F9svBERERFpCBtvRERERBqiycab0WjEa6+9hqysLMTExCApKQljx47Frl27giK/YOTN95iXlwedTtfi3+jRo33yHvzpwIEDePXVVzFx4kR07drV/N6OHDnS4mscnWNXrwHrpWtYL21jvXQd66Vr3KmXjoRCPfOIojHXrl1TcnNzFQBKUlKS8stf/lIZMWKEotfrFb1er6xZsyag+QUjb7/HqVOnKgCUO+64Q5k6darV32uvveaz9+Iv999/vwLA6u+bb76xub+jc7xq1SqXrgHrJeulLayXvsd66TpX66UjoVDPPKW5xtv//M//KACU7Oxspbq62vz8tm3blLCwMCUyMlIpKysLWH7ByNvvsfnDaPXq1T4qceC9/PLLyrPPPqt8/PHHyqlTp5SMjAy7H0aOznFYWJhL14D1kvXSFtZL32O9dJ2r9dKRUKhnntJU462xsVFJSkpSAChffPGF1fZHHnlEAaDMnTs3IPkFI1+8x1D4MJLZ+zBydI4ffvhh8y9RZ64B6yXrpbNYL72L9dI7PGm8hUI98wZNxbzt2bMHVVVV6NKlCwYPHmy1fdKkSQCAjRs3BiS/YBQK7zHQHJ3jPn36ADeWgHPmGoTCNQuF9xhorJeuC4X3GOx4DZwTFugCuOLgwYMAgAEDBtjc3vz8iRMnUFtbi/j4eL/mF4x8+R537NiBr7/+GvX19UhNTcXw4cMxYsQIL5VcOxyd46tXrwI3AnBtnWP5GrBesl56A+ul61gvAy8U6pk3aKrxVlZWBgDo3Lmzze3x8fFISEhATU0NysrK0LdvX7/mF4x8+R7Xrl0rpBcvXoxBgwbhvffeQ0ZGhocl1w5H5/jcuXPCvvI5lq8B6yXrpTewXrqO9TLwQqGeeYOmbptevnwZABAbG9viPnFxcQCA2tpav+cXjHzxHvv374+VK1fi6NGjqKurw6lTp/D+++8jMzMT+/btQ25urvm4ocDROVafi5bOsfoasF5ex3rpGdZL17FeBl4o1DNv0FTPGwWHOXPmCOmYmBhMmDAB9957L3JyclBSUoI333wT+fn5ASsjhR7WSwpGrJfkC5rqeWtubdfV1bW4T3Or3Zn74N7OLxj58z0mJCRg9uzZAICtW7d6lJeWODrHzdth5xyrrwHr5XWsl55hvXQd62XghUI98wZNNd6a4wJOnTplc3ttbS1qamqEff2ZXzDy93vs0aMHAODMmTMe56UVjs5x+/btrfZVk68B6yXrpTewXrqO9TLwQqGeeYOmGm/Z2dkAgP3799vc3vx8t27dnGqRezu/YOTv91hdXQ1Iv+pbO0fnOCoqCrgxJYOtcyxfA9ZL1ktvYL10Hetl4IVCPfMGTTXebr/9diQlJaG0tBR79+612r5u3ToAwLhx4wKSXzDy93tcv349AGDgwIFeyU8LHJ3jb7/9FrgxJYMz14D1kvXSG1gvXcd6GXihUM+8ItCzBLuqedmMnJwc5cKFC+bn7S2b8dBDDyk9evRQ/vCHP3glP63x5jk7ePCgsnnzZsVoNArP19fXK/Pnz1cAKAaDQfn66699/K78y9lliCIjI5VXXnnF/Ly8DJH6Gjz00ENK586dFb1e3+IyRKyXrJf2sF56H+ul55xZYSHUv5c9pbnGm60Fa0eOHKno9XpFp9Mpf/nLX6xeM2zYMAWAsnDhQq/kpzXePGcff/yxOZ977rlHefDBB5V77rlHSU5OVgAoERERrWIZmM2bNyuDBg0y/0VERCgAlKysLPNzv//97837X7t2TWnbtq0CQImOjrY6x2+//bbVNWjeH4DVNWC9ZL20hfXS91gvXedqvVT4vewxzTXelBsX9tVXX1X69OmjREVFKW3btlVGjx6t7Ny50+b+9iqJO/lpkbfO2YkTJ5QnnnhCGTJkiNKhQwclIiJCiY6OVnr06KE8+uijyrfffuund+Rbq1evNn+BtfQ3bNgw4TV33XWXAkBJTk62eY7la2AwGBQASl5ens0ysF5aY71kvfQH1kvXuFMv+b3sGZ2iKEqgb90SERERkXM0NWCBiIiIKNSx8UZERESkIWy8EREREWkIG29EREREGsLGGxEREZGGsPFGREREpCFsvBERERFpCBtvRERERBrCxhsRERGRhrDxRkRERKQhbLwRERERaQgbb0REREQawsYbERERkYaw8UZERESkIf9/uYB9JVObZP8AAAAASUVORK5CYII=", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "fig = plt.figure(figsize=(7, 5))\n", + "grid = ImageGrid(fig, 111, nrows_ncols=(1, len(runs)), axes_pad=0.1,\n", + " share_all=True, cbar_mode=\"single\", cbar_location=\"right\")\n", + "\n", + "for ax, s, p in zip(grid, solvers, runs):\n", + " rho = p.get_var(\"density\")\n", + " g = p.get_grid()\n", + " im = ax.imshow(rho.T,\n", + " extent=[g.xmin, g.xmax, g.ymin, g.ymax],\n", + " vmin=0.9, vmax=2.1)\n", + " ax.set_title(s, fontsize=\"small\")\n", + "grid.cbar_axes[0].colorbar(im)" + ] + }, + { + "cell_type": "markdown", + "id": "c1843160-9953-4f58-835c-ba5254e12afe", + "metadata": {}, + "source": [ + "We see that the 4th order solver has more structure at fine scales (this becomes even more apparently at higher resolutions)." + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.12.5" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/docs/source/index.rst b/docs/source/index.rst index 32224e92d..ff23e0c99 100644 --- a/docs/source/index.rst +++ b/docs/source/index.rst @@ -62,6 +62,7 @@ new ideas. :caption: Examples :hidden: + compressible-rt-compare.ipynb advection-error.ipynb compressible-convergence.ipynb diff --git a/pyro/pyro_sim.py b/pyro/pyro_sim.py index 7dde5ae4f..1d1481c4a 100755 --- a/pyro/pyro_sim.py +++ b/pyro/pyro_sim.py @@ -239,7 +239,10 @@ def single_step(self): tm_vis.end() def __repr__(self): - """ Return a representation of the Pyro object """ + s = f"Pyro('{self.solver_name}')" + return s + + def __str__(self): s = f"Solver = {self.solver_name}\n" if self.is_initialized: s += f"Problem = {self.sim.problem_name}\n"