{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Отчет по лабораторной работе №6. Задача Коши и решение систем ДУ.\n",
"## Выполнил Сарафанов Ф.Г., 430 группа/2018, 19 вариант\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Постановка задачи\n",
"\n",
"Рассмотрим численное решение уравнения \n",
"\n",
"$$\n",
"y''+y'-4y=t\\cdot e^{-t},\\quad y(0)=1, y'(0)=0,\\quad t\\in[0,2]\n",
"$$\n",
"\n",
"запишем его в виде системы ДУ 1-го порядка:\n",
"\n",
"$$\n",
"\\begin{cases}\n",
"y'=z\\\\\n",
"z'=-z+4y+t\\cdot e^{-t}\n",
"\\end{cases}%,\\qquad i=2,3,\\ldots, n\n",
"$$\n",
"\n",
"В программе используется векторный подход: массив переменных Y=[y z], y=Y[0], z=Y[1]
, и система имеет вид\n",
"\n",
"$$\n",
"\\begin{cases}\n",
"Y[0]'=Y[1]\\\\\n",
"Y[1]'=-Y[1]+4\\cdot Y[0]+t\\cdot e^{-t}\n",
"\\end{cases}\n",
"$$\n",
"\n",
"Здесь сопоставляется вектор производных и вектор правых частей\n",
"\n",
"$$\n",
"Y'=f(t,Y),\\; \\text{где}\\;f(t,Y)=\n",
"\\begin{bmatrix}\n",
"Y[1]\\\\\n",
"-Y[1]+4\\cdot Y[0]+t\\cdot e^{-t}\n",
"\\end{bmatrix}\n",
"$$\n",
"\n",
"## Метод Тейлора 3-го порядка\n",
"\n",
"Можно показать, что в векторном виде итерационный процесс решения методом Тейлора 3-го порядка будет записываться в виде\n",
"\n",
"$$\n",
"Y^{n+1}=Y^n+hf+\\frac{h^2}{2}\\left(\\frac{\\partial f}{\\partial t}+Jf\\right)+\\frac{h^3}{6}\\left(\\frac{\\partial ^2f}{\\partial t^2}+\\frac{\\partial J}{\\partial t}f +J\\frac{\\partial f}{\\partial t}\\right)\n",
"$$\n",
"\n",
"Где все функции (вектора, матрицы) берутся в точке $Y^n, t^n$."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Импорт библиотек"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"from sympy import Symbol,diff,exp,cos,sin,Matrix,var,symbols\n",
"import sympy as sym\n",
"import numpy as np\n",
"from matplotlib import pyplot as plt\n",
"from scipy.integrate import solve_ivp\n",
"from scipy import interpolate\n",
"from numpy.linalg import norm"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Зададим условия расчетов:"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"n=2 # количество переменных\n",
"y00=np.array([1,0]) # начальные условия y(0)=1,y'(0)=0\n",
"Tn0,tmax0=0,2 # временной интервал\n",
"h=0.01 # начальный шаг\n",
"eps=0.01 # точность"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Задаем символьные переменные для расчета всех необходимых величин:"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"t=Symbol('t')\n",
"y = sym.Matrix(n, 1, lambda i,j:var('y[%d]' % (i)))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Задаем правую часть уравнения"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"def f(T,y):\n",
" return Matrix([\\\n",
" y[1],\\\n",
" -y[1]+4*y[0]+T*exp(-T),\\\n",
" ])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Чтобы не загружаться излишними символьными вычислениями, вычислим символьно все нужные функции и преобразуем в численный вид"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"A='def pfpt(t,y):\\n\\treturn np.array(%s)'%(str(np.array2string(np.array(f(t,y).diff(t)).T[0],separator=', ')).replace('Matrix(',''))+'\\n'+\\\n",
" 'def p2fpt2(t,y):\\n\\treturn np.array(%s)'%(str(np.array2string(np.array(f(t,y).diff(t).diff(t)).T[0],separator=', ')).replace('Matrix(',''))+'\\n'+\\\n",
" 'def J(t,y):\\n\\treturn np.array(%s)'%(str(f(t,y).jacobian(y)).replace('Matrix(',''))+'\\n'+\\\n",
" 'def model(t,y):\\n\\treturn np.array(%s)'%(str(np.array2string(np.array(f(t,y)).T[0],separator=', ')).replace('Matrix(',''))+'\\n'+\\\n",
" 'def pJpt(t,y):\\n\\treturn np.array(%s)'%(str(f(t,y).jacobian(y).diff(t)).replace('Matrix(',''))\n",
"\n",
"exec(A.replace('))',')',100))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Зададим массив времени и массив координат решения:"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
"Time=np.array([])\n",
"Y = np.empty((0,n),np.float32)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Итерационный процесс:"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"итерация 1 401\n",
"итерация 2 801\n",
"итерация 3 1601\n",
"погрешность 0.02829945\n",
"итерация 4 3200\n",
"погрешность 0.020543108\n",
"итерация 5 6401\n",
"погрешность 0.014714718\n",
"итерация 6 12800\n",
"погрешность 0.010471105\n",
"итерация 7 25601\n",
"погрешность 0.0074286712\n"
]
}
],
"source": [
"N=0\n",
"while True:\n",
" h=h/2\n",
"\n",
" yn=y00\n",
" Tn=Tn0\n",
" tmax=tmax0 \n",
"\n",
" N+=1\n",
"\n",
" OldTime=Time\n",
" OldY=Y\n",
"\n",
" Time=np.array([])\n",
" Y = np.empty((0,n),np.float32)\n",
"\n",
" while Tn<=tmax:\n",
"\n",
" Tn=Tn+h\n",
" yo=yn\n",
" Time=np.append(Time,Tn)\n",
"\n",
" Y=np.append(Y,[yo],axis=0)\n",
"\n",
" const_pfpt=pfpt(Tn,yn)\n",
" const_J=J(Tn,yn)\n",
" const_F=model(Tn,yn)\n",
" yn=np.array(yn+h*(const_F+h/2*(const_pfpt+np.dot(const_J,const_F))+h**2/6*(p2fpt2(Tn,yn)+\\\n",
" np.dot(pJpt(Tn,yn),const_F)+np.dot(const_J,const_pfpt))))\n",
" print('итерация',N,len(Y))\n",
" if N>=3:\n",
" YY=Y\n",
" YO=OldY\n",
" if len(YY)%2!=0:\n",
" YY=YY[:-1]\n",
" if len(YO)%2!=0:\n",
" YO=YO[:-1]\n",
" YY=np.array(YY,dtype=np.float32)\n",
" YO=np.array(YO,dtype=np.float32)\n",
" delta=np.max(norm(YY[0::2]-YO,axis=0))\n",
" print('погрешность',delta)\n",
" if N>10*2:\n",
" break\n",
" if delta"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"f, (ax1, ax2,ax3) = plt.subplots(1,3, gridspec_kw = {'width_ratios':[3, 3,3]})\n",
"f.set_size_inches(11,3)\n",
"\n",
"ax1.set_title('Графики решений')\n",
"ax1.plot(Time,Y1,'*', markersize=1)\n",
"ax1.plot(Time2,Z1)\n",
"ax1.set_ylabel(\"y\")\n",
"ax1.set_xlabel(\"t\")\n",
"\n",
"ax2.plot(Time,Y2,'*', markersize=1)\n",
"ax2.plot(Time2,Z2)\n",
"ax2.set_ylabel(\"y'\")\n",
"ax2.set_xlabel(\"t\")\n",
"\n",
"ax3.plot(Y1,Y2,'*', markersize=1)\n",
"ax3.plot(Z1,Z2)\n",
"ax3.set_xlabel(\"y\")\n",
"ax3.set_ylabel(\"y'\")\n",
"\n",
"f.tight_layout()\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Чтобы не возиться с вычитанием массивов разной длины, можно решения разными способами простейшим способом интерполировать на нужном отрезке (соединяя точки) и построить разность интерполирующих функций на одном и том же разбиении:"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [],
"source": [
"x = np.linspace(np.min(Time), np.max(Time),1000)\n",
"f = interpolate.interp1d(Time, Y1,axis=0, fill_value=\"extrapolate\")\n",
"Y1f=f(x)\n",
"f = interpolate.interp1d(Time2, Z1,axis=0, fill_value=\"extrapolate\")\n",
"Z1f=f(x)\n",
"\n",
"x = np.linspace(np.min(Time), np.max(Time),1000)\n",
"f = interpolate.interp1d(Time, Y2,axis=0, fill_value=\"extrapolate\")\n",
"Y2f=f(x)\n",
"f = interpolate.interp1d(Time2, Z2,axis=0, fill_value=\"extrapolate\")\n",
"Z2f=f(x)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Разностные графики:"
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {
"scrolled": false
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAxAAAADQCAYAAACX+YfUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAIABJREFUeJzs3Xl8XHW9//HXJ2vTfUm6t3QvdGEtZd8KSAG1KiAFF/CCuIDoxQ30XkQWBUQQf4AKgiAXLb2AUhBlR7a2tKULdN/bdEvTpGmSZpuZz++POeEOIWnS5EwySd7PxyOPnDnnez7ne6ZpTj7z3czdERERERERaYq0tq6AiIiIiIi0H0ogRERERESkyZRAiIiIiIhIkymBEBERERGRJlMCISIiIiIiTaYEQkREREREmkwJhIiIiIiINJkSCGl1ZrbJzCrMrMzMdpnZn8yse1vXS8DMLjezt9u6HiIiIpK6lEBIW/mMu3cHjgaOBf6rjesjIiIiIk2gBELalLtvA/4JTAIws6+Z2UozKzWzDWb2jdqyZnaSmW0MWi62mtnVCcfeMLMrE16fZWabEl4fFpTZa2bLzeyzCcdyzOzXZrbZzErM7O1g39LgWhVmFgu2y8zsJ8F5bmZjEuLcamaPBtsjguMZde/ZzPLN7PRgO83Mrjez9Wa2x8xmm1nf+t4rMzu9Tj3KzOzbwbFNZnaDma0ws+KgVadLcKyPmT1vZruDY8+b2dCEuEcE520ATgKyzexZMys0s9sTyt1kZv+T8PqBxPfAzB41s1uD7X5BzG8l1D0/4dwvBud+9G8mIiIi7YMSCGlTZjYMOA9YHOwqAD4N9AS+BtxjZkcHx9YApwQtF58FfmVmvZpwjUzgOeAloD/wHeAJMxsfFLkLOAY4EegL/AiIufsRwbXOBba7e/fg6xctve8E1wKfA04DBgPFwP0HKJ9Yj+7u/kDCsS8B5wCjgXH8X6tOGvAn4BBgOFAB3Jdw3uPA74EjgBHAYcAdwOHARYnJVi0zG0v8ffmEoDvaP4G/uPvv6jmeCdwC7DjAfYqIiEiKUgIhbeXvZrYXeBv4N/ALAHf/h7uv97h/E/+j/5Tg2G53r/0U24DVQHkTrnU80B243d2r3f014HngEjNLA/4D+K67b3P3qLu/6+5VId7rgXwD+Km75wfXvAm4sL6Wiya4z923unsRcBtwCYC773H3p919v7uXBsdOAzCzkcB44IHg2BPAkuA92E48ufhCPdf6JfEkoK5s4O/AKne/9QD3PJ94QigiIiLtjBIIaSufc/fe7n6Iu3/b3SsAzOxcM5tnZkVBgnEekFt7kpmdbGalxP8A/bu7RxJi/jboorSX+B+xtQYDW909lrBvMzAkiN0FWN/M+3g/4Zo/qOd4YdBtaKWZfbme44cAf0uIsRKIAgOaUZetCdubid83ZtbVzP4QdNHaB7wJ9Daz9OA6e+u8j4kKgIGJO8zsOOBQ4LF6yl8NdAVOMLOcugfNrAfxFp7/Pqg7ExERkZShBEJShpllA08T71I0wN17Ay8Qb20AwN3fdvcewATgW2aW2I3m2iAp6U28W1Ct7cCwoLWh1nBgG1AIVBLv9tMcRydc8656jue6ex/gGuBR++RsU1uBc2tjBF9dgrEhB2tYwvZw4vcN8H3irQzHuXtP4NRgvwG7iScTDbV49Ad21dl3J3C9u0frKf9uEH8B8ZaOun4IzHb3zY3ci4iIiKQoJRCSSrKId4HZDUSC5OBTtQfNbFTCH+DZxH9+K5oQdz7xrk4/MrPMYADzZ4BZQavEI8DdZjbYzNLN7IQgmQlTMfE/2K3O/t8Dt5nZIQBmlmdmM5p5javNbGgwCPsnwJPB/h7E36e9wbGfJZyzIfj6ZtA6cClwhJmdaGaDgC8THz9Saxrg7v58A3WYF7RmfId4F7ETEo71ID6upb7EQkRERNoJJRCSMoI++NcCs4n/wX0pMCehyOnAGjMrIz6G4R53f6MJcauJD7o+l3iLwwPAV919VVDkB8AHxD81LyI+gDis/xubgtmHZgNXBfeY6F7i9/hS0DVrHnBcM6/1F+JjRmqTgtoxCL8Bcojf+zzgX7UnuLsDXyX+B/9SYBOwCvgx8GFQt6cTrjGIeBekA3L3PUHMR2pngyI+MP637l7crLsTERGRlGDxvx9EpD2z+JS1V7r7Ky2Mc3kQ5+Qw6iUiIiIdj1ogRERERESkyZRAiIiIiIhIkymBEOkA3H1ES7svBXEeVfclSTYzm25mq81snZldX8/xbDN7Mjg+38xGJBy7Idi/2szOCfZ1MbP3LL56/HIz+3lC+UctvoL9kuDryNa4RxGRjqw5i1WJiIg0S7D+yP3A2UA+sMDM5rj7ioRiVwDF7j7GzGYSn9jgYjObAMwEJhJf5+QVMxsHVAHT3L0sWOn8bTP7p7vPC+L90N2fap07FBHp+Dp1ApGbm+sjRoxo62qIiLS6RYsW7QPmuvv0Vr70VGCdu28AMLNZwAwgMYGYQXxVdoCngPvMzIL9s4JV2zea2TpgqrvPBcqC8pnBV7NnCNGzQUQ6q0WLFhW6e15j5Tp1AjFixAgWLlzY1tUQEWl1Zra2DZIHiK8An7hqej6fnLr4ozLuHjGzEqBfsH9enXOHwEctG4uAMcD97j4/odxtZnYj8CrxRRCr6lbKzK4CrgIYPny4ng0i0imZWZMWetUYCBERaU11F1OET7YWNFSmwXPdPeruRwJDgalmNik4fgNwKHAs0Jf4GiefDOL+oLtPcfcpeXmNfvgmItKpKYEQEZHWlA8MS3g9FNjeUBkzywB6EV/ksdFz3X0v8AYwPXi9w+OqgD8R70IlIiItoARCRERa0wJgrJmNNLMs4oOi59QpMwe4LNi+EHgtWDV9DjAzmKVpJDAWeM/M8sysN4CZ5QBnEV9RHTMbFHw34HPEV1gXEZEW6NRjIEREpHUFYxquAV4E0oFH3H25md0MLHT3OcDDwOPBIOki4kkGQbnZxAdcR4Cr3T0aJAmPBeMg0oDZ7v58cMknzCyPePenJcA3W+9uRUQ6JiUQIiLt0NtrC5m7oZAffGo88Q/X2w93fwF4oc6+GxO2K4GLGjj3NuC2OvuWAUc1UH5aS+srItIeuDu/eWUtRw3vzenj+yf1WurCJCLSDl32p/e4//X17NxX2dZVERGRFFBQWsW9r67lu7OWJP1aaoEQEWlH3J17Xl5DNBafuKiqJtbGNRIRkbZWVF7NQ29uAKCkoibp11MCISLSjqzfXc5vX1v30evKSLQNayMiIqngR08t5ZWVBa12PSUQIiLtxNpdpdz03PKP7VMLhIhI57ajpIJFm4tb9ZpKIERE2okr/7yQzXv2f2xfVUQJhIhIZxWLOSf88rVWv64GUYuItAOPvL3xE8kDQJW6MImIdFq7y6ra5LpKIEREUlzJ/hpufn5Fvccq1YVJRKRTenbJNo77xattcu2kJhBmNt3MVpvZOjO7vp7j2Wb2ZHB8vpmNSDh2Q7B/tZmdE+wbZmavm9lKM1tuZt9NKN/XzF42s7XB9z7JvDcRkdawtWg/1/z1/QaPqwVCRKRzas1B03UlLYEIVgS9HzgXmABcYmYT6hS7Aih29zHAPcAdwbkTiK88OhGYDjwQxIsA33f3w4DjgasTYl4PvOruY4FXg9ciIu3aT//+IW+tLWzwuAZRi4h0Plc+tpDnlm5vs+snswViKrDO3Te4ezUwC5hRp8wM4LFg+yngTIsvqToDmOXuVe6+EVgHTHX3He7+PoC7lwIrgSH1xHoM+FyS7ktEpFU89OYG3lyz+4BlqqNKIEREOptXVu5q0+snM4EYAmxNeJ3P//2x/4ky7h4BSoB+TTk36O50FDA/2DXA3XcEsXYA9a7hbWZXmdlCM1u4e/eBH8wiIm3F3bnthZWNlosEC8qJiEjHt7VoP+fd+1aj5dyT+2xIZgJh9eyrezcNlTnguWbWHXga+J677zuYSrn7g+4+xd2n5OXlHcypIiKtYlNhOZ+9750mlY2qBUJEpNOYs3Q7K3Y0/qdvsj9bSmYCkQ8MS3g9FKjbWeujMmaWAfQCig50rpllEk8ennD3ZxLK7DKzQUGZQUDbjSwREWmBh97awAfbSppUVi0QIiKdw+wFW/nVi6ubVDYSS+6HS8lMIBYAY81spJllER8UPadOmTnAZcH2hcBrHm9zmQPMDGZpGgmMBd4Lxkc8DKx097sPEOsy4NnQ70hEJMl+8L9LeWL+liaXj7bDBCIJM/R1MbP3zGxpMEPfzxPKjwxirA1iZrXGPYqIhO1HTy9rctlkPxuSlkAEYxquAV4kPth5trsvN7ObzeyzQbGHgX5mtg64jmDmJHdfDswGVgD/Aq529yhwEvAVYJqZLQm+zgti3Q6cbWZrgbOD1yIi7UZZVYSnFuUf1DntrQUiSTP0VQHT3P0I4EhgupkdH8S6A7gnmKGvOIgtItJuVNZEufaviw/qnGQnEBnJDO7uLwAv1Nl3Y8J2JXBRA+feBtxWZ9/b1D8+AnffA5zZwiqLiLSJFdv3cd5vGx8YV1c7bIH4aIY+ADOrnaEvcaW8GcBNwfZTwH11Z+gDNgYfPk1197lAWVA+M/jy4JxpwKXBsceCuL9Lzq2JiITvnXWFzDnIKVvbbQuEiIg0TVlVhMfnbWrWue2tBYIkzdBnZulmtoT4+LeX3X1+cM7eIEZD1yI4XzP0iUjKeXttId96ouHFRBuS7GdDUlsgRESkcafc8RrF+2uadW6k/c3ClJQZ+oJurkeaWW/gb2Y2CahvovR6n6ru/iDwIMCUKVPaXVYmIh3Tlx+e33iheqgFQkSkA3vm/fxmJw/QLrswJWWGvlruvhd4g/gYiUKgdxCjoWuJiKScWMz58VNNHzRdlxIIEZEOaln+Xq6bvbRFMdphF6ZkzNCXF7Q8YGY5wFnAquCc14MYoBn6RKSd+GBbCU8u3Np4wQa060HUIiJSv0Wbi7jisYUtjtPeWiDcPWJmtTP0pQOP1M7QByx09znEZ+h7PBgkXUQ8ySAoVztDX4Rghr5g7Z/HghmZ0ojP+vd8cMkfA7PM7FZgcRBbRCRlzV2/hy/9cV6LYmgMhIhIB+PuXPC7uaHESvZiQcmQhBn6lgFHNVB+A/GZn0RE2oVLHmpZ8gAQbccLyYmISB3VkRiH//yl0OK1txYIERGpn7vz5T82b9B0XclugVACISLSSqIx538XbaW0MtJ44SaqiSqBEBHpCLYU7eftdYWhxNIYCBGRDuJ7Ty7huYNcDKgxaoEQEWn/nlu6ne8c5GrTB6JZmEREOoCnF+WHnjxAu5yFSUREEpRXRfjZnOWhxtQgahGRdq5kfw3f/9+WTdfakGQPlBMRkeRxdyb+7MXQ46oFQkSkHZu7fg9H3BzeoOlE/Xtk06drVlJii4hIcrk7L6/YlZTYGgMhItJObd9bwTceb/laD/XZdPv5SYkrIiKt46G3NvCLF1YlJbYSCBGRdsjdOfH215IS+7tnjk1KXBERaR2LNhfzy38mJ3kAjYEQEWl3SvbX8JO/f5CU2F8/ZST/efa4pMQWEZHWccHv3g095qjcbtx98ZFkphvD+3YNPX4iJRAiIiGqica49R8r+MeyHaHH/ul5h/Hl4w8JPa6IiLSOsqoIP3kmOR8wTZ80kCOH9U5K7LqUQIiIhOjce99iXUFZUmJ//dRRSYkrIiKt4+G3NjInCVN6P/+dkzl0YI/Q4zZECYSISAiiMedP72xMSvIwpn93nrjyuNDjiohI67n+6WXMWrA1KbEnDemVlLgN0TSuIiIt5O78bfE2bv3HyqTEv376oQzo2SUpsVubmU03s9Vmts7Mrq/neLaZPRkcn29mIxKO3RDsX21m5wT7hpnZ62a20syWm9l3E8rfZGbbzGxJ8HVea9yjiEhd6wrKkpI89OmayZIbzw49bmPUAiEi0kI3PPNB0j5V6kjTtZpZOnA/cDaQDywwsznuviKh2BVAsbuPMbOZwB3AxWY2AZgJTAQGA6+Y2TggAnzf3d83sx7AIjN7OSHmPe5+V+vcoYjIJ60rKOWsu99MSuxvnDaa3m2wHpASCBGRZqpteUhW8vCrCw9PStw2NBVY5+4bAMxsFjADSEwgZgA3BdtPAfeZmQX7Z7l7FbDRzNYBU919LrADwN1LzWwlMKROTBGRNrF2Vyk/enpZUmJv/OV5xH89tj51YRIRaaYPt+3jutlLkxL72jPHctGUYUmJ3YaGAInZVn6wr94y7h4BSoB+TTk36O50FDA/Yfc1ZrbMzB4xsz4tvwURkaY7+543Wbxlb+hxxw3o3mbJAyQ5gQi7r2uw/xEzKzCzD+vEUl9XEWk197y8hs/c93ZSYq+97Vyu65hrPdT3tKu72lFDZQ54rpl1B54Gvufu+4LdvwNGA0cSb6X4dYMVM7vKzBaa2cLdu3c3fAciIk0QjTnn//atpMV/6T9PS1rspkhaApHQ1/VcYAJwSdCHNdFHfV2Be4j3daVOX9fpwANBPIBHg331ucfdjwy+XgjzfkREaj3zfj73vro2KbEf+uoUMtM7bONwPpDYrDIUqDuf4UdlzCwD6AUUHehcM8sknjw84e7P1BZw913uHnX3GPAQ8S5U9XL3B919irtPycvLa+btiYhAJBrj6UX5LN++r/HCB+niKcOYc81Jocc9WMl8Sn3U19Xdq4Havq6JZgCPBdtPAWfW7evq7huBdUE83P1N4g8TEZFWt2F3WVK6LfXpmsmdFxzO2RMGhB47hSwAxprZSDPLIv5B0Zw6ZeYAlwXbFwKvubsH+2cGLdcjgbHAe8Ez42FgpbvfnRjIzAYlvPw88LGWaxGRZLjmL4uTNu7hR9PHc/jQ1lks7kCSOYi6vv6qdScy/1hfVzNL7Os6r865dfvJ1ucaM/sqsJD4rBzFdQuY2VXAVQDDhw9v2p2IiAD/79W1/PrlNUmJ/a/vndphpmptSPB7/hrgRSAdeMTdl5vZzcBCd59DPBl4PBgkXUQ8ySAoN5v44OgIcLW7R83sZOArwAdmtiS41E+CVug7zexI4l2dNgHfaLWbFZFO6bF3N/Gv5TuTEvuV606jX/fspMQ+WMlMIJLW17UBvwNuCcrdQryv6398Ioj7g8CDAFOmTGkspogIAA+/vTFpycM710/r8MlDreAP+xfq7LsxYbsSuKiBc28Dbquz723qf2bg7l9paX1FRJqqsKyKn81ZHnrcAT2z+eE5hzKmf/fQYzdXMhOIg+nrmt/Uvq4Ncfddtdtm9hDwfLNrLiIScHfmbtjDLc+HPyvo+ZMH8Y3TRjGkd07osUVEpPU8835+0mbl++NXj2Xy0NZdaboxyRwDEXpf1wNdTH1dRSQZfvXiai59aH7jBQ9Sj+wMfj5jYkr0ZRURkeZ7b2NR0pKHOdeclHLJAySxBSIZfV0BzOyvwOlArpnlAz9z94dRX1cRCdn3Zy/l6ffzQ4+blZ7GBz8/p/GCIiKS0iLRGF/8w9zQ42ZlpPFf5x+Wsh8yJXUl6rD7ugb7L2mgvPq6ikgoKmuiPPTmhqQkD989cyxfP3VU6HFFRKR1fbithE//v+SsB/T7Lx/NtENTd1a+pCYQIiLtTVlVhFueW8GTC7c2Xvgg3XHBZC44eigZHXedBxGRTmFdQRmX/+mAveub7ZXrTkupAdP1UQIhIpLg6FtepjoSCz3uf51/GBcfq6mjRUQ6grPu/nfoMbMz0rj9gskpnzyAEggREQAKSiu5/JEFSUkeHrl8Sko3RYuISNPkF+/n4j/Ma7xgM9z9xSM5//BBjRdMAUogRKTTW7ylmOtmL2VjYXnosTfdfn7oMUVEpPVtKizne08uYdveitBjr7plOl0y00OPmyxKIESkU1uxfR+ff+Dd0OMeO6IP9196dOhxRUSkbZx+1xuhx8xMNx65/Nh2lTyAEggR6aTcnbUFZZz327dCj/37Lx/DCaP70SsnM/TYIiLSugrLqjjuF6+GHndAz2zuuugIThmbF3rsZFMCISKd0s+fW8Gj724KPe5PzzuM6ZMGhh5XRERa37qCUm58djnRmIcee/5Pzgo9ZmtRAiEinUos5lw3ewl/X7I91Lij87rxw3PGM31S+xgAJyIiB1ZeFeGsu98MPW6vnExe/N6pocdtTUogRKTT2Fq0n3teXhN68gDwj2tPaXd9WEVEpH5lVREm/ezF0OPecO6hnDtpEAN7dQk9dmvSakYi0iksy9/LKXe+zjOLt4Ua94RR/dh0+/lKHg6CmU03s9Vmts7Mrq/neLaZPRkcn29mIxKO3RDsX21m5wT7hpnZ62a20syWm9l3E8r3NbOXzWxt8L1Pa9yjiLRfi7cUc34SxscN6JnNlaeMYni/rqHHbm1KIESkQ3N3NhWW89n73gk99rNXn8Sfr5gaetyOzMzSgfuBc4EJwCVmNqFOsSuAYncfA9wD3BGcOwGYCUwEpgMPBPEiwPfd/TDgeODqhJjXA6+6+1jg1eC1iEi9luXv5fMPvMvmPftDjTt1ZF/m/+Qs0tMs1LhtRQmEiHRY7s49r6xNytR7L37vVI4Y1pvMdP0aPUhTgXXuvsHdq4FZwIw6ZWYAjwXbTwFnmpkF+2e5e5W7bwTWAVPdfYe7vw/g7qXASmBIPbEeAz6XpPsSkXauojqalA+b/nLlcfz5PzrWh00aAyEiHdYJv3yNnfsqQ4151amj+MrxhzCsb/tvgm4jQ4CtCa/zgeMaKuPuETMrAfoF++fVOXdI4olBd6ejgPnBrgHuviOItcPM+tdXKTO7CrgKYPjw4Qd7TyLSzs1esJUfPb0s9Li/ufhIThyTG3rctqYEQkQ6nFU79/G1Py0IPXn45Rcmc/GUYaR1kCboNlLfm1d3fsSGyhzwXDPrDjwNfM/d9x1Mpdz9QeBBgClTpoQ/X6OIpKw/vrWBW/+xMtSYI3O78YNPjef8wzvmzHxKIESkQ3l87ib++9nlocYcndeNn31mIqeOa3+L/aSgfGBYwuuhQN1psWrL5JtZBtALKDrQuWaWSTx5eMLdn0kos8vMBgWtD4OAgjBvRkTar0g0xr+W7ww9eQB49GvHcki/bqHHTRVKIESkQ9hfHeHul9bwx7c3hhr366eM5JppY7WqdHgWAGPNbCSwjfig6EvrlJkDXAbMBS4EXnN3N7M5wF/M7G5gMDAWeC8YH/EwsNLd724g1u3B92eTc1si0p5UR2L84oWVSVlQdNPt54ceM9UogRCRdm/xlmJuem4FS7fuDTXuLTMm8pUTRoQas7MLxjRcA7wIpAOPuPtyM7sZWOjuc4gnA4+b2TriLQ8zg3OXm9lsYAXxmZeudveomZ0MfAX4wMyWBJf6ibu/QDxxmG1mVwBbgIta725FJFVNufVl9lVGQo157bQxXHHKqFBjpiolECLSbrk7c5Zu57uzljRe+CB8ccpQrj5jTIdufm5LwR/2L9TZd2PCdiUN/KHv7rcBt9XZ9zb1j4/A3fcAZ7awyiLSQWzfW8E3/2dR6MnDfZcexbmTBnWYaVobowRCRNqlXfsq+et7W/jNK2tDjfvN00Zz/bmHhhpTRETa3jvrCvn2E+9TUlETWsz0NOPemUfy6cMHhxazPWgwgTCz14nPblHk7he2XpVERA5s8ZZiPv/Au6HGnHHkYC47cQRHD9dCxU0RdEN6wt2L27ouIiKNeWddIV/64/zGCx6EE0b14+czJjJuQI9Q47YHB2qBuDz4fpmZ9dFDQkTaWnlVhHkb9nDFYwtDjXvpccO57XOTiI/FlSYaCCwws/eBR4AX3V3Tn4pISqmOxHh15S6+9cT7ocYd0juHv151fKgx25MGEwh33wxgZlnoISEibWxdQSnf/J/3WVdQFlrMaYf25+unjOKE0f1Ci9lZuPt/mdl/A58CvgbcFwxwftjd17dt7URE4mY+OJf3t4Q7wcb3zhrL984aF2rM9iatsQLu/l/Ep8p7mHirxFoz+4WZjW7sXDObbmarzWydmV1fz/FsM3syOD4/WEG09tgNwf7VZnZOwv5HzKzAzD6sE6uvmb1sZmuD7+qHINIBVEdiLNhUxFl3vxlq8nD5iSN46KtTlDy0QPBh0s7gKwL0AZ4yszvbtGIi0umVVUUY99N/hp48rLn1XL575thQY7ZHjSYQ0LyHhJmlA/cD5wITgEvMbEKdYlcAxe4+BrgHuCM4dwLxafsmAtOBB4J4AI8G++q6HnjV3ccCrwavRaQdW7+7jCseW8BFv58bWszTxuXxlyuP46bPTuw0s2Ukg5lda2aLgDuBd4DJ7v4t4BjggjatnIh0av9es5tpd71BdTQWWszxA3ow74YzycpIU3dXmjALk5ldS3zxnULgj8AP3b3GzNKAtcCPGjh1KrDO3TcEcWYBM4jP311rBnBTsP0U8SZwC/bPcvcqYGMwF/hUYK67v5nYUlEn1unB9mPAG8CPG7s/EUk9kWiM+RuLQh/w9p9njeOaaWOUOIQjF/hCbXfXWu4eM7NPt1GdRKST+90b67njX6tCjfno147lqOF9tKBogqZM49rch8QQYGvC63zguIbKBIsLlQD9gv3z6pw7pJF6DnD3HUGsHWbWv75CZnYVcBXA8OHDGwkpIq1t1c593PXiGl5ZuSu0mOcfPogrTx7JUZphKTSJ6zbUc2xla9ZFRKSwrIq/zt/Cr19eE2rcP11+LKePr/dPyk6t0QSiBQ+J+j7iqzv4uqEyTTm3Wdz9QeBBgClTpmgwuEiKiERj/H3Jdn7wv0tDjXvfpUdx3qRBpKnVQUSkQ1qWv5evPvIee/eHt77DFSeP5JKpwxjTv/NN0doUyVxILh8YlvB6KLC9gTL5ZpYB9AKKmnhuXbvMbFDQ+jAIKGhJ5UWk9SzcVMTdL6/h3fV7Qov5jVNH8aXjDmF4v66hxRQRkdQRizkfbCthxv3vhBr3hnMP5eunjNIHTweQzARiATDWzEYC24gPir60Tpk5xMdXzAUuBF5zdzezOcBfzOxuYDDxWaDea+R6tbFuD74/G9aNiEhyVEdi3PvqGu5/PdxZP9/60RkM7ZOjgW4iIh3U3v3V3PXSav5n3pbQYh5zSB+uO3scJ43JDS1mR5W0BCIY03AN8CKQDjzi7svN7GZgobtP7v4xAAAgAElEQVTPIT417OPBIOki4kkGQbnZxAdcR4Cr3T0KYGZ/JT5YOtfM8oGfufvDxBOH2WZ2BbAFuChZ9yYiLffS8p3c+o+VbCnaH1rMe2ceyfmTB5GR3qQJ5kREpB3aV1nDkTe/HGrMK08eyQ/OGU+XzPTGC0tSWyBw9xeAF+rsuzFhu5IG/tB399uA2+rZf0kD5fcAZ7akviKSfJU1Ub7x+CL+vWZ3aDG/cNQQ7r74yNDiiYhI6olEY7y8ItxVpXt3zeTmGZP47BGDQ4vZGSQ1gRARqVVZE+WPb23grpfCmyHj26eP5nNHDWHcAA1yay/MbDpwL/GW6T+6++11jmcDfya+nsQe4GJ33xQcu4H4+kFR4Fp3fzHY/wjwaaDA3SclxLoJ+DpQm63+JPhgS0TaGXfnU/e8yYbC8tBifu2kEVw7bSx9umWFFrOzUAIhIklVFYmytaiCs+7+d2gx+3bL4vYvTOZTEweGFlOSL2GB0bOJT5axwMzmuHvi+kAfLTBqZjOJLzB6cZ0FRgcDr5jZuKB766PAfcQTj7rucfe7knZTIpJ0q3buY/pv3go15p0XHM4Xjx3WeEGplxIIEUmauev38NtX1zJ3Q3izKz39rRMYN6AHPbpoQZ92qLUXGBWRdu4bjy/kxeXhrQv0n2eN4+JjhzGwV5fQYnZGSiBEJHS79lXy8NsbefDNDaHFvO7scVxzxhhNq9e+tfYCowDXmNlXgYXA9929uL5CWmRUJHW4O7tLq5j6i1dDjXvXRUdw4TFDQ43ZWSmBEJFQ3fmvVTzwRnjTsl516ii+/6lxZGdoZowOoLUXGP0dcEtQ7hbg18B/1FdQi4yKpIZozLnx2Q95Yn5407P+4vOTOX18HoN754QWs7NTAiEiLVZcXs3cDXv4dogzY3zvrLGcP3kQYzVAuiNp1QVG3f2jfg9m9hDwfLNrLiJJFYs5JRU1HHVLeNOzZqYbv/vSMZw1YUBoMSVOCYSINFs05jzzfj4/fGpZaDFPG5fHPRcfSV/NitERteoCo2Y2yN13BC8/D3wY2p2ISGhK9tdw98ureWzu5tBi3n/p0Rw3qi+53bNDiyn/RwmEiDTLCx/s4DevrGHNrrJQ4n379NFcduIIBvTUwLaOqg0WGL3TzI4k3oVpE/CN1rtbEWlMNOZs3lPOtF+HN0vfiH5d+cNXpjB+oFqvk0kJhIg02f7qCBsLy/ni7+dSXh0NJeaXjx/OzGOHM2lIr1DiSWpr5QVGv9KiyopI0mwsLOcP/17PrAVbGy/cRM9efRLjB/bQatKtQAmEiDTJM+/n86sXV7OjpDKUeBceM5Q7LzgcM4jP0ikiIh1ddSTGSyt2cs1fFocW88JjhnLXRUeEFk8apwRCRBpUHYmxYFMRd7+8hkWb65398qDdddERnDC6H0M0G4aISKfy1trd/OrF1SzLLwkl3qcmDOAXX5iscQ5tQAmEiHxCeVWE9zYW8bVHF4QSr3+PbL4zbQyfOWIwvbtqcLSISGdSWlnDr15czZ9DGiQ9Oq8bP/jUeM6dPCiUeHLwlECIyEfcnQfeWM89L68hEmv5VPg9umTwk/MOY+axw9RNSUSkE/rz3E3c+Ozy0OL97ktHc8ah/TXOoY0pgRAR1uwq5c01u7n1HytDiXfK2Fx+8fnJ5PXI1i95EZFOpiYaI7+4gnPueZPqaCyUmNeeOZZLpg5jUC91f00FSiBEOqlINEZReTW3/mMlc5YecD2uJvvBp8Zx+vj+mlFJRKST+nBbCT+bszy0cXMnju7HX75+fCixJDxKIEQ6mVjM+dvibfz+3+tZW9DyNRyG9+3KFSeP5KsnHKJuSiIindTWov288MEOfvnPVaHE+9bpo5lx5GAOHdgzlHgSLiUQIp3Elj37eXPtbu59dS27S6taHO/IYb2544LDGZnbjayMtBBqKCIi7Y2788t/ruLBNzeEEu+wQT357/MP48QxuaHEk+RQAiHSwS3L38vfFm/jT+9sCiXevTOP5PhR/bRitIhIJ1ZSUcPrqwq4bvYSQphzgx5dMnjh2lMY0LOLPpRqB5RAiHRAReXVvLu+MLSFev7jpJF8+ohBHDWst7opiYh0YpU1Ud5aW8jX/7wwlHhnHtqfr544gtPG5YUST1qHEgiRDiISjbF8+z4eemsDzy/b0eJ4Uw7pw5ePP4TjRvXVrBciIsKTC7bw65fWUBBCN9jc7lncc/GRnDwmVx9MtUNKIETaue17K/jb4m089u6mFv9ST08zrjx5JOdMGsjRw/uEVEMREWmv9lXWsHL7Pq7+y2IKy1qeOJx1WH9u+uxEBvfKIS1NiUN7pQRCpB0qrazhtVUFPLtkO6+tKmhxvNPH5/GV4w/h9PH9SdcvdEkyM5sO3AukA39099vrHM8G/gwcA+wBLnb3TcGxG4ArgChwrbu/GOx/BPg0UODukxJi9QWeBEYAm4Avuns480uKdGDRmDNrwRZ++rcPQ4n3+aOGcPGxwzh+VL9Q4knbSmoCkaSHRL0xzexR4DSgJAh/ubsvSeb9ibSm6kiMhZuKeHnlrlAGRI/K7cYPzhnPGeP7k5Olxd6kdZhZOnA/cDaQDywwsznuviKh2BVAsbuPMbOZwB3AxWY2AZgJTAQGA6+Y2Th3jwKPAvcRf6Ykuh541d1vN7Prg9c/Tt4dirRv+6sjPL9sB3e/tIad+ypbHO/8yYO49XOT6N01U12VOpCkJRDJeEgE5xwo5g/d/alk3ZNIa6uJxli1o5TXVxdw98trWhyvV04md1wwmakj+9G3W1YINRQ5aFOBde6+AcDMZgEzgMRnwwzgpmD7KeA+i//lMQOY5e5VwEYzWxfEm+vub5rZiHquNwM4Pdh+DHgDJRAin7B9bwVz1+/h+/+7NJR491x8BCeNzqW/ZuzrkJLZApGMhwRNiCnSrlVHYmwpKuelFbu481+rWxwvt3sWt8yYxFHD+zCwl36RS5sbAmxNeJ0PHNdQGXePmFkJ0C/YP6/OuUMaud4Ad98RxNphZv1bUHeRDiUWcyojUW6as5zZC/NbHO+4kX352kkjOHFMLj27ZIZQQ0lVyUwgkvWQOFDM28zsRuBV4PogAfkYM7sKuApg+PDhB3lLIslRXhWhoLSKvy3exm9fXdvieIcO7MEPzxnP+IE9GNqnawg1FAlNfX0Y6s4i31CZppzbLHo2SGfzxuoCnn5/G88t3d7iWMeN7MuVp4zi7AkDQqiZtAfJTCCS8ZCob2WR2pg3ADuBLOBB4k3UN3+isPuDwXGmTJkSyoNHpDmKy6spLKvi8Xmb+fPczS2Od/r4PL512miG9u3KkN6adlVSVj4wLOH1UKDuXzC1ZfLNLAPoBRQ18dy6dpnZoKD1YRBQ76wDejZIZ1BcXs3qXaXc/dIa3ttU1OJ4l584guvPPZQumRpH19kkM4FI1kOi3v21TdRAlZn9CfhBCPcgEqodJRXsKKnkwX9v4F/Ld7Y43sxjh/HVE0aQ1yObvB7ZIdRQJOkWAGPNbCSwjfh4t0vrlJkDXAbMBS4EXnN3N7M5wF/M7G7i4+PGAu81cr3aWLcH358N60ZE2ouS/TU8NndTKGPpDunXlTsvOJzR/buT213Pnc4qmQlEMh4S1lDMhE+YDPgcEM68YyIttGF3Gat2lvKHNzewdOveFsXqnp3BNdPG8OnDB9G3WxZdszQTs7QvQXfVa4AXic+m94i7Lzezm4GF7j4HeBh4PBj/VkT8dz1BudnEx71FgKuDGZgws78SHyyda2b5wM/c/WHiicNsM7sC2AJc1Iq3K9Jm3J03Vu/mifmbeWVly6f7/u6ZYzl1XB7HHKI1giSJCUQSHxKfiBlc8gkzyyOeZCwBvpmsexM5kILSSrbvreSN1QU88MZ6qiOxFsU767D+XDJ1OEcN76OZk6RDcPcXgBfq7LsxYbuSBv7Qd/fbgNvq2X9JA+X3AGe2pL4i7cnWov38e81uHnxzA1uK9rco1rEj+nDRlGF8ccqwxgtLp5LUjy+T9JD4RMxg/7SW1lekuRZvKWbb3gqeeX9bixd2G5XbjbMnDuCSY4czoGcXrdEgIiIHVB2JsXBzEXOWbGfWgq2Nn3AAGWnGFaeM5KpTRtGna5ZWi5Z6qf+DyEFyd7YWVfD+lmJW7Szlf+Ztpqwq0ux4OZnpfPn44Uwe2ptzJw0kI8202I6IiDRqY2E5b6wu4OfPtXw2+3MmDuCH54xnSO+u+uBKGqUEQqQRkWiMmqjz3LLtbCws58UPd7KhsLzZ8TLSjC8cPYTjR/Vj6si+mmZVRESarKi8mnfWFXLL8ysoKP3EbPUHZdKQnvzsMxMZlduNfhoQLQdBCYRIHe5OQWkVG3aXM3vhVj7YVsK6grJmx+ualc55kwdx1mEDOKRfVw4b1DPE2oqISEdXWRPljdW7+et7W/j3mt0tijU6rxs/Pf8wxvbvwbC++gBLmkcJhAjwQX4JheVVvLpyF/9YtoPi/TUtiveZIwYzfeJADh/aS7+gRUTkoEWiMV5ZuYvXVhW0eJXoiYN78p1pYxk7oDuj87qHVEPpzJRASKdTE42xYGMRH24vYe2uMuYs3U5VC2ZK6tklgzMPG8ApY3M5Y3x/enTJICO9vjUPRUREGhaLOc9/sIP3Nxfz6LubWhTrhFH9uPykEYzK7cbYAT3CqaBIQAmEdGiVNVHyiyvYWVLJM4vzWburjA+2lbQo5lHDe3Pi6H6MG9CDMw7tT4/sDA16FhGRZnF3/ndhPsu3l/DY3M0tinXOxAHMnDqcQ/p2ZZRaGiSJlEBIh1EViVJWGeG1VQWUVNTw1trCFvcV7dElg8OH9uKzRwxmwqBeDOmTo7UYRESk2WIxp7w6wuyF+SzZupfnlm5vUbwvHTecTx8+mJG53RjYq0tItRQ5MCUQ0i7FYs6aglLKKiO8va6Qd9YVsmL7Psqroy2KO2lITyYP6cWEQT05//DBdMlM02rPIiLSIpFojI2F5Ty3bAdvr93N+1v2NjtW325ZXHXqqI9awrtkaspVaX36y0hS3r7KGmoiMV74cCf5xftZtaOUN9fuxr1lcbMz0pg+aSCnjcujT7csjjmkj7ojiYhIKKoiUV5ZUcC76wt5ftkOSiqaPznHKWNz+dpJIxiZ252Rud1CrKVI8yiBkJRRE41RFYkPcC7eX82y/BLmrt/D6l2lLY7dt1sWZ4zvT99umZw3eRDD+nald06mBjuLiEgo9lXWsK+ihifmb2HxlmLmbShqdqzBvbrwmSMHc9Exw8jtnkXvruo6K6lFCYS0Ondn7/4a1uwqpawqwjvr9rB4azGrd5ayv4VdkACG9slh3IAenD95EMP7dSW3e7Y+sRERkVDVRGPkF1eweuc+Hp+3mXfW7WlRvC8cPYTzJg1i6qi+ag2XlKcEQpLC3TEzthbtp6C0it2lVby5djcrd+xj5Y59VNY0f9rURFNH9OXksbnk9chmTP/uTB7Si6z0NNLS9ItXRETCtSx/L7tLq3hp+S5mL9ra7K60WelpXHzsMCYO7smnJg6kd06mnlvSriiBkBaJxZyaWIxVO0opr46wakcpizYXs66gLJSuR7UOHdiDE0b3IyczndPG5TEyrxtdMtPp2SUztGuISOsws+nAvUA68Ed3v73O8Wzgz8AxwB7gYnffFBy7AbgCiALXuvuLB4ppZo8CpwG18zdf7u5Lknl/0jFEojFW7Sxl4aYi1u8u5+9LtlFaGWlWrD5dMzl9fH9OHN2PYw7pw8jcbmphkHZNCYQ0KhaLf8SybncZ2/dWULy/mgWbilmxfR9rdoXT7ajW5CG9mDSkJ92yMpg+aSD9umfTOyeTPpo6VaRDMLN04H7gbCAfWGBmc9x9RUKxK4Bidx9jZjOBO4CLzWwCMBOYCAwGXjGzccE5B4r5Q3d/Kuk3J+3arn2V7N1fw98Wb2PVzn28sbr504D36ZrJmYcN4OwJAxiZ241xWshNOhglEALEuxztq4hQVh1hwcYiSqsibNhdxgf5JSzbVkJ1C1Zqrs9JY/oxuFcOo/K6c9KYfmRnpDO0Tw7dsvUjKdLBTQXWufsGADObBcwAEhOIGcBNwfZTwH0W/7h2BjDL3auAjWa2LohHE2KKfCQWc97bVERxeTULNxfzzw92sL2kstnxDh3Yg6MP6cOnDx/E6LzuDOip9RikY9Nfa52Eu1MTdbbtraAqEmX5tn2s213GrpJKPtxewvrd5URjLZwXtY7hfbty/Ki+ZKSncdq4PEbldiMtzRilpluRzmwIsDXhdT5wXENl3D1iZiVAv2D/vDrnDgm2DxTzNjO7EXgVuD5IQD7GzK4CrgIYPnz4Qd6SpLKyqgg1kRiLNhfz9rpCVu8sZe6G5g947tstizMP7c/hQ3sxKq87J4zqp/EL0ukogehAKqqj1MRiLN26l8qaGNv3VrB4SzEb9+xn3a7SFi+yVp+MNOMzRwxmWN+uZGekMe3Q/nTLyqBX10x65Wh8goh8Qn1/adX99KKhMg3tr28+5tqYNwA7gSzgQeDHwM2fKOz+YHCcKVOmhPtpirSKaMwxYNm2EjYWlpFfVMErqwpYurX5i7alGYwb0IMzDu3PtEP7a1Y/kYASiHbC3amOxthWXEFN1FlbUMqanaUU7a9m5Y5SNuwuo3h/8xepOZCsjDTOmTiQw4f0IisjjTH9uzNuQA8y0kxjE0TkYOUDwxJeDwW2N1Am38wygF5AUSPn1rvf3XcE+6rM7E/AD0K4B2lD7o47rCkopbi8hrUFpby1tpAP8kvYua/53ZAATh6Ty3Ej+9K7aybTDhtAv25ZWulZpB5KIFJELOaUVkaoqImyeEsxpVURCsuq+HBbCRt2l7OuoIxIyF2MamWkGUcN702vnEzOOLQ/Q3rn4A5HDutNTla6fnmKSJgWAGPNbCSwjfig6EvrlJkDXAbMBS4EXnN3N7M5wF/M7G7ig6jHAu8Rb5moN6aZDXL3HcEYis8BHyb7BiUc7k5lTYz1u8vYXx1l4eYilm7dy5pdZWwsLG9R7BH9ujJ5aG8G9Mjm3MmDyO2exYCeXfS8E2kiJRCtJBKNsbusivziCsoqI6zeFW81yC+uYM2uUgrLqpN6/RH9unLC6Fx6dslgYK8uHDuiL5npafTrnkVu9+ykXltEpFYwpuEa4EXiU64+4u7LzexmYKG7zwEeBh4PBkkXEU8ICMrNJj44OgJc7e5RgPpiBpd8wszyiCcZS4Bvtta9SuOiMSfNYNXOUor3V7O7tIp5G4pYvXMfK0JYM2h4364M7ZPD5CG9OGlMLt2yMxiV243eXTM1Fk+kBZRAhGRPWRVlVREKSqv4IL+EPeVVbCwsZ1txBUvzSxoP0EIDemYzpHcOp47LY0S/blRHYxw9vA+53bPIykija5b+qUUkNbj7C8ALdfbdmLBdCVzUwLm3Abc1JWawf1pL6ystU1pZQySYxGPF9n3sq6xhxY59rN9d3qLxCYmOHdGH7tkZnDI2j0MH9SDdjIlDetE1M10DnEWSQH9VNsOba3bz98XbyC+uYP3uMvaUJ7f1AOLjEKaN78/AXl3o2SWDk8fmkZOZTtfsdM1qJCIibaK8KkLUnVU7StlfHWFPWTXLt+9jS1E5GwrL2bxnfygz/GWlp3HC6H4M7t2FblkZnDY+jz5ds+iala5F2UTaQFITiFZebXQkMAvoC7wPfMXdk/KX/XWzl1JY9olZAFskPc0Y0juHsw4bwNgB3YnEnCOG9mJw7xzSTYOVRURS2fa9FWwp2s+Y/t3pnZOJA5np9U0OlfoqqqOUVsVbDbYW7Wd3WRVF5dVsLCxnZ0klGwvL2VhYTlWI6wONH9CDnKx0JgzuyWGDemLAqNxujBnQHYC87tlKEkRSSNISiDZYbfQO4B53n2Vmvw9i/y7s+6qOxJqdPBw+tBdj+/egW3Y6Z4zvT8+cTDLSjMMG9STNIKOdPmxERDq7215YyT+W7fjYPjM4engfBvfOoToSZWz/Hgzo1YWyyghj+ndnaJ8cSisjdM1KZ3DvHCpqovHZ7bpmEYnFSDM7qEG90ZgTc2d/VZTKSHza7p0llZRVRaiojlJQWkVZVQ17yqsp2FdFdSTGuoIytu+tIC3NKKkIdya/NIPJQ3rRMyeTnMx0jj6kD71zMqmOxpg6si/dsjIwg0G9ckhXNyORdiWZLRCtttqoma0EpvF/M3k8FsQNPYHYWrz/gMfPGJ/H4UN7k55mTBzck3EDemAGQ3rn6NMTEZEOatr4/p9IINxh0eZiFm0uBuDF5buaFCvNIDsjnapIlJ45maSZEY05OZnppKcZVZEo3bIzqKqJUR2NUR2JURWJUhMNf6a+rIw0cjLT6dElg8G9cujTLZNI1Bk3sAd53bPZV1nD+AE9GNQ7h/KqCIN6dSGvRzaxGPTMydBzT6SDSmYC0ZqrjfYD9rp7pJ7yH9PS1Ub3VdQwol9Xzpk4kCOG9aYmGuO4kf3IyUrXwmkiIp3Up48YxPyNez6aUS/NjDSDkooa3KFofzUV1VF6dMlgW3EFVdEYed2z2ba3AoA+XTMpqagh5tC/RxfS04xeOZnBVNppdMlIJyPdSE8zcjIziMZiZGWkxb/S08nOTCM7I42K6ijds//v0d6jSwbp6Wnsr4rQr3s23bPTKa2MMKBnF/p0zWJfZQ29cjLp3zOb8qooXTLTyOueTUVNlOyMdLIy0nB3JQIi8jHJTCBac7XRplwrvrOFq40eNbwPb/zwjIM9TUREOrDsjHTuvPCItq5Gy/RI2EzoUqvkQUTqSman+4NZbZQmrjba0P5CoHcQo6FriYiIiIhICyUzgfhotVEzyyI+KHpOnTK1q41Cwmqjwf6ZZpYdzK5Uu9povTGDc14PYhDEfDaJ9yYiIiIi0iklrQtTG6w2+mNglpndCiwOYouIiIiISIgs/uF95zRlyhRfuHBhW1dDRKTVmdkid5/S1vVIRXo2iEhn1dRngxYeEBERERGRJuvULRBmthvY3MzTc4kP3k4Vqk/jUq1OqVYfSL06qT6Na26dxgJz3X16yPVp91rwbOhIPx/Jovo0LtXqlGr1gdSrU6rVB5pfp0PcPa+xQp06gWgJM1uYSs3/qk/jUq1OqVYfSL06qT6NS8U6dVap+G+RanVSfRqXanVKtfpA6tUp1eoDya+TujCJiIiIiEiTKYEQEREREZEmUwLRfA+2dQXqUH0al2p1SrX6QOrVSfVpXCrWqbNKxX+LVKuT6tO4VKtTqtUHUq9OqVYfSHKdNAZCRERERESaTC0QIiIiIiLSZEogRERERESkyZRA1GFm081stZmtM7Pr6zmebWZPBsfnm9mIhGM3BPtXm9k5rVSf68xshZktM7NXzeyQhGNRM1sSfM0Joz5NrNPlZrY74dpXJhy7zMzWBl+XtVJ97kmoyxoz25twLPT3yMweMbMCM/uwgeNmZr8N6rvMzI5OOBb6+9PEOn0pqMsyM3vXzI5IOLbJzD4I3qNQludtQn1ON7OShH+bGxOOHfDfO0n1+WFCXT4Mfm76BsdCf3+CuMPM7HUzW2lmy83su/WUafWfpc5Kz4ZQ6qRnQwo9G1LtudDEOnXqZ0NKPRfcXV/BF5AOrAdGAVnAUmBCnTLfBn4fbM8Engy2JwTls4GRQZz0VqjPGUDXYPtbtfUJXpe10Xt0OXBfPef2BTYE3/sE232SXZ865b8DPJLk9+hU4GjgwwaOnwf8EzDgeGB+st6fg6jTibXXAs6trVPwehOQ28rv0enA8y399w6rPnXKfgZ4LZnvTxB3EHB0sN0DWFPP/7VW/1nqjF9N/L2nZ4OeDY3VKaWeDU2oT6s+F5pYp9PpxM8GUui5oBaIj5sKrHP3De5eDcwCZtQpMwN4LNh+CjjTzCzYP8vdq9x9I7AuiJfU+rj76+6+P3g5Dxjawmu2uE4HcA7wsrsXuXsx8DLQ0lVwD7Y+lwB/beE1D8jd3wSKDlBkBvBnj5sH9DazQSTn/WlSndz93eCa0Ao/R014jxrSkp+/sOqT9J8hAHff4e7vB9ulwEpgSJ1irf6z1Enp2RBCnQ5Az4a4Vv3/nGrPhabU6QA6xbMhlZ4LSiA+bgiwNeF1Pp/8h/mojLtHgBKgXxPPTUZ9El1BPOus1cXMFprZPDP7XAvrcrB1uiBoOnvKzIYd5LnJqA9BE/5I4LWE3cl4jxrTUJ2T8f40R92fIwdeMrNFZnZVK9bjBDNbamb/NLOJwb42fY/MrCvxX7hPJ+xO+vtj8e4wRwHz6xxK9Z+ljkLPhvDqpGdDw1L5/3OqPBdAz4baa46gDZ8LGc09sYOyevbVnee2oTJNOTcZ9YkXNPsyMAU4LWH3cHffbmajgNfM7AN3X98KdXoO+Ku7V5nZN4l/Kjetiecmoz61ZgJPuXs0YV8y3qPGtObP0EExszOIPyhOTth9UvAe9QdeNrNVwacyyfQ+cIi7l5nZecDfgbG0/Xv0GeAdd0/8RCqp74+ZdSf+UPqeu++re7ieU1LiZ6mD0bMhnDrp2XBgKfn/OYWeC6BnA5AazwW1QHxcPjAs4fVQYHtDZcwsA+hFvHmrKecmoz6Y2VnAT4HPuntV7X533x583wC8QTxTbalG6+TuexLq8RBwTFPPTUZ9EsykTvNikt6jxjRU52S8P01mZocDfwRmuPue2v0J71EB8Dda3v2iUe6+z93Lgu0XgEwzy6WN3yMO/DMU+vtjZpnEHxJPuPsz9RRJyZ+lDkjPhhDqpGdDo1Lu/3MqPReC63X6Z0PKPBc85AEw7fmLeIvMBuJNmbWDcCbWKXM1Hx8oNzvYnsjHB8ptoOUD5ZpSn6OIDxwaW2d/HyA72M4F1hLOgKKm1GlQwvbngXnBdl9gY1C3PsF232TXJyg3nviAJkv2exTEG0HDg8DO5+MDnN5L1vtzEHUaTrxv9ol19ncDevz/9u6YNYogDAPwO50gQkpLiY2FpZX4A9RfoZoeG+gAAAIjSURBVKYJ4j+4RmwM2KSxs7K1E0vBShGxiEYLxcpGLMTC0uIsdg4WkWSK3M6d9zywkOzdJB8zu/vmu2XJ6OtXSa5OUM/ZxVpluOh+rfPVtN4nXU99ffEH4emJ5qckeZxk/4j3dDmWNm1rvO7JBtnQUtdR173Jz+dj6pk8Fxpq2uhsyArlwoks9v+0ZXh6/XOGC++s7ruX4ROcJDmV5Ek9qd4k2R6NndVxn5Jcm6ie50m+Jzmo29O6/3KSw3oSHSbZmXCO7if5WH/3iyQXRmNv1bn7kuTmFPXU7+8m2ftr3FLmKMOnEN+S/M7Q8e8k2U2yW18vSR7Weg+TXFrm/DTW9CjJz9Fx9Lbu367z866u6Wyieu6MjqHXGQXYv9Z72fXU99zI8DDseNxS5qf+7CsZbi+/H63L9d7H0qZux11nIhtaapINK5QNDfVMmguNNW10NmSFcmHRxQEAABzLMxAAAEAzDQQAANBMAwEAADTTQAAAAM00EAAAQDMNBHRQStkqpdzuXQcAq0EusE40ENDHVhJBAcCCXGBtaCCgj70k50spB6WUB72LAaA7ucDa8I/koINSyrkkz+bz+cXOpQCwAuQC68QdCAAAoJkGAgAAaKaBgD5+JTnTuwgAVoZcYG1oIKCD+Xz+I8nLUsoHD8sBIBdYJx6iBgAAmrkDAQAANNNAAAAAzTQQAABAMw0EAADQTAMBAAA000AAAADNNBAAAECzPw6kVO/IoCigAAAAAElFTkSuQmCC\n",
"text/plain": [
"