{ "cells": [ { "cell_type": "code", "execution_count": 23, "metadata": { "vscode": { "languageId": "python" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 103 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'unable to convert 2.27685316369070 + 1.49035927873741*I to float; use abs() or real_part() as desired'\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "<>:4: DeprecationWarning: invalid escape sequence \\d\n", "<>:5: DeprecationWarning: invalid escape sequence \\l\n", "<>:6: DeprecationWarning: invalid escape sequence \\d\n", "<>:4: DeprecationWarning: invalid escape sequence \\d\n", "<>:5: DeprecationWarning: invalid escape sequence \\l\n", "<>:6: DeprecationWarning: invalid escape sequence \\d\n", ":4: DeprecationWarning: invalid escape sequence \\d\n", " d=plot(delta,(lamb,RealNumber('0.01'),Integer(3)),title=\"$\\delta$ vs $\\lambda$\",hue=RealNumber('0.8'),alpha=RealNumber('0.9'),transparent=True)\n", ":5: DeprecationWarning: invalid escape sequence \\l\n", " tl=text(\"$\\lambda$\",(RealNumber('2.4'),RealNumber('0.1')))\n", ":6: DeprecationWarning: invalid escape sequence \\d\n", " td=text(\"$\\delta$\",(RealNumber('1.01'),Integer(2)))\n" ] } ], "source": [ "var('gamma lamb')\n", "delta=sqrt((gamma+1)/(gamma-1))*arctan(sqrt((gamma-1)/(gamma+1)*(lamb^2-1)/(1-(gamma-1)/(gamma+1)*lamb^2))) -arctan(sqrt((lamb^2-1)/(1-(gamma-1)/(gamma+1)*lamb^2)))\n", "delta=delta.subs(gamma=1.4)\n", "d=plot(delta,(lamb,0.01,3),title=\"$\\delta$ vs $\\lambda$\",hue=0.8,alpha=0.9,transparent=True)\n", "tl=text(\"$\\lambda$\",(2.4,0.1))\n", "td=text(\"$\\delta$\",(1.01,2))\n", "(d+td+tl).save(\"aerodyn_3_delta.png\")" ] }, { "cell_type": "code", "execution_count": 64, "metadata": { "vscode": { "languageId": "python" } }, "outputs": [ { "data": { "text/html": [ "\n", "\n" ], "text/plain": [ "Graphics3d Object" ] }, "execution_count": 64, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from sage.plot.colors import rainbow\n", "clr=rainbow(8);\n", "\n", "# Shock Wave\n", "delta,Ma1,gamma,beta=SR.var('delta Ma1 gamma beta')\n", "eq1=tan(delta)==(Ma1^2*sin(beta)^2-1)/(tan(beta)*(Ma1^2*((gamma+1)/2-sin(beta)^2)+1))\n", "eq1=eq1.subs(gamma=1.4)\n", "sol1=solve(eq1,Ma1)[1]#Critical: Solve for Ma1 instead of beta\n", "sol1m=solve(sol1.subs(beta==1.2),Ma1,solution_dict=True)\n", "sol1m=sol1m[0][Ma1]# Find Ma1 for specific delta\n", "#delta==0\n", "fg1=implicit_plot(sol1.subs(delta=0),(Ma1,1,4),(beta,0.2,pi/2),color=clr[0],\n", " axes_labels=(r'$Ma_1$',r'$\\beta$'),title=r'$\\beta\\ vs\\ Ma_1$')\n", "fg1+=text(r'$\\delta=0^\\circ$',(1.2,0.8),color=clr[0])\n", "# plot3d(sol1,(Ma1,1,4),(beta,0,pi/2))\n", "# contour_plot(sol1,(Ma1,1,4),(beta,0,pi/2))\n", " # sol1=solve([eq1.subs(delta==3*pi/36)],Ma1)\n", " # # to_poly_solve=True)\n", " # implicit_plot(sol1[1],(Ma1,1,4),(beta,0,pi/2))\n", "# implicit_plot(sol1[1],(Ma1,1,4),(beta,pi/15,pi/2))\n", "\n", "# f11+=text3d(r'$\\delta$',(0.2,0.2,0.2))\n", "# dmdb=derivative(sol1[1],beta);\n", "# solve(dmdb==0,beta)\n", "for i in range(1,8):\n", " fg1+=implicit_plot(sol1.subs(delta=i*pi/36),(Ma1,1,4),(beta,0.2,pi/2),color=clr[i])\n", " delt='$'+str(i*5)+'^\\\\circ$'\n", " fg1+=text(r'%s'%delt,(sol1m.subs(delta=i*pi/36)+0.1,1.1),color=clr[i])\n", "end\n", "fg1.save_image('1_2d.png')\n", "fg1\n", "fg1_3d=Graphics()\n", "fg1_3d=implicit_plot3d(eq1,(Ma1,1,4),(beta,pi/15,pi/2),(delta,0,pi/3.5),\n", " color='#c0ebd7',axes_labels=['Ma1','beta','delta'])\n", "fg1_3d\n", "# fg1==implicit_plot([sol1[1].subs(delta=k*pi/18)) for k in range(0,7)],(Ma1,1,3),(beta,pi/15,pi/2))" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "vscode": { "languageId": "python" } }, "outputs": [ { "data": { "text/plain": [ "[delta == -arctan(35/251)]" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "#find text pos (deprecated)\n", "delta,Ma1,gamma,beta=SR.var('delta Ma1 gamma beta')\n", "eq1=tan(delta)==(Ma1^2*sin(beta)^2-1)/(tan(beta)*(Ma1^2*((gamma+1)/2-sin(beta)^2)+1))\n", "eq1=eq1.subs(gamma=1.4)\n", "# sol1=solve(eq1,Ma1)\n", "# sol1=sol1[1]\n", "\n", "sol1=eq1\n", "sol1d=solve(sol1.subs(Ma1==1.2,beta==pi/4),delta)\n", "sol1d\n", "# find_root(sol1d.subs(beta=pi/4),0,1)\n", "# sol1m[0][Ma1].subs(delta=0.2)\n", "# dmdb=derivative(sol1,beta)\n", "# sol1d=solve(dmdb.subs(delta=0.4)==0,beta)\n", "# find_root(sol1,1,pi/2-0.05)" ] }, { "cell_type": "code", "execution_count": 63, "metadata": { "vscode": { "languageId": "python" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 16 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.20379237810420833)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 15 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.19771342522750943)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 32 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 31 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 50 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 50 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 69 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 68 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 88 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 89 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 110 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n", "verbose 0 (3839: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 110 points.\n", "verbose 0 (3839: plot.py, generate_plot_points) Last error message: 'Unable to compute f(1.5207963267948965)'\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlIAAAGXCAYAAAB4Er7pAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAA9hAAAPYQGoP6dpAACGpklEQVR4nOzddZxVZf4H8PehuztEwMRCsbu7u3NNdFfdVVd/1rrmWmu7dndjJ2Bii6igiHRJd8zM+f1xBhhgkHBm7gzzffu6r3vnzjnnPhcnPvPE90nSNBVCCCGEEJZflVw3IIQQQgihooogFUIIIYSwgiJIhRBCCCGsoAhSIYQQQggrKIJUCCGEEMIKiiAVQgghhLCCIkiFEEIIIaygCFIhhBBCCCsoglQIIYQQwgqKIBVCCCGEsIIiSIUQQgghrKBquW7AkiRJ0hpnoSo2xhu4MY3NAUMIIYRQTpTLIJUkSYIrcVaapjOSJGmD72XtvTanjQshhBBCKFReh/ZWw2ZYFdI0HYnX0D2HbQohhBBCWEh5DVLT0QJrFHluHJrlpjkhhBBCCItLcjXlKEmSf6Ir7seasrlQm+DkNE1nFnN8L8xO03TXsmxnCCGEEMKS5CRIJUmyg6zXac/C2w5pmk5PkuQGTErT9MpFjt8EvdEtTdMfy7zBIYQQQgjFyNXQ3mx8IZsHdU+aptMLn6+F9YoemCRJHdyFwyJEhRBCCKE8yUmQStP0k8LX3go9i3xqcwyZ90Hh6r27cVGapq+UZRtDCCGEEJYml5PNN8KUNE0HQmGJg/XxWJFjLsEjaZq+XXjMyWXeyhBCCCGEJchlHantMbHIx5fjpjRN+0KSJMehCaolSbI7EmxYxm0MIYQQQliiXK7aexWD8YssME1J0/TGws+tISvAWWOR0x5L0/SYsmxnCCGEEMKS5GrVXhVMwKZpmv5c5g0IIYQQQigBuZoj1RUzI0SFEEIIoSIr8yCVJMnWuFM29+mysn79EEIIIYSSkrM5UiGEEEIIFV153WsvhBBCCKHciyAVQgghhLCCIkiFEEIIIaygnBTkTJKkRCZmpWmalMR1QgghhBBWREw2DyGEEEJYQeV6aC/JHJ0kyZ25bksIIYQQwqJyudfeH0qS5AhsjG3xQ46bE0IIIYSwmHI/tJckyUOQpunxuW1JCCGEEMLCyvXQXgghhBBCeZazob0kSf4p23PvfqyJqtgEJ6dpOjNX7QohhBBCWFa5Kn+wA95HLVyFHdI0nZ4kyQ34O67MRbtCCEuWJEkL7IoCdMMFaZrm5bZVIYSQW7ka2puNL7AZ7knTdHrh87WwXo7aFEL4YzuiUZqmT6A6ds5xe0IIIedy0iOVpuknSZJUxVY4q8inNpf1VIUQcihJkm2wN77ELDRM0/SxIoe0Qv9ctC2EEMqTXE423whT0jQdCEmStMH6eOwPzwohlIVhaIHX8SkOn/eJJEn2wrNpmg7OTdNCCKH8yGWQ2h4Ti3x8OW5K07RvTloTQpivMCTVLBx23x1vQJIkm2BMmqbPJknSJYdNDCGEciGXQWo79EqS5G9JkvwLA9I0/ee8TyZJsn+SJA9gP+yfJMkDSZLsn6O2hlCpFIakGkmS7ItWaZrekSTJFrJVttcmSdIbZbLXZZIkDZIkuSRJkvFJkjyzhGNOTZIkL0mS65Ik2bws2hVCCOSoIGeSJFUwAZumafpzmTcghPCHkiQ5E33TNO2d67ZAkiTN8B90S9N0g0U+tzZORdc0TbfPQfNCCJVYrnqkumJmhKgQyp8kSdriKNTJdVuK2BH3YLWiTyZJUg0boDHey0G7QgiVXJkHqSRJtsadqJYkyWVl/fohhD+WpumINE23SNP0zRU5P0mSqkmSdE+S5NHCOVWSJHk8SZJjixyzR5IkpxTezlyGy3ZM0/Qz5BUuTJnnQLyEnfBu4bWrJ0lyZZIkpxVOHXg7SZImK/JeQghhaco8SKVp+lGappunado8TdN/lfXrhxBK3X54WtajtWrhcz3QAJIkqY9r0jS9B59gr2W45rw5CL9i9cLrrIuB6ID6stp0cC9GpWl6t2ySfFcLL2wJIYQSE3vthRAWkyRJujy3RU5/X1b9fBu8Wvjcr/i88HEeGidJ8rWsrMLBS2lLB/xW+OFArJ4kSXWsl6bp17LCoL3SNM1LkqSrrJfqvsLj18f7aeFk0CRJ6iZJsleSJO+s0D9MCCEsIoJUCGExaZomy3Nb5NxJsl6mnkX2zVxPVtxT4XNdcA32wdFLac4O+KDw8S+yeVIH44XC5+YP68nmUn2UpunsIh+/lyRJo8LXnp6m6WuyyuwhhPCnRZAKIZSGlhgKSZI0xrQ0TQuSJOmQJMnvmJum6bOyArwjlnKtFmmajit8PBB7oH+aprMLVwBvb8FE80kYU+R190ZvHFFSbyyEEIqKIBVCWCbLOSz2FFZJkuQY7I/nCp8fjxtxWJIkZ2BymqavFneBJEnWT5LkXpydJMl5hU/3xytpmn6TJMmBeACNCq/XtvB1JUlyBPbFk7I5Wz8t9xsOIYRlkJM6UiGEiitJkp4VvV7TyvAeQgjlQ042LQ4hlG/LsGlxhZQkSU0chFULe61eKDKfKoQQllsM7YUQirPETYsrsjRNZ6dp+kSapqumafpkhKgQwp8VQSqEsJglbVocQghhYRGkQgiLWcKmxTWTJDlS4bBY4TBZCCFUajHZPISwmPK2aXEIIZRX0SMVQlhIOd20OIQQyqVy1yOVJEki2zdralreGhdCCCGEUER57JGqj8mF98srjVvc4ha3lfn2ipPTI7VLE9unl/gwrUF6QeHnjrlK2nAv6ZDRC58zY7z0htbSx/aUpmkJtOP9iamkZ+qN8Tn/9yi127sTUu0+STX5KPXM2Ny3p2LcKqXyGKRCCCEUo59n9PSI13Sxn509bWvr4994+n0ee4c7/kaHVguf9/pZ5M1i3/tIkuKuvJyuHsKG9di9SQlcrJyZlc85A9mlL2vX5ftNOLRFrlsVyrEoyBlCCBXARIO97BTv2E59DdR1vpESr2HsWE67iUN34KhdFj7vx+fp9yQHPkaDNiXQkD5TeG8Sz3YpoVRWjnw3jaN+ZOBMbu7MX9tRZSV7j6HELXePVJIk2yZJ0iNJkpFJkqRJkuy/yOeTJEkuL/z8zCRJeiZJsk6JtTiEECqZfHlecJTPrOIHc53gek+o4XZ0LuD4a6lbi7vOWTjbTBvLq6ex1gGsd2QJNeaaIaxZmwOal9AFy4GClBuGsulXVE34shtnt48QFZbJigzt1cV3OHMJnz8f5xZ+fhOMxjtJkqzInKcQQqj0erlCHz94UzMn+4vbdHQYjsMtz/Pe1zz0T5o0WHBOmvLa6dnjve8uoc6jftN4eTz/XCULHCuDobPY6TvOH8Tf2vF5N9atl+tWhQpkuYf20jR9Q2GV42SR78zCFXdn46o0TV8ofO44jMGR+N+i1yss6le0sF8ErhBCKDRYL++61mt2soFVfe8IjXA3fhjEhfdw9sHsvPHC5333KD+9wKHPUa+kpvjcMIz2NTmqZQldMMeeGMMZP9OgGu9vwPaNc92iUAGV9GTzjmiFt+c9UbiXVS9suYRzLpSt0pt3G17CbQohhApphglecLSetjQFW7pOH1U8gVqzOeoqVmvLNScvfN6kobxxFusfQ5eDSqgxY+bw5FjOakv1Cr5O6fc5HPYDR/3E3k3pu3GEqLDCSnqy+by1ImMWeX4MOizhnGtwU5GP64swFUKo5FKpV/xFHzV8iotd4WoNXCb7q/TsexgwlM/vplaRPv2CAl46npoN2ePWEmzQ/0ZSLeEvrUvwojnw/O+c/nM2L+qpLhwWK/LCn1Naq/YWrSeRFPNcdmDWYzV/B/ZFhwtDCKEy+tL/fOxNr9rKYfb0sE1shf/DG32yuVG3nMX6nRc+r8+tDP6AY9+ldqMSasycAu4aybGtaFy9hC5axsbNofsvPPM7BzbjzjVoWSPXrQorgZIOUqML71thVJHnW1i8lyqEEEIxxvrBa871pm211tRsZ5mKxzBuQrZKb4/NOOvARc77kXf/yWZ/o9NOJdigZ8Yyek42rFcRzeuFyk95cu2sFyr+aA8lpKQHun+Than5lUySJKmB7fBJCb9WCCGsdOaa6TmH+8i6hshzhBu9pJp70T7lhOuy4x68YOEskD+XF4+lcUd2vqYEG5Sm3DKCXRrTpW4JXrgMjJvDET9y8A9s3ZAfN+XwlhGiQola7h6pJEnqYbUiT3VMkqQrJqRpOjRJkv/ioiRJfsEvuAgz8MSfb24IIazc3naePn73nrVc4O9u0spfcDBufSEb1nv9OlouUlS8178Z8x0nfUr12iXYoK+m8uVUXlm3BC9aBl4o7IXKS3l8bY6IXqhQOlZkaG9jfFDk43kTxR/G8fgPauNONEYf7Jqm6dQVb2YIIaz8+nvZB+7Rww52tZG37a49/ou+v3Le3fztoGxYr6jhffjwara7lLYbF3PhP+Oh0bSuwR4VZDuYEbM56xdeHMd+Tbl7DVrVXPp5IaygJE3L1z6DSZI0kJVBaJim6ZTlPL18vZkQQlhGU4xwh/U9aUNj1be/Jzygtj5YazYbn0q1qvS5c+FVenNm8L8NqdWIEz+maknOfJ1dQOtPOLk113Ve+vG5lJ9y90guHETdqty6Ggc3j16oslUp/7Fjr70QQsixAvlecLRPtPa9PFe7zkVquxEbovtdDBrJl/9bOETBO+czeRhHvFLCIQp6jGNiHse3WvqxudR3GqcMoM9UTm3NtZ1oVEFXF4YKJ4JUCCHkWG9X+dhX3rSx7o5zizXsJtsm4pWPufMl7jibdToufN7At/niDva4jWZrlkLDHhrNZvVZu5xOMp+ZzxVDsorrq9fmw65s3WjhY/pMyYYm29bMVh4Om8XmDRe/1pNjmJafXSeKc4blUMHL04YQQsX2m57ecqXXbamrNQ10jFQ26XT0OE68jn225PT9Fj5v5kRePpFOu7DJGaXQsCl5vD2RI8rpdjDvTmC9L7hpGJd24L+rMWIOD4zi0B/4unBa7l0j6PAZ1Xqxfz+aFNNT9co4tm/EyW2YkJfNswphGUWPVAgh5Mg0Yz3vSL1saaoqdneNK1XxOpoXsNs1VK/G/ecvPtXnte7Mnc5+D1ClNP4kfnMCc1P2b1YKF/8Ths/ivEE8NZbtGvLa+qxZh6YfcdNqnNg6G4489AcGbs6qtRi+RVbGoV2tJV+3SuE/cM0qlXSmT1hREaRCCCEHChR40bE+V9en+LdL/VsjZ2MP3Pgs737F2zfQvNHC5/Z7mn5PcuDjNGxXSg18aRxd69HhD8JHWZpTwM3D+ffgbDL5g2tyXKsFCbNnVzoVqfswLX/B47ZLWbW3bzOeHsusAjrUpE2s8gvLLoJUCCHkwMf+4zMfetWWjrGXh21uXVyLr3/mwnv5+6Hsskg5gykjefV01jmU9Y4opcYVpLw9gdPLSSXztydkJQ1+nclZ7bh8VRou8utrvXoLHvcYt2CV4fQCHhyV9TS9O5Hz2hc/5yv23AsrKIJUCCGUsaE+9qZLvW4nHbU0y1nG4A3kzeTIf7NuR676y8LnpWk2L6paLfa6sxRX9v80g/F52byhXPp5BhcMynrHtmvI8+uwbr0lH//5FF4dzyYNOKx59tzBzdmsQfa4RXUO6JdVOK8S43ehZMRk8xBCKEMzjPecw/W2ubFSh7vBs6q6R7ZlxDl3MHQsT1xMzUX21P3ybn59i/3up07TUmzkR5Opis0blOKL/IFxc/jrL6zzRTZp/Mm1+aDrH4co2LQBV3Skcy22+ZZpeXQrcs5qtRkwMyuXEEIJKTdBKkmS7kmS/IjPc92WEEIoDanUS473uRo+UsX/udg1mvkLDsezPbn3Vf57Jmt1WPjc8b/w9j/Y+DRW36OUG9pvejaBu27VUn6hRczK5/qhrNaHh0dzVUcGLMP+eJ9NpuXH/DYz+3j7Rtm2Nv8dTqOPsuvC1ML7GuXmV99KY3CvXLcgd8rN0F6apnfgjiKVzUMIYaXyqZt94j09bOUou3rC1jrhFgwexcnXc8j2nLz3wufl52UbEtdrzS7Xl0FD+89grTpl8EKFCtJsFd5Fg7ISBqe1yUoaNK+x9HOhWpLVf5o3qXzQLKonbNeIAtQqDIQfT2bLBhVv8+VybO5M3ruIz/7L5ZV0b5FyE6RCCGFlNtzn3nCh1+yogxaqO8cgfInqeRzxbxrV456/L9758vF1jPicEz+i5lJGt0qmsbPLZm+9gjTbXPhfQ7JesP2b8VanrDdseWzcgHPbc+twqibZ0OSr67FNo2zT4huGki+bc/ViBdt8uRwb3ocXj2PSYHa9MdetyZ0IUiGEUMpmmuQ5h+llE2MUuNx/XKCa+9EFFz3IF/358DYa1V/43FHf0PNytv4n7bcoowZPzaNBKQ7rFaTZBPJ/DabvdHZpzD1rsEUxFceX1YHNFzw+p/2Cxzs0zm6hxOTNode/+OhaWm/Ead/QfO1ctyp3IkiFEEIpSqVecZLP8KFqrvJP/9bSkTgB737JtU9w9V/YYp2Fz507kxeOocW6bHdZGTa6WpL15JS0/JQXf+eqoXw7jZ0a8dGGbPUnAlQoU6O/y4aZf/+R7f+VBfwS3+Oxgqnkbz+EEErX5+7wkTf0sJXD7eh522uFuzF2AsdczU4bcX4xNaHevZAJAznlS6ot43ShEtGiRja8V1Km52e1nG4ens1f2qERvbqybaOSe41QqvLz+Pg/We9os7U4+Qtad811q8qHCFIhhFBKRvraa/7hNdtrp7nG/qGfxKeoW8Ah12SjXI9etPg2LwPfps8t7P5fWpb1tJ5N6/PGhKxw1Z8pVjVgBveNyva/m5zHIS14uks2pylUGL/356XjGPklW13A9pdRLYq/zxdBKoQQSsEsUzzrUL1sZKR81/qPc1R3KzbCDc/w1he8+R9aLVITasZ4Xjo+25B407Ny0PiDm3PHSF4Zz37LudfetLxs/tO9o+g9mabVOL4VZ7Vl1dpLPz+UGwUF9LmV9y6k4Sqc+DHtN891q8qfCFIhhFDCUqkeTvaZOXqr4Wrnu0Jr++NM2cTyC+/lvMPZbdNFzk3pcQr5s9n/oVLakHhptmvE7k04qT/N12PLpcxhGjU7237lhXHZZsezCrLhuyfX5oDm2fYsoUKZ+BsvncCQXmz2V3a6hhplWBGjIknStHwVfihSR6phmqZTlvP08vVmQgiVUh+3e9x5HrC1PW1rhIsNl/gWVaez4ck0bcBHt1Gj+sLnfvMQL5/AIc+yzsE5aPw8E+ayz/d8MiULVbs1zjYFrpZkhS1/m5nVm/pkCr8UFsLcvAEHNeOg5nSM3qeKKE356l7e/ju1m7L/g3TcYZlPr5T77kSQCiGEEjRcH/fY1pN2RBP7etAtavgQm6fZPnqv9+Gbe+nUZuFzJwzi7g3ocnD2Cyzn8gp4eAxPjOHDycwt8iO2YdWs3tPG9bMerG0b0iomzlRkk4dlvaED32Sjv2S1oWot33S2ShmkYmgvhBBKyAzjPetQH9vEUHlu9B9nqeFabIEH3uCp93nyksVDVH4eLx5DnebsfksuWl+MalU4qXV2yytg3NysUnidKjSsVoq7JoeylKZ8fX/WC1WjHke+yhp75bpVFUcEqRBCKAEFCrzgGF+r6n3VXOkfLtfWnjgPPw3hrFs5aU8O32nx8z+6luGfcULv5e4FKBvVqkSP00po0lBe+QuD3qHrCex2E7Ub5bpVFUsEqRBCKAEfucbnenrJNg60pdftqg4ewZzZHH4Fq7TglmJW4Q3/PKvPs81FrLJVGTc8VEppylf3ZBth12rEUa+XwWbYK6kIUiGE8CcN8r63Xe51u2qtsbYu9IpEbzRF97sYMJTP76buInOw50znhaNpvSHbXZqL1ofKZuLgrBfqt/cK50LdQK0oLr/Cyk2QSpKkO7oj1smGECqMKUZ63hE+tLWR8lzren9Tw02yeVHP9uTOl7jzHNbvvPj5b/2dqSOyeSlVqy/++RBKSkEBX97NO+dTuwlHv8Vqu+a6VRVfrNoLIYQVlC/Pw3b0lt89o6Xr/ctVtrM9XsCvw9noFPbYjKcuXXxu9oAePLkve9/Nxqfm4A2ESmPCIF45icE96XYqu/ynVObiVcrVB+WmRyqEECqa9/2fL3znZZs4wZ6esZ3GeACzZ3Pov2jRiHv/sXiImjaGl09ijb3pdkoOGh8qhYICvriDd/+ZrQg99l06FbPYIay4CFIhhLAC+nvZe270qt2srbVa/uo7fIzGOPMufhjMZ3fSoO7C56ZpFqKShH3vjyoCoXSMH5j1Qg3pzSZnsPO11Kyf61atfCJIhRDCcppgkBcc5wM7mIrzXeevqrkdG+O5ntzxEneczYarL37+l//jl9c4ogf1WpRly0NlUJBPn9t47yLqteK495erOnlYThGkQghhOcw1yzMO9rkOPjfXXa5xvqYOxRn4dQQnXc8h23P6foufP24Ab53Lxqex5t5l3Piw0hvTL1uRN+JzNj2Tna+hRt2lnxdWXEw2DyGE5dDDaV72vIes5yxH6Olk0/Elas5hyzOZPC2r0dOw3sLn5s/lvi2YM5VTv45fcKHk5M3mw2v48GqadM6GjFfZssybUSkHqaNHKoQQltF3HtXbA16xky2tbpqT9MdnaICz7qLfb3x6x+Ihiqzo5pjvOOnTCFGh5Az7LJsLNf5ntr6Qbf+PalGEvsxEkAohhGUw1g9ecZq37aiKWg51lbNUcS82kM2Luv3FbF7URmssfv6QD7Megx3/TduNy7btYeU0exrvX0yfW2mzMad8Rav1c92qyieG9kIIYSlmmexem3pdQ2+r5z53OMvaDsTDGDQiqxe12yY8fdniq/BmTeauDWjQjhN6UaVqLt5FWJkMfJsepzB9LDteyeZ/KxdfVzG0F0IIYWEFCrzkeF+Z4U2tXeZU11nbqrgLc+Zk9aKaNyy+XhS88VdmTuD4D8rFL7tQgc0Yny1W+O4ROu7Ice9lc6JC7kSQCiGEP/Cx63zmDS/Z3p66+sGhRuIL1LX0eVH9nsl+6e3/MI07lnHjw0ojTfnhWd44K5tYvu/9bHhC1CArDyJIhRDCEvzqHW+7xJt211g9W7rcxRLPYU0L5kXd/rfi50VNHs6rp9HlEDY4powbH1YaU0bw2hkMeIW1D2TP26nfOtetCvNEkAohhGJMMsRzjvCx7fxmttvd4XS1/AMHYdDIrF7Uwdtxxv6Ln19QwEvHU712tpde9ByE5VVQwNf38c55VK/Doc/T5cBctyosKoJUCCEsYq5ZnnaQb7XygXw3+KfLdLAVrsHsORx6Oc0act95xYekPrfw23sc8w51mpTxGwgV3u8/Zb2ZQ3qz4Unsej21G+e6VaE4EaRCCGERrzvTtwZ7UVcn2l0Pu0jxlOyH5rl38/1vfHJ78fOixnyfbRK7+Tl03rmMGx8qtLmzsqKaH11Lo1U59j067ZjrVoU/Um6CVJIk3dEdVXLdlhBC5fWV+3zsYS/Z1fraa+hsH+EDtMbzvbjthWxeVLc1Fz9/7iyeP4qma7DT1WXc+FChDXo/64WaNJit/8k2F1G9Vq5bFZam3ASpNE3vwB1F6kiFEEKZGuELr+rubbvIU81fXOdU1dyIbfDLcE78z5LnRZFtFDt+ACd/Eb8Ew7KZPo63/56t7lxlG454meZr57pVYVmVmyAVQgi5NN3vnnaQPjbxvVkecL2zNHQQzsGMWRx8GS0bc//5xc+L+vUdPruZXW+MCtNh6dKUbx/m7X+QFrDvfXQ9gSoxLlOhRJAKIVR6+fI85wjfqe4N1fzbaW6wjjZ4QFau+cxbsh6pPnfRoJh98maMz1bpddqZzc8u0+aHCmjcgGwYb3BP1juK3W6iXotctyqsiAhSIYRK7wOX+NJnXrClA21ugEMMwueyzYgfeJ0H3+Chf7Jep8XPT9Nsu468Wez/UPQohCXLm81H1/HhVTRozzFv03mXXLcq/BkRpEIIldpPXvSe671qd+01tbWLnSvxONbBt7/Q/b+cvDfH7V78Nb55kJ9e4NDnaNC2DBsfKpTBvehxKhN/Zavz2fbirM5YqNhi0+IQQqX1u5/cYzM9bO5H1dznfkdr7VTchklT2fjUbCjvk9upVXPxa4wfyN1dWedQ9n+gjN9AqBBmjOft8/j2Qdpvyd7/o+W6uW5VqaiUZWejRyqEUCnNNMmT9vOVNfUxx4Mu9w+tdcONsuG6E65j3GTevqH4EJU/lxeOpl4r9rilrN9BKO/SlL6PZZsM58/NKtxvdHIM/a5sIkiFECqdAvmed6TvzfCyts53pEdtaRaeRQ3c8DQvfcTLV9GpTfHX6X0lI7/kxI+oWb8M30Ao937vz+vd+e191j2c3W6mfqtctyqUhghSIYRK5wOX+toHXrCDHa1trpP0wntohw/78s97OP8I9t2q+GsM/SQLUttdSvvNy7DxoVybM53eV/HJDTRchaPeYPUlzK0LK4eYIxVCqFR+8KwnHeZ5+5ihlos94BR13YyzMWYCG57M6u1470aqFfPn5qwp2byoeq04oTdV40/SSi9NGfAKb/yVaWPY5kK2uqDSFWWNOVIlIUmSargcR6EVRuEhXJmmaUFJv14IISyr0fp60fE+spshZnnQLU5U15H4G/LyOOLfFBTw1KXFhyiyX5YzfufYdyNEBSb+ln1N/Pwqq+2e7Y/XdLVctyqUldL4EXABTsNx+AEb40FZL1NMxwwh5MQMEzxlf32tq5eZ/udyF1nV6rhX9qf0ZQ/R67usJ6p10+Kv88OzfPcw+z1Ik2JqSoXKI282H1+f1YSq04xDn2ftA4qveh9WXqURpLbAy2mavlb48eAkSY6QBaoQQihzWeXyw/0g3wvqOtehXrS9SXgXdfDap1z9GNeewvYbFn+dycOzOkBdDqbrcWXX/lD+/PoOr3Vn0m9scS7bXkLNerluVciF0ghSH+G0JEnWSNP05yRJNsDWsukHi0mSpCaKLiyOtS8hhBL1ngt97WPP29aO1lHTKd7GG+iIwaM45mr22ZLzDi/+GgUFvHQc1etkdYCi16FymjIiK2fwwzN02I7DX6TFOrluVcil0ghS16Eh+idJko+q+L80TZ9cwvEX4rJSaEcIIejrcT3d5HV7a6yuY1zhWFVcg10xazaHXE6jejx84ZJr/Hx6U7aU/dh3qdOkDN9AKBfy59LnNnpeloXpAx5l/aMiUIfSCVKH4WgcKZsj1RX/TZJkZJqmDxdz/DW4qcjH9TG8FNoVQqhkRvray/6it92MNNsTbne0Og6STeaEs27l+0F8cgeNl9AfPupb3ruILf5Op53KqPGh3Bjci9fP4vcf2OQMdvg3tRvlulWhvCjx8gdJkgzDtWma3lHkuYtxdJqmay3D+VH+IITwp00z1j029rF2XlDDI650la1VxWeyv9jue5WTb+CBCzhhj+KvM3cm/+tG1Rqc3IdqxVQ4DyunycN55zz6PUXbzdjrDtp0y3WryrVK2T9XGj1SdbBomYN8RFH8EEKZyDfXsw71o2peVNNFjvK8rY3CF7IQ9WV/zryFU/dZcoiCd87PJhSf8lWEqMoib3Y2lNv7SmrUy1ZobnBsbO0SilcaQaoH/i9JkqGyob0NcS5iO88QQpl4y9996yvP2sqeuqrpRC/LfjitgXGTOOgyNliNW85a8nV+eYPPb2eP22jRpWzaHnLr59d482wmDWbTs9j+Mmo1zHWrQnlWGkN79fFvHIAWGIkncUWapnOW4fwY2gshrLCvPeBZp3janqpr6Cr3OFRtl8oqBefns/v5fDuQr++lfYvirzNtLHetT+uNOOq1mFS8shs/MAtQv7xGx53Y49YIzyugUn6XlHiPVJqmU2WlDs4u6WuHEMIfGeIjrzhNT7sbJ88zrnWE2vbCpYXHXPIA73/DOzcsOUSlKa/8hTSf/R6IELUym7c33qc3Zlv+HPocax8Y/8/DsovNDUIIK4WJBnvagfra0memecp1ztVWCzwqm6T50odc8zjXncqOGy35Wl/dw889OOIV6rcqozcQylSa0u9p3v4HM8ax9T+zvfFq1Ml1y0JFE0EqhFDhzTbNk/Y1QHMv4jIneMJmRqAPGuHnYRx7DQduu+SimzBuAG+eQ7dTWXOfMml+KGOj+2Z74w3pxVr7s9tNNO6Y61aFiiqCVAihQitQ4AVH62+0p3R1sC3kO1YP2eTytTBtBgdcQttmPHjBkodt8ubw/FE0bM9uN5bhmwhlYuZEPriML+6gyeoc/Rar7ZrrVoWKLoJUCKFC+8ClvvaG5+1mTS3t5yLHSFyDvWRDOCddz9AxfH43Deou+Vo9L2fMd/zlM2r8wXGhYsnP4+v7+OAS8max83Vs9leq1ch1y8LKIIJUCKHC6usJH7ja2/aTj2td6wA1HWZB5fL/PsczH/Ds5azdYcnXGvIhH13LTldF0cWVya/v8tY5jO1H1+PZ6Wrqt851q8LKJIJUCKFCGuELLzvJp/bQ3zQvuc2pmlpTVrQuQa9vOe8u/nEYB2+/5GvNmsILx7DK1mx1fpk0P5Sy8b/w1t+zRQOrbM0pX0ZADqUjglQIocKZYoQn7ecH3bxvuvtd5HprmY5esu0VRvzOof9i2w245uQ/vt6bZzNzPMd/QJWqpd/+UHpmTqL3v7MNhuu34ZBn6HJwlDMIpSeCVAihQplrpqfs72f1PKu68x3uW7v6EO9hFcyZyyGXU70qT11KtT/4Sdf/Zb59kH3vj5VbFVl+XlZg9YNLs/0Rt7+cLc6heu1ctyys7MpNkEqSpDu6iz35QghLkEq97EQ/+NVTNrWHrlZzslNwF7YtPO7cO/hyAL1voUXjJV9v2hheOZk192XDE8rgDYRSEfOgQi6VmyCVpukduKPIFjEhhLCQD13jC895yZ7aa+xvLrOXKk7FaYXHPPQGd7zEXeew+TpLvlaaZiEK9rk3hn4qopgHFcqDchOkQgjhj/zkJe+42Lv2NV2BJ1zrILVthlsLj/miP6fdxEl7cuq+f3y9bx7IfgEf/jL1lrBVTCifZoyn95V8fkfMgwq5V+KbFv9ZsWlxCGFRo3zjAVv70DbeM9erbnSJrkbjC7TEmAlsfGpWdLPXLdT8gxpBEwZx9wascyj73V9GbyL8aXmz+fz2LEQV5GfbusQ8qHKlUkbZ6JEKIZRrU4zwhH38aH1vmeUu53pIV/3wkSxEzc3LVujNzeP5K/44RBXk89Jx1GnGbjeX0ZsIf0qa8sMzvPtPJg/LAvN2l0VPYigfIkiFEMqtOaZ7wj5+Vdsz6jjTvsbZ11N4FvP2HT73Dj7pxwc307b5H1/zkxsY+jEn9KJWg1J+A+FPG/JRtrHwiD7ZooCj3qD5WrluVQgLRJAKIZRLBfI97yi/GOIpm9ne2rZxpsPwLxxceNxDb3D7i9x5Dluv/8fXHP0d71/CVufRYZtSfgPhTxn/S9YD9dMLtO7GcR/QcftctyqExcUcqRBCufS287zvFs/aW1V13ONue6lnHzwpm4zxRX+2OYujd+He8/54snHebO7ZJDvm5M+pVrOM3khYLtPHZQU1v7iTeq3Z+RrWPYIqURinIog5UiGEUB585V4futF7DjLRLD38x+Hq6YIHZT+tx0zggIvpuhp3nL30FVvvX8L4AZz8RYSo8mjuLD6/jd5XIWXHK7ONhWMiecXwzS9suHquW5EbEaRCCOXKIO951Rm+sJ/vTPSqG/xdGwV4GbVlk8oPuZz8gqVPLofBvbK5UTtfS6ulDP+FslWQz3eP8MFlTB3JJqez3aXUXcpct1A+fP0z/3qYVz4m7Znr1uRGBKkQQrnxu5887SADbOctE93vAg/rqi96o03hcefewWc/Ltvk8tnTeOn4rGDjln8v3faHZZemDOjBexfx+w90OSTrhWq2Rq5bFpbFVwOyANXjE1Zvx8MX5rpFuRNBKoRQLkw3zhP2NlgnT8n3T0cabQ9P4GlsXHjcvMnld53DVust/brvnM/03zn2vdiQuLwY8lE2kXzYx3Tckf0fpO0muW5VWBZf9OdfD/HaZ6zRnkcv4vAd/3g/y5VdJX7rIYTyIs9sT9nfb+Z6TDMH2tzG/uJgXIpDC4+bV7n8L3stvXI5DHqfL+9iz9tp0qkU30BYJmP6ZT1QP/egVVeOfovOu0RF8org85+yHqjXP2PN9jx+MYftQNX44yRW7YUQciuVetGxPvWSx+ygndZudaud1bSnrDeqimxyebdTaN+Cnv9d+ryo2VO5cz0adyzsjYpVXzkzaSg9L+Pbh7P/HzteyTqHxf+TiuDDvlz5KG9/wVqrcOlxHLr9EgNUpYzE0SMVQsip3q70pSe85gDVVXGfa+yjprXwsCxEzZm7fJPLyYb0Zozj+A/iF3auzBjPh9dk27rUasget2ZhuNoy/P8LuZOmvPMlVz1G7+9YvzNPXcrB20UPVHEiSIUQcuZ7T3nfpT50sGGmec/tTtNEnmyFXp3C44pOLm/TbOnX/fVdvrybPe/IekBC2Zoznc9u4ePrSAvY5kK2OJea9XPdsvBH0jSbPH7lo9kw+iZr8crV7L1FDL/+kXITpJIk6Y7usj9AQwgruaE+9pLjfWdvnxjvZVe5TWffohfaFR734Bvc8dKyTy6fNYVXTmLVHdj4tFJrfijG3Fl89T8+vJqZE7N//20vjj3xyrv8fJ7rxdWP0/dXtt2At29g524RoJZFzJEKIZS5cX52vy30t56H8F9nmuZgF+MJHFF43Kc/sP3ZHLsr9/xj2X6o9ziVvo9zxvfRG1VW8ufyzQP0+jfTRtH1eLa9hMar5rpl4Y/MzeOJd7MA9fMwdtuE/zuGbVa81lqljF3lpkcqhFA5TDPW4/YwXDuPqeoM+2jpIGfjMgtC1PCxWeXyTddatsrlZKv0vrqHve6MEFUWCvKz0NrrX0z8jXUPZ/vLoxZUeTd7TtbTe92TDB7Nflvx2P9lQ3lh+UWQCiGUmTlmeNK+hprrEW3sZF2HO9MuEkfJghTMmMV+F1Ojeja5vEb1pV977kx6nEKHbel2amm+i1BQwE/P88GljOvPWvtz+Eu0XIah15A702dy76tc/zSjxmflC165mvWiNMifEkEqhFAmCuR7wVF+1d/TdtBRE9e7wk6q2QT3y8YF0pQTr6P/UD6+jRaNl+36va5gyjCOei1W6ZWWNOXn1/jgEkZ/y2q7c8CjtN14qaeGHJoyPZtneNMzTJzKMbvyzyNZc5Vct2zlEEEqhFAm3vJ3fb2qhwOR7wnXOVhtDfAi5u0jfM3jPP0Bz/2Lrsu4Cero7/j4era/jGZrlk77K7M05bf3ef9ihn9Gh+044UM6bJ3rloU/MnYitzzPHS8ycw4n7sEFR7Bq61y3bOUSQSqEUOo+9V+fukVvhxtmsg/c5q+aGYXPMK+iwcsf8X/3cdlxHLTdsl27IJ9XTqbZWmx1QSm9gUpqXoDq+S+GfkjbTTnmHTrtFKu5yrPfRnHj09z/OlWrcNq+nHPI0velDCsmglQIoVT96HlvOtc3DvGZMV5znbt11gtvY14H0veDOPoqDto2q568rD6/nZFfctLHUeixpKQpg97NAtSwj2mzMUf0YI29IkCVZ98P4roneOp9GtXjoqPpvj9NGuS6ZSu3CFIhhFIzzKdecLRB9vKqse5znr42cR8exA6Fx42bxL4X0blNtov8ss5xmjSE9/6PTc6g/Ral8x4qkzTl17ezADX806wH6sjXWH2PCFDl2cffZ0Pir33GKi25+cxsGK9u7Vy3rHKIIBVCKBXjDfSkfY2wicdM9X+O1theTsZFOL7wuLl5HHwZ02dle+gt6w//NOW1M6jViJ2uLo13UHmkKQPfzALUiD6025yj3mC13SJAlVdpmm0gfO0TfPQ9XVblkYs4fEeqx2/2MhX/3CGEEjfdOI/bwyitPKimo2xnXyfZHgfj30WO/eutfPID799Eh1bL/hr9X+aX1znsRWrF0MUKSdPs37Dnvxj5Be235Oi36LxLBKjyKi+PZ3pmAer7QWzehZevyrZxidWquRFBKoRQouaa6Un7Gm6Wh61hc6u5zPm2lVjfgo2I4a6XufsV7juPrZejmvLcmbx1DqvtwVr7lcKbWMmlKT+/mpWMGPklq2wdk8jLu5mzsyKa1z+VFdHcYzNu/1tWhTz+n+VWBKkQQonJl+c5R/hVP8/YRUt1PexK+6ihlmwj4nkjdx98w1m38NeDOGmv5Xudj//D1JEc83b8ElkeBQUMeJneVzLq66x46bHv0XGH+HcsryZN5a5X+O9zjJucFdF88d/LXhoklL4IUiGEEpFKva67H7zuLQebaaZ3XOd09Q3Bp2hZeOygkdm8qB025MbTl+91Jg7mo2vZ4lyaxi+TZZI/l++f4KPrGPcTq27PcR/QcfscNyws0YjfufWFrNd2zlxO2IN/HEbntrluWVhUBKkQQono5d++cI9PHK2/sXq6xS1aextvoEvhcVNnZCv0mtTn6cuotpw/hd46l9pN2eb/SvgNrITmzuTr+/nkeiYPZY192O/+WOFYnvUbxA1P88R71K7JGftx9sG0aprrloUlKTdBKkmS7uhuwfSJEEIF8ZV7feAy3zlKbyO87CofWcsd+B92KTyuoCCrFTXsdz67c/nr2/z6Lv1f5OCnqFmvhN/ESmTWZL64k8/+y4xx2WbCW/8z9sIrr9KUnt9m85/e6EO75lx7Cn/ZiwZ1c926sDRJmqa5bsNCkiRpgMlomKbplOU8vXy9mRAqgf5e8bQD/Opgjxrjfuerb0+H4p+4psix/3dfVu/m1WvYc/Ple52CAu7dhGq1OPGjmNNTnGljs/D0xR3kzaLrCWx1Hk0657ploTh5eTzXK9tE+OufWb8z5x3GYRW3hEGl/K6smP+rQgjlwlCfeM5hRtvDo8a4wonWsKedcQSuKnLsk+9x9WP857TlD1HwwzPZBOkTPowQtagJv/LpzXxzP1WqsfFpbH4ODdrkumWhONNm8MAb3PxstgJv5268dT27bBxf2xVR9EiFEFbI735yv62Mtr7bVfEXeznbubaSWA9vWbAR8Zf92eavHLJ9Vrl8eX9Z5M3hjrVpvg5HvlKy76MiG/oJn97ITy9Spxmbnpnd6jTJdctCcUaP57YXswnkU6ZnxTP/cdhKtQKvUsbA6JEKISy3KUZ41G4m6OBeNe1lY5f4m20kWuJFC0LUqPHsfzEbdOaev6/YX9xf3cOkwRz+csm9h4qqID8rRvrJDdk2Lk3XYO+72eAYqseWIOVS/yHc+AyPvE2Napy8dzaBfJWWSz83lH8RpEIIy2WmSR6zu99V94B2NrCK+1xiL9XMxAdoXHjsrNkccHHWVfzildSqueTrLsncmVndow2OpeW6Jfc+Kpo50/n2oWwIb+KvWQ2oI15h9b2ionV5lKb0+pabnqXHJ7RuyhUncOo+NKqf69aFkhRBKoSwzOaa5Sn7GWmsJ2ytubpecrW/qOVH9EKHwmPTlFNu5Ltf6X1r9otkRXx9PzN+Z9tLSuhNVDDTxvD57dkqvFmT6HIwBz1Bu01z3bJQnDlzeer9bP7TtwNZZ1UeuIAjd6JmjVy3LpSGCFIhhGVSIN8LjjbIV162r7lme8N/XKGhV9ED3Yocf+PTPPo2j1/MJmut2GvmzcmqmK97BE06lcCbqEBG96XPrfR9lCrV2egvbP43GnfMdctCccZN4n89uP1FRk9g9015+4ZsInlMIF+5RZAKISxVKvWGv+rnJe87ylDj9HaL57SaXytqzyLHv/4Z5/+PC4/iyJ1X/HX7PcmUYVkNpMogP4/+L/H5bQzpTf027HAF3U6hduOlnh5y4Kch2fYtj7yVfXzsbtn8p7U7/PF5YeURq/ZCCEvVy5Xed4kvHO9dw7zpP8bZyGG4EFcXOfanIWx+BtttwEtX/rn5O/dskq1GO/qNP/kGyrnpv/PVvXx5F1OGs8o2bHYWa+1P1eq5bl1YVJry7lfc9Axvfk6rJpx5QDb/qVmjXLcupypl31v0SIUQ/tAX7vKBS/ziWG/4zZMuUd1GjsGRFq4VNXFqtv1L++Y89n9/LkSN+IKRX3JEjz/5BsqxEV9m85/6PZUN/6x3VFa+oHXXXLcsFGfWbB5/N+uB6vcbXVfLynkctkPMf6rMIkiFEJaon6e9prsxjvSkIW7S3YZ2sgW2xAMW/Amal8ehlzNhKl/c/ee3tvj6fhquwup7/LnrlDd5c/jxuWz4bvhnNOyQDd9tdBJ1Yj+1cmnMhKz2050vM24y+2zJbX9lu64x/ymUUpBKkqQtrsMeqI2fcVKapl+VxuuFEEreQG95wTGmOsD/jPR3hzrCIbZEa7xgQa0o+MddfPAN79xIpz9ZUTt/bhY2NjqJKlX/3LXKi/ED+fpevnkwW4XYcScOf4k19l553uPK5vtB2eq7x9+lWlVO2J2/Hczq7XLdslCelHiQSpKkMT6WlZPZA2PRGZNK+rVCCKVjqE887UAz7Ow2kx1lZxc7zY6YZeFaUXD/a9zyPHeczQ4b/vnX/+0DZo5nncP+/LVyKW92Nnn8q3v47X1qNc7qYXU7hRZdct26UJyCgmyxxC3PZ/Og2jbj3ydmRTQbR/2nUIzS6JG6AMPSND2hyHODS+F1QgilYIzvPWEvM23sNvl2tbG7XWB/VfyKDy2oFQUf9eX0mzltX87Yv2Ta8Ovb1G9L6xIIZbnwyU38+lY2dDd7SlY884BH6XJQVB8vryZP48E3svIFv47MSnY8cQkHb1dhNxAOZaQ0vjz2xVtJkjyL7TACd6Zpem9xBydJUtPCIwSR+UPIkQkGedRuZujkLg1sZFVPudwpqumNN7F+keOHjuHAS9lyHW45q+TaMaQXHXfI3fyT33oyfUxWTXzgm1n5hTYbZRPg+z6eBbyhH7PNhQvqOs2anG2s/OnNjPuJ2k3osB1bnEvH7XPzPsLS/TQkC08Pv8nsuRy6Q7ZQYrMuMf8pLJvSCFKdcDpukq2K3hS3JkkyO03TR4o5/kJcVgrtCCEsh6lGe9SupmrgPqtop4kerna5Wp7EU9ihyPHTZ2Yr9OrW4rl/UaOElukXFDC2H+sdWTLXWxHPHMRuN7HRicyayLOH0v0HnjmEk/tQr2W2x91Lx7PVBXz3CANeJn8OLTdgzzvZ8AR+ey87LpQv84bvbnuRt7+gRWP+flhWvqBNs1y3LlQ0pRGkquDLNE0vKvz4myRJ1pGFq+KC1DWy0DVPfQwvhXaFEJYg2z9vNxPN9bhuaqvqTf/xgPpuxK04tMjxBQUcdw0DR/DpnSVbO2fqSPJm0WT1krvm8jq+J42LVFKfMy0rkFmzQRaiRvfNep+G9M5uzdfJVt6tdxQN2tDvaX54OluR1+BPTrwPJWfR4buN1+TRizhk+yhfEFZcaQSpUfhxked+wkHFHZym6WzMnvdxEn2pIZSpOWZ4wt7GGulFu5tuuo/drKem/i6b9LjoqN2/H+H53rz4b9Yr4a1b5kzN7nNZybvlegseD+jBLtcx4vMsUN3VlTHfZYVCq9Vit5vZ+NSFh4HWreCT5Fc2MXwXSlNpBKmPseYiz62BIaXwWiGEPyHfXM86xDDfedshhhijt1v9po3jcKysy7io53tx+UNceRL7b1PybUoKSwHkzf7j40rb8M/p+xhSPruF0d+QVGHtg7Lep9X34JZO1G0Rv4zLoxi+C2WlNILUzfgkSZKL8IxsjtQphbcQQjlRIN+LjvWzd3zqeN8Y7B03KNDZAdgR91l4z4dvf+HYa7JKzhcdXTrtarRqFqbG9afTjqXzGkuSpoz+NuuF+rlHVlm9SjVq1GeD47N9/w59ZsHxc6ZlPVOh/Bg/mYfezApoxvBdKAslHqTSNP0iSZIDZH/IXorfcHaapo+X9GuFEFZMgQI9nKqfZ/Vzivf197KrtLW+LWVdys+h6PzxsRPZ72LWWoUHLii9XphqNWi/BT89z6ZnlM5rFDV3ZlbjaUAPfn6VqSOQZD1OBz2RTRa/Z+Ps44+KdM/lzWb2VBrF5rTlwhf9ufMlnnqfgjQrWxDDd6EslEp1jDRNX8WrpXHtEMKfk0q95RzfeMAQZ3pBX4+6yKa2sBXq4jXUK3LOnLkcdCmz52QbEdepVbpt3OxvPHsIfZ9g/RJevVeQz6hvGPxBVvhzcE/yZtK4M+sckk1y7/s4h72YhbqfX6dK9azkwfSxTB5Ow3YM7kXbTSJI5dKMWTz9QRagvhxAh5Zcdhwn7pkN5YVQFqLMWAiVzPsu0cetJvqbh3zrv850gF3tiMn4BC2KHJ+mnHEzn/en539p36LYy5aoLgex/tFZeYGZE9j4NKqu4E+rWZOz4DTyy6w+1ZDeWZHM6nVovxXbX86a+9BsrQU9F/Va0efWbOuWoR9x5Ks0XZ39H+bDq2i3eRakDox+9pz4ZTh3v5KtwJs0jd03pcfV7LEZVWO7nVDGkjRNc92GhSRJ0kD287xhmqZTlvP08vVmQihnPnSt91xolrNcp6//c7TL/MX+6I2e6LbIObe9wF9v5aF/ctzuZdfW/Lm8eTZf3JnNm1rvqKxIZ/N1shIERYdr5s5k5kQmD2Xir0wcxNgfGPUVEwZmx1SrTfsts2usuj1tNsl6nELFkJfHa59lGwe//QVNGnDSntnk8c5tc926UKhSDqJGkAqhkujjdm84S77TXeUn3e3vFn/1F4lHZcN5uy5yzrtfsvv52UatN5bBfKXijPqWz/7LL68xY1z2XJVqWemBKtWZO4P8RVb41WlOszVp3Y3WG9GmG03XXPFerZA7YyZw32v8rwfDxrLZ2tlWRIdsT+2aSzs7lLEIUuVBBKkQSt43HvSyE3GSKw1ytF084AL/VMX1eAxHLXLOwOFsejqbrsVr1+Z+yKSgIFvJN2FgNiE8b3ZWSbx6HWo1ym4N2mWFNGs1yG1bw5+Tpnz0fTb36fneVKvKkTtx+n50W7S4TihPIkiVBxGkQihZ/TzjeUeo6ij/Nso+tvCUS92kmgtwC/66yDlTprNFd+bm0eeu2PU+lI0p03n83SxA9fuN1dtxxn7ZkHJ8DVYIlTJIRUd3CCuxAV71gqNUd7CrjLWjDT3hEo8UhqiLLR6iCgo4+iqG/06fO+MXWChdaZqVLrinB0++z6w57LslN3dnx42oUiXXLQzhj0WQCmElNch7nnGwGvZ0rWk2tqbnXeF11Z2M03BFMedd+gCvfsqr17BWLO0PpWTyNJ54LwtQ3w7MVoP+80hO2J12ZbAyNISSEkEqhJXQUJ940r5q2M4NWEN7PVztc7UcjoNxu8X74Z/5gKse49pT2HPzMm92WMmlKZ//xD2vZoUzZ89h7y246i/stknu5+GFsCLKzRypJEm6ozuqyAorxxypEFbASF972I6q6upWTTXWQE//NVRD22EzWbXcRRc8ffMLW53JAdtkFaGjGnQoKZOnZXOf/teDvr+ySktO3osT9qBt81y3LpSgSvlTo9wEqXlisnkIK26Ubz1sR1Ws4S4dVFPNh241VVNbowPet3DVcrIl5pucRotGfHhbLCsPf16a0ufHBb1Pc+ayz5acsg+7bhy9TyupShmkYmgvhJXEGN97xM6q6+w+q8k3R083KtDUrmiC1y0eouZt/zJnbrb9S4So8GdMmpr1Pt3zatb71KEl/3d01vvUJjZ4DiuhCFIhrATG+tHDdlLDKh62rgkm+dCtGmhlW+ThLSz6eyxNOfMWvhiQbf8Sk3zDikhTPv6e+1/P9r6bM5d9t+K6U9glep/CSi6CVAgV3O/6e9iOqmvtWZsaaoSe/qud9nbBKHyIVYo5986XuPdVHriALdYp02aHlcCo8TzyFg+8wc/D6Niai4/Jep9aN81160JZKpBNcK6MIkiFUIGN94uH7aiaZl60tX5+844brK2z/fCdbE7U2sWc+8E3/O02/nZQ9osvhGUxN4/XP8t6n17/jOrVOGhb7j6X7TaIuk+VSb7sj7Tn8QJG5LY5ORNBKoQKaoJBhavzGnrVjr72szf9xya6OBbvyuZEbVrMub+N4pDL2GFDbji9bNsdKqb+Q7Kep0feYsxEuq3BbX/jiB1pFEVbK425ss3Nn8NLGIv2OCR3Tcq5CFIhVEATDfawHVDHm3bzqZ+87lpbWd85eAJPY+dizp0xi/0vplE9nr6MavFTICzBtBk80zPrffqkH00acPQunLgHG6yW69aFsjJb9ofZ83gZE9AJx8lq0m2iki7XKxQ/QkOoYCYZ6mE7SFX3vr310lcP19jehv4t2zvvLsX/hZimnHwDA0fw2Z3ZL8YQikpTPv2hcOL4+8yYnU0Yf+pS9tuKWrGqs1KYKVug8hx6YIqswOPpOAhdVe7wVFQEqRAqkClGeNiOhXMTDvC2b7zkSrvY2C24FFfKtn8pzq3P88S7PHkJ63Uqs2aHCmDoGB59m0feziaOd2jJeYdz/O50aJXr1oWyME02HeB5vIbpWA/nynqeuojwVJwoyBlCBTHVKA/azlyzfeYIL/nCc/5lP1t7ECfiPFyn+B92vb5lp3P528HceEZZtjyUV9Nn8nxvHn4rW3xQuyYHbsNxu8WGwZXFZNlOB8/hTczCRrLgdBDWWL7LVcqcFUEqhApgmjEesr2ZpvrGcZ7xqadc6mDbex6H4mTZkF5xP8lG/M5Gp9ClA+/cEPOiKrOCAnp/x0Nv8lwvps/KVtsdtzsHb0f9OrluYShtE2RznZ7HO5iDzWXB6SB0XPFLR5AqDyJIhbCwqUZ72I5mmKSfkzzmI4/5P0fa2VvYRzYf6hEUV/dw9hy2P5thY/n6Xlo0LsvWh/Ji4PBs2O6Rtxgyhk5tsp6nY3bN6j+FldsIWXh6UbbqLh9by3qeDpCtvCsBlTJIxd+lIZRjU4z0sB3NNs3PTvWoDzzoAkfa2UeyH4C74SHFhyg4+3a+/oUPb40QVdlMnpatunv4TT7uR4O6HLp9FqC2Wi82pl7ZDZAFpxfxuewX/g64VfazI6a+lYwIUiGUU1OM8JAdzDXTr053v3f8z98dbw9fYy9Zd/wzqL6EazzwOne/wr3/YNPiqnKGlU5eHu9+lc17eukj5uSxSzeeuIT9t469FFdmKb6U1Xd6ET+hDnbHWbKfGfG3VMkrN0N7SZJ0R3dZlfk1xdBeqMQmG+YhO8g311Dd3eJ1t/ub7g7wE7aV1XF5F0uqhfjtL2x+RjZ0c+95Zdb0kANpypcDePwdnno/K5jZZdWs5+noXWKz4JVZHnrLgtNLGC7boHxfWa/TLqhdds2plH2c5SZIzRNzpEJlN8kQD9lBKjXCmW7Uw43OcK5DDZbNa2iMXrIfmMWZMp1up2QThz+5PWr/rKwGDufxd7PbL8Np1YTDd+Sonem2Zgzdraxm4G1ZeHpVNnm8vSw47Y9t5Gy4qVJ+xcXQXgjlyLyK5alkfoi63mnOdahRskrltWQ/RJcUouYV3RwzkTeuixC1shkzgac/yMLT5z9lYfnAbbjjbHbckKpLmiwXKrSJstD0oqxQ5gxZXafTZQFqI5U0xZQDEaRCKCeyvfN2kKhupDPc6BXXO80/HG4CdpXVePkIf7TI6u5XeOYDnrmc1dqVRctDaZs6I5vv9Pg72fynKlXYY7Nsi599tox5Tyur4bKVdi/JVtrlyeZFXibreVrOGk+hlMTQXgjlwHgDCzcgrmWY093oZTc43d8dZqqsJ2qQbKf1tf7gOl//zBbdOXkvbj+7LFoeSsvcPN76POt5evljZs5mm/WzYbuDt6Npw1y3MJS0FN/ilcLb1xastDsA+6FNrhq3bCplp1gEqRByafxrxtdr4KGaR6iunqFOXShEzcSesh+oH8i675dk8rRsXlTDetm8qJo1yuQdhBJUUJDtc/f4u1mv4vgprLMqR+3CkTvFVi0ro9my+Y4vy/a0G4aGsu/7fWUr7hrlqnHLr1IGqRjaCyFXxjxh3OBjPNSttppWMdQpbvTS/BA1V1axvI9sTtQfhSg48xbGTuKt6yNEVSRpms11eqYnz/bMCqe2a86Je2a9T+t3jknjK5sJsj3tXpFtyzIVq8p6nfaVTRaPb+GKI4JUCLkw8j6/Dz/Zw91qq1V1VYP9xY1emr86Lx/HySaV9pCt1PsjT73HY+/w6EV0blvqrQ9/Uppmw7BPf5D1PA0ZQ8vG2ZDdoTuw9Xqxz93K5ldZr9MrsnmO+dgEF8jC07oqaXfOSiCCVAhlbfgtRo0+26Mb1VK36mp+S05cKESlOANPy4pt7raUyw0by+k3c9gO2RBQKJ/SlL6/LghPv46kWUMO2pbDdmTb9WPF3cokX1ZNfN58px9RUzbf8U7srdzPdwrLKOZIhVCWhlxt2Pj/8/gGNTWusr6BydFu8IKbdHeOQ6Syv1Cvx4M4fimXKyhg579nNYT6PkDjJVXnDDnTb1A2bPf0B/w8LPt/dOA2WXjaoWtsIL0ymSHbBPgVWamCsWguC037yopj1s1Z68pEpexUi2/hEMpCmvLb//ltyjWe6FpDq2QTA5LD3OC5+SEKrpaFqFssPUTBzc/S81veuylCVHnSf0hhz1NPfhxMw7ocsA3/PZOdu1E9fvKuNIbiNVlwel9WomQt2ffvvrJyBdHRuHKLHqkQSluaMvBsv8y81dPrVtO+yvZ+sL/rPbtQiLoNf8UVuGQZLvvLcNY/kdP346bupdj+sFRpmgWmFz7kuV7ZEF79Ouy3VTbkusvGsQBgZZGHzywIT/1kPRLbyPay20elru9UKXukIkiFUJrSfH4+zY9z7/PcOlWtluzlczu5xQtu1t3ZhSHqEdnk8r/LeqSW9tOooIAdz2HY73z/AHVqle7bCItLU74akIWnF3ozYBj1arP3Fll42n3TqCq/spggW133Gt6QVRlvLitRsJesWG6U9UIlDVLRwRxCaSmYw0/H+K7Ks15ap4q1k4P1tIm7vOAOZzvD/uB5nIC/WLYQBfe+Sq/vsiG9CFFlJz+fT37IgtMLHzJ0DE0aZD1PN57BThtFeFoZpPhB1uP0Gj5BATbEmbLwtAliYWUgglQIpSN/Ov0O8mXtd726Ohskx3nNWh7wqnv9w1/sjay8wRGyelF3W7YQNXYiF/yPk/Zkx6UVlwp/2tw8PvgmC08vfZTtYdi6aTbn6aBts9V2MWG84pspK3o7LzwNlU0M3xn/k/U+xSq7UJz49g+hpM2dyPd7+bThV97qnK+bMz2rtce96WEXOsauyLZ7OUBWufgRyz4h9ZIHshpD151aOs0PTJvB219mW7O88jGTptGxNUfvwkHbsdnaUedpZTDMwhPFZ6KTbCuWvbGdrGRBCH+k3ASpJEm6o7voLQ0V2exR0r676t3iVx90mGNz53tYXc961xMudpgdwZey4YEtZLWiqi/j5b/5JRvWu+Ws2GutpI34nR6f8MonvPc1c+bSZVXOPCDredpgtagwXtHlyXYKeK3w1lf2S3Br2SKPvbGmSjrRJ6ywmGweQkmZOUj63c7ebfe7j9tNs40r3GGOV33qKZc60LbI5l5sh9VlNWfqLePl05Ttz2bcZL69L5bQ/1nzCmS+8knW6/TlAKpWYdsN2HdL9tkyqsSvDEbKJoq/Kft+m4RmFp4o3ihHbVsJVcoMGj+KQygJ075X0HdXb3Sa6YtW02zvBjca521fesEV9rYlsm0idkFb2V5byxqi4P2v6f0dr14TIWpFzZmbTdJ/5eMsQA0dk5Up2GMzzj6YPTePelwV3RzZ5PA3ZSvs+sp+u2+Gs7EHuonaTqHkRI9UCH/W5E/l99vDS2tV8X2TSXZJ7nCFQXr7zsuusqtNwHALNiPtjZbL8RJpynZ/Y+ZsPr87hpiWx+jxvPk5b3ye3U+Zziots16nfbdiuw2osaxjq6FcGmpBcHpPtglwC9n8wz1kf7w0zVnrKpVK+ZMp/q4N4c+Y8LY5P+7v2fVq+bXBNHslj7nAF/r4yWuutaNsWd1Y2Q/zArxr+UIUfP4TH/bl5asiRC1Nfj59fuL1z7Lw9PXP2b/Zxmty7iFZqYKY71SxzZYt1nhDFqB+lPUwbSHbYmkPdBUTbkPZiCAVwooa+6yZPx/pyQ3rG1Vntn2S5/3V27410Jv+YxvrI5uTsVvh/YdovwIv9b8edGiZFXvMlT4/Zsv+2zZj9IRss+TN11n8uCffY9pMVm/L9huWTdvGTizsderDW18wcWpW32m3TTjnYHbblOaNyqYtoXQMsiA4vS/b166NrNfpX7IyBY1y1bhQqUWQCmFFjLzX1MGneGzjRibXZN+kh5M9q7+h3nGDzWUJY5psUutQ9MJqK/BS02bw1PtcdFTZLLnv8yMffc/UGXz6Y/a623Xlrpd5+K3smI3X5PGLFz/3lY/ZvmsWuF7ona2Ea9u85NuYn88X/Xm9TxaevhywoF1nHpDNedp0LarGRJgKa7rse+YtWYD6xYIVdpfJAtR6KulYUihXIkiFsDzSlCFXmTj6Eo9s3FBe9Tr2Tp51pPuMNM57btLNmsg2L90f38v+gl53BV/y/W+yuVGH7Vgi7+APzZiVFZ285pTs4+d6sscF/PIYq7Zi+LPZP0G7Fku+RpXC32w1q5fc8FmaMmgk736V3d7/hglTsonhu27MWQdmvU8tm5TM64WyV4Bv8Hbh7WPMlfXg7oH/YEc0yFUDQ1iCCFIhLKs0n1/OMmbyXR7duL4aVZvbOXnC/m420xy93aqLVZH9Ajhc9svgTQqnm6+Yd7+iUxtWb/en38FSDRzBtU9kVdNXa5cNic2czcf9ss8vrXdp3614+n1mzcmGIts0W/G2/D4pW6k4LzwNHp2VJ9isS9brtOvGWWHMqCpecQ2TlSR4WzZ3cLxsJev2uFFWmmAN0esUyrf4ERTCssifyU9HGTbnZY93q6NRlc42dp89XaO2Gj5ym45aI/vL+gRZeYOXZTWj/oyfhtB1RcYEV8B6nfj49gX1k4aOye5Xb5cNpT34RtbT9O5XnHc4a3dY/Bor2nM2fWYW2N79ine+5NuB2fNrd8gC2s7dshV2Dequ2PVD7k1DTwvCU39ZSNoEp8mC0+ayla0hVBQRpEJYmrkT+H5fA6t/4ekNq2udbGQNN9vdv7TRzFuu17pwcXUqK8//JJ6SDUn8WUPHlt0k8yRhyyJjkNc+kdVX2nD1rAbTZl2y51s04oCL+fHhFZ+3NXlaNherd196fctXP5OXn82v2rkb5xySbQJcGnOsQtnIx9ey0PSOrL7TXHSQhaYrsBNiRDZUZKUepJIkuRBX45Y0Tc8u7dcLoUTNGkbf3fVrMNQLaxbonOympUvs6VLrWNVrrtWkcNZGin/KNh9+AIeUUBOqVc0mV5e1B16nVRP+c1r2cbc1FnxutbYMGJZVBu+6+rJdb/zkrIRDr++y8PTtQAoKsuC03QYcv3tWVXztDlGaoCIbYkGP03uYgPrYATfLAtRqYrgurDxKNUglSbIJTpEVlw2hYpnWj76769N6hjdWnW695EjVnGofl9jKul5ypXrqzD/8P4W3/8qG9kpKu+b8MqIEL7gMXv+M/AKuP51Zs7PQs/PfGfcytWoydWZ23JIKWRYU0H9oVs/psx/5pB/9fss+16Fltgqw+/5su342jBjBqeKaLFtdNy88/Syr37SprHd2V1lV8ah5GlZWpRakkiSph8dxMopZKB1COTapt/T7fby3Wi0ftZ5oC383zm6Od7l9bOlJl6hZZCbHfbLeqEvxtxJuyu6bcsH/GD72j1fLlZTe3zFyHPtskVUF/+zHbIPk84/IQhR8/H02BNhl1ezj3ydlZRPmBafP+2cVxJOEdVZl8y7Z+duuT4dWpf8eQumZJVtE8b6sx+lL2RDeqrJ6aVfLVtc1zlH7QihrpbZFTJIkD2NCmqbnJEnSE98WN7SXJElN1CzyVH3ZbhqxRUzIjd+fl//TkV5Zt4nvmoy2qxv11dkZ/utYu7rPeaoV+RvkBdkw3mm4XckPWUyZzmpHsX4nXrmaOrVK+AWKGDSSrn/JakgVNfm1bA7TVwOYPJ0vBmQr5n4ZnoWnQSOz41o0zkLTZmtn9xuvGZPDK7o8WViaF5w+llUWbyYLTDsV3ncWw3Whcn4JlEqQSpLkcPwfNknTdNZSgtTlsvpqi4ogFcre8NvN/u0sz27Y2qC64xyQPOw1BS50r786yM26q1Jk44n3ZRPKD8ATSm9Lip7fsNeFrNeR+89nnY6l9EJFzJydBaWfhvDdr9nw3rcDGTU++3zdWqzfOSt8uXmXbCL6qq1imK6iS9FPFprelw3bTZH9hbudBeFpXbEFS1hMpfzuL/EglSRJe9kfMLumafpd4XM9RY9UKM/SAgZdZNro6zzRraVxNWc4LHnBvX51nSdd5jiXOV5S5OfEl7IJtFuih9Jfsv35Txx1Jb+OzCp3H7xdVoSyddMVCy9pyqRp2VYvw8ZmdZp+Hs6Aodn8pqFjs2PI5mlt0Dkrw7DBatl95zZlU2k9lL5BFgSn92V7Q9bAVhYEp43FPKewVBGkSuSCSbI/XpQNm89TVRZyClAzTdMlrkFKkqSBbP5iBKlQNvJn0f84E6Y+49FuTc2tXt3hXnOlD/xPDzfr7uxF1uANkG1VsZqskGBZjV7NnsOjb/PQm3zyQxZ0WjRm7VWyFXYtGtO0wYKAk8gmjU+atuA2cSpjJmbhafqsBdeuVjULR2utwpqrFN63z25NG5bRGwxlYrQFQ3XvY7Csd2ljC4bqtkLtHLUvVFgRpErkgklSX1YmpKgHZbXXrkvTtN9Szo8gFcrOnHH0299IX3i8a201q7RwuNf8zZOe0dN9znPCItWghsl+ydSXbUKcqxo4o8dne+F9OzDrRfp9crZ574SpWcCa961dpQoN62bbqTSql92aN6R9iyK35lnPVuxNt3IaIxui61l4+6nw+XUsCE7biU1/w58WQarUXuQPhvaKOTaCVCgbM37h+z0NrPe7p7vM0SJZ3/6ed7zbvONLT7nUgbZd6JQJsp6oGbJJt21z0OwQlmZJwWkN2fYr8+Y6xQLKUMIqZZCKyuahcpr0Ef32913rGl7uNF3nZDe7edBBrvSlAV51jV0X2SFv3ibEY2UVmiNEhfJitIWDU//C59eUhaZLCu/b5KBtIazsyqRHanlEj1QodWOekvY/1kert/Nem99s6ESbuc5eLjLISK+51pbWXeiUAhyBV2RzSspox5YQivVHwWl7C3qdWpd5y0Jl9btJmmsUPVIhrNTSlKHXyh98kdfWX83XjQfazqVW090OzjXBVD391wYW3yH4AjyL50WICmVvlAXBqZcFwWktWWC6TASnUHZGGe9rP/vaL74ywNd+McxYqZ65blpORJAKlUPBXH45w+yx93l2k84G1R5sPw+oZxfbOEuKj9xmde0WO/U23IBbZPWiQihNqWy/ug8Lb71lq0TJgtP2uFwWnGKOUyhNqdRwvy8WmkbJisk1Vl83azjCTjayjJturoRiaC+s/PKm8MMhpsx4zxMbtzWx+iSHet50He3qHxqr7x03aGfx/VdexEE4BzeWdbtDpVCAHywITh/JCumRrarbRlavbFsRnELpSaUGG+1rP/vKz/PD0+8mgeYa6WYNG1nDRlbXzRo6aLVQbT2VdLJ5BKmwcps5mH77GFNlsMc3rE2VWo7ymkGq2NMFOmnjTf/RvJiF39/IyhzshadFFedQMubIirnOC00fY6JseGBj2arQbWRfe01z1MawcitQ4FcjFwtNE00FbTSzkdVtZI3C8LS6tpovGpqKE0GqPIggFUrM5I/pd4Bfm1T19FpTNUlWd6TX9DHS/i7WzRp6uFpD9RY7dbRs9/qWsqGVKEwYVtRUfGpBj1Mf2QrQurL5dtsU3jZDnRy1May8ZpvjB4N94xffGugbv/jOr6aZCVbRcn4P07zeplYrHuEjSJUHEaRCiRj9KAP+4utVO3h1ld90SnZxiKe95htH+LedbOQ5/1LH4jsAz5YNpQzGF6LMQVg+Y2Q9TfOC07ey4btmssA0r8epq9hyJZSsiab6zkDfGOjbwtuPBsuTL5FYU3tdrWZDq8+/L643/k+IIFUeRJAKf0pawG8XS4de44N1N9C72Xe6OdWebveod5zkeofYziMuUqOYX2MpTsBTsp6oTcu4+aFiycePsrpin8h6nn4p/NyqFvQ2bSMrTVApf8uEEjdvEnjRXqZvDTTYaFBLDevrXBiWVtPVatbTSd3S71uvlF/iEaTCyiN/Oj8dI2/8i17ZeCN9635tZ9fZynlu9byz3e4U+7jT2aoqfi+UO9Edj+Gosmx7qBAmy4bm5oWmzzBFtpnoBrINrLeUBafF13+GsPzy5BlgWJHQlPU0TZD9emyigQ0X6WVaQzvVcrMoP4JUeRBBKqyQWcPot68Zc3/x9MYdDa/+swM8Yh2HusLDLveQ8x3hWqcsccLkl7IJvqfi1rJseyiXUgy0IDR9gn6FzzeRzW+aF5w2UXYbV4eV13Qz9TVofmj61kDfG2SWOaCj1osMza22rJPAy0q5aUhZKjdBKkmS7rLOgCqyXvAIUmHZTPmcfvv5vXYVT2xQ1ewqMx3uJe1s4Vx3uMXzrnGyf/5BH9MEbCSbXP4hapRR00P5MUMWpueFpk8wrvBz61g4OK2hkv7GCCVi3tBcX7/6zq/6+tW3BvrZcKlUNVV1ser8YbkNrW4DnTVSP9dNX5pK+W1RboLUPNEjFZbLmKcYcIJBrTt5ZrUR6idtHaGHBlZxshs87C13ONvp9lviJVIcItv65Rt0KKOmh9xJMUg2TNdHFp6+QR7qyVbQzQtNm6FxbpoZVgIzzPKDwQuFpr4GzS810FBd6+tsA53n9zStY1U1K+afcxGkyoMIUmGZpPn8dglDr/Hlmpt6rdVXOiU7O8TTErUd4d9e8bFHXORIO//hpZ6QzYd6FgeXRdtDmRuPz2Wh6fPC2/jCz3W2cG/TuixhBl0IS5ZKDTN2obD0nYF+MUKBAonE6trZQGfr62QDq1lfJ6toWZ6G5v6sleaNLI8IUqHiyZvMj0cqmPCGtzfa2mcNPrSpM+3mZjPNcaBLfaivZ11ub1v+4aVGyH5x7onHy6LtodTNkpUcmBec+uDXws81la3E3KzwflNR9DIsvxlm6ee3Ir1Mg/T1q0mmgUbqze9lmhea1rFqseVWVjIRpMqDCFLhD03vT7/9zC4Y47luaxpY40u7u8VmzjTRVHu6wA8G6+Fq2+m61MvtL/uF2082gThULAWycgNFe5u+xVzUxIYWhKbN0Ekl/UkfVkgqNdSY+b1M80LTL4VzmaqoYnVt5/cuZcGps/ZarEy9TMujcr7pCFKhwhj3Kj8dZVL9Fp5Yv6rJVUY52NNWt7vRxtvN+Ub43Zv+Y2NrLfVyb2IP2fYvh5Z228OfNm8z36+K3D6ncCewbIXKvMC0GdYXiwbCsptiun5+m9/TNK+XabLpyHqZ5gWlefeVpJdpeUSQKg8iSIXFpClDr+G3iw1rt7WnOvdXPanrSK9qYR2DjbKLf5hhtnfcoItVl3rJAtkv2mb4QPn97n8S07A6ts9tU8pUiqGyVXRFg9O8eU1t0M2C4LSxmBAels1sc/Q3VD+/+d6gwvvfDDUGFuplKjo01658lRkoryrlP1BOKnaFsMzyp9P/BH5/1vddDvJS81e1STZ2uBfV1Vx/Q+zs72qp4SO36aj1Ml32BfwgW+Jemt/5fWTbhczbb+0ibCfbeuZx2dDTx7gQHRc59xVZeGpd2N4RVs7tauaFpnlhaV54WjQ0nSULTN3QquybGSqYAgUGGTk/KM0LTj8bJl8BaK+FdXV0uB2tq6P1dLSWVdRSM8etDxVJBKlQfs38jX77S2cO1GuTw/Ss+7T1HWNf96qmpm/8Yjfnaamxt92g9XJMG74eO8pWa5WWGXgJ1xR+/JxsKPEnWbmFPrK6VWvItqXpWcw1qhTe17Ry/KlXICs78B2+tiA8zavX1NqC0NSt8LZs0ThUVqnUGBMWCkv9/OYHg80wCzRW33o62dGG/upA6+lkHatWhLpMoQKIob1QPk18nx8ONadGfS9vuKYfqr9lR1fZxoUSiU/9YA/nW0N7b7hOUw2X+dL9sbYs2BxUWu1HX9m2Ib9gNVmvVAPZPn5XFX6ebL+2urL5Py0XucbTslVoHVS8ob2Zskn838kmgH9b+Hha4efnhaZuFvQ0RWgKf6ToPKYFw3KDjC/cLqW2mrroYD2dCnuYsvvWmsawXNmolP/I0SMVypc0Zdj1DLrQpBZbemrtycYnHznUc7oUxp73fGU//6ebNfVwtQbLuTnH66iNvUq+9QtZTzZs17nw46FF7ouuEKyK+rKhxkWD1GGl2cASNNaCoPRt4a2/rAdq3lYFXbFv4f0GFn+vIcwz3Uw/GeJHQ/xgsB8WmcdUVRVraG9dHf3NwdbV0bo66qT1EvfRDKG0RJAK5UfeFPofz7gXDV7jaM+0flONpJ6TfKKV9cErPnaIy+1oQ8+7YoVWzHwh6/0o7bU2CQtVsboWZ8uqZy/62rUsWH1Wns3GAFno+96C0DSq8PP1ZJP4t5e9166y7VXqlGUjQ4UxzQw/GeoHv80PTT8abLDR849ZVStdrOqI+fOYOllT+5jHFMqNCFKhfJj+A/0OZM5oX2x8mjfq3WcV2zjEM+pqBp70nmNcZX9be8Ilaqi+Qi/1u2wCc1l6QDZB+j+4y+Jj0NMofJflQ55sw94fZMNz/Qof/ywbiiSb+N4VJxbed5XVaaoihIVNNcNPRYLSvPshhT1MZBvydtHBoXbQRQfrWNVaVlEvYngo5yJIhdwb+zT9T5JXp4M3Nt/TV9XvLqxUfpOqhWHpXq861Y2Otav7nKfan/jSbWzB5Oay8LosfFwvm+9USzYUNs9s2fypXOzxV4DBFg9MP1G433wW8NbDTvibrBJ8F1FuICxuiul+NHih3qUfDDas8Cs+kcwPTIfbURerzg9MddXOcetDWDERpELuFMxl0PkM/69prQ/0zBqjDU+et497dfOX+Yf917POcYczHeAWZ6nyJ/s8NsUlmKj0w0BvjMQ+GI3PZDWhxmI42qEXNlG6QWqarDdpgGzu0oDC28+y1YXQUBaSNsNJhY/XQYtSbFeomCab5kdDFupd+sFgw/2OLDB10to6OjrKztaxqi6FgSkKWIaVTblZtZckSXd0t2BuaqzaW5nNHsWPhzLlM6PWOsdTLZ6Wl8x2mBesUmRm0U2e8Xd3usARrnFKiay8GS1bRbcPHlN6G9QOkg13TV3k+cmyitzPY3NZkLpYNiz2Z+RhmCwcFQ1LA2Q1qOZpKfsGm3dbRxaa2qqkS25CseaVFehvqP6G+snQwgngg40o7NNNJDprMz8oFQ1MtWMOU2VUKX+ElJsgNU+UP6gEJn3Ej4cg0a/raV6qc63mujjcixpqP/+wGzzlPHe7yNGudFKJLl9+FofLVu79T8VZdj9JFtCKuw2RhSmyulOrWxCW1iq8XwONyrLBodzLk+c3o/1kSJHQlD2etwlvNVWtpq21rFIkMHWwZgSmsLAIUuVBBKmVWJoy4lZ+/YeCBlt4f72uPqp2m/UcaV/3qV5kjsT1nnK+u13sGFc4sVRqwLyCkzEFR+EUWT2jXE2WniMbBhwhG/YrevtNFpYmFjm+vqy0QqdFbmtgFaXX0xYqpmlmGGDYQkGpv6F+McIcc0F9daxlFWvrYC2rFD5eRWdtVY+ZIGHpIkiVBxGkVlJ5k+l/EuOeN2uV7l7o+Jufkzfs4jpb+sdCQek/nnSB/7nEsf7lhFItpDdBtoruLlmAaS6rdr6BbBn/qrIVfi0tXzBJZZPIp8m+mMfJVgsWvZ/3eEThbcwi16iL9rIht44WD0xNVNKfWmGJig7H/VQYlOYFp2FFlji00cza84PSgtDURrMoXBn+jEr5xRNBKpS+qV/zwyHMHWf8Old5sskdphrlYE9a3R4LHTovRF3qOJc7vsx+qOfJ9t17W1Znqi9FKtlkvVT1itxqWvDFNu9+jiw4TS+8L1jCazWQBbZmhbe2hbd2hbd5jxuopD+VwlLNNscgowww1M+GL3U4bl5YWlsHa2q/3EVsQ1hGlfJHVgSpUHrSlJF3MfAc6q5n4HrdPVfzXHW1dISXNbPmQoff6SXd/deljvMvJ+So0Qv8LhtWGyEbcpsqC0jTZGUMkkVu1WUhq27hbd7jecGpOZqiRlm+iVBhFSgwwjg/G+Znw+eHpgGGGWy0gsKoXl+dxYbi1orhuJAbEaTKgwhSK4m8KQw4md+fkbbt7tPO7b1T5SKr2d1BnlBrkb3xHveOo13lHIe40RkxvBAqjUmmGlAYln42rPDxML8YMX/T3Wqq6qyNNa1iDe2sob01tbeGdlpqEt8vobyolF+IEaRCyZv6bbYqb84Yc9e8S48Wb+rrMVv7px1dqcois416+MQBLnas3dzv/PilEFY6s83xq5FFgtKCHqbfi2wO1FYzaxQGpKKhqaNWf6oIbQhlpFL+8I4gFUpOmjLqHn75G3W7mLLOrZ6qfa6xvrefB6zniMVO+dxPtvM3e9rc0y6NXxahwsqXb6ixBhrhl8IhuHnDckWH4hqoWxiU2i8UmlbXNrZDCRVdBKnyIIJUBZU3lZ9PZeyTtDndsM6HerrqEaqo5nAvaaPbYqcMMdpmTtdJG++7KTYhDeVenjxDjDHQiPm3Xww30AiDjDK3sJJXddV01qZIUGo/fziuhcbR6xpWVpXyCzuCVPjzpn7Lj4cxZxRr3uubFjO86jRtbOIwz6un5WKnTDHdlrqbYbY+7tI8ykSGcmKuPEOM9ksxYek3o+QVbttcXTWdtLaatlbXzmrazr910DJ6V0NlVCmDVHynhxWXpgy/hUEXULeL/G59vF3nf/q4xUZOtqfbVStmjVoqdYobDDU2QlTIiTnmGmz0/IBUNCwNNlp+4TBcDdV11sZq2trHlguFpVW0UDXKnoZQ6UWQCitmzlj6n8CE12l3thmdzvNslWMN0cue7rCJ05c4fPE/r3jaB552mbVLdavehfXxo9aaaquZ0SYYZqzNrbPYcU96zzQzra6t7W1YZu0LJWuiqQYZaZBRC93/aqQhxsyfs1RLjflh6QDbLBSW2mkeYSmE8IciSIXlN+Fd+h9Dms96rxnbtIMnbWOWyY7xjo62X+KpPxvmbLc73X4OtUOpNXGaGY5zrZt1t0rh0OJdXvawt8DG1vS4ixc77xUf215XrTX1gt5G+F1bzUutnWHF5ckz1Nhiw9Igo0wssl10Q3V11lYnrR1su4WG4tpqpkrONgYKoeLJM8c0o0010jSjTDXSVKPs5MpcNy0nIkiFZVcwl98uYdh/aLwTaz1iQM0vPW9zjXR0rHc11nGJp6dSp7pRO83d6IxSa+YDXjfEGC/o7Uanz39+Va0M96xUqp0WSzy/SmFPWk3VY1Jwji2pV2mQUYYUGYKroopVtNBZGxtZw8G200kbnbTWSRuN1Y//lyEsRZ7Zphq1UDhaNCxNM8oM4xY6r4rq6mkVQSrXkiTpju5yt2ds+CMzf+XHI5j2DZ2ulbb/uw+T67zvYmvZ3wEeUVO9P7zEaz7V07fecF2p7hh/oj3BFR5e7HNL613a11ae9r5Z5uigpTaalUobQ2am2YYYbYgxBhvtN6OW2KvUQF2dC8NRFpRazw9Lq2gZVbxDWIK5Zv1hOJr38UwTFjqvqhrqaa2+1upro4Nt1ddm/sf1Cu9ra1Kpe3Vj1V5YujGP8/PpVG9BlyfMabCul53oB0/bzmW2c+kyfRNt4QzVVdPLLWXSO5DY3m+etKrW4Dx36WJVNVX3rq+c5/AynaNVGU03c35Imnc/7zbEaGNMnH9sVVW003z+EFwnbeYHp+hVCmFhqdRME00zusht1PzHRXuWZhUp+gpV1VwsDBX3ce3lr5pfKb9B40+4sGR5U/mlO2MepeXRrH6HydUme8o2xunvEM9ax8HLdKlfDPeZHz3j8pz9MjzYdjbTBbTQyAEu9qOHK/VfUn/WVDMMWSgcjVnocdGq3dVU1V4Lq2plHavay+Y6aGlVrayqlbaaRcmAUOnNNdM0YxYLRkVvU40y3Rj55ix0bg311NNKPa3V01IL6xQblmppFH+UlKD4qRWKN/lTfjqGuWNY6xFaHWOYTz3lANXUdKKPtdZ1mS/Xy7eqqGJPm5Vem5eimzXmP15NWwMM09evulo9Z20qzwoUzF/dmN1+N9QYQ42dH5YmWNBpXF01HbTUQUsb6Gw/W1lVKx0Kg1IbTWMFXKiUChSYYVyxPUcL9yCNNtvkhc5NVC0MR9mthfV0sot6Wqmv9fzn62q51OkVoXREkAoLK5jLkH8z5Crqb8L6b1JnNd96WA+naGtTh3pevT+YrF2ckcZrrqG6apdSw//YZ36ws78b52W11DTVTGR1giqjVGqiqYYZa+j8oDR2oY9HGDe/UjfUUUt7zXXQSjdrOMi283uTOmiptabRuxcqldmmLVM4mm6stLCQ6zy1NC4SkFprZcPFwlE9rdSO76tyL4JUWGDGz/x0NFO/ZtVLWeX/pFWq+sAlervShk6ylzuLLbK5NM00NMFU083MSZhqp7nzHTF/G5qPfW9L6+pi1TJvS2lLpcabbKTxRhpnhHHze5SyoDTGML+bYdb8c6qpqp3m2muhvRa2su78x+01t4qWMUcpVApzTDfNGNONMd3Y+Y/n3U8tEpjmmr7QufMmZ88LQW1tulgwmtd7VF2tHL3DUNJisnnIKpSP/B+/nkvNdqz9KA02M9csLztRP0/a2XW2ct4K/yIdbJTVHe1ix7jM8SXb/kU86T29fedurzjMDra1gTPs7wPf+MoA+Qr8bLhrnKyFxqXalpI23UwjjJsfkIqGpaLPzTF3/jmJRCtN5oei9lpYRcv5j9troaXGMewWVkqp1GxTFwpDf3S/aDhKVFFHc3W1UE/L+avYFg1H9bSKuUeVdLJ5BKnKbvZoBvyFCa/R+lRWu5GqdU03ztMOMNKXDvDoMk8q/yOXesC/PeIaJ/uHw2JicaF5w2xjTDTGhML7iUYZPz8gjSx8PGWRH/IN1NVWM2001Vbzwvtm2mg2/76VJlEaIKxU5q1YW5ZwNN1YeUV6X6GKaupqoa6W6mn5h/d1NFUl/shYVhGkyoMIUmVo3MtZiFKFNe+n2d5gvF88bk+zTHaEV7S3eYm8XCp1sftd43FrWcVfHegIO2m4Ek6QzJdvvCmLhaOxhfcLbhOMNWmhuUhkxUBbFxOKigalNpqqp06O3mEIJatAvhnGLzUUTSu8LyjS60q2pH9poWjefS2NYt5R6YggVR5EkCoDedMYeDaj76fpvqx5LzWyyeOjfONRu6mtiaO8rolOJf7yXxngCo941acSbKaLnXXTzRrW1kFHrcpNb1W+fNPMNNl0E0013pTC2+Qijxf/eKKp0kW+HOuqpaUmWmpc5NZEC43mP573fAN1K/sQQajgUqk5ppturOnGmuH3wse/z38uu2VBaYbfpYWV6uepru5SQtGCXqWaGsT3TO5Vyv8BEaQqm3llDeaMZrX/0vokkuxrf4iPPGEvTa3hKG+oW8pVvYcb61Wfes/X3vfN/KX0NVS3urbaaFYYMrKwUU9ttdVUW011Cu9rq6mKRFoYW9Ii/xVIzTHXLHPMMsfMwvvs8WyzzDHVDFPMMMV0k003xfSFPp5WuLpvUVVV0UQDTeffGhZ5nH3cXMOFwlGuViyGUFLmmrlQGJqxUCj6fbHQtOiQGtTWRF0t1NF8fiiaN/9o0aBUQ90cvMvwJ0SQKg8iSJWSomUNGmzKWo9SZ7X5n/7Fm552oHY2c7iX1dKgTJuXSo00zk+G6G+oAYYZXTgkNtZEY00yzczFhsCWV1VV1FZTLTXUUkN9dTRQRwN1NVRXA3U1UGeRx/XUV1tj9TXVQDMNNVA3hgZChZdnTrFhaEnPzTFtsWvU1LCwZ6h54X2L+ZOzF32ujqaqVtKSI5VEBKkSuWCSXIgDsRZm4hNckKbpgGU8P4JUSZven/7HMPWbwrIGF1FlwdDZQG950r46280hnla9HPec5Ms302wzC3uVZpqtQCqRrU5b8F/2cQ3V1VJD7cLgVF6GDEMoDfnyzDBuiUNpiz63aPFHsurYdYoEoHlhqPjnmqlWivtmhgonglSJXDBJ3sRT+EJWp+oqrIcuaZpO/6NzC8+PIFVS0gJG3M6gC6i5SmFZg00XOmSw3h6zu452dJgXVqhGVGn73pPmmKaJ1XW0fa6bE0KZyJbtTzHd72YYV+S24ONFPzeryN6F81RTq5ieoiUHpPL8h1Qo9yJIlcoLJElzjMV2aZr2XobjI0iVhFlD6X8Ck96n7V/pdA1VF17hNdYP7rOFtjZxpNdKtUDcbNO85Di7uVkjq4ARvtDX41rb0FAf28aFGuu40Hn9vaKtTdTX2o9e0M5mGmhbau0MobTMNWuxIFRcGCr6+YJihrJraqiOZupqro5mRW7Z0NmioalGLFwIZadSfqGVxThHw8L7CcV9MkmSmizUN1y/1Fu0MktTxjzGL2dSrQHrv0OTnRc7bIbxnrSvRlZ1uJdLNUR97QGTDfGTF+zqRpBntmcc4mR91NNSU2t4yQlO0HOx85PCuUjZEEKl/D4N5UyBfDNNWKYwNO/zixZ6JFuyXzQQZRvNrjv/40XDUm1Ny2WvcQiVWakGqSRJEtyEj9I07beEwy7EZaXZjkpjzu/8fBrjXqDl0ax2G9UbFXtoD6eaZbJjvVvqG11u5ETQyxXznxuit5oaqKclaGtTw31mmjHzn4O17Kufp+WZpaEOGmhTqm0NlU+eOWYab4bxZhpfGJDGL/TcDOMXCkYzTbB4B3iijqZFeoeaaWXDYgLRgsfRWxRCxVfaPVK3Y31s/QfHXCMLW/PUx/DSbNRKaVwPBpxMmkeXZ2mx5Erk/b3sJ8872FOLDaWVlUkGq63J/I+rqKqm+sb6YaEgBes6rKybFyqgAgVmm7xICJqwWCBa9HPFrUSDWhqprak6mqqtqcY6aWezYgNR1lvUOCpgh1AJlVqQSpLkNuyLbdM0XWIwStN0NmYXOa+0mrRyypvCwHMLi2vuzRr3UrPVEg9PpXq7Ukc7WsehZdjQhc0wTrVFhhOrqWWWSblpUChX5pq5WPhZUk/RvEA004TFCjqSfV0VDUR1NNVY58KPmyz2udqaRigKISyzEg9ShcN5t+EAbJ+m6W8l/Rqh0KTe9D+OueOyAFWkuOaSTDbUSF861HM5HVKoqaFFh0bmmKZOKRcBDWWnQL5ZJptlopkmmmlCsY9nFX487zbDeHnFFEJNVFFL44VCT1NrFAk/TRYLRHU0jVVoIYRSVRo9UnfgSOyHqUmSzOsemZymafFlosPyyZ/Fbxcz/CYabsUG71F72bZymSjLtc2tU5otXKpm1vK1e+d/nGe22aZqpEMOWxUWlS3Bn7oMIWjxx7NMVtxC2nmBqLbGamuilsbqaaW5LkVC0OI9RTU1jCKoIYRypzSC1OmF9z0Xef4EPFQKr1e5TPuen45kxs90uo7255Is+xDEvKAy3GeaW6u0WrlUHWxrurEmG66hdgbrpa1NIkiVgjxzzDa5sHdoUuHjSQv1Fi0eiOYFpklS+cVet5ZGiwWixjrNf1xb48L7Jgs9rql+TLAOIaw0SjxIpWkaPyFLQ5oy4g5+/Qd1Vqfbl9Rbb7kv01hHXRzsDX9VX2ur2a0UGruw7z1piKyE2Lv+qYNtbeoM+3vYh67SzuYG6+VAj5d6WyqabOPXaUsIQYt/XNxzxQ2TzVNDvYWCT21NNNBuiSFo3uNaGsYcohBCEHvtVQxzfmfAiYx/lbZn0uk/VF3xeR+zTfOMg/zqbZ3tZmOnWd0esdVDCSpQYK7pZptitqnmmDr/8WxTCj8u+niKWSYXG5SKm0ANVVRTSyM1NSzsHWq4XB/X1EDV2DInhFByKmVHSgSp8m7Cu/Q/lnQuaz5Is71L5LIFCvzgGR+7zmjfqq6udja3iq210U1Ta2qsY6XYYDSVyjfHXDPMMd1c04s8nrFQIFpSCFo0LGVL6pf85VhFdTU1UFN9NTVQQ/3lDkLV1Y4hshBCeVIpfyBFkCqv0gIGX8GQf9F4Z9Z6hJqtS+WlxvrRz141zMfmTulpfI0pptSkwZyq2sxqaVrDjupppb7W6mmlzZiBauQn0tod5TXeQk311VBPdXVUUU0V1QvvF9yqqi5RRSqVKliGW758c+SZLd/sFbiftVAQWtLjeYFpSfOAiqqq5kLBZ+HHDQr/HeovFpCKOy56/0IIK6EIUuVBBCnMnchPRzPhDVa9gg4XkZTNaqX0p+MlYx4G0+uv4tu1dzW2zhxTjTLNaK3GDTGk/kyTa8619u8Mb8DUcpAJqqiuqhqqqamqmqqrrbq6qqujRuF9dXWX+HjBccWfU0P92JojhLDySVMKZlMwM7vlz1jwuGAm+TMpmFHk8R8c1+WxShmkYoJEeTNjAH33Im8C671O093L9OWTWquyxXDSVN1a7Wy12BGvkGwmT2N5VZ43Syez1THbVHlmKpCnQJ58c+c/Lpj/OF+iylJvVVSVqKKqGqqqOT8cLem+qhqxLD6EUPGlKQWzCm/zQsqswgCzyHMLfW6R5+aHnnnPzVg4BC0UhGZZrj6IKrWoUpsqdbK5ulWK3CqpCFLlyZTP6bsnNVpkq/KWsTZUiavZdsmfa7YvY59WrWCWajU7qFVzs7JrVwghlIU0JZ2zSBgp6XBTzPXS2Utv20KqFIaYWkVCTa0i4abw+WoNiw8+VesUeVx7yQFp/nE1y2x0pCKJIFVeTPmcb3ek3gas14PqTZZ+TmkomM6oB7NvmInv0v486q698DEtYu+7EEIZSNNs/9DlCSMlEm6Ws5eGJYeYos9Vb/oHwadWkUBTa9mul1Rb6m4WofRFkCoPZg3h+72otz4bvE3VurlrS/ODaVDYy1S9Bf0OYNMf46+QECqreb0z84ecSvC2LOFmCeU/liipWXzoKBpOqjWkSqulB5Wi5/xR8ElqRKCpxMpNkEqSpDu6Uwknu/zcPfuGXK9HbkMU1Ou24HHt1Zg5gGl9qd81Z00KoVKbPxm4FIJMcbe0mNdabsmCoaAqtYq51SwMIf/f3v3H1lXedxx/f32T2PllA2VACiFLxkiWqUACGWOhjTOVjq1qKHQMWrVbq0kdELbCoGiZJmX7Y426FkTLMiEBHQh1qGpVkrVdNFpIWkoNCwk/VlDS4gQINMDiNLZJiBM73/3xnGsfn9zre8+x7w/7fF7S0b3nuc+557lfP/b9+nnOj9kw/fTxJTGjytr0T5/UXdMkUu6+EdgYO2svH46+AQd/AEseCsO+jdT7NLzwYVh5AAptMNQfylt0tprkmA+VSGTqmdgcy9DoltHJRbmlMGckkSmb9GRYNOUkOdI0iVRuDUVXeJjW0dh2ALSeA+feEZIogN6noP0PYPbSxrZL8mV4KmkglrAMRCMlA6PLs5Z57H0rlflg+s9g00uMwJRIOKa1Q8sZGRKVSkmP/rSL1IuuI9Vo7vDcB+HIy7BoA5z5mXCGRKP8eiv07wCGohsjbwhnEcrUdmIwXXJxIlGeuqxcgpR1BCZiM2JJS2uUcLROQFlrdSM8w9vqPoSSS7kchlQi1QyOH4Rf/jW880j4Y91xOZyyGuYuh1lLwkiR5v0np+JZR8MjLMei58dGykatlynzMttmfr9EedoDeocVJihRmYgyHfAr0mC5/AVUItVM3tsL//dtOLQVDj0JJw6H8pZZ0PabMOOsaJkXRommdUChHQpzwxRB/LH4xdIyY3L+d+wejk3x49EyGB5PJNaLz4fLk3Xj9ZJ1Y69VVbdMMlIpmRm3QvRzLI62zIj9bGckftbRY6mycnWH3zdL8jIJ+5aI1IoSqWaQ60QqzofCZRGO7IYju+A7j0HfAXj/cbjgKBx/Gwb7qG4koaXEF2f0aAVC37fov3krsU75OpwI9wVkKDz6UFQWexx+faj6ujXTEh0IOx1apo88H16mReWJ14brttYmmSlZtzXar5IVEZkUcplI6YjEtPr74f77YedOeOkleOIJOOWU7O+3fTt885uwbBk89RSsWwcLF4Yvz5mLwvLUcbj+GzBvHnz3u3DOpXD22dFp0UdCQjXUHz1Gz08aJSnz6IPhfYgvnFw2nHAnyqwl+qJPPFoLUKjwejXbTC+d2IxKeGKvVayrKVIREZk4SqTS6u6GJUvg1luhp2d8SdTAAFx7LTzzDJx5Jpx/Pnzuc7Bt28l1W6IEoLV15DgQs3AdlsJsYF72doiIiEgmmtpL69gxWLkSNm0Ko0Lj8cMfwm23wYsvhvWhIZg9G157LSRWcd/6Fhw9CgsWQGfn+PYrIiIy8TS1J1XYsydMva1ZE0aSpkUhPHgQvvKV2BRYCTNnwvr1I+uvvgqnxe6pVyjA3LlhyjCZSF2n+9uJiIg0GyVSaTz8MOzdGx4XLYJHHw1TcxASog0b0r3fgQPQ1ja6rK0NDh2akOaKiIhIbSmRqtbevXDPPdDVFUaOVq+G118f33t2dJw8gvXuu3D66eN7XxEREakLJVLV2rwZrrkmJFEQRo0uuGDk9Z4e+OpX003tLVkC9903sj4wEM4KXLBgQpsuIiIitZG/ROrIEdi1K/12hw6F0aKdO8PxUN3dcOqpYb2oOM03lnj9OXPgzTdhy5ZwTFRXFyxdGpKynp70bRQREWmU5csb3YKGaJqz9sxsLbAWaAEWU6uz9nbuhIsvTt9AERERKc89l2ftNU0iVVTzyx9kHZESERGR8pYvVyLVDJr+OlIiIiJSSi4TKd0vQ0RERCQjJVIiIiIiGSmREhEREclIiZSIiIhIRkqkRERERDJSIiUiIiKSkRIpERERkYyUSImIiIhkpERKREREJCMlUiIiIiIZKZESERERyUiJlIiIiEhGSqREREREMmqaRMrM1prZy8D/NLotIiIiItVomkTK3Te6+1Lg9xrdFhEREakvM/tTM3vBzPrNzM3sPTObVcV215jZsWibfWb2vJmdV482QxMlUiIiIpJf7v4dd78Q2ABsB9qARWNtY2bvA64GpgPfd/f57n6Ru79S8wZHlEiJiIhIM1kBPBQ9/60KdW8HdkbP/6tmLRqDEikRERFpCmZWIIxEPR8VlZ2iM7M/Ap4GlkdFW2vauDKUSImIiEizuJgwwtQdrZcckTKz2cCH3X0z0Ansd/dddWlhwrRG7FRERESkhE5gq7u/ZWaHKT+1dztwZ3RQ+TnAfyQrmNlcYD1wPvAWcBDYA1zj7ldOVIM1IiUiIiLNYiXws+j5HkpM7ZnZpcCr7v4WsDoq3pqocxrwE+As4Cp3/zzwC+BOYN5ENliJlIiIiDScmU0Dprn7kaioGzg3Ki/WmQF8yt2LB6MXE6knEm/3baAD+Ct396jse8Ac4Mcl9m1m9mkz+7e07dbUnoiIiDSDFcCzsfVuQp6ygJFjpr4AfD1WpxN43d33FAvM7HrgD4EvuvvhWN0Lo8dt8Z2a2SeBS4APAS+lbbRGpERERKQZdDJ6iq6YPJ0HYGZLgQF3747WFxOm6ZJn690QPW5KlK8CnMSIlLs/4u63kSGJAiVSIiIi0hwuI1zOoGj4zD0zawFuBDbGXj9pWi+aBrwc2FfiopydwP+6e89ENrpmiZSZ3WRme83sqJntMLMP1mpfIiIiMnmZ2XTA3P1orDh+CYQbgAfcfSj2eqkDzd8HFIDnEu8/k3ALum3R+rVmNmci2l6TRMrMrgPuBv4ZWAY8CWwxs3NrsT8RERGZ1C4l3BYm7jVgkDAld6q7P594fRXQ7e77YmXvAIeB5KjTdcAM4KfR+lXu/u4EtLtmB5v/LSFzvD9avyW6AumNwLp4RTNrBVpjRXOz7NDMrLe3N8umIiIiMk4dHR3tQH/sLLmqRBfXXEfiFi/uPmhmrwHthMsWxLdZDpwJbEls42Z2P9BpZhatXwH8WVTlV9GxVXuYIJby81Z+w3Bq4hHgWnd/NFb+NeAid1+VqP+PhAtmJXW4e1+K/bYDyqREREQap+rvbjNrI1wzajEwCzhBGIW60t1/EdV5BLjX3X8crf8L8FFgITCTkG/sAR509zujOrOAe4HTgP3AbuAu4J8Ix2G9Ddzs7r9OtOdBAHf/bJoPXItE6v3Am8BKd/9ZrPzvgb9w98WJ+qVGpN4gfSJlvb29JyrV6+vrY/78+ezbt4/29vZq337YihUr2L49OfpY+20buW/FLD3FLD3FLD3FLD3FLL1qY9bR0dFBhhGpZpE1karldaSSgbQSZbj7ADAwXMks285S/uDa29sz/RIVCoVM241320bvGxSzLBSz9BSz9BSz9BSz9CrFLM3gx1RSi0TqADBEuCx73BmE4bRJbe3atQ3ZttH7Hg/FrL77Vszqu22j9z0eill9953XmDU7M/s4sAa4KqzaN4D/dPdNVW1fixE4M3sG2OHuN8XKXgY2u/u68luOOtYp1dRepOKH6evro6Ojg97e3nFn53mhmKWnmKWnmKWnmKWnmKWXImbZppQmuVpN7d0FPGxmzwJdwOeBcwkHfzVUa2sr69evp7W1tXJlARSzLBSz9BSz9BSz9BSz9BSzsdVkRArCBTmBOwiXb/85cKu7/6SK7cYzIiUiIiJSNzVLpLKycLT5XCbxkf8iIiKSD02XSImIiIhMFrppsYiIiEhGSqREREREMlIiJSIiIpLRlEqkzOxDZvY9M/uVmXl0ka1K26wysx1mdtTM9pjZDXVoatNIGzMz64zqJZcldWpyw5nZOjPbbmb9ZvaOmW2KboJZabtc9rUs8VI/AzO70cxeNLO+aOkysz+usE0u+1hR2pipn40W/a66md1doV6u+1nSlEqkgNnAC8DN1VQ2s4WEu00/CSwDvgR83cw+UbMWNp9UMYtZTLi0RXH55QS3q5mtAjYCvw9cQbge22PRHcxLynlfSx2vmDz3szeAvwMuiZYngM1m9rulKue8jxWlillMnvsZAGa2gnDNxxcr1FM/S5iyZ+2ZmQNXj3WJdzP7MrDG3X8nVnYvcKG7X1b7VjaXKmPWCWwFTnX3Q3VpWJMzs98A3gFWlbtWmvraiCrj1Yn62UnM7CDwRXd/oMRr6mMlVIhZJ+pnmNkcYCdwE/APwPPufkuZuupnCVNtRCqty4DHEmX/DVxiZtMb0J7J5Dkz229mj5vZ6kY3psE6oseDY9RRXxtRTbyK1M8AMyuY2fWEEeSuMtXUx2KqjFlR3vvZRuAH7v6jKuqqnyXU6hYxk8VZnHwj5bcJcTkd2F/3FjW//YTh3x1AK/AZ4HEz66zmyvVTjZkZ4ZZIP3X3n49RVX2NVPFSPwPM7AOEJKANeJcwYvxymerqY6SOWe77WZRsLgdWVLmJ+llC3hMpOPlGx1amXAB33w3sjhV1mdl84HYgF394Ev4VuAC4vIq66mtVxkv9bNhu4CLgFOATwENmtmqMxEB9LEXM8t7Pos/6NeAj7n40xabqZzF5n9p7i5Bdx50BDAI99W/OpPU08NuNbkS9mdk9wBpgtbu/UaF67vtayniVkrt+5u7H3P0Vd3/W3dcRTgz5Qpnque9jkDpmpeSpn11M6CM7zGzQzAYJJ4f8TbReKLGN+llC3kekuoCPJco+Ajzr7scb0J7Jahk5Gs6NpqfuAa4GOt19bxWb5bavZYxXKbnqZ2UYYQqqlNz2sQrGilkpeepnjwMfSJT9O7AL+LK7D5XYRv0sYUolUtGZB+fFihaa2UXAQXd/3cw2AGe7+59Hr98L3GxmdwH3EQ6i+0vgk3VsdkOljZmZ3QK8CrwEzAA+TRg+z9OprxuBTwFXAf1mVvzvrNfd3wNQXxsldbzUz8DMvgRsAfYRbuR+PdAJXBm9rj6WkDZmee9n7t4PjDpW0cwOAz3FYxjVz6rg7lNmIfzCeInlwej1B4FtiW1WEU77HAD2Ajc0+nM0c8yAO4BXgPcIZ109CfxJoz9HnWNWKl4OfDZWR31tHPFSP3OABwhf8gOEy0X8CLhCfWziYqZ+VjKG24C7y8UsKst1P0suU/Y6UiIiIiK1lveDzUVEREQyUyIlIiIikpESKREREZGMlEiJiIiIZKRESkRERCQjJVIiIiIiGSmREhEREclIiZSIiIhIRkqkRERERDJSIiUiIiKSkRIpERERkYz+H8rFFHBe2xp8AAAAAElFTkSuQmCC", "text/plain": [ "Graphics object consisting of 24 graphics primitives" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#Shock Wave 2'\n", "from sage.plot.colors import rainbow\n", "clr=rainbow(8);\n", "delta,Ma1,gamma,beta,pp=SR.var('delta Ma1 gamma beta pp')\n", "eq1=tan(delta)==(Ma1^2*sin(beta)^2-1)/(tan(beta)*(Ma1^2*((gamma+1)/2-sin(beta)^2)+1))\n", "eq1=eq1.subs(gamma=1.4)\n", "eq2=pp==2*gamma/(gamma+1)*Ma1^2*sin(beta)^2-(gamma-1)/(gamma+1)\n", "eq2=eq2.subs(gamma=1.4)\n", "# sol2sb=solve(eq1,sin(beta),solution_dict=True)#sin(beta)\n", "# eq22=sin(beta)==sol2sb[0][sin(beta)]\n", "# eq2.subs(sin(beta)=sol2sb)\n", "sol2=solve([eq1,eq2],[pp,Ma1],solution_dict=True)\n", "# sol2=solve([eq2.subs(eq22)],[Ma1,pp])#poor ability to solve equations. only use substitution can help.\n", "# sol2=solve(sol2,Ma1)\n", "sud=0;bb=0.05\n", "fg2=parametric_plot((sol2[1][Ma1].subs(delta=sud),sol2[1][pp].subs(delta=sud)),(beta,0.1,pi/2-bb),color=clr[0],\n", " axes_labels=(r'$Ma_1$',r'$\\frac{p2}{p1}$'),title=r'$\\frac{p2}{p1}\\ vs\\ Ma_1$')\n", "fg2+=parametric_plot((sol2[0][Ma1].subs(delta=sud),sol2[0][pp].subs(delta=sud)),(beta,0.1,pi/2-bb),color=clr[0])\n", "fg2+=text(r'$\\delta=0^\\circ$',(1.2,1.1),color=clr[0])\n", "for i in range(1,8):\n", " sud=i*pi/36\n", " fg2+=parametric_plot((sol2[1][Ma1].subs(delta=sud),sol2[1][pp].subs(delta=sud)),(beta,0.1,pi/2-bb),color=clr[i])\n", " fg2+=parametric_plot((sol2[0][Ma1].subs(delta=sud),sol2[0][pp].subs(delta=sud)),(beta,0.1,pi/2-bb),color=clr[i])\n", " tx=sol2[1][Ma1].subs(delta=sud,beta=pi/3)+0.2;\n", " ty=sol2[1][pp].subs(delta=sud,beta=pi/3)+0.2;\n", " delt='$'+str(i*5)+'^\\\\circ$'\n", " fg2+=text(r'%s'%delt,(tx,ty),color=clr[i])\n", "end\n", "fg2.show(xmin=1,xmax=4,ymin=0,ymax=10,aspect_ratio=0.2)\n", "fg2.save_image('2_2d.png',xmin=1,xmax=4,ymin=0,ymax=10,aspect_ratio=0.2)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "vscode": { "languageId": "python" } }, "outputs": [ { "ename": "AttributeError", "evalue": "'Sequence_generic' object has no attribute 'subs'", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mAttributeError\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 7\u001b[0m \u001b[0msol2\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0msolve\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0meq2\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mbeta\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 8\u001b[0m \u001b[0;31m# plot3d(pp,(Ma1,1,4),(beta,pi/15,pi/2),xlabel='Ma1',ylabel='beta',zlabel='p2/p1')\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 9\u001b[0;31m \u001b[0mfg1\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mplot\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0msol2\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0msubs\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mpp\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mInteger\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m5\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mMa1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mInteger\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mInteger\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m4\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0maxes_labels\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'Ma1'\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m'p2/p1'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mcolor\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mclr\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mInteger\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mtitle\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34mr'$\\frac{p_2}{p_1}$ vs $Ma_1$'\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 10\u001b[0m \u001b[0;31m# for i in range(1,8):\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 11\u001b[0m \u001b[0;31m# print(pp.subs(beta=i*pi/36))\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/usr/local/SageMath/local/lib/python3.9/site-packages/sage/structure/sequence.py\u001b[0m in \u001b[0;36m__getattr__\u001b[0;34m(self, name)\u001b[0m\n\u001b[1;32m 882\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m__hash\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 883\u001b[0m \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 884\u001b[0;31m \u001b[0;32mraise\u001b[0m \u001b[0mAttributeError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"'Sequence_generic' object has no attribute '%s'\"\u001b[0m\u001b[0;34m%\u001b[0m\u001b[0mname\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 885\u001b[0m \u001b[0mseq\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mSequence\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 886\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;31mAttributeError\u001b[0m: 'Sequence_generic' object has no attribute 'subs'" ] } ], "source": [ "#Shock wave 2 p/p\n", "var('pp gamma Ma1 beta')\n", "pp=2*gamma/(gamma+1)*Ma1^2*sin(beta)^2-(gamma-1)/(gamma+1)\n", "pp=pp.subs(gamma=1.4)\n", "# plot3d(pp,(Ma1,1,4),(beta,pi/15,pi/2),xlabel='Ma1',ylabel='beta',zlabel='p2/p1')\n", "# fg1=plot(pp.subs(beta=0),(Ma1,1,4),axes_labels=['Ma1','p2/p1'],color=clr[0],title=r'$\\frac{p_2}{p_1}$ vs $Ma_1$')\n", "# for i in range(1,8):\n", "# print(pp.subs(beta=i*pi/36))\n", "# fg1+=plot(pp.subs(beta=i*pi/36),(Ma1,1,4),color=clr[i])\n", "# end\n", "# fg1\n", "# f=s.function(x)" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "vscode": { "languageId": "python" } }, "outputs": [ { "ename": "NameError", "evalue": "name 'clr' is not defined", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mNameError\u001b[0m Traceback (most recent call last)", "\u001b[1;32m/media/kai/Data/Program Files/SHU/TEX/Aerodyn.ipynb Cell 6\u001b[0m in \u001b[0;36m1\n\u001b[1;32m 8\u001b[0m sol32\u001b[39m=\u001b[39msolve(eq3\u001b[39m.\u001b[39msubs(beta\u001b[39m=\u001b[39mpi\u001b[39m/\u001b[39mInteger(\u001b[39m2\u001b[39m)),[Ma2],solution_dict\u001b[39m=\u001b[39m\u001b[39mTrue\u001b[39;00m)\u001b[39m#another solution\u001b[39;00m\n\u001b[1;32m 9\u001b[0m sud\u001b[39m=\u001b[39mInteger(\u001b[39m0\u001b[39m);bb\u001b[39m=\u001b[39mRealNumber(\u001b[39m'\u001b[39m\u001b[39m0.05\u001b[39m\u001b[39m'\u001b[39m)\n\u001b[0;32m---> 10\u001b[0m fg3\u001b[39m=\u001b[39mparametric_plot((sol3[Integer(\u001b[39m3\u001b[39m)][Ma1]\u001b[39m.\u001b[39msubs(delta\u001b[39m=\u001b[39msud),sol3[Integer(\u001b[39m3\u001b[39m)][Ma2]\u001b[39m.\u001b[39msubs(delta\u001b[39m=\u001b[39msud)),(beta,RealNumber(\u001b[39m'\u001b[39m\u001b[39m0.1\u001b[39m\u001b[39m'\u001b[39m),pi\u001b[39m/\u001b[39mInteger(\u001b[39m3\u001b[39m)\u001b[39m-\u001b[39mbb),color\u001b[39m=\u001b[39mclr[Integer(\u001b[39m0\u001b[39m)],\n\u001b[1;32m 11\u001b[0m axes_labels\u001b[39m=\u001b[39m(\u001b[39mr\u001b[39m\u001b[39m'\u001b[39m\u001b[39m$Ma_1$\u001b[39m\u001b[39m'\u001b[39m,\u001b[39mr\u001b[39m\u001b[39m'\u001b[39m\u001b[39m$Ma_2$\u001b[39m\u001b[39m'\u001b[39m),title\u001b[39m=\u001b[39m\u001b[39mr\u001b[39m\u001b[39m'\u001b[39m\u001b[39m$Ma_2\u001b[39m\u001b[39m\\\u001b[39m\u001b[39m vs\u001b[39m\u001b[39m\\\u001b[39m\u001b[39m Ma_1$\u001b[39m\u001b[39m'\u001b[39m)\n\u001b[1;32m 12\u001b[0m fg3\u001b[39m+\u001b[39m\u001b[39m=\u001b[39mplot(sol32[Integer(\u001b[39m1\u001b[39m)][Ma2],(Ma1,Integer(\u001b[39m1\u001b[39m),Integer(\u001b[39m4\u001b[39m)),color\u001b[39m=\u001b[39mclr[Integer(\u001b[39m0\u001b[39m)])\n\u001b[1;32m 13\u001b[0m fg3\u001b[39m+\u001b[39m\u001b[39m=\u001b[39mtext(\u001b[39mr\u001b[39m\u001b[39m'\u001b[39m\u001b[39m$\u001b[39m\u001b[39m\\\u001b[39m\u001b[39mdelta=0^\u001b[39m\u001b[39m\\\u001b[39m\u001b[39mcirc$\u001b[39m\u001b[39m'\u001b[39m,(RealNumber(\u001b[39m'\u001b[39m\u001b[39m1.15\u001b[39m\u001b[39m'\u001b[39m),RealNumber(\u001b[39m'\u001b[39m\u001b[39m1.35\u001b[39m\u001b[39m'\u001b[39m)),color\u001b[39m=\u001b[39mclr[Integer(\u001b[39m0\u001b[39m)])\n", "\u001b[0;31mNameError\u001b[0m: name 'clr' is not defined" ] } ], "source": [ "from sage.plot.colors import rainbow\n", "clr=rainbow(8);\n", "#Shock Wave 3\n", "delta,Ma1,Ma2,gamma,beta,pp=SR.var('delta Ma1 Ma2 gamma beta pp')\n", "eq1=tan(delta)==(Ma1^2*sin(beta)^2-1)/(tan(beta)*(Ma1^2*((gamma+1)/2-sin(beta)^2)+1))\n", "eq1=eq1.subs(gamma=1.4)\n", "eq3=Ma2^2==(Ma1^2+2/(gamma-1))/(2*gamma/(gamma-1)*Ma1^2*sin(beta)^2-1)+2/(gamma-1)*Ma1^2*cos(beta)^2/(Ma1^2*sin(beta)^2+2/(gamma-1))\n", "eq3=eq3.subs(gamma=1.4)\n", "sol3=solve([eq1,eq3],[Ma1,Ma2],solution_dict=True)#this gives 4 sets of solution\n", "sol32=solve(eq3.subs(beta=pi/2),[Ma2],solution_dict=True)#another solution\n", "sud=0;bb=0.05\n", "fg3=parametric_plot((sol3[3][Ma1].subs(delta=sud),sol3[3][Ma2].subs(delta=sud)),(beta,0.1,pi/3-bb),color=clr[0],\n", " axes_labels=(r'$Ma_1$',r'$Ma_2$'),title=r'$Ma_2\\ vs\\ Ma_1$')\n", "fg3+=plot(sol32[1][Ma2],(Ma1,1,4),color=clr[0])\n", "fg3+=text(r'$\\delta=0^\\circ$',(1.15,1.35),color=clr[0])\n", "for i in range(1,8):\n", " sud=i*pi/36\n", " fg3+=parametric_plot((sol3[3][Ma1].subs(delta=sud),sol3[3][Ma2].subs(delta=sud)),(beta,0.1,pi/3-bb),color=clr[i])\n", " fg3+=parametric_plot((sol3[1][Ma1].subs(delta=sud),sol3[1][Ma2].subs(delta=sud)),(beta,0.1,pi/3-bb),color=clr[i])\n", " tx=sol3[1][Ma1].subs(delta=sud,beta=pi/2.4);\n", " ty=-sol3[1][Ma2].subs(delta=sud,beta=pi/2.4)+0.2;\n", " delt='$'+str(i*5)+'^\\\\circ$'\n", " fg3+=text(r'%s'%delt,(tx,ty),color=clr[i])\n", "end\n", "fg3.show(xmin=1,xmax=4,ymin=0,ymax=4,aspect_ratio=0.5)\n", "fg3.save_image('3_2d.png',xmin=1,xmax=4,ymin=0,ymax=4,aspect_ratio=0.5)\n", "fg3.save_image('3_2d.eps',xmin=1,xmax=4,ymin=0,ymax=4,aspect_ratio=0.5)\n", "# implicit_plot3d(eq3,(Ma1,1,4),(Ma2,0.3,4),(beta,0,pi/3.5))" ] }, { "cell_type": "code", "execution_count": 26, "metadata": { "vscode": { "languageId": "python" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 35 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.20987392654089496)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 35 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.20975028620814265)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 68 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.31258607163312907)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 68 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.3132840686322348)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 102 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.4217456296240393)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 102 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.4221174568757481)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 136 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.5325893459440663)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 137 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.532942461154469)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 172 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.6454503562143244)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 172 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.6470570860352423)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 200 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.7353981633974482)'\n", "verbose 0 (3838: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 200 points.\n", "verbose 0 (3838: plot.py, generate_plot_points) Last error message: 'Unable to compute f(0.7353981633974482)'\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAnYAAAGeCAYAAAAZq0yLAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAABhgElEQVR4nO3dd3hUZd7G8e8kIQkJJPQQWuhNVCCAgCAdBASxIKjLYgHltZe1RF3ruiiyCmqodgVEBQQFlCCEDgLSBKR3ApGWBAKp5/3jSSUJZCDJmUzuz3Wda2ZOzsz8Mjsud57qsCzLQkRERESKPQ+7CxARERGRgqFgJyIiIuImFOxERERE3ISCnYiIiIibULATERERcRNedheQH6mpqRw/fpyUlJQcP6tatSpeXsXi1xAREREpVC7dYnfixAmGDx9O2bJlqVatGjVr1sxxfPDBB3aXKSIiIuISHK66jt2uXbvo0qUL/v7+PPLII1SpUoUvv/ySX3/9lbfeeovmzZvj6elJp06d8PPzs7tcEREREdu5ZLA7d+4cLVu2pFKlSixcuJDSpUsDcOHCBUJCQrjzzjsJDw+3uUoRERER1+KSg9NGjx7N3r17+fXXXzNCHYCvry916tRh7969NlYnIiIi4ppcboxdamoqEyZMoH///tSuXTvHz6OioggMDCz6wkRERERcnMsFux07dnDs2DG6deuW42f79u3j4MGDdOrUyYbKRERERFybywW72NhYAOrUqZPjZ5MmTaJixYr84x//KOqyRERERFyey42xa9y4Mb6+vhkBL92OHTsIDw/n008/pWzZsjZVJyIiIuK6XK7FLjAwkBdffJEPP/yQ+Ph4AP744w9uv/12Ro8ezcCBA22uUETcXWxsLG+99RYVK1bkrrvuyvWaiRMn4uXlxQsvvMDq1auLuEIRkdy55HInAN988w3z5s3Dx8eH8uXLM2LECBo2bGh3WSJSQpw4cYLnn3+e9evXs2nTpmw/2759OxMnTmTjxo1ERkbaU6CISC5cNtiJiNjpu+++o1atWnTr1o1z585lnE9OTuaHH35g/vz51K9fn3//+982Vikikp3LdcWKiGQ1efJkWrZsSXBwMMHBwTRu3Jjbb7892zUpKSmEh4czZMgQ1q5dC8C9997LV199lXHN/PnzmTRpEpMmTeLjjz++7Pvu27ePtm3b4uXlxdGjRzPOz5w5kwEDBvDbb7/RvXt3AJKSknjllVeYMGECY8eOpWfPnpw6daogfn0REaco2ImIy3r++ec5evQoq1ev5s8//yQgIIBNmzYxc+bMbNfNnj2bQYMGER8fz/79+wHo169fxiSsuLg4wsLCeOihh2jfvj1z58697Hs7HA4A6tWrx65duwD4888/qV+/PgcOHCAuLo7WrVsDMHz4cIKDgxkxYgS9e/dm48aNlC9fvqA+BhGRfFOwExGXtH79ev744w9ee+01vL29qVixIn5+fjlmzAN07doVDw8Pli1bxi233AKYQNamTRsAvLy8OH36NC1btuTbb7/lhx9+uOR7HzhwIGPJpfr167Nr1y6SkpLYsmULLVu2ZOHChXTq1AkvLy82btzIzJkzGTZsGACbN2+ma9euGcHw3LlzzJ07lx49ehTYZyMikhcFOxFxSQsXLswIaWDClr+/P5UrV85xbbly5Zg7dy6dO3fO2IZwy5YttGrVCoDSpUuzbds2wsLC+Omnn/jmm28u+d6LFy+mS5cuADRo0IDdu3fzww8/ZHQBZ+2GXbRoER06dMDHxyfjcbdu3Thz5gwA/v7+9O3bl6SkpKv4NERE8kfBTkRcUvPmzTNa586fP8+///1vJk+enOf1x48fp1atWgCcPn2aMmXK4OHhwYEDB6hcuTKlSpVi4MCB/OMf/6B69eqXfO/o6GgqVaoEmBa7+fPn07hxY3x8fEhNTSUyMjJjd5xy5coRFBSU8b4///wzN910E9OmTbvqz0BExFkut0CxiAhAr169iI+P54svvuDs2bOMHj2aKlWq5Hn94MGDeeaZZ/j6669JTk5m6NChAFSsWJFnn32W6dOnExcXR2BgYLaWwKw2b97MRx99xNy5c3E4HDz33HM0btyY/v3706JFC2bOnMmcOXM4c+YM06dP5+GHH2bw4MEsW7aMadOmkZiYyN13383s2bMzuoFFRIqSljsRESkCnTt31pp3IlLo1BUrIlKIEhISmDp1Kvv372fatGkkJCTYXZKIuDG12ImIiIi4CbXYiYiIiLgJBTsRERERN6FZscVd7BrwDgaf6pB4DC4cgsC22a85Pg1SzkLpBlC+sx1VioiISBFQi11xd2Q8rA6BJV7w5wAoVSH7z0/MgXKdodpwSD4FCUfsqFJERESKgFrsijvf2tDuMFgW+NbI/RpHWn738AEcRVWZiIiIFLFiEewsyyIuLo6yZctm7L8oWfhcYhX9Sv0hejqkXgCfEPCpVnR1iYiISJEqFsudxMbGEhgYSExMDAEBAXaX41r2PAd+TU1r3OmFUPM58G9id1UiIiJig2LRYieXUPlOCLjB3C9VBf68Ddpsy+x+FRERkRLD6X/9ly5dSr9+/ahWrRoOh4Mff/zxss9ZsmQJoaGh+Pr6UrduXSZMmHAltUpuyoRm3i9dH87vgLOb7atHREREbON0sDt37hzXX389H3/8cb6u37dvH3369KFjx45s2LCBl156iSeeeIIZM2Y4XaxcJGY1LC8HKRfM45Q4c+vhbVtJIiIiYh+nu2J79+5N79698339hAkTqFWrFmPGjAGgSZMmrFu3jtGjR3PHHXc4+/aSlU8NqPU8ePqaxzErIKA9+De1ty4RERGxRaGPsVu1ahU9e/bMdq5Xr158+umnJCUlUapUqRzPSUhIyLZRdmxsbOEWuX49nD4Nnp7g4WGOS9339QU/P3P4+0Muv0OR8K0BgR3h4GggBeJ3QrNZ9tSSVXw8/PWX3VWIiIgUjJYt7a4g3wo92B07doygoKBs54KCgkhOTubEiRMEBwfneM7IkSN54403Cru0TC+8AL/9duXP9/KCsmWhYkWoVMncph9VqkBIiDlq14aqVU1ALCjlu5jDlfz1F4SGXv66gvQG4AdkneNtZTkuPkcet/m5tjDeI69zF/88NY/7Vtrjy92/1HGp18/vexXEUVR1pF7mNj/XiUjJ4PoLiGQoklmxF689l77CSl5r0oWFhfHMM89kPI6NjaVmzZqFV+CUKaaVKTXVHCkped9PSYGEBHP9uXPmNj4eYmLg5MnMY/duWLMGjh+HM2cy38vbG2rVMkGvbl1o0gSuucYc1aqBO6zT17ixaQUtSuffACueSyesSyWqtFvr4uvz8Zwc73GJGqz8vl4+Dis/yYUs97P+jrmlpdwe5/Y7SXYemIW/PS66n+XWcamfp19zuddJ+4PQcZnXwSMf71dQr+PIUnseR8bPPNOu98zlZ7n8PMdrZvl5vt7P2dd0g//vFaEIgl3VqlU5duxYtnPR0dF4eXlRsWLFXJ/j4+ODj49PYZeW6aIWxQIXGwsHDsD+/dlv166Fr7+GC2mTHypVMi1d6Ufr1lCYgbaw+PnZ0Gw9u4jfrwSyLg6VqXncT3t82fsWWKm537/Uz7Ldv9o6UtNeL6/b/FyT31ur4F7nql4jqeBrSj9IASsly/1U8zj9vks3c6aH2fTg6HmJ+2nX5Xk//XXyuJ/Xzx2e4PDK8n5e2W/Jeq4Iz+da0yXOZ/xRIHYo9GDXrl07fvrpp2znFixYQKtWrXIdX+eWAgLg2mvNcbGUFBP0tm6FDRtg3Tr47DP473/Nz+vUgS5dMo/ql9hlQqQwOdJaZ8i8EXFK1jBtpaTdT8n7/uV+nmuYvESwzHHfFd4zyewMlPE6KWAlZ3le8pWftzVIOxEEnTp/cQjNcniUynmuoI5ynWz8LJ3jdLA7e/Ysu3fvzni8b98+Nm7cSIUKFahVqxZhYWEcOXKEr776CoARI0bw8ccf88wzzzB8+HBWrVrFp59+yrRp0wrutyjOPD2hXj1z9O+fef7oUVi9GiIjYfFiE/YAGjQwAa9rV+jRAypUsKXswrCFaSRylgo0oA6d7S7nqizmNZJJwBNvPCmFJ954pN164o0XPnjic0X3PZxfpUjENTjSuz09gRLyh72dLKtgAqKrnU9NzPJzZ4+k7I/zG347F5/hKE5vKRYZGUmXLjkH6w8dOpQvvviC++67j/379xMZGZnxsyVLlvD000+zdetWqlWrxgsvvMCIESPy/Z7aUgyIjoYlS0zIW7zYTFDw9ISOHU0g7N/fhEObJXCWHxlKLz6gHLUAOMJaNjOFYFpwkBV0JIzy1Mn2vL+YQ3VaU5ZgtjGTGtxAAMW3dXIybYnnb1JIIpUkUkgkJeM2Ea5izJoHXjkCn1fa48vd98I37SiNF76USrt19rGn/lEWEXeQ3pJ6uTDo19DuSvNNe8UWV4cPw9y5MGeOmdGbkGDG5d19NwwaBDVqFHlJf/AZMRxgCW/yJPsoT22SSeAjGjGcNZQhiEOsYiFh3E9ktuf+xRxqcANlCGInc6lKCwKoVuS/Q1GwsEglmWQS0oJeQi73zePc7jvzvIt/nsyFtON8xv0kzpNKklO/gwPPjJB4JeGwFH5ph3+W+9kP7yw/88CzkP7XEBFxLwp27uDsWZg/H6ZNg3nzIDHRtOTdfTcMHGiWXSlCr+PICHZ7iOBXnuURzDZnqaTwNv48zQHKkH3Syp9MJ5kLBBJS7Ltii5tUUrIFvazhLylLCMx5LvfHl7omiXiSOU8i57BIyVd9nng7HQZzHjl/Zq73xxt/hUcRcQtFstyJFLIyZUyAGzjQLLvy448m5D32GDz5JNx2GwwbZsblFeQaevlwhv2UJnMcoAee+FCWaLbmCHbNGFSktUkmDzzxTgs4RSmFJJKIJ5FzJBGfx5H3z9KfF8+JPJ9jur4vz4vSeFPmosM/l3Nl0sJgzvMXP8cLXxyaaSIiRUjBzt0EBsLQoeaIjoZvvoFPPjETLWrXhoceguHDzdIqRSCeE3jhm+2cF75c4EyRvL+4NjOxJBBfAgvtPVJIJpnzOcJgegBM5ByJnM1xJGU5H8PBi35+jkTisC4z8NqBhxOBsAw+BOBDWXwIwDvtNutjb8po8oyIXJKCnTurUgWeeQaefhpWrYLJk+HNN+GNN2DwYHjiiUJfb86HQC6eKJDIWfwommAp4okXnpTFh7IF+roWFskkXDIQ5gyDmY8vcIZYDmc5F0cCcaSQcMn3TQ+AuQW//D5WSBRxXwp2JYHDAe3bm2P0aNOCN24cfPkldOsGYWGmm7YQFpSsRGP+YHLG42QSSCCOcoQU+HuJFCUHDkrhSyl88S/AP1SSSUwLebEkpN3m5/EZ9uX4+eW6ofMKiT5prajZb8vlci5Q3c0iLkbBrqSpWNHsjfuvf8GMGTByJHTvbna5CAuDW28t0HF4IdzEOaKJ4TCB1GA/S6hOawU7kTx44Y0XFfHj6ic9pf8hlTUIXi4kJhDLOXZxgRgSiOECMSQSl+d7eOKdRxA0YTA/IdGLItxpSMTNaVZsSWdZ8OuvJuAtXWr2rn3hBbjnHnByZ5AtTOMAS1nHBK5hECHcRBseYQ8L2c4MatCW/SzhJl6hAnUL6RcSkYKWSgoJxGYLe+m3FziT41xC2vms55I4l+fre+KTRxAsR2nK40v5tNtyWe5nntO6iiKZFOwk08qV8M478NNPZo/al16CBx90OuCJiFwsheS0cJi/IJh+7jynucBpznM6z+VxSuF/Udi7VBAsny0wlrpocpdIcadgJzlt2WJa8L791uxm8Z//mKVUinipFBGRdBZWxqSTrGEv8/bMJc/lNSnFC988Q1/W29JUoDQVKU0F/KiIL+Xx1GgmcUEKdpK3TZvg5ZfNDhctW5qw16NHoUyyEBEpTEmcz1cAzHqbHiLz6kb2ISBb2Ls4/OX22JdyWgxbCpVLB7vw8HDCw8NJSUlh586dCnZ2WbYMXnzRdNV27Wq6a1u3trsqKSTnOE8KqZTFT7MdRTAzlS9wmnhOcp5TnE+7vdzjJOJzeTVHWqtg/sNgaSrgQ6CWp5F8celgl04tdi7AsszYu5degq1bzeSKUaOgenW7K5MCNo4feZQx+FCKICpQhXJUoTxBlKcK5TLOmcflqUoFKukfHZEc0lsJ8xsEzeOTuS5T48AjLexVSjsqZ7mfefhnOe9NGf1xVgIp2IlzUlLM+ndhYXDunAl6zzwDvhqA7C72EcXvbCeaMxznVNrtaaI5nXF7jgvZnuOFJ0FpIS+YilmOCtnOVaUC3prBKJInC4sk4vMIfieI5yTx/J12P/NI5GyO1/LE26kg6EclLT3jBhTs5MrExMBbb8HYsVCrFvzvf2YNPI2/KxHOcT4j+B3jFFGcIoqTafdPph2nOM4pUi7adqsCAQRfFABN+DPnqlOJ6lSmtP6BEcm3JC5wnpOcyyX0ZR7Zf5Zby6A3ZfIdBP2pgi/l1VrvYhTs5Or89ZfZsuyXX8zEivBwaNDA7qrERaSQwklis4W99PuZIdDcnr9o1mIFAqhBZWpQmepUynabfj9QXU0iVyR9lnF6yMtfIDzJxVtEOvDMCHmZR+WLHlfBL+2cD2X132whU7CTq2dZZubsE0/A0aPwyivw/PPg7W13ZVJMWFjEco4oTnKUkxzmb45wgsP8nXbf3B7ndLbn+eNL9bSgl1vwq0FlKlNOLQoiBSCVFC5wJksQ/JtzROdymPPxnODiIOiJT74CYPrPS1Hanl+2GFOwk4ITH2+6Z0ePhoYNYeJE6NDB7qrEjSSSRNRFwe/iAHiEEyRnWci2FF5UpxK1CCKEIGpRhVpptyFUpSaVKYOfjb+ViHtKJYV4TqaFvEuHwHNEk0BMjtfwpmyOAOiXSyAsQ1X8qKilZFCwk8KwZQs89BCsXm1u33kHype3uyopIVJJ5W/OZIS+Q0RzmL85SDQHOMZBojnCCVKzjP2rQEBG4Ls4/NUiiKpUUKufSCFLJiFLa19uQTDrueMkXzSJy4FHWtALogxVKZN2ax5nvV+V0m7837SCnRSOlBTTYhcWBqVLm/u33mp3VSIAJJPMUU5ygOMc5DgHic64PcAxDnCcs5zPuL4UXtSkykWtfUHUpip1CKYWQZTSLgQiRcaMETyXEfLOcjzt9hhn024zHx/LsaagB16XCIHZH5emfLEaF6hgJ4Xr6FEYMcKsgfePf5hZtBUq2F2VyCVZWMRwNi3o5R7+jnISK238kAce1KAyddKCXnrgM0dVqlHJbVsHRIqDBM5mBMDsoS/rY3M/OcsfdQAelOLVXGYQuyoFOyl8lgXffGMmV5QuDZMmwS232F2VyFVJJImDHGcfx9jPMfYRlXHs51i2iR7elKIWVbKFvfT7talKZcoVqxYBEXeVPlv44tDXhkftLi3fFOxcwbRpcPasWSakc2e7qyk8R46YMXfz5sHQoTBmDJQrZ3dVIoUingvZAp+5n/n4TJYFZf3xzdbKV49qaUd16lAVX63pJyL5pGCXVVwcfPIJ/PGH2TZr0aKrCx5r18KUKdCiBaxYYcab1amT/Zo5c8y+q8HBMHMm3HCDe2/TZVnwxRfw1FNQtqy53727zUWJFL0zxOUIe+Y4xl6OciGt68eBgxpUzgh79alOPapnPA6kjM2/iYi4Eo32zWrPHmjc2Cy4e/Lk1YW6hAQYOBDWrIGgILP8x/33Q2Rkzms90sbe+Pi4/84NDof5HLp3hwceMIsaP/ssvP22+f1FSohylKU5ZWlOzgW9U0klipPs4Si7OcIejrKHI2xkNzNYmq21rxKBGa179bMEvvpUp0oxG/QtIldPLXZZJSbCjTfCjz9efatZRIQJLJs3m8cpKeDvDwcOmKCX1fTpcOEChIS4d1fsxVJT4YMPTEvmNdfA1KnQpIndVYm4vFPEZgt8mfePEsXJjOv88c1o3aufFvwaUoOG1CSYigp9Im5ILXZZ7d1rukr79zctbV5pH8+pU/Dee6YbMS+lS8Nrr2U+3r8/++xPT0/T9bh1a85gN2hQgf0KxYqHhwm/3brB3XdDy5bw/vtmFq27t1yKXIUKBNCGANqQ8w+hc5xnL1HsSQt76aHvB5ZwgOMZ6/f540uDtJDXkBpp983jCrjhWGaREkLBLt3XX8O+fea2bl2YNct0pYIJaCNHOvd6J06Ar2/2c76+cOZMgZTrVpo3h/Xr4V//gkcegfnz4dNPoXLlIi9lDdsyNqI/xikOEU1brslx3TR+4yznaUB1OtOiyOsUyYs/pbmWulxL3Rw/SySJvRxlF0fYySF2cpidHGIFWzjCiYzrKhKQLfQ1pCYNqE4DauCvLZ5EXJqCHZhA99FHsGqVaVnr0gUOHry61wwMzNnCd/YsVKp0da/rrvz8YNw46N3bjL277jozW7gQu6bPEs9Q3uEDHqUWphV1PLP5kl8BaEUjpvBKjufNYQWdaU4wFZnJUo7wN9Up+hAq4ixvStGYEBoTkuNnZ4lnN0cywl56+PuJlZwmLuO66lRKC3qZLXwNqUEdgvGmVFH+OiKSC5cOduHh4YSHh5OSknL5i6/G7Nlw++0m1IFpVbvuusyfnzxp9j91piu2cWOYPDnzcUKCmXUbkvP/UCWLfv3MuMR77zVdtG++acbgeRTs4q6fMY8DHGcmS/kf/5dxvjZVOcz3WFjUoEqez/dIG5vkQymNUxK3UAY/mtMg18kcJ4lhF4czQt9ODrOGbXxDBPFp2zp54kEdgmlMrYyjCSE0ppa6dkWKkCZPgFlyIyYGnnwSoqNNq9HatVcXJpKTTYhbswZq1IAFC0zwW7WqwMp2aykp8MYb8J//QM+epou8ELpmHXRmH9OoTTAAr/M5r3P/ZZ83nUVcIJEQgtQVKyWWhcVRTmQEvh1px3YOsJ9jGTtzVKZcjrDXmFqEEKQdOUQKmIIdmBD22mtQr55pLXrsMahf/+pfd+FCmDED2raFJUvglVfM+D3JvwULzFZk3t7w7bfQoUOBvvzFwe45xtOU2vhQioWs5zkG0ySXbisRubTzJLCLw/zFQbZzgL84yF8cZAeHOE8CAL5404ia2Vr5GlOLhtTED9/LvIOI5EbBTlzfkSNm1uzKlWa9u+eeK7Cu2YuD3Rq2cQNNAVjIOh5jLNv4Uq0KIgUklVQOEZ0t7KUf6duwOXAQQlCurXzafk3k0hTspHhIToZXXzWzk/v2ha++yr6czBW6ONglk4xX2tDT/URRh7vZwORcxx2JSME6RSw7OJQt7G3nAHs5SkraMi0VCaAptbmG2hm311BbizGLpHHpyRMiGby84L//NV2xQ4aYbdhmzco+yeUqrWYr3XmWE8zGFx/iOA+gmX4iRaQCAbTjGtpdtMRQAons4SjbOcA2DrCN/azgTz5jPokkZTzXhL0QrqFOxv0gKijwSYmiYCfFS58+sG6dmcXcrh189lmBLfBcg8o8z90ZG66vYAvtaUZTahfI64vIlfHBm6ZpLXR3ZDmfTDJ7OMpW9rON/WxlPyvZyuf8ki3wmbCXvYVPgU/clbpipXiKj4fhw802ZM89Z1rzvPL/d8o0fmMpm5jAHAbRhZu4nkcYwGI2sJ4dpJDKTg4zkuFUoXwh/iIiUtDSA982DrCVfWnBz4zpSw985SmbI+w1pTZVFfikmFOwk+LLsmDMGBPsunQxs2YrVrS7KhFxUckks5cotrKfrexLC377cwS+ZtRJ273D3DajDoGUsbl6kfxRsJPib/FiuOsuKFPGjLtr3tzuikSkGLk48P3JPrawjx0czJi0UZMq2cLetdSlMbU0BldcjoKduIeDB+G222D7dvjkE7jnHrsrEpFiLoFE/uIgW9jLFval3e7lMH8D4IUnjaiZEfTSg18IVdWdK7ZRsBP3cf48jBhhlkJ5+mkYNcqpcXciIvlxmji2prXqpYe9LewlhnMAlMUv1+7cigTaXLmUBAp24l4sC8LDTbDr2hWmT4dy5eyuSkTcnIXFYf7OFvS2sI/tHCCJZACCqZijO7cpIRkz8UUKgoKduKfFi+GOO6BKFfjpJ2igBYZFpOglkcxODrGFvRlj97awl31EAeCJB42oxfXU43rq0Zz6XE89qqKJYHJlFOzEfe3aBf36QXS02bO3Sxe7KxIRASCOeLayj83sZRN72MhuNrOHs2kLo1ehfI6w14halNLys3IZCnbi3s6cMTNmFy+Gjz+Ghx+2uyIRkVylkso+ojKC3ib2sIndHOA4YHbBuYba2cLe9dSnPGVtrlxciYKduL/kZDPm7uOP4Ykn4H//06QKccouYCDQFrgh7bYR4GFnUVJinCGOzezNFvb+ZB8JaWvv1aTKRa179alHNTz0DS2RXDrYhYeHEx4eTkpKCjt37lSwk6szfjw8/jh0724WM9akCsmnXcA7wBpgG2ABgUAboD3QDhP4ytlUn5Q8ySSzk8PZwt4m9nCMUwD448u11M1o1WtOfa6lDmXws7lyKWwuHezSqcVOCszChTBwIFStCj//DPXq2V2RFDOxwFpgddqxCjgJOIBrMCEvPew1TDsvUlSOcyot6GWGve0cIIVUHDhoSA1a0ICWNKQF9WlBAy3D4mYU7KTk2bkT+vY14+9mz4b27e2uSIoxC9OitxIT8lYCW9POV8R027ZPO1oD/vaUKSXYBRLYxgE2spsN7OIPdrGJ3ZzjAgC1CMoIeemBrzqVtchyMaVgJyXTyZMwYACsXWsWNL7rLrsrEjcSg+m2TQ96qzEtfZ7A9WRv1auNWvWk6KWQwm6O8Ae72JDlOEksAJUpRwsaZAt8GrdXPCjYScmVkAAPPABTp8LIkfDCC+DQP7FS8FKA7WRv1duZ9rOqmJDXAegINActaCG2sLA4RHRGyEsPfelbqJXFj+upl60btym1tQSLi1Gwk5LNsuC11+Ctt+DBB80Ei1La1FsK3wkyx+itwLTwXcB01bbDhLyOmEkZGu4udvqbM2xkN3+wMy307WYXh7Gw8KYU11InWzfuddTDD1+7yy6xFOxEAL78EoYPh06dzGLG+p5JEUsA1gPLgWVpt2cwrXehZAa9G0F7Eojt4ohnE7vZkDZubwO7+JN9JJOCBx40phYtaUAoDWlFI5pTXzNyi4iCnUi6xYvhttugdm2YNw+qVbO7IinBUjFLqyzLchxO+1lTTMhL774NsaNAkYskkMhW9md04/7BTjaymwsk4sBBY2rRikaE0pBQGirsFRIFO5Gs/vwTevcGDw+YPx+aNrW7IhHAzLI9SGbIW44JfgA1yQx6N2GCn0aLiitIIpntHGA9O1nHDtazg03s4QKJGS176a166WHPn9J2l12sKdiJXOzwYejTBw4dgjlzoGNHuysSydUJzPi89KC3HkgGKgOdgM5pt03RLhniOpJIZhv7Wc/OjMC3id0kkJQR9i5u2VPYyz8FO5HcxMSYbtmVK+Gbb+DOO+2uSOSyzmEmZESmHWuAJKAS2YPeNSjoiWtJD3umVc8Evqxhrwm1CKVRtpY9TdDInYKdSF4SEuC++2D6dHj/fXjqKbsrEnFKPGbW7RJM0FuNCXoVMQEvPew1Q0FPXE8SyWxlX5Zu3J1sYg+JaWGvKSFprXom8F2v2biAgp3IpaWmwosvwnvvwbPPwqhRZvydSDEUT2aL3pK0+4lABbIHvWtR0BPXlEgSW9nP+rSgt44dbGZvtrDXika0pjGtacx11MUHb7vLLlJXFOzGjRvHe++9R1RUFNdccw1jxoyh4yXGIU2ZMoVRo0axa9cuAgMDufnmmxk9ejQVK+Zv0r6Cndjuww9Ni92gQfDFF+DjY3dFIlftPNm7btODXnmyB73rUNAT15Ue9tInZ6xlB5vZQzIpeFOK66mXFvRM4GtMLTzxtLvsQuN0sJs+fTpDhgxh3Lhx3HjjjUycOJFPPvmEbdu2UatWrRzXL1++nE6dOvHBBx/Qr18/jhw5wogRI2jQoAGzZs3K13sq2IlLmDED7r0X2raFH3+EcuXsrkikQJ3HjMuLJDPoJWCCXmegG9AVaIxm3Ypru0ACm9jDWv5iLTtYy1/8xUEsLMpQmlAaZrTqtaYxtanqNnvjOh3sbrjhBlq2bMn48eMzzjVp0oQBAwYwcuTIHNePHj2a8ePHs2fPnoxzH330EaNGjeLQoUP5ek8FO3EZy5dD//5mjbv586FmTbsrEik0FzBBbzHwGyboJQPBmICXHvS0jp4UB7GcYz0708KeOQ5wHIBKBNKKRrShSUbLXhAVbK74yjgV7BITE/Hz8+P777/ntttuyzj/5JNPsnHjRpYsWZLjOStXrqRLly7MmjWL3r17Ex0dzV133UWTJk2YMGFCru+TkJBAQkJCxuPY2Fhq1qypYCeuYft2s9ZdcrIJd9dea3dFIkXiLGZZlUWYoLcBs75ePTKDXhegil0FijgpmtNZgp5p2fubMwDUpAqtaUwbGvMC99hbqBOcCnZHjx6levXqrFixgvbt22ec/+9//8uXX37Jjh07cn3eDz/8wP3338+FCxdITk6mf//+/PDDD5TKY0/O119/nTfeeCPHeQU7cRlRUWatu3374KeftNadlEinMF226UHvr7Tz15LZmtcJ0P9rS3FhYXGQ4/yepVVvHTuIY77dpeXbFY2HdTiy90NblpXjXLpt27bxxBNP8Oqrr7J+/Xp++eUX9u3bx4gRI/J8/bCwMGJiYjKO/HbZihSZ4GCIjIQWLaBnTxPuREqYCsDtwMfAduAI8DXQCpgJ9E+7pi3wMib8nbelUpH8ceAghKoMpDOjGMFixnCGn+0uyymF3hU7ZMgQLly4wPfff59xbvny5XTs2JGjR48SHBx82ffVGDtxWRcuwD33mB0qPvnErHsnIljAXkyYW5R2/A34AO0xLXrdMCHQy6YaRdyRUy123t7ehIaGEhERke18REREtq7ZrOLj4/G4aN0vT08zzbgYLKEncmm+vvD99/DAA3D//Wa9OxHBgRl79xDwLXAM2AS8A/gD7wLtMLti3AFMAPbk+koi4gyn/1B65plnGDJkCK1ataJdu3ZMmjSJgwcPZnSthoWFceTIEb766isA+vXrx/Dhwxk/fjy9evUiKiqKp556ijZt2lCtWrWC/W1E7ODpCRMnQlAQPP88REebhYzzGJ4gUhJ5YNbDuw54CjO79ncgIu14DEgB6gA9gJ6YMXrlbahVpDhzOtgNGjSIkydP8uabbxIVFUWzZs2YN28eISFmwntUVBQHDx7MuP6+++4jLi6Ojz/+mGeffZZy5crRtWtX3n333YL7LUTs5nDAW29B5crw5JPw998weTLkMUFIpKTzwnTJtgdeA2IwEzHSg94kTBhshQl6PTAtfCVrDwER52lLMZGCNnUqDB0KN99s9pn187O7IpFi5yCZIW8hcBLThduZzKDXBC2ULHIxBTuRwvDLL3DHHWbW7E8/QXl1KIlcqVTMmnnpQW85Zuuz6mSGvO5o/TwRULATKTxr1pi17qpVg19/NbcuaBpm4dkGmNYQEVcXDywlM+htSTt/PWZsXg+gA1DalupE7KVgJ1KYtm8369x5esKCBdCwYYG/xRpMC0YcsAp4CbMo7FpgCtACWAGEYQamZzUHaI3ZImomcAOmFUSkODmG6a5dgAl6xzDLqnQEeqUdzVC3rZQMV7RAsYjkU5MmsHKlGWfXoQOsX1+gLx8P/Ag8C7wODAd6AweAgZgwNxS4P+3ITfr/Cfigf/ikeKoK/AP4CjiKacEbiZmg8SpmJm4N4AHgO8yOGSLuSi12IkXh5Eno2xe2boUff4Ru3QrkZTdjup92AfUxrXYBmHXD3k77OZhlJPwxgS/ooteYjtnsPQR1xYr7uQAsA34BfgW2Yv6YaQPcjGnNaw142lWgSAFTi51IUahYEX77DW680Yy7mzWrQF72Wkw3a720xwez3FbIcp0nUBbzj9rFBmFa9ToXSEUirsUXM+buf8CfmP82JmFa8D7ALKFSBRgMfI5p8RMpzhTsRIqKv7/ZemzAALjzTvjii6t+SQdmHbD0LtR3yFz81feia32BM1f9jiLFW03gQeB74ATmD6NHgX1p56tjum6fx2yDlmBPmSJXTMFOpCh5e5t17oYNM1uQjRlTYC/9GWas0ftAIGavzqzOYrZvEhEjfZHkNzGTkKIxs8RbAl9j9rKtAPQDPgZ221OmiFO097JIUfP0hAkTzNp2Tz8Np07BG29c1RZk8zDj6N7DjCnyxfwjlS4BM/4u5CrKFnF3lTBdsoMxfxhtJnNs3jPA40BdzNi8m4EuQBlbKhXJm4KdiB0cDnjnHRPuXnwRTp+GsWPBw/lG9KWYcUH9MMs8rMasSRcNHMaMJVqCGSCuYCeSPw7MxKTrgRcwfxhFYoLeL8A4oBRmvbxemNno16KZ5WI/l54VGx4eTnh4OCkpKezcuVOzYsU9TZ4MDz8M994Ln33m1P6ye4HmmH90sorBbLA+A2iLCXavYFobROTq7SazNW8RZumh6kAfTMjrjpmwJFLUXDrYpdNyJ+L2vvsO/vGPzP1lS2vNfJHiIgGzpMq8tGMHpjWvIybk9UH72krRUbATcRW//AK33w5t2pjZs/quixRLe4H5mJC3GDiPGQaR3prXFbOupEhhULATcSUrVpiFjOvXh/nzoXJluysSkatwHjMUIr01bw/gjVk3Mr01rwFqzZOCo2An4mo2bTL7y1aoYPaXrVnT7opEpIDsIjPkRQKJmAXG01vzOgMaiCFXQ8FOxBXt2gU9eoBlQUQENGxod0UiUsDOYbpq04PeAcxSRV3JDHqa8CTOUrATcVWHD5uWu5Mn4ddfoXlzuysSkUJiAX+RGfKWAUlAIzJD3k2Aj10FSrGhYCfiyk6cgN69YedOmDsXOnSwuyIRKQJxwG9kBr0jmAkX3cgMerVsq05cmYKdiKuLjYVbb4U1a2DGDBP0RKTEsIA/yQx5KzA7zVxDZsjrgFliRUTBTqQ4uHABBg2CefPgm2/MfREpkc4ACzEhbz5mx5kAzDZn/TBBr6JdxYntFOxEioukJHjgAZgyBcaPN7tViEiJlgpsBH4GfgLWAR5Ae0zIuwUtjlzSKNiJFCepqfDUU/DRRzBypNlnVkQkTRQwFxPyIjDr6NXFBLx+mAkY3rZVJ0XBy+4CRMQJHh4wdqxZ4y4sDE6dgnffBYf+HhcRCAaGpR3nMcup/AzMAj7E7F/bCxP0+gBaAt39qMVOpLgaO9a03g0bBhMmgKen3RWJiIuygM2Ylryfgd/Tzrcls8u2GeqydQcKdiLF2VdfmXF3t91mJlX4aJUrEbm845jJFz8BCzCLJYeQ2WXbGa2ZV1wp2IkUd7Nnm1myN90EM2dCmTJ2VyQixUgCZnuz9AkYBzBr5vXEhLw+QJBdxYnTFOxE3MHixdC/PzRrZhYyrlDB7opEpBiygK1kdtmuSjvfhszWvOtQl60rc+lgFx4eTnh4OCkpKezcuVPBTuRS1q2Dm2+G4GBYsMDciohchb8xXbY/A79idsSoiQl5t2D2tfW1rTrJjUsHu3RqsRPJp+3bzf6ypUpBRATUq2d3RSLiJhKBpZjWvJ+AfYAf0B3TktcXMytX7KVgJ+JuDhyAHj0gLs603F17rd0ViYibsYDtZI7LW4lZLLkVmbNsW6AuWzso2Im4o+ho6NUL9u8325C1a2d3RSLixk5itjf7GfgFiAGqkTkuryumdU8Kn4KdiLuKiYF+/WD9epg1y3TRiogUsiRgGZmtebsx4/C6kzk2r7pt1bk/BTsRdxYfD3fdZbpkp0yBgQPtrkhEShAL2EnmLNvlQAqmmza9yzYUs7+tFAwFOxF3l5QE990H06bBxIkwfLjdFYlICXUKM7v2J0zX7RmgKmbiRT9Mq56/XcW5Ce0VK+LuSpWCr7+G8uXhoYfM/rIvvGB3VSJSAlUA7k47koEVZHbZforZ7aIrmbNsa9lTZrGmFjuRksKy4LXX4K234Pnn4Z13wKE5ayLiGnZhQt7PmGVVkoHryZyA0Rp12eaHgp1ISTNmDDz9NAwbBhMmgKen3RWJiGRzBtNl+zNmgeRTQBVMK94tQA+grF3FuTgFO5GS6Msv4cEH4fbbTTetj7b7FhHXlAysJnMCxjbAG+hM5gSM2jbV5ooU7ERKqtmzYdAg6NQJZs4Efw1ZFhHXtweYiwl6SzDLqzQjs8v2BqAk90Mo2ImUZIsXQ//+ZneKn3+GChXsrkhEJN9igQWYlry5wAmgEtAHE/R6ASUtNSjYiZR069bBzTdDcLBZ7y5Yuz2KSPGTAqwhc5btn0ApoBOZrXl1bauu6CjYiQhs3272l/XxgYgIqFsS/u9PRNzZfjJn2S4GEoEmZI7La4d7rvmmYCcixoEDJtydPQu//mq6Z0VE3EAcsBDTkjcXiMasqdcbE/R6AeXsKq6AKdiJSKboaOjVC/bvh3nzoF07uysSESlQqcBaMrtsN2Fa7jqS2WXbwLbqrp5Lr/UXHh5O06ZNad26td2liJQMVapAZKRpreve3Yy5ExFxIx6YmbNvARuBg8CHQGngJaAh0Bj4F5mzbosTtdiJSE7x8XDXXSbYTZ0Kd95pd0UiIoXuHPAbmWvmHcN00Z62sSZnuXSLnYjYxM8PZs2CgQPNWncTJ9pdkYhIofMH+gOTgSOYLtsnbK3Iee44IURECkKpUmZXiooVYcQIOH4c/v1v7S8rIiWCB9Aq7ShOFOxEJG8eHjB2LFStCi+/bMLdhx9qf1lxK3/uhRV/wsDOUEGjfaSYU1esiFyawwEvvQSTJ8OECTB4MCQk2F2VSIFZuhkeGQNVb4dbX4bvFsN5fcWlmNLkCRHJv9mzTbBr1w5+/BH036O4ieOnYPpimLIQft8OZf3g9o5wbw/o0hy81L8lxYSCnYg4Z9ky6NfP7E4xfz4EBdldkUiB2nUYpi40IW/XYQgqD4O7mpDXqpGGmYprU7ATEedt2WIWMvbzM7tU1Ktnd0UiBc6yYN0OmBIB3y6C46ehYU24pxvc2x3q17C7QpGcFOxE5Mrs32/CXUyMablr0cLuikQKTXIyLNpgWvJmLoO4eGjTxAS8QV0gqILdFYoYCnYicuX+/hv69IEdO8z4uy5d7K5IpNCdT4CfVpqu2vlrICUVuoeakHdbRzM+T8QuCnYicnXOnoU77jBbkU2Zol0qpEQ5GQM/LDEhb9lmKO0D/dub8Xi9WoN3KbsrlJJGwU5Erl5iItx3H3z7LYSHw//9n90ViRS5A8fMWLwpC2HLXrMm3l2dTUte+2ZmWUiRwnZFX7Nx48ZRp04dfH19CQ0NZdmyZZe8PiEhgZdffpmQkBB8fHyoV68en3322RUVLCIuyNsbvvkGnngCHnkEXnvNjDwXKUFCqsIL98Dmz8wxvC/MXQ0dn4C6d8NLk81iyCKFyekWu+nTpzNkyBDGjRvHjTfeyMSJE/nkk0/Ytm0btWrVyvU5t956K8ePH+c///kP9evXJzo6muTkZNq3b5+v91SLnUgxYVnw7rsQFgYPPWRa77QAmJRgqamwfItpxfs+Ek7HwXX1TCve3d2gZhW7KxR343Swu+GGG2jZsiXjx4/PONekSRMGDBjAyJEjc1z/yy+/MHjwYPbu3UuFClc2bUjBTqSY+eILGDYMevc23bP+/raVsmYbBFeE6pXg2Ck4FA1tr8l53bTf4Ox5aFAdOmuCrxSCxCT45XcT8uasgIQkuOk6Mx7vzk5QvqzdFYo7cKorNjExkfXr19OzZ89s53v27MnKlStzfc6cOXNo1aoVo0aNonr16jRs2JB//etfnD9/Ps/3SUhIIDY2NtshIsXIfffB3LmweDF07QrR0YX+lmu2wf+mw+ufQ6/nYMlGc378bAgZBF7dYMArue8FOmcFdG4Ow2+BU3Fw5O9CL1dKIO9S0P9GmP4aHJ8Fn79gzo1432xndtsr8EMkXNB2ZnIVnOojOXHiBCkpKQRdtNJ8UFAQx44dy/U5e/fuZfny5fj6+jJr1ixOnDjBI488wqlTp/IcZzdy5EjeeOMNZ0oTEVfTqxcsXWqWQ2nfHn75BerXL5S3ir8APy6HkQ+Zxz9EQu8XYNc3ULsqHP7e9BLXuES3l0fabgI+pbSzgBS+AH8YerM5ok7C9LRJFwNfNz+7vaPpru3SAjw97a5WipMrmjzhuOj/9SzLynEuXWpqKg6HgylTptCmTRv69OnD+++/zxdffJFnq11YWBgxMTEZx6FDh66kTBGxW8uWsGqVGWfXrh2sWVMob7P7CLwzFXYfNo97tTFrja340zyuXvnSoa7/jRC5Eb78Bfx9oVqlQilTJFfBFeGpgbB2Iuz4Gp6+04zL6/EvqHkXPBMO63doPpLkj1PBrlKlSnh6euZonYuOjs7RipcuODiY6tWrExgYmHGuSZMmWJbF4cOHc32Oj48PAQEB2Q4RKabq1IEVK6BhQ7OA8U8/FfhbXFsXVnwM9aqbxwePm9sGNeDcBfh8vtkx4IF3YfuB3F9jUFfTeqLxdWKnhjXh9fth5zewZjwM7Gxa8lo9DE3+CW99BXuO2F2luDKngp23tzehoaFERERkOx8REZHnDNcbb7yRo0ePcvbs2YxzO3fuxMPDgxo1tNGeSIlQsSIsXAg33wwDBsDEiQX68g6HWScsvePgnanw1J3QooEZlH5/b7inu9nj87ZXzExFEVfmcJgty8Y+Dke+h19GmcejpkH9e6HdI/DRTIg+bXel4mqueLmTCRMm0K5dOyZNmsTkyZPZunUrISEhhIWFceTIEb766isAzp49S5MmTWjbti1vvPEGJ06cYNiwYXTq1InJkyfn6z01K1bETaSkwNNPw0cfwcsvw1tvFfiAts/mmVa5USPMSycnZ664sj8K6twNGyZD8wYF+rYiRSL+gpnsM2WhmWFrWdCjlRmPN6ADlNF2ZiWe0wtMDRo0iJMnT/Lmm28SFRVFs2bNmDdvHiEhIQBERUVx8ODBjOvLlClDREQEjz/+OK1ataJixYrcdddd/Oc//ym430JEigdPTxg7FmrWhOefh8OHYfJkKFUw+y7NW2327Xzv/8zMwo27ofuzcGI2+PpAXNqwXm3zJMWVny8M7maOkzFmbbwpC2HIf83Pbr3RhLyeraGUlpAskbSlmIjYY+pUsyxK587www9wlf9tL90EOw/BLe3M49XboGIgLN4Arw415ybMhq8jzHg8EXeyPwqmLYIpEbB1P1QKhLu6mJDX7hrN9C5JFOxExD6LF5sxd3XrmnXvqlW7opfZexSaD4O4+OznY+bC+p1mRmFKKuw8DCOHQ5XyV1+6iCuyLNi8x7TiTfsNDv9tlvy5p7sJeU1r212hFDYFOxGx15YtZocKT0+YPx+aNrW7IhG3kJpqWrKn/ma6bM+cheb1M7czq17Z7gqlMCjYiYj9Dh824e7QIZg50+xWISIFJiER5q8xLXk/rYTEZLPbyr3d4Y6boJy2M3MbCnYi4hpiYmDgQNM9O3myGX8nIgUu5izMXGbG4y3aYCZZ9G1rQl7ftmaikRRfCnYi4jqSkuDRR02we+UVePNNjfoWKURHT8C3i8wC3ut3mu3M7rjJrPmo7cyKJwU7EXEtlgWjRsGLL8I998Bnn4GPmhBECttfB8zM2qkLzTZ9QeXNjiz3dDOLI+tvrOJBwU5EXNP338OQIdC6Nfz4o9m9QkQKnWXBuh1mVu23iyDqJNStBnd3NbNrNbPWtSnYiYjrWrUK+veH8uVh3jyoX9/uikRKlJQUWLLJhLwflpiZtdfVM614g7tCSFW7K5SLKdiJiGvbswf69oUTJ0zLXYcOdlckUiIlJMKva01X7ZyVcD4BOlxrlk4Z2Bkql7O7QgEFOxEpDk6dgttvNy14X3wBd99td0UiJVpcvNmzdupv8Ovv5lzP1ibkDegAZbVnrW1cOtiFh4cTHh5OSkoKO3fuVLATKckSEmD4cPj6a3j9dXj1VY3mFnEBf58x3bTTfoNlm8HXG/rfaMbk9b4BfLztrrBkcelgl04tdiICmFHdb78N//43DB5sZsyWLm13VSKS5uBxM+Fi2m+wcTcE+sMdncyYvM7NtXxKUVCwE5Hi54cf4J//hGuvNePugoPtrkhELrJtvwl4036DPUehagUY1MXMrG3dWA3uhUXBTkSKp3Xr4NZbwcMD5syBFi3srkhEcmFZsPYvM+ni20Vw/DTUq2YC3t3doEmI3RW6FwU7ESm+jhwx4W77dpgyBQYMsLsiEbmElBSI3GgmXcxYAjHnoHl9E/IGd4WaVeyusPhTsBOR4i0+HoYONd2zI0fCCy+oj0ekGEhIhPlrTMj7aSVcSISO15nxeHd2gkrl7K6weFKwE5HiLzXVzJR96y0z9m7SJG1DJlKMxJ6D2StMd23EOvO3WY9WphVvQAezh63kj4KdiLiPqVPhgQegVSuYOROqqF9HpLj5+wx8t9iMx1u+BXxKQZ+2JuT1bQv+mgh/SQp2IuJeVq82Y+1Kl4affoJmzeyuSESu0KFo+D7ShLy1f4GfL/Rvb0LezW20Rl5uFOxExP0cPAj9+sHevTBtGtxyi90VichV2nMEvos0IW/zHtM9e1sHE/K6hUIpL7srdA0KdiLins6ehX/8wyyF8p//QFiYJlWIuIntB2D6Ipi2CHYegooBZiHkwV3hputK9kLICnYi4r5SU+GNN+DNN+Guu8xOFf4ahS3iLiwLNu02rXjTF8P+Y2Yh5IGdTchr29QsdVmSKNiJiPubOdPMlq1f3+xUUbu23RWJSAGzLPh9uwl530XC0RNQKwju6mxCXsuGJaPRXsFOREqGLVvMYsaxsWbNu86d7a5IRApJaqqZUfvtIvhhiZlpW7+62dJscFdoVtfuCguPgp2IlBwnT8KgQRAZCWPGwKOPlow/4UVKsORkWLzRhLyZS+HMWbimNgzqaoJew5p2V1iwFOxEpGRJTobnn4cPPoAHH4TwcC1mLFJCJCbBgrUm5M1eAWfPmy7awV1Nl21IVbsrvHoKdiJSMn31FTz0ELRsCTNmQHCw3RWJSBGKvwDzVptJFz+vMluatbvGtOLd2QmqV7a7wivj0sEuPDyc8PBwUlJS2Llzp4KdiBSs33+H224z92fOhBtusLceEbFFXLzZr3bab/DrWkhKhhubmdm1xS3kuXSwS6cWOxEpNFFRcMcdsH49fPihacXTuDuREutMHMxZaXa8SA95VqTdVeWfgp2ISGIiPPOMGW93330wbpzZkkxESrT0kPfPXnZXkn8KdiIi6b7+Gh5+GBo3NuPu6tSxuyIREaeUsPWYRUQuYcgQWLUKYmIgNBTmz7e7IhERpyjYiYhkdf31sG4dtG8Pffua7chSU+2uSkQkXxTsREQuVr48zJlj9pl9/XXo1w9On7a7KhGRy1KwExHJjYcH/PvfMG+e6Z5t1Qo2brS7KhGRS1KwExG5lJtvNkuhBAZCu3ZmYWMRERelYCcicjl16sCKFTB4MAwdCo88AgkJdlclIpKDgp2ISH6ULg2ffQYTJsCnn0KHDrBvn91ViYhko2AnIpJfDodZ527lSjh1Clq0gB9/tLsqEZEMCnYiIs4KDTXj7rp2NXvNPvOM2b1CRMRmCnYiIleiXDmzO8WYMfDxx3DTTXDwoN1ViUgJp2AnInKlHA548klYvhyOHTNds3Pn2l2ViJRgCnYiIlerTRv44w+48Ua45RZ44QVISrK7KhEpgRTsREQKQoUKMHs2vPcevP8+dOyoWbMiUuRcOtiFh4fTtGlTWrdubXcpIiKX53DAv/5lumajo6F5c5g+3e6qRKQEcViWZdldxOXExsYSGBhITEwMAQEBdpcjInJ5MTFmaZTp0+HBB2HsWPD3t7sqEXFzLt1iJyJSbAUGwrRpZjHjadPMXrObN9tdlYi4OQU7EZHC4nDAAw/AunXg7W0mWYSHg+t3lIhIMaVgJyJS2Jo0gTVr4KGH4LHHzKLGJ0/aXZWIuCEFOxGRouDrCx9+aGbOLltmJlYsXWp3VSLiZhTsRESKUv/+sGkT1K0LXbrAyy9rOzIRKTAKdiIiRa1GDVi0CN58E0aNgvbt4a+/7K5KRNyAgp2IiB08PU1r3apVcPYstGwJ48ZpYoWIXBUFOxERO7VqZbYju/9+ePRR6NvX7DsrInIFFOxEROzm52eWQZk714S8Zs1g1iy7qxKRYkjBTkTEVfTpA1u2mH1mb78dhg2DuDi7qxKRYkTBTkTElVSuDDNnwiefwLffmmVRVq2yuyoRKSauKNiNGzeOOnXq4OvrS2hoKMuWLcvX81asWIGXlxfNmze/krcVESkZHA6zv+ymTRAUBB06wCuvaFkUEbksp4Pd9OnTeeqpp3j55ZfZsGEDHTt2pHfv3hw8ePCSz4uJieGf//wn3bp1u+JiRURKlHr1zCLGr78O776bOdFCRCQPDstybm79DTfcQMuWLRk/fnzGuSZNmjBgwABGjhyZ5/MGDx5MgwYN8PT05Mcff2Tjxo35fs/Y2FgCAwOJiYkhICDAmXJFRNzDxo1w332wdSu89JJZKsXb2+6qrtqWaZB4Fio0gDqd7a5GpPhzqsUuMTGR9evX07Nnz2zne/bsycqVK/N83ueff86ePXt47bXX8vU+CQkJxMbGZjtEREq05s3h999NoPvvf6FNGxP2bLYvEv6cDn98Bt/dBUfTGhSPrIX5T8HGL2HOQ3B6X87n/jUHaneG0OFw/hTEHinCwkXclFPB7sSJE6SkpBAUFJTtfFBQEMfyWHdp165dvPjii0yZMgUvL698vc/IkSMJDAzMOGrWrOlMmSIi7snb23TL/v67Wci4dWvz2Maxd9/dAckXoOUDUOMG+P4uSE6A7wZCxzBoPhRa3A8/3p/78x1p/wp5+QCOIitbxG1d0eQJhyP7f32WZeU4B5CSksI999zDG2+8QcOGDfP9+mFhYcTExGQchw4dupIyRUTcU4sWsHYthIXBf/5jWu82bbKllPsioemdmY8Tz8KBpeATAGXS2gCqt4HDq+Hs8ezPbdwf9keaVr1S/hBQraiqFnFf+WtCS1OpUiU8PT1ztM5FR0fnaMUDiIuLY926dWzYsIHHHnsMgNTUVCzLwsvLiwULFtC1a9ccz/Px8cHHx8eZ0kREShZvb7PX7IABZuxdq1Zm5uxLL0GpUkVWRtC1mfd3/AQ93oUz+6F0hczzHp7gUxait2aGvXTNBhVJmSIlhlMtdt7e3oSGhhIREZHtfEREBO3bt89xfUBAAFu2bGHjxo0Zx4gRI2jUqBEbN27khhtuuLrqRURKupYtYd06ePFFeOstuOEG2Ly5SEs4/DssehWqt4ZrBkH8CfDyzX6Nly9cOFOkZYmUSE53xT7zzDN88sknfPbZZ2zfvp2nn36agwcPMmLECMB0o/7zn/80L+7hQbNmzbIdVapUwdfXl2bNmuHv71+wv42ISEnk7W1C3Zo1kJRkWu/efLPIxt7VaANd34Ty9eDzjqZblYvWW0g8C36ViqQckRLN6WA3aNAgxowZw5tvvknz5s1ZunQp8+bNIyQkBICoqKjLrmknIiKFIDTUtN49/7wJdi1bwurVhfZ2h1bDe0GZM15rd4aj66BMVTgXnXldcgIkxEG5kEIrRUTSOL2OnR20jp2IiJM2bTJ7za5fD489Bm+/DWXLFuhbHFkHvzwFQxeBlzfsnAffDoBHt8IXnWHYGgisAbsXQORrMEw7o4kUOgU7ERF3lZICH35oJlVUrAjjx0PfvgX6Fttmwum9ZoLEweUQ+jDU7wl7FsL2GVCjLexfAje9AhXqFuhbi0guFOxERNzd/v0wYgT8+isMHgxjxpg9aEXE7VzROnYiIlKM1K4N8+fDN9/AwoXQpAl8/rlZ5FhE3IqCnYhISeBwwL33wvbtcMst8MAD0KMH7Npld2UiUoAU7ERESpJKleCrr0y37N690KwZvPYanD9vd2UiUgAU7ERESqKePeHPP83SKCNHwrXXmu5aESnWFOxEREoqPz+zsPGWLWYcXp8+cOedoP25RYotBTsRkZKuUSOIiIBp02DFCjO5YvRos4uFiBQrCnYiImImVwweDH/9BQ8+CC+8YHauWL7c7spExAkKdiIikikwEMaONVuT+ftDx45w//3w9992VyYi+eDSwS48PJymTZvSunVru0sRESlZWrSAlSth4kSYPRsaNoSPPoLkZLsrE5FL0M4TIiJyaX//DS+/DJ98Ak2bmha9bt3srkpEcuHSLXYiIuICKleGSZNM92xgIHTvDnfcYbYqExGXomAnIiL5kz6ZYsoUWL3azJ599VWIj7e7MhFJo2AnIiL553DAPffAjh3wzDPw7rvQuDFMn669Z0VcgIKdiIg4r0wZePtts/dsaKhZKqVzZ/jjD7srEynRFOxEROTK1a0Ls2bBggVw4gS0agVDh8Lhw3ZXJlIiKdiJiMjV69EDNm2C8ePhl1/M8ij//jecPWt3ZSIlioKdiIgUDC8vePhh2LULnnoK3nsPGjQwy6SkpNhdnUiJoGAnIiIFKyAA/vtf2LnTrHc3fLhZ8Dgiwu7KRNyegp2IiBSOWrXgm29gzRqz/l3PnnDzzbBhg92VibgtBTsRESlcbdrA0qUwY4ZZ1LhlS7j7btizx+7KRNyOgp2IiBQ+hwNuvx3+/BMmT4Zly8z6d48+CseO2V2diNtQsBMRkaLj5QXDhpkJFm+/DVOnQr16ZgZtTIzd1YkUewp2IiJS9EqXhuefh7174fHHYfRoE/Defx8uXLC7OpFiS8FORETsU748vPMO7N4Nd9xhwl7DhvD551oiReQKKNiJiIj9qleHiRNh61a44QZ44AG49lr4/ntITbW7OpFiw6WDXXh4OE2bNqV169Z2lyIiIkWhUSMT5n7/HWrWhLvugubNYeZMBTyRfHBYlmXZXcTlxMbGEhgYSExMDAEBAXaXIyIiRWXFCnjtNfjtNxPwXn8d+vc3s2xFJAeXbrETEZES7sYbYeFCWLIEypWDAQOgdWv4+Wdw/XYJkSKnYCciIq7vpptg8WJYtAj8/KBfPzMW75dfFPBEslCwExGR4qNLF9N6FxEBpUpB797Qvj0sWKCAJ4KCnYiIFDcOB3TvDsuXZ7bY9eoFbdvC7NmaZCElmoKdiIgUTw6HCXSrVsGvv4KvrxmDd/31ZkeL5GS7KxQpcgp2IiJSvDkc0LOn6aJdtswsk3LvvWYv2smTISHB7gpFioyCnYiIuI8OHWDePFi/Hlq0gIcfNluVjRkD587ZXZ1IoVOwExER99OypVnoeOtWMx7vX/+C2rXh7bfhzBm7qxMpNAp2IiLivpo0gS++MHvR3nUXvPUWhISYPWmPHLG7OpECp2AnIiLur3ZtCA+HfftgxAiYNMmcGzoUtmyxuzqRAqNgJyIiJUdwMLz7Lhw8CKNGmUWPr7sObr7Z7HChtfCkmFOwExGRkicgAJ5+GvbsgSlT4Phx6NHDjM2bMgWSkuyuUOSKKNiJiEjJVaoU3HMP/PGHabELCoJ//MPMpP3f/zTRQoodBTsRERGHA7p1MztZbNpkti4LC4Pq1eGRR+Cvv+yuUCRfFOxERESyuu46+PJLMw7vuedgxgwzu/bmm80aedqyTFyYSwe78PBwmjZtSuvWre0uRURESpqqVeH1103A++orOHEC+vY1O1p89BHExdldoUgODsty/SlAsbGxBAYGEhMTQ0BAgN3liIhISWRZZl/asWNNK56/PzzwADz2mBmTJ+ICXLrFTkRExGU4HNC+PUyfbtbDe/RR+PpraNAAeveGOXMgJcXuKqWEU7ATERFxVs2a8N//wqFD8OmncPIk3Hor1KljdreIirK7QimhFOxERESuVOnScP/98PvvsG4d9OwJI0dCrVowcCAsWqRFj6VIKdiJiIgUhNBQ+OQTOHoU3n8ftm41S6g0bgwffACnTtldoZQACnYiIiIFqVw5ePxxE+wiI6FFC3jhBbMm3v33w4oVasWTQqNZsSIiIoXt+HH47DOYNAn274dGjcyM2n/+0yyrIlJAFOxERESKSmoqLF5sJlzMnAnJyWZtvAcfhD59wMvL7gqlmFNXrIiISFHx8DDj7qZONTNnx46Fw4fNjNqaNU2X7Y4ddlcpxZha7EREROy2caPpqv3mGzh9Gm680bTiDRwIZcrYXZ0UI2qxExERsVvz5vDhh2ZG7bRp4Odngl1wsJlwsWiR9qiVfFGLnYiIiCs6cAC++MLsbrFnD9SoAffeC0OGwDXX2F2duKgrarEbN24cderUwdfXl9DQUJYtW5bntTNnzqRHjx5UrlyZgIAA2rVrx6+//nrFBYuIiJQIISHw2muwaxesXAn9+sHkydCsGbRsadbGO3bM7irFxTgd7KZPn85TTz3Fyy+/zIYNG+jYsSO9e/fm4MGDuV6/dOlSevTowbx581i/fj1dunShX79+bNiw4aqLFxERcXsOB7RrB+PGmQkXs2aZrctefNG04vXubSZjxMfbXam4AKe7Ym+44QZatmzJ+PHjM841adKEAQMGMHLkyHy9xjXXXMOgQYN49dVX83W9umJFREQucuoUfP+96apdsQL8/c3s2sGDoVcv8Pa2u0L7rYmFYG+o7gPHEuHQBWgbmPO6acfhbAo0KA2dyxd9nQXIqRa7xMRE1q9fT8+ePbOd79mzJytXrszXa6SmphIXF0eFChXyvCYhIYHY2Nhsh4iIiGRRoQI8/DAsX27G4IWFwebN0L8/BAXBsGGwcKFZK8/dRZ6G6dHwWRTctRX+iDPnxx+BkNXgtQQG/AkVSuV87pwT0LkcDK8Gp5LhSEKRll7QnAp2J06cICUlhaCgoGzng4KCOJbPfv7//e9/nDt3jrvuuivPa0aOHElgYGDGUbNmTWfKFBERKVnq1oWXX4YtW8zx2GNmO7MePcxWZo89ZgKgu86svWMrXEiFB4LhhgAT7gBq+8LhdnCoLawNhYZ+uT/fw2FufTzAUTQlF5YrmjzhcGT/rS3LynEuN9OmTeP1119n+vTpVKlSJc/rwsLCiImJyTgOHTp0JWWKiIiUPM2awVtvmUkXa9eaWbSzZ0PHjmZCxrPPmskY7hTyIpvDnZUzH59Nybxf3Qdq+Ob93P6VIPIMfHkM/D2gmk9hVVkknNq7pFKlSnh6euZonYuOjs7Rinex6dOn8+CDD/L999/TvXv3S17r4+ODj0/x/mBFRERs5XBAq1bmGDXKhLlp08xEi/ffN2vk3XYb3HEH3HRT8d7O7Nosizj/dALerWfun0uFz6NMS9zC0/BcTWjin/P5g/JubCpunGqx8/b2JjQ0lIiIiGznIyIiaN++fZ7PmzZtGvfddx9Tp06lb9++V1apiIiIXBkPD+jQAcLDzRZmy5bBoEHw889mi7PgYDMmb/58SEy0u9or83ssvLoPWgfAoLTWuzsrw/3BcE8Q3FMFbvsTUl1++d6r4vSs2OnTpzNkyBAmTJhAu3btmDRpEpMnT2br1q2EhIQQFhbGkSNH+OqrrwAT6v75z38yduxYbr/99ozXKV26NIGBucxMyYVmxYqIiBQCy4J162DGDHPs3g2BgWbNvNtvh5tvhtKl7a7SOROOwKfHYPH14OsBXmltWPvPQ501sCEUmpe1t8ZC5PQYu0GDBjFmzBjefPNNmjdvztKlS5k3bx4hISEAREVFZVvTbuLEiSQnJ/Poo48SHByccTz55JMF91uIiIiI8xwOaN0a3nkHdu40s2qfesrsXXv77VCpktmv9ttvISbG7mpztzoGglbAvvPmcedysC4OxhyGcsvhQtp4u7i0W2/33k1VW4qJiIhITjt3mla8mTNNq16pUtCpk2nN69fPLJLsCtbFwlO7YVFzE9rmnTRLm/x2PSw+A6/WNtdNOAJfH4cVLW0stvAp2ImIiMilHTgAP/1kjsWLISnJzL7t39+EvDZtzDg+u8z8G/aeB08HLI+Bh6tBzwqw+DSsj4MUYGc8jKwLVdx74WYFOxEREcm/2FhYsMCEvLlz4eRJqFIFbrnFhLwePcwuGGILBTsRERG5MikpsGqVCXlz5sBff4GPj5lp26eP2ce2bl27qyxRFOxERESkYOzendllu2yZ2c6sYUMzu7Z3bzNGr7jNsi1mFOxERESk4MXFwW+/mbXx5s+HQ4fA1xc6d84Meg0amJm5UmAU7ERERKRwWRZs3w6//GJC3tKlZiHkunVNyLv5ZhP4yrrv+nJFRcFOREREitbZsxAZmdmat2+f2dKsbVvo3t0cbdqYJVbEKQp2IiIiYh/LMmPzFi40x6JFcOYMlCljWvF69DBBr0kTddvmg4KdiIiIuI6UFPjjDxPyIiJgxQrTbRscnNma1707VKtmd6UuScFOREREXFd8PCxfntmit2GDOd+okZll27mzuVXQA1w82IWHhxMeHk5KSgo7d+5UsBMRESnp/v7b7H4RGWmO7dvN+QYNsge9GjVsLNI+Lh3s0qnFTkRERHJ1/LiZZbtkiQl6W7ea8/XqZQ96tWrZWWWRUbATERER9/H33yboRUaasLdlizlfp44JeB07mqN+fbecjKFgJyIiIu7rxAmzC0Z60Nu82czEDQqCm26CqVPNUituQsFORERESo6YGFi50oS9w4fhq6/srqhAKdiJiIiIuAkPuwsQERERkYKhYCciIiLiJhTsRERERNyEgp2IiIiIm1CwExEREXETCnYiIiIibkLBTkRERMRNKNiJiIiIuAkFOxERERE3oWAnIiIi4iYU7ERERETchEsHu/DwcJo2bUrr1q3tLkVERETE5Tksy7LsLuJyYmNjCQwMJCYmhoCAALvLEREREXFJLt1iJyIiIiL5p2AnIiIi4iYU7ERERETchIKdiIiIiJtQsBMRERFxEwp2IiIiIm5CwU5ERETETSjYiYiIiLgJBTsRERERN6FgJyIiIuImFOxERERE3ISCnYiIiIibULATERERcRMKdiIiIiJuQsFORERExE24dLALDw+nadOmtG7d2u5SRERERFyew7Isy+4iLic2NpbAwEBiYmIICAiwuxwRERERl+TSLXYiIiIikn8KdiIiIiKF7IcffuD666+nbNmyOBwOSpcuTXx8/GWfN3PmTLy9vXE4HNSsWZPmzZuze/fuPK9XsBMREREpZHfeeSebNm0iLCyM1q1bc+HCBfbu3XvJ55w8eZJZs2aRlJTELbfcwqFDh9i4cSP169fP8zkKdiIiIiJFZO3atQwdOhSAPXv2XPLa0aNH07JlSwD69OmTr9dXsBMREREpAikpKVy4cIHmzZsDXLJL9ddff6Vt27b88ccfAHTp0iVf76FgJyIiIlIE1q9fT8uWLalXrx6Qd4vduXPnWLhwIbfeeiuRkZEEBwfTuHHjfL2HV4FVKyIiIiJ5ioyMpEuXLlStWhV/f/88g93o0aN59tln2b17N4cPH+aee+7JcU1cXBxvvPEGO3fupGrVqlSoUIG6deuqxU5ERESkKKxYsYL27dsDULdu3Vy7YtesWUPt2rWpWrUqixcvBnJ2w546dYqbbrqJY8eOMXv2bCZNmkTDhg159tlnFexERERECltycjLJycn4+fkBUK9ePQ4ePEhycnLGNYmJiUydOjVjckV6sOvatWu21xo4cCAxMTFMnDgRh8MBQL9+/Th79qyCnYiIiEhhW7t2La1atcp4XK9ePZKTkzlw4EDGubFjx/LEE09kPI6MjKRWrVrUrVs349y3337LokWLeOSRR/D39884v2nTJkCTJ0REREQKXfr4unTpEyjSu2O3bduGj49PxvkdO3YQFRWVoxt2woQJAAwYMCDb+SVLluBwOK4s2I0bN446derg6+tLaGgoy5Ytu+T1S5YsITQ0FF9fX+rWrZtRlIiIiEhJsGrVKtq2bZvxOOvM2NTUVMaPH8+jjz6a8fPcumGTk5NZvnw5NWvWzLFIcWRkJNdee63zwW769Ok89dRTvPzyy2zYsIGOHTvSu3dvDh48mOv1+/bto0+fPnTs2JENGzbw0ksv8cQTTzBjxgxn31pERESk2ElKSsKyLHx9fTPOZQ12EyZM4MEHH8TT0zPj57lNnDh58iQpKSm0aNEi2+ufP3+e33//nc6dOzsf7N5//30efPBBhg0bRpMmTRgzZgw1a9Zk/PjxuV4/YcIEatWqxZgxY2jSpAnDhg3jgQceYPTo0c6+tYiIiEixs2bNGlq3bp3tXEhICF5eXixZsoTTp09nLFqcbsmSJdSrV4+aNWtmnKtSpQr+/v5UrFgx27XTp08nMTGRDh06OLeOXWJiIuvXr+fFF1/Mdr5nz56sXLky1+esWrWKnj17ZjvXq1cvPv30U5KSkihVqlSO5yQkJJCQkJDxODY21pkyM1iWRVxc3BU9V0RERCRd2bJlM2agOuPcuXOMHDkyx5ZgXl5ehISEEBsby7PPPpvtZ3/88QfHjx+nd+/e2c47HA6GDRtGZGQklmXhcDiIiIjgu+++A6BatWpgOeHIkSMWYK1YsSLb+bfffttq2LBhrs9p0KCB9fbbb2c7t2LFCguwjh49mutzXnvtNQvIccTExDhTrhUTE5Pr6+jQoUOHDh06dDhzOJtBzp8/b7Vo0cLy8/OzAMvDw8OqU6eOtWPHjoxrBg8ebEVGRmY8fu6556ymTZtapUuXtgDLz8/PatasmTV69OiMa86dO2cNGTLE6tu3rzVs2DDrvffes1JSUqxXXnnF6tatm3VFO09cnFittNTozPW5nU8XFhbGM888k/E4NjY2W1NkfpUtW5aYmJjLXpf++ocOHSIgIMDp90nXunVr1q5de8XPd6fX0Gda8K/hbp+pq9Tibp+rK7yGPtOCfw19pgX/Gs5+pmXLlnXq9X19fTP2ec3LtGnTsj0eNWoUo0aNuuRz/Pz8+Oqrr3Kcf+uttwAntxSrVKkSnp6eHDt2LNv56OhogoKCcn1O1apVc73ey8srRx9xOh8fH3x8fJwpLVcOh8Op/wACAgKu6j8YT0/Pq3q+u70G6DMt6NcA9/lMXa0Wd/lcXeU1QJ9pQb8G6DMt6NeAq/9MXY1Tkye8vb0JDQ0lIiIi2/mIiIiMLTIu1q5duxzXL1iwgFatWuU6vq44yzpNWa9RMFzld3GV1ygIrvS7uFItV8tVfhdXeY2C4Cq/i6u8RkFwld/FVV7DHTms9H7RfJo+fTpDhgxhwoQJtGvXjkmTJjF58mS2bt1KSEgIYWFhHDlyJKOZcN++fTRr1oyHH36Y4cOHs2rVKkaMGMG0adO444478vWesbGxBAYGEhMTUyipurBfvyTSZ1rw9JkWDn2uBU+facHTZ1rw3PUzdXqM3aBBgzh58iRvvvkmUVFRNGvWjHnz5hESEgJAVFRUtjXt6tSpw7x583j66acJDw+nWrVqfPjhh/kOdUXBx8eH1157rUC6f8XQZ1rw9JkWDn2uBU+facHTZ1rw3PUzdbrFzg7umqpFRERECpL2ihURERFxE8Wixc5KW2j4ShcHFBERESkJikWwExEREZHLU1esiIiIiJtQsBMRERFxE24f7JYuXUq/fv2oVq0aDoeDH3/88bLPWbJkCaGhofj6+lK3bl0mTJhQ+IUWM85+rpGRkTgcjhzHX3/9VTQFu7iRI0fSunVrypYtS5UqVRgwYAA7duy47PP0Xb20K/lc9V29tPHjx3PddddlrNbfrl075s+ff8nn6Ht6ac5+pvqOOm/kyJE4HA6eeuqpS17nDt9Vtw92586d4/rrr+fjjz/O1/X79u2jT58+dOzYkQ0bNvDSSy/xxBNPMGPGjEKutHhx9nNNt2PHDqKiojKOBg0aFFKFxcuSJUt49NFHWb16NRERESQnJ9OzZ0/OnTuX53P0Xb28K/lc0+m7mrsaNWrwzjvvsG7dOtatW0fXrl259dZb2bp1a67X63t6ec5+pun0Hc2ftWvXMmnSJK677rpLXuc231WrBAGsWbNmXfKa559/3mrcuHG2cw8//LDVtm3bQqyseMvP57p48WILsE6fPl0kNRV30dHRFmAtWbIkz2v0XXVefj5XfVedV758eeuTTz7J9Wf6nl6ZS32m+o7mX1xcnNWgQQMrIiLC6tSpk/Xkk0/mea27fFfdvsXOWatWraJnz57ZzvXq1Yt169aRlJRkU1Xuo0WLFgQHB9OtWzcWL15sdzkuKyYmBoAKFSrkeY2+q87Lz+eaTt/Vy0tJSeHbb7/l3LlztGvXLtdr9D11Tn4+03T6jl7eo48+St++fenevftlr3WX76rTW4q5u2PHjhEUFJTtXFBQEMnJyZw4cYLg4GCbKivegoODmTRpEqGhoSQkJPD111/TrVs3IiMjuemmm+wuz6VYlsUzzzxDhw4daNasWZ7X6bvqnPx+rvquXt6WLVto164dFy5coEyZMsyaNYumTZvmeq2+p/njzGeq72j+fPvtt/zxxx+sXbs2X9e7y3dVwS4XFy+CbKUt9afFka9co0aNaNSoUcbjdu3acejQIUaPHq3/I7rIY489xubNm1m+fPllr9V3Nf/y+7nqu3p5jRo1YuPGjZw5c4YZM2YwdOhQlixZkmcQ0ff08pz5TPUdvbxDhw7x5JNPsmDBAnx9ffP9PHf4rqor9iJVq1bl2LFj2c5FR0fj5eVFxYoVbarKPbVt25Zdu3bZXYZLefzxx5kzZw6LFy+mRo0al7xW39X8c+ZzzY2+q9l5e3tTv359WrVqxciRI7n++usZO3Zsrtfqe5o/znymudF3NLv169cTHR1NaGgoXl5eeHl5sWTJEj788EO8vLxISUnJ8Rx3+a6qxe4i7dq146effsp2bsGCBbRq1YpSpUrZVJV72rBhQ7Fp2i5slmXx+OOPM2vWLCIjI6lTp85ln6Pv6uVdyeeaG31XL82yLBISEnL9mb6nV+ZSn2lu9B3Nrlu3bmzZsiXbufvvv5/GjRvzwgsv4OnpmeM5bvNdtWvWRlGJi4uzNmzYYG3YsMECrPfff9/asGGDdeDAAcuyLOvFF1+0hgwZknH93r17LT8/P+vpp5+2tm3bZn366adWqVKlrB9++MGuX8ElOfu5fvDBB9asWbOsnTt3Wn/++af14osvWoA1Y8YMu34Fl/J///d/VmBgoBUZGWlFRUVlHPHx8RnX6LvqvCv5XPVdvbSwsDBr6dKl1r59+6zNmzdbL730kuXh4WEtWLDAsix9T6+Es5+pvqNX5uJZse76XXX7YJc+LfziY+jQoZZlWdbQoUOtTp06ZXtOZGSk1aJFC8vb29uqXbu2NX78+KIv3MU5+7m+++67Vr169SxfX1+rfPnyVocOHay5c+faU7wLyu2zBKzPP/884xp9V513JZ+rvquX9sADD1ghISGWt7e3VblyZatbt24ZAcSy9D29Es5+pvqOXpmLg527flcdlpU2MlBEREREijVNnhARERFxEwp2IiIiIm5CwU5ERETETSjYiYiIiLgJBTsRERERN6FgJyIiIuImFOxERERE3ISCnYiIiIibULATERERcRMKdiIiIiJuQsFORERExE38Pw+F5EvBhQ3vAAAAAElFTkSuQmCC", "text/plain": [ "Graphics object consisting of 22 graphics primitives" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "from sage.plot.colors import rainbow\n", "clr=rainbow(8);\n", "#Shock Wave 4\n", "delta,Ma1,Ma2,gamma,beta,pp,sigma=SR.var('delta Ma1 Ma2 gamma beta pp sigma')\n", "eq1=tan(delta)==(Ma1^2*sin(beta)^2-1)/(tan(beta)*(Ma1^2*((gamma+1)/2-sin(beta)^2)+1))\n", "eq1=eq1.subs(gamma=1.4)\n", "eq4=sigma==(2*gamma/(gamma+1)*Ma1^2*sin(beta)^2-(gamma-1)/(gamma+1))^(-1/(gamma-1))*((gamma+1)*Ma1^2*sin(beta)^2/((gamma-1)*Ma1^2*sin(beta)^2+2))^(gamma/(gamma-1))\n", "eq4=eq4.subs(gamma=1.4)\n", "sol4=solve([eq1,eq4],[Ma1,sigma],solution_dict=True)\n", "sol42=solve(eq4.subs(beta=pi/2),sigma,solution_dict=True)#another solution\n", "\n", "sud=0;bb=0.05\n", "fg4=parametric_plot((sol4[1][Ma1].subs(delta=sud),sol4[1][sigma].subs(delta=sud)),(beta,0.1,pi/4-bb),color=clr[0],\n", " axes_labels=(r'$Ma_1$',r'$\\sigma$'),title=r'$\\sigma\\ vs\\ Ma_1$')\n", "# fg4+=parametric_plot((sol4[4][Ma1].subs(delta=sud),sol4[4][Ma2].subs(delta=sud)),(beta,0.1,pi/4-bb),color=clr[0])\n", "fg4+=plot(sol42[0][sigma],(Ma1,1,4),color=clr[0])\n", "fg4+=text(r'$\\delta=0^\\circ$',(1.1,0.9),color=clr[0])\n", "for i in range(1,8):\n", " sud=i*pi/36\n", " fg4+=parametric_plot((sol4[1][Ma1].subs(delta=sud),sol4[1][sigma].subs(delta=sud)),(beta,0.1,pi/4-bb),color=clr[i])\n", " fg4+=parametric_plot((sol4[0][Ma1].subs(delta=sud),sol4[0][sigma].subs(delta=sud)),(beta,0.1,pi/4-bb),color=clr[i])\n", " tx=sol4[1][Ma1].subs(delta=sud,beta=pi/2.4)+0.05;\n", " ty=sol4[1][sigma].subs(delta=sud,beta=pi/2.4)+0.02;\n", " delt='$'+str(i*5)+'^\\\\circ$'\n", " fg4+=text(r'%s'%delt,(tx,ty),color=clr[i])\n", "end\n", "fg4.show(xmin=1,xmax=4,ymin=0,ymax=1,aspect_ratio=2)\n", "# fg4.save_image('4_2d.png',xmin=1,xmax=4,ymin=0,ymax=1,aspect_ratio=2)\n" ] }, { "cell_type": "code", "execution_count": 56, "metadata": { "vscode": { "languageId": "python" } }, "outputs": [], "source": [ "import numpy as np\n", "delta,Ma1,Ma2,gamma,beta,pp,sigma=SR.var('delta Ma1 Ma2 gamma beta pp sigma')\n", "eq1=tan(delta)==(Ma1^2*sin(beta)^2-1)/(tan(beta)*(Ma1^2*((gamma+1)/2-sin(beta)^2)+1))\n", "eq1=eq1.subs(gamma=1.4)\n", "eq4=sigma==(2*gamma/(gamma+1)*Ma1^2*sin(beta)^2-(gamma-1)/(gamma+1))^(-1/(gamma-1))*((gamma+1)*Ma1^2*sin(beta)^2/((gamma-1)*Ma1^2*sin(beta)^2+2))^(gamma/(gamma-1))\n", "eq4=eq4.subs(gamma=1.4)\n", "sol4=solve([eq1,eq4],[Ma1,sigma],\n", " #solution_dict=True\n", " )\n", "sol4=sol4[1]\n", "sol4_2=solve([sol4[1],eq1],[Ma1,delta],\n", " # solution_dict=True\n", " )\n", "\n", "\n", "# iterate Delta\n", "# sol4_2=sol4_2\n", "# sol4_tdel=sol4_2[0].subs(delta==pi/8)\n", "# sol4_sig=sol4_2[1].subs(delta==pi/8)\n", "# betas=np.arange(0,pi/2,0.1)\n", "# #iterate beta\n", "# #ma1 -> sigma\n", "# print(sol4_tdel)\n", "# sol4_bet_itr=find_root(sol4_tdel.subs(beta=betas[4]),1,4)\n", "\n", "# print(sol4_bet_itr)\n", "\n", "\n", "# solsig=sol4_2[1][sigma].subs(sol4_2[0]) #sigma\n", "\n", "# solsig_itr=solsig.subs(beta==pi/8)\n", "# solsig_itr=solsig_itr.subs(Ma1==Ma1s[1]).n()\n", "# sol1m=solve(sol1.subs(beta==1.2),Ma1,solution_dict=True)\n", "# sol1m[0][Ma1].subs(delta=0.2)\n", "# dmdb=derivative(sol1,beta)\n", "# sol1d=solve(dmdb.subs(delta=0.4)==0,beta)\n", "# find_root(sol1,1,pi/2-0.05)" ] }, { "cell_type": "code", "execution_count": 34, "metadata": { "vscode": { "languageId": "python" } }, "outputs": [ { "data": { "text/plain": [ "array([1. , 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2. , 2.1, 2.2,\n", " 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9, 3. , 3.1, 3.2, 3.3, 3.4, 3.5,\n", " 3.6, 3.7, 3.8, 3.9])" ] }, "execution_count": 34, "metadata": {}, "output_type": "execute_result" } ], "source": [] } ], "metadata": { "kernelspec": { "display_name": "SageMath 9.5", "language": "sage", "name": "sagemath" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "sage", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.2" }, "orig_nbformat": 4, "vscode": { "interpreter": { "hash": "7f31cd35f582712ee55286e4e453a72bdfcbe905b252282b8b04a3c94a3089b6" } } }, "nbformat": 4, "nbformat_minor": 2 }