diff options
Add fixed velocity boundaries to generated LBM kernel
Interestingly this increased performance to ~750 MLUPS compared to ~665 MLUPS.
| -rw-r--r-- | codegen_lbm.py | 84 | ||||
| -rw-r--r-- | lbm_codegen.ipynb | 1776 |
2 files changed, 419 insertions, 1441 deletions
diff --git a/codegen_lbm.py b/codegen_lbm.py index cd93649..8e7b0a7 100644 --- a/codegen_lbm.py +++ b/codegen_lbm.py @@ -51,50 +51,48 @@ __kernel void collide_and_stream(__global __write_only float* f_a, const float f_curr_7 = f_i(f_b, cell.x , cell.y+1, 0,-1); const float f_curr_8 = f_i(f_b, cell.x-1, cell.y+1, 1,-1); + const float ux0 = f_curr_3 + f_curr_6; + const float ux1 = f_curr_1 + f_curr_2; + const float ux2 = 1.0/(f_curr_0 + f_curr_4 + f_curr_5 + f_curr_7 + f_curr_8 + ux0 + ux1); + const float ux3 = f_curr_0 - f_curr_8; + + float u_x = -ux2*(-f_curr_2 - f_curr_5 + ux0 + ux3); + float u_y = ux2*(-f_curr_6 - f_curr_7 + ux1 + ux3); + + if ( m == 2 ) { + u_x = 0.0; + u_y = 0.0; + } + const float x0 = f_curr_0 + f_curr_1 + f_curr_2 + f_curr_3 + f_curr_4 + f_curr_5 + f_curr_6 + f_curr_7 + f_curr_8; - const float x1 = 2*f_curr_0; - const float x2 = 2*f_curr_8; - const float x3 = -f_curr_3 + f_curr_5; - const float x4 = pow(x0, -2); - const float x5 = 9*x4; - const float x6 = f_curr_0 - f_curr_8; - const float x7 = f_curr_1 - f_curr_7; - const float x8 = f_curr_2 - f_curr_6; - const float x9 = x6 + x7 + x8; - const float x10 = 6/x0; - const float x11 = x10*x9; - const float x12 = f_curr_3 - f_curr_5; - const float x13 = -f_curr_2 + f_curr_6 + x12 + x6; - const float x14 = pow(x13, 2); - const float x15 = 3*x4; - const float x16 = -x14*x15 + 2; - const float x17 = x11 + x16; - const float x18 = pow(x9, 2); - const float x19 = x15*x18; - const float x20 = -x19; - const float x21 = x10*x13; - const float x22 = x20 + x21; - const float x23 = 1.0/$tau; - const float x24 = (1.0/72.0)*x23; - const float x25 = 6*x4; - const float x26 = x18*x25; - const float x27 = (1.0/18.0)*x23; - const float x28 = x5*pow(2*f_curr_2 - 2*f_curr_6 + x3 + x7, 2); - const float x29 = x20 - x21; - const float x30 = x14*x25 + 2; - const float x31 = -f_curr_0 + f_curr_8 + x3 + x8; - const float x32 = x15*pow(x31, 2) - 2; - const float x33 = x19 + x32; - - f_a[0*$nCells + gid] = f_curr_0 - x24*(72*f_curr_0 - x0*(x17 + x22 + x5*pow(-f_curr_1 + f_curr_7 - x1 + x2 + x3, 2))); - f_a[1*$nCells + gid] = f_curr_1 - x27*(18*f_curr_1 - x0*(x17 + x26)); - f_a[2*$nCells + gid] = f_curr_2 - x24*(72*f_curr_2 - x0*(x17 + x28 + x29)); - f_a[3*$nCells + gid] = f_curr_3 - x27*(18*f_curr_3 - x0*(x22 + x30)); - f_a[4*$nCells + gid] = f_curr_4 - 1.0/9.0*x23*(9*f_curr_4 + 2*x0*x33); - f_a[5*$nCells + gid] = f_curr_5 - x27*(18*f_curr_5 - x0*(x29 + x30)); - f_a[6*$nCells + gid] = f_curr_6 - x24*(72*f_curr_6 + x0*(x10*x31 + x11 - x28 + x33)); - f_a[7*$nCells + gid] = f_curr_7 - x27*(18*f_curr_7 + x0*(x11 - x26 + x32)); - f_a[8*$nCells + gid] = f_curr_8 - x24*(72*f_curr_8 - x0*(-x11 + x16 + x29 + x5*pow(x1 + x12 - x2 + x7, 2))); + const float x1 = 6*u_y; + const float x2 = 6*u_x; + const float x3 = pow(u_y, 2); + const float x4 = 3*x3; + const float x5 = pow(u_x, 2); + const float x6 = 3*x5; + const float x7 = x6 - 2; + const float x8 = x4 + x7; + const float x9 = x2 + x8; + const float x10 = 1.0/$tau; + const float x11 = (1.0/72.0)*x10; + const float x12 = 6*x3; + const float x13 = x1 - x6 + 2; + const float x14 = (1.0/18.0)*x10; + const float x15 = -x4; + const float x16 = 9*pow(u_x + u_y, 2); + const float x17 = -x2; + const float x18 = x15 + 6*x5 + 2; + + f_a[0*$nCells + gid] = f_curr_0 - x11*(72*f_curr_0 + x0*(-x1 + x9 - 9*pow(-u_x + u_y, 2))); + f_a[1*$nCells + gid] = f_curr_1 - x14*(18*f_curr_1 - x0*(x12 + x13)); + f_a[2*$nCells + gid] = f_curr_2 - x11*(72*f_curr_2 - x0*(x13 + x15 + x16 + x2)); + f_a[3*$nCells + gid] = f_curr_3 - x14*(18*f_curr_3 - x0*(x17 + x18)); + f_a[4*$nCells + gid] = f_curr_4 - 1.0/9.0*x10*(9*f_curr_4 + 2*x0*x8); + f_a[5*$nCells + gid] = f_curr_5 - x14*(18*f_curr_5 - x0*(x18 + x2)); + f_a[6*$nCells + gid] = f_curr_6 - x11*(72*f_curr_6 + x0*(x1 - x16 + x9)); + f_a[7*$nCells + gid] = f_curr_7 - x14*(18*f_curr_7 + x0*(x1 - x12 + x7)); + f_a[8*$nCells + gid] = f_curr_8 - x11*(72*f_curr_8 + x0*(x1 + x17 + x8 - 9*pow(u_x - u_y, 2))); moments[gid] = x0; }""" diff --git a/lbm_codegen.ipynb b/lbm_codegen.ipynb index 7d13f70..34e9e6b 100644 --- a/lbm_codegen.ipynb +++ b/lbm_codegen.ipynb @@ -14,28 +14,30 @@ }, { "cell_type": "code", - "execution_count": 30, + "execution_count": 2, "metadata": {}, "outputs": [], "source": [ - "c = [Matrix((x,)) for x in [(-1, 1), ( 0, 1), ( 1, 1), (-1, 0), ( 0, 0), ( 1, 0), (-1,-1), ( 0, -1), ( 1, -1)]]" + "c = [Matrix(x) for x in [(-1, 1), ( 0, 1), ( 1, 1), (-1, 0), ( 0, 0), ( 1, 0), (-1,-1), ( 0, -1), ( 1, -1)]]" ] }, { "cell_type": "code", - "execution_count": 31, + "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/latex": [ - "$$\\left [ \\left[\\begin{matrix}-1 & 1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}0 & 1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}1 & 1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}-1 & 0\\end{matrix}\\right], \\quad \\left[\\begin{matrix}0 & 0\\end{matrix}\\right], \\quad \\left[\\begin{matrix}1 & 0\\end{matrix}\\right], \\quad \\left[\\begin{matrix}-1 & -1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}0 & -1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}1 & -1\\end{matrix}\\right]\\right ]$$" + "$$\\left [ \\left[\\begin{matrix}-1\\\\1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}0\\\\1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}1\\\\1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}-1\\\\0\\end{matrix}\\right], \\quad \\left[\\begin{matrix}0\\\\0\\end{matrix}\\right], \\quad \\left[\\begin{matrix}1\\\\0\\end{matrix}\\right], \\quad \\left[\\begin{matrix}-1\\\\-1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}0\\\\-1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}1\\\\-1\\end{matrix}\\right]\\right ]$$" ], "text/plain": [ - "[[-1 1], [0 1], [1 1], [-1 0], [0 0], [1 0], [-1 -1], [0 -1], [1 -1]]" + "⎡⎡-1⎤ ⎡0⎤ ⎡1⎤ ⎡-1⎤ ⎡0⎤ ⎡1⎤ ⎡-1⎤ ⎡0 ⎤ ⎡1 ⎤⎤\n", + "⎢⎢ ⎥, ⎢ ⎥, ⎢ ⎥, ⎢ ⎥, ⎢ ⎥, ⎢ ⎥, ⎢ ⎥, ⎢ ⎥, ⎢ ⎥⎥\n", + "⎣⎣1 ⎦ ⎣1⎦ ⎣1⎦ ⎣0 ⎦ ⎣0⎦ ⎣0⎦ ⎣-1⎦ ⎣-1⎦ ⎣-1⎦⎦" ] }, - "execution_count": 31, + "execution_count": 3, "metadata": {}, "output_type": "execute_result" } @@ -46,7 +48,7 @@ }, { "cell_type": "code", - "execution_count": 32, + "execution_count": 4, "metadata": {}, "outputs": [], "source": [ @@ -55,7 +57,7 @@ }, { "cell_type": "code", - "execution_count": 33, + "execution_count": 5, "metadata": {}, "outputs": [ { @@ -68,7 +70,7 @@ "[1/36, 1/9, 1/36, 1/9, 4/9, 1/9, 1/36, 1/9, 1/36]" ] }, - "execution_count": 33, + "execution_count": 5, "metadata": {}, "output_type": "execute_result" } @@ -79,7 +81,7 @@ }, { "cell_type": "code", - "execution_count": 34, + "execution_count": 6, "metadata": {}, "outputs": [ { @@ -92,7 +94,7 @@ "1" ] }, - "execution_count": 34, + "execution_count": 6, "metadata": {}, "output_type": "execute_result" } @@ -103,7 +105,7 @@ }, { "cell_type": "code", - "execution_count": 35, + "execution_count": 7, "metadata": {}, "outputs": [ { @@ -118,7 +120,7 @@ "3 " ] }, - "execution_count": 35, + "execution_count": 7, "metadata": {}, "output_type": "execute_result" } @@ -130,7 +132,7 @@ }, { "cell_type": "code", - "execution_count": 36, + "execution_count": 22, "metadata": {}, "outputs": [], "source": [ @@ -139,7 +141,7 @@ }, { "cell_type": "code", - "execution_count": 37, + "execution_count": 23, "metadata": {}, "outputs": [ { @@ -149,7 +151,7 @@ " f_next_6, f_next_7, f_next_8], dtype=object)" ] }, - "execution_count": 37, + "execution_count": 23, "metadata": {}, "output_type": "execute_result" } @@ -161,7 +163,7 @@ }, { "cell_type": "code", - "execution_count": 38, + "execution_count": 24, "metadata": {}, "outputs": [ { @@ -171,7 +173,7 @@ " f_curr_6, f_curr_7, f_curr_8], dtype=object)" ] }, - "execution_count": 38, + "execution_count": 24, "metadata": {}, "output_type": "execute_result" } @@ -183,7 +185,7 @@ }, { "cell_type": "code", - "execution_count": 39, + "execution_count": 25, "metadata": {}, "outputs": [ { @@ -197,7 +199,7 @@ "_curr_7 + f_curr_8" ] }, - "execution_count": 39, + "execution_count": 25, "metadata": {}, "output_type": "execute_result" } @@ -209,95 +211,41 @@ }, { "cell_type": "code", - "execution_count": 40, + "execution_count": 55, "metadata": {}, - "outputs": [ - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaUAAAAlCAYAAADiFfbCAAAABHNCSVQICAgIfAhkiAAACCJJREFUeJztnVuM1FQYx3/LIoKg3BRCZA1e42VXN4KioogEifeAGHnwZY0YEjEYHwzGROXFqC9qDPEaSdUoRAWioJiIKBqVi2YBkUQQXQRUBEHAIAq4PnynbLd0mE57Oj2d/X5Js3O6nfb/9fs6p+2Z6R8URVEURak5pgA/AoeA53PWkpSix6D686foMRRdv6IAcC5SxJOAIUDvfOUkougxqP78KXoMRdevKEd4APgmbxEhPGBmBcu7GEMluKjfI34OXNRfKS7G4NG1cqAobADaA9O8fOUcwSP+wehqDHFxVb9HvBy4qr8SXI3Bo+vkoCbonreAGuBK4HPgNeAVYH++chJR9BhUf/4UPYai668ZuuUtoAbYC5wBfAH8Ztq2uQH4HtgI3JPB+rOMoQH4FFgPrAFutbhun2rk4EtE/zrgEcvrrob+NmAtsBr4JIP1Zx3DMGApUkffASdbXn81cjANqZ91SMdXn8E2lBplJp0v5aOmMWbZEabdz9K2w4V6PNIZNQAnIAflkIj3PQT8FZgOAv+G5l1VYps2YwjrPxW4yLweBGxB4ohiJvH3e5Csc1APnBR4/RXQHPG+pDmohv42oE+Z980k2f6HbGuoHlgGjDbtvshxEYWrORgMbAJ6AnXAQuAmS9uqKfT2XTSzgLlllvnZ/G0GNgN/Bv7XYNZxGtADKb6fgFZgPLADGI4cQJOQAl2PHPCzzfJ++wPzeotZ9wLz/5dDel4A3gq0nwS2Ac8G5m0rEUvcGOYn0D8beNGs83dgN3KW6++/IJXsd5v6qSCGHmaKImkOsq6h2SX0hkm6/yuJIUkNvYl0Lp+Z9e45hj5XczAf+bztCRwGegHbjxFHl0U7pWh2mikOzcgtEZ8eSEcyHblN0hc4gOzr/kghAzQhl/EAjcgtrpGm/WCgfRtwSmD9W5GrjzC7zOSzz7R/sBhDEv1BRgDH0dHBhqlkv9vUD/FiWAGcDzwX2p5P0hxkXUN+exnwH/AM8EaEjqT7v5IYkuifgOzLd5EOYiGlb6G6nIP9SKd+EJgDrCqjqUuiY0rpCRfzRGAlHfft9wD/AOfQ+cC4ECnmE5ED4CkzP9yui9hmuw3hAeLEcCbJ9PsMRAaR78It/RA/hpHICUEz8gFki6xrCGAUclZ/C3Jm32RRP2RbQ92Rq43pwKVIHBNz0J8mB/2BG5GxsaFI/YyxG0JtoJ1SOuqQogwWcxNSzGEa6XxmPsK0G5FB9PbAcsH2NuQ2gs9Q4Je0wgPEjSGpfpD7/wuAx83/bJJWv/+/cjH47EUG3K9LpbqDatQQdNTMr8gVwPC0wgNkXUNbkd8PbUauMhYRPaaXlGrkYBwyprQL+Bt4H+lglRB6+y4d7XQMgPtsp+MsuhtyhvQHMAApRoDLkA+FjcCdyLeifBpD7ZXABUjHtBM5gxsXQ1uL5RiS6q9DfiuyFHg9pqZKSKsfjtYcbvdFznp3IGMC44GnY2hrsag/TQ31NuvZh3zZYSydx13SknUNrUKutAciH+pXIyc5cWixqD9NDrYAlyP1cxC5SnopXghdC71Sso+H3KZYh5zdnWfmL0YGPucgX/HeiDzSpInOxRtuHwLuBz4GvkWex2XzSikKj6NjSKp/FDAZGRdYbSbbt47CeMTXH6U53O4PfGjmfY2MzSzKTL39GhqM/AZnDbAcuY2a9XiGh70aOgzMQG6lrUVOzmx2qlF42M3BcvPeVjN/E/BeZuoVRVEURVEURVEURVEURVEURVEUpYbxfwNj+3cjiqIoiqIoiqIotUXR7YBVf/4UPQbVnz9Fj6Ho+sGRGIpuB6z686foMaj+/Cl6DEXXDw7F4JodsEfxLb09im3H7KE5yBuPrqMfih9D0fVDjjEEHzO0ATjbvG5HHrU+6ah3uIvqz5+ix6D686foMRRdP+QcQ7BTKrodsOrPn6LHoPrzp+gxFF0/5BxD8Nl3WdsBF9nSG+Q5VbuBdyyv1ydr/bVgS56lJTlUxxIb5Lhbhf1aUlv18gyj2LbqTXQ8Q3I18nDYCZa3kXUMsW3hs7QzLrqldz1wDXAz5T9IXLVjjmtL7nIO4liSg7s58NvTEIfXUrXksv42ytuqu1xDRbdVD7b7IA+njfoSgqs5qMgWfgpScD4NiNNjK3JGcbqZ30qHE+pwYJ55vRCxHl4BTA21HzXr8nkMuDtCwwDgrMA0D7ExDs7rFVN/qRiS6J9qlhlD+U4paQxx9ZMyBp+1iIunLf2VxJBWfy9kELZUp5R1DkrpjxPDIOSJ72MpXUsu11Ab5TslV2voPmBJGe1pY6jmcXwHpe3rXc3BDMQXqx/iQr0EuCS4oeCYUtB5US29S1t6l8NlO2afY9mSu56DcpbkaWKohq36q8DDZt2lcLmG4tiqu1pDtWKr7nM78q06m/oriSETW/jgmFJQiFp6R9thZ0k1LLHBXVvyvC3JIXtb9dHIPrftvuujtuq1b6vu0xcxDVxsUbtPljkoawvvd0phO2C19C5th50F1bLEdtmWPE9LcqiOrfoVwLXIrZG5wPXIQK8N1Fa9a9iq+0xAjCcPpFbdmaxzUNYW3r99F7YDVkvvzu00tMRYphp2zEltyVtiLpd1DpJakoPdHKSxVX/CTCBnh/ciV6zlaImxjMu26nH0VxJDV7ZV95kMzIqhyacl5nJZ56CsLXwpO3QPtfQOF8IS4G3zvq3Ijs0SD7s5cN2W3DVLcrBvq15tPNRWvdZs1UE6hIuBj7ISHcLDXg7UFl5RFEVRFEVRFEVRFEVRFEVRFEVJzf8pYN+WylHVPgAAAABJRU5ErkJggg==\n", - "text/latex": [ - "$$\\frac{- f_{curr 0} + f_{curr 2} - f_{curr 3} + f_{curr 5} - f_{curr 6} + f_{curr 8}}{f_{curr 0} + f_{curr 1} + f_{curr 2} + f_{curr 3} + f_{curr 4} + f_{curr 5} + f_{curr 6} + f_{curr 7} + f_{curr 8}}$$" - ], - "text/plain": [ - " -f_curr_0 + f_curr_2 - f_curr_3 + f_curr_5 - f_curr_6 + f_curr\n", - "──────────────────────────────────────────────────────────────────────────────\n", - "f_curr_0 + f_curr_1 + f_curr_2 + f_curr_3 + f_curr_4 + f_curr_5 + f_curr_6 + f\n", - "\n", - "_8 \n", - "──────────────────\n", - "_curr_7 + f_curr_8" - ] - }, - "execution_count": 40, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ - "u_x = sum([ (c_i*f_curr[i])[0] for i, c_i in enumerate(c) ]) / rho\n", - "u_x" + "#u_x = sum([ (c_i*f_curr[i])[0] for i, c_i in enumerate(c) ]) / rho\n", + "#u_x" ] }, { "cell_type": "code", - "execution_count": 41, + "execution_count": 56, "metadata": {}, - "outputs": [ - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaUAAAAlCAYAAADiFfbCAAAABHNCSVQICAgIfAhkiAAAB3pJREFUeJztnWuIFWUYx3+blpqlrZomaXSny25JWlaWbWLRnbWiPvRlI0PIMPoQSVHtl6i+VEQUBMLJCKWyMC2DSovIUgsvWZBmrK3dcDPSyC7W9uF5hzM7zpwzl/edi/v8YOC8s+fM+T/7PDPvO/POmT8oiqIoilJp5gLfAgeA5wvWkgbVXzxVj0H1K0pJOAMp5JuAicDIYuUkRvUXT9VjUP2KUiLuAz4vWoSPGtCd4P1l0w/JYiij/qRUPQbVryglYRvQ71uWFSsHSHZAL6N+iB9DWfUnoeoxqH5FKRHjga+BB4HjgFHFygGSdUpl1A/xYyir/iRUPQbVrzjlsKIFVIy9wMnAx8BPpm2Ta5AdZjtwl+Vtg3v9bwK/Aq9Z3q6Ha/2TgQ+Ar4DNwI2Wtw/uYzgRWI3E8CUwzvL2XetvBzb5lv1Ap8Xtu9YPMB/YapZFwBAH36EoAExDTvmPsbCtYKEOQzqjycCRyEFlYsjnHgB+9y3/AH8H1l0a8Z0u9Q8BLgeup3mnlDYGm/rh4BiOB841r8cDvUgugnQz8BJQ2NIR8Z2uc/AhMNO0RyN1FUY36WJwnQN/+yigj/AbEbopp/4JwA5gONACrACus/RdinIQc4EeX3sysBzYiIxKTzLrNwLHmtdTqV+3XgE8AawD5gXaj5hteTwK3BmiYQxwqm9ZBjwTWDcipv6oGNLon2fe00HzTiltDHH1kzEGjy3ACSE6xiF3cDVawjqzJDGk0X8P8F7E99qKIav+ZjH4c3AbsLQg/TSIoZH++4GdSKd3OJKP8yNiUEIYWrSAijEFuaQAcATwNrAAWIOMSv9E/qetwG7zvnbkNB6gDbk8NN20F/raN1PfAQB2ISP3IHvM4rHPtL9JqL9RDGn0JyFtDHH1p82Bn2nIQaU3REefWdLgMgedyP9yOdKZrgAejtCRNoas+pvF4OcWZL4xDNf6s9TQH8B3yBWAJcCGFDoHLTqnlAx/Qc8B1iOFDPAb8BdwOgMPrucgxXw0sgM8adYH2y0h39dvS7ghuEOGxXAK6fTnQRz9WXLgMRZYDNxBtXIwFDlTXQBcgIzu55RIf5jmqByMBi4CVtkSbnBdQ63Atcjc3iSkA+uwG8KhjXZK8WlBCtMr6HakmIO0MXBUOM2024C11A9ywfb3yGUEj0nADzaEG4L6ITyGtPpdE1c/ZIthGPAG8Jj5m01c52AX8vubncgofSVyELZFVv1hmqPqqBN4BzlrsUUeNTQbmVPag9yk8RYyQFBiopfv4tPPwNtHf0YKEqRzbwV+QeZL9pv1FyKj1e3A7cgchUdboL0eOBvpmPqQEdzsGLq6UuqPiiGt/ix0xXhPXP1ZctCCXC5aDbwUV3wCXOdgA3KWNxY5KF6GdLC2yKo/THNUHd0KPJtd8gDyqKFe5AxvODIw6ABesBXAYEDPlNJTQy5TbEVGp2ea9auQu22WILd4b0ceadLOwOINtg8A9wLvA18gz+OyeaYURo2DY0irH2RS91XzuV3IzumSGnZzMAM5GHZSvyW53Zl6oYa9HPyLTLSvMev7gFecqk+mP0xzWB2NAc4D3nUl2kcNuzX0qfnsRrN+B/JTCUVRFEVRFEVRFEVRFEVRFEVRFEVRLOH9Niav23oVRVEURVEURVGqSdVtglV/8VQ9BtVfPFWPoer6oSQxVN0mWPUXT9VjUP3FU/UYqq4fShRD2WyCa6jVd9HU0BwUTY3Box+qH0PV9UOBMfgfM7QNOM287gdeR3rJqqD6i6fqMaj+4ql6DFXXDwXH4O+ULgE+Qp6OvAh5/HqVUP3FU/UYVH/xVD2GquuHgmPwP/tOrb4bo1bfzXEdw1pE+1aifYKykIdVNsh+twH7tZSH/h7kmW6bqNs92ETt4hvj2i4e3McQ2y5erb7rDGar7zLnYJTv9SdE2zKUNQdeez7iqBpVS2XW34PYlDeizDUU1y6+zDnwaGQXX9YcJLKLV6tvtfrOoj9JDFn1j0AmYaM6Jdc5yGL1PR55EvwsomupzDXUQ/NOqaw1lMQuvsw58GhkF1/WHDS1i/fPKanVt1p9Z9GfJIYs+tcBZwHPBb7LRgx5WH2/CDxkth1FmWtoIXK28R/wNPCyRf1JYnBtF1/mHHg0sosvaw6giV28f05Jrb7V6jsreeRgOjKgmULdnM0Wrq2+ZyL/c9uOth551NAMZFR8A3KJyLbf1GCwi7exH7uyiwe3OWhqF+91Smr1rVbfWckzB3uRyeqrMquuk4fV98XAFcilkaXA1chErw3yqiFvv/0RGT1PzSrcx2Cxi/e0ZdkPXNjFg/scNLWL9y7fqdW3Wn1H0RXzfa5zMBoZce1GJkmvBJ6Kqa0rB/1hmoPtx80CMjq8GzljbUZXjPfkUUMjzXb2IfNKs4jnbBtHf5IYirCL77KoP0sOPJLaxXfFfJ/rHDS1i4+yQ6+hVt9q9Z2dGvZy0IqMDLcAnyHzGiudqndj9Z0nNezW0ATk9yubEdvvxQTmAxxQ49Czi4dsx6I87eLBbg7ULl5RFEVRFEVRFEVRFEVRFEVRFEXJzP/mgKdK6IgqEQAAAABJRU5ErkJggg==\n", - "text/latex": [ - "$$\\frac{f_{curr 0} + f_{curr 1} + f_{curr 2} - f_{curr 6} - f_{curr 7} - f_{curr 8}}{f_{curr 0} + f_{curr 1} + f_{curr 2} + f_{curr 3} + f_{curr 4} + f_{curr 5} + f_{curr 6} + f_{curr 7} + f_{curr 8}}$$" - ], - "text/plain": [ - " f_curr_0 + f_curr_1 + f_curr_2 - f_curr_6 - f_curr_7 - f_curr_\n", - "──────────────────────────────────────────────────────────────────────────────\n", - "f_curr_0 + f_curr_1 + f_curr_2 + f_curr_3 + f_curr_4 + f_curr_5 + f_curr_6 + f\n", - "\n", - "8 \n", - "──────────────────\n", - "_curr_7 + f_curr_8" - ] - }, - "execution_count": 41, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ - "u_y = sum([ (c_i*f_curr[i])[1] for i, c_i in enumerate(c) ]) / rho\n", - "u_y" + "#u_y = sum([ (c_i*f_curr[i])[1] for i, c_i in enumerate(c) ]) / rho\n", + "#u_y" ] }, { "cell_type": "code", - "execution_count": 42, + "execution_count": 26, "metadata": {}, "outputs": [ { "data": { "text/latex": [ - "$$\\left[\\begin{matrix}\\frac{- f_{curr 0} + f_{curr 2} - f_{curr 3} + f_{curr 5} - f_{curr 6} + f_{curr 8}}{f_{curr 0} + f_{curr 1} + f_{curr 2} + f_{curr 3} + f_{curr 4} + f_{curr 5} + f_{curr 6} + f_{curr 7} + f_{curr 8}}\\\\\\frac{f_{curr 0} + f_{curr 1} + f_{curr 2} - f_{curr 6} - f_{curr 7} - f_{curr 8}}{f_{curr 0} + f_{curr 1} + f_{curr 2} + f_{curr 3} + f_{curr 4} + f_{curr 5} + f_{curr 6} + f_{curr 7} + f_{curr 8}}\\end{matrix}\\right]$$" + "$$\\left[\\begin{matrix}u_{x}\\\\u_{y}\\end{matrix}\\right]$$" ], "text/plain": [ - "⎡ -f_curr_0 + f_curr_2 - f_curr_3 + f_curr_5 - f_curr_6 + f_cur\n", - "⎢─────────────────────────────────────────────────────────────────────────────\n", - "⎢f_curr_0 + f_curr_1 + f_curr_2 + f_curr_3 + f_curr_4 + f_curr_5 + f_curr_6 + \n", - "⎢ \n", - "⎢ f_curr_0 + f_curr_1 + f_curr_2 - f_curr_6 - f_curr_7 - f_curr\n", - "⎢─────────────────────────────────────────────────────────────────────────────\n", - "⎣f_curr_0 + f_curr_1 + f_curr_2 + f_curr_3 + f_curr_4 + f_curr_5 + f_curr_6 + \n", - "\n", - "r_8 ⎤\n", - "───────────────────⎥\n", - "f_curr_7 + f_curr_8⎥\n", - " ⎥\n", - "_8 ⎥\n", - "───────────────────⎥\n", - "f_curr_7 + f_curr_8⎦" + "⎡uₓ ⎤\n", + "⎢ ⎥\n", + "⎣u_y⎦" ] }, - "execution_count": 42, + "execution_count": 26, "metadata": {}, "output_type": "execute_result" } @@ -309,514 +257,103 @@ }, { "cell_type": "code", - "execution_count": 43, + "execution_count": 27, "metadata": {}, "outputs": [ { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAZU8AAAAqCAYAAADTjjwVAAAABHNCSVQICAgIfAhkiAAAIABJREFUeJzs3Xu0LVtd2PnvhXsVBEVAedhCrrb4aC/qULAFL2YrMUNJ2/Ic2sbHSsTWiJp02rSRNs0xIwTQmA4+SNqYpEQdRtpHG1REwLeCMhRBhKEgbuUtCgZUlNfpP2qtnHXWWY+qWjVr/n5V388YZ5xz9lp779+sWat+q+b8zblAkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJ6uEhwOuA3wRetP7z4A7fd4eSQUmSJEmSJElSRfcDfh54BfBS4NFVo5EkSZIkReS9oyRJmoubagcgSZIkSdKZXO8qSdrlfK4kSd2YMyVJUlbWPEmSJElaKutkJEmSJEmSJEVkTaokSZIk6RTvHSVJ0ly4vlGSJEmSlJ1rVSUpL+ddJUmSJEn7eL84b9YrSZIkSVoy61wkSdscC5UkaTjz6Pw5ryhJkiRJys75YUnKyzFoSZIkSdIu7xXnz3olSZIkSUtljYskSZIkSZKkObLeR5IkSZLUhfePkiRpLlwjKUmSJEnKzvWukqRtzuVKktSdeVOSJGVlzZMkSZKkpbJORpIkSZIkSVI01qNKkiRJkrrw/lGSJM2F6xslSZIkSV18GvCi9Z+XAK8DHnzoyRfAVeD2Hr/gkcBnDI9PkiRJkiRJkkK7L/CJ63/fC3gt8AH1wpEkSZIkBeS9oyRJmosnATfXDkKSJEmSpIFc7ypJ2sf5XEmSujFnSpKkrKx5kiRJkrRE1slIkiRJkiRJisqaVEmSJEnSKd47SpKkuXB9oyRJkiQpM9eqSlJuzrtKkiRJkvbxfnHerFeSJEmStFTWuUiSdjkWKknScObR+XNeUZIkSZKUmfPDkpSbY9CSJEmSpF3eK86f9UqSJEmSlsgaF0mSJEmSJElzZb2PJEmSJKkL7x8lSdJcuEZSkiRJkpSZ610lSbucy5UkqTvzpiRJysqaJ0mSJElLZJ2MJEmSJEmSpIisR5UkSZIkdeH9oyRJmgvXN0qSJEmS+roArgK3D37Cjs8FnnJuVCfcqfDPVz+PB14DvAf4t5VjGSp7G7LHD/nbYPzS/MzhdZG9DcZfX/Y2ZI9/7m4BbqodhFK5O/Bm4L8/8+dkvzZkjx/yt8H466vRhh8G/vGJ57wMuP8EsUhdmDdb2eOH/G0w/vqi5k0pimM5M/s1wPjry96G7PFDzDzovaMyM2/GlT1+yN+G7PFDzLypWB4EfE/tICRJkiRJGmCK9a7qxpqmVvb4IX8bsscP+dsQdUza+VxldSjHZr9WQP42GH992dtgzpTOZ56MK3v8kL8Nxj+MNU+5WPMkSZIkaWncFz4fx2jqyh4/5G9D9vhhHm2QxjSH10T2Nhh/fdnbkD3+uXNPd/XlOo5W9vghfxuyxw/522BNqnSaebOVPX7I34bs8UP+NkTNm1IU7k0bV/b4IX8bsscPMfOg947KzLwZV/b4IX8bsscPMfOmYnF9oyRJkiQpK/d078/ao2uytyF7/JC/DcY/jPOuisw8eU32NmSPH/K3wfiHcZ5VmbifelzGX1/2NkStNfJ+cV6sV5IkSZK0RNa5xOCc4DXZ25A9fsjfBuMfxtoZLYl7GMSVPX7I34bs8UPMeUXz6Pw4ryhJkiRJysr54f6cy21ljx/yt8H464s4/gyOQasu82Qre/yQvw3GX1/UPClF4R4IcWWPH/K3wfiHsV5peaxXkiRJkrQ0U9S43Knwz1c/jhPVl70N2eOH/G3IHr9UQvbXhfHXl70N2eOHebRhzm4BbqodhNKw5vWa7G3IHj/kb4PxD2O9jzIxb16TvQ3Z44f8bTD+YVwvokxcLxKX8deXvQ1R10x6/6jM5po3s8cP+duQPX7I34aoeVOxuEZSkiRJkpSVe7rHMFZNEzgeV1v2+CF/G7LHDzHHpZ3LVVZ+vnds2dtg/PVFzJlg3lQe5snYsrfB+OuLmicVizVPkiRJkpbGfeFzcYymPuOvL3sbsscvlZD9dZE9fsjfBuOvbw5tmDP3dFcfruO4xvjry96G7PFDzFob61EViXnzGuOvL3sbsscPMfOmFMlc96aF/G0w/vqytyFqDvT+UZmZN+My/vqytyFq3lQsrm+UJEmSJPV1AVwFbh/8hC0fCfwK8H4jBHbIo4HPLvjz1c/H0g5YPQa4L3CXuuEMkr0N2eOH/G0wfml+5vC6yN4G468vexuyx78EHwg8qXYQCuEJtAV7b1//eSHwd/Y871uB7z3zd2W/NmSPH/K3wfjrq9WGTwDeCtztwOMPAl6JmzKoPPNmd9njh/xtMP76ouZNaQrn5szs1wDjry97G7LHDzHzoPeOisq8afy1ZW9D9vghZt5UTFeAf1A7CEmSJEmSephivausaeoje/yQvw3Z44f8bYg6Ju18riI6J8dmv1ZA/jYYf33Z22DOlI4zT+ZuQ/b4IX8bjH84a57yuYI1T5IkSZKWwX3h83GMpq7s8UP+NmSPH+bRBmlMc3hNZG+D8deXvQ3Z418C93TXhus4usseP+RvQ/b4IX8brEnV0pk3u8seP+RvQ/b4IX8bouZNaQruTWv8tWVvQ/b4IWYe9N5RUZk3jb+27G3IHj/EzJuK6Qqub5QkSZIk5eKe7tez9qif7G3IHj/kb4PxD+e8q2owT/aTvQ3Z44f8bTD+4ZxnVQTup567DcZfX/Y2RK018n5xnq5gvZIkSZKk5bDOpTznBPvJ3obs8UP+Nhj/cNbOaA7cw8D4a8vehuzxQ8x5RfPofF3BeUVJkiRJUi7OD1/PudzusscP+dtg/PVFHH8Gx6BVjnmyu+zxQ/42GH99UfOkNAX3QMjdhuzxQ/42GP9w1ist0xWsV5IkSZK0DFPUuDwa+OyCP1/9OE5UX/Y2ZI8f8rche/xSCdlfF8ZfX/Y2ZI8f5tGGuftA4Em1g1B11rz2k70N2eOH/G0w/uGs91EE5s1+srche/yQvw3GP5zrRRSB60Vyt8H468vehqhrJr1/VFRLzpvZ44f8bcgeP+RvQ9S8qZiu4BpJSZIkSVIu7uk+jS7jzGPUNIHjcbVljx/ytyF7/BBzXNq5XEXk53vnjh/yt8H464uYM8G8qTjOuRfNfo3IHj/kb4Px1xc1TyqmK1jzJEmSJGkZ3Bc+F8do6jP++rK3IXv8UgnZXxfZ44f8bTD++ubQhrlzT3fBtHvTQv5rg/HXl70N2eOHmLU21qNqKubNfoy/vuxtyB4/xMyb0lSWvDct5G+D8deXvQ1Rc6D3j4rKvJm7DcZfX/Y2RM2biukKrm+UJEmSJHV3AVwFbh/8hC0vAD5njKgO+JvANxf8+ervnwC/UTuIHQ3tAElXEdvQR8T4G+yD2hq690HE+KXaIr4uGnJfWxtyx99X9vghfxuyx78UjwC+oXYQqu7zac+FBwAfDTwZeDdtwcLGBwB/BjzszN8V8drQsKz37hHb0LCcPogYf0Oe94m/QVvIueuewCuAh04bjhbKvLmcazbEbEPDcvogYvwN+fOmNJVzc2b2a4Dx15e9Ddnjh3h50HtHRTa3vNmQ533zIQ3mzUyyxw/x8qbiujPwu7Q5Q5IkSZKkDEqvd1XLmibHdGtrsA9qa8jRB9YCK5tzcmz2awXkb0PE+PvIHj/kb4M5UzrOPJm7DRHj7yt7G4z/PNY85WLNkyRJkqSlcF/4fGqPcezTsJxxvojxN+QfK+4jYvwNy+oDaWwRXxMN+V/XDebnTLK3IXv8S+Ge7gLXcTQsJz9CzDY02Ae1NeToA2tSFYF5M8f1YiwR29BgH9TWkKMPrE9VbXPbmxbyvP4Pacgdf1/Z25A9foiXB713VGTmTeOvLXsbsscP8fKm4nJ9oyRJkiQpG/d0v561R/nX5fQRMf4G+6C2hjxzEM67amrmSa/RtTXYB7U15M6T0pTcTz1WGxpyx99X9vghfxsi1hp5vzhf1itJkiRJWhLrXMpzTtCxxNoa7IPaGnLPCToWqkzcw8D4a8vehuzxQ7x5RfPovDmvKEmSJEnKxvnh6zmX6/htTQ3O49bWkKcPXNeiGsyT5snaGpbTBxHjb8ifJ6WpuAdC7jZEjL+v7G0w/vNYr7Q81itJkiRJWorSNS5/E/jmgj9f/dUeZ9mnwbHG2hrsg9oaltUH0tgivi4acr+uG3LH31f2NmSPH+bRhiV4BPANtYNQVda85q/36SNi/A32QW0Ned4nWu+j2sybXrNra7APamvInTelKbleJFYbGnLH31f2+CF/GyKumfT+UZEtOW9mjx9itqGP7PFD/jZEzJuKyzWSkiRJkqRs3NN9GqfGmceqaYKY43ENucel+4gYf4Pj6rU15OkDa4GViZ/vnTt+iNmGPoy/vmg5E8ybiuWce9GI14gG82Qmxl9fxDypuKx5kiRJkrQU7gufS8QxmgbH+WpryD1W3FfENjQsqw+ksUV8XTTkfl03mJ8zyR4/zKMNS+Ce7ppyb1qIeW1oyJ3j+4gYf4PvUWpryNMHruNQbUvPmw15rhdjiBh/g31QW0OePrBGVbUteW9ayN+GiPH3kT1+yN+GiDnQ+0dFZt7M3YaI8feRPX7I34aIeVNxub5RkiRJktTHBXAVuH3wE9a+EPj1saLa40OBXwLev+DvUD+/R3tubP78SN1w/puG7oO3UdvQVdT4G+yD2hq69UHU+KWaor4uGnJfWxtyx99H9vghfxuyx78038/p+00tz1uBr9z6/+PWX7vpjJ8Z9drQsJz37lHb0LCMPogaf0Oe94lPAn5l52vvD/wi8CUTxyJtM2/eKGr8fURtQ8My+iBq/A2586ZUW9ecWfv1c0hD7mtwQ+74+8jehuzxQ/027OZB7x2VUea82ZDnffMhDebNLLLHD/Xb4P1jPl8FPKd2EJIkSZIkdVB6vauOs6bpRlHj7yNqGxrsg9oacvSBtcCagy45tvZr7ZAG53KzyB4/5G9D7fjNmcrKPFlXwzLyJORvg/Gfz5qnfKx5kiRJkjR37gufT4Qxjn0aljHOFzX+hvxjxV1Fjb9hOX0gjS3qa6Ih/+u6wfycRfY2ZI9/adzTXfu4juNGUePvI2obGuyD2hpy9IE1qYrKvHmjqPH3EbUNDfZBbQ05+sD6VEWUeW9ayPP6P6Qhd/x9ZG9D9vihfht286D3jsrIvFlXQ+74+8jehuzxQ/02eP+Yj+sbJUmSJElZuKd7N9Ye7Re1DV1Fjb/BPqitIc8chPOuisA8uV/UNnQVNf4G+6C2hrx5UorA/dTracgdfx/Z44f8bagdv2s8l8l6JUmSJElLYJ1LPc4J7he1DV1Fjb/BPqitIe+coGOhmgP3MKirIXf8fWRvQ/b4oX4bzKPL5LyiJEmSJCkL54e7cS73RlHj7yNiGxqcx62tIU8fuK5FUZgnbxQ1/j6itqFhGX0QNf6G3HlSqs09EOpqWEYOgfxtMP7zWa+0TNYrSZIkSZq70jUuHwr8Eu19tGKIMM6yT4NjjbU12Ae1NSynD6SxRX1dNOR+XTfkjr+P7G3IHj/Mow1L8v3A7bWDUCjWvO4XtQ1dRY2/wT6orSHP+0TrfRSReXO/qG3oKmr8DfZBbQ1586YUgetF6mnIHX8f2eOH/G2oHf++HOj9ozJaQt7MHj/EbUNX2eOH/G2oHb/3jjm5RlKSJEmSlIV7ute1Pc48Rk0T1B/POqQh97h0V1Hjb3BcvbaGPH1gLbCy8/O962rIc707l/HXV7sN1kApqy73orVfX4c0mCezMP76arfBmqecrHmSJEmSNHfuC59L7fGNQxoc56utIfdYcR9R29CwnD6Qxhb1ddGQ+3XdYH7OInv8MI82LIl7umtXib1pIe61oSF3ju8qavwNvkeprSFPH7iOQxEtKW825LlenCtq/A32QW0NefrAGlVFtIS9aSF/G6LG31X2+CF/G2rH73pGzYV5s64G82YW2dtQO37vHXNyfaMkSZIkqasL2jGHg3XCJ59AOyj5SuCLRwxs1/fSLmxVHPcCfhf4P4H7AB9UN5z/pqH7AHTUNnQVNf4G+6C2hm59EDV+qaaor4uG3NfWhtzx95E9fsjfhuzxL839gZcCtxx4/BuBFwNvB94CPBu4bZrQVMEdgf8FeBfwwK2vPx143pk/O+q1oWE5792jtqFhGX0QNf6GPO8TP4f2+nTn9f9vAn6Q7vFLYzNvHhY1/j6itqFhGX0QNf6GvHlTqqlvzqz9+jmkIfc1uCF3/H1kb0P2+KF+G7bzoPeOymYOebMhz/vmQxrMm1lkjx/qt8H7x3zuBPwx8PDagUiSJEmSdMQU6121nzVNh0WNv4+obWiwD2pryNEH1gIrsz45tvZr7ZAG53KzyB4/5G9D7fjNmcrGPBlDwzLyJORvg/Gfz5qnfKx5kiRJkjRn7gufU4Qxjn0aljHOFzX+hvxjxV1Fjb9hOX0gjS3qa6Ih/+u6wfycRfY2ZI9/adzTXdtcx3FY1Pj7iNqGBvugtoYcfWBNqqIxbx4WNf4+orahwT6orSFHH1ifqkjmsDct5Hn9H9KQO/4+srche/xQvw3u6a7MzJsxNOSOv4/sbcgeP9Rvg/eP+bi+UZIkSZKUgXu6n2bt0XFR29BV1Pgb7IPaGvLMQTjvqprMk8dFbUNXUeNvsA9qa8iZJ6Xa3E+9vobc8feRPX7I34ba8bvGc5msV5IkSZI0d9a51OGc4HFR29BV1Pgb7IPaGnLOCToWquzcwyCGhtzx95G9Ddnjh/ptMI8uk/OKkiRJkqQMnB8+zbncw6LG30fENjQ4j1tbQ54+cF2LajNPHhY1/j6itqFhGX0QNf6GvHlSqsk9EGJoWEYOgfxtMP7zWa+0TNYrSZIkSZqzKWpcvhf4woI/X/1FGGfZp8Gxxtoa7IPaGpbTB9LYor4uGnK/rhtyx99H9jZkjx/m0YYluT/wUuCWA49/I/Bi4O3AW4BnA7dNE5omZs3rcVHb0FXU+Bvsg9oa8rxPtN5HkZg3j4vahq6ixt9gH9TWkDNvSrW5XqS+htzx95E9fsjfhtrxu7eAsltS3sweP8RtQ1fZ44f8bagdv/eOOblGUpIkSZKUgXu617NvnHmMmiaoP551SEPucemuosbf4Lh6bQ15+sBaYGXl53vH0JDnencu46+vdhusgVI2fe5Fa7++DmkwT2Zh/PXVboM1TzlZ8yRJkiRpztwXPp/a4xuHNDjOV1tD7rHiPqK2oWE5fSCNLerroiH367rB/JxF9vhhHm1YEvd010bJvWkh7rWhIXeO7ypq/A2+R6mtIU8fuI5DkSwxbzbkuV6cK2r8DfZBbQ15+sAaVUWypL1pIX8bosbfVfb4IX8basfvekZlZ96MocG8mUX2NtSO33vHnFzfKEmSJEnq6gK4CtwOcIeBP+TzgQ8DfnScmG7wSbQBPqvQz9cwbwc+EvgV4E3r/4/tEbSDY68CvrrAzy/ZhvsBPw+8grYQ/9Ej/uyNKfrgV2njfznwfxX4+aXbcAm8DPgt4OdG/tkwTR/cCvws7bn0O8CHjPizp4j/CbTnz8uB/0A70SVFZn477b8AbwN+eOSfu1E6/tI5OntugLJteCBtXt78eSfwyBF/Ppjfsvkj2tfi1xx4/AJ4BvBQ4LOA9wDPB+4xRXCazAOBPwf+mra/Hwn89tbjtwJvOPN3ZM/x2fMjlG9D9vtXKHsP63uU007l9zfQbrTwYev/fzrwBbTXrM17mwcilTeHvJn9vti8eZp587RbWVbelGoYmjPNg6eVHCM2j592K7lzyBzGiPvkQe8dlYV507y5j3nzNPPmad4/zs9f0W7M8Y21A5EkSZIk6YjS6111oznUNIG1wKdY03TaJbnXs95K7nF1OD4ubS2wMhqSY7PnTHCN7CmukT3NNbKnmTM1B+ZJ8+Q+c7g/vpW8eRLK58raeRKsecrImidJkiRJc+a+8Dm5NvO4OcyHZa+DuSR3DQzkH2t1z1tl4zzcac7DHWduOG4J9Srqxz3dBa7j6CJ7fgTXcXRxSe572FvJ+x5lw5pUZWDePM28eZp58zTz5ml98qZUg3vTujftPnPI47eSO4fMYYzYPd01R+ZN8+Y+5s3TzJunub5xflzfKEmSJEnKwD3dD7P2qJuSbZjD2Cvkrz+6JHftEZQdA69dewTOu6oO82Q35snTzJPHmSdPc55VGbifuvs47OM+Dqe5j8NprvHULuuVJEmSJM2ddS7Tck6wG+cET3NO8Ljs47m1x0LB2hnNg3sYuIfBPnN4L3AreecUYR7ziuZR7eO8oiRJkiQpA+eHD5vDXG728ec5jN9mn8eFsnO52edxIdb4MzgGremYJ08zT55mnjzOPHmaeyAoA/dAcA+EfeaQx28ldw4pXa9k3a9qsV5JkiRJ0pyVrnH5JOB24FmFfr6GyT5e6ljjac4ZnuZ472mnxkulaKyJOc29C44zNxy3hL0L1M8f0b4ev+bA4xe0NR0PBT4LeA/wfOAeUwSnScyh5hVy5/g55EfIfw97Se77V3CfP+t9NAXzZjfmzdPMm8eZN09zvYgycL2I60V2zSWP30reMeI5jA+Dn5mpeVpi3sw+rwrmzS5uxbx5TLS8qRxcIylJkiRJysA93ad3bJz5Vs6vaQLX+5xSelzXmqbTLrGm6ZTa632sBVY2fr63c7n7OJd7mmtkT/PzVTQXQ+5FrXk6zTx5nHnyuCXmSeVgzZMkSZKkOXNf+Hyyj1OC43ynZK+Bgfx1MLeSuwYG3PdW+VivcprzcMc5D3fcEubh1I97umuKvWnBHH+K6zhOcx3HabeS+x7WdRzKwLzZjXnzNPPmcebN06xRVQZL3JsW3NP9lDnk8VvJm0PmMD4M7umueTJvmjf3MW+eZt48zT3d58f1jZIkSZKk0VwAV2kXfR7yU8B/LhjDDwLfVPDn65ortP197M/F+rkPWv//g0f63btFDe9PO/B8P+ADaAcQ77vn+55IO3i++fNu4F07X3vYgd85Zht24//vgE9c//tewGtp27HPFbof922l++COwAdt/fuFtIvAd0XtgzvSFjDdtcP3XiFuH/wC8Bnr/9+N9rWxa2gflI7/3sDvA3cCbgKeDfxPI/0uqRTz2/X2XZc+E/g8jk+wRo7/vnTL0VfInRuuMCx+KN8HG3cF/gS4y4HvvULMPjC/je/TaCem953Lu+4KvJf2OqT5eD/go2hfv0+hvTbctvX4c4H/Z+d7rtDvGlH62nYzZXN86Wtb13vYK8TNL13uXyFuH3S9h71CzPxY+v4Vyp5DXfL7A9a//zakurLnza73xRD3mm3erN8H5k3zptTFkJwJ5sFt+zZsKz1GbB6/3tAcArHz4MaxMeIrxOwD86DmKnPejPq+efN/86Z5sw/zprL4H2n77eM7PPeZwB+z/9z9IODbacd53rX+md8wTog6g/0Sj30Sj30Sk/0Sj30Sj30Sk/0Sj33SzafQHpsvrx3IAaXXu+pG2WqaYPpxacd0r2dNU8w+6DquHrUPTo1LOyatjIbk2Ow5E2KskYW417vSeT9KHe0V4ubMjbmukTVnKgvz5PnMkzHzTOY8ufu1EjVPtfMkmCuz6lrzdKzeCZxHmILHuCyPb3ke47I8vuV5jMvzGJfl8S1vicd46XUy7gvfzRXq1cHMcU+DOcyHuad7/T6Y+94+UjTOw13Pebi8ueEKMefhllCvov7c012u45h/fgTXcexyHUe8PrAmVVmYN+tfL8yb9fvAvFl//sO8qQwy700LMV7/MP3nfprHr+fetPH6wHUamivz5jjMm7Gu2ebN+n1g3pyvsfZ0h2XWc0dnn8Rjn8Rjn8Rkv8Rjn8Rjn8Rkv8Rjn3Sz9LWqmQ0dD75C9zGo7LVHUHbcrOvYK8Qeu3Tdaf0+KLnutHbNLjh+rDrMk/Wv0ebJ+n1gnjRPSl25n/r5pq41gvLx3xf3cahdL7Ux130czIHLZb1SHh7f8jzGZXl8y/MYl+XxLc9jXJ7HuKwlHt/oNS5gncvUnBOsP5bonGD9PnBO0DlBaSruYTAO9zCINyfnHgZ1+8A8qmPG+pznIZY47hKdfRKPfRKT/RKPfRKPfRKT/RKPfRLPEvvE+eHcsn0GyxzHn+cwfpt9HhfKfgZL9nlccPxZy2WeNE/25WeVmSfPZZ5UVu6BcL6htUoQN4fMIY9br1Q3fvOgjhmrXmmJY/pT8viW5zEuz2Nclse3PI9xWR7f8jzGZS3x+FrjAj8IfFPBn6/WFerNGcL046WONV7POcOYfdB1vDdqH3QZL5WisSbmeu5dkDc3XCHmXNoS9i5Qf58GvIHD5/O2uwLvpb0OaR7c569+vY/7/NXvA/f5q/8+0XofZWHerH/NNm/W7wPzpnlT6sr1Iueb22dmziWP+5mZdfvAz8zUXC0xb2afVwXz5i7zZrw+MG/Ol2sk58s+icl+icc+icc+icl+icc+icc+icc+6cb1rtrn2DjzGDVN4HqfXVOvl7Gm6cZ4d/9/iTVNNfvAmibNkZ/vPQ7ncvPmnCvEjX9jrmtkzZnKYsi9aOn7gpux5sk82Y950jyp6VjzlIPHtzyPcVke3/I8xuV5jMvy+JbnMS5ricfXOhn3he/iCvVqYGB+exrMYc9Y93Sv3wfu6S5Nz3qV6y11Hu4KMXPD5mvOw8Weh1N/7um+bFPsTQvm+G378ovrOFzH0ZfrOKQ6lpQ3o14vwLxZuw/Mm/X7wLypLJa4Ny3E2NMdYtz/wjz3MXBv2rp94N60mivz5vnMm/Gu2ebN+n1g3pwv1zfOl30Sk/0Sj30Sj30Sk/0Sj30Sj30Sk/1yWte1qhfr590+9An3pO2E1YlfdBPw94G7mrEZAAAgAElEQVQXAe8A/hJ4CfB17C/82LgX8B7awVGV9yHAx574sxmAfDzti3Db/YAfp+3b3wE+Yv31lwAfuv73pwA/sv73s4GnAb8GfOXO/5+0/lkbTwa+Yk/M96A9PzZ/foT2ArH9tTsfaG/XNgyJ/yt3fu7LgPsfiKPPcR8SPyO04c7Ab7C/gClyH1zSrYCpdB8civ9UG/4h8PwO8Q/tg9Ln0DcAf0g7wH3Lui0P7tAeqabSr+sa+a1Ebrjg+ARr5Nyw61COjpqfu+aGofH3acO5ffB3Ob55TdT8bH4r4zWcvq+EtijlKvDpR55zb9rFrd++/v89ac+LHwNeDbwT+K/AL9PeNN9hUMQq6fnA92z9/weAZ+08p+81ovS17aGUzfFT3v/B+PmxTxtK3r9C+T4Yml+63sNmf48S9X1il/y+KcS5d4d2SlPKlje73heDebNLG8ybx5k3T7fBvKkl6ZIzoez4aqQ82DX+Q1+7oMwYsXl8nBxyThumyOMbx8aIx+yDqec/zYOag0x5M+r75u3r3QXmTfOmedO8Gc9jge8Afgl4O20/fH/H772JdjH1vzzxvAcB7wP+8YHHn73+vT8J/AvaRdgf1zGGWnbnoqE9519Au8D/ncBbaV8zT6Kdqz7kYbSvpTcCf73++2eAR4wedT/Z+mXMPoGY/ZKtT+DGfrFP6vP6FY/Xr5i8fsXj9Sser18x1bp+/dj68S7rtqZ0T7qtd1VZEWuaoO64tGO6eWqaOKMNl5SraerThtLj6lHryk6NSzsmrTnokmNLX+/mMpd7welNiKNe70rn/Sh1tFFzZqQ1sgxsgzlTc2WeNE/Wuj9eep489DUoU/NUO0+CubKWc+qdoFvN06l6J5jHPAIM3xNyCtmOsce3PI9xebvH2OM7Ls/h8jKdw5DvGHsOl+c5XNZY53D2Ohn3hS/PvX3KjTFBvvmwrnUwkceKL3FP97594J7uWjrn4ZyH24ianzPOwy2tXkXDuKe7trmO40au43Adx0bUPnAdh1SPefNG5k3z5kbUPsiQN0/Nf5g3ldHQvWlhWfOfp9rQdYw4cg7Zli2PR9qbloFtKL1Oo08bXKchHWbeNG+aN82b5s35i7CnO8yjnvtpuH9lTfZJTO6JHI97Isfj9Ssmr1/xeP2Kx+tXTO7pfr174p7ufXQdD+4zBuW6nHHGXiHH2OWx+qPIfXCJ605rzkE4fqwszJP7mSfNkxvmyePxc6QN5knN0Zj7qUddq7KvDSVyyAVlao32xX+oDWPkcXAfh32m6oPa+zi4xlP7WK80jHuoluceqmV5DpfnMS7L41uex7g8c11ZnsPljXEOR61xAetcInBOcD/nBJ0T3Ig6luicoJSTexi4h0Gf9wJR5+Tcw6B+H5hH5y3K5zxvfMk6hqu05/Ux2cZdso1rbcy5TyDfeO5G136ZQ59Ajn6Z82vFPonJ61c8vlbisU/isU9icn5Y2yJ+BstSxp83phi/dR73xvgvKfcZLFHW5UTtA8eflYl58kbOc46TJyH2e5VLzJPmSamfMfdAgGXX/V5wvFbpnDaYx/PUK40x5lCiXsm6X50jSr3SUsf0pzKH4wse47FlOoch3zH2HC7Pc7gsz+HyPMZleXzLG+MYz6HGZeh5dS/gPbTjRSpr6jlDqDte6lhjnjlDzmjDJeXmDPu0ofR4b9R52y7jpVI0S6mJce+CVrb8HGkukIFtWMLeBRrmNXRbO3Ff2nGBTz/w+O49+D1pz4kfA15N+xmd/xX4ZeDLgTsMjliluM/fflPd/4H7/NXqg0vc56/m+0TrfZSVeXM/86Z5c8O8eTx+jrTBvKk58jMz/czMOeTxrjkk6vhq1/Hhc9pQe/7THKi5WELenOre6wLzpnnTvGnejMc1ksPtzkc/DXgB8Frauei30r5mnkQ7V33Mw2hfT28E/nr9988Ajxg96u7m0CcwvF8i9gnk65cx+wRi9ot9Yp+MwetXPL5W4rFPYvI9cTy1XitzWO+qsrbHmceoaQLX+5xqQ+n1MtY0WdO0EbWmwZomLYGf7+1cbrS53NI5J8P7lrmukTVnKqsu96Klr9EPxZon86R5csM8qbFZ8zSM+4yW5zEuy+Nbnse4PPeFL8tzuDzP4bI8h8tzX/iW+8KXtbQ93fe1wT1j3dM9Uh/02bc3ai2Y+94qo6XUq5xb43nBdPNwS8/Pu19zHi7+PJyGcU93bZTYmxbM8afyi+s46tcKXeI6DtdxSP1NlTeXvqf7oa9tmDeHt8G8eZx5UxrXmHvTwnLy4JDx4XPaMNUY90a2PO7etPXzuHvTainMm+ZN86Z507w5b65vHM79K+NxT+R43H83Jq9f8Xj9isfrV0xev+Lx+hVP5D3dL2iP5e1Dn/Bl68c/4sDjG89cP+/NtAObT6cd2LhKOzh304HvewLwqhM/W3V8J/D/bf3//YDfBj5z/f+7Ae8P3Mz1A1wr4JvX//4D4H/femz7/49d/46Nf0D7Yj+l6fg86N6Gy63nrOgW/7YHAa/k8Hk+VOk+2Pg12sXiT+sYV0OcPvgD2sKrF9MO4o7t3Pg3MR5qwyNpL+abAeF/3jGuhm59MMU59LW0A5x/yvWvaSmq0q/rqfNbqdxwwekJ1m0NcXLDthI5uvS1dWhu6GOqPvhx4FFnR3uj0q9jML+V8F10Kx74Idpz/9jmRf8r7b3mps+/av3/N9AW1zwF+I/An3H6nlR1/CzXT45/PfDyM39m6Wtb6Rw/1f0f1L+Hvdx6zoqy968wfh9cbj1nRb/8UvIeNup7lIY47xNP5ffHA6/vGKs0pWx5c2jOBPOmedO82WDelM7RNWeWHF+NlAcvt56z4nD8h752QZkxYvN4/vFhiDVGXGP+0zyoOcicN5uOz4Pp5sQuMG8eYt40b5o36/kt2rnSd9C+vvsuXv0R2o2Dj/kZ2jnZfYvvP3b9O3+6x++MYHcuGtoP7nkR7Rz0U2kXBb94/bzX0y4I3/VN68ffAvwn2kXA373+vm8pFHsXGftlrD6BmP2SsU/gxn6xT+rz+hWP16+YvH7F4/UrHq9fMdW6fn3q+vlPHKMRI/oyuq13VVkRa5qg7ri0Y7p5aprO6YPo61lLjas3xKkrOzYu7Zi05qBLjnUPiDJrZCHW9W7DNbL7LWWNbKm5XHOmsjJPtsyT10xxf2yePPw1KJMra+dJMFfWcm69E5yueTpW7wTzmUeA4XtClpbxGHt8y/MYl7d7jD2+4/IcLi/LOQw5j7HncHmew2WNdQ5nr5NxX/h4lrK3z+XWc1YMG2PKOB+24Z7uh2Xf0x3qr82UonEeznm4c1mvYr2KhnFPd21zHceNXMfhOo4xuI7DmlTNk3nzRuZN8+YYlp43u8x/mDeVzZC9aWF5859QZoy4IU4O2ciYxyOND4/RB7XHiC+3nrPCdRrSNvNmy7x5jXlzP/OmeTOz2nu6w3zqud2/si77JCb3RI7HPZHj8foVk9eveLx+xeP1Kyb3dL/el+Ge7n3UrD1a+rocqD/26rrTuLVHmxhrrTuNsC7H8WNFYJ7czzxpnjyXedI8qXkacz/1y63nrIizVgXK788C5WqNYLo8Du7jcMhS9nG43HrOCtd4qmW90jDuoVqee6iW5Tlcnse4LI9veR7j8sx1ZXkOlzfGORy1xgWsc4nAOcH9nBN0TvBczgk6Jyjt4x4GLfcwuKb2e4HLreesyDen6B4G5tE5ivA5zxv3Wz/vHes4Hn/kuRnHXTKNa23MvU8g13juRtd+mUufQPx+mftrxT6JyetXPL5W4rFP4rFPYnJ+WNsifgbLUsafYZrxW+dxD38eSKm53Kjrchri9IHjz8rCPHkj5znzfFbZ5dZzVvR7r2KePM48Kd1ozD0Qll73e0G/WiUwj8+xXuly6zkrho05lKhXsu5X54hQr7TkMf0pzOX4gsd4bFnOYch5jD2Hy/McLstzuDyPcVke3/LGOMZzqHEZel49AXjVKJFqbNn3jHWsMc+c4Tl9UHLOEOKO9zbEmbc9NV4qRbOUmpjLreescO+CMUXNDV25d4G5rZbvAp7f4Xk/RHv+3/HA47v34F+1/v8bgB8AnkL7WZ1/Rv1xDu3nPn/7uc/f/NeHRL9/3cRYq+7Veh9pP/PmfuZN8+a5zJvmTc2Tn5npZ2Zuy5rHI40RX249Z8V8xofBz8yUYBl5c6p7rwvMm4eYN82b5s16XCM53O589LuAF9HOQT8V+A7aMaurtOf3/Q78nG9aP+ctwH8C/iXw3evv/ZZCsZ8ylz6BYf0SsU8gZ7+M1ScQs1/sE/tkLF6/4vG1Eo99EpPvieOp9VqZw3pXlbU9zjxGTRO43gfqjuta02RN0yENcWoarGnS3Pn53i3ncq+pPZd7ufWcFeVzTleukTVnarm63IuWfn1Z83SaedI8eS7z5HJZ8zSM+4yW5zEuy+Nbnse4PPeFL8tzuDzP4bI8h8tzX/iW+8LHkn1tJkyzpwHEHecrUQMDscaKo9fBlBpnbYhTC+a+t8pmKfUql1vPWRFnT3cwP++L13m4fmrPw2kY93TXRom9acEcfyq/bHMdx36u43Adh7U2imiKvOme7oe/BubNQ8yb5k3zpiIac2/aJeXBMcaHId4YN+TM4+5Ne1rt+U9zoObCvNkyb15j3tzPvGnezMr1jcO5f2U87okcj/vvxuT1Kx6vX/F4/YrJ61c8Xr/iibyn+8X6ObcPfcJ/oL25POaR65/xGuBDtr5+C+3g0FXaAY19nkd7kBTPL3P9IOkX0J4Pu/4Hri9A/9fA44APpD3ZNwOau/9/HDcOQD+pQ1wN3Qegu7RhaPwb9wReATy0Y0x9lO6DbR8EPBe4rUNcDXH64MPWf9+XdqH7AzvG1dU58cPpNjwWeBvwN2ivmT9Jt4Hohm59UPocujttwr4H7SDkz9PmFZ1nRZs7L+qGMVulX9dT57dSueGCcgtXp8jPUC5Hl762Ds0NfUzRB3cD/hi40ygRX6/069j8VsZjgPfSvjYP+VbgjcBHnfhZzwH+hGuLWz8L+DzgDjvPuw/wR7R59TE949V4ngo8DLiV9v36U4D3AZ+79ZwH0p4fH7L7zT2UvraVzvFT3f/Vvoed+v4Vxu2Dc/ILlL2HjfoepSHG+8Qu+f2Ze36fNLU55M2hORPMm+ZN82aDeVPq6pycWXJ8NUoe7Bo/B74G5caIzeO5x4ch3hhxjflP82BdK5xL7WtuebMhxvvmbReYNw8xb5o3zZvDrTgv530m8ADa/rqg/+LVb1x/z4cfePyjafPJdx94/KvX3/8VPX5nBLtz0XD4tflk2jY+Y+frj1t//Xm0r5tdtwyMbcX574My9ssYfQJl+mXFMvsEbuyXKH0C5/fLXPoEvH7V5vUrJq9f8Xj9isfrV0w1r1+vBP6Qwx9IU0OX9a4aV4aaJqg/Lu2Ybo6aJpjvetaS4+oNMerKTo1LOyY9jRXO2Y5laI51D4gya2QhzvVuwzWyhy1hjSyUm8s1Z05rhblzCPPkNebJuvfHS8+THPgalMuVtfMkmCvPsWJ43ju33gmO1zydqneC+cwjnLMn5DErljlXM9XxBeeNPYfL2T7GnsPj8xwuL8s5DDmPsbmuPM/hssY8h7PWybgvfExL2NtnjDGmrPNh29zTfb/SY8VT7+3jnu55rHBurhTn4ZyHO5f1KtaraBj3dF8u13HEWUMAruM4xnUcruOwzmYaK7zfPca8GeeaDebNY8ybufNml/kP86aiG2tvWljW/CeUGyNuiJNDIG8ejzI+DPnHiF2nMV8rvK/sy7x5jXnTvNmXedO8WdOK3Hu6wzzqucH9K2uzT2JyT+R43BM5Hq9fMXn9isfrVzxev2JyT/fruaf7YZFqj8B1ObXHXsF1pxCz9gjqrzuNsC7H8WNNzTwZ5xoN5slTzJPmSfOkaiu5n3rktSpQdn+WjQvK1BrBdGtV3MfhsCXs4+Aaz/laYb1SDe4DXJ57qJblPsDleZ0oy3O4PM/h8sx1ZXkOlzfWORyxxgWsc5mac4KxxhKdEzzOOUHnBB0PVXbuYXCNexjEfC+QfU4R5juvaB7Nb0Xuz3lm/fufD/w+7Z74V4HHH3l+xnGXTOPfsIw+gVzjudCvX+bQJxB7DBiW8Vrx+hWT1694vH7FY5/EY06Jyfnh5Yo0l+v483Tjt87j7o+/5Fxu1HU5DTH6wPHnela4p/wx5knnObdl/6yyc9+rmCePM09q6UrugQDW/V7Qr1YJzONzq1caY8yhVL2Sdb/LtqLuHghwfr3S0sf0T1nh2luIPb84h2Mc+RyGnMc403xsxuMLnsOleQ6XZ64ry3O4vLGOceYal3POq+cB//G8EFVI9j1jHWvMMWcI890noeR4b0OMedsu46VSNEuoiXHvgrz5OcpcILh3gbltfI+hrdG455HnfCvwRuCjjjxn9x78s4DPA+6w87z7AH9Eez/6mAHxahyRal7Bep/a+RGWvT4k8v0r1K97td5HMm82xLlmg3nzFPOmedO8qdr8zMyWn5k5vzweZYx4ruPD4GdmapmWmjenuve6wLx5iHnTvGneHG6FayRr2Z2PPvTafDJt+56x57HHrR97Hu3rZtctA2NbYZ3+Rt9+idonkLNfxugTKNMvK+yTjSh9Al6/vH7FEfm1ssI+2YjSJ+D1y/fEcdR8rWRe76pxnRpnHqOmCVzvU3tc15oma5oOaYhR02BNUxwr3Dt+DH6+9zXO5cacy62Rc7pyjaw5U8sw9F609DXamqfTzJPmyXOZJ/NaYc1TDe4zWp575ZblOVye53B57gtflteJ8jyHy/IcLs994d0XPqLsazNhmj0Noo7zlaqBgVhjxZHrYEqOszbEqAVz31tltIR6lch7uoP5eTde5+H6qz0Pp2Hc032ZptqbFpad47vklw3XcRzmOg7XcVhrM40VruM4pFbeXPqe7hz4Gpg3jzFvmjfNm6qt5N60sJw8ONb4MMQb486ax92b9rTa85/mQGVk3rzGvGne7Mu8ad6sacXwsVTXNw7n/pXxuCdyPJH3312xzD4Br18Ref2Kx+tXTF6/4vH6FU/kPd0v1j/zdrixkLeLTwdefuI5j17//W20B2Lj3cA/W//7a/d8383AQ4DfGhCXyroJ+ASu75sHAr++57m3cf058qD1/28DfpX2BGTP/18P3G/r+z4ceMO5gW/p2oah8QO8P+3C6qesHxvTFH2w7e3AzwKfc1bU15uiDzbnzBuBnwI+ZYzA186Nf/PYsTa8DvgN2gv5u4GfAD5pnPAnOYf+FvD7wFuBd9JOAnzqOOFLRUzxup46v5XMDSVMkRugXI6e4tpaMjfAdH3wSNpJyr8aJeprpngdm9/KeAntPentBx7/NuBLgYcDrz7yc+5Gu1j12bRFNNC+j3027eTotjcB/27974veEWss96GdCP9d4AXAg2kLop6z9ZzfBn4N+MKBv2OKa1vJHD/V/V+Ee9i53r9uHqt1D+t7lPPz+51p37/8+zPbIZ1rDnkzyrjvJjbzpnlzl3nTvKl5GJozS4+vRsmDXePnwNdKMY/nzyH7Yq49Rjz1/Kd5UBmZN68pde9VgnnTvDkG8+ay/RzwKoZftzbnxMMOPP73ac+xH9r5+mPWv/O71v//7vX/rwIfOzCWqeybi4bDr81nrf9+wNbX7gA8DfhL4IuAd+z5vnefF+YgWftljD6BmP2StU9gf7/YJ3V5/YrH61dMXr/i8foVj9evmGpfv/4zcH/acYAouqx31bgy1DRtHqs1Lu2Ybp6apk1sc1zPuoRx9WPj0o5JK6MhOdY9IPLM5Z47j+ga2cOWskZ2E9vYc7nmTGVhnrzGPFn3/njpeZIDX4MyubJ2ngRzZU3n1jvB8ZqnQ/VOML95hKF7QpaU9Rh7fMvzGJe3e4w9vuPyHC4vwzkMeY+x53B5nsNljX0OZ62TcV/4eJayt8+5Y0yZ58O2RdzbJ3sNzOaxSGOt7umupXMeznm4c1mvYr2KhnNP9+VyHcdpruPIc/+6iW2O97DZ36Psi9maVGVk3jzNvGneHMPS82aX+Q/zpqIba29aWNb85+axWmPE5vE8OWQTW+YxYtdpSNeYN68xb5o3+zBvmjezO3eN47H1jTDPNY7uXxmPfRJT7T1FS5lbn4B7Itfk9Ssmr1/xeP2Kx+tXTLWvX1nXqi5VpNojmO/4cZax101smeuP5lp7tHms1hh4hHU5jh+rBvNkN+ZJ8+S5zJPmSc1Dyf3Uo65VgfL7s5Q21VoV93E4bCn7OLjGU4dYr9Sfe6iW5x6qZXkOl+cxLsvjW57HuDxzXVmew+WNeQ5HrHEB61ym5pxgN84JOid4LucEnROUwD0M3MMg/nuB7HOKm9jmOK9oHl22c+cU4fi84rE5xW1fRzsm8feAvzjyvKzjLlnGtbbNvU8gz3juti79Mqc+gfj9MvfXin0Sk9eveHytxGOfxGOfxOT88LJFmstd+vjzlOO3zuPuj7/UXK7rchx/Vl7mydOc58yfJzeP1ap5Mk+aJzUPJfdAgPx5MPseCLDcPB49h+yLecp6Jet+da6a9UqO6Zc1p+MLHuOxZTiHIe8x9hwuz3O4LM/h8jzGZXl8yxvzGGeucRna5puBh3D9vbxiyL5nrGONeeYMN7HNcZ+EJYz3nhovlaJZSk2Mexfkzc9RcsMmNvcuMLeN6SW0n415+4HHvw34UuDhwKsPPGffPfjPrv//vp3nvgn4d+t/XwyKWGOIVPMK883xWfLjJrbM97BzvX/dPOY+f9b7qC7zZjfmTfPmucyb5k3Ng5+Z2fIzM+eXx6OMEc91fHjz2NA8aA5UVkvMm1Pfe5Vg3jRvnsu8uWyukRxm33z0odfms9Z/P2Dn63cAngb8JfBFwDv2fO+7zwuzt7n1CfTrl4h9Ann7ZYw+gZj9Yp/YJ2Py+hWPr5V47JOYfE8cT+3XSub1rhrXqXHmc2uawPU+m8dqjeta02RN0xhK1zRY06S58fO9r3EuN+ZcbtSc4xpZc6aWY8i96BTXaGuejjNPmifPZZ5cNmue+nOf0fI8xmV5fMvzGJfnvvBleQ6X5zlcludwee4L33Jf+Fiyr82EafY0iDzOl6UGZhPbHOtgotXAgHu6S0upVzm3xrMk8/ON8ToP10/teTgN557uyzTF3rRgju+SX8B1HMe4jiPePazrOLREtfLm0vd058DXzJuHmTfNm+ZNRVByb1pYTh6sPT4M7mNQO4eAe9P2yYPmQGVl3rzGvGne7MO8ad7MzPWNw7h/ZTzuiRxP7T1FS8raJ+D1KyKvX/F4/YrJ61c8Xr/iqX396rVW9eYuT9pyC/BRtAXGx9xn/fdr9jy2+donAx8M/NnWY58A3AV4Zc+4VN5V4IN2vvZm2sEnaE/YuwN/CtyDduAJ4NNoC0heBfw94GVb33/bzv9/Hfh42kHoPwEeRbcTeTVyG4bGfxPQ0Bb9fF/HmPqYog/uBrwf8BbgTsDfBv7vDrGtRm7D0Pjvsv4Z7wDuSnsxfhbjOTf+fTHv/v/FwD3Xf94K/E3aCYFTViPGf8459FrazQPuRJukLmgTuBTVFK/rqfPb2K/roVYdn1c6N0DZHD3FtXVobuhqij4A+ALgO8cJ+TpTvI7Nb2X8AfAXwEOBH9957Ntpz5nPBF5x4uf8Hdr3sD/a8fdubqTf0/H5Gt+q4/O+GfgO4N9y/QBLF1Nc20rm+CnyS5R72KnvX2HcPjgnv5S8h438HmU1chtK3b9+OfCi9R+pplXH50XOm0NzJpg3zZvXmDfHaYN5U3O26vi83ZxZenw1Sh7sGj9MO0ZsHs8/Prwv5tpjxFPPf5oHldGq4/Oy5M2u7ZlqTmyIVYfnmDfNm2Mwb+ocv7f++7YDj/8t2nyx279vpM0pX017jj15/fWrXHttRNV3Lvrz1n9vvy4eCnwE8MPA29Y/8zbaovJfB144SqT9Ze2XMfoEYvZL1j6Bfv1in0zD61c8Xr9i8voVj9eveLx+xVT7+vUr678/G3hut5CL6rreVeNadXxezZomqLvexzHdPDVN5/RB9PWsJcfVVxPEvy/mPuPSjkkro1XH5+3mWPeAyDGXe871zjWyxy1ljWypuVxzprJYdXyeedI8ObbS9TbZ8iRMmytr50kwV2Z3rObpUL0TzG8eYeiekCVlPcYe3/I8xuXtHmOP77g8h8vLcA5D3mPsOVye53BZY5/DWetk3Bc+nqXs7XPOGFP2+TD3dD8u+57uUH9tphSNn1niPNy5rFexXkXDuaf7cq06Ps91HKfjdx2H6ziOcR2H6zg0D6uOzzNvno7fvGnePGbpebPL/Id5U9GtOj7v1N60sKz5z30xjzVGvOrwHPN4nvHhMT47tvYYses0pGtWHZ9n3jRvjs28ad40b+Y2dE93mFc99zHuX1mefRJT7T1FS1lKn4B7Ik/B61dMXr/i8foVj9evmGpfv7KuVV2qVcfnTVF7tOR1OVHGXl13Grf2aF/MU647jbAux/Fj1bDq+DzzpHnSPHke86R5UvOw6vi8IfupR12rAuX3Zxlq1fF5U+wh4D4Oxy1lHwfXeKoU65VaEff4zHp8wT1US/McLs9jXJbHtzyPcXnmurI8h8sb8xyOVuMC1rnUsOr4POcEnRN0TvA8zgk6JyiBexgMjX9fzO5h0I97GOSfVzSP6lxDP+d54+OApwJPB36R9v3yIVnHXbKMa20soU8gz3juRtd+mVOfQOx+WcJrxT6JyetXPL5W4rFP4rFPYnJ+eNlWHZ9X8zNYljD+PPX4rfO4N47dlpzLjbwuZzVyGxx/1tysOj7PPHk6/iXPc0b5rLJz3quYJ8dpg3lSc7bq+LwheyAsue73HKsOzzGP56lXGqPmrUS9knW/imBovZJj+mXN6fiCx3hsGc5hyHuMPYfL8xwuy3O4PI9xWR7f8sY8xplrXIa2+RNox/xfeUaMKkqAJgsAACAASURBVGOKcaKS46WONeaZMzynD6Lvk1ByvHc1Qfz7Yu47XipFs5SaGPcuyJufo8wFuneBua2EPwD+gvbzMX9857Fvpz1nPhN4xZGf0fdzU9+9/vs93cPUyFYdn+c+f+7zN/f1IdHvX/fF7D5/0vRWHZ9n3jRvmjfPY940b2oeVh2f52dmLuszM+eQx6OMEc91fBjOm/80ByqrVcfnzSlvTjkn1teq4/PMm+bNc5k3dS7XSB73eeu/d1/XDwU+Avhh4G3rn3kb8Fe0ue+Fo0Taz1L6BPb3S8Q+gbz9MkafQMx+sU/skzF5/YrH10o89klMvieOp/ZrJfN6V41r1eE559Q0get99sU85biuNU3WNB2z6tIAytc0WNOkuVl1fN5SP9/budz6c7k1ck4XrpE1Z2o5Vh2eM3WeBGueTjFPmifPZZ7UuZZW8+Q+o+V5jMvy+JbnMS7PfeHL8hwuz3O4LM/h8twXvuW+8LFk39Mdys//Rx/nc0/380Xet3c1Qfz7Yu47VilFs5R6lXNe10OtOj7P/Hzj8XYerp/a83Aazj3dl2nV8Xmu4yifX1zHMU4b5nr/ui9m13FI01t1fN7YeXPpe7qDeXNjNXIbzJvmTamkVcfnDdmbdkl5cKzxYYgzxp09j7s37XG15z/Ngcpq1fF55k3z5pjMm+ZN82Z+S1vfCO5fGZF7IsdTe0/RkrL2CXj9isjrVzxev2Ly+hWP1694al+/eq1VvbljkBu3AncE3nDieX+y/vsj9jz2kVv//liuv1H96PXfr+8Zl+pogB8CXk470PS1wC8DzwF+gvZ8edX6z3uABwIv2Pr+3f+/B/jf1l+7A/B0Tp9r52q4sQ1D4/902oHDlwGPXH/tS4DfLhU84/fB3YEfoV1gfgfa4p+fKBg/jNsH96a9+N5Ee63697SD0iU1dI9/X8y7/38v8A3Az9G24xcYtwhrV8O459CL1t/7EuB968f+S8H4pRIaxn1dT53fGsZ9XQM8H/hE2kLR1wGPo+xNQsN4uQGmz9EN4/bB1LkBxu+De9BuOvO8gjFvaxj3dWx+K+Mq8EfAx+x8/RnAF9O+Xt/KtU2N/nz9Z9ejaBfAdjm/bga+dP3vn+4Zr6b3XNrCjA8H/nCEn9cw7rUte46Pcg871/vXfTHXvodt8D1Kn/z+rvXPl7KInDejjPuaN82bfTSYN82bmqtTObNh2fOfHPjalGPEDebxTDlkX8y1x4gbpp3/NA9qzsyb/ceNwbzZV4N507ypLN64/vsj9zx2F+CTaD9w5i92HvtV4NeA/2P9+JVC8ZVwai7662k3Wrsb8CDgdtpr6FO3nvPg9d9vBn6T9nWz7ReBx9JuOjelrP0yRp9AzH7J2idwvF/skzq8fsXj9Ssmr1/xeP2Kx+tXTLWvX5u5+s/oG3ght9JtvavqqFnTBPHGpRsc041Y03ROH0Rfz7qEcfVj49KOSWvOnMud1xrZfTHXzvsN1tFGXCNbai7XnKm5MU+aJ8E8eUrD+PXm1jwpqkM1T8fqnWB+8whD94QsKesx9viW5zEub/cYe3zH5TlcXoZzGPIeY8/h8jyHyxr7HM5aJ+O+8Dk0zG9vn3PGmLKP80Wpg5lrDcy+mGuPtTa4p7u0q8F5OOfhztOw7Hm42nNwYH7Lwj3ddYrrOK7X4DqOiPevruOI+x5lX8zWpGrOzJvXazBvmjfP17CcvNll/sO8qbnokjMbljP/uS/m2mPEDebxiDnk3M+OjTBG7DoNqT/zpnkTzJunNJg3zZvLMXRPd5hXPfc296+cnn0SU+09RUuZa5+AeyLX4PUrJq9f8Xj9isfrV0y1r19Z16rquClqj5a8LifK2KvrTuPWHu2LufYYeIP740kb5skbNZgnzZPnaTBPmic1V0P2U1/yWhWIvY9Dl/inzuUN7uOQqdYIXOOp81iv1Iq4x2fW4wvuoVqa53B5HuOyPL7leYzLM9eV5Tlc3pjncLQaF7DOJTLnBG/U4Jygc4LnaXBO0DlBaT/3MHAPA3APg1Ma3MPAPKo+hn7OM7T71X8f7T75T+zwu7KOu2QZ14Ll9AnkGc+Ffv0ypz6BuP2ylNeKfRKT1694fK3EY5/EY5/E5Pywuqj5GSxLGH/OviYk+zwuTD+X2+C6HMefNSfmyes1OM+ZKU/ui7l2zVODedI8qbkasgfCkut+IcYeCEvO49FzyL6Ya9crNVj3q2kNrVdyTL+sOR1f8BiPLcM5DHmPsedweZ7DZXkOl+cxLsvjW96YxzhzjcvQNn/0+u/XDwlOk2twz9g5jTVGmTM8pw+i75OwhPHeU+OlUgYN86uJce+CvPk5Sm5w7wJzWwlXadfCfMzO158BfDHta/atwH3WX//z9Z9tpz43ddvNwJeu//3TA+LVtNzn70YN7vMX8R52rvev+2Ku/T6lwXof6RDz5o0azJvmzfM0mDfNm5orPzMz53qRJefxKGPEcx0fBj8zUzrGvGneNG+e1mDeNG8uh2skr/f1wF2BuwEPAm6nvYY+ded5D17//WbgN2lfN9t+EXgs8JZxQu5krn0C3folYp9A3n4Zo08gZr/YJy37ZBxev+LxtRKPfRKT74njqf1aybzeVdMbu6YJXO9Te1y3wZoma5rO02BNkzSGpX6+t3O59edyl5xzIPZcrjlTumbqPAnWPA3RYJ40T56nwTyp7pZW8+Q+o+V5jMvy+JbnMS7PfeHL8hwuz3O4LM/h8twXvuW+8PE1uDYz0zhflhoY93SPu7Z0X8x9xyqlDBrmV6/inu558jN7vuY8XH8N7umegXu66xjXcdyoYdz84jqO/hpcxzGn9yjW2mhO3NP9Rg3mTfPmeRrMm+ZNzdWQvWmXlAdrjw+D+xjUziHg3rR98qA5UHNn3jRvgnnzlAbzpnlzOZa2vhHcvzIi90SOp/aeoiVl7RPw+hWR1694vH7F5PUrHq9f8dS+fp29VvWCtoD49j2PPXz92D898TO+aP28V9MOpmzcTFsccXX953N3vu+J66/fvW/QkiTN2Io2P17UDUOSpNG9AHjpzteuHvhzZc/33wl4B/D/dvx9/2r9s35yQKySJEmSJEmSpHhWOJcqSVqGFePlvIv1z/r+nt/3l8Cv7Pn6R69/3s8c+L4Hrh9vev6+mrrMRb+J6+e0n0O7Cdu2p6wfew/tQvCH0xaRfzzthstXgZ8fGOOK886JbP0yVp9AuX5Zsaw+gdP9UrtP4Lx+mWOfgNevqXn9isnrVzxev+Lx+hVTlOvXO9e/K4Ku610lSVI5K5yzlSSpjxXmTknScqwYJ+9dMKzeCfbXPJ2qd4J5zSMM3RPylBXLmquZ+viC88YbnsPj2XeMPYfH4zlcXqZzGPIdY3NdeZ7DZZU6hzPWybgvvCRJ16xwbk6SND/u6S5JklZ4vytJkiRJGm6F95WSpGVYMV7Ou2DYGsehe7rDfOq5t7l/5bTsk5ii7Cl6zArX++xyT+Rpef2KyetXPF6/4vH6FVOU61fGtaqSJEmSJEmSpH5WWK80JfcBLs89VMtyH+DyvE6U5Tlcnudweea6sjyHyytxDkeqcQHrXCRJkiRJ87RinHnFC4bNKcLwz3n+58B7gYdsfe3K+vsef+B7so27ZBv/XkKfQL7x3L79Mpc+gbhjwEt4rXj9isnrVzxev+KxT+Ixp8Tk/LAkScu2Yrz1KpIkSZKkZVlRfw8EGF6v5Jj+aStcewtx5xfncoyjnsOQ7xhnm4/NdnzBc7g0z+HyzHVleQ6XV+IYZ61xGdrmJ66/fvcxgpUkSZIkhfYC4KU7X7t64M+Vned1+dzUbf9q/XN+cmCskiRJkiRJkiRJklTDCtdITu3UfPSbuH4++znAvfc87ynrx98DvIq2DvuuwMcDP71+7OcHxrjCOv1dXfolap9Avn4Zq0+gXL+ssE921e4T8Pq1j9evaWV4raywT3bV7hPw+rXLPplelNdK1vWukiSpnBXuHS9JkiRJutEKa56m5D6j5blXblmew+V5DpfnvvBleZ0oz3O4LM/h8twX/hr3hZckSZKkeXNPd0mStMJ1HJIkSZIkSZJ0yopxxlIvcH1jF+5fGY97IscTZU/RY1Ysq0/A61dEXr/i8foVk9eveLx+xRPl+nVsrerF+ntvPxTgsSd8/vqxf3Tom9fuAPzU+rlvAr4b+DfAy9fB/d76sb+9833ftv76LSd+/iWHi5n3/Rlygy1JUg2X9MtxTY0gJUkayY8Crzvj+/9n2nz4RR2e+3Xr576S6zdHOuQS7zslSRrbJeZXSdJyXGLekyRpbJc4lypJWoZLyua8C4bdi/4p8LI9X3/I+uf90IHv+7L141/X8/fV1Gcu+t7Ao4DfBd4AfPLWY9+y/jnvBT5x5/vuDLx2/fhDTvyOS8Y/J7L1y1h9AuP0yyX2CXTvlyn6BMbvlzn3CXj9morXr5i8fsXj9Sser18xRbl+vZ520WsEXde7gnO4kiSN4ZLx35dJkjRnl5g7JUnLcUm5vHfB8HHbfTVPp+qdYF7zCEP3hNx2iXM1JY8vOG8MnsNT2HeMPYfH4zlcXqZzGPIdY3NdeZ7DZZU6hzPWyZTeF16SpKguGf89kSRJEbmnuyRJy3KJ97uSpOW4xPtKSZLGdon3lZKkZbikbM67YNi96NA93WE+9dz7uH/lNOyTmKLsKbpxybj9Muc+AfdEnorXr5i8fsXj9Sser18xRbl+ZVyrKknSElzS7z2T9UqSpKwuMedJkjS2S8Yfh992wbC8bL2S+wCPyT1Uy3If4PK8TpTlOVye53B55rqyPIfLK3EOR6pxgX51Lpf0Oy8cB5Uk6UaXmE8lSRrbJeOPY2xcMDwnD/mc50+lHTf4lp2vX1l/3+MPfF+2cZdM499L6RPINZ47pF/m0icQcwx4Ka8Vr18xef2Kx+tXPPZJPOaUmJwfliRpOS4Z932sJElRXdIv51mrJEnSaZeUvae8YHheHlKvBI7p77pk3D6ey/GFmPOLMJ9jHPUchnzHONN8LOQ7vuA5XJrncHnmurI8h8srcYyz1rgMbfO3rb9+y4mff0m/880xVklSFpeY4yRJy/GjwOsGfm+fz039uvVzXwnco+PPv8ScLEnS2C4xv0qSluESc54kSWO7xPwqSVqGS8aved92wfBcufQ1kvvcG3gU8LvAG4BP3nn8W9Y/573AJ+48dmfgtevHH3Li91xinX6fGoFj/RK1TyBfv4zVJzBOv1xin0TrE/D6BV6/Ior2WrnEPonWJ+D1C3xPHFGU10rW9a7gOLMkSWO4ZPx7GEmS5uwS70UlSctwSdn7xQuG58ql1Dy5z2h57pVbludweZ7D5bkvfFleJ8rzHC7Lc7g894W/pvS+8JIkRXRJv/dDzstJkjJzT3dJkpblkvHnWSRJiugS7yklSRrbJeZXSdJyXFJuLPWC4blyKesbwf0rI3JP5Hii7Cm6cYl9Al6/IvL6FY/Xr5i8fsXj9SueKNevY2tVL9bfdzvAzR0C3XaX9d9/deJ576M9GP8Q+JL1n3cDv0rbqd8JPAD44z0//33r5x7z+x1i2PaGHs+VJKmmfwN88M7XPol244fvpX3Dte23JohJkqRS/opr95lDPAp4F/CTJ573BODpwCuAhwNv7fCzve+UJGl85ldJ0pKY9yRJGp9zqZKkpYia8w7N775z/fedDnzfpkD6JaNHNMxXA/8EuC/wO7QfxPNLO8/pOhcN8Gbgx4DfpP1QnmcCt60fe9v679cAL935vncCzwW+HPhU4IVHfkeJcyJKv3wGbX98Cm2fPA744T3PG6tPYJx+mXOffCPwaOBjgL8GXrT+2sv3PLdrv0zRJzB+v0TpkycAXwncuv7/7wD/gv3H3etXHU8Engx8F/A1O495/ZrOFeBJO197M3CfPc/1+jWd+wJPBR4B3BV4Ne3r5Bd2nuf1azqXwN/Y8/Vn0OacDa9f07kj7TXsi2lfM28EfmD9td0FolGuX3fm2r1gbV3Xu4JzuJIkjSHq+LUkSVGZOyVJSxI17+2reTpV7wRx5hGgW70THJ5HGLon5LY5z9WcW780xvGF+c4bj1GL5Dl83Lm1RZ7D/QypE/IcPu4K59X8RD2HIc4xhm41POa6YS7pVosDnsNDda2tKXUOZ6yTKb0vvCRJUUUdo5QkaWzu6S5J0rJ4vytJWhLvKyVJGp/3lZKkpYia84bu6Q5x6rlh2j3dwf0rT3FP92ui9Il7ul8TpU+g+9pL90Suwz3dW7X75Aru6b4RpU82zlkTuo/Xr/Nc4p7uG1H6xD3dz9NnT3dJkubOeiVJ0lKY8yRJGp/1SmWdU6/kPsDHda09AvdQHaprLZH7AA9z7n7s4HWij1r7scN8z2HoXrPjOTxcl9obMNcNdcl5tTSew8edWxczp/3YoV+di+OgkiSdz3wqSdL4Ms0rHptTvBn4Pto1N/+s5++KMu4C9echYbxzYg59Mse5y6H9EqVPIMZ8Z4Q+gTj9cu6+CdFeK3Pok23H5lTB69dUrnD+vglev8qo/Vna9sn1Lqk/1wvmlG3OD1/PfRAkSWpFHVeWJGls1ipJkjS+yPeUfeuVNqKM6UO32iPIM08MsY7vuZ/BEm1+cSPKMT7381SinsMQ5xjPrcZnI8rx3TaXep2NKMf4CvX3QIB5n8Pn7IEA5rpTLplXXcxGlOML3etcShzjrDUuQ9t8l/X3vvvEz3eMVZI0V+Y4SdKSHNpLvIuun8/5BODpwCuAhwNv7fjzzcmSJI3P/CpJWgpzniRJ4zO/SpKWwjWSZZ2zRnKfNwM/Bvwm7Z6QzwRu23r8beu/XwO8dOd73wk8F/hy4FOBFx75PXOt0z93f/ZDjvVL1D6BOP1y7v7s+0zxWplzn5y7dnMfr1/jGbLe8xCvX+e5wnlrRPfx+nW+c/Zn38fr13ku6bYOFXxPPKVz167uU7Jfsq53BceZJUkaQ+QxbEmSIvJeVJK0FJHvF5dS8+Q+o8OMUb/kXrnHnVuL5Dl8WoR94ed8Dm9zX/gyruC+8KWdsy981OMLcY7xJe4LX1LXuhpwX/htpfeFlyQpIuflJElL4p7ukiQtS+S6HEmSxuQ9pSRJ4zO/SpKWJOpY6lLWN4L7V07FPd2vF6Vf3NP9mih90nXdJXj9qsE93eP0yxXc030jSp/AeWtC9/H6dZ5Lxt/THbx+nWtRe7pfAFeB2/c89nfXj33l0B/OteD+ |