{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"\n",
"# Standalone Fishbone-Moncrief C Code\n",
"\n",
"We start with the NRPy+ expressions generated in the [Tutorial-FishboneMoncriefID](Tutorial-FishboneMoncriefID.ipynb), and output them to the C file \"FishboneMoncriefID/FMstandalone.h\".\n",
"\n",
"Further, $\\Gamma = \\alpha u^0$ is given by (as shown [here](Tutorial-u0_smallb_Poynting-Cartesian.ipynb)):\n",
"$$\n",
"\\Gamma = \\alpha u^0 = \\sqrt{\\frac{1}{1 - \\gamma_{ij}v^i_{(n)}v^j_{(n)}}}.\n",
"$$"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"execution": {
"iopub.execute_input": "2020-12-17T02:14:18.841310Z",
"iopub.status.busy": "2020-12-17T02:14:18.837158Z",
"iopub.status.idle": "2020-12-17T02:14:24.365260Z",
"shell.execute_reply": "2020-12-17T02:14:24.364668Z"
},
"scrolled": true
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Wrote to file \"FishboneMoncriefID/FM_standalone.h\"\n"
]
}
],
"source": [
"import sympy as sp # SymPy: The Python computer algebra package upon which NRPy+ depends\n",
"from outputC import lhrh # NRPy+: Core C code output module\n",
"import indexedexp as ixp # NRPy+: Symbolic indexed expression (e.g., tensors, vectors, etc.) support\n",
"import finite_difference as fin # NRPy+: Finite difference C code generation module\n",
"import FishboneMoncriefID.FishboneMoncriefID as fmid\n",
"\n",
"# Step 1: Set up the Fishbone-Moncrief initial data. This sets all the ID gridfunctions.\n",
"fmid.FishboneMoncriefID(\"Spherical\")\n",
"\n",
"gammaDD = ixp.zerorank2()\n",
"\n",
"DIM = 3\n",
"for i in range(DIM):\n",
" for j in range(DIM):\n",
" if i<=j:\n",
" gammaDD[i][j] = fmid.IDgammaDD[i][j]\n",
" else:\n",
" gammaDD[i][j] = fmid.IDgammaDD[j][i]\n",
"\n",
"# gamma_{ij} v^i_{(n)} v^j_{(n)}\n",
"Gammacontraction = sp.sympify(0)\n",
"for i in range(DIM):\n",
" for j in range(DIM):\n",
" Gammacontraction += gammaDD[i][j] * fmid.IDValencia3velocityU[i] * fmid.IDValencia3velocityU[j]\n",
"\n",
"Gammafactor = sp.sqrt(1 / (1 - Gammacontraction))\n",
"\n",
"# -={ F-M quantities: Generate C code from expressions and output to file }=-\n",
"FishboneMoncrief_to_print = [\\\n",
" lhrh(lhs=\"alpha\",rhs=fmid.IDalpha),\\\n",
" lhrh(lhs=\"betaU0\",rhs=fmid.IDbetaU[0]),\\\n",
" lhrh(lhs=\"betaU1\",rhs=fmid.IDbetaU[1]),\\\n",
" lhrh(lhs=\"betaU2\",rhs=fmid.IDbetaU[2]),\\\n",
" lhrh(lhs=\"Gammafactor\",rhs=Gammafactor),\\\n",
" lhrh(lhs=\"Gamma_times_ValenciavU0\",rhs=Gammafactor*fmid.IDValencia3velocityU[0]),\\\n",
" lhrh(lhs=\"Gamma_times_ValenciavU1\",rhs=Gammafactor*fmid.IDValencia3velocityU[1]),\\\n",
" lhrh(lhs=\"Gamma_times_ValenciavU2\",rhs=Gammafactor*fmid.IDValencia3velocityU[2]),\\\n",
" lhrh(lhs=\"uKS4U1\",rhs=fmid.uKS4U[1]),\\\n",
" lhrh(lhs=\"uKS4U2\",rhs=fmid.uKS4U[2]),\\\n",
" lhrh(lhs=\"uKS4U3\",rhs=fmid.uKS4U[3]),\\\n",
" lhrh(lhs=\"uBL4U1\",rhs=fmid.uBL4U[1]),\\\n",
" lhrh(lhs=\"uBL4U2\",rhs=fmid.uBL4U[2]),\\\n",
" lhrh(lhs=\"uBL4U3\",rhs=fmid.uBL4U[3])\n",
" ]\n",
"# print(fmid.uKS4U[3])\n",
"fin.FD_outputC(\"FishboneMoncriefID/FM_standalone.h\",FishboneMoncrief_to_print,params=\"outCverbose=False,CSE_enable=False\")"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"execution": {
"iopub.execute_input": "2020-12-17T02:14:24.369070Z",
"iopub.status.busy": "2020-12-17T02:14:24.368440Z",
"iopub.status.idle": "2020-12-17T02:14:24.370477Z",
"shell.execute_reply": "2020-12-17T02:14:24.370973Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Writing FishboneMoncriefID/FM_standalone.c\n"
]
}
],
"source": [
"%%writefile FishboneMoncriefID/FM_standalone.c\n",
"\n",
"#include \"stdio.h\"\n",
"#include \"stdlib.h\"\n",
"#include \"math.h\"\n",
"\n",
"const double a = 0.9375;\n",
"const double M = 1.0;\n",
"const double r_at_max_density = 12.0;\n",
"const double r_in = 6.0;\n",
"\n",
"int main(int argc, const char *argv[]) {\n",
"\n",
" // Step 0a: Read command-line input, error out if nonconformant\n",
" double xx0,xx1,xx2;\n",
"/*\n",
" if(argc != 4) {\n",
" printf(\"Error: Expected three command-line arguments: ./FM_standalone r theta phi\\n\");\n",
" exit(1);\n",
" }\n",
" xx0 = strtod(argv[1],NULL);\n",
" xx1 = strtod(argv[2],NULL);\n",
" xx2 = strtod(argv[3],NULL);\n",
"*/\n",
"\n",
"// printf(\"# Output: r,th,ph, alpha, betaU0, betaU1, betaU2, Gamma, Gamma*vValenciaU0, Gamma*vValenciaU1, Gamma*vValenciaU2\\n\");\n",
" for(double xx0=1.6;xx0<50.0;xx0+=0.2) {\n",
" xx1 = 1.56463634120e0; //M_PI/2.0;\n",
" xx2 = 0.0;\n",
" double alpha,betaU0,betaU1,betaU2,Gammafactor,Gamma_times_ValenciavU0,Gamma_times_ValenciavU1,Gamma_times_ValenciavU2;\n",
" double uKS4U1,uKS4U2,uKS4U3,uBL4U1,uBL4U2,uBL4U3;\n",
"#include \"FM_standalone.h\"\n",
" if(xx0 < r_in) {\n",
" Gammafactor = 1.0;\n",
" Gamma_times_ValenciavU0 = Gamma_times_ValenciavU1 = Gamma_times_ValenciavU2 = 0.0;\n",
" uKS4U1 = uKS4U2 = uKS4U3 = 0.0;\n",
" uBL4U1 = uBL4U2 = uBL4U3 = 0.0;\n",
" }\n",
" printf(\"%e %e %e %.15e %.15e %.15e %.15e %.15e %.15e %.15e %.15e %.15e %.15e %.15e %.15e %.15e %.15e\\n\",\n",
" xx0,xx1,xx2,\n",
" alpha,betaU0,betaU1,betaU2,\n",
" Gammafactor,\n",
" Gamma_times_ValenciavU0, // util1(1) in FMtorus.f90; util(1,i,j,k) near the write statement\n",
" Gamma_times_ValenciavU1, // util1(3) in FMtorus.f90.\n",
" Gamma_times_ValenciavU2, // util1(2) in FMtorus.f90.\n",
" uKS4U1,uKS4U2,uKS4U3,\n",
" uBL4U1,uBL4U2,uBL4U3);\n",
" }\n",
" return 0;\n",
"}"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"execution": {
"iopub.execute_input": "2020-12-17T02:14:24.376108Z",
"iopub.status.busy": "2020-12-17T02:14:24.375824Z",
"iopub.status.idle": "2020-12-17T02:14:24.733826Z",
"shell.execute_reply": "2020-12-17T02:14:24.734257Z"
}
},
"outputs": [],
"source": [
"!gcc -O2 FishboneMoncriefID/FM_standalone.c -o FM_standalone -lm"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"execution": {
"iopub.execute_input": "2020-12-17T02:14:24.739653Z",
"iopub.status.busy": "2020-12-17T02:14:24.737603Z",
"iopub.status.idle": "2020-12-17T02:14:24.849288Z",
"shell.execute_reply": "2020-12-17T02:14:24.849011Z"
}
},
"outputs": [],
"source": [
"!./FM_standalone > out.txt"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"execution": {
"iopub.execute_input": "2020-12-17T02:14:24.864880Z",
"iopub.status.busy": "2020-12-17T02:14:24.861157Z",
"iopub.status.idle": "2020-12-17T02:14:25.608280Z",
"shell.execute_reply": "2020-12-17T02:14:25.607793Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"--2020-12-16 21:14:25-- http://astro.phys.wvu.edu/zetienne/torus_cuts.csv\n",
"Resolving astro.phys.wvu.edu (astro.phys.wvu.edu)... 157.182.3.45\n",
"Connecting to astro.phys.wvu.edu (astro.phys.wvu.edu)|157.182.3.45|:80... connected.\n",
"HTTP request sent, awaiting response... 200 OK\n",
"Length: 16445 (16K) [text/csv]\n",
"Saving to: ‘torus_cuts.csv’\n",
"\n",
"torus_cuts.csv 100%[===================>] 16.06K --.-KB/s in 0.01s \n",
"\n",
"2020-12-16 21:14:25 (1.40 MB/s) - ‘torus_cuts.csv’ saved [16445/16445]\n",
"\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAZ0AAAEWCAYAAAC9qEq5AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABPVklEQVR4nO3dd3xV9fnA8c+TnZCEQAhhB2QP2SBWUMEFto66N/ZnReuqda+2bm3raJ1VW0WsojhQFBBRZKgoMkQQlI2AbAgJZOc+vz/OueEmZFzgruQ+79frvpJ7vmc85+bmPOf7Pd/zPaKqGGOMMaEQE+4AjDHGRA9LOsYYY0LGko4xxpiQsaRjjDEmZCzpGGOMCRlLOsYYY0LGkk6UEZF7ReR/Ydz+DyJyfC3lM0Xk96GLKHKI4xUR2S0i88IdjwkvEbleRH4RkcUBXOdiEdksIn8K1DoPVr1MOiKyTkS2iUgjn2m/F5GZPu9VRPaJyF4R2SQiT4hIrE/5TBEpcst3iMh7ItLSj23fKyL3HkLMvvHsFZFcd/rxbtnEKvP3cafPrG599ZWq9lTVmXD4CdD97Dw+n+leEfnQZ90qIn+ssswf3en3+qxj46HvUUANBU4C2qjqYKj4Xq9y9+1jEWnlndlNUn8TkZ3u628iIv5sKFL22/1fPvEQlxURuU5EvheRAhHZ4v5fXxDoOMPkXuAaVe3jnVD1f0ZEWovIjyLylPt59BSRT0Rkl4jkisgCETnVO7+7riuBvx5OYCLSWUTeFJHtIpInIitF5GkRaVPXsvUy6bhigT/WMU8fVU0FjgPOB/6vSvl1bnkXIAN4MtBBVheP+8rwmb4dOFpEMn2mjQZWBDmehuAXn880VVVP8ylbAVxWZf6Afq4iEheodQE5wDpV3eeu+3jgYeAMoCmwFhjvM/8Y4EygD9AbOA24KoDx1Mn3RC4MngJuBG4GMoHWwD3AyDDGFEhNgaU1FYpIDjAbmKSqN6hzp/+HwHSgBdAcuAHIq7LoUqDxof7tRKQT8A3wC9BPVdOBY4DVOCdOtarPSecfwC0iklHXjKq6CvgS6FtD+S7gXaCXiAwSka1VakVnVVfFFZFmIvKRe0axS0TmiMihfKYlwPvABe56Y3GS5Os1LSAiSSLyP/cMN1dEvhWRbLeslYhMcmNaJSJX1rCOqSJyXZVpi0XkLPf3biIy3V3PTyJyXg3rGS4iS3zeTxeRb33ezxGRM93f14nIiSIyErgLON89i/f9fHNE5EsRyXfP2prV+MnV7lsgRUR6utvuCSS50xGnpjwVaOVTU2olIoki8k9xmjZ+cX9PdJc5XkQ2isjtIrIFeEVELheRL6p8Jur+cyIip4rIMnd/NonILdV8hlcA/8E5+dgrIvcBvwHeVtUfVLUEeAA4VkQ6uouNBh5X1Y2qugl4HLj8ED8r31i6uzWGXHGaQ0/3KRsrIs+LyBQR2QcMdz+zd92z3rUicoPP/PeKyAQRGefu/w8iMtAtew1oB3zo7vNtIvKMVK65lkk1LQsi0gW4BrhAVaeraqGqlqvqF6p6uc98vxOR5e6214jIVT5l3r/lbeK0nGwWkTPdv9cK93t/V5V9edv9v8sXkSUi0kVE7nSX3yAiJ/uzbT/+Bt7jj6eG8o44Ced1Vb3NndYM6AC8pKol7utLVf2iyuLedR7qCdO9wJeqepOqbgRQ1W2q+k9VfbOuhetz0pkPzAQO+AeuSkS6AcOAVTWUNwPOBhap6rfATuBkn1kuBcYBqOq9qnqvO/1mYCOQBWTjHEQPdVyhcew/Kz8F52zkl1rmHw00BtrinOVdDRS6ZW+6cbUCzgEeFpER1axjPHCh942I9MA5257sHpCnA2/gnDFdADznzlPV10BncZJwPM5ZdysRSRORZGAgMMd3AVX9GOcs/i23htLHp/gi4HfudhPw429ci9fY/7mOdt97Y9gHjKJybekX4G5gCM5JSh9gMM4ZtFcLnLPQHJzaRl3+C1ylqmlAL2BG1RlU9b84f8O5bhze5g/f5jLv773cnz0B32S92J3mzOw0O13kR3z7N+D8/T4EPsH5/K8HXheRrj6zXQQ8BKQBX7nzL8apaZwA3Cgip/jMfzrOdzIDmAQ84+7zpcDPwGnuPv9dVa/z/i1wzpp3Ax9UE+oIYIOqzq9jl7bhJO90nO/UkyLS36e8Bc6JSGvgL8BLwCXAAJxjxp9FpIPP/KfhfIeaAIuAaTjH0dbA/cALB7Ht2pwEFANbqyk7AifhvKCqf/GZvhPnGPc/N3lm17Du7UCRu40KInKRe6JR06udO+uJOCfph6Q+Jx1wviTXi0hWDeUL3bOx5TgJ6rkq5U+Jc21lMbAZuMmd/irOFw8RaYqTBN6oZv2lQEsgR1VLVXWO1j6Y3UKfP+BTvgWq+hXQ1P3nvgw3ydWiFCfZdHLP8Baoap6ItMWp6t6uqkWq+h3OGXTVZiaAiUBfcarpABcD76lqMc4/yzpVfUVVy1R1Ec4X7dyqK1HVQpzaw7E4/6yLcWqWx+AcvFeq6s469sfXK6q6wl3vBGqoobpaVfnHqFob+x9woXswvcB9X5eLgfvds7ftwH04Jx5eHuCvqlrsxliXUqCHiKSr6m5VXejHMgAfA+eJSG83ef8F56QmxS1PBfb4zL8HSBVxruuoam9Vre57W5sh7nofdc+UZwAf4XNyAnzgnkF7gCOBLFW9351/Dc6B2/e6yheqOkVVy3EO2L4nGNVy/6ffB653v3tVNQO2VFlmo/sdKPJ+p1V1sqquVscsnGQ6zGexUuAhVS3FSYzNgH+par6q/gAsqxLvHFWdpqplwNs4J5yP+izfXtzWFz+2XdO+z8epgd/mbWqtohfQCHjLd6J77BkOrMOp9W4Wkdki0rnKfIXAbcAkEfnOZ/obqppRy+tnd9ZKn70419Vy3ZrpS3XtX71OOqq6FOcf4o4aZumP8w90PnAUzh/K1w3uh9laVS92DzDgHJhOc8/2z8P5om2uZv3/wDmz+MStPtcUR0U8Pn/AG6opfw24DueLM7Ga8qrzTgPeFKcJ6O/ugbUVsEtV833mXY9zJlaJO89k9h8gLmR/k14OcJTvAR3nYNyihnhmAcfjJJ5ZOEn+OPc1q459qcr3YFKA8zesyS9V/jEm+Ba6/yircGpVK1V1gx/bb4XzmXmtd6d5bVfVIj/W43U2cCqwXkRmicjR/iykqp/iXPB9F+dAsg7Ix6nFAuzFOYv2Sgf21nHiU5dWODUI32adqt8f388whyqJH6fG73uWXfXvmSS1XAtzv8fvAG/U0lyzE+eEr4KqtsE5ICbi1gpFZJSIfO02leXi/B18m2t3uskQ9rcU+NYuCqn8/atatqOa5VP93HZNBuH8T97rfhZVTQJeBmb4nDACoE5T63Wq2hHnb7OPKiew7mfvPZHq50c8VVX67FX1GXWuUf8TqC7eSup10nH9Fac3xgEHVXCyv3sgmotzplgnddrH5wJn4fxhXqthvnxVvVlVj8BpQrhJRE44+F2o8BpOO/UUVS2oI8ZSVb1PVXsAv8KpmVyG0yTXVETSfGZvB2yqYVXjcWoCR+M0M3zuTt8AzKpyQE9V1T/UsJ6qSWcWdSedUA1xPg6nKbS62mN1MfyC8w/r1Y7KTZ1Vl9nH/toHIlIpMavqt6p6Bk5z1fs4tTe/qOqzqtpZVbNxkk8c+y8u/0Dls/A+7rTD8QvQVipfm6z6/fHd/w3A2irfkzRVPRX/VPf5P41z8fueasq8ZgBtxL0+VB1xrsO9CzwGZLsHxilUbrIMisPZtnvS8D5OE161PWpV9SacE+4ZIlLTsW8D8Cz7m2O9st11v+97giIiF0vl62lVX97mtc9wjo2HpN4nHXU6CbyF00ujNo8CV1Y9INRiHE4V9EjgvepmEJHfiEgntzljD1BODRf+/KGqa3EO0nfXNa84F++PdC845uE0E3jcL9pXwCPidDboDVxBzc1KU3AOsPfjXF/xxv8R0EVELhWRePc1SES617Cer4CuONc/5rlNEzk4NczZNSyzFac5Itjfw7dwrtFVd7DfCmSKSGOfaeOBe0Qky73e9xdqb5ZbDPQUkb4ikoRzoRUAEUlw/5kbu00wefj5HXH/fr3E0Q54EafpZ7c7yzicE53W4nSlvhkY68+6q2yj4gXMw6mN3Ob+zY/HuY5RU41jHpAvTseKZBGJdWMe5GcIW3GuUXjjuQrnf+DiKrWtSlT1J5zrJ2+KyEnebeOcgHkl4NR6tgNlIjKKytdqg+mwtq1OE7d3PTW5Duck8TMRyRaRJiJyn3tMinG/u/+Hc83Vl7c2Uuw7UVVf18o9Qau+vM1r9wLDxLkNpTVUXBev6dhQSb1POq77ObDprBJVXYJz8LvVz3VOxDloTqyl1tEZ+BSnmWMu8Jyqfl7DvH5Rp/dNbR0IvFrgNEHk4VyzmsX+GtmFQHucs9aJONcfPq1he8U4SfVEfK5buU1vJ+NU83/BaSL5G84/UnXr2QcsBLw9rcD5TNar6rYa9uFt9+dOEfH3OsdBU6dn06dazfUXVf0RJ8mscZuHWgEP4nRU+R5YgrNfD9ay/hU438FPgZVA1d5ClwLrRCQPp7PAxQAi0q7KGWRVSTh/k704B/e5wJ99yl/AuYi/BKf2MxmfC9ni9BS7uKa4cVoHCqu82uIkmVHADpzroJe5n1N1+16OU8vui9OlewfONcTG1c1fjUdwEnyuOL36LsRJQr/4nGHfVcOy1+J0m34C2IXT7PgATnP6z+53+Aack43dOB0gJvkZ12EJ0LaVWo7Rbi1lDM5341Oc5tX27u95ON+JYg7s0Vhrz7g6g3K+70cBbYDFIpKPcw33Fyp/P6slh9f827CJyGqcXkfVHrCNMSZYROQXnM4EAR1BRJzONs+oavNArtdfDaWmE3AicjbOmcYB3VuNMSYE7gD+KiILArVCd10PA7cHap0HHYPVdA4kztAzPYBLVXVamMMxxjQgbnNhdU2Gc1R1VKjjCTVLOsYYY0LGmteMMcaETCAHK4xozZo10/bt24c7DGOMqVcWLFiwQ1VrGvXloEVN0mnfvj3z59c1TJMxxhhfIrK+7rn8Z81rxhhjQsaSjjHGmJCxpGOMMSZkLOkYY4wJGUs6xhhjQsaSjjHGmJCJmi7TxjQEHo+HzRvXs2P7Vso0DiToj4YxDZUqKfHQsXtvEhJqe4JCYFnSiTATJ05kyJAhtGxZ7bObTJRbvfx7ZOO3dFv8OAmFW5GaHzljTK08EsfmjuezaPfptGzTjnadeoRku9a8FkHKy8s566yzGDTI3+dfmWiTV1jKEd/cRWLBZks45rDEaBktV79FbGom7//3SX5etSw02w3JVsxB2bSppidLm6gnMcR4Suqezxg/xGgZIjEkpaYyZ8rbdS8QiG2GZCvGL7GxsXTp0oXMzMxwh2KMiSLJKWns2VnTA34Dy5JOhCkrK2PXrl3YIyeMMSEjgsdTHpJNWdKJIEVFRaxZswZVJTc3N9zhGGPC6PL3Czlx3L5whxFw1nstgpSX7z/TWLduHU2aNAljNMYE3uXvF/Lq4lJu/VUCfz8pqWL6xjwPbZ/cy+ejUzi+fRxyX15FWXIctM+I4ff947np6MSK6etyPXT4196K9+mJ0K1ZDHcNTeSMbvF+x+S7La/EWCi6Jx2A9v/MZ/0e5V8jE7nhqMRK8/3p4yL++U0JJ3SI5dPLGgHgUeWJuSWM/a6Utbke4mMgJyOG07rE8eCIpAO2FW0s6UQQj2d/b6SZM2fSr1+/MEZjTHAkxcFT35Rw7aAEcjJqbmx5ZlQSZ/eIo7AUPlldxnVTi0hNEMYMqHxPyQcXJDO4dSy7C5W/fVnC2RMK+eL/hCFt/D+8ebflVfXup3aNhf8sLK2UdIrKlHHfl5LTuPLc988q5l/flPD0qCSObhNHUZmydJuHrzeGpvkq0lnzWgTxvY4zderUMEZiTPD8qm0sfVrEcNeMolrna5wELVJj6NAkhqsGJtA7O4Zpq8sOmK9pstAiNYbuWbG8dFoSCbHwwY9lzFxXRuz9eWzYU7lr+bjFJTR+NI99Jfv/37zb8r6yUysfGi/oGc+a3R6+2bh/++8sK6VJEhzXvnJye//HMq7ol8AlvRPo2DSGns1jOb9XPE+OrFzLefW7Eno8u5eEB/Jo80Q+98woosxT/bXc6audfdmYV3lf3lpaSspDeeQVO8tt3evh8vcLyfpHPmmP5HHMy/uYvf7AzyycrKYTQXxrOuvWrQtfIKbeOX7sgW3/5/WM55pBCRSUKqe+XnBA+eV947m8bwI7CjycM6HwgPI/DEzg/F7xbNjj4dKJB5bPvLzRIcUqwGMnJXHc2AL+NKScga1ia51fVZmxtpzl2z10yaz9PDkuBuJjodQDx7ePo3PTGF5eVMpfj99fQ3lpYSkX9YqnUYL/ozmkJQoX9IrnpYWlHOXWoF5cUMrv+yfw447KiaBlmjBrfRmb8jy0Tq8+3skrSvm/SUU8ODyRs3vEsWizh6snFyLAA9U0wZ1wRCwtU4XXvy/l9qH79+XVxaWc2S2O9EShsFQZ/moB3bNimHpxChlJwltLSznptQK+u6oR3bNq/5xDxWo6ESQuLo7s7GwAtm0LTfdFY8JhWE4cZ3SL45ZPaq7t/H5SEakP55HwYD4nvlaACNx4VM3DtRSVKffNKiavGE48wkkMYwbE8/J3JXjcVoQfd5Tzxc/lBzTRebflfT0wq/iA9Y8ZkMCbS0vJL1Z+3FHO1xvL+V3fA68dPXlKEoWl0PbJvXR9Zi+j3y/k9e9LK9ViHv2yhLO7x3HnsES6ZDo1oXuPS+SxuSWUlB9Y24kR4ZLe8bz2fWnFtK17PXyyuozRfZx9eeuHUvKKlbfOSWZgq1g6NY3h7mMTOaZtLC8sKD1gneFiNZ0IkpqayoABA5gyZQp79+6tewFjXLXVOlLipdbyZikxtZa3bVx7+aH624mJ9HxuH5N+KqV/ywPPwh8a4XQI2LLXwx2fFnNuj7iKWoavk18rIEagsAyaJAlPnpLIyE7OfKP7xHP3jGKmrSpjVOd4/rOwlAEtY+hXZXvebXk1TT6wFjS4tXMgH7+0lB93eDita9wBzXAA3ZrFsuQPjfhui4cvfi7nq41l/P7DQp78OoY5v2tEcrzww7Zyzu9ZuUZzXPs4isqKWb3LU22tZHSfeP72ZQkLN5fTv2Usry8ppXkj4cQjnHm/3VTOlr1KxqP5lZYrLofk+MgZo8+SToTxXtcpLy9nz549NG7cOMwRGRMcXTJjuWpAPLd/WszUi1MOKM9OFTo1jaFT0xjevyCGLk/vpV/LWI7NqXzYeuWMZAa0iiUjyUmgvjJTYjinh9MsdsIRcYxbXMqDIyr3QPPdVl3GDEjg+fklbNijvH5Wco3ziQj9WsbSr2Us1x+VwBc/lzHslQIm/FDK6L6HNrhm96xYBraKYdxiJ0mPW1zKJb3jiY1xEopHoXtWDBPPP/CzTPG/M1/QBbV5TUReFpFtIrK0hvJuIjJXRIpF5Baf6V1F5DufV56I3OiW3Ssim3zKTg3mPoTS9u3b+fTTTyveb9myJYzRGBN8fz0ukV/yPby4oPahfZqlxHDtoARumFp0wI3TrdOdhFE14XhdNSCeD1eU8cL8UgrLlAt7HfoR+JLe8azc6SEtEU7q6P81ku7NnNi27XNi79k89oAL/LPWlZEcBx1rSX6j+yQwfmkpCzeXs3irh8v67N+Xga1iWbPbQ3oiFcna+2qVFjlXUoIdyVhgZC3lu4AbgMd8J6rqT6raV1X7AgOAAmCizyxPestVdUpgQw6f8vJySktLEXe4+rVr14Y5ImOCK6tRDHcck8g/v657PLnrBjsX7d9cenC9sYa2i6NrZgy3TC/igp7xpCUeelNTeqKw6aY0vr86lZgaHitx9oQCHv+qmLkbylif6+GrDWVcOrGI+Bj4dRenlnbn0ATeXV7Go18Us2JnORN+KOXeWcXcfHQCCbE1x3dhrzh2FypXTCqkf8sYejXfn/gu7h1Ph4wYfv1GAZ+sLmNdrtPb7pE5xbz/Y+Rc0wlq0lHV2TiJpabybar6LVDbJ3ICsFpV1wc6vkjj7b2WkZEBwJQpDSafGlOjPx2dQLOUuhNBdmoMl/WJ5y8zi2vsWlyTK/vHU1LOAR0IDkXjJKk1cY3sGMfHq8s4a0IhXZ7Zy7lvF5IQC7MuT6GHe63m1M7xvHx6Eq8uLqXXc/v407QirhmYUKmXXXUyU2L4dZc4vtvi4bLelWtsSXHCrMtTGNgylt99UEiXp/dy1oRC5v1STk7jyKnpSLDH+BKR9sBHqtqrlnnuBfaq6mPVlL0MLFTVZ3zmvRzIA+YDN6vq7hrWOwYYA9CuXbsB69dHdt7atGkTbdq0oUOHDqxdu5a+ffuyaNGicIdlIsiCBQsY8OGIcIdR79w2vYjpa8pYdFVquEOJOAtOm8GyOZPYs3Mr1z3w7wPKRWSBqg4M1PYiJ/1VQ0QSgNMB3zG3nwc6An2BzcDjNS2vqi+q6kBVHZiVlRXMUAPCW9MREZKSkvjll1/CHJEx9dueIuXbTeW8uKCEPw0J3dMxTc0ivffaKJxazlbvBN/fReQl4KNwBBYMKSkptGzZkuTkZDIzM+1eHWMO0xlvFvDNpnIu6BXPJb0jqAtXFIv0pHMhMN53goi0VNXN7tvfAtX2jKuPMjMz6dmzJwUFBaSnp7Np0yaKiopISrJBAo05FMG4v8gcnmB3mR4PzAW6ishGEblCRK4Wkavd8hYishG4CbjHnSfdLWsEnAS8V2W1fxeRJSLyPTAc+FMw9yHUVBURoXv37gB89dVXYY7IGGMCJ6g1HVW9sI7yLUCbGsr2AQc8QlNVLw1MdJFn9erVzJgxg86dO/PQQw/x8ssvs29fw3uehjEmekV0R4Jo4/F4Km58O/LIIwH46aefwhmSMcYElCWdCOLbe61p06akpKTw7rvvhjkqY4wJHEs6EcRby/GOSBAXF8eyZcvCGZIxxgSUJZ0I4q3pxMQ4f5ZWrVqRn59/wFhTxhhTX1nSiSBNmjShVatWJCc7o9d27twZVWX16tVhjswYczju/LSI7MfykfvyGPtd3ePMNWSRfp9OVGnZsiWdOnWqqOn069ePDz/8kM8++4xOnTqFOTpjDt/l7xfy6uIDh1ocf3YyP+4o575ZJfRrEcPCKsPVLN5STt8XnJ6cG/6USpv0mIrHBaz9YyrtMyL3/PmbjWU8+mUJ75+fzFFtYmmcKJR5lCfmlvDfRaWsz/XQtnEMNx6VwLWDK4+asGJnOddPLWLO+nJS4oVzesTx+MlJB/XU00hjSSeCeDwePB4PsbHOoIDHHnssYD3YTMMyrF0sE86t/CyajCTh0S/KyUoRftzhqXhQmdcLC0rIaSys3xOcpuaScq11dOfDsXKXhxih0kPi7v6siBcXlvLib5Lo0yKWuRvKGfORMzDole6gpHtLlBPGFdA7O5avrmjErkLl/z4oJLeokDfPOfCZOfVF5J4eRKHFixfzxRdfsGPHDgCGDRtGXFycjUhgGpSEWGiRGlPplRTnHPDTE3EeuubzfJ2CUuX1JaVc0W9/LWBdrodhrxQA0OFfe5H78jh+rFMTUlUe+6qYI/6VT8IDeXR8Kp9/fl358dPt/5nPPTOKuGZyIZl/z2fYK86ycl8e//u+cvPXieP2cfn7hRXvP/ixlH4v7CXloTwyHs1j8Et7WbS5vNp9vfz9Qi6dWIRHnXXLfXkAvLq4lJuPTuC33eM5okkMF/eO5/f9Enhozv4431hSyo4C5Y2zkunbIpYRHeJ49tQk3vqhjLW7PQf3oUcQq+lEkKq91xISEujatStLlzaYkX5MENz4cRHfban+oBdMfVvE8s+RgT8hGjMgnl+/UcDjpySREi+8ubSUVmkxDMvZX/Npmy58cEEyZ7xZyLzfN6JtY6moqTz3bSl//ryYf41MYnj7WD5bW86NHxeRliBc0X9/4nrqmxJuOjqBuVekUObnMXzLXg/nvl3IgyMSObdHPEVlyqItHuJqOH3/18gk+rWI4eZPitl40/4mw6IySKpy9E2Oh/V7lPW5HnIyYvhyQzlHt4mlcdL+GtjJHeOIEfhyQxkdmtTPAUwt6UQQ3/t0vBISEio9TdSY+m7munJSH86reN86PYafrtt/QB7aLo426TG8tbSU3/VL4MUFpVzZv/JgnbExQtNk5/8kq5HQInX/Uf/RL4u5fnBCxbNzOmfG8tMODw/NKa6UdAa1juXe4w8uaW7OV0o9cF7P+IrrSN2zan6CaOMkqUgavjGO6hzHU9+UcEKHOHo1j2HepnJeXuRc6/ol30k6m/M9tEit3OQXH+vs9+b8+tuj1ZJOBKla0wGnR1thYSHbt2+nPjyewYReMGobwXRUm1hePXP/NZ3qaglX9o/npYWlDGgVy3dbyvnoomSWbqu7OpJXrGzMU47NqZwIjmsfy7++KaGgVEmJd/6/Brfy/3HTXr2zYzilYyy9ntvLSR3jOD4nlrO6x9P2IB+S9q+RSVz9USF9X9iHAK3ShCv6xfPolyXE1N8+An6xazoRpLqaTv/+/QH45JNPwhKTMYGWHAedmsZUvKrreXZZnwQWbSnnpmlF/LZ7HM1SAn+oqq4HmABVb4sr9cl1sTHC1ItTmDG6EYNaxfLu8jK6PLOXj1Yc3OOgmyYLE85NofDuNNbdmMq6G1MrEtcRTZyfLdNi2LK3cjCl5cquQqVlWv3NTJZ0IkjLli1p1aoVKSn7e6Ycd9xxAMyZMydcYRkTck2ThXN6xPPZ2nKu7F/9tQvvNZxyn6SQnii0SRdmr698jWvWunI6NJGKWk5NmjcSfvFpuiouU5Ztr1zDEhEGt47lrmGJzP5dI47LieWV7w4u6fjuQ5v0GGJEGL+0lGNzYslq5ByWj2kby9yN5eQV749n+poyPArHtK2/jVSWdCJIu3btaNu2LY0a7X8GyIgRzqOJFy9eHK6wjAmLl05LYvutqYzoUP0BNqexECMwZWUZ2/Z52FPkHJzvHJrI0/NKeGlBCSt3lvPC/BKen1/CXUMT69zmiUfE8e8FJczdUMbSbeVc/kEhJeX7D/pfbSjjgVnFfLOxjJ/3ePhsTRnfb/XQo9n+Q2m3Z/byzLzabwD9dlM5b/9QyupdHuZuKOOcCQV8t6Wcp3yaSi86Mp5mKcJF7xayeEs5n68t49opRZzfM44OTervobv+pssGqLS0lLKyskrTUlJSaNy4MXv27AlTVMaER1KcVHSlrk52agyPnJDIo18Wc+O0Ioa1i2Xm5Y34w8B49pUoD39RzDVTlLbpwqMnJlbqRFCTx05O5MoPlVP+V0DjJOGuoYls37c/6TROFOZuLOfZb0vYXaS0SBUuPjKePx+3P6H9tNPDjoLarz8Vlyv3zSpm9W4PCbFwbE4cX/1fI47M3n+dKTVB+PTSFK6fWsTR/91HcrxwTvc4njilfl3Dq0qiZVyvgQMH6vz588MdRq1mz57Ncccdx6BBg5g3b17F9EsuuYTZs2fz888/hzE6EwkWLFjAgA9HhDsM04AsOG0Gy+ZMYs/OrVz3wL8PKBeRBao6MFDbq791tAaout5rAD179mTDhg3k5uaGISpjjAmcYD+u+mUR2SYi1d7dKCLdRGSuiBSLyC1Vyta5j6X+TkTm+0xvKiLTRWSl+7NJMPchlKrrvQb7R52eMGFCyGMyxphAqjXpiMgVInKrz/tNIpInIvkicrUf6x8LjKylfBdwA/BYDeXDVbVvlardHcBnqtoZ+Mx93yDUVNM55phjAKf5zRhj6rO6ajpXAy/7vN+mqulAFnBhXStX1dk4iaWm8m2q+i1wMP0NzwBedX9/FTjzIJaNaDXVdI4++mgAvv/++5DHZIwxgVRX0hFV3enz/m0AVS0CkqtfJGAU+EREFojIGJ/p2aq62f19C5Bd0wpEZIyIzBeR+du3bw9mrAHRsWPHA+7TAYiNjSU1NdU6EhhQxSPW6dQEhkoMaGgHD60r6WT4vlHVhwFEJAZoFqSYvIaqan9gFHCtiBxbdQZ12qNq7H6nqi+q6kBVHVgfhpDp0KED2dnZByQdcO7h2bNnD+XloR/Y0USOBCmjIKNruMMwDURJcgu0pCCk26wr6XwiIg9WM/1+IKjjsqjqJvfnNmAiMNgt2ioiLQHcn9uCGUcoFRQUUFxcXNHM5mv48OEArFmzJtRhmQjSun1nVh31CHub9LQajzksnpgE1vW8gd3bNqEeD7Fxofk+1bWVW4H/iMgqwHtLfB9gPvD7YAUlIo2AGFXNd38/GSfRAUwCRgOPuj8/CFYcofb555+zbNkyMjMzDyi7+OKLefbZZ1m+fDmdO3cOQ3QmEjTNzCQvtyU/9L2XuEZNQOyuB3OIPOUU5W5hx9rv2btnF1kt24Zks7UmHVXdB1woIkcAPd3Jy1R1tT8rF5HxwPFAMxHZCPwViHfX/W8RaYGTwNIBj4jcCPTAabqb6F5QjwPeUNWP3dU+CkwQkSuA9cB5/u1q5PP2XvN2kfbVu3dvRISZM2dy+umnhzo0E0FyjujMlrXLmfXW66SkZRAbF3dA5xNj/KGqFO3bS3xiIiN+e2lItulXfUpV14jIDlXNE5F0f1euqrX2cFPVLUCbaorycGpU1S2zEzjB3xjqk5p6rwE0atSI+Ph4xo8fzxNPPBHq0EwEERGOOuF0GqU1Yd2K7ynclx/ukEw9FRMTQ0bX3vQbehKZ2a1Css2DacSbCfT3+WkCrLakA5Cdnc0vv/wSypBMhBIRjjzqWI486oD+NcZEtENpELZ6fJDUdHOoV48ePSgvL2ft2rWhDMsYYwLGrkJGkN69e9OyZUuSk6u/Ber4448H4J133glhVMYYEziWdCJIx44dadq0KUlJ1Q9dfvbZZwNOLzdjjKmPDiXpRMezEMJgz549FBUVVTSzVdW5c2eaNQv2PbnGGBM8B5N0pMpPE2BTp05l9erV7N27t8Z5hg0bxqpVq0IYlTHGBM7BJJ3zq/w0AVbbfTpe3bp1Y+XKlTYOmzGmXvI76ajqCt+fJvDq6jINkJ7u3CY1fvz4kMRkjDGBdNDXdETk8iDEYdifdGqr6VxwwQUAfPrppyGJyRhjAumgko6IdASeFpEaHydgDl1d9+kAtG/fnoSEBJYsWRKqsIwxJmDqenLon0TkZxHZIyL5wAyc8dN+cJ8gullE/hSSSKPAr371K7Kzs0lMTKx1vjZt2rB9+/ZqR6M2xphIVldN5ypgkKo2VtU0Vc1R1SdUtZn7BNE+OCM9mwDo1KkT6enpJCQk1DrfgAED8Hg8fPPNNyGKzBhjAqOupHOfqm6tqdB91s3fAhtS9Nq6dSuFhYV11mCuuuoqANavXx+KsIwxJmBqTTqqWmcXKX/mMf6ZNGkSGzdupKSkpNb5jj32WJKTk62mY4ypd+q6pjPIfeaN9/1lIvKBiDwlIk2DH1508acjAUB8fDzt27e3MdiMMfVOXc1rLwAlACJyLM4D1MYBe4AXgxta9PHnPh2vlJQUNm7cSF5eXrDDMsaYgKkr6cSq6i739/OBF1X1XVX9M9ApuKFFn4NJOiNGjABgwoQJQY3JGGMCqc6kIyLeB72dgNNl2qvOB8CJyMsisk1EltZQ3k1E5opIsYjc4jO9rYh8LiLLROQHEfmjT9m9IrJJRL5zX6fWFUd94c8wOF4XX3wx4FwHMsaY+qKuxDEemCUiO4BCYA6AiHTCaWKry1jgGZwmuersAm4AzqwyvQy4WVUXikgasEBEpqvqMrf8SVV9zI/t1yunnHIKWVlZdXaZBujTpw9xcXHMnz8/BJEZY0xg1NV77SHgZpzkMVT3j7kfA9xS03I+y8/GSSw1lW9T1W+B0irTN6vqQvf3fGA50Lqu7dV3nTp1Ijk5mdjYWL/mb9++PTt27KjxUQjGGBNp6uq99hdV/VpVJ6rqPp+ircAjwQ2tIob2QD/At3/wdSLyvdt816SWZceIyHwRmb99+/Zgh3rY1q9fX+vzdKq69dZbKS0t5aeffgpyZMYYExh1XTwYKiIP+U5wx12bReXrO0EhIqnAu8CNqurtpvU80BHoC2wGHq9peVV9UVUHqurArKysYId72N577z22bdvm9/A2w4cPB2DWrFnBDMsYYwKmrqRzOtBHRJ4AEJHOwJfAv1X1/mAGJiLxOAnndVV9zztdVbeqarmqeoCXgMHBjCOUDqb3Guxvjnv88RrzrjHGRJS6rukUAb8F2ovIeOBT4FZV/XcwgxLnqPtfYLmqPlGlrKXP298C1faMq4/8ebSBLxEhMzOTNWvW2HUdY0y9UNc1nZuA63Gup5wMLAI6iMhNblmt3EQ1F+gqIhtF5AoRuVpErnbLW4jIRuAm4B53nnTgGOBSYEQ1XaP/LiJLROR7YDjQYEa5PpTEMWTIEMrLy/niiy+CEJExxgRWXV2m03x+f6qaabVS1QvrKN8CtKmm6Aug2jYmVb3U3+3XNwdb0wHnfp133nmHV155hWHDhgUrNGOMCYhak46q3heqQAycd955/OMf//C7yzTA6aefjogwc+bM4AVmjDEBUueAn3WtwJ95jH+OOOIIEhISDqqmExMTQ7du3cjNzbXrOsaYiFfX0e2V2gpFJBanB5kJgOXLlx/UfTpe1113Hbt372bVqlVBiswYYwKjrqTzlTvO2ffuxfyPRGSUiMxyL+SvA74KfpjR4e233yY3N9fvLtNeJ598MgBvvfVWMMIyxpiAqeuazhi3N1kGEIszMsA7OF2VfwIKVDXyb/WvJw72Ph2vTp06ER8fzzPPPMM999wTjNCMMSYg6hwp2h0JwDsawFoRGayqnwQ3rOh0MKNMV9W5c2eWLVvGvn37aNSoUaBDM8aYgDjoo5uq3lF1moj0DEw40e1QazoAp57q3Mb06quvBjQmY4wJpIM/pa7eawFaT1Tzd8y16lxzzTWAXdcxxkS2QCWdgz81NwcYM2YMjRs3PqSaTocOHUhNTWXBggVBiMwYYwIjUEnHbhAJgJycHGJjYw8p6QCce+65FBQUsGPHjgBHZowxgRGopGMCYP78+RQXFx9y0hkzZgyqyiefWD8PY0xkClTSKQnQeqLaW2+9xb59++qesQaDBg0iIyODJ598MoBRGWNM4NTZZdrLfUJnZyDJO819HDWqOiTwoUWfw+m9BhAbG0vz5s2ZP38+e/fuJTU1NZDhGWPMYfOrpiMivwdmA9OA+9yf9wYvrOjkvU/nUJMOwDnnnAPA888/H5CYjDEmkPxtXvsjMAhYr6rDcUYmyA1WUNHqcGs6AH/84x8BeOONNwISkzHGBJK/SafIfYooIpKoqj8CXYMXVnQKRNJp3rw5WVlZLF261EadNsZEHH+TzkYRyQDeB6aLyAfA+mAFFa3uuOMOUlJSDns9J510EmVlZcyZMycAURljTOD4lXRU9beqmquq9wJ/Bv4LnFnXciLysohsE5GlNZR3E5G5IlIsIrdUKRspIj+JyCoRucNnegcR+cad/paIJPizD/VBq1atiImJOayaDsDDDz8MwOzZswMRljHGBIzfXaZFpImI9AbygY1ALz8WGwuMrKV8F3AD8FiVbcUCzwKjgB7AhSLSwy3+G/CkqnYCdgNX+LsPkW7mzJmUlpYedtLJyclhyJAhvPfeewGKzBhjAsPf3msPAN8DTwOPu6/Hal2Iii7Vu2op36aq3wKlVYoGA6tUdY2qlgBvAmeIczQegfN4BYBX8aPGVV+8+eabh3VzqK/+/fuzaNEiPv744wBEZowxgeFvTec8oKOqHqeqw93XiCDG1RrY4PN+ozstE8hV1bIq06slImNEZL6IzN++PfIf+xOIjgRel19+OQCPPVbnuYExxoSMv0lnKc6D3OoVVX1RVQeq6sCsrKxwh1OnQPY2GzRoEGlpaXz55ZcBW6cxxhwuf5POI8AiEZkmIpO8ryDGtQlo6/O+jTttJ5AhInFVpjcIgazpAJx44okUFRUxderUgKzPGGMOl79J51WcC/iPsv+azuPBCgr4Fujs9lRLAC4AJqlTFfgcOMedbzTwQRDjCKlAJ5277roLgH/84x8BWZ8xxhwuf8deK1DVpw525SIyHjgeaCYiG4G/AvEAqvpvEWkBzAfSAY+I3Aj0UNU8EbkOZ7idWOBlVf3BXe3twJsi8iCwCKf7doPw97//nTfeeCNgSWfgwIG0b9+etWvXoqoBW68xxhwqf5POHBF5BJgEFHsnqurC2hZS1QvrKN+C00RWXdkUYEo109fg9G5rcLKysgI+isDtt9/OH/7wBxYtWkT//v0Dum5jjDlY/jav9QOGAA9zEF2mzcGZNGkS5eXlAa2RnHfeecTHx/PnP/85YOs0xphD5VdNxx3k0wTZW2+9hcfjCWjSadq0KS1atGDq1Kn2uANjTNj5e3NohojcICJPiMhT3lewg4s2ge5I4PX73/8eVeWhhx4K6HqNMeZg+du8NgVoDywBFvi8TAAFK+ncdtttxMTE8OqrrwZ0vcYYc7D87UiQpKo3BTUSE7RHESQlJdGvXz8WLFjA8uXL6d69e1C2Y4wxdfG3pvOaiFwpIi1FpKn3FdTIolCwajpARUeCJ598MuDrNsYYf/mbdEqAfwBz2d+0Nj9YQUWrl19+GQhO0jnjjDM46qijmDlzZkVyM8aYUPM36dwMdFLV9qrawX0dEczAolF6ejoQnKQDcN1117Fy5UomT54clPUbY0xd/E06q4CCYAZiYNy4cUFd/znnnEN8fDxXX311ULdjjDE18bcjwT7gOxH5nMojEtwQlKii1Ntvvw0Er6bj7VAwb948vv32WwYNGhSU7RhjTE38rem8DzwEfIV1mQ6aYHYk8Hr00UcBuPXWW4O2DWOMqYm/IxLYDR4hEIqkM3z4cLKyspgzZw579uyhcePGQduWMcZU5e+IBJ1F5B0RWSYia7yvYAcXbUKRdABuueUWPB4PN91kt14ZY0LL3+a1V4DngTJgODAO+F+wgopWoerKfMstt5CTk8PcuXODdkOqMcZUx9+kk6yqnwGiqutV9V7g18ELKzp98IHzPLpg13RiYmJ44IEHWL58OR9//HFQt2WMMb78TTrFIhIDrBSR60Tkt4ANVxxgcXHOJbZQPGzt/PPPJyMjgyuuuCLo2zLGGC9/k84fgRTgBmAAcCnOo6JNAD3zzDNAaJJOQkICRx99NJs3bw76/UHGGOPlV9JR1W9Vda+qblTV36nqWar6dV3LicjLIrJNRJbWUC7uYxJWicj3ItLfnT5cRL7zeRWJyJlu2VgRWetT1tf/3Y1s77//fki398ILLyAi3HbbbSHdrjEmetWadERkqIhc5vP+HRGZ4b5G+LH+scDIWspHAZ3d1xiczgqo6ueq2ldV+wIjcEZD+MRnuVu95ar6nR9x1Auh6r3m1bZtW4YPH87WrVt5/fXXQ7JNY0x0q6umcx+VB/bsCtwK3AvUeXqsqrOBXbXMcgYwTh1fAxki0rLKPOcAU1W1wQ/DE+qkA/DKK68gItx8880h26YxJnrVlXTSVXWZz/uVqrrATSZpAdh+a2CDz/uN7jRfFwDjq0x7yG2Oe1JEEmtauYiMEZH5IjJ/+/btAQg3uLzdl0OZdNq1a8dvf/tbtm7dyvz5NnC4MSa46ko6Gb5vVPUsn7fZAY+mCrfWcyQwzWfynUA3YBDQFLi9puVV9UVVHaiqA7OysoIaayDExDh/jlAmHXBqO02aNOG+++4L6XaNMdGnrqTzo4gccD+OiPwG+CkA298EtPV538ad5nUeMFFVS70TVHWz2xxXjHPT6uAAxBERpk2bVvdMQZCens7111/PRx99xAsvvBCWGIwx0aGupPMn4AkReUVErndfY4En3LLDNQm4zO3FNgTYo6qbfcovpErTmveajzjVgTOBanvG1UfhaF7zuvbaa4mJiakYIscYY4Kh1qSjqquA3sAcoL37mg30VtUVda1cRMbjPG20q4hsFJErRORqEfE+0GUKsAbneT0vAdf4LNsepxY0q8pqXxeRJcASoBnwYF1x1BePPPIIEJ6k07x5cy655BL27t1rXaiNMUEj0TL21sCBAzXSL5T37duXxYsX8/jjj4dlMM7i4mIaN26Mx+Nh165dpKbaoBPGRDsRWaCqAwO1Pn9HJDAhEM7mNYDExETuvvtuSktLueiii8ISgzGmYbOkE0Ei4VrKPffcQ6dOnZg5cyb1oZu5MaZ+8fd5OpeKSFqVab8JTkjRK9w1He+2P/jgAwoLC+3ajjEm4Pyt6TwNzBGR7j7T7g9CPFHNew0lnEkHoEePHlxzzTWMHTuWZ599NqyxGGMaFn+Tzlrg/4B3RORcd1p4j4wN0JQpU4DwJx1wmtliY2O5+eab2bt3b7jDMcY0EP4mHVXVhcBxwBgReQyIDV5Y0SmSehJmZWVx1113UVxczJlnnhnucIwxDYS/SWczgKruAE4BFOgZrKCi1R133AFERk0H4P7776dt27Z89tlnIX/sgjGmYfI36Yz1/qKqHlW9FWe0ABNAM2bMACIn6YAzNI+IcPnll1NaWlr3AsYYUwt/k86d1Uy7I5CBmMjovVZV9+7deeCBB9izZw+PPvpouMMxxtRzcbUVisgo4FSgtYg85VOUDpQFM7BoFI7n6fjj7rvv5ocffuD+++9n4MCBjBo1KtwhGWPqqbpqOr/gPMStCFjg85qEc23HBFAkdSSo6umnnyYuLo4zzzyTrVu3hjscY0w9VdeAn4tV9VWgk6q+6vN6T1V3hyjGqJGd7TyiKNJqOgCZmZk8+OCDlJSUMGzYsHCHY4ypp/y9pjNYRKaLyAoRWSMia0VkTVAji0KTJk0CIjPpANx8880cc8wxrFy5kjFjxoQ7HGNMPeRv0vkvzjN0huI8sXOg+9MEUCR2JKjqs88+IyMjg5deeomJEyeGOxxjTD3jb9LZo6pTVXWbqu70voIaWRT6wx/+EO4Q6pSYmMjMmTOJj4/nxhtvZPdua2U1xvjP36TzuYj8Q0SOFpH+3ldQI4tCX331FRDZNR2APn36MGPGDDZv3swFF1xASUlJuEMyxtQTtXaZ9nGU+9P3QT4KjAhsONGtPjSveQ0dOpRnnnmGq666il/96ldE+gPyjDGRwa+ajqoOr+ZVZ8IRkZdFZJuILK2hXETkKRFZJSLf+9aeRKRcRL5zX5N8pncQkW/cZd4SkQR/9qE+iNT7dGoyZswYevfuzYIFCxg9enS4wzHG1AP+Pk8nW0T+KyJT3fc9ROQKPxYdC4yspXwU0Nl9jQGe9ykrVNW+7ut0n+l/A55U1U7AbsCfOOqF+pZ0AL7++msyMzMZN24cf/7zn8MdjjEmwh3M2GvTgFbu+xXAjXUtpKqzgV21zHIGME4dXwMZItKyppnFORqPAN5xJ70KnFlXHPVFp06dwh3CQUtOTmbJkiU0atSIBx98kOeeey7cIRljIpi/SaeZqk4APACqWgaUB2D7rYENPu83utMAkkRkvoh8LSJnutMygVx3+1Xnr/feffddoH7VdABatmzJ/PnzSUxM5Pbbb2fp0mpbU40xxu+ks09EMnE6DyAiQ4A9QYvKkaOqA4GLgH+KSMeDXYGIjHET1/zt27cHPsIAq08dCarq1q0by5cvJz09nVNOOcUSjzGmWv4mnZtwxlvrKCJfAuOA6wOw/U1AW5/3bdxpqKr35xpgJtAP2InTBBdXdf7qqOqLqjpQVQdmZWUFINzguuiii4D6mXQAOnTowLRp08jLy6Nv375Mnz493CEZYyJMnUlHRGJxnhh6HPAr4Cqgp6p+H4DtTwIuc3uxDcG5CXWziDQRkUR3+82AY4Bl6lQFPgfOcZcfDXwQgDgiwqJFi8IdwmHr1asXY8eOxePxMHLkSKZOnRrukIwxEaTOpKOq5cCFqlqmqj+o6lJV9etpXiIyHpgLdBWRjSJyhYhcLSJXu7NMAdYAq4CXgGvc6d2B+SKyGCfJPKqqy9yy24GbRGQVzjWe//q3q5GvPjev+Tr77LN57733UFV+85vfVFyrMsYYf28O/VJEngHeAvZ5J6rqwtoWUtVany7q1lyurWb6V8CRNSyzBhjsR8z1Tn3sMl2TM888k0mTJnHGGWdw7rnnMnnyZHsOjzHG72s6fYGewP3A4+7rsSDFFLUaSk3H6ze/+Q0zZsygWbNmnH322UyePDncIRljwiyoIxKYg9OnTx+g4SQdgOOOO46lS5fSvXv3ilpPJD+szhgTXP6OSNBYRJ7wdj8WkcdFpHGwg4s2b7zxRrhDCIrmzZszc+ZMWrduzTvvvEPv3r0pKCgId1jGmDDwt3ntZSAfOM995QGvBCuoaNXQmtd8paWlsWrVKgYNGsTSpUtp3bo1P/74Y7jDMsaEmL9Jp6Oq/lVV17iv+4AjghlYNDrttNOAhpl0AOLj45k3bx5jxowhNzeXXr162XUeY6KMv0mnUESGet+IyDFAYXBCik6qyooVK4CGm3S8XnjhBcaNG0dSUhJnnXUWTz/9tF3nMSZK+Jt0rgaeFZF1IrIOeAbnJlETIL4H3YaedAAuvfRS1q9fz8knn8wNN9xATk4Oq1atCndYxpgg87f32mJV7QP0Bnqraj/sAW4B5b1HJ5pkZmYyadIkRo8ezYYNG+jatSsPP/xwuMMyxgSRvzUdAFQ1T1Xz3Lc3BSGeqOWbdKKhpuMlIowdO5bXX3+d+Ph47r77bnr27Mkvv/wS7tCMMUFwUEmniug5MoaAiHD00UdX/B5tLrroIjZv3syAAQNYtmwZ3bt3Z8qUKeEOyxgTYIeTdOzKbwDFx8fz0ksvAdGZdACaNGnC/PnzefHFF2nRogW//vWvOfHEE/n222/DHZoxJkBqTToiki8iedW88tn/FFETINaDy3HllVeyZMkSHn74YT7//HMGDx7MqFGjyM/PD3doxpjDVGvSUdU0VU2v5pWmqv4OFmr8sG/fPk4++WQgems6vhISErjzzjv57LPPyM7O5uOPPyYzM5N77rknKjtdGNNQHE7zmgmg8vJyNm/eDFjS8XX88cezZcsW7rvvPkSEhx56iO7du/PFF1+EOzRjzCGwpBMhorX3mr/+8pe/kJuby6WXXkp+fj7Dhg1jyJAhTJgwIdyhGWMOgiWdCBFtN4ceiuTkZMaNG8eqVat4+OGH+fbbbzn//PNp06YN48ePD3d4xhg/WNKJEHadwn8pKSnceeedrF69mmHDhrFp0yYuuugisrKy+M9//hPu8IwxtbCkEyHi4+MZPnw4YDUdf7Vv357Zs2ezbt06Tj75ZHbu3MmVV17JUUcdxdixYykstOEBjYk0QU06IvKyiGwTkaU1lIuIPCUiq0TkexHp707vKyJzReQHd/r5PsuMFZG1IvKd++obzH0IlfT0dB5//HHAks7BysnJYdq0aWzfvp2nn36a3Nxcfve735GamspJJ53E0qXVfv2MMWEQ7JrOWGBkLeWjgM7uawzwvDu9ALhMVXu6y/9TRDJ8lrtVVfu6r+8CHXS4NOTn6YRCZmYm1113HcuXL+e+++6jSZMmfPrppxx55JG0a9eOxx9/3JoxjQmzoCYdVZ0N7KplljOAcer4GsgQkZaqukJVV7rr+AXYBmQFM9Zw27JlS0Xzmjk8MTEx/OUvf2HHjh1MnjyZfv36sXHjRm655RY6dOjAvffeyzfffBPuMI2JSuG+ptMa2ODzfqM7rYKIDAYSgNU+kx9ym92eFJHEmlYuImO8j9jevn17IOMOuPLycvLynLFUraYTOKeeeioLFy5k586dPP/883Tr1o3777+fIUOGkJGRwWWXXcbKlSvDHaYxUSPcSadWItISeA34nap620XuBLoBg4CmwO01La+qL6rqQFUdmJUV2RUlu08nuJo0acLVV1/NtGnTmD9/PiNGjKCgoIDXXnuNLl260KxZMx5++GF27twZ7lCNadDCnXQ2AW193rdxpyEi6cBk4G636Q0AVd3sNscVA68Ag0MYb9DYfTqh079/fz777DOKior43//+x5AhQ8jLy+Puu+8mOzubfv36cckll7Bw4cJwh2pMgxPupDMJuMztxTYE2KOqm0UkAZiIc73nHd8F3NoP4hyZzwQaRNckq+mEXkxMDBdffDFz586luLiYhQsXcuutt7Jy5Upef/11BgwYQEpKCscccwzPP/88paWl4Q7ZmHov2F2mxwNzga4islFErhCRq0XkaneWKcAaYBXwEnCNO/084Fjg8mq6Rr8uIkuAJUAz4MFg7kOoNGrUiFGjRoU7jKglIvTr149HHnmEvLw83n33XUaOHElCQgJfffUV11xzDc2aNePss8/m+uuvt8ctGHOIJFqG0x84cKDOnz8/3GHU6uuvv+boo49mypQploAiyPr16/noo49YvHgxkydPrniqaXx8PB07duSEE07gyiuvpE+fPmGO1JjAE5EFqjowUOuzxxNEIGteiyw5OTlce+21gNPLcOLEiYwbN4558+bx448/8uOPP/Lss8/SsWNHBg8eTEpKCueccw4nnXQSsbGxYY7emMhiNZ0IsWLFCnr27ElZWRkff/wxp5xySrhDMn4oLi5mwoQJrF27loULF/Lpp5+yb98+wDl5aNasGb169eLKK6/klFNOoWnTpmGO2JiDYzWdBsrj8VBWVhbuMMxBSkxM5NJLL614X1ZWxuTJk3nvvff4+uuvWb9+PZ9//jmff/45ANnZ2SQlJdGrVy+OPfZYzjzzTLp06RKu8I0JOavpRIhly5bRs2dPAKZNm1bxFFFT//3888+sXr2ab775hjfeeINly5ZRXl5eUR4bG8uJJ55I3759ycjIoEePHpx00kkkJyeHMWpjHFbTaaDsPp2Gq127drRr147hw4dzxx13AE5z6vvvv8/s2bNZu3YtW7Zs4YknnqjULTspKYns7Gx69uzJZZddRpcuXejUqRNpaWnh2hVjDpvVdCLEkiVL6N27NwDTp0/nxBNPDHNEJtRKSkr48MMPmT59OgsXLmTt2rXs3r27Uq0InJpRRkYGrVq1okuXLgwdOpTTTjuNnJwc4uLsPNIEVqBrOpZ0IsSmTZu4/vrrmThxIp9++iknnHBCuEMyESI/P5+1a9eyYsUKnn/+eVavXs327dspKCioNF9CQgIxMTFkZGTQokULcnJy6N69O8OGDWPo0KGkp6eHaQ9MfWZJ5xBFetIBmD17Nscdd5wlHeMXj8fDsmXLWLFiBXv27GHp0qWMHz+eXbt2UVxcfMD8GRkZFBUV0bhxY5o3b067du3o2LEjxx57LAMHDqRVq1bEx8eHYU9MJLNrOg2Ux+OhpKQEsGs6xj8xMTH06tWLXr16VUzzPgjQ4/GwatUq5s2bR15eHvv27eOHH35g6tSp5ObmsnXrVpYsWQLAU089BTjfOxEhKSmJ9PR0mjZtSosWLTjhhBMqOjmkp6fTtWtXS07mkFnSiRALFy7kpJNOAizpmMMXExNDly5dau2O/fPPPzN//nxKS0vJz89n8eLFfPLJJ+zatYvdu3ezdetWli1bxowZM6pdf2JiIo0aNWLo0KF07doVVWXXrl20adOG9u3b06lTJ7p27UpmZqZ9p00FSzoRwnqvmVDz9qqrTV5eHlu3bmXHjh189913zJgxg02bNrFz505yc3PZu3cvX331FZMnT65xQNTExESaN29OXFwc+fn5pKamkpGRQdOmTWnevDm//vWvadGiBeXl5SQkJNC+fXvatGljtakGypJOhLBRpk0kSk9PJz09nc6dO3P00Ufzhz/8odr5VJVt27axbNkyVq9ezfr169m4cSObN2+ma9eu7Nmzh0WLFrF161Zyc3NZt25dxbJvvvlmtesUEWJjY+nZsydNmjSpaCZMS0ujcePGNGnShObNm3PqqafSpEkTCgsLSU1NpWXLlrRs2dKSVoSypBMhoqVDh2mYRITs7Gyys7P9fux6fn4+a9asQUTIzc1l0aJFLF26lB07dlTUpPbt20e7du3Izc3l559/ZteuXZVO0ACee+65GrfhHZQ1LS2N7du3U1JSQnJyMikpKTRq1IjmzZtzwgknkJ6ezrp164iLi6uogfXs2ZPMzEyaN29+WJ+NqcySToSwmo6JNmlpaZVG5j722GP9XjY/P58NGzawefNm0tLS2L17N3PnzmXNmjXk5uayZ88e8vPz8Xg8dOnShfz8fH7++Wdyc3MpKyvD4/FUnOhNmjSpxu0cd9xxzJw5E4B+/fqxadMm0tLSSE1NJTU1leHDh/Pgg87TVW666aZKSS05OZk+ffpUjBj/4YcfEh8fX1GWkpJCVlZWRVIrLS2NitqZJZ0I0bp1ay644ALefPNNSzrG1CEtLY0ePXrQo0ePimmHMkjuvn37KCwsJD8/n6VLl7J582Z27tzJY489Rrdu3bj99tsr5j377LPZtGkTe/fuJT8/n/z8/EotFNOnT2fz5s0UFBRQWFgIwOjRoyuSztlnn33Ada9rrrmGZ599ltLSUhISEoiNja2UlK677jpuvvlm8vLyOOOMM0hKSqr0Ovfcczn11FPJzc3lySefPKD89NNPJzMz86A/l2CypBMhcnJyuOKKKyzpGBNCjRo1olGjRjRr1owOHTpUTJ88eTIJCQmVnmt1zz331Loubxd0cJrLi4qKKrVgzJs3j8LCwoqkVFBQwBFHHFEx/4MPPliprKCggLZt2wJU1M527dpFUVERRUVFFBcXM2jQIAB27NjB/ffff0BM33//vSUdU73S0lL27NkT7jCMMUCLFi1YtmzZIS8vIgcM2Nq3b98a509ISODuu++usbxp06bMmjWrxvJOnTpRXl5OSUlJRUIqKiqiZcuWBx17sAX1cdUAIvKyiGwTkaU1lIuIPCUiq0TkexHp71M2WkRWuq/RPtMHiMgSd5mnpAFUDb744gvOOeccwK7pGBNuLVq0YMuWLeEO46DExMSQlJRERkYG2dnZ5OTkkJCQEO6wDhD0pAOMBUbWUj4K6Oy+xgDPA4hIU+CvwFHAYOCvItLEXeZ54Eqf5Wpbf71g9+kYEzlatGjB7t27qx1OyByeoDevqepsEWlfyyxnAOPUOep+LSIZItISOB6Yrqq7AERkOjBSRGYC6ar6tTt9HHAmMDUY8d944428+OKLB3TTTEtLq+h1snr16gOWa9y4Mc2aNUNVWbNmzQHlGRkZZGZm4vF4WLt2rfVeMyaCeJulevfu3SBG7n733Xfp1q1buMMAIuOaTmtgg8/7je602qZvrGb6AURkDE7tqc47r2vSpk0bmjdvfsBTPZs3b07Hjh0BKrpm+mrRogUdOnTA4/GQl5d3wHpbt25Nu3btKCsrIz8/H3DuKRg2bBhHHnnkIcVqjAmMkSNHcskll1BUVBTuUAIiMTEx3CFUiISkEzSq+iLwIjijTB/KOm655RZuueWWgMZljIlsrVu35rXXXgt3GA1SKK7p1GUT0NbnfRt3Wm3T21Qz3RhjTISLhKQzCbjM7cU2BNijqpuBacDJItLE7UBwMjDNLcsTkSFur7XLgA/CFr0xxhi/Bb15TUTG43QKaCYiG3F6pMUDqOq/gSnAqcAqoAD4nVu2S0QeAL51V3W/t1MBcA1Or7hknA4EQelEYIwxJrDsyaHGGGNqFOgnh0ZC85oxxpgoYUnHGGNMyFjSMcYYEzKWdIwxxoRM1HQkEJHtwPoqk5sBO8IQTiSwfY9Otu/R53D3O0dVswIVTNQkneqIyPxA9sqoT2zfbd+jTbTue6TttzWvGWOMCRlLOsYYY0Im2pPOi+EOIIxs36OT7Xv0iaj9juprOsYYY0Ir2ms6xhhjQsiSjjHGmJCJ2qQjIiNF5CcRWSUid4Q7nmASkZdFZJuILPWZ1lREpovISvdnk3DGGAwi0lZEPheRZSLyg4j80Z0eDfueJCLzRGSxu+/3udM7iMg37vf+LRFJCHeswSIisSKySEQ+ct9Hxb6LyDoRWSIi34nIfHdaxHznozLpiEgs8CwwCugBXCgiPcIbVVCNBUZWmXYH8JmqdgY+c983NGXAzaraAxgCXOv+naNh34uBEaraB+gLjHSfV/U34ElV7QTsBq4IX4hB90dguc/7aNr34ara1+f+nIj5zkdl0gEGA6tUdY2qlgBvAmeEOaagUdXZwK4qk88AXnV/fxU4M5QxhYKqblbVhe7v+TgHoNZEx76rqu5138a7LwVGAO+40xvkvgOISBvg18B/3PdClOx7DSLmOx+tSac1sMHn/UZ3WjTJdp/CCrAFyA5nMMEmIu2BfsA3RMm+u81L3wHbgOnAaiBXVcvcWRry9/6fwG2Ax32fSfTsuwKfiMgCERnjTouY73zQnxxqIp+qqog02L7zIpIKvAvcqKp5zkmvoyHvu6qWA31FJAOYCHQLb0ShISK/Abap6gIROT7M4YTDUFXdJCLNgeki8qNvYbi/89Fa09kEtPV538adFk22ikhLAPfntjDHExQiEo+TcF5X1ffcyVGx716qmgt8DhwNZIiI92SzoX7vjwFOF5F1OE3nI4B/ER37jqpucn9uwznZGEwEfeejNel8C3R2e7MkABcAk8IcU6hNAka7v48GPghjLEHhtuP/F1iuqk/4FEXDvme5NRxEJBk4Ceea1ufAOe5sDXLfVfVOVW2jqu1x/rdnqOrFRMG+i0gjEUnz/g6cDCwlgr7zUTsigYicitPuGwu8rKoPhTei4BGR8cDxOEOcbwX+CrwPTADa4Tzy4TxVrdrZoF4TkaHAHGAJ+9v278K5rtPQ9703zgXjWJyTywmqer+IHIFz9t8UWARcoqrF4Ys0uNzmtVtU9TfRsO/uPk5038YBb6jqQyKSSYR856M26RhjjAm9aG1eM8YYEwaWdIwxxoSMJR1jjDEhY0nHGGNMyFjSMcYYEzKWdIwJExG5Q0QuFpF7RURFpJNP2Y3utIG1rcOY+saSjjEhJo4Y4BTgE3fyEpwbGb3OBX4IdWzGBJslHWNCQETau89vGodzh3hbIEFVt7uzvI870rmIdAT2ADvCEasxwWRJx5jQ6Qw8p6o9gQE4zzXxygM2iEgvnBrPW2GIz5igs6RjTOisV9Wv3d9HAlOrlL+Jk3DOZP9QJsY0KJZ0jAmdfT6/DwbmVSn/CLgU+FlV80IWlTEhZM/TMSbERKQn8KP7vJsKqlogIrcDK8ITmTHBZ0nHmNAbBXxcXYGqvhniWIwJKRtl2pgQE5HpwGU+jw82JmpY0jHGGBMy1pHAGGNMyFjSMcYYEzKWdIwxxoSMJR1jjDEhY0nHGGNMyFjSMcYYEzL/D29d44NuHfb/AAAAAElFTkSuQmCC\n",
"text/plain": [
"