{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# 1D steady-state diffusion\n", "\n", "Diffusion problems study the *spreading of substances* (e.g., heat, particles, or pollutants) over time due to random motion, governed mathematically by the diffusion equation. At steady state, the solution depends on boundary conditions (e.g., fixed values or fluxes at edges), with applications spanning heat transfer, environmental pollutant dispersion, biological transport, and material science. The problems are analytically solvable in simple geometries but often require numerical methods for complex scenarios." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Problem setup\n", "\n", "The adiabatic rod has a length $L=0.5$ m, corss-sectional area $A=10\\times10^{-3}$ m$^2$, left-end temperature $T_A=100$ °C, right-end temperature $T_B=500$ °C, and thermal conductivity $k=1000$ W/(m$\\cdot$K). Determine the steady-state temperature distribution.\n", "\n", "![Adiabatic rod](images/diffusion1d.jpg)\n", "\n", "The mathematical model for one-dimensional steady-state diffusion problem is\n", "\n", "$$\\frac{\\mathrm{d}}{\\mathrm{d}x} \\left( k\\frac{\\mathrm{d}T}{\\mathrm{d}x} \\right) = 0.$$\n", "\n", "The boundary conditions are\n", "\n", "$$T|_{x=0} = T_A, \\quad T|_{x=L} = T_B.$$\n", "\n", "The analytical solution for the temperature distribution at steady state is\n", "\n", "$$T = T_A - \\frac{T_A - T_B}{L}x.$$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Solve problem\n", "\n", "### Define 1D grid\n", "\n", "Divide the rod evenly into 5 control volumes, as a result, $\\delta x=0.1$ m.\n", "\n", "\n", "\n", "![Grid](images/diffusion1d_grid.jpg)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Discretize 1D diffusion equation\n", "\n", "Make $k = k_w$ and $A = A_w = A_e$, derive the discrete equations for the internal and boundary nodes as following.\n", "\n", "**Internal Nodes 2, 3 and 4**\n", "\n", "The discrete equation satisfied by the internal nodes in the rod is\n", "\n", "$$a_P T_P = a_W T_W + a_E T_E,$$\n", "\n", "where\n", "\n", "$$a_W = a_E = \\frac{k}{\\delta x}A, \\quad a_P = a_W + a_E = \\frac{2k}{\\delta x}A.$$\n", "\n", "**Boundary Node 1**\n", "\n", "The discrete equation for the boundary node 1 is\n", "\n", "$$a_P T_P = a_W T_W + a_E T_E + S_u,$$\n", "\n", "where\n", "\n", "$$a_W = 0, \\quad a_E = \\frac{k}{\\delta x}A, \\quad a_P = a_W + a_E - S_p,$$\n", "\n", "$$S_P = -\\frac{2k}{\\delta x}A, \\quad S_u = \\frac{2k}{\\delta x}A \\cdot T_A$$\n", "\n", "**Boundary Node 5**\n", "\n", "The discrete equation for the boundary node 5 is\n", "\n", "$$a_P T_P = a_W T_W + a_E T_E + S_u,$$\n", "\n", "where\n", "\n", "$$a_W = \\frac{k}{\\delta x}A, \\quad a_E = 0, \\quad a_P = a_W + a_E - S_p,$$\n", "\n", "$$S_P = -\\frac{2k}{\\delta x}A, \\quad S_u = \\frac{2k}{\\delta x}A \\cdot T_B$$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Solve algebraic equations" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "# Step 0: import required libraries\n", "import numpy as np # for array operation\n", "from matplotlib import pyplot as plt # for plotting figures" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "# Step 1: parameter declarations\n", "nx = 5 # number of spatial grid points\n", "k = 1000 # thermal conductivity\n", "L = 0.5 # length of the rod\n", "A = 0.01 # cross-sectional area of the rod\n", "TA = 100 # temperature at x = 0\n", "TB = 500 # temperature at x = L\n", "nt = 30 # number of iterations\n", "dx = L / nx # spatial grid size\n", "coeff = k*A/dx # coefficient for west terms\n", "x = np.linspace(0.5*dx, L-0.5*dx, nx) # spatial grid points" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkYAAAG2CAYAAACap0noAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjAsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvlHJYcgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAS79JREFUeJzt3X18j/X////7a+c2NtvMThimFjFqOcsIZegE6Uy1Tij56S1lreUkCr1r3rwVlVK8FSnxfifVtyQjfAw5V84iImQzaTY7aZvt+P3htdell432mtfJNrfr5fK6tNfzeL6O1+PYgePe83gex2EyDMMQAAAA5ObqAgAAAKoLghEAAIAZwQgAAMCMYAQAAGBGMAIAADAjGAEAAJgRjAAAAMwIRgAAAGYEIwAAADOCEQAAgJnLg9Fvv/2mhx9+WMHBwfL19dX111+vbdu2WZYbhqGJEycqIiJCderUUY8ePbRnzx6rdRQWFurpp59WgwYN5Ofnp/79++v48ePO3hQAAFDDuTQYZWVlqUuXLvL09NQ333yjvXv36rXXXlP9+vUtfaZOnarXX39dM2fO1JYtWxQWFqZevXrp7Nmzlj6JiYlaunSpFi1apLS0NOXm5qpv374qKSlxwVYBAICayuTKh8iOGTNG69ev17p16ypcbhiGIiIilJiYqNGjR0s6PzoUGhqqKVOmaNiwYcrOzlZISIgWLFig+++/X5J04sQJRUZGatmyZerTp4/TtgcAANRsHq788i+//FJ9+vTRfffdp7Vr16pRo0YaPny4hg4dKkk6fPiwMjIy1Lt3b8tnvL291b17d23YsEHDhg3Ttm3bVFxcbNUnIiJCMTEx2rBhQ4XBqLCwUIWFhZb3paWl+uOPPxQcHCyTyeTALQYAAPZiGIbOnj2riIgIubnZ5ySYS4PRL7/8olmzZikpKUkvvPCCNm/erGeeeUbe3t569NFHlZGRIUkKDQ21+lxoaKh+/fVXSVJGRoa8vLwUGBhYrk/Z5y80efJkTZo0yQFbBAAAnO3YsWNq3LixXdbl0mBUWlqq9u3bKyUlRZIUGxurPXv2aNasWXr00Uct/S4cxTEM429Hdi7VZ+zYsUpKSrK8z87OVpMmTXTs2DH5+/tXdXMAAIAT5eTkKDIyUvXq1bPbOl0ajMLDw9WqVSurtmuvvVZLliyRJIWFhUk6PyoUHh5u6ZOZmWkZRQoLC1NRUZGysrKsRo0yMzMVFxdX4fd6e3vL29u7XLu/vz/BCACAGsae02BcelValy5dtH//fqu2AwcOqGnTppKkqKgohYWFKTU11bK8qKhIa9eutYSedu3aydPT06pPenq6du/efdFgBAAAUBGXjhg9++yziouLU0pKigYOHKjNmzdr9uzZmj17tqTzCTAxMVEpKSmKjo5WdHS0UlJS5Ovrq4SEBElSQECAhgwZoueee07BwcEKCgpScnKy2rRpo/j4eFduHgAAqGFcGow6dOigpUuXauzYsXr55ZcVFRWlGTNm6KGHHrL0GTVqlAoKCjR8+HBlZWWpU6dOWrFihdX5xOnTp8vDw0MDBw5UQUGBevbsqXnz5snd3d0VmwUAAGool97HqLrIyclRQECAsrOzLznHqKSkRMXFxU6sDLh8np6e/E8CgFqpssdvW7h0xKimMAxDGRkZOnPmjKtLAaqkfv36CgsL4z5dAPA3CEaVUBaKGjZsKF9fXw4uqDEMw1B+fr4yMzMlyerqTgBAeQSjv1FSUmIJRcHBwa4uB7BZnTp1JJ2/hUXDhg05rQYAl+DSy/VrgrI5Rb6+vi6uBKi6sj+/zJEDgEsjGFUSp89Qk/HnFwAqh2AEAABgRjCChclk0ueff37JPoMHD9aAAQMqvc4jR47IZDJp586dl1WbI/x1eytbZ48ePZSYmOjw2gAArsHk61pq8ODBOnPmzN8Gnb9KT0+3PG/uyJEjioqK0o4dO3T99ddb+rzxxhuqjbe+ioyMVHp6uho0aCBJWrNmjW6++WZlZWWpfv36ln6fffaZPD09XVQlAMDRCEawKHto76UEBAQ4oRLnc3d3r9T2BwUFOaEaAICrcCrtCtGjRw8988wzGjVqlIKCghQWFqaJEyda9fnrqaWoqChJUmxsrEwmk3r06CGp/Km05cuXq2vXrqpfv76Cg4PVt29fHTp0yKbaCgsLNWrUKEVGRsrb21vR0dGaO3euZfnatWvVsWNHeXt7Kzw8XGPGjNG5c+ds2raff/5Z3bp1k4+Pj1q1amX10GHJ+lTakSNHdPPNN0uSAgMDZTKZNHjwYMt3/fVUWlZWlh599FEFBgbK19dXt912m37++WfL8nnz5ql+/fr69ttvde2116pu3bq69dZblZ6ebtPvCADgHIwY2ajshnmucLk3l5w/f76SkpK0adMmbdy4UYMHD1aXLl3Uq1evcn03b96sjh07auXKlWrdurW8vLwqXGdeXp6SkpLUpk0b5eXl6aWXXtJdd92lnTt3ys2tcrn70Ucf1caNG/Xmm2/quuuu0+HDh/X7779Lkn777TfdfvvtGjx4sD788EP99NNPGjp0qHx8fKzCz6W2rbS0VHfffbcaNGig77//Xjk5OZecJxQZGaklS5bonnvu0f79++Xv72+5F9CFBg8erJ9//llffvml/P39NXr0aN1+++3au3ev5ZRbfn6+pk2bpgULFsjNzU0PP/ywkpOT9fHHH1fq9wMAcB6CkY3y8/NVt25dl3x3bm6u/Pz8qvz5tm3basKECZKk6OhozZw5U6tWraowGIWEhEiSgoODL3mK6Z577rF6P3fuXDVs2FB79+5VTEzM39Z04MAB/fe//1Vqaqri4+MlSc2bN7csf+eddxQZGamZM2fKZDKpZcuWOnHihEaPHq2XXnrJEr4utW0rV67Uvn37dOTIETVu3FiSlJKSottuu63Cmtzd3S2nzBo2bGg1x+ivygLR+vXrFRcXJ0n6+OOPFRkZqc8//1z33XefpPP3Dnr33Xd11VVXSZJGjBihl19++W9/NwAA5+NU2hWkbdu2Vu/Dw8Mtj4qoqkOHDikhIUHNmzeXv7+/5RTc0aNHK/X5nTt3yt3dXd27d69w+b59+9S5c2erkbIuXbooNzdXx48ft7Rdatv27dunJk2aWEKRJHXu3LlyG3gJ+/btk4eHhzp16mRpCw4OVosWLbRv3z5Lm6+vryUUXVgbAKB6YcTIRr6+vsrNzXXZd1+OC6+mMplMKi0tvax19uvXT5GRkZozZ44iIiJUWlqqmJgYFRUVVerzFztFVcYwjHKnD8uuivtr+6W2raKr6Oxxw8OLXZ13Yc0V1VYbr+wDgNqAYGQjk8l0WaezaoqyOUUlJSUX7XP69Gnt27dP7733nm666SZJUlpamk3f06ZNG5WWlmrt2rWWU2l/1apVKy1ZssQqbGzYsEH16tVTo0aNKvUdrVq10tGjR3XixAlFRERIkjZu3HjJz1Rm+1u1aqVz585p06ZNllNpp0+f1oEDB3TttddWqjYAQPXCqTRUqGHDhqpTp46WL1+ukydPKjs7u1yfwMBABQcHa/bs2Tp48KC+++47JSUl2fQ9zZo106BBg/T444/r888/1+HDh7VmzRr997//lSQNHz5cx44d09NPP62ffvpJX3zxhSZMmKCkpKRKT+6Oj49XixYt9Oijj+qHH37QunXrNG7cuEt+pmnTpjKZTPrqq6906tSpCkcJo6Ojdeedd2ro0KFKS0vTDz/8oIcffliNGjXSnXfeadPvAQBQPRCMUCEPDw+9+eabeu+99xQREVHhgd7NzU2LFi3Stm3bFBMTo2effVb//ve/bf6uWbNm6d5779Xw4cPVsmVLDR06VHl5eZKkRo0aadmyZdq8ebOuu+46PfnkkxoyZIjGjx9f6fW7ublp6dKlKiwsVMeOHfXEE0/o1VdfveRnGjVqpEmTJmnMmDEKDQ3ViBEjKuz3wQcfqF27durbt686d+4swzC0bNkybgIJADWUyWCyg3JychQQEKDs7Gz5+/tbLfvzzz91+PBhRUVFycfHx0UVApeHP8cAaqNLHb+rihEjAAAAM4IRAACAGcEIAADAjGAEAABgRjACAAAwIxgBAACYEYwAAADMCEYAAABmBCMAAAAzghEuW7NmzTRjxozLWseaNWtkMpl05swZu9R05MgRmUwm7dy50y7ru5C96nV0nQAA2xCMarkNGzbI3d1dt956q6tLsejRo4cSExOt2uLi4pSenq6AgADXFOUEgwcP1oABA6zaIiMjlZ6erpiYGNcUBQCwQjBykry8PJlMJplMJssDUp3h/fff19NPP620tDQdPXrUad9rKy8vL4WFhclkMrm6FKdyd3dXWFiYPDw8XF0KAEAEo1otLy9P//3vf/WPf/xDffv21bx586yWl50OWrVqldq3by9fX1/FxcVp//79lj6HDh3SnXfeqdDQUNWtW1cdOnTQypUrL/qdjz/+uPr27WvVdu7cOYWFhen999/X4MGDtXbtWr3xxhuWoHjkyJEKT02tX79e3bt3l6+vrwIDA9WnTx9lZWVJkpYvX66uXbuqfv36Cg4OVt++fXXo0CGbfj/vvPOOoqOj5ePjo9DQUN17772WZYWFhXrmmWfUsGFD+fj4qGvXrtqyZctF1zVx4kRdf/31Vm0zZsxQs2bNLMvnz5+vL774wrLda9asqfBU2tq1a9WxY0d5e3srPDxcY8aM0blz5yzLe/TooWeeeUajRo1SUFCQwsLCNHHiRJu2HQBQMYJRLbZ48WK1aNFCLVq00MMPP6wPPvhAhmGU6zdu3Di99tpr2rp1qzw8PPT4449bluXm5ur222/XypUrtWPHDvXp00f9+vW76OjTE088oeXLlys9Pd3StmzZMuXm5mrgwIF644031LlzZw0dOlTp6elKT09XZGRkufXs3LlTPXv2VOvWrbVx40alpaWpX79+KikpkXQ+9CUlJWnLli1atWqV3NzcdNddd6m0tLRSv5utW7fqmWee0csvv6z9+/dr+fLl6tatm2X5qFGjtGTJEs2fP1/bt2/X1VdfrT59+uiPP/6o1PovlJycrIEDB+rWW2+1bHdcXFy5fr/99ptuv/12dejQQT/88INmzZqluXPn6pVXXrHqN3/+fPn5+WnTpk2aOnWqXn75ZaWmplapNgDAXxgwsrOzDUlGdnZ2uWUFBQXG3r17jYKCApvXm5uba3mdPHnSkGRIMk6ePGm1zFHi4uKMGTNmGIZhGMXFxUaDBg2M1NRUy/LVq1cbkoyVK1da2r7++mtD0iW3t1WrVsZbb71led+0aVNj+vTpVsunTJlieT9gwABj8ODBlvfdu3c3Ro4cabXOslqysrIMwzCMBx980OjSpUultzUzM9OQZOzatcswDMM4fPiwIcnYsWNHhf2XLFli+Pv7Gzk5OeWW5ebmGp6ensbHH39saSsqKjIiIiKMqVOnVljvhAkTjOuuu85qPdOnTzeaNm1qeT9o0CDjzjvvtOpzYZ0vvPCC0aJFC6O0tNTS5+233zbq1q1rlJSUGIZx/vfXtWtXq/V06NDBGD16dIXbahiX9+cYAKqrSx2/q4oRIweqW7eu5RUaGmppLzstVfZyhP3792vz5s164IEHJEkeHh66//779f7775fr27ZtW8vP4eHhkqTMzExJ50dmRo0apVatWql+/fqqW7eufvrpp0vOV3riiSf0wQcfWNbz9ddfW41CVUbZiNHFHDp0SAkJCWrevLn8/f0VFRUlSZWeR9WrVy81bdpUzZs31yOPPKKPP/5Y+fn5lnUXFxerS5culv6enp7q2LGj9u3bZ9N22Grfvn3q3Lmz1VyrLl26KDc3V8ePH7e0/XWfSef3W9k+AwBUHTM+a6m5c+fq3LlzatSokaXNMAx5enoqKytLgYGBlnZPT0/Lz2UH5LJTUs8//7y+/fZbTZs2TVdffbXq1Kmje++9V0VFRRf97kcffVRjxozRxo0btXHjRjVr1kw33XSTTfXXqVPnksv79eunyMhIzZkzRxERESotLVVMTMwl6/qrevXqafv27VqzZo1WrFihl156SRMnTtSWLVsspxsvnAhuGMZFJ4e7ubmVO01ZXFxcqVr+7jsqquev+6xsWWVPIwIALo4RIwfKzc21vE6ePGlpP3nypNUyezt37pw+/PBDvfbaa9q5c6fl9cMPP6hp06b6+OOPK72udevWafDgwbrrrrvUpk0bhYWF6ciRI5f8THBwsAYMGKAPPvhAH3zwgR577DGr5V5eXpa5QhfTtm1brVq1qsJlp0+f1r59+zR+/Hj17NlT1157rWVSti08PDwUHx+vqVOn6scff9SRI0f03Xff6eqrr5aXl5fS0tIsfYuLi7V161Zde+21Fa4rJCREGRkZVuHownsTVWa7W7VqpQ0bNlitZ8OGDapXr55VyAUAOAYjRg7k5+d30faLLbOHr776SllZWRoyZEi5+wLde++9mjt3rkaMGFGpdV199dX67LPP1K9fP5lMJr344ouVGpl44okn1LdvX5WUlGjQoEFWy5o1a6ZNmzbpyJEjqlu3roKCgsp9fuzYsWrTpo2GDx+uJ598Ul5eXlq9erXuu+8+BQUFKTg4WLNnz1Z4eLiOHj2qMWPGVGp7ynz11Vf65Zdf1K1bNwUGBmrZsmUqLS1VixYt5Ofnp3/84x96/vnnFRQUpCZNmmjq1KnKz8/XkCFDKlxfjx49dOrUKU2dOlX33nuvli9frm+++Ub+/v5W2/3tt99q//79Cg4OrvCeTcOHD9eMGTP09NNPa8SIEdq/f78mTJigpKQkubnx/zEA4Gj8S1sLzZ07V/Hx8RUeeO+55x7t3LlT27dvr9S6pk+frsDAQMXFxalfv37q06ePbrjhhr/9XHx8vMLDw9WnTx9FRERYLUtOTpa7u7tatWqlkJCQCucFXXPNNVqxYoV++OEHdezYUZ07d9YXX3whDw8Pubm5adGiRdq2bZtiYmL07LPP6t///neltqdM/fr19dlnn+mWW27Rtddeq3fffVeffPKJWrduLUn617/+pXvuuUePPPKIbrjhBh08eFDffvut1SnIv7r22mv1zjvv6O2339Z1112nzZs3Kzk52arP0KFD1aJFC7Vv314hISFav359ufU0atRIy5Yt0+bNm3XdddfpySef1JAhQzR+/Hibtg8AUDUm48KJEVegnJwcBQQEKDs72+r/8CXpzz//1OHDhxUVFSUfH58qf0deXp5lonVubq5DR4yqg/z8fEVEROj999/X3Xff7epyrnj2+nMMANXJpY7fVcWpNCfx8/Or8B5CtU1paakyMjL02muvKSAgQP3793d1SQAAVBrBCHZ19OhRRUVFqXHjxpo3bx6PugAA1CgctWBXzZo1uyJGxgAAtROTrwEAAMwIRpXEKAhqMv78AkDlEIz+RtkdhsseFwHURGV/fi+8YzYAwBpzjP6Gu7u76tevb3kOla+v70UfCwFUN4ZhKD8/X5mZmapfv77c3d1dXRIAVGsEo0oICwuTJB7SiRqrfv36lj/HAICLIxhVgslkUnh4uBo2bFilB4MCruTp6clIEQBUEsHIBu7u7hxgAACoxZh8DQAAYEYwAgAAMHNpMJo4caJMJpPV668TRA3D0MSJExUREaE6deqoR48e2rNnj9U6CgsL9fTTT6tBgwby8/NT//79dfz4cWdvCgAAqAVcPmLUunVrpaenW167du2yLJs6dapef/11zZw5U1u2bFFYWJh69eqls2fPWvokJiZq6dKlWrRokdLS0pSbm6u+ffuqpKTEFZsDAABqMJdPvvbw8KjwMmLDMDRjxgyNGzdOd999tyRp/vz5Cg0N1cKFCzVs2DBlZ2dr7ty5WrBggeLj4yVJH330kSIjI7Vy5Ur16dPHqdsCAABqNpePGP3888+KiIhQVFSUHnjgAf3yyy+SpMOHDysjI0O9e/e29PX29lb37t21YcMGSdK2bdtUXFxs1SciIkIxMTGWPhUpLCxUTk6O1QsAAMClwahTp0768MMP9e2332rOnDnKyMhQXFycTp8+rYyMDElSaGio1WdCQ0MtyzIyMuTl5aXAwMCL9qnI5MmTFRAQYHlFRkbaecsAAEBN5NJgdNttt+mee+5RmzZtFB8fr6+//lrS+VNmZS58/IZhGH/7SI6/6zN27FhlZ2dbXseOHbuMrQAAALWFy0+l/ZWfn5/atGmjn3/+2TLv6MKRn8zMTMsoUlhYmIqKipSVlXXRPhXx9vaWv7+/1QsAAKBaBaPCwkLt27dP4eHhioqKUlhYmFJTUy3Li4qKtHbtWsXFxUmS2rVrJ09PT6s+6enp2r17t6UPAABAZbn0qrTk5GT169dPTZo0UWZmpl555RXl5ORo0KBBMplMSkxMVEpKiqKjoxUdHa2UlBT5+voqISFBkhQQEKAhQ4boueeeU3BwsIKCgpScnGw5NQcAAGALlwaj48eP68EHH9Tvv/+ukJAQ3Xjjjfr+++/VtGlTSdKoUaNUUFCg4cOHKysrS506ddKKFStUr149yzqmT58uDw8PDRw4UAUFBerZs6fmzZvHM80AAIDNTIZhGK4uwtVycnIUEBCg7Oxs5hsBAFBDOOL4Xa3mGAEAALgSwQgAAMCMYAQAAGBGMAIAADAjGAEAAJgRjAAAAMwIRgAAAGYEIwAAADOCEQAAgBnBCAAAwIxgBAAAYEYwAgAAMCMYAQBwBcnLy5PJZJLJZFJeXp6ry6l2CEYAAABmBCMAAAAzD1s/UFhYqM2bN+vIkSPKz89XSEiIYmNjFRUV5Yj6AADAZfrrKbOL/SxJfn5+Tqupuqp0MNqwYYPeeustff755yoqKlL9+vVVp04d/fHHHyosLFTz5s31//1//5+efPJJ1atXz5E1AwAAG9StW7fC9tDQUKv3hmE4o5xqrVKn0u68807de++9atSokb799ludPXtWp0+f1vHjx5Wfn6+ff/5Z48eP16pVq3TNNdcoNTXV0XUDAADYXaVGjHr37q3//e9/8vLyqnB58+bN1bx5cw0aNEh79uzRiRMn7FokAACoutzcXMvPeXl5lpGikydPcvrsAiaDcTPl5OQoICBA2dnZ8vf3d3U5AAA4TF5enuXUWm5ubo0ORo44ftt8VdqxY8d0/Phxy/vNmzcrMTFRs2fPtktBAAAArmJzMEpISNDq1aslSRkZGerVq5c2b96sF154QS+//LLdCwQAAHAWm4PR7t271bFjR0nSf//7X8XExGjDhg1auHCh5s2bZ+/6AACAHfn5+ckwDBmGUaNPozmKzcGouLhY3t7ekqSVK1eqf//+kqSWLVsqPT3dvtUBAAA4kc3BqHXr1nr33Xe1bt06paam6tZbb5UknThxQsHBwXYvEAAAwFlsDkZTpkzRe++9px49eujBBx/UddddJ0n68ssvLafYAAAAaqIqXa5fUlKinJwcBQYGWtqOHDkiX19fNWzY0K4FOgOX6wMAUPM44vht87PSJMnd3d0qFElSs2bN7FEPAACAy9gcjKKiomQymS66/JdffrmsggAAAFzF5mCUmJho9b64uFg7duzQ8uXL9fzzz9urLgAAAKezORiNHDmywva3335bW7duveyCAAAAXMXmq9Iu5rbbbtOSJUvstToAAACns1sw+vTTTxUUFGSv1QEAADidzafSYmNjrSZfG4ahjIwMnTp1Su+8845diwMAAHAmm4PRgAEDrN67ubkpJCREPXr0UMuWLe1VFwAAgNNV6QaPtQ03eAQAoOZx2Q0ec3JyLF+Yk5Nzyb4ECwAAUFNVKhgFBgYqPT1dDRs2VP369Su8waNhGDKZTCopKbF7kQAAAM5QqWD03XffWa44W716tUMLAgAAcBXmGIk5RgAA1ESOOH7bfB+j5cuXKy0tzfL+7bff1vXXX6+EhARlZWXZpSgAAABXsDkYPf/885YJ2Lt27VJSUpJuv/12/fLLL0pKSrJ7gQAAAM5i832MDh8+rFatWkmSlixZon79+iklJUXbt2/X7bffbvcCAQAAnMXmESMvLy/l5+dLklauXKnevXtLkoKCgv72Un4AAIDqzOYRo65duyopKUldunTR5s2btXjxYknSgQMH1LhxY7sXCAAA4Cw2jxjNnDlTHh4e+vTTTzVr1iw1atRIkvTNN9/o1ltvtXuBAAAAzsLl+uJyfQAAaiKXPhKksggWAACgpqpUMLrYY0D+ikeCAACAmq5SwYjHgAAAgCtBpYJR9+7dHV2HJk+erBdeeEEjR47UjBkzJJ0fhZo0aZJmz56trKwsderUSW+//bZat25t+VxhYaGSk5P1ySefqKCgQD179tQ777zDFXIAAMBmNl+VJknr1q3Tww8/rLi4OP3222+SpAULFlg9KsQWW7Zs0ezZs9W2bVur9qlTp+r111/XzJkztWXLFoWFhalXr146e/aspU9iYqKWLl2qRYsWKS0tTbm5uerbty+n9AAAgM1sDkZLlixRnz59VKdOHW3fvl2FhYWSpLNnzyolJcXmAnJzc/XQQw9pzpw5CgwMtLQbhqEZM2Zo3LhxuvvuuxUTE6P58+crPz9fCxculCRlZ2dr7ty5eu211xQfH6/Y2Fh99NFH2rVrl1auXGlzLQAA4MpmczB65ZVX9O6772rOnDny9PS0tMfFxWn79u02F/DUU0/pjjvuUHx8vFX74cOHlZGRYbmztiR5e3ure/fu2rBhgyRp27ZtKi4utuoTERGhmJgYS5+KFBYWKicnx+oFAABg852v9+/fr27dupVr9/f315kzZ2xa16JFi7R9+3Zt2bKl3LKMjAxJUmhoqFV7aGiofv31V0sfLy8vq5Gmsj5ln6/I5MmTNWnSJJtqBQAAtZ/NI0bh4eE6ePBgufa0tDQ1b9680us5duyYRo4cqY8++kg+Pj4X7XfhbQLKbgtwKX/XZ+zYscrOzra8jh07Vum6AQBA7WVzMBo2bJhGjhypTZs2yWQy6cSJE/r444+VnJys4cOHV3o927ZtU2Zmptq1aycPDw95eHho7dq1evPNN+Xh4WEZKbpw5CczM9OyLCwsTEVFRcrKyrpon4p4e3vL39/f6gUAAGBzMBo1apQGDBigm2++Wbm5uerWrZueeOIJDRs2TCNGjKj0enr27Kldu3Zp586dllf79u310EMPaefOnWrevLnCwsKUmppq+UxRUZHWrl2ruLg4SVK7du3k6elp1Sc9PV27d++29AEAAKgsm+cYSdKrr76qcePGae/evSotLVWrVq1Ut25dm9ZRr149xcTEWLX5+fkpODjY0p6YmKiUlBRFR0crOjpaKSkp8vX1VUJCgiQpICBAQ4YM0XPPPafg4GAFBQUpOTlZbdq0KTeZGwAA4O9UOhiVlJRoz549io6OVp06deTr66v27dtLkgoKCvTjjz8qJiZGbm5VujVShUaNGqWCggINHz7ccoPHFStWqF69epY+06dPl4eHhwYOHGi5weO8efPk7u5utzoAAMCVwWQYhlGZjvPmzdPMmTO1adOmcqGjpKREnTp1UmJioh5++GGHFOpIjng6LwAAcCxHHL8rPbwzd+5cJScnVzgS4+7urlGjRmn27Nl2KQoAAMAVKh2M9u/frxtvvPGiyzt06KB9+/bZpSgAQO2Sl5cnk8kkk8mkvLw8V5cDXFSlg1FeXt4l7xB99uxZ5efn26UoAAAAV6h0MIqOjr7kYzbS0tIUHR1tl6IAAABcodLBKCEhQePHj9ePP/5YbtkPP/ygl156yXIZPQAAeXl5Vq+/aweqg0pflVb2sNa0tDTFx8erZcuWMplM2rdvn1auXKkuXbooNTXV6sGyNQVXpQGA/f3d45vKVPIwBJTjiON3pYORdD4cTZ8+XQsXLtTPP/8swzB0zTXXKCEhQYmJifLy8rJLUc5GMAIA+yMYwdFcHoxqK4IRANjfhafPyp5hefLkSfn5+VmW/fVnwBaOOH5X6ZEgAAD8nYsFHj8/P8IQqi37Pb8DAACghiMYAQAAmHEqDQDgcH5+fkyyRo1Q5RGjoqIi7d+/X+fOnbNnPQAAAC5jczDKz8/XkCFD5Ovrq9atW+vo0aOSpGeeeUb/+te/7F4gAACAs9gcjMaOHasffvhBa9askY+Pj6U9Pj5eixcvtmtxAAAAzmTzHKPPP/9cixcv1o033mh1865WrVrp0KFDdi0OAADAmWweMTp16pQaNmxYrj0vL6/SdzkFAACojmwORh06dNDXX39teV8WhubMmaPOnTvbrzIAAAAns/lU2uTJk3Xrrbdq7969OnfunN544w3t2bNHGzdu1Nq1ax1RIwAAgFPYPGIUFxen9evXKz8/X1dddZVWrFih0NBQbdy4Ue3atXNEjQAAAE7BQ2TFQ2QBAKiJHHH8rtINHg8dOqTx48crISFBmZmZkqTly5drz549dikKAADAFWwORmvXrlWbNm20adMmLVmyRLm5uZKkH3/8URMmTLB7gQAAAM5iczAaM2aMXnnlFaWmpsrLy8vSfvPNN2vjxo12LQ4AAMCZbA5Gu3bt0l133VWuPSQkRKdPn7ZLUQAAAK5gczCqX7++0tPTy7Xv2LFDjRo1sktRAAAArmBzMEpISNDo0aOVkZEhk8mk0tJSrV+/XsnJyXr00UcdUSMAAIBT2ByMXn31VTVp0kSNGjVSbm6uWrVqpW7duikuLk7jx493RI0AAABOYdN9jAzD0NGjRxUSEqKMjAxt375dpaWlio2NVXR0tCPrdCjuYwQAQM3jiOO3TY8EMQxD0dHR2rNnj6Kjo9W8eXO7FAEAAFAd2HQqzc3NTdHR0Vx9BgAAaiWb5xhNnTpVzz//vHbv3u2IegAAAFzG5melBQYGKj8/X+fOnZOXl5fq1KljtfyPP/6wa4HOwBwjAABqHpfPMZKkGTNm2OWLAQAAqhubg9GgQYMcUQcAAIDL2RyMcnJyKmw3mUzy9va2en4aAABATWJzMKpfv75MJtNFlzdu3FiDBw/WhAkT5OZm89xuAAAAl7E5GM2bN0/jxo3T4MGD1bFjRxmGoS1btmj+/PkaP368Tp06pWnTpsnb21svvPCCI2oGAABwCJuD0fz58/Xaa69p4MCBlrb+/furTZs2eu+997Rq1So1adJEr776KsEIAADUKDaf69q4caNiY2PLtcfGxmrjxo2SpK5du+ro0aOXXx0AAIAT2RyMGjdurLlz55Zrnzt3riIjIyVJp0+fVmBg4OVXBwAA4EQ2n0qbNm2a7rvvPn3zzTfq0KGDTCaTtmzZop9++kmffvqpJGnLli26//777V4sAACAI9l852tJOnLkiN59910dOHBAhmGoZcuWGjZsmJo1a+aAEh2PO18DAFDzOOL4XaVgVNsQjAAAqHkccfyu0o2G1q1bp4cfflhxcXH67bffJEkLFixQWlqaXYoCAABwBZuD0ZIlS9SnTx/VqVNH27dvV2FhoSTp7NmzSklJsXuBAAAAzmJzMHrllVf07rvvas6cOfL09LS0x8XFafv27XYtDgAAwJlsDkb79+9Xt27dyrX7+/vrzJkz9qgJAADAJWwORuHh4Tp48GC59rS0NDVv3twuRQEAALiCzcFo2LBhGjlypDZt2iSTyaQTJ07o448/VnJysoYPH+6IGgEAAJzC5hs8jho1StnZ2br55pv1559/qlu3bvL29lZycrJGjBjhiBoBoEry8vJUt25dSVJubq78/PxcXBGA6q5Kl+u/+uqr+v3337V582Z9//33OnXqlP75z3/avJ5Zs2apbdu28vf3l7+/vzp37qxvvvnGstwwDE2cOFERERGqU6eOevTooT179lito7CwUE8//bQaNGggPz8/9e/fX8ePH6/KZgEAgCtclYKRJPn6+qp9+/bq2LGj5f/IbNW4cWP961//0tatW7V161bdcsstuvPOOy3hZ+rUqXr99dc1c+ZMbdmyRWFhYerVq5fOnj1rWUdiYqKWLl2qRYsWKS0tTbm5uerbt69KSkqqumkAAOAKVak7X999992VXuFnn312WQUFBQXp3//+tx5//HFFREQoMTFRo0ePlnR+dCg0NFRTpkzRsGHDlJ2drZCQEC1YsMDybLYTJ04oMjJSy5YtU58+fSr1ndz5Gqg98vLyrH4ODQ2VJJ08edLqVBqn1YCaz2V3vg4ICLC8/P39tWrVKm3dutWyfNu2bVq1apUCAgKqXEhJSYkWLVqkvLw8de7cWYcPH1ZGRoZ69+5t6ePt7a3u3btrw4YNlu8tLi626hMREaGYmBhLn4oUFhYqJyfH6gWgdqhbt67lVRaKJCk0NNRqGQBUpFKTrz/44APLz6NHj9bAgQP17rvvyt3dXdL5UDN8+PAqpbVdu3apc+fO+vPPP1W3bl0tXbpUrVq1sgSbv/7DVvb+119/lSRlZGTIy8tLgYGB5fpkZGRc9DsnT56sSZMm2VwrAACo3WyeY/T+++8rOTnZEookyd3dXUlJSXr//fdtLqBFixbauXOnvv/+e/3jH//QoEGDtHfvXstyk8lk1d8wjHJtF/q7PmPHjlV2drbldezYMZvrBlA95ebmWl4nT560tJ88edJqGQBUxObL9c+dO6d9+/apRYsWVu379u1TaWmpzQV4eXnp6quvliS1b99eW7Zs0RtvvGGZV5SRkaHw8HBL/8zMTMsoUlhYmIqKipSVlWU1apSZmam4uLiLfqe3t7e8vb1trhVA9XexuUN+fn7MKwLwt2weMXrsscf0+OOPa9q0aUpLS1NaWpqmTZumJ554Qo899thlF2QYhgoLCxUVFaWwsDClpqZalhUVFWnt2rWW0NOuXTt5enpa9UlPT9fu3bsvGYwAAAAqYvOI0bRp0xQWFqbp06crPT1d0vnHhIwaNUrPPfecTet64YUXdNtttykyMlJnz57VokWLtGbNGi1fvlwmk0mJiYlKSUlRdHS0oqOjlZKSIl9fXyUkJEg6Pyl8yJAheu655xQcHKygoCAlJyerTZs2io+Pt3XTAADAFc7mYOTm5qZRo0Zp1KhRlqu5qnqJ3MmTJ/XII48oPT1dAQEBatu2rZYvX65evXpJOn+X7YKCAg0fPlxZWVnq1KmTVqxYoXr16lnWMX36dHl4eGjgwIEqKChQz549NW/ePKs5UACuTH5+fqrEHUkAwKJS9zGq7biPEQAANY/L7mN06623XvK+QGXOnj2rKVOm6O23377swgAAAJytUqfS7rvvPg0cOFD16tVT//791b59e0VERMjHx0dZWVnau3ev0tLStGzZMvXt21f//ve/HV03AACA3VX6VFpRUZE+/fRTLV68WOvWrdOZM2fOr8BkUqtWrdSnTx8NHTq03GX8NQGn0gAAqHkccfyu8hyj7OxsFRQUKDg4WJ6ennYpxlUIRgAA1DyOOH7bfFVambJnpwEAANQWNt/gEQAAoLYiGAEAAJgRjAAAAMwIRgAAAGZVCkZnzpzRf/7zH40dO1Z//PGHJGn79u367bff7FocAACAM9l8VdqPP/6o+Ph4BQQE6MiRIxo6dKiCgoK0dOlS/frrr/rwww8dUScAAIDD2TxilJSUpMGDB+vnn3+Wj4+Ppf22227T//3f/9m1OAAAAGeyORht2bJFw4YNK9feqFEjZWRk2KUoAAAAV7A5GPn4+CgnJ6dc+/79+xUSEmKXogAAAFzB5mB055136uWXX1ZxcbGk889KO3r0qMaMGaN77rnH7gUCAAA4i83BaNq0aTp16pQaNmyogoICde/eXVdffbXq1aunV1991RE1AgAAOIXNV6X5+/srLS1N3333nbZv367S0lLdcMMNio+Pd0R9AAAATmNTMDp37px8fHy0c+dO3XLLLbrlllscVRcAAIDT2XQqzcPDQ02bNlVJSYmj6gEAAHAZm+cYjR8/3uqO1wAAALWFzXOM3nzzTR08eFARERFq2rSp/Pz8rJZv377dbsUBAAA4k83BaMCAAQ4oAwAAwPVMhmEYri7C1XJychQQEKDs7Gz5+/u7uhwAAFAJjjh+2zzHCAAAoLay+VSam5ubTCbTRZdzxRoAAKipbA5GS5cutXpfXFysHTt2aP78+Zo0aZLdCgMAAHA2u80xWrhwoRYvXqwvvvjCHqtzKuYYAQBQ81TrOUadOnXSypUr7bU6AAAAp7NLMCooKNBbb72lxo0b22N1AAAALmHzHKPAwECrydeGYejs2bPy9fXVRx99ZNfiAAAAnMnmYDR9+nSrYOTm5qaQkBB16tRJgYGBdi0OAADAmWwORrfccosiIyMrvGT/6NGjatKkiV0KAwAAcDab5xhFRUXp1KlT5dpPnz6tqKgouxQFAADgCjYHo4td3Z+bmysfH5/LLggAAMBVKn0qLSkpSZJkMpn00ksvydfX17KspKREmzZt0vXXX2/3AgEAAJyl0sFox44dks6PGO3atUteXl6WZV5eXrruuuuUnJxs/woBAACcpNLBaPXq1ZKkxx57TG+88QZ3iAYAALWOzXOMPvjgA0IRcAl5eXkymUwymUzKy8tzdTkAABvYfLm+JG3ZskX/+9//dPToURUVFVkt++yzz+xSGAAAgLPZPGK0aNEidenSRXv37tXSpUtVXFysvXv36rvvvlNAQIAjagQAAHAKm0eMUlJSNH36dD311FOqV6+e3njjDUVFRWnYsGEKDw93RI1AtffXU2YX+1mS/Pz8nFYTAMB2No8YHTp0SHfccYckydvb2zKf4tlnn9Xs2bPtXiBQE9StW9fyCg0NtbSHhoZaLQMAVG82B6OgoCCdPXtWktSoUSPt3r1bknTmzBnl5+fbtzoAAAAnsvlU2k033aTU1FS1adNGAwcO1MiRI/Xdd98pNTVVPXv2dESNQLWXm5tr+TkvL88yanTy5ElOnwFADWJzMJo5c6b+/PNPSdLYsWPl6emptLQ03X333XrxxRftXiBQE1ws/Pj5+RGMAKAGMRkXe/hZBc6dO6ePP/5Yffr0UVhYmCPrcqqcnBwFBAQoOzubezThsuXl5VnmE+Xm5hKMAMBBHHH8tmmOkYeHh/7xj3+osLDQLl8OAABQndg8+bpTp06W56YBKM/Pz0+GYcgwDEaLAKCGsXmO0fDhw/Xcc8/p+PHjateuXbl/+Nu2bWu34gAAAJzJpjlGkuTmVn6QyWQyyTAMmUwmlZSU2K04Z2GOEQAANY/L5xhJ0uHDh8u9fvnlF8t/bTF58mR16NBB9erVU8OGDTVgwADt37/fqo9hGJo4caIiIiJUp04d9ejRQ3v27LHqU1hYqKeffloNGjSQn5+f+vfvr+PHj9u6aQAA4ApnczBq2rTpJV+2WLt2rZ566il9//33Sk1N1blz59S7d2+rxyhMnTpVr7/+umbOnKktW7YoLCxMvXr1stxkUpISExO1dOlSLVq0SGlpacrNzVXfvn1r5OgVAABwHZtPpUnSggUL9O677+rw4cPauHGjmjZtqhkzZigqKkp33nlnlYs5deqUGjZsqLVr16pbt24yDEMRERFKTEzU6NGjJZ0fHQoNDdWUKVM0bNgwZWdnKyQkRAsWLND9998vSTpx4oQiIyO1bNky9enT52+/l1NpAADUPNXiVNqsWbOUlJSk22+/XWfOnLGMytSvX18zZsy4rGKys7MlnX/siHT+tF1GRoZ69+5t6ePt7a3u3btrw4YNkqRt27apuLjYqk9ERIRiYmIsfQAAACrD5mD01ltvac6cORo3bpzc3d0t7e3bt9euXbuqXIhhGEpKSlLXrl0VExMjScrIyJAkq4dylr0vW5aRkSEvLy8FBgZetM+FCgsLlZOTY/UCAACo0uTr2NjYcu3e3t5Wc4NsNWLECP3444/65JNPyi0zmUxW78uugLuUS/WZPHmyAgICLK/IyMgq1w0AAGoPm4NRVFSUdu7cWa79m2++UatWrapUxNNPP60vv/xSq1evVuPGjS3tZY8duXDkJzMz0zKKFBYWpqKiImVlZV20z4XGjh2r7Oxsy+vYsWNVqhsAANQuNgej559/Xk899ZQWL14swzC0efNmvfrqq3rhhRf0/PPP27QuwzA0YsQIffbZZ/ruu+8UFRVltTwqKkphYWFKTU21tBUVFWnt2rWKi4uTJLVr106enp5WfdLT07V7925Lnwt5e3vL39/f6gUAAGDzna8fe+wxnTt3TqNGjVJ+fr4SEhLUqFEjvfHGG3rggQdsWtdTTz2lhQsX6osvvlC9evUsI0MBAQGqU6eOTCaTEhMTlZKSoujoaEVHRyslJUW+vr5KSEiw9B0yZIiee+45BQcHKygoSMnJyWrTpo3i4+Nt3TwAAHAFq9Ll+mV+//13lZaWqmHDhlX78ovMAfrggw80ePBgSedHlSZNmqT33ntPWVlZ6tSpk95++23LBG1J+vPPP/X8889r4cKFKigoUM+ePfXOO+9Ueu4Ql+sDAFDzOOL4XeVglJmZqf3798tkMqlFixYKCQmxS0GuQDACAKDmqRb3McrJydEjjzyiiIgIde/eXd26dVNERIQefvhhy32IAAAAaiKbg9ETTzyhTZs26euvv9aZM2eUnZ2tr776Slu3btXQoUMdUSMAAIBT2Hwqzc/PT99++626du1q1b5u3Trdeuutl3UvI1fhVBoAADVPtTiVFhwcrICAgHLtAQEB5e4+DQAAUJPYHIzGjx+vpKQkpaenW9oyMjL0/PPP68UXX7RrcQAAAM5k86m02NhYHTx4UIWFhWrSpIkk6ejRo/L29lZ0dLRV3+3bt9uvUgfiVBoAADWPI47fNt/gccCAAXb5YgAAgOrmsm7wWFswYgQAQM1TLUaM/io3N1elpaVWbQQLAABQU9k8+frw4cO644475OfnZ7kSLTAwUPXr1+eqNAAAUKPZPGL00EMPSZLef/99hYaGXvR5ZwAAADWNzcHoxx9/1LZt29SiRQtH1AMAAOAyNp9K69Chg44dO+aIWgAAAFzK5hGj//znP3ryySf122+/KSYmRp6enlbL27Zta7fiAAAAnMnmYHTq1CkdOnRIjz32mKXNZDLJMAyZTCaVlJTYtUAAAABnsTkYPf7444qNjdUnn3zC5GsAAFCr2ByMfv31V3355Ze6+uqrHVEPAACAy9g8+fqWW27RDz/84IhaAAAAXMrmEaN+/frp2Wef1a5du9SmTZtyk6/79+9vt+IAAACcyeZnpbm5XXyQqaZOvuZZaQAA1DzV4llpFz4bDQAAoLaweY7RX/3555/2qgMAAMDlbA5GJSUl+uc//6lGjRqpbt26+uWXXyRJL774oubOnWv3AgEAAJzF5mD06quvat68eZo6daq8vLws7W3atNF//vMfuxYHAADgTDYHow8//FCzZ8/WQw89JHd3d0t727Zt9dNPP9m1OAAAAGeyORj99ttvFd7csbS0VMXFxXYpCrbLy8uTyWSSyWRSXl6eq8sBAKBGsjkYtW7dWuvWrSvX/r///U+xsbF2KQoAAMAVKn25/uOPP6433nhDEyZM0COPPKLffvtNpaWl+uyzz7R//359+OGH+uqrrxxZKwAAgENVesRo/vz5KigoUL9+/bR48WItW7ZMJpNJL730kvbt26f/9//+n3r16uXIWnGBvLw8q9fftQMAgEur9J2v3dzclJGRoYYNGzq6JqerqXe+NplMlepn483NAQCoERxx/LZpjlFlD8QAAAA1kU2PBLnmmmv+Nhz98ccfl1UQKi83N9fyc15enkJDQyVJJ0+elJ+fn6vKAgCgxrIpGE2aNEkBAQGOqgU2ulj48fPzIxgBAFAFNgWjBx54oFbOMQIAAJBsmGPE/CIAAFDbVXrEiCubqjc/Pz/2EQAAl6nSwai0tNSRdQAAALiczY8EAQAAqK0IRgAAAGYEIwAAADOCEQAAgBnBCAAAwIxgBAAAYEYwAgAAMCMYAQAAmBGMAAAAzAhGAAAAZgQjAAAAM4IRAACAGcEIAADAjGAEAABg5tJg9H//93/q16+fIiIiZDKZ9Pnnn1stNwxDEydOVEREhOrUqaMePXpoz549Vn0KCwv19NNPq0GDBvLz81P//v11/PhxJ24FAACoLVwajPLy8nTddddp5syZFS6fOnWqXn/9dc2cOVNbtmxRWFiYevXqpbNnz1r6JCYmaunSpVq0aJHS0tKUm5urvn37qqSkxFmbAQAAagmTYRiGq4uQJJPJpKVLl2rAgAGSzo8WRUREKDExUaNHj5Z0fnQoNDRUU6ZM0bBhw5Sdna2QkBAtWLBA999/vyTpxIkTioyM1LJly9SnT59KfXdOTo4CAgKUnZ0tf39/h2yfoxmGofz8fFeXAQCoQXx9fWUymVxdRpU54vjtYZe1OMDhw4eVkZGh3r17W9q8vb3VvXt3bdiwQcOGDdO2bdtUXFxs1SciIkIxMTHasGHDRYNRYWGhCgsLLe9zcnIctyFOkp+fr7p167q6DABADZKbmys/Pz9Xl1GtVNvJ1xkZGZKk0NBQq/bQ0FDLsoyMDHl5eSkwMPCifSoyefJkBQQEWF6RkZF2rh4AANRE1XbEqMyFQ3yGYfztsN/f9Rk7dqySkpIs73Nycmp8OPL19VVubq6rywAA1CC+vr6uLqHaqbbBKCwsTNL5UaHw8HBLe2ZmpmUUKSwsTEVFRcrKyrIaNcrMzFRcXNxF1+3t7S1vb28HVe4aJpOJ4VAAAC5TtT2VFhUVpbCwMKWmplraioqKtHbtWkvoadeunTw9Pa36pKena/fu3ZcMRgAAABVx6YhRbm6uDh48aHl/+PBh7dy5U0FBQWrSpIkSExOVkpKi6OhoRUdHKyUlRb6+vkpISJAkBQQEaMiQIXruuecUHBysoKAgJScnq02bNoqPj3fVZgEAgBrKpcFo69atuvnmmy3vy+b9DBo0SPPmzdOoUaNUUFCg4cOHKysrS506ddKKFStUr149y2emT58uDw8PDRw4UAUFBerZs6fmzZsnd3d3p28PAACo2arNfYxcqTbcxwgAgCuNI47f1XaOEQAAgLMRjAAAAMwIRgAAAGYEIwAAADOCEQAAgBnBCAAAwIxgBAAAYEYwAgAAMCMYAQAAmBGMAAAAzAhGAAAAZgQjAAAAM4IRAACAGcEIAADAjGAEAABgRjACAAAwIxgBAACYEYwAAADMCEYAAABmBCMAAAAzghEAAIAZwQgAAMCMYAQAAGBGMAIAADAjGAEAAJgRjAAAAMwIRgAAAGYEIwAAADOCEQAAgBnBCAAAwIxgBAAAYEYwAgAAMCMYAQAAmBGMAAAAzAhGAAAAZgQjAAAAM4IRAACAGcEIAADAjGAEAABgRjACAAAwIxgBAACYEYwAAADMCEYAAABmBCMAAAAzghEAAIAZwQgAAMCMYAQAAGBGMAIAADAjGAEAAJgRjAAAAMwIRgAAAGa1Jhi98847ioqKko+Pj9q1a6d169a5uiQAAFDD1IpgtHjxYiUmJmrcuHHasWOHbrrpJt122206evSoq0sDAAA1iMkwDMPVRVyuTp066YYbbtCsWbMsbddee60GDBigyZMn/+3nc3JyFBAQoOzsbPn7+zuyVAAAYCeOOH572GUtLlRUVKRt27ZpzJgxVu29e/fWhg0bKvxMYWGhCgsLLe+zs7Mlnf8FAwCAmqHsuG3PMZ4aH4x+//13lZSUKDQ01Ko9NDRUGRkZFX5m8uTJmjRpUrn2yMhIh9QIAAAc5/Tp0woICLDLump8MCpjMpms3huGUa6tzNixY5WUlGR5f+bMGTVt2lRHjx612y8WVZOTk6PIyEgdO3aM05ouxr6oXtgf1Qf7ovrIzs5WkyZNFBQUZLd11vhg1KBBA7m7u5cbHcrMzCw3ilTG29tb3t7e5doDAgL4Q15N+Pv7sy+qCfZF9cL+qD7YF9WHm5v9riWr8VeleXl5qV27dkpNTbVqT01NVVxcnIuqAgAANVGNHzGSpKSkJD3yyCNq3769OnfurNmzZ+vo0aN68sknXV0aAACoQWpFMLr//vt1+vRpvfzyy0pPT1dMTIyWLVumpk2bVurz3t7emjBhQoWn1+Bc7Ivqg31RvbA/qg/2RfXhiH1RK+5jBAAAYA81fo4RAACAvRCMAAAAzAhGAAAAZgQjAAAAsysmGL3zzjuKioqSj4+P2rVrp3Xr1l2y/9q1a9WuXTv5+PioefPmevfdd51Uae1ny75IT09XQkKCWrRoITc3NyUmJjqv0CuALfvis88+U69evRQSEiJ/f3917txZ3377rROrrd1s2RdpaWnq0qWLgoODVadOHbVs2VLTp093YrW1n63HjDLr16+Xh4eHrr/+escWeAWxZV+sWbNGJpOp3Ounn36q/BcaV4BFixYZnp6expw5c4y9e/caI0eONPz8/Ixff/21wv6//PKL4evra4wcOdLYu3evMWfOHMPT09P49NNPnVx57WPrvjh8+LDxzDPPGPPnzzeuv/56Y+TIkc4tuBazdV+MHDnSmDJlirF582bjwIEDxtixYw1PT09j+/btTq689rF1X2zfvt1YuHChsXv3buPw4cPGggULDF9fX+O9995zcuW1k637o8yZM2eM5s2bG7179zauu+465xRby9m6L1avXm1IMvbv32+kp6dbXufOnav0d14Rwahjx47Gk08+adXWsmVLY8yYMRX2HzVqlNGyZUurtmHDhhk33nijw2q8Uti6L/6qe/fuBCM7upx9UaZVq1bGpEmT7F3aFcce++Kuu+4yHn74YXuXdkWq6v64//77jfHjxxsTJkwgGNmJrfuiLBhlZWVV+Ttr/am0oqIibdu2Tb1797Zq7927tzZs2FDhZzZu3Fiuf58+fbR161YVFxc7rNbarir7Ao5hj31RWlqqs2fP2vXhjVcie+yLHTt2aMOGDerevbsjSryiVHV/fPDBBzp06JAmTJjg6BKvGJfzdyM2Nlbh4eHq2bOnVq9ebdP31oo7X1/K77//rpKSknIPlA0NDS334NkyGRkZFfY/d+6cfv/9d4WHhzus3tqsKvsCjmGPffHaa68pLy9PAwcOdESJV4zL2ReNGzfWqVOndO7cOU2cOFFPPPGEI0u9IlRlf/z8888aM2aM1q1bJw+PWn9YdZqq7Ivw8HDNnj1b7dq1U2FhoRYsWKCePXtqzZo16tatW6W+94rZgyaTyeq9YRjl2v6uf0XtsJ2t+wKOU9V98cknn2jixIn64osv1LBhQ0eVd0Wpyr5Yt26dcnNz9f3332vMmDG6+uqr9eCDDzqyzCtGZfdHSUmJEhISNGnSJF1zzTXOKu+KYsvfjRYtWqhFixaW9507d9axY8c0bdo0glGZBg0ayN3dvVy6zMzMLJdCy4SFhVXY38PDQ8HBwQ6rtbaryr6AY1zOvli8eLGGDBmi//3vf4qPj3dkmVeEy9kXUVFRkqQ2bdro5MmTmjhxIsHoMtm6P86ePautW7dqx44dGjFihKTzp5kNw5CHh4dWrFihW265xSm11zb2OmbceOON+uijjyrdv9bPMfLy8lK7du2Umppq1Z6amqq4uLgKP9O5c+dy/VesWKH27dvL09PTYbXWdlXZF3CMqu6LTz75RIMHD9bChQt1xx13OLrMK4K9/l4YhqHCwkJ7l3fFsXV/+Pv7a9euXdq5c6fl9eSTT6pFixbauXOnOnXq5KzSax17/d3YsWOHbVNgqjxtuwYpu9xv7ty5xt69e43ExETDz8/POHLkiGEYhjFmzBjjkUcesfQvu1z/2WefNfbu3WvMnTuXy/XtxNZ9YRiGsWPHDmPHjh1Gu3btjISEBGPHjh3Gnj17XFF+rWLrvli4cKHh4eFhvP3221aXwZ45c8ZVm1Br2LovZs6caXz55ZfGgQMHjAMHDhjvv/++4e/vb4wbN85Vm1CrVOXfqb/iqjT7sXVfTJ8+3Vi6dKlx4MABY/fu3caYMWMMScaSJUsq/Z1XRDAyDMN4++23jaZNmxpeXl7GDTfcYKxdu9aybNCgQUb37t2t+q9Zs8aIjY01vLy8jGbNmhmzZs1ycsW1l637QlK5V9OmTZ1bdC1ly77o3r17hfti0KBBzi+8FrJlX7z55ptG69atDV9fX8Pf39+IjY013nnnHaOkpMQFlddOtv479VcEI/uyZV9MmTLFuOqqqwwfHx8jMDDQ6Nq1q/H111/b9H0mwzDPKgYAALjC1fo5RgAAAJVFMAIAADAjGAEAAJgRjAAAAMwIRgAAAGYEIwAAADOCEQAAgBnBCMAVY/DgwRowYECVPtutWzctXLjwsr5/5syZ6t+//2WtA4BjEYwA2NXlhA97OXLkiEwmk3bu3GmX9X311VfKyMjQAw88cFnrGTp0qLZs2aK0tDS71AXA/ghGAPA33nzzTT322GNyc7u8fzK9vb2VkJCgt956y06VAbA3ghEAp9q7d69uv/121a1bV6GhoXrkkUf0+++/W5b36NFDzzzzjEaNGqWgoCCFhYVp4sSJVuv46aef1LVrV/n4+KhVq1ZauXKlTCaTPv/8c0lSVFSUJCk2NlYmk0k9evSw+vy0adMUHh6u4OBgPfXUUyouLr5ovb///rtWrlxZ7hSYyWTSe++9p759+8rX11fXXnutNm7cqIMHD6pHjx7y8/NT586ddejQIavP9e/fX59//rkKCgps/M0BcAaCEQCnSU9PV/fu3XX99ddr69atWr58uU6ePKmBAwda9Zs/f778/Py0adMmTZ06VS+//LJSU1MlSaWlpRowYIB8fX21adMmzZ49W+PGjbP6/ObNmyVJK1euVHp6uj777DPLstWrV+vQoUNavXq15s+fr3nz5mnevHkXrTktLc0SfC70z3/+U48++qh27typli1bKiEhQcOGDdPYsWO1detWSdKIESOsPtO+fXsVFxdbagRQvXi4ugAAV45Zs2bphhtuUEpKiqXt/fffV2RkpA4cOKBrrrlGktS2bVtNmDBBkhQdHa2ZM2dq1apV6tWrl1asWKFDhw5pzZo1CgsLkyS9+uqr6tWrl2WdISEhkqTg4GBLnzKBgYGaOXOm3N3d1bJlS91xxx1atWqVhg4dWmHNR44cUWhoaIWn0R577DFLqBs9erQ6d+6sF198UX369JEkjRw5Uo899pjVZ/z8/FS/fn0dOXJE3bt3r/wvD4BTMGIEwGm2bdum1atXq27dupZXy5YtJcnqlFPbtm2tPhceHq7MzExJ0v79+xUZGWkVeDp27FjpGlq3bi13d/cK112RgoIC+fj4VLjsr3WGhoZKktq0aWPV9ueffyonJ8fqc3Xq1FF+fn6lawbgPIwYAXCa0tJS9evXT1OmTCm3LDw83PKzp6en1TKTyaTS0lJJkmEYMplMVa7hUuuuSIMGDZSVlfW36yqrqaK2C9f/xx9/WEa1AFQvBCMATnPDDTdoyZIlatasmTw8qvbPT8uWLXX06FGdPHnSMkqzZcsWqz5eXl6SpJKSkssrWOcncGdkZCgrK0uBgYGXvb5Dhw7pzz//VGxs7GWvC4D9cSoNgN1lZ2dr586dVq+jR4/qqaee0h9//KEHH3xQmzdv1i+//KIVK1bo8ccfr3SI6dWrl6666ioNGjRIP/74o9avX2+ZfF02QtOwYUPVqVPHMrk7Ozu7ytsSGxurkJAQrV+/vsrr+Kt169apefPmuuqqq+yyPgD2RTACYHdr1qxRbGys1eull15SRESE1q9fr5KSEvXp00cxMTEaOXKkAgICKn2PIHd3d33++efKzc1Vhw4d9MQTT2j8+PGSZJkL5OHhoTfffFPvvfeeIiIidOedd1Z5W9zd3fX444/r448/rvI6/uqTTz656ERvAK5nMgzDcHURAHA51q9fr65du+rgwYMOGYk5efKkWrdurW3btqlp06ZVXs/u3bvVs2dPHThwQAEBAXasEIC9EIwA1DhLly5V3bp1FR0drYMHD2rkyJEKDAx06KM2vvjiCwUFBemmm26q8jpWrFghwzAsl/MDqH4IRgBqnA8//FD//Oc/dezYMTVo0EDx8fF67bXXFBwc7OrSANRwBCMAAAAzJl8DAACYEYwAAADMCEYAAABmBCMAAAAzghEAAIAZwQgAAMCMYAQAAGBGMAIAADAjGAEAAJj9/6xw3GJ7osIMAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Step 2.1: set initial condition\n", "T = np.ones(nx) * TA # a numpy array with nx elements all equal to TA\n", "\n", "# Step 2.2: visualize initial condition and analytical solution\n", "plt.figure() # create a new figure\n", "plt.plot(x, T, '-k', label='Initial condition')\n", "plt.scatter(x, TA-(TA-TB)/L*x, s=50, c='k', marker='+', label='Analytical solution')\n", "plt.xlabel('Length (m)')\n", "plt.ylabel('Temperature (degree Celsius)')\n", "ax = plt.gca()\n", "ax.set_xlim(0, 0.5)\n", "ax.set_ylim(0, 600)\n", "plt.legend(loc='upper left')\n", "plt.show() # show the figure" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Iteration 10: Tdiff = 26.337\n", "Iteration 20: Tdiff = 3.468\n", "Iteration 30: Tdiff = 0.457\n" ] } ], "source": [ "# Step 3: finite volume calculations\n", "Tn = np.ones(nx) # placeholder array to advance the solution\n", "for n in range(nt):\n", " Tn = T.copy() # copy the existing (old) values of T into Tn\n", " \n", " # left boundary\n", " T[0] = (coeff*Tn[1] + 2 * coeff*TA) / (coeff + 2*coeff)\n", "\n", " # internal nodes\n", " for i in range(1, nx-1): # skip the first and the last elements\n", " T[i] = (coeff*Tn[i-1] + coeff*Tn[i+1]) / (2*coeff)\n", "\n", " # right boundary\n", " T[-1] = (coeff*Tn[-2] + 2*coeff*TB) / (coeff + 2*coeff)\n", " \n", " # calculate the temperature difference for convergence\n", " Tdiff = np.sum(np.abs(T - Tn)) # total temperature difference\n", " if (n+1) % 10 == 0: # print every 100 iterations\n", " print('Iteration {}: Tdiff = {:.3f}'.format(n+1, Tdiff))" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkYAAAG2CAYAAACap0noAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjAsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvlHJYcgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAZMhJREFUeJzt3XdYFNf7NvB76U2qyIKioGIFFbFENHZF87WmaCSxayxBQUQUC00DigVjj71ExRhL/MUS0agRNYq9YBfEAlakd+b9A5g3K5iwuMtS7s917RX3zNnZZ5zA3J45MyMRBEEAEREREUFN1QUQERERlRcMRkREREQFGIyIiIiICjAYERERERVgMCIiIiIqwGBEREREVIDBiIiIiKgAgxERERFRAQYjIiIiogIMRkREREQFVB6Mnj17hm+//RZmZmbQ09NDixYtcOnSJXG5IAjw9/eHlZUVdHV10blzZ9y6dUtmHZmZmZg0aRKqV68OfX199OvXD0+fPi3rTSEiIqIKTqXBKCEhAe3bt4empiYOHz6MqKgoLF68GMbGxmKfkJAQLFmyBCtWrEBkZCSkUil69OiB5ORksY+Hhwf27duHsLAwREREICUlBX369EFubq4KtoqIiIgqKokqHyI7Y8YMnDlzBqdPny52uSAIsLKygoeHB6ZPnw4gf3TIwsICCxYswLhx45CYmAhzc3Ns27YNgwcPBgA8f/4c1tbWOHToEFxcXMpse4iIiKhi01Dllx84cAAuLi746quvcOrUKdSsWRMTJ07E2LFjAQDR0dGIj49Hz549xc9oa2ujU6dOOHv2LMaNG4dLly4hOztbpo+VlRXs7e1x9uzZYoNRZmYmMjMzxfd5eXl4+/YtzMzMIJFIlLjFREREpCiCICA5ORlWVlZQU1PMSTCVBqNHjx5h9erV8PT0xMyZM3HhwgVMnjwZ2traGDZsGOLj4wEAFhYWMp+zsLDA48ePAQDx8fHQ0tKCiYlJkT6Fn39fcHAwAgIClLBFREREVNaePHmCWrVqKWRdKg1GeXl5aNWqFYKCggAAjo6OuHXrFlavXo1hw4aJ/d4fxREE4T9Hdv6tj4+PDzw9PcX3iYmJqF27Np48eQJDQ8PSbg4RERGVoaSkJFhbW6NatWoKW6dKg5GlpSWaNGki09a4cWPs2bMHACCVSgHkjwpZWlqKfV6+fCmOIkmlUmRlZSEhIUFm1Ojly5dwdnYu9nu1tbWhra1dpN3Q0JDBiIiIqIJR5DQYlV6V1r59e9y9e1em7d69e6hTpw4AwNbWFlKpFOHh4eLyrKwsnDp1Sgw9Tk5O0NTUlOkTFxeHmzdvfjAYERERERVHpSNGU6ZMgbOzM4KCgjBo0CBcuHABa9euxdq1awHkJ0APDw8EBQXBzs4OdnZ2CAoKgp6eHlxdXQEARkZGGD16NKZOnQozMzOYmprCy8sLDg4O6N69uyo3j4iIiCoYlQaj1q1bY9++ffDx8UFgYCBsbW2xdOlSfPPNN2Ifb29vpKenY+LEiUhISEDbtm1x9OhRmfOJoaGh0NDQwKBBg5Ceno5u3bph8+bNUFdXV8VmERERUQWl0vsYlRdJSUkwMjJCYmLiv84xys3NRXZ2dhlWRlT+aWpq8h8hRKQSJT1+y0OlI0YVhSAIiI+Px7t371RdClG5ZGxsDKlUyvuAEVGFx2BUAoWhqEaNGtDT0+Mvf6ICgiAgLS0NL1++BACZq0eJiCoiBqP/kJubK4YiMzMzVZdDVO7o6uoCyL9FRo0aNXhajYgqNJVerl8RFM4p0tPTU3ElROVX4c8H5+ARUUXHYFRCPH1G9GH8+SCiyoLBiIiIiKgAgxGVOzY2Nli6dKnC1te5c2d4eHgobH3vU1S9yq6TiIj+G4NRJTVixAhIJBLMnz9fpn3//v3l/rRHZGQkvvvuO1WXoTQnT56ERCIpcvuHvXv3Yu7cuUr97jdv3qBXr16wsrKCtrY2rK2t4ebmhqSkJJl+N27cQKdOnaCrq4uaNWsiMDAQvOUZEVUFDEaVmI6ODhYsWICEhARVl1IiWVlZAABzc/MqOdnd1NRUoU+ILo6amhr69++PAwcO4N69e9i8eTOOHTuG8ePHi32SkpLQo0cPWFlZITIyEsuXL8eiRYuwZMkSpdZGRFQeMBhVYt27d4dUKkVwcPAH+/j7+6NFixYybUuXLoWNjY34fsSIERgwYACCgoJgYWEBY2NjBAQEICcnB9OmTYOpqSlq1aqFjRs3yqzn2bNnGDx4MExMTGBmZob+/fsjJiamyHqDg4NhZWWFBg0aACh6aurdu3f47rvvYGFhAR0dHdjb2+P3338HkD8CMmTIENSqVQt6enpwcHDAzp075fp7unbtGrp06YJq1arB0NAQTk5OuHjxorh8z549aNq0KbS1tWFjY4PFixd/cF0xMTGQSCS4evWqTP0SiQQnT55ETEwMunTpAgAwMTGBRCLBiBEjABQ9lZaQkIBhw4bBxMQEenp66N27N+7fvy8u37x5M4yNjfHHH3+gcePGMDAwQK9evRAXF/fB+kxMTDBhwgS0atUKderUQbdu3TBx4kScPn1a7LN9+3ZkZGRg8+bNsLe3x+eff46ZM2diyZIlHDUiokqP9zGSU+EN7VRB3ptLqqurIygoCK6urpg8eTJq1apV6u/+888/UatWLfz11184c+YMRo8ejXPnzqFjx444f/48du3ahfHjx6NHjx6wtrZGWloaunTpgk8//RR//fUXNDQ0MG/ePPTq1QvXr1+HlpYWAOD48eMwNDREeHh4sQfdvLw89O7dG8nJyfj5559Rr149REVFiffKycjIgJOTE6ZPnw5DQ0McPHgQQ4cORd26ddG2bdsSbds333wDR0dHrF69Gurq6rh69So0NTUBAJcuXcKgQYPg7++PwYMH4+zZs5g4cSLMzMzEQCMPa2tr7NmzB1988QXu3r0LQ0ND8T5A7xsxYgTu37+PAwcOwNDQENOnT8dnn32GqKgosb60tDQsWrQI27Ztg5qaGr799lt4eXlh+/btJarn+fPn2Lt3Lzp16iS2nTt3Dp06dYK2trbY5uLiAh8fH8TExMDW1lbu7SYiqigYjOSUlpYGAwMDlXx3SkoK9PX15frMwIED0aJFC/j5+WHDhg2l/m5TU1MsW7YMampqaNiwIUJCQpCWloaZM2cCAHx8fDB//nycOXMGX3/9NcLCwqCmpob169eLYW7Tpk0wNjbGyZMn0bNnTwCAvr4+1q9fLwal9x07dgwXLlzA7du3xRGlunXristr1qwJLy8v8f2kSZNw5MgR7N69u8TBKDY2FtOmTUOjRo0AAHZ2duKyJUuWoFu3bpgzZw4AoEGDBoiKisLChQtLFYzU1dVhamoKAKhRowaMjY2L7VcYiM6cOQNnZ2cA+SM51tbW2L9/P7766isA+fcNWrNmDerVqwcAcHNzQ2Bg4H/WMWTIEPz2229IT09H3759sX79enFZfHy8zIghAFhYWIjLGIyIqDLjqbQqYMGCBdiyZQuioqJKvY6mTZtCTe3//+9iYWEBBwcH8b26ujrMzMzER0NcunQJDx48QLVq1WBgYAADAwOYmpoiIyMDDx8+FD/n4ODwwVAEAFevXkWtWrXEUPS+3Nxc/PDDD2jWrBnMzMxgYGCAo0ePIjY2tsTb5unpiTFjxqB79+6YP3++TH23b99G+/btZfq3b98e9+/fR25ubom/Q163b9+GhoaGTLgzMzNDw4YNcfv2bbFNT09PDEVA/iM5CvfBvwkNDcXly5exf/9+PHz4EJ6enjLL3x+ZLBzNK+8T94mIPhZHjOSkp6eHlJQUlX13aXTs2BEuLi6YOXNmkVEONTW1Iqewirt7ceGpm0ISiaTYtry8PAD5p8CcnJyKPaVjbm4u/vm/RsA+dJqp0OLFixEaGoqlS5fCwcEB+vr68PDwECdyl4S/vz9cXV1x8OBBHD58GH5+fggLC8PAgQMhCMIHQ0JxCsPjP/uU5m7QH/qO9+spbh+UZB6QVCqFVCpFo0aNYGZmhk8//RRz5syBpaUlpFIp4uPjZfoXhq3CkSMiosqKwUhOEolE7tNZ5cH8+fPRokWLIiMv5ubmiI+Plzng/nPicGm1bNkSu3btQo0aNWBoaFjq9TRr1gxPnz7FvXv3ih01On36NPr3749vv/0WQH4gu3//Pho3bizX9zRo0AANGjTAlClTMGTIEGzatAkDBw5EkyZNEBERIdP37NmzaNCgQbHPBCsMfXFxcXB0dARQ9O+zcITs30acmjRpgpycHJw/f148lfbmzRvcu3dP7m37L4VBKjMzEwDQrl07zJw5E1lZWWKtR48ehZWVVZFTbERElQ1PpVURDg4O+Oabb7B8+XKZ9s6dO+PVq1cICQnBw4cPsXLlShw+fPijv++bb75B9erV0b9/f5w+fRrR0dE4deoU3N3d8fTp0xKvp1OnTujYsSO++OILhIeHIzo6GocPH8aRI0cAAPXr10d4eDjOnj2L27dvY9y4cUVGO/5Neno63NzccPLkSTx+/BhnzpxBZGSkGD6mTp2K48ePY+7cubh37x62bNmCFStWyMxr+iddXV188sknmD9/PqKiovDXX39h9uzZMn3q1KkDiUSC33//Ha9evSp2BNLOzg79+/fH2LFjERERgWvXruHbb79FzZo10b9//xJv3/sOHTqETZs24ebNm4iJicGhQ4cwYcIEtG/fXgw9rq6u0NbWxogRI3Dz5k3s27cPQUFB8PT05Kk0Iqr0GIyqkLlz5xY5zdK4cWOsWrUKK1euRPPmzXHhwoUPHvTloaenh7/++gu1a9fG559/jsaNG2PUqFFIT0+XewRpz549aN26NYYMGYImTZrA29tbHG2ZM2cOWrZsCRcXF3Tu3BlSqRQDBgwo8brV1dXx5s0bDBs2DA0aNMCgQYPQu3dvBAQEAMgf+frll18QFhYGe3t7+Pr6IjAw8F8nXm/cuBHZ2dlo1aoV3N3dMW/ePJnlNWvWREBAAGbMmAELCwu4ubkVu55NmzbByckJffr0Qbt27SAIAg4dOlTk9Jk8dHV1sW7dOnTo0AGNGzeGh4cH+vTpI97+AACMjIwQHh6Op0+folWrVpg4cSI8PT2LzEMiIqqMJAJvTIKkpCQYGRkhMTGxyEE7IyMD0dHRsLW1hY6OjooqJCrf+HNCRKrwb8fv0uKIEREREVEBBiMiIiKiAgxGRERERAUYjIiIiIgKMBgRERERFWAwIiIiIirAYERERERUgMGIiIiIqACDEREREVEBBiP6aDY2Nli6dOlHrePkyZOQSCR49+6dQmqKiYmBRCJRyANxi6OoepVdJxERyYfBqJI7e/Ys1NXV0atXL1WXIurcuTM8PDxk2pydnREXFwcjIyPVFFUGRowYUeQ5btbW1oiLi4O9vb1qiiIiIhkMRmUkNTUVEokEEokEqampZfa9GzduxKRJkxAREYHY2Ngy+155aWlpQSqVVrmnt6urq0MqlUJDQ0PVpRARERiMKrXU1FT88ssvmDBhAvr06YPNmzfLLC88HXT8+HG0atUKenp6cHZ2xt27d8U+Dx8+RP/+/WFhYQEDAwO0bt0ax44d++B3jho1Cn369JFpy8nJgVQqxcaNGzFixAicOnUKP/74oxgUY2Jiij01debMGXTq1Al6enowMTGBi4sLEhISAABHjhxBhw4dYGxsDDMzM/Tp0wcPHz6U6+9n1apVsLOzg46ODiwsLPDll1+KyzIzMzF58mTUqFEDOjo66NChAyIjIz+4Ln9/f7Ro0UKmbenSpbCxsRGXb9myBb/99pu43SdPniz2VNqpU6fQpk0baGtrw9LSEjNmzEBOTo64vHPnzpg8eTK8vb1hamoKqVQKf39/ubadiIiKx2BUie3atQsNGzZEw4YN8e2332LTpk0QBKFIv1mzZmHx4sW4ePEiNDQ0MGrUKHFZSkoKPvvsMxw7dgxXrlyBi4sL+vbt+8HRpzFjxuDIkSOIi4sT2w4dOoSUlBQMGjQIP/74I9q1a4exY8ciLi4OcXFxsLa2LrKeq1evolu3bmjatCnOnTuHiIgI9O3bF7m5uQDyQ5+npyciIyNx/PhxqKmpYeDAgcjLyyvR383FixcxefJkBAYG4u7duzhy5Ag6duwoLvf29saePXuwZcsWXL58GfXr14eLiwvevn1bovW/z8vLC4MGDUKvXr3E7XZ2di7S79mzZ/jss8/QunVrXLt2DatXr8aGDRswb948mX5btmyBvr4+zp8/j5CQEAQGBiI8PLxUtRER0T8IJCQmJgoAhMTExCLL0tPThaioKCE9PV3u9aakpIivFy9eCAAEAMKLFy9klimLs7OzsHTpUkEQBCE7O1uoXr26EB4eLi4/ceKEAEA4duyY2Hbw4EEBwL9ub5MmTYTly5eL7+vUqSOEhobKLF+wYIH4fsCAAcKIESPE9506dRLc3d1l1llYS0JCgiAIgjBkyBChffv2Jd7Wly9fCgCEGzduCIIgCNHR0QIA4cqVK8X237Nnj2BoaCgkJSUVWZaSkiJoamoK27dvF9uysrIEKysrISQkpNh6/fz8hObNm8usJzQ0VKhTp474fvjw4UL//v1l+rxf58yZM4WGDRsKeXl5Yp+VK1cKBgYGQm5uriAI+X9/HTp0kFlP69athenTpxe7rWXhY35OiIhK69+O36XFESMlMjAwEF8WFhZie+FpqcKXMty9excXLlzA119/DQDQ0NDA4MGDsXHjxiJ9mzVrJv7Z0tISAPDy5UsA+SMz3t7eaNKkCYyNjWFgYIA7d+7863ylMWPGYNOmTeJ6Dh48KDMKVRKFI0Yf8vDhQ7i6uqJu3bowNDSEra0tAJR4HlWPHj1Qp04d1K1bF0OHDsX27duRlpYmrjs7Oxvt27cX+2tqaqJNmza4ffu2XNshr9u3b6Ndu3Yyc63at2+PlJQUPH36VGz75z4D8vdb4T4jIqLS44zPSmrDhg3IyclBzZo1xTZBEKCpqYmEhASYmJiI7ZqamuKfCw/Ihaekpk2bhj/++AOLFi1C/fr1oauriy+//BJZWVkf/O5hw4ZhxowZOHfuHM6dOwcbGxt8+umnctWvq6v7r8v79u0La2trrFu3DlZWVsjLy4O9vf2/1vVP1apVw+XLl3Hy5EkcPXoUvr6+8Pf3R2RkpHi68f2J4IIgfHByuJqaWpHTlNnZ2SWq5b++o7h6/rnPCpeV9DQiERF9GEeMlCglJUV8vXjxQmx/8eKFzDJFy8nJwdatW7F48WJcvXpVfF27dg116tTB9u3bS7yu06dPY8SIERg4cCAcHBwglUoRExPzr58xMzPDgAEDsGnTJmzatAkjR46UWa6lpSXOFfqQZs2a4fjx48Uue/PmDW7fvo3Zs2ejW7duaNy4sTgpWx4aGhro3r07QkJCcP36dcTExODPP/9E/fr1oaWlhYiICLFvdnY2Ll68iMaNGxe7LnNzc8THx8uEo/fvTVSS7W7SpAnOnj0rs56zZ8+iWrVqMiGXiIiUgyNGSqSvr//B9g8tU4Tff/8dCQkJGD16dJH7An355ZfYsGED3NzcSrSu+vXrY+/evejbty8kEgnmzJlTopGJMWPGoE+fPsjNzcXw4cNlltnY2OD8+fOIiYmBgYEBTE1Ni3zex8cHDg4OmDhxIsaPHw8tLS2cOHECX331FUxNTWFmZoa1a9fC0tISsbGxmDFjRom2p9Dvv/+OR48eoWPHjjAxMcGhQ4eQl5eHhg0bQl9fHxMmTMC0adNgamqK2rVrIyQkBGlpaRg9enSx6+vcuTNevXqFkJAQfPnllzhy5AgOHz4MQ0NDme3+448/cPfuXZiZmRV7z6aJEydi6dKlmDRpEtzc3HD37l34+fnB09MTamr8dwwRkbLxN20ltGHDBnTv3r3YA+8XX3yBq1ev4vLlyyVaV2hoKExMTODs7Iy+ffvCxcUFLVu2/M/Pde/eHZaWlnBxcYGVlZXMMi8vL6irq6NJkyYwNzcvdl5QgwYNcPToUVy7dg1t2rRBu3bt8Ntvv0FDQwNqamoICwvDpUuXYG9vjylTpmDhwoUl2p5CxsbG2Lt3L7p27YrGjRtjzZo12LlzJ5o2bQoAmD9/Pr744gsMHToULVu2xIMHD/DHH3/InIL8p8aNG2PVqlVYuXIlmjdvjgsXLsDLy0umz9ixY9GwYUO0atUK5ubmOHPmTJH11KxZE4cOHcKFCxfQvHlzjB8/HqNHj8bs2bPl2j4iIiodifD+xIgqKCkpCUZGRkhMTJT5Fz4AZGRkIDo6Gra2ttDR0Sn1d6SmpooTrVNSUpQ6YlQepKWlwcrKChs3bsTnn3+u6nJIyRT1c0JEJI9/O36XFk+llRF9ff1i7yFU2eTl5SE+Ph6LFy+GkZER+vXrp+qSiIiISozBiBQqNjYWtra2qFWrFjZv3sxHXRARUYXCoxYplI2NTZUYGSMiosqJk6+JiIiICjAYlRBHQYg+jD8fRFRZMBj9h8I7DBc+LoKIiir8+Xj/jtxERBUN5xj9B3V1dRgbG4vPodLT0/vgYyGIqhpBEJCWloaXL1/C2NgY6urqqi6JiOijMBiVgFQqBQA+pJPoA4yNjcWfEyKiiozBqAQkEgksLS1Ro0aNUj0YlKgy09TU5EgREVUaDEZyUFdX5wGAiIioEuPkayIiIqICDEZEREREBVQajPz9/SGRSGRe/5zAKQgC/P39YWVlBV1dXXTu3Bm3bt2SWUdmZiYmTZqE6tWrQ19fH/369cPTp0/LelOIiIioElD5iFHTpk0RFxcnvm7cuCEuCwkJwZIlS7BixQpERkZCKpWiR48eSE5OFvt4eHhg3759CAsLQ0REBFJSUtCnTx/k5uaqYnOIiIioAlP55GsNDY1iL/MVBAFLly7FrFmz8PnnnwMAtmzZAgsLC+zYsQPjxo1DYmIiNmzYgG3btqF79+4AgJ9//hnW1tY4duwYXFxcynRbiIiIqGJT+YjR/fv3YWVlBVtbW3z99dd49OgRACA6Ohrx8fHo2bOn2FdbWxudOnXC2bNnAQCXLl1Cdna2TB8rKyvY29uLfYqTmZmJpKQkmRcRERGRSoNR27ZtsXXrVvzxxx9Yt24d4uPj4ezsjDdv3iA+Ph4AYGFhIfMZCwsLcVl8fDy0tLRgYmLywT7FCQ4OhpGRkfiytrZW8JYRERFRRaTSYNS7d2988cUXcHBwQPfu3XHw4EEA+afMCr3/+A1BEP7zkRz/1cfHxweJiYni68mTJx+xFURERFRZqPxU2j/p6+vDwcEB9+/fF+cdvT/y8/LlS3EUSSqVIisrCwkJCR/sUxxtbW0YGhrKvIiIiIjKVTDKzMzE7du3YWlpCVtbW0ilUoSHh4vLs7KycOrUKTg7OwMAnJycoKmpKdMnLi4ON2/eFPsQERERlZRKr0rz8vJC3759Ubt2bbx8+RLz5s1DUlIShg8fDolEAg8PDwQFBcHOzg52dnYICgqCnp4eXF1dAQBGRkYYPXo0pk6dCjMzM5iamsLLy0s8NUdEREQkD5UGo6dPn2LIkCF4/fo1zM3N8cknn+Dvv/9GnTp1AADe3t5IT0/HxIkTkZCQgLZt2+Lo0aOoVq2auI7Q0FBoaGhg0KBBSE9PR7du3bB582Y+04yIiIjkJhEEQVB1EaqWlJQEIyMjJCYmcr4RERFRBaGM43e5mmNEREREpEoMRkREREQFGIyIiIiICjAYERERERVgMCIiIiIqwGBEREREVIDBiIiIiKgAgxERERFRAQYjIiIiogIMRkREREQFGIyIiIiICjAYERERERVgMCIiIqpCUlNTIZFIIJFIkJqaqupyyh0GIyIiIqICDEZEREREBTTk/UBmZiYuXLiAmJgYpKWlwdzcHI6OjrC1tVVGfURERPSR/nnK7EN/BgB9ff0yq6m8KnEwOnv2LJYvX479+/cjKysLxsbG0NXVxdu3b5GZmYm6deviu+++w/jx41GtWjVl1kxERERyMDAwKLbdwsJC5r0gCGVRTrlWolNp/fv3x5dffomaNWvijz/+QHJyMt68eYOnT58iLS0N9+/fx+zZs3H8+HE0aNAA4eHhyq6biIiISOFKNGLUs2dP7N69G1paWsUur1u3LurWrYvhw4fj1q1beP78uUKLJCIiotIRBAG//fYbgoOD8ffff8ssu3btGurVq6eiysonicBxMyQlJcHIyAiJiYkwNDRUdTlEREQfTRAEHD16FAEBATh37hwAQFtbG6NGjcLq1asBACkpKRV6XpEyjt9yX5X25MkTPH36VHx/4cIFeHh4YO3atQopiIiIiEpPEAQcOXIEzs7O6NWrF86dOwcdHR24u7vj0aNHWLhwoapLLNfkDkaurq44ceIEACA+Ph49evTAhQsXMHPmTAQGBiq8QCIiIvpvhYGoXbt26N27N/7++2/o6OjAw8MDjx49wtKlS2FlZaXqMss9uYPRzZs30aZNGwDAL7/8Ant7e5w9exY7duzA5s2bFV0fERER/QtBEHD48GF88skn6N27N86fPw9dXV1MmTIF0dHRCA0NhaWlpdhfX18fgiBAEIQKfRpNWeS+j1F2dja0tbUBAMeOHUO/fv0AAI0aNUJcXJxiqyMiIqJiFQYif39/REZGAgB0dXUxYcIETJs2DVKpVMUVVkxyjxg1bdoUa9aswenTpxEeHo5evXoBAJ4/fw4zMzOFF0hERET/nyAIOHjwINq0aYP//e9/iIyMhK6uLqZOnYro6GgsXryYoegjyD1itGDBAgwcOBALFy7E8OHD0bx5cwDAgQMHxFNsREREpFiFgSggIAAXL14EAOjp6WHixInw8vIqcrNGKp1SXa6fm5uLpKQkmJiYiG0xMTHQ09NDjRo1FFpgWeDl+kREVF4JgoDff/8dAQEBuHTpEoD8QPT999/Dy8urQh53FUUZx2+5R4wAQF1dXSYUAYCNjY0i6iEiIiLkB6L/+7//Q0BAAC5fvgwgPxC5ubnBy8sL5ubmKq6wcpI7GNna2kIikXxw+aNHjz6qICIioqpMEAQcOHAAAQEBuHLlCoD8K8nc3NwwdepUBiIlkzsYeXh4yLzPzs7GlStXcOTIEUybNk1RdREREVUphY/uCAwMlAlEkyZNwtSpU1G9enUVV1g1yB2M3N3di21fuXKlOBmMiIiISiYvL08MRFevXgUAGBgYYNKkSfD09GQgKmMKe1bao0eP0KJFCyQlJSlidWWKk6+JiKis5eXlYf/+/QgMDMS1a9cA5AeiyZMnw9PTk7fAKYFyM/m6OL/++itMTU0VtToiIqJKKS8vD/v27UNgYCCuX78OAKhWrRomT56MKVOmMBCpmNzByNHRUWbytSAIiI+Px6tXr7Bq1SqFFkdERFRZ5OXlYe/evQgMDMSNGzcA5Acid3d3TJkyhYML5YTcwWjAgAEy79XU1GBubo7OnTujUaNGiqqLiIioUigMRAEBAbh58yYAwNDQEO7u7vDw8GAgKmcUNseoIuMcIyIiUrS8vDzs2bMHgYGBMoHIw8MDHh4eRe4HSPJT2RyjpKQk8Qv/a3I1gwUREVVleXl5+PXXXxEYGIhbt24BYCCqSEoUjExMTBAXF4caNWrA2Ni42Bs8CoIAiUSC3NxchRdJRERU3uXm5oqBKCoqCgBgZGQEDw8PuLu7MxBVECUKRn/++ad4DvTEiRNKLYiIiKgiyc3Nxe7duxEYGIjbt28DyA9EU6ZMgbu7O4yNjVVbIMmFc4zAOUZERCS/3Nxc/PLLLwgMDMSdO3cAAMbGxpgyZQomT57MQFQGlHH8VpP3A0eOHEFERIT4fuXKlWjRogVcXV2RkJCgkKKIiIjKq9zcXOzYsQP29vZwdXXFnTt3YGxsjMDAQMTExMDX15ehqAKTOxhNmzZNnIB948YNeHp64rPPPsOjR4/g6emp8AKJiIjKg9zcXGzfvh1NmzbFN998gzt37sDExARz585FTEwM5syZAyMjI1WXSR9J7vsYRUdHo0mTJgCAPXv2oG/fvggKCsLly5fx2WefKbxAIiIiVcrJyUFYWBjmzZuHu3fvAgBMTU3h6emJSZMmcQpGJSN3MNLS0kJaWhoA4NixYxg2bBiA/P9JKuJz0oiIiIqTk5ODnTt3Yt68ebh37x6A/GPd1KlT4ebmxkBUSckdjDp06ABPT0+0b98eFy5cwK5duwAA9+7dQ61atRReIBERUVnKycnBjh07MG/ePNy/fx9AfiDy8vKCm5sbqlWrpuIKSZnknmO0YsUKaGho4Ndff8Xq1atRs2ZNAMDhw4fRq1cvhRdIRERUFnJycrBlyxY0btwYw4cPx/3792FmZobg4GDExMTAx8eHoagK4OX64OX6RERVWU5ODrZv3465c+fi4cOHAAAzMzNMmzYNEydOZBgqx1T6SJCSYrAgIqKKICcnBz///DPmzZsnBqLq1auLgcjAwEDFFZIqlCgYfegxIP/ER4IQEVFFkJ2dLQaiR48eAcgPRN7e3pgwYQIDURVXomDEx4AQEVFFl52djW3btmHevHmIjo4GAJibm4uBSF9fX8UVUrkglBNBQUECAMHd3V1sy8vLE/z8/ARLS0tBR0dH6NSpk3Dz5k2Zz2VkZAhubm6CmZmZoKenJ/Tt21d48uSJXN+dmJgoABASExMVsSlERFSOZGVlCevXrxdsbW0FAAIAoUaNGsKiRYuElJQUVZdHH0EZx2+5r0oDgNOnT+Pbb7+Fs7Mznj17BgDYtm2bzKNC5BEZGYm1a9eiWbNmMu0hISFYsmQJVqxYgcjISEilUvTo0QPJycliHw8PD+zbtw9hYWGIiIhASkoK+vTpw1N6RERVXFZWFtavX48GDRpgzJgxiI6ORo0aNbB48WJER0dj6tSpHCWiIuQORnv27IGLiwt0dXVx+fJlZGZmAgCSk5MRFBQkdwEpKSn45ptvsG7dOpiYmIjtgiBg6dKlmDVrFj7//HPY29tjy5YtSEtLw44dOwAAiYmJ2LBhAxYvXozu3bvD0dERP//8M27cuIFjx47JXQsREVV8WVlZWLduHRo0aICxY8ciJiYGFhYWWLJkCaKjo+Hp6Qk9PT1Vl0nllNzBaN68eVizZg3WrVsHTU1Nsd3Z2RmXL1+Wu4Dvv/8e//vf/9C9e3eZ9ujoaMTHx6Nnz55im7a2Njp16oSzZ88CAC5duoTs7GyZPlZWVrC3txf7FCczMxNJSUkyLyIiqtiysrKwdu1a2NnZ4bvvvsPjx49hYWGB0NBQPHr0CFOmTGEgov8k952v7969i44dOxZpNzQ0xLt37+RaV1hYGC5fvozIyMgiy+Lj4wEAFhYWMu0WFhZ4/Pix2EdLS0tmpKmwT+HnixMcHIyAgAC5aiUiovIpKysLmzZtQlBQEGJjYwEAUqkU06dPx3fffccwRHKRe8TI0tISDx48KNIeERGBunXrlng9T548gbu7O37++Wfo6Oh8sN/7twkQCm4L8G/+q4+Pjw8SExPF15MnT0pcNxERlQ+ZmZlYs2YN7OzsMH78eMTGxsLS0hJLly7Fo0eP4OHhwVBEcpN7xGjcuHFwd3fHxo0bIZFI8Pz5c5w7dw5eXl7w9fUt8XouXbqEly9fwsnJSWzLzc3FX3/9hRUrVohPMI6Pj4elpaXY5+XLl+IoklQqRVZWFhISEmRGjV6+fAlnZ+cPfre2tja0tbVLXCsREZUfmZmZ2LhxI4KDg8V/2FpaWmLGjBkYO3YsdHV1VVwhVWRyByNvb28kJiaiS5cuyMjIQMeOHaGtrS0+XK+kunXrhhs3bsi0jRw5Eo0aNcL06dNRt25dSKVShIeHw9HREUD+cOmpU6ewYMECAICTkxM0NTURHh6OQYMGAQDi4uJw8+ZNhISEyLtpRERUjmVmZmLDhg0IDg7G06dPAeTPKy0MRP929oGopOQORgDwww8/YNasWYiKikJeXh6aNGki951Cq1WrBnt7e5k2fX19mJmZie0eHh4ICgqCnZ0d7OzsEBQUBD09Pbi6ugIAjIyMMHr0aEydOhVmZmbi048dHByKTOYmIqKKKSMjAxs2bMD8+fPFQFSzZk3MmDEDY8aMYSAihSpxMMrNzcWtW7dgZ2cHXV1d6OnpoVWrVgCA9PR0XL9+Hfb29lBTK9WtkYrl7e2N9PR0TJw4EQkJCWjbti2OHj0q80C/0NBQaGhoYNCgQUhPT0e3bt2wefNmqKurK6wOIiIqe4WBKDg4WLxnXs2aNeHj44PRo0czEJFSSARBEErScfPmzVixYgXOnz9fJHTk5uaibdu28PDwwLfffquUQpVJGU/nJSKi0snIyMD69esRHByM58+fAwBq1aolBiLOEaVCyjh+l3h4Z8OGDfDy8ip2JEZdXR3e3t5Yu3atQooiIqKqJyMjA8uXL0e9evUwadIkPH/+HLVq1cKqVavw4MEDTJw4kaGIlK7Ep9Lu3r2LTz755IPLW7dujdu3byukKCIiqlxSU1PFuagpKSkyj+JIT0/HunXrMH/+fMTFxQEArK2tMXPmTIwcOZJhiMpUiYNRamrqv94hOjk5GWlpaQopioiIKr/09HSsXbsW8+fPF2/Ka21tjVmzZmHEiBEMRKQSJQ5GdnZ2OHv2bJEHvRaKiIiAnZ2dwgojIqLKqXCEaMGCBWIgql27thiItLS0VFwhVWUlDkaurq6YPXs2nJ2di4Sja9euwdfXF97e3govkIiIKqbU1NRi/9ykSRO8evUKQH4gmj17NoYPH85AROVCia9KK3xYa0REBLp3745GjRpBIpHg9u3bOHbsGNq3b4/w8HCZB8tWFLwqjYhI8f7r8U2FSngYIipCGcfvEgcjID8chYaGYseOHbh//z4EQUCDBg3g6uoKDw+PCpv2GYyIiBSPwYiUTeXBqLJiMCIiUpzU1FSsXr0aISEhMqfMYmNjAQAvXryQuSrtn38mkodK72NERET0b1JTU7Fw4ULY2tpi2rRpePXqFerWrYuNGzfi2rVrYj99fX2ZF1F5UqpnpRERERVKSUnBqlWrsGjRInGEqG7dupg9eza+/fZbaGpqyky+JirPGIyIiKhUUlJSsHLlSixatAivX78GANSrVw+zZ8/GN998UyEvxiFiMCIiIrkkJyeLgejNmzcAgPr164uBSEOj6KFFX1+fk6ypQih1MMrKykJ0dDTq1atX7A8BERFVLsnJyVixYgUWL14sBiI7OzvMnj0brq6uPBZQpSD35Ou0tDSMHj0aenp6aNq0qXiVweTJkzF//nyFF0hERKqVnJyM4OBg2NjYYObMmXjz5g3s7OywdetWREVFYdiwYQxFVGnIHYx8fHxw7do1nDx5Ejo6OmJ79+7dsWvXLoUWR0REqpOUlISgoCAxEL19+xYNGjTAtm3bEBUVhaFDhzIQUaUj9//R+/fvx65du/DJJ5/I3LyrSZMmePjwoUKLIyKispeUlITly5dj8eLFSEhIAAA0aNAAvr6++Prrr6Gurq7iComUR+5g9OrVK9SoUaNIe2pqaonvckpEROVPYmIili9fjiVLloiBqGHDhvD19cXgwYMZiKhKkPtUWuvWrXHw4EHxfWEYWrduHdq1a6e4yoiIqEwkJiZi7ty5sLGxwZw5c5CQkIBGjRphx44duHXrFlxdXRmKqMqQe8QoODgYvXr1QlRUFHJycvDjjz/i1q1bOHfuHE6dOqWMGomISAnevXuHZcuWITQ0FO/evQMANG7cGL6+vvjqq68YhqhKknvEyNnZGWfOnEFaWhrq1auHo0ePwsLCAufOnYOTk5MyaiQiIgV69+4dAgICYGNjAz8/P7x79w5NmjRBWFgYbty4wXlEVKXxIbLgQ2SJqGp49+4dli5diqVLlyIxMRFA/oUzfn5++PLLL6GmxsdnUsWijON3qa6zfPjwITZt2oRHjx5h6dKlqFGjBo4cOQJra2s0bdpUIYUREZFiJCQkiIEoKSkJANC0aVP4+fnhiy++YCAi+ge5fxpOnToFBwcHnD9/Hnv27EFKSgoA4Pr16/Dz81N4gUREVDoJCQnw9fWFjY0NAgMDkZSUBHt7e/zyyy+4fv06vvrqK4YiovfI/RMxY8YMzJs3D+Hh4dDS0hLbu3TpgnPnzim0OCIikt/bt28xZ84c2NjYYO7cuUhKSoKDgwN2796Na9euMRAR/Qu5T6XduHEDO3bsKNJubm4uPjuHiIjK3tu3b7FkyRIsW7YMycnJAIBmzZrB19cXAwcOZBgiKgG5g5GxsTHi4uJga2sr037lyhXUrFlTYYUREVHJvHnzBkuWLMHy5ctlApGfnx8GDBjAQEQkB7l/WlxdXTF9+nTEx8dDIpEgLy8PZ86cgZeXF4YNG6aMGomIqBivX7/GzJkzYWNjg6CgICQnJ6N58+bYu3cvrly5gs8//5yhiEhOcl+un52djREjRiAsLAyCIEBDQwO5ublwdXXF5s2bK+S9L3i5PhFVJK9fv8bixYuxYsUK8QKYFi1awM/PD/369WMYoipDGcdvuYKRIAiIjY2Fubk54uPjcfnyZeTl5cHR0RF2dnYKKUgVGIyIqCJ49eqVGIhSU1MB5Acif39/9OvXj8+rpCpH5fcxEgQBdnZ2uHXrFuzs7FC3bl2FFEFERB/26tUrLFq0CCtXrhQDkaOjI/z9/dG3b18GIiIFkisYqampwc7ODm/evKnQI0RERBVBYSBasWIF0tLSAAAtW7aEv78/+vTpw0BEpARyn4gOCQnBtGnTcPPmTWXUQ0RU5b18+RLe3t6wsbFBSEgI0tLS4OTkhP/7v//DxYsXOUpEpERyT742MTFBWloacnJyoKWlBV1dXZnlb9++VWiBZYFzjIioPHj58iUWLlyIVatWiSNErVq1gr+/Pz777DOGIaL3qHyOEQAsXbpUIV9MRET5Xrx4IQai9PR0AEDr1q3h7++P3r17MxARlSG5g9Hw4cOVUQcRUZUTHx+PhQsXYvXq1WIgatOmDfz8/BiIiFRE7mBU+GTm90kkEmhra8s8P42IiIqKj49HSEgIVq9ejYyMDABA27Zt4efnh169ejEQEalQqR4J8m8/tLVq1cKIESPg5+fHm4wREf1DXFwcQkJCsGbNGplA5O/vDxcXFwYionJA7mC0efNmzJo1CyNGjECbNm0gCAIiIyOxZcsWzJ49W7y8VFtbGzNnzlRGzUREFUpcXBwWLFiAn376SQxEn3zyCfz9/dGzZ08GIqJyRO5gtGXLFixevBiDBg0S2/r16wcHBwf89NNPOH78OGrXro0ffviBwYiIqrTnz5+LgSgzMxMA0K5dO/j7+6NHjx4MRETlkNznus6dOwdHR8ci7Y6Ojjh37hwAoEOHDoiNjf346oiIKqBnz55h8uTJqFu3LpYtW4bMzEw4Ozvj6NGjOHPmDEeJiMoxuYNRrVq1sGHDhiLtGzZsgLW1NQDgzZs3MDEx+fjqiIgqkGfPnmHSpEmoV68eli9fjszMTLRv3x7h4eGIiIjgKBFRBSD3qbRFixbhq6++wuHDh9G6dWtIJBJERkbizp07+PXXXwEAkZGRGDx4sMKLJSIqj54+fYr58+dj3bp1yMrKApA/cu7v74+uXbsyDBFVIHLf+RoAYmJisGbNGty7dw+CIKBRo0YYN24cbGxslFCi8vHO10RUGk+ePMH8+fOxfv16MRB9+umn8Pf3R5cuXRiIiJRMGcfvUgWjyobBiIjk8eTJEwQHB2PDhg1iIOrYsSP8/f3RuXNnBiKiMqKM43epbjR0+vRpfPvtt3B2dsazZ88AANu2bUNERIRCiiIiKo9iY2MxYcIE1KtXD6tXr0ZWVhY6deqEEydO4NSpUxwlIqoE5A5Ge/bsgYuLC3R1dXH58mXxEtTk5GQEBQUpvEAiIlV7/Pgxxo8fj/r162PNmjXIzs5G586dceLECZw8eRKdO3dWdYlEpCByB6N58+ZhzZo1WLduHTQ1NcV2Z2dnXL58WaHFERGp0uPHjzFu3DjY2dnhp59+QnZ2Nrp06YKTJ0/ixIkTDERElZDcV6XdvXsXHTt2LNJuaGiId+/eKaImIiKVevz4MYKCgrBp0yZkZ2cDALp06QI/Pz906tRJxdURkTLJPWJkaWmJBw8eFGmPiIhA3bp1FVIUEZEqxMTE4LvvvkP9+vWxdu1aZGdno2vXrjh16hT+/PNPhiKiKkDuYDRu3Di4u7vj/PnzkEgkeP78ObZv3w4vLy9MnDhRGTUSESlVdHQ0xo4dCzs7O6xbtw45OTno1q0b/vrrLxw/frzYUXIiqpzkPpXm7e2NxMREdOnSBRkZGejYsSO0tbXh5eUFNzc3ZdRIRFQqqampMDAwAACkpKRAX19fZvmjR48QFBSELVu2ICcnBwDQvXt3+Pn5oUOHDmVeLxGpXqku1//hhx/w+vVrXLhwAX///TdevXqFuXPnyr2e1atXo1mzZjA0NIShoSHatWuHw4cPi8sFQYC/vz+srKygq6uLzp0749atWzLryMzMxKRJk1C9enXo6+ujX79+ePr0aWk2i4iqiEePHmH06NFo0KABNmzYgJycHPTo0QMREREIDw9nKCKqwkoVjABAT08PrVq1Qps2bcR/kcmrVq1amD9/Pi5evIiLFy+ia9eu6N+/vxh+QkJCsGTJEqxYsQKRkZGQSqXo0aMHkpOTxXV4eHhg3759CAsLQ0REBFJSUtCnTx/k5uaWdtOIqJJ6+PAhRo0ahQYNGmDjxo3Izc1Fz549cebMGRw9ehTt27dXdYlEpGIluvP1559/XuIV7t2796MKMjU1xcKFCzFq1ChYWVnBw8MD06dPB5A/OmRhYYEFCxZg3LhxSExMhLm5ObZt2yY+m+358+ewtrbGoUOH4OLiUqLv5J2viSqP1NRUmT9bWFgAAL7++mvs3r1b/EeTi4sL/Pz80K5dO5XUSUQfT2V3vjYyMhJfhoaGOH78OC5evCguv3TpEo4fPw4jI6NSF5Kbm4uwsDCkpqaiXbt2iI6ORnx8PHr27Cn20dbWRqdOnXD27Fnxe7Ozs2X6WFlZwd7eXuxTnMzMTCQlJcm8iKhyMDAwEF+FoQgAwsLCZEaSjxw5wlBEREWUaPL1pk2bxD9Pnz4dgwYNwpo1a6Curg4gP9RMnDixVGntxo0baNeuHTIyMmBgYIB9+/ahSZMmYrD55y+2wvePHz8GAMTHx0NLSwsmJiZF+sTHx3/wO4ODgxEQECB3rURERFS5yT3HaOPGjfDy8hJDEQCoq6vD09MTGzdulLuAhg0b4urVq/j7778xYcIEDB8+HFFRUeLy9587JAjCfz6L6L/6+Pj4IDExUXw9efJE7rqJqPy5ffs2Bg0aBDW1or/aXrx4gZSUFPFFRFQcuYNRTk4Obt++XaT99u3byMvLk7sALS0t1K9fH61atUJwcDCaN2+OH3/8EVKpFACKjPy8fPlSHEWSSqXIyspCQkLCB/sUR1tbW7wSrvBFRBVXVFQUhgwZgqZNm+KXX35BXl4e+vbti7/++kvso6+vL/MiIiqO3MFo5MiRGDVqFBYtWoSIiAhERERg0aJFGDNmDEaOHPnRBQmCgMzMTNja2kIqlSI8PFxclpWVhVOnTsHZ2RkA4OTkBE1NTZk+cXFxuHnzptiHiCqvmzdvYtCgQbC3t0dYWBgEQUD//v1x6dIlHDhwAC1btlR1iURUwch9g8dFixZBKpUiNDQUcXFxAPIfE+Lt7Y2pU6fKta6ZM2eid+/esLa2RnJyMsLCwnDy5EkcOXIEEokEHh4eCAoKgp2dHezs7BAUFAQ9PT24uroCyJ8UPnr0aEydOhVmZmYwNTWFl5cXHBwc0L17d3k3jYgqiOvXryMwMBB79uwR2wYOHAhfX1+0aNFCdYURUYUndzBSU1ODt7c3vL29xau5Snsq6sWLFxg6dCji4uJgZGSEZs2a4ciRI+jRoweA/Ltsp6enY+LEiUhISEDbtm1x9OhRVKtWTVxHaGgoNDQ0MGjQIKSnp6Nbt27YvHmzzBwoIqocrl69isDAQOzbt09s+/LLLzFnzhw0a9asSH99fX2U4I4kRESiEt3HqLLjfYyIyrdLly4hMDAQBw4cAJB/UcZXX32FOXPmwN7eXsXVEZGqqOw+Rr169frX+wIVSk5OxoIFC7By5cqPLoyI6OLFi+jbty9atWqFAwcOQCKRYMiQIbh58yZ27drFUERECleiU2lfffUVBg0ahGrVqqFfv35o1aoVrKysoKOjg4SEBERFRSEiIgKHDh1Cnz59sHDhQmXXTUSV2IULFxAQEIBDhw4ByD+FP2TIEMyePRuNGjVScXVEVJmV+FRaVlYWfv31V+zatQunT5/Gu3fv8lcgkaBJkyZwcXHB2LFj0bBhQ2XWqxQ8lUZUPpw7dw4BAQH4448/AOQHom+++QazZ89GgwYNVFwdEZU3yjh+l3qOUWJiItLT02FmZgZNTU2FFKMqDEZEqnXmzBkEBASIt95QV1fH0KFDMXPmTNjZ2am4OiIqr5Rx/Jb7qrRChc9OIyIqrb/++gsBAQH4888/AQAaGhoYNmwYZs6ciXr16qm4OiKqikodjIiISuvkyZMICAjAyZMnAeQHopEjR8LHxwe2traqLY6IqjQGIyIqE4Ig4MSJEwgICBAf1aGpqYlRo0bBx8cHderUUXGFREQMRkSkZIIg4NixYwgMDERERASA/Gckjh49GjNmzEDt2rVVXCER0f/HYERESiEIAo4ePYqAgACcO3cOQP4DnMeOHYvp06ejVq1aKq6QiKioUgWjd+/e4ddff8XDhw8xbdo0mJqa4vLly7CwsEDNmjUVXSMRVSCCIODIkSMICAjA+fPnAQA6Ojr47rvv4O3tzd8RRFSuyR2Mrl+/ju7du8PIyAgxMTEYO3YsTE1NsW/fPjx+/Bhbt25VRp1EVM4JgoCDBw8iMDAQkZGRAABdXV2MHz8e06ZNg6WlpYorJCL6byV6JMg/eXp6YsSIEbh//z50dHTE9t69e4sTKomo6hAEAQcOHEDr1q3Rt29fREZGQldXF1OnTsWjR4+wZMkShiIiqjDkHjGKjIzETz/9VKS9Zs2aiI+PV0hRRFT+5eXl4bfffkNgYCCuXr0KANDT08P3338PLy8v1KhRQ7UFEhGVgtzBSEdHB0lJSUXa7969C3Nzc4UURUTlV15eHvbt24fAwEBcv34dAKCvrw83NzdMnTqVvweIqEKT+1Ra//79ERgYiOzsbAD5z0qLjY3FjBkz8MUXXyi8QCIqH/Ly8rB79260aNECX375Ja5fv45q1aph5syZiImJwfz58xmKiKjCk/tZaUlJSfjss89w69YtJCcnw8rKCvHx8WjXrh0OHToEfX19ZdWqNHxWGtGH5ebmYvfu3Zg7dy6ioqIAAIaGhpg8eTKmTJkCU1NTFVdIRFVVuXhWmqGhISIiIvDnn3/i8uXLyMvLQ8uWLdG9e3eFFERE5UNubi527dqFefPm4fbt2wDyn5Ho4eEBd3d3mJiYqLhCIiLFkysY5eTkQEdHB1evXkXXrl3RtWtXZdVFRCqSk5ODsLAwzJs3D3fv3gUAGBsbY8qUKZg8eTKMjY1VWyARkRLJFYw0NDRQp04d5ObmKqseIlKRnJwc7NixA/PmzcP9+/cBACYmJvD09MSkSZNgZGSk4gqJiJRP7snXs2fPho+PD96+fauMeoiojGVnZ2PTpk1o1KgRhg8fjvv378PMzAxBQUGIiYnB7NmzGYqIqMqQe47RsmXL8ODBA1hZWaFOnTpFJltfvnxZYcURkfJkZ2dj69at+OGHHxAdHQ0AqF69Ory8vDBx4kRUq1ZNxRUSEZU9uYPRgAEDlFAGEZWVrKwsbN68GUFBQXj8+DEAoEaNGpg2bRomTJhQIa8sJSJSFLkv16+MeLk+VQWZmZnYtGkTgoODERsbCwCwsLCAt7c3xo8fDz09PRVXSEQkn3JxuT4RVSwZGRnYuHEjgoOD8fTpUwCAVCrF9OnT8d133zEQERH9g9zBSE1NDRKJ5IPLecUaUfmQkZGBdevWYcGCBXj27BkAwMrKCjNmzMCYMWOgq6ur4gqJiMofuYPRvn37ZN5nZ2fjypUr2LJlCwICAhRWGBGVTnp6OtauXYsFCxYgLi4OAFCrVi3MmDEDo0ePho6OjoorJCIqvxQ2x2jHjh3YtWsXfvvtN0WsrkxxjhFVBmlpafjpp58QEhKC+Ph4AIC1tTVmzpyJkSNHQltbW8UVEhEpVrmeY9S2bVuMHTtWUasjohJKTU3F6tWrsXDhQrx8+RIAUKdOHcycORMjRoyAlpaWiiskIqo4FBKM0tPTsXz5ctSqVUsRqyOiEkhJScGqVauwaNEivHr1CgBgY2ODWbNmYdiwYQxERESlIHcwMjExkZl8LQgCkpOToaenh59//lmhxRFRUcnJyVixYgUWL16MN2/eAADq1auHWbNm4dtvv4WmpqaKKyQiqrjkDkahoaEywUhNTQ3m5uZo27Ytn7ZNpERJSUlYvnw5lixZIj6Sp379+pg9eza++eYbaGjw7htERB9L7t+kXbt2hbW1dbGX7MfGxqJ27doKKYyI8iUmJmLZsmUIDQ1FQkICAKBBgwaYM2cOvv76awYiIiIFkvs3qq2tLeLi4lCjRg2Z9jdv3sDW1pb3MSJSkHfv3mHp0qVYunQpEhMTAQCNGjXCnDlzMHjwYKirq6u4QiKiykfuYPShq/tTUlJ4fxQiBXj79i2WLl2KH3/8EUlJSQCAJk2aYM6cOfjqq68YiIiIlKjEwcjT0xMAIJFI4OvrK/MYgdzcXJw/fx4tWrRQeIFEVcWbN28QGhqKZcuWITk5GQBgb2+POXPm4Msvv4SampqKKyQiqvxKHIyuXLkCIH/E6MaNGzKXAmtpaaF58+bw8vJSfIVEldzr16+xePFirFixAikpKQCAZs2awdfXFwMHDmQgIiIqQyUORidOnAAAjBw5Ej/++CPvEE30kV69eoVFixZh5cqVSE1NBQC0aNECvr6+6N+/PwMREZEKyD3HaNOmTcqog6jSSE1NhYGBAYD8uXf6+voyy1+8eIGFCxdi9erVSEtLAwC0bNkSfn5+6Nu3778+pJmIiJSrVNf5RkZGYvfu3YiNjUVWVpbMsr179yqkMKLKJj4+HiEhIVizZg3S09MBAK1atYKfnx/+97//MRAREZUDco/Vh4WFoX379oiKisK+ffuQnZ2NqKgo/PnnnzAyMlJGjUQV2vPnz+Hh4QFbW1uEhoYiPT0dbdq0wcGDB3HhwgX06dOHoYiIqJyQe8QoKCgIoaGh+P7771GtWjX8+OOPsLW1xbhx42BpaamMGonKvcI5Qu//2d3dHT///DMyMzMBAO3atYOfnx969uzJMEREVA5JhA/dmOgD9PX1cevWLdjY2KB69eo4ceIEHBwccPv2bXTt2hVxcXHKqlVpkpKSYGRkhMTERE4qp1IpacjJy8tjICIiUhBlHL/lPpVmamoq3mOlZs2auHnzJoD8u/QWTiQlouIxFBERlW9yB6NPP/0U4eHhAIBBgwbB3d0dY8eOxZAhQ9CtWzeFF0hU3sXExGDkyJHFPtX+xYsXSElJEV9ERFS+yT3HaMWKFcjIyAAA+Pj4QFNTExEREfj8888xZ84chRdIVF49evQIQUFB2LJlC3JycgDkP2TZ29sbvXr1ApB/6vn9y/WJiKj8kmuOUU5ODrZv3w4XFxdIpVJl1lWmOMeI5PHgwQMEBQVh69at4kOTe/ToAV9fX3To0OE/72NERESKofI5RhoaGpgwYYJ4hQ1RVXL//n0MHz4cjRo1wqZNm5CbmwsXFxecPXsWR48eRYcOHVRdIhERfSS55xi1bdtWfG4aUVVw584dDB06FI0aNRJHiXr37o2///4bR44cQbt27WT66+vrQxAECILA0SIiogpG7jlGEydOxNSpU/H06VM4OTkV+cXfrFkzhRVHpEq3b9/G3LlzERYWhsIzzn369IGvry9at26t4uqIiEgZ5L6PUXEPtpRIJBAEARKJRJxzUZFwjhH9061btzB37lz88ssvYiDq168ffH194eTkpOLqiIiokMrnGAFAdHR0kdejR4/E/8ojODgYrVu3RrVq1VCjRg0MGDAAd+/elekjCAL8/f1hZWUFXV1ddO7cGbdu3ZLpk5mZiUmTJqF69erQ19dHv3798PTpU3k3jaq4GzduYNCgQXBwcMCuXbsgCAIGDhyIy5cv47fffmMoIiKqAuQORnXq1PnXlzxOnTqF77//Hn///TfCw8ORk5ODnj17yjxSISQkBEuWLMGKFSsQGRkJqVSKHj16iDeZBAAPDw/s27cPYWFhiIiIQEpKCvr06VMhR6+o7F27dg1ffPEFmjVrht27d0MQBHzxxRe4evUq9u7dC0dHR1WXSEREZUUoha1btwrOzs6CpaWlEBMTIwiCIISGhgr79+8vzepEL1++FAAIp06dEgRBEPLy8gSpVCrMnz9f7JORkSEYGRkJa9asEQRBEN69eydoamoKYWFhYp9nz54JampqwpEjR0r0vYmJiQIAITEx8aPqp4rl8uXLwoABAwQAAgBBIpEIX331lXD9+nVVl0ZERCWgjOO33CNGq1evhqenJz777DO8e/dOHJUxNjbG0qVLPyqkJSYmAsh/7AiQf9ouPj4ePXv2FPtoa2ujU6dOOHv2LADg0qVLyM7OluljZWUFe3t7sQ/RP126dAn9+vVDy5YtsX//fkgkEnz99de4ceMGfvnlFzg4OKi6RCIiUhG5g9Hy5cuxbt06zJo1C+rq6mJ7q1atcOPGjVIXIggCPD090aFDB9jb2wMA4uPjAQAWFhYyfS0sLMRl8fHx0NLSgomJyQf7vC8zMxNJSUkyL6r8IiMj0adPH7Rq1Qr/93//BzU1Nbi6uuLWrVvYuXMnmjZtquoSiYhIxUo1+bq4ORfa2toyc4Pk5ebmhuvXr2Pnzp1Flr3/4E2h4Aq4f/NvfYKDg2FkZCS+rK2tS103lX/nz5/HZ599hjZt2uDgwYNQU1PD0KFDERUVhe3bt6Nx48aqLpGIiMoJuYORra0trl69WqT98OHDaNKkSamKmDRpEg4cOIATJ06gVq1aYnvhY0feH/l5+fKlOIoklUqRlZWFhISED/Z5n4+PDxITE8XXkydPSlU3lW9nz56Fi4sLPvnkExw+fBjq6uoYPnw47ty5g61bt6Jhw4aqLpGIiMoZuYPRtGnT8P3334uXM1+4cAE//PADZs6ciWnTpsm1LkEQ4Obmhr179+LPP/+Era2tzHJbW1tIpVKEh4eLbVlZWTh16hScnZ0BAE5OTtDU1JTpExcXh5s3b4p93qetrQ1DQ0OZF1UeERER6NGjB9q3b4+jR49CXV0dI0eOxN27d7F582bY2dmpukQiIiqvSjNje+3atULt2rUFiUQiSCQSoVatWsL69evlXs+ECRMEIyMj4eTJk0JcXJz4SktLE/vMnz9fMDIyEvbu3SvcuHFDGDJkiGBpaSkkJSWJfcaPHy/UqlVLOHbsmHD58mWha9euQvPmzYWcnJwS1cGr0iqHkydPCl26dBGvMtPQ0BDGjBkjPHz4UNWlERGREijj+F2qYFTo1atXwosXL0r/5QUHsPdfmzZtEvvk5eUJfn5+glQqFbS1tYWOHTsKN27ckFlPenq64ObmJpiamgq6urpCnz59hNjY2BLXwWBUceXl5Ql//vmn0KlTJ/H/H01NTeG7774ToqOjVV0eEREpkTKO33I/EqTQy5cvcffuXUgkEjRs2BDm5uYfP3ylInwkSMUjCAL+/PNPBAQE4PTp0wAALS0tjB49GjNmzEDt2rVVXCERESmbMo7fcj9ENikpCd9//z127tyJvLw8AIC6ujoGDx6MlStXwsjISCGFERVHEASEh4cjMDAQZ86cAZAfiMaOHYvp06fzCkMiIvoock++HjNmDM6fP4+DBw/i3bt3SExMxO+//46LFy9i7NixyqiRCIIg4MiRI3B2doaLiwvOnDkDbW1tTJo0CY8ePcKKFSsYioiI6KPJfSpNX18ff/zxBzp06CDTfvr0afTq1euj7mWkKjyVVn4JgoBDhw4hMDAQFy5cAADo6Ohg/PjxmDZtGqysrFRcIRERqUq5OJVmZmZW7OkyIyOjInefJiotQRDw+++/IzAwEBcvXgQA6OrqYsKECZg2bZp4jysiIiJFkvtU2uzZs+Hp6Ym4uDixLT4+HtOmTcOcOXMUWhxVPYIg4LfffkOrVq3Qr18/XLx4EXp6evDy8kJ0dDQWL17MUEREREoj96k0R0dHPHjwAJmZmeKVP7GxsdDW1i5y47zLly8rrlIl4qk01cvLy8P+/fsRGBiIa9euAcg/bevm5oapU6dW6KseiYhIOcrFqbQBAwYo5IuJgPxAtHfvXsydOxfXr18HABgYGGDSpEnw9PRE9erVVVwhERFVJaW+j1FlwhGjspebm4tff/0Vc+fOxa1btwAA1apVw+TJkzFlyhSYmZmpuEIiIirvysWI0T+lpKSI9zIqxGBB/yY3Nxe//PIL5s6di9u3bwPIn7jv7u4ODw8PTuAnIiKVkjsYRUdHw83NDSdPnkRGRobYLggCJBIJcnNzFVogVQ65ubkICwvDvHnzcOfOHQCAsbExPDw84O7uDmNjY9UWSEREhFIEo2+++QYAsHHjRlhYWEAikSi8KKo8cnJysHPnTsybNw/37t0DAJiYmMDT0xOTJk3indKJiKhckTsYXb9+HZcuXULDhg2VUQ9VEjk5Ofj555/xww8/4MGDBwAAU1NTTJ06FW5ubjzlSkRE5ZLcwah169Z48uQJgxEVKzs7G9u2bcMPP/yAR48eAci/KaiXlxe+//57VKtWTcUVEhERfZjcwWj9+vUYP348nj17Bnt7e2hqasosb9asmcKKo4ojKysLW7ZsQVBQEGJiYgAA5ubmmDZtGiZMmAADAwPVFkhERFQCcgejV69e4eHDhxg5cqTYJpFIOPm6isrKysKmTZsQFBSE2NhYAECNGjXg7e2N8ePHQ19fX8UVEhERlZzcwWjUqFFwdHTEzp07Ofm6CsvMzMTGjRsRHByMJ0+eAACkUimmT5+O7777Dnp6eiqukIiISH5yB6PHjx/jwIEDqF+/vjLqoXIuIyMD69evx/z58/Hs2TMAgKWlJWbMmIGxY8dCV1dXxRUSERGVntzBqGvXrrh27RqDURWTnp6OdevWYcGCBXj+/DkAoGbNmpgxYwbGjBkDHR0dFVdIRET08eQORn379sWUKVNw48YNODg4FJl83a9fP4UVR6qXlpaGn376CSEhIYiPjwcAWFtbw8fHB6NGjYK2traKKyQiIlIcuZ+Vpqam9uGVVdDJ13xWWlGpqalYs2YNFi5ciBcvXgAAateujZkzZ2LEiBEMREREpHLl4llp7z8bjSqXlJQUrFq1CosWLcKrV68AADY2Npg1axaGDRsGLS0tFVdIRESkPB/1ENmMjAzOLakkkpOTsXLlSixevBivX78GANStWxezZs3C0KFDi5wyJSIiqow+fF7sA3JzczF37lzUrFkTBgYG4t2N58yZgw0bNii8QFKupKQkBAUFwcbGBj4+Pnj9+jXq1auHTZs24c6dOxg1ahRDERERVRlyB6MffvgBmzdvRkhIiMxpFQcHB6xfv16hxZHyJCYmYu7cueJpsrdv36JBgwbYunUr7ty5gxEjRjAQERFRlSP3qbStW7di7dq16NatG8aPHy+2N2vWDHfu3FFocaR47969w48//oilS5fi3bt3AICGDRtizpw5+Prrr6Gurq7aAomIiFRI7mD07NmzYu9hlJeXh+zsbIUURfJLTU0Vn0eWkpJS5FEcCQkJWLp0KX788UckJiYCABo3bow5c+Zg0KBBDEREREQoRTBq2rQpTp8+jTp16si07969G46OjgorjBTjzZs3CA0NxbJly5CcnAwgfx/6+vriyy+//NfbLxAREVU1JQ5Go0aNwo8//gg/Pz8MHToUz549Q15eHvbu3Yu7d+9i69at+P3335VZK8nh9evXWLJkCZYvX46UlBQA+fPAfH198fnnnzMQERERFaPEN3hUV1dHXFwcatSogT/++ANBQUG4dOkS8vLy0LJlS/j6+qJnz57KrlcpKuoNHlNTU2X+bGFhAQBwc3PDpk2bxOXNmzeHr68vBgwYwEBERESVhjKO3yUORmpqaoiPj0eNGjUU8sXlSUUNRhKJpET98vLyStyXiIioolDG8Vuu4QMeXCsm7jciIqKSkWvydYMGDf7zIPv27duPKohKrnDuEACcOHECffv2BQDEx8eLV6gRERFRyckVjAICAmBkZKSsWkhO/7wkv0uXLuKfDQwMilyuT0RERP9NrmD09ddfV8o5RkRERESAHHOMOE+FiIiIKrsSjxiV8OI1UhF9fX3uIyIioo9U4mCUl5enzDqIiIiIVI53+yMiIiIqwGBEREREVIDBiIiIiKgAgxERERFRAQYjIiIiogIMRkREREQFGIyIiIiICjAYERERERVgMCIiIiIqwGBEREREVIDBiIiIiKgAgxERERFRAQYjIiIiogIMRkREREQFVBqM/vrrL/Tt2xdWVlaQSCTYv3+/zHJBEODv7w8rKyvo6uqic+fOuHXrlkyfzMxMTJo0CdWrV4e+vj769euHp0+fluFWEBERUWWh0mCUmpqK5s2bY8WKFcUuDwkJwZIlS7BixQpERkZCKpWiR48eSE5OFvt4eHhg3759CAsLQ0REBFJSUtCnTx/k5uaW1WYQERFRJSERBEFQdREAIJFIsG/fPgwYMABA/miRlZUVPDw8MH36dAD5o0MWFhZYsGABxo0bh8TERJibm2Pbtm0YPHgwAOD58+ewtrbGoUOH4OLiUqLvTkpKgpGRERITE2FoaKiU7SMiIiLFUsbxu9zOMYqOjkZ8fDx69uwptmlra6NTp044e/YsAODSpUvIzs6W6WNlZQV7e3uxT3EyMzORlJQk8yIiIiIqt8EoPj4eAGBhYSHTbmFhIS6Lj4+HlpYWTExMPtinOMHBwTAyMhJf1tbWCq6eiIiIKqJyG4wKSSQSmfeCIBRpe99/9fHx8UFiYqL4evLkiUJqJSIiooqt3AYjqVQKAEVGfl6+fCmOIkmlUmRlZSEhIeGDfYqjra0NQ0NDmRcRERFRuQ1Gtra2kEqlCA8PF9uysrJw6tQpODs7AwCcnJygqakp0ycuLg43b94U+xARERGVlIYqvzwlJQUPHjwQ30dHR+Pq1aswNTVF7dq14eHhgaCgINjZ2cHOzg5BQUHQ09ODq6srAMDIyAijR4/G1KlTYWZmBlNTU3h5ecHBwQHdu3dX1WYRERFRBaXSYHTx4kV06dJFfO/p6QkAGD58ODZv3gxvb2+kp6dj4sSJSEhIQNu2bXH06FFUq1ZN/ExoaCg0NDQwaNAgpKeno1u3bti8eTPU1dXLfHuIiIioYis39zFSJd7HiIiIqOKpUvcxIiIiIiprDEZEREREBRiMiIiIiAowGBEREREVYDAiIiIiKsBgRERERFSAwYiIiIioAIMRERERUQEGIyIiIqICDEZEREREBRiMiIiIiAowGBEREREVYDAiIiIiKsBgRERERFSAwYiIiIioAIMRERERUQEGIyIiIqICDEZEREREBRiMiIiIiAowGBEREREVYDAiIiIiKsBgRERERFSAwYiIiIioAIMRERERUQEGIyIiIqICDEZEREREBRiMiIiIiAowGBEREREVYDAiIiIiKsBgRERERFSAwYiIiIioAIMRERERUQEGIyIiIqICDEZEREREBRiMiIiIiAowGBEREREVYDAiIiIiKsBgRERERFSAwYiIiIioAIMRERERUQEGIyIiIqICDEZEREREBRiMiIiIiAowGBEREREVYDAiIiIiKsBgRERERFSAwYiIiIioAIMRERERUQEGIyIiIqICDEZEREREBRiMiIiIiApUmmC0atUq2NraQkdHB05OTjh9+rSqSyIiIqIKplIEo127dsHDwwOzZs3ClStX8Omnn6J3796IjY1VdWlERERUgUgEQRBUXcTHatu2LVq2bInVq1eLbY0bN8aAAQMQHBz8n59PSkqCkZEREhMTYWhoqMxSiYiISEGUcfzWUMhaVCgrKwuXLl3CjBkzZNp79uyJs2fPFvuZzMxMZGZmiu8TExMB5P8FExERUcVQeNxW5BhPhQ9Gr1+/Rm5uLiwsLGTaLSwsEB8fX+xngoODERAQUKTd2tpaKTUSERGR8rx58wZGRkYKWVeFD0aFJBKJzHtBEIq0FfLx8YGnp6f4/t27d6hTpw5iY2MV9hdLpZOUlARra2s8efKEpzVVjPuifOH+KD+4L8qPxMRE1K5dG6ampgpbZ4UPRtWrV4e6unqR0aGXL18WGUUqpK2tDW1t7SLtRkZG/J+8nDA0NOS+KCe4L8oX7o/yg/ui/FBTU9y1ZBX+qjQtLS04OTkhPDxcpj08PBzOzs4qqoqIiIgqogo/YgQAnp6eGDp0KFq1aoV27dph7dq1iI2Nxfjx41VdGhEREVUglSIYDR48GG/evEFgYCDi4uJgb2+PQ4cOoU6dOiX6vLa2Nvz8/Io9vUZli/ui/OC+KF+4P8oP7ovyQxn7olLcx4iIiIhIESr8HCMiIiIiRWEwIiIiIirAYERERERUgMGIiIiIqECVCUarVq2Cra0tdHR04OTkhNOnT/9r/1OnTsHJyQk6OjqoW7cu1qxZU0aVVn7y7Iu4uDi4urqiYcOGUFNTg4eHR9kVWgXIsy/27t2LHj16wNzcHIaGhmjXrh3++OOPMqy2cpNnX0RERKB9+/YwMzODrq4uGjVqhNDQ0DKstvKT95hR6MyZM9DQ0ECLFi2UW2AVIs++OHnyJCQSSZHXnTt3Sv6FQhUQFhYmaGpqCuvWrROioqIEd3d3QV9fX3j8+HGx/R89eiTo6ekJ7u7uQlRUlLBu3TpBU1NT+PXXX8u48spH3n0RHR0tTJ48WdiyZYvQokULwd3dvWwLrsTk3Rfu7u7CggULhAsXLgj37t0TfHx8BE1NTeHy5ctlXHnlI+++uHz5srBjxw7h5s2bQnR0tLBt2zZBT09P+Omnn8q48spJ3v1R6N27d0LdunWFnj17Cs2bNy+bYis5effFiRMnBADC3bt3hbi4OPGVk5NT4u+sEsGoTZs2wvjx42XaGjVqJMyYMaPY/t7e3kKjRo1k2saNGyd88sknSquxqpB3X/xTp06dGIwU6GP2RaEmTZoIAQEBii6tylHEvhg4cKDw7bffKrq0Kqm0+2Pw4MHC7NmzBT8/PwYjBZF3XxQGo4SEhFJ/Z6U/lZaVlYVLly6hZ8+eMu09e/bE2bNni/3MuXPnivR3cXHBxYsXkZ2drbRaK7vS7AtSDkXsi7y8PCQnJyv04Y1VkSL2xZUrV3D27Fl06tRJGSVWKaXdH5s2bcLDhw/h5+en7BKrjI/52XB0dISlpSW6deuGEydOyPW9leLO1//m9evXyM3NLfJAWQsLiyIPni0UHx9fbP+cnBy8fv0alpaWSqu3MivNviDlUMS+WLx4MVJTUzFo0CBllFhlfMy+qFWrFl69eoWcnBz4+/tjzJgxyiy1SijN/rh//z5mzJiB06dPQ0Oj0h9Wy0xp9oWlpSXWrl0LJycnZGZmYtu2bejWrRtOnjyJjh07luh7q8welEgkMu8FQSjS9l/9i2sn+cm7L0h5Srsvdu7cCX9/f/z222+oUaOGssqrUkqzL06fPo2UlBT8/fffmDFjBurXr48hQ4Yos8wqo6T7Izc3F66urggICECDBg3KqrwqRZ6fjYYNG6Jhw4bi+3bt2uHJkydYtGgRg1Gh6tWrQ11dvUi6fPnyZZEUWkgqlRbbX0NDA2ZmZkqrtbIrzb4g5fiYfbFr1y6MHj0au3fvRvfu3ZVZZpXwMfvC1tYWAODg4IAXL17A39+fwegjybs/kpOTcfHiRVy5cgVubm4A8k8zC4IADQ0NHD16FF27di2T2isbRR0zPvnkE/z8888l7l/p5xhpaWnByckJ4eHhMu3h4eFwdnYu9jPt2rUr0v/o0aNo1aoVNDU1lVZrZVeafUHKUdp9sXPnTowYMQI7duzA//73P2WXWSUo6udCEARkZmYqurwqR979YWhoiBs3buDq1avia/z48WjYsCGuXr2Ktm3bllXplY6ifjauXLki3xSYUk/brkAKL/fbsGGDEBUVJXh4eAj6+vpCTEyMIAiCMGPGDGHo0KFi/8LL9adMmSJERUUJGzZs4OX6CiLvvhAEQbhy5Ypw5coVwcnJSXB1dRWuXLki3Lp1SxXlVyry7osdO3YIGhoawsqVK2Uug3337p2qNqHSkHdfrFixQjhw4IBw79494d69e8LGjRsFQ0NDYdasWarahEqlNL+n/olXpSmOvPsiNDRU2Ldvn3Dv3j3h5s2bwowZMwQAwp49e0r8nVUiGAmCIKxcuVKoU6eOoKWlJbRs2VI4deqUuGz48OFCp06dZPqfPHlScHR0FLS0tAQbGxth9erVZVxx5SXvvgBQ5FWnTp2yLbqSkmdfdOrUqdh9MXz48LIvvBKSZ18sW7ZMaNq0qaCnpycYGhoKjo6OwqpVq4Tc3FwVVF45yft76p8YjBRLnn2xYMECoV69eoKOjo5gYmIidOjQQTh48KBc3ycRhIJZxURERERVXKWfY0RERERUUgxGRERERAUYjIiIiIgKMBgRERERFWAwIiIiIirAYERERERUgMGIiIiIqACDERFVGSNGjMCAAQNK9dmOHTtix44dH/X9K1asQL9+/T5qHUSkXAxGRKRQHxM+FCUmJgYSiQRXr15VyPp+//13xMfH4+uvv/6o9YwdOxaRkZGIiIhQSF1EpHgMRkRE/2HZsmUYOXIk1NQ+7lemtrY2XF1dsXz5cgVVRkSKxmBERGUqKioKn332GQwMDGBhYYGhQ4fi9evX4vLOnTtj8uTJ8Pb2hqmpKaRSKfz9/WXWcefOHXTo0AE6Ojpo0qQJjh07BolEgv379wMAbG1tAQCOjo6QSCTo3LmzzOcXLVoES0tLmJmZ4fvvv0d2dvYH6339+jWOHTtW5BSYRCLBTz/9hD59+kBPTw+NGzfGuXPn8ODBA3Tu3Bn6+vpo164dHj58KPO5fv36Yf/+/UhPT5fzb46IygKDERGVmbi4OHTq1AktWrTAxYsXceTIEbx48QKDBg2S6bdlyxbo6+vj/PnzCAkJQWBgIMLDwwEAeXl5GDBgAPT09HD+/HmsXbsWs2bNkvn8hQsXAADHjh1DXFwc9u7dKy47ceIEHj58iBMnTmDLli3YvHkzNm/e/MGaIyIixODzvrlz52LYsGG4evUqGjVqBFdXV4wbNw4+Pj64ePEiAMDNzU3mM61atUJ2drZYIxGVLxqqLoCIqo7Vq1ejZcuWCAoKEts2btwIa2tr3Lt3Dw0aNAAANGvWDH5+fgAAOzs7rFixAsePH0ePHj1w9OhRPHz4ECdPnoRUKgUA/PDDD+jRo4e4TnNzcwCAmZmZ2KeQiYkJVqxYAXV1dTRq1Aj/+9//cPz4cYwdO7bYmmNiYmBhYVHsabSRI0eKoW769Olo164d5syZAxcXFwCAu7s7Ro4cKfMZfX19GBsbIyYmBp06dSr5Xx4RlQmOGBFRmbl06RJOnDgBAwMD8dWoUSMAkDnl1KxZM5nPWVpa4uXLlwCAu3fvwtraWibwtGnTpsQ1NG3aFOrq6sWuuzjp6enQ0dEpdtk/67SwsAAAODg4yLRlZGQgKSlJ5nO6urpIS0srcc1EVHY4YkREZSYvLw99+/bFggULiiyztLQU/6ypqSmzTCKRIC8vDwAgCAIkEkmpa/i3dRenevXqSEhI+M91FdZUXNv763/79q04qkVE5QuDERGVmZYtW2LPnj2wsbGBhkbpfv00atQIsbGxePHihThKExkZKdNHS0sLAJCbm/txBSN/And8fDwSEhJgYmLy0et7+PAhMjIy4Ojo+NHrIiLF46k0IlK4xMREXL16VeYVGxuL77//Hm/fvsWQIUNw4cIFPHr0CEePHsWoUaNKHGJ69OiBevXqYfjw4bh+/TrOnDkjTr4uHKGpUaMGdHV1xcndiYmJpd4WR0dHmJub48yZM6Vexz+dPn0adevWRb169RSyPiJSLAYjIlK4kydPwtHRUebl6+sLKysrnDlzBrm5uXBxcYG9vT3c3d1hZGRU4nsEqaurY//+/UhJSUHr1q0xZswYzJ49GwDEuUAaGhpYtmwZfvrpJ1hZWaF///6l3hZ1dXWMGjUK27dvL/U6/mnnzp0fnOhNRKonEQRBUHURREQf48yZM+jQoQMePHiglJGYFy9eoGnTprh06RLq1KlT6vXcvHkT3bp1w71792BkZKTAColIURiMiKjC2bdvHwwMDGBnZ4cHDx7A3d0dJiYmSn3Uxm+//QZTU1N8+umnpV7H0aNHIQiCeDk/EZU/DEZEVOFs3boVc+fOxZMnT1C9enV0794dixcvhpmZmapLI6IKjsGIiIiIqAAnXxMREREVYDAiIiIiKsBgRERERFSAwYiIiIioAIMRERERUQEGIyIiIqICDEZEREREBRiMiIiIiAowGBEREREV+H9HI+muHps/IgAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Step 4: visualize results after advancing in time\n", "fig = plt.figure()\n", "plt.plot(x, T, '-k', label='Numerical solution {}'.format(n+1))\n", "plt.scatter(x, TA-(TA-TB)/L*x, s=50, c='k', marker='+', label='Analytical solution')\n", "plt.xlabel('Length (m)')\n", "plt.ylabel('Temperature (degree Celsius)')\n", "ax = plt.gca()\n", "ax.set_xlim(0, 0.5)\n", "ax.set_ylim(0, 600)\n", "plt.legend(loc='upper left')\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "```{note}\n", "The variable $n$ here represents a *pseudo-time* that advances the simulation to converge towards the analytical solution.\n", "```" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Exercise\n", "\n", "1. Try to do more iterations and observe the final temperature difference.\n", "\n", "2. Modify the code so that a small tolerance for the error (e.g., 10$^{-5}$) can be specified to obtain a convergent result." ] } ], "metadata": { "kernelspec": { "display_name": "cfd-python", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.10.14" } }, "nbformat": 4, "nbformat_minor": 2 }