{ "cells": [ { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "# Modelos Metapopulacionais\n", "\n", "O conceito de metapopulações pode ser visto como uma maneira relativamente simples de tratar a heterogeneidade em populações sem abandonar o paradigma de populações homogêneas modeladas com Equações diferenciais ordinárias.\n", "\n", "![metapops](https://upload.wikimedia.org/wikipedia/commons/thumb/9/94/Metapopulation_%281%29.svg/2560px-Metapopulation_%281%29.svg.png)\n", "\n", "#### Exercício 1:\n", "\n", "Desenvolva um modelo metapopulacional com dois patches (duas sub-populações) A e B nos quais se propaga uma doença de transmissão direta modelada como SIR. Considere que existe um fluxo migratório entre as cidades.\n", "\n", "1. Escreva as equações das duas populações representando as migrações $A \\rightarrow B$ como $k_1$, e $B \\rightarrow A$ como $k_2$.\n", "1. Utilizando uma parametrização que não permita equilíbrio endêmico, encontre os equilíbrios $S_A^* + S_B^*$ e $I_A^* + I_B^*$ .\n", "1. Calcule a Matriz Jacobiana do sistema e determine a estabilidade dos equilíbrios.\n", "1. Crie uma função interativa que nos permita alterar os valores de $k_1$ e $k_2$. Realize simulações com $k_1 = k_2$, $k_1 < k_2$ e $k_1 > k_2$, com $\\beta_A \\neq \\beta_B$" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:06:36.606257Z", "start_time": "2021-10-27T12:06:36.597247Z" }, "collapsed": false }, "outputs": [ ], "source": [ "%display typeset" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "#### Resposta ao item 1:\n", "\n", "$$ \\dot{S}_A = -\\beta_A S_A I_A + k_2 S_B - k_1 S_A$$\n", "$$\\dot{I}_A = \\beta_A S_A I_A + k_2 I_B - k_1 I_A -\\gamma I_A$$\n", "$$\\dot{S}_B = -\\beta_B S_B I_B + k_1 S_A - k_2 S_B $$\n", "$$\\dot{I}_B = \\beta_B S_B I_B + k_1 I_A - k_2 I_B -\\gamma I_B$$" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "2) Vamos agora encontrar os equilíbrios do sistema usando o solver do sage:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:46:57.626884Z", "start_time": "2021-10-27T12:46:57.583435Z" }, "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\left[\\left[S_{A} = r_{1}, I_{A} = 0, S_{B} = \\frac{k_{1} r_{1}}{k_{2}}, I_{B} = 0, R_{A} = r_{2}, R_{B} = r_{3}\\right]\\right]\\)" ], "text/latex": [ "$\\displaystyle \\left[\\left[S_{A} = r_{1}, I_{A} = 0, S_{B} = \\frac{k_{1} r_{1}}{k_{2}}, I_{B} = 0, R_{A} = r_{2}, R_{B} = r_{3}\\right]\\right]$" ], "text/plain": [ "[[S_A == r1, I_A == 0, S_B == k_1*r1/k_2, I_B == 0, R_A == r2, R_B == r3]]" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" } ], "source": [ "var('S_A I_A S_B I_B beta_A beta_B R_A R_B gamma k_1 k_2')\n", "#beta_A = 1.2\n", "#beta_B = 1.2\n", "#gamma = 0.5\n", "#k_1 = 0.5\n", "#k_2 = 0.5\n", "solve([-beta_A*S_A*I_A + k_2*S_B - k_1*S_A, \n", " beta_A*S_A*I_A+ k_2*I_B - k_1*I_A - gamma*I_A,\n", " -beta_B*S_B*I_B - k_2*S_B + k_1*S_A,\n", " beta_B*S_B*I_B- k_2*I_B + k_1*I_A - gamma*I_B,\n", " gamma*I_A,\n", " gamma*I_B],[S_A, I_A, S_B, I_B, R_A, R_B])" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "3) Calculamos a Matriz jacobiana para podermos realizar a análise de estabilidade do equilíbrio encontrado acima" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:50:25.456043Z", "start_time": "2021-10-27T12:50:25.444750Z" }, "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\left(\\begin{array}{rrrr}\n", "-I_{A} \\beta_{A} - k_{1} & -S_{A} \\beta_{A} & k_{2} & 0 \\\\\n", "I_{A} \\beta_{A} & S_{A} \\beta_{A} - \\gamma - k_{1} & 0 & k_{2} \\\\\n", "k_{1} & 0 & -I_{B} \\beta_{B} - k_{2} & -S_{B} \\beta_{B} \\\\\n", "0 & k_{1} & I_{B} \\beta_{B} & S_{B} \\beta_{B} - \\gamma - k_{2}\n", "\\end{array}\\right)\\)" ], "text/latex": [ "$\\displaystyle \\left(\\begin{array}{rrrr}\n", "-I_{A} \\beta_{A} - k_{1} & -S_{A} \\beta_{A} & k_{2} & 0 \\\\\n", "I_{A} \\beta_{A} & S_{A} \\beta_{A} - \\gamma - k_{1} & 0 & k_{2} \\\\\n", "k_{1} & 0 & -I_{B} \\beta_{B} - k_{2} & -S_{B} \\beta_{B} \\\\\n", "0 & k_{1} & I_{B} \\beta_{B} & S_{B} \\beta_{B} - \\gamma - k_{2}\n", "\\end{array}\\right)$" ], "text/plain": [ "[ -I_A*beta_A - k_1 -S_A*beta_A k_2 0]\n", "[ I_A*beta_A S_A*beta_A - gamma - k_1 0 k_2]\n", "[ k_1 0 -I_B*beta_B - k_2 -S_B*beta_B]\n", "[ 0 k_1 I_B*beta_B S_B*beta_B - gamma - k_2]" ] }, "execution_count": 3, "metadata": { }, "output_type": "execute_result" } ], "source": [ "J = jacobian([-beta_A*S_A*I_A + k_2*S_B - k_1*S_A, \n", " beta_A*S_A*I_A+ k_2*I_B - k_1*I_A - gamma*I_A,\n", " -beta_B*S_B*I_B - k_2*S_B + k_1*S_A,\n", " beta_B*S_B*I_B- k_2*I_B + k_1*I_A - gamma*I_B],[S_A, I_A, S_B, I_B])\n", "J" ] }, { "cell_type": "markdown", "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:08:00.845976Z", "start_time": "2021-10-27T12:08:00.841739Z" }, "collapsed": false }, "source": [ "Agora Substituímos o valor do equilíbrio na matriz jacobiana:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:51:04.898288Z", "start_time": "2021-10-27T12:51:04.888269Z" }, "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\left(\\begin{array}{rrrr}\n", "-k_{1} & -\\beta_{A} & k_{2} & 0 \\\\\n", "0 & \\beta_{A} - \\gamma - k_{1} & 0 & k_{2} \\\\\n", "k_{1} & 0 & -k_{2} & -\\frac{\\beta_{B} k_{1}}{k_{2}} \\\\\n", "0 & k_{1} & 0 & -\\gamma + \\frac{\\beta_{B} k_{1}}{k_{2}} - k_{2}\n", "\\end{array}\\right)\\)" ], "text/latex": [ "$\\displaystyle \\left(\\begin{array}{rrrr}\n", "-k_{1} & -\\beta_{A} & k_{2} & 0 \\\\\n", "0 & \\beta_{A} - \\gamma - k_{1} & 0 & k_{2} \\\\\n", "k_{1} & 0 & -k_{2} & -\\frac{\\beta_{B} k_{1}}{k_{2}} \\\\\n", "0 & k_{1} & 0 & -\\gamma + \\frac{\\beta_{B} k_{1}}{k_{2}} - k_{2}\n", "\\end{array}\\right)$" ], "text/plain": [ "[ -k_1 -beta_A k_2 0]\n", "[ 0 beta_A - gamma - k_1 0 k_2]\n", "[ k_1 0 -k_2 -beta_B*k_1/k_2]\n", "[ 0 k_1 0 -gamma + beta_B*k_1/k_2 - k_2]" ] }, "execution_count": 4, "metadata": { }, "output_type": "execute_result" } ], "source": [ "J_eq = J(S_A=1,I_A=0,S_B=k_1/k_2,I_B=0)\n", "J_eq" ] }, { "cell_type": "markdown", "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:09:00.121973Z", "start_time": "2021-10-27T12:09:00.117775Z" }, "collapsed": false }, "source": [ "e também os valores dos parâmetros" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:53:36.962542Z", "start_time": "2021-10-27T12:53:36.947112Z" }, "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\left[\\frac{3}{5}, -\\frac{2}{5}, -1, 0\\right]\\)" ], "text/latex": [ "$\\displaystyle \\left[\\frac{3}{5}, -\\frac{2}{5}, -1, 0\\right]$" ], "text/plain": [ "[3/5, -2/5, -1, 0]" ] }, "execution_count": 5, "metadata": { }, "output_type": "execute_result" } ], "source": [ "J_eq(gamma=1.0,k_1=.5,k_2=.5,beta_A=1.6,beta_B=1.6).eigenvalues()" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Como o maior autovalor é real e positivo o Equilíbrio é instável" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "4) Agora vamos criar uma simulação interativa para estudar o efeito da existência de patches na dinâmica do sistema." ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:09:45.258939Z", "start_time": "2021-10-27T12:09:45.242162Z" }, "collapsed": false }, "outputs": [ ], "source": [ "def model_2_pop(t, y, params):\n", " beta_A, beta_B, gamma, k_1, k_2, mu = params\n", " S_A, I_A, S_B, I_B, R_A, R_B = y\n", " return [-beta_A*S_A*I_A + k_2*S_B - k_1*S_A, \n", " beta_A*S_A*I_A+ k_2*I_B - k_1*I_A - gamma*I_A,\n", " -beta_B*S_B*I_B - k_2*S_B + k_1*S_A,\n", " beta_B*S_B*I_B - k_2*I_B + k_1*I_A - gamma*I_B,\n", " gamma*I_A,\n", " gamma*I_B]" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Agora vamos criar uma versão do modelo anterior, com nascimentos e mortes, ou seja demografia:\n", "\n", "$$ \\dot{S}_A = -\\beta_A S_A I_A + k_2 S_B - k_1 S_A -\\mu S_A +\\mu(S_A+I_A+R_A)$$\n", "$$\\dot{I}_A = \\beta_A S_A I_A + k_2 I_B - k_1 I_A -\\gamma I_A -\\mu I_A$$\n", "$$\\dot{R}_A = \\gamma I_A-\\mu R_A$$\n", "$$\\dot{S}_B = -\\beta_B S_B I_B + k_1 S_A - k_2 S_B -\\mu S_B +\\mu(S_B+I_B+R_B)$$\n", "$$\\dot{I}_B = \\beta_B S_B I_B + k_1 I_A - k_2 I_B -\\gamma I_B -\\mu I_B$$\n", "$$\\dot{R}_B= \\gamma I_B-\\mu R_B$$" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T13:02:05.137728Z", "start_time": "2021-10-27T13:02:05.124770Z" }, "collapsed": false }, "outputs": [ ], "source": [ "def model_2_pop_demog(t, y, params):\n", " beta_A, beta_B, gamma, k_1, k_2, mu = params\n", " S_A, I_A, S_B, I_B, R_A, R_B = y\n", " return [-beta_A*S_A*I_A + k_2*S_B - k_1*S_A -mu*S_A +mu*(S_A+I_A+R_A), \n", " beta_A*S_A*I_A+ k_2*I_B - k_1*I_A - gamma*I_A - mu*I_A,\n", " -beta_B*S_B*I_B - k_2*S_B + k_1*S_A - mu*S_B+mu*(S_B+I_B+R_B),\n", " beta_B*S_B*I_B - k_2*I_B + k_1*I_A - gamma*I_B - mu*I_B,\n", " gamma*I_A-mu*R_A,\n", " gamma*I_B-mu*R_B]" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Agora passamos à instanciação do nosso solver e à implementação de uma função de plotagem." ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:10:00.599063Z", "start_time": "2021-10-27T12:10:00.587040Z" }, "collapsed": false }, "outputs": [ ], "source": [ "T=ode_solver()\n", "T.function = model_2_pop\n", "T.algorithm = \"rk8pd\"#\"rk4imp\"\n", "from itertools import cycle\n", "labels = [r\"$S_A$\", \"$I_A$\", \"$S_B$\", \"$I_B$\",\"$R_A$\", \"$R_B$\"]\n", "\n", "def plot_sol(sol):\n", " plots=None\n", " c = cycle(['red','blue','green', 'black', 'yellow', 'orange', 'magenta', 'gray', 'pink', 'brown'])\n", " for i,co in zip(range(len(sol[0][1])),c):\n", "# co = c.next()\n", " if plots is None:\n", " plots = list_plot([(j[0],j[1][i]) for j in sol],color=co, plotjoined=True, legend_label='%s'%labels[i], alpha=.8, gridlines=true)\n", " else:\n", " plots += list_plot([(j[0],j[1][i]) for j in sol],color=co, plotjoined=True, legend_label='%s'%labels[i], alpha=.8, gridlines=true)\n", " show(plots)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T13:09:10.707213Z", "start_time": "2021-10-27T13:09:09.635603Z" }, "collapsed": false }, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "5d89c8df33e04dcfb38141bf261b31c7", "version_major": 2, "version_minor": 0 }, "text/plain": [ "Interactive function with 7 widgets\n", " beta_A: FloatSlider(value=4.6, descr…" ] }, "execution_count": 9, "metadata": { }, "output_type": "execute_result" } ], "source": [ "@interact()\n", "def simula(beta_A=(4.6,(.2,5,.1)), beta_B=(3.6,(.2,5,.1)), gamma=(1.2, 2,.1), k_1=(0.5,(0,1,.1)), k_2=(0.5,(0,1,.1)), mu=(.04,(0, 1,.1)), model=selector(['no demog', 'demog'], buttons=True)):\n", " if model == \"demog\":\n", " T.function = model_2_pop_demog\n", " else:\n", " T.function = model_2_pop\n", " inits = (.9,0.1, 1, 0, 0, 0)\n", " t_span = [0,30]\n", " T.ode_solve(t_span,inits,num_points=300, params=(beta_A, beta_B, gamma, k_1, k_2, mu))\n", " show(inits[0]*k_1/k_2)\n", " plot_sol(T.solution)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "#### Exercício 2:\n", "\n", "Modelos evolutivos em 5 patches. Adapte o modelo de seleção natural sem mutação, para 4 subpopulações, mantendo a generalidade para multiplos tipos.\n", "\n", "1. Escreva as equações em forma matricial\n", "2. Encontre os equilíbrios do sistema\n", "3. Escreva uma função interativa para estudar os parâmetros\n", "\n" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T13:16:27.024486Z", "start_time": "2021-10-27T13:16:26.971079Z" }, "collapsed": false }, "outputs": [ { "data": { "text/html": [ "a) Vamos escrever o modelo evolutivo, começando com 4 tipos e 5 localidades, denotadas pelos vetores \\(X\\) e \\(L\\), respectivamente:" ], "text/plain": [ "a) Vamos escrever o modelo evolutivo, começando com 4 tipos e 5 localidades, denotadas pelos vetores \\(X\\) e \\(L\\), respectivamente:" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "\\(\\displaystyle \\verb|X=| \\left(\\begin{array}{rrrr}\n", "x_{1} & x_{2} & x_{3} & x_{4}\n", "\\end{array}\\right) \\verb|,|\\verb| |\\verb|L=| \\left(\\begin{array}{rrrrr}\n", "A & B & C & D & E\n", "\\end{array}\\right)\\)" ], "text/latex": [ "$\\displaystyle \\verb|X=| \\left(\\begin{array}{rrrr}\n", "x_{1} & x_{2} & x_{3} & x_{4}\n", "\\end{array}\\right) \\verb|,|\\verb| |\\verb|L=| \\left(\\begin{array}{rrrrr}\n", "A & B & C & D & E\n", "\\end{array}\\right)$" ], "text/plain": [ "'X=' [x_1 x_2 x_3 x_4] ', L=' [A B C D E]" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "Seja \\(D\\) a matriz com as densidades em cada uma das localidades:" ], "text/plain": [ "Seja \\(D\\) a matriz com as densidades em cada uma das localidades:" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "\\(\\displaystyle \\verb|D=| \\left(\\begin{array}{rrrrr}\n", "x_{\\mathit{1A}} & x_{\\mathit{1B}} & x_{\\mathit{1C}} & x_{\\mathit{1D}} & x_{\\mathit{1E}} \\\\\n", "x_{\\mathit{2A}} & x_{\\mathit{2B}} & x_{\\mathit{2C}} & x_{\\mathit{2D}} & x_{\\mathit{2E}} \\\\\n", "x_{\\mathit{3A}} & x_{\\mathit{3B}} & x_{\\mathit{3C}} & x_{\\mathit{3D}} & x_{\\mathit{3E}} \\\\\n", "x_{\\mathit{4A}} & x_{\\mathit{4B}} & x_{\\mathit{4C}} & x_{\\mathit{4D}} & x_{\\mathit{4E}}\n", "\\end{array}\\right)\\)" ], "text/latex": [ "$\\displaystyle \\verb|D=| \\left(\\begin{array}{rrrrr}\n", "x_{\\mathit{1A}} & x_{\\mathit{1B}} & x_{\\mathit{1C}} & x_{\\mathit{1D}} & x_{\\mathit{1E}} \\\\\n", "x_{\\mathit{2A}} & x_{\\mathit{2B}} & x_{\\mathit{2C}} & x_{\\mathit{2D}} & x_{\\mathit{2E}} \\\\\n", "x_{\\mathit{3A}} & x_{\\mathit{3B}} & x_{\\mathit{3C}} & x_{\\mathit{3D}} & x_{\\mathit{3E}} \\\\\n", "x_{\\mathit{4A}} & x_{\\mathit{4B}} & x_{\\mathit{4C}} & x_{\\mathit{4D}} & x_{\\mathit{4E}}\n", "\\end{array}\\right)$" ], "text/plain": [ "'D=' [x_1A x_1B x_1C x_1D x_1E]\n", "[x_2A x_2B x_2C x_2D x_2E]\n", "[x_3A x_3B x_3C x_3D x_3E]\n", "[x_4A x_4B x_4C x_4D x_4E]" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "Seja \\(M\\) a matriz com as taxas de migração \\(k_{ij}\\) entre cada uma das localidades. como não faz sentido fala de migração de uma localidade " ], "text/plain": [ "Seja \\(M\\) a matriz com as taxas de migração \\(k_{ij}\\) entre cada uma das localidades. como não faz sentido fala de migração de uma localidade " ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "para si mesma, sua diagonal é zero. " ], "text/plain": [ "para si mesma, sua diagonal é zero. " ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "\\(\\displaystyle \\verb|M=| \\left(\\begin{array}{rrrrr}\n", "0 & k_{\\mathit{AB}} & k_{\\mathit{AC}} & k_{\\mathit{AD}} & k_{\\mathit{AE}} \\\\\n", "k_{\\mathit{BA}} & 0 & k_{\\mathit{BC}} & k_{\\mathit{BD}} & k_{\\mathit{BE}} \\\\\n", "k_{\\mathit{CA}} & k_{\\mathit{CB}} & 0 & k_{\\mathit{CD}} & k_{\\mathit{CE}} \\\\\n", "k_{\\mathit{DA}} & k_{\\mathit{DB}} & k_{\\mathit{DC}} & 0 & k_{\\mathit{DE}} \\\\\n", "k_{\\mathit{EA}} & k_{\\mathit{EB}} & k_{\\mathit{EC}} & k_{\\mathit{ED}} & 0\n", "\\end{array}\\right)\\)" ], "text/latex": [ "$\\displaystyle \\verb|M=| \\left(\\begin{array}{rrrrr}\n", "0 & k_{\\mathit{AB}} & k_{\\mathit{AC}} & k_{\\mathit{AD}} & k_{\\mathit{AE}} \\\\\n", "k_{\\mathit{BA}} & 0 & k_{\\mathit{BC}} & k_{\\mathit{BD}} & k_{\\mathit{BE}} \\\\\n", "k_{\\mathit{CA}} & k_{\\mathit{CB}} & 0 & k_{\\mathit{CD}} & k_{\\mathit{CE}} \\\\\n", "k_{\\mathit{DA}} & k_{\\mathit{DB}} & k_{\\mathit{DC}} & 0 & k_{\\mathit{DE}} \\\\\n", "k_{\\mathit{EA}} & k_{\\mathit{EB}} & k_{\\mathit{EC}} & k_{\\mathit{ED}} & 0\n", "\\end{array}\\right)$" ], "text/plain": [ "'M=' [ 0 k_AB k_AC k_AD k_AE]\n", "[k_BA 0 k_BC k_BD k_BE]\n", "[k_CA k_CB 0 k_CD k_CE]\n", "[k_DA k_DB k_DC 0 k_DE]\n", "[k_EA k_EB k_EC k_ED 0]" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "Logo, \\(F=D \\times M\\) é uma matriz em que cada elemento \\(f_{ij}\\) representa a IMIGRAÇÃO do tipo \\(i\\) em cada patch \\(j\\) vindo de todos os outros patches:" ], "text/plain": [ "Logo, \\(F=D \\times M\\) é uma matriz em que cada elemento \\(f_{ij}\\) representa a IMIGRAÇÃO do tipo \\(i\\) em cada patch \\(j\\) vindo de todos os outros patches:" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "\\(\\displaystyle \\verb|Fi=| \\left(\\begin{array}{rrrrr}\n", "k_{\\mathit{BA}} x_{\\mathit{1B}} + k_{\\mathit{CA}} x_{\\mathit{1C}} + k_{\\mathit{DA}} x_{\\mathit{1D}} + k_{\\mathit{EA}} x_{\\mathit{1E}} & k_{\\mathit{AB}} x_{\\mathit{1A}} + k_{\\mathit{CB}} x_{\\mathit{1C}} + k_{\\mathit{DB}} x_{\\mathit{1D}} + k_{\\mathit{EB}} x_{\\mathit{1E}} & k_{\\mathit{AC}} x_{\\mathit{1A}} + k_{\\mathit{BC}} x_{\\mathit{1B}} + k_{\\mathit{DC}} x_{\\mathit{1D}} + k_{\\mathit{EC}} x_{\\mathit{1E}} & k_{\\mathit{AD}} x_{\\mathit{1A}} + k_{\\mathit{BD}} x_{\\mathit{1B}} + k_{\\mathit{CD}} x_{\\mathit{1C}} + k_{\\mathit{ED}} x_{\\mathit{1E}} & k_{\\mathit{AE}} x_{\\mathit{1A}} + k_{\\mathit{BE}} x_{\\mathit{1B}} + k_{\\mathit{CE}} x_{\\mathit{1C}} + k_{\\mathit{DE}} x_{\\mathit{1D}} \\\\\n", "k_{\\mathit{BA}} x_{\\mathit{2B}} + k_{\\mathit{CA}} x_{\\mathit{2C}} + k_{\\mathit{DA}} x_{\\mathit{2D}} + k_{\\mathit{EA}} x_{\\mathit{2E}} & k_{\\mathit{AB}} x_{\\mathit{2A}} + k_{\\mathit{CB}} x_{\\mathit{2C}} + k_{\\mathit{DB}} x_{\\mathit{2D}} + k_{\\mathit{EB}} x_{\\mathit{2E}} & k_{\\mathit{AC}} x_{\\mathit{2A}} + k_{\\mathit{BC}} x_{\\mathit{2B}} + k_{\\mathit{DC}} x_{\\mathit{2D}} + k_{\\mathit{EC}} x_{\\mathit{2E}} & k_{\\mathit{AD}} x_{\\mathit{2A}} + k_{\\mathit{BD}} x_{\\mathit{2B}} + k_{\\mathit{CD}} x_{\\mathit{2C}} + k_{\\mathit{ED}} x_{\\mathit{2E}} & k_{\\mathit{AE}} x_{\\mathit{2A}} + k_{\\mathit{BE}} x_{\\mathit{2B}} + k_{\\mathit{CE}} x_{\\mathit{2C}} + k_{\\mathit{DE}} x_{\\mathit{2D}} \\\\\n", "k_{\\mathit{BA}} x_{\\mathit{3B}} + k_{\\mathit{CA}} x_{\\mathit{3C}} + k_{\\mathit{DA}} x_{\\mathit{3D}} + k_{\\mathit{EA}} x_{\\mathit{3E}} & k_{\\mathit{AB}} x_{\\mathit{3A}} + k_{\\mathit{CB}} x_{\\mathit{3C}} + k_{\\mathit{DB}} x_{\\mathit{3D}} + k_{\\mathit{EB}} x_{\\mathit{3E}} & k_{\\mathit{AC}} x_{\\mathit{3A}} + k_{\\mathit{BC}} x_{\\mathit{3B}} + k_{\\mathit{DC}} x_{\\mathit{3D}} + k_{\\mathit{EC}} x_{\\mathit{3E}} & k_{\\mathit{AD}} x_{\\mathit{3A}} + k_{\\mathit{BD}} x_{\\mathit{3B}} + k_{\\mathit{CD}} x_{\\mathit{3C}} + k_{\\mathit{ED}} x_{\\mathit{3E}} & k_{\\mathit{AE}} x_{\\mathit{3A}} + k_{\\mathit{BE}} x_{\\mathit{3B}} + k_{\\mathit{CE}} x_{\\mathit{3C}} + k_{\\mathit{DE}} x_{\\mathit{3D}} \\\\\n", "k_{\\mathit{BA}} x_{\\mathit{4B}} + k_{\\mathit{CA}} x_{\\mathit{4C}} + k_{\\mathit{DA}} x_{\\mathit{4D}} + k_{\\mathit{EA}} x_{\\mathit{4E}} & k_{\\mathit{AB}} x_{\\mathit{4A}} + k_{\\mathit{CB}} x_{\\mathit{4C}} + k_{\\mathit{DB}} x_{\\mathit{4D}} + k_{\\mathit{EB}} x_{\\mathit{4E}} & k_{\\mathit{AC}} x_{\\mathit{4A}} + k_{\\mathit{BC}} x_{\\mathit{4B}} + k_{\\mathit{DC}} x_{\\mathit{4D}} + k_{\\mathit{EC}} x_{\\mathit{4E}} & k_{\\mathit{AD}} x_{\\mathit{4A}} + k_{\\mathit{BD}} x_{\\mathit{4B}} + k_{\\mathit{CD}} x_{\\mathit{4C}} + k_{\\mathit{ED}} x_{\\mathit{4E}} & k_{\\mathit{AE}} x_{\\mathit{4A}} + k_{\\mathit{BE}} x_{\\mathit{4B}} + k_{\\mathit{CE}} x_{\\mathit{4C}} + k_{\\mathit{DE}} x_{\\mathit{4D}}\n", "\\end{array}\\right)\\)" ], "text/latex": [ "$\\displaystyle \\verb|Fi=| \\left(\\begin{array}{rrrrr}\n", "k_{\\mathit{BA}} x_{\\mathit{1B}} + k_{\\mathit{CA}} x_{\\mathit{1C}} + k_{\\mathit{DA}} x_{\\mathit{1D}} + k_{\\mathit{EA}} x_{\\mathit{1E}} & k_{\\mathit{AB}} x_{\\mathit{1A}} + k_{\\mathit{CB}} x_{\\mathit{1C}} + k_{\\mathit{DB}} x_{\\mathit{1D}} + k_{\\mathit{EB}} x_{\\mathit{1E}} & k_{\\mathit{AC}} x_{\\mathit{1A}} + k_{\\mathit{BC}} x_{\\mathit{1B}} + k_{\\mathit{DC}} x_{\\mathit{1D}} + k_{\\mathit{EC}} x_{\\mathit{1E}} & k_{\\mathit{AD}} x_{\\mathit{1A}} + k_{\\mathit{BD}} x_{\\mathit{1B}} + k_{\\mathit{CD}} x_{\\mathit{1C}} + k_{\\mathit{ED}} x_{\\mathit{1E}} & k_{\\mathit{AE}} x_{\\mathit{1A}} + k_{\\mathit{BE}} x_{\\mathit{1B}} + k_{\\mathit{CE}} x_{\\mathit{1C}} + k_{\\mathit{DE}} x_{\\mathit{1D}} \\\\\n", "k_{\\mathit{BA}} x_{\\mathit{2B}} + k_{\\mathit{CA}} x_{\\mathit{2C}} + k_{\\mathit{DA}} x_{\\mathit{2D}} + k_{\\mathit{EA}} x_{\\mathit{2E}} & k_{\\mathit{AB}} x_{\\mathit{2A}} + k_{\\mathit{CB}} x_{\\mathit{2C}} + k_{\\mathit{DB}} x_{\\mathit{2D}} + k_{\\mathit{EB}} x_{\\mathit{2E}} & k_{\\mathit{AC}} x_{\\mathit{2A}} + k_{\\mathit{BC}} x_{\\mathit{2B}} + k_{\\mathit{DC}} x_{\\mathit{2D}} + k_{\\mathit{EC}} x_{\\mathit{2E}} & k_{\\mathit{AD}} x_{\\mathit{2A}} + k_{\\mathit{BD}} x_{\\mathit{2B}} + k_{\\mathit{CD}} x_{\\mathit{2C}} + k_{\\mathit{ED}} x_{\\mathit{2E}} & k_{\\mathit{AE}} x_{\\mathit{2A}} + k_{\\mathit{BE}} x_{\\mathit{2B}} + k_{\\mathit{CE}} x_{\\mathit{2C}} + k_{\\mathit{DE}} x_{\\mathit{2D}} \\\\\n", "k_{\\mathit{BA}} x_{\\mathit{3B}} + k_{\\mathit{CA}} x_{\\mathit{3C}} + k_{\\mathit{DA}} x_{\\mathit{3D}} + k_{\\mathit{EA}} x_{\\mathit{3E}} & k_{\\mathit{AB}} x_{\\mathit{3A}} + k_{\\mathit{CB}} x_{\\mathit{3C}} + k_{\\mathit{DB}} x_{\\mathit{3D}} + k_{\\mathit{EB}} x_{\\mathit{3E}} & k_{\\mathit{AC}} x_{\\mathit{3A}} + k_{\\mathit{BC}} x_{\\mathit{3B}} + k_{\\mathit{DC}} x_{\\mathit{3D}} + k_{\\mathit{EC}} x_{\\mathit{3E}} & k_{\\mathit{AD}} x_{\\mathit{3A}} + k_{\\mathit{BD}} x_{\\mathit{3B}} + k_{\\mathit{CD}} x_{\\mathit{3C}} + k_{\\mathit{ED}} x_{\\mathit{3E}} & k_{\\mathit{AE}} x_{\\mathit{3A}} + k_{\\mathit{BE}} x_{\\mathit{3B}} + k_{\\mathit{CE}} x_{\\mathit{3C}} + k_{\\mathit{DE}} x_{\\mathit{3D}} \\\\\n", "k_{\\mathit{BA}} x_{\\mathit{4B}} + k_{\\mathit{CA}} x_{\\mathit{4C}} + k_{\\mathit{DA}} x_{\\mathit{4D}} + k_{\\mathit{EA}} x_{\\mathit{4E}} & k_{\\mathit{AB}} x_{\\mathit{4A}} + k_{\\mathit{CB}} x_{\\mathit{4C}} + k_{\\mathit{DB}} x_{\\mathit{4D}} + k_{\\mathit{EB}} x_{\\mathit{4E}} & k_{\\mathit{AC}} x_{\\mathit{4A}} + k_{\\mathit{BC}} x_{\\mathit{4B}} + k_{\\mathit{DC}} x_{\\mathit{4D}} + k_{\\mathit{EC}} x_{\\mathit{4E}} & k_{\\mathit{AD}} x_{\\mathit{4A}} + k_{\\mathit{BD}} x_{\\mathit{4B}} + k_{\\mathit{CD}} x_{\\mathit{4C}} + k_{\\mathit{ED}} x_{\\mathit{4E}} & k_{\\mathit{AE}} x_{\\mathit{4A}} + k_{\\mathit{BE}} x_{\\mathit{4B}} + k_{\\mathit{CE}} x_{\\mathit{4C}} + k_{\\mathit{DE}} x_{\\mathit{4D}}\n", "\\end{array}\\right)$" ], "text/plain": [ "'Fi=' [k_BA*x_1B + k_CA*x_1C + k_DA*x_1D + k_EA*x_1E k_AB*x_1A + k_CB*x_1C + k_DB*x_1D + k_EB*x_1E k_AC*x_1A + k_BC*x_1B + k_DC*x_1D + k_EC*x_1E k_AD*x_1A + k_BD*x_1B + k_CD*x_1C + k_ED*x_1E k_AE*x_1A + k_BE*x_1B + k_CE*x_1C + k_DE*x_1D]\n", "[k_BA*x_2B + k_CA*x_2C + k_DA*x_2D + k_EA*x_2E k_AB*x_2A + k_CB*x_2C + k_DB*x_2D + k_EB*x_2E k_AC*x_2A + k_BC*x_2B + k_DC*x_2D + k_EC*x_2E k_AD*x_2A + k_BD*x_2B + k_CD*x_2C + k_ED*x_2E k_AE*x_2A + k_BE*x_2B + k_CE*x_2C + k_DE*x_2D]\n", "[k_BA*x_3B + k_CA*x_3C + k_DA*x_3D + k_EA*x_3E k_AB*x_3A + k_CB*x_3C + k_DB*x_3D + k_EB*x_3E k_AC*x_3A + k_BC*x_3B + k_DC*x_3D + k_EC*x_3E k_AD*x_3A + k_BD*x_3B + k_CD*x_3C + k_ED*x_3E k_AE*x_3A + k_BE*x_3B + k_CE*x_3C + k_DE*x_3D]\n", "[k_BA*x_4B + k_CA*x_4C + k_DA*x_4D + k_EA*x_4E k_AB*x_4A + k_CB*x_4C + k_DB*x_4D + k_EB*x_4E k_AC*x_4A + k_BC*x_4B + k_DC*x_4D + k_EC*x_4E k_AD*x_4A + k_BD*x_4B + k_CD*x_4C + k_ED*x_4E k_AE*x_4A + k_BE*x_4B + k_CE*x_4C + k_DE*x_4D]" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "Seja \\(f_i\\) o vetor de fitness de cada tipo \\(i\\):" ], "text/plain": [ "Seja \\(f_i\\) o vetor de fitness de cada tipo \\(i\\):" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "\\(\\displaystyle \\verb|f=| \\left(\\begin{array}{rrrr}\n", "\\frac{1}{5} & \\frac{1}{3} & \\frac{2}{5} & \\frac{1}{2}\n", "\\end{array}\\right)\\)" ], "text/latex": [ "$\\displaystyle \\verb|f=| \\left(\\begin{array}{rrrr}\n", "\\frac{1}{5} & \\frac{1}{3} & \\frac{2}{5} & \\frac{1}{2}\n", "\\end{array}\\right)$" ], "text/plain": [ "'f=' [1/5 1/3 2/5 1/2]" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "Nosso modelo fica então \\(\\dot{x_{ij}}=f_i \\times x_{ij} + F_{ij}- E_{ij} -\\phi x_i\\):" ], "text/plain": [ "Nosso modelo fica então \\(\\dot{x_{ij}}=f_i \\times x_{ij} + F_{ij}- E_{ij} -\\phi x_i\\):" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" } ], "source": [ "pretty_print(html('a) Vamos escrever o modelo evolutivo, começando com 4 tipos e 5 localidades, denotadas pelos vetores $X$ e $L$, respectivamente:'))\n", "var('x_1 x_2 x_3 x_4 A B C D E')\n", "X = Matrix(SR,[x_1, x_2, x_3, x_4])\n", "L = Matrix(SR,[A,B,C,D,E])\n", "show(\"X=\", X,\", L=\",L)\n", "\n", "D = Matrix(SR, X.ncols(),L.ncols(), lambda i,j: var('{}{}'.format(X[0,i],L[0,j])))\n", "pretty_print(html(\"Seja $D$ a matriz com as densidades em cada uma das localidades:\"))\n", "show(\"D=\",D)\n", "pretty_print(html(\"Seja $M$ a matriz com as taxas de migração $k_{ij}$ entre cada uma das localidades. como não faz sentido fala de migração de uma localidade \"))\n", "pretty_print(html(\"para si mesma, sua diagonal é zero. \"))\n", "M = Matrix(SR, L.ncols(),L.ncols(), lambda i,j: var(f'k_{L[0,i]}{L[0,j]}'))\n", "MD = Matrix(SR, L.ncols(),L.ncols(), lambda i,j: 0 + int(i==j)*sum(M[i])) # Diagonal de M\n", "MLT = Matrix(SR, L.ncols(),L.ncols(), lambda i,j: (1-int(j>=i))*var(f'k_{L[0,i]}{L[0,j]}')) # Triangulo inferior\n", "MUT = Matrix(SR, L.ncols(),L.ncols(), lambda i,j: (1-int(j<=i))*var(f'k_{L[0,i]}{L[0,j]}')) # Triangulo superior\n", "M = M-(identity_matrix(5).elementwise_product(M))#+MD # MLT+MD+MUT\n", "show(\"M=\",M)\n", "pretty_print(html(r\"Logo, $F=D \\times M$ é uma matriz em que cada elemento $f_{ij}$ representa a IMIGRAÇÃO do tipo $i$ em cada patch $j$ vindo de todos os outros patches:\"))\n", "Fi = D*M\n", "show(\"Fi=\", Fi)\n", "\n", "pretty_print(html(\"Seja $f_i$ o vetor de fitness de cada tipo $i$:\"))\n", "f_i = Matrix(SR, [1/5, 1/3, 2/5, 1/2])\n", "show('f=', f_i)\n", "pretty_print(html(r\"Nosso modelo fica então $\\dot{x_{ij}}=f_i \\times x_{ij} + F_{ij}- E_{ij} -\\phi x_i$:\"))\n" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Onde a emigração $E_{ij}$de cada tipo $i$ a partir de cada local $j$ é dada pelo efluxo total de cada localidade:\n", "\n" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\left(\\begin{array}{rrrrr}\n", "k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}} & k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}} & k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}} & k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}} & k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\n", "\\end{array}\\right)\\)" ], "text/latex": [ "$\\displaystyle \\left(\\begin{array}{rrrrr}\n", "k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}} & k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}} & k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}} & k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}} & k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\n", "\\end{array}\\right)$" ], "text/plain": [ "[k_AB + k_AC + k_AD + k_AE k_BA + k_BC + k_BD + k_BE k_CA + k_CB + k_CD + k_CE k_DA + k_DB + k_DC + k_DE k_EA + k_EB + k_EC + k_ED]" ] }, "execution_count": 5, "metadata": { }, "output_type": "execute_result" } ], "source": [ "# Vamos somar as linhas de M para gerar os efluxo de cada localidade\n", "of=matrix(sum(M.columns()))\n", "of" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Multiplicado pelas densidades relativas de cada tipo em cada localidade" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\\(\\displaystyle \\verb|Ei=| \\left(\\begin{array}{rrrrr}\n", "{\\left(k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}}\\right)} x_{\\mathit{1A}} & {\\left(k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}}\\right)} x_{\\mathit{1B}} & {\\left(k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}}\\right)} x_{\\mathit{1C}} & {\\left(k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}}\\right)} x_{\\mathit{1D}} & {\\left(k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\\right)} x_{\\mathit{1E}} \\\\\n", "{\\left(k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}}\\right)} x_{\\mathit{2A}} & {\\left(k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}}\\right)} x_{\\mathit{2B}} & {\\left(k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}}\\right)} x_{\\mathit{2C}} & {\\left(k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}}\\right)} x_{\\mathit{2D}} & {\\left(k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\\right)} x_{\\mathit{2E}} \\\\\n", "{\\left(k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}}\\right)} x_{\\mathit{3A}} & {\\left(k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}}\\right)} x_{\\mathit{3B}} & {\\left(k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}}\\right)} x_{\\mathit{3C}} & {\\left(k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}}\\right)} x_{\\mathit{3D}} & {\\left(k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\\right)} x_{\\mathit{3E}} \\\\\n", "{\\left(k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}}\\right)} x_{\\mathit{4A}} & {\\left(k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}}\\right)} x_{\\mathit{4B}} & {\\left(k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}}\\right)} x_{\\mathit{4C}} & {\\left(k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}}\\right)} x_{\\mathit{4D}} & {\\left(k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\\right)} x_{\\mathit{4E}}\n", "\\end{array}\\right)\\)" ], "text/latex": [ "$\\displaystyle \\verb|Ei=| \\left(\\begin{array}{rrrrr}\n", "{\\left(k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}}\\right)} x_{\\mathit{1A}} & {\\left(k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}}\\right)} x_{\\mathit{1B}} & {\\left(k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}}\\right)} x_{\\mathit{1C}} & {\\left(k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}}\\right)} x_{\\mathit{1D}} & {\\left(k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\\right)} x_{\\mathit{1E}} \\\\\n", "{\\left(k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}}\\right)} x_{\\mathit{2A}} & {\\left(k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}}\\right)} x_{\\mathit{2B}} & {\\left(k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}}\\right)} x_{\\mathit{2C}} & {\\left(k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}}\\right)} x_{\\mathit{2D}} & {\\left(k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\\right)} x_{\\mathit{2E}} \\\\\n", "{\\left(k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}}\\right)} x_{\\mathit{3A}} & {\\left(k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}}\\right)} x_{\\mathit{3B}} & {\\left(k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}}\\right)} x_{\\mathit{3C}} & {\\left(k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}}\\right)} x_{\\mathit{3D}} & {\\left(k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\\right)} x_{\\mathit{3E}} \\\\\n", "{\\left(k_{\\mathit{AB}} + k_{\\mathit{AC}} + k_{\\mathit{AD}} + k_{\\mathit{AE}}\\right)} x_{\\mathit{4A}} & {\\left(k_{\\mathit{BA}} + k_{\\mathit{BC}} + k_{\\mathit{BD}} + k_{\\mathit{BE}}\\right)} x_{\\mathit{4B}} & {\\left(k_{\\mathit{CA}} + k_{\\mathit{CB}} + k_{\\mathit{CD}} + k_{\\mathit{CE}}\\right)} x_{\\mathit{4C}} & {\\left(k_{\\mathit{DA}} + k_{\\mathit{DB}} + k_{\\mathit{DC}} + k_{\\mathit{DE}}\\right)} x_{\\mathit{4D}} & {\\left(k_{\\mathit{EA}} + k_{\\mathit{EB}} + k_{\\mathit{EC}} + k_{\\mathit{ED}}\\right)} x_{\\mathit{4E}}\n", "\\end{array}\\right)$" ], "text/plain": [ "'Ei=' [(k_AB + k_AC + k_AD + k_AE)*x_1A (k_BA + k_BC + k_BD + k_BE)*x_1B (k_CA + k_CB + k_CD + k_CE)*x_1C (k_DA + k_DB + k_DC + k_DE)*x_1D (k_EA + k_EB + k_EC + k_ED)*x_1E]\n", "[(k_AB + k_AC + k_AD + k_AE)*x_2A (k_BA + k_BC + k_BD + k_BE)*x_2B (k_CA + k_CB + k_CD + k_CE)*x_2C (k_DA + k_DB + k_DC + k_DE)*x_2D (k_EA + k_EB + k_EC + k_ED)*x_2E]\n", "[(k_AB + k_AC + k_AD + k_AE)*x_3A (k_BA + k_BC + k_BD + k_BE)*x_3B (k_CA + k_CB + k_CD + k_CE)*x_3C (k_DA + k_DB + k_DC + k_DE)*x_3D (k_EA + k_EB + k_EC + k_ED)*x_3E]\n", "[(k_AB + k_AC + k_AD + k_AE)*x_4A (k_BA + k_BC + k_BD + k_BE)*x_4B (k_CA + k_CB + k_CD + k_CE)*x_4C (k_DA + k_DB + k_DC + k_DE)*x_4D (k_EA + k_EB + k_EC + k_ED)*x_4E]" ] }, "execution_count": 6, "metadata": { }, "output_type": "execute_result" } ], "source": [ "for i in range(2):\n", " of=of.stack(of)\n", "E = of.elementwise_product(D)\n", "show(\"Ei=\",E)" ] }, { "cell_type": "code", "execution_count": 40, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:14:53.755001Z", "start_time": "2021-10-27T12:14:53.740582Z" }, "collapsed": false }, "outputs": [ ], "source": [ "def sel_model(t, D, params):\n", " D = np.array(D).reshape((4,5))\n", "# print(D)\n", " f, M= np.array(params[0]).reshape(4,1),np.array(params[1])\n", " F=D@M\n", " of = (M.sum(axis=1)).reshape(1,5)\n", " E = np.vstack([of,of,of,of])*D\n", "\n", " phi = sum(f*D, axis=0)\n", " res = f*D + F - E - phi*D \n", " return list(res.ravel())" ] }, { "cell_type": "code", "execution_count": 41, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:15:02.959817Z", "start_time": "2021-10-27T12:15:02.945998Z" }, "collapsed": false }, "outputs": [ ], "source": [ "T2=ode_solver()\n", "T2.function = sel_model\n", "T2.algorithm = \"rk8pd\"" ] }, { "cell_type": "code", "execution_count": 42, "metadata": { "ExecuteTime": { "end_time": "2021-10-27T12:15:04.360267Z", "start_time": "2021-10-27T12:15:04.093535Z" }, "collapsed": false, "scrolled": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[1. 1. 1. 1. 1.]\n" ] } ], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "inits = np.random.random((4,5))\n", "# Normalizing population on every locality\n", "s = inits.sum(axis=0)\n", "inits /= s.reshape(1,5)\n", "print(inits.sum(axis=0))\n", "inits = inits.ravel()\n", "t_span = [0,70]\n", "m = np.random.random((5,5)) # random migration rates\n", "np.fill_diagonal(m,0)\n", "m = np.around(m,2)\n", "f = np.array(f_i)\n", "T2.ode_solve(t_span,list(inits),num_points=300, params=(f,m))" ] }, { "cell_type": "code", "execution_count": 47, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/av/WaAAAACXBIWXMAAAsTAAALEwEAmpwYAABnDElEQVR4nO2dd3yT1f7H3yfpLtCyNxQQGSJLcbFREQT3RBEVJyp6vb/rQL3XfXFvxb2Fy3AgVsGBiCwRFQEFRIECsgttadORJt/fH+dpSdMkTdukTdrzfr3yKuQ5z3lO2iSf55zz/X6+SkQwGAwGgyHSsNX2AAwGg8Fg8IURKIPBYDBEJEagDAaDwRCRGIEyGAwGQ0RiBMpgMBgMEYkRKIPBYDBEJEagDAaDwRCRGIEyGAwGQ0RiBMpgMBgMEYkRKIPBYDBEJEagDAaDwRCRGIEyGAwGQ0RiBMpgMBgMEYkRKIPBYDBEJEagDAaDwRCRGIEyGAwGQ0RiBMpgMBgMEYkRKIPBYDBEJEagDAaDwRCRGIEyGAwGQ0QSU9sDMBgMBkP1SZuS3hjoBMQDhcCWrVPHHKzdUVUPJSK1PQaDwWAwVJK0KekKGALcDAwCGgMOjyZJwEFgCfAcsHjr1DFR9YVvBMpgMBiijLQp6aOBaUBTIBlQAZoLkAdkApO2Th3zRfhHGBqMQBkMBkOUkDYlPQV4GTgTPUOqLA5gLlqoskM5tnBgBMpgMBiigLQp6S3Qy3XtgYRqdFUAbAcGbZ06Zm8oxhYujEAZDAZDhGPNnH4E0oDYEHTpBLYCAyJ5JmWi+AxRS12MWjIY/DAN6EBoxAmrn/ZWv5eEqM+QY2ZQhqihPkQtGQzeWAERcwiw5+TM3M6hn9Ip2LaG4kP7EWcR9qRGxLXsTOKRJ9HgqOGoGJ/a5gDOj9TACSNQhqigvkQtRSgKaAu0RM9WE4L4mYj+OyVZjwTgb2A6sLxmhx+9WDdlW4CO/tpkLZlB9tIZIG7i23YnrlVXVFwCrrwsCretpThrN3GtjqD1Fc/46yID6BSJN3Nmic8Q0VQxakkBDazHnLQp6VETtRSBdASWoWerRV7HlMfDBtg9Hv5uICYCFwLp4RhsHWQI+qbMJ9nLZpG95APsjZrT/Ow7iW/TrVwbx58ryVn5caBrNAUGA4urO9hQY2ZQhoilPkYtRRjNgF/RMyd7CPstANqgl2MNAUibkv4hcA4+BL84aw9/v3odAK2vfIa45ml++5Fip78lPtArDh9tnTrm/OqON9SYGZQhIrFmTksITdRSgtXPkrQp6REdtRRhjAZSCK04gQ5oOQr9960XKKXsBL88WvIzof3/zTnVFpvgczaau/ZrcBeT1GNIQHECAokTaPEbVJnXU1MYgTJEKj6jljIeGVu+pT0Ge3ITEjr0IuWEC4ht1t5Xf1ERtRRhtMTHzYHb7SYnJwcRoWHDhsTElP0aycvLY/ny5WRnZ5OSksKIESOw2cr5UqeGbdQeKKViqLwwhONnpb9rbfHJKJv/e4PCHb8BkJDWp7Jd+6JJ2pT01K1Tx2SForNQYQTKEHFYARFnoT/YPkkZOK703+5CB0W7/iBv3UIcG5fRavxjxLXs7Ou0BOCstCnpo03gRFA40cs/pRQWFvLCCy/w1FNP0a5dO6666iquvfZaRASlFMXFxcyZM4d33nmH1q1bU1BQwKZNm5g0aZJnm5g5c+YcM27cuGLCLwyhnv3VGDGprRBnEcrue/bjytUrpDENm4XicvlAZ+DnUHQWKoxAGSIKK2ppGhUERKQOvrTccwe+fJlDP39Gzo9zaTb2Vn+nJgHT0qakR2TUUnVQSin0jCckX/CPPvpo/1tuucUeH3/4PmHt2rW89957/P333/zxxx9ceeWVjB07ljZt2gCwbds25syZw+OPP84xxxzDV199xeTJk8sIlMPhSJ43b959NfE7iWrssYQu7alChAA3hLWFEShDpBEwaikQCZ36cejnz3DnV7jFFNKoJUsY4oiMpaRA4feV4vfff6eoqAhPgdq0aRO9e/cGoEmTJowcOZKPP/6YG2+8EYCmTZvSt29fpk+fzoYNG1i5ciW33XYbQOkyn81mIyGhOjEv9QSXk0B/TnuDxjgzt1N8KDMUV1PovcGIwgiUIdK4GZ0/U2kKtq4GIK5V14DtRCTZuS/jeaXUDEInDHWOgoICPKN8RYTs7GxatWoFQFJSEo0bN2b37t2lx1NSUjjrrLO45ppr2LBhAwUFBZx66qll+rXZbGVEr/6Rgn6L7wzYqjhrNyrWf03Z+HZHUZCxhoKMX2nYZ2R1B5UIbK5uJ6HGCJQh0hhEELOArO8/KP23u8hB0a5NFO5YT2KXATQ67pyA5yqllD0ppTfQu7qDrcsUFhbinYYSExOD0+kEQCmFzWYrbaOUIjs7m1deeYWXXnqJE088kXXr1jF27FjGjj0c3GK32+vZDKopOp2sUxF0ckNjF+wthGf3o2ctBb5+ugvzCsTpvETZ7T5/WQ2OPoXs5bNxbFxK0f6LiGvWwe8IKggzBzgQaQESYATKEEFY3nqNg2mbvXRGuedim3UguedQbPEV5/PaEhug4pORwrxKj7O+4D2DUkrRrl075s6dC2ih2bJlC0cffXRpG6fTyXfffccTTzwBQJs2bUoFrQS73V6TMyjBjwCE52dyEVyZDKOaQ7cW0Kw1JMVDrMvrvkvgmcdEyA80+I63f95URM5UqvxNW0xqS1IHXULW4nfZO/t+najbuvzqQf7mn8he8SGtLvlvoN9RRIb8G4EyRBKd0N5gKRU17HjnZ6X/dhcV4NyfwcFF77B/3hMU7d9G46ETAp4vziJiU1tRtOev6o65zuJrBnX88cdzyy23sGXLFgoLC1m0aBF33nknO3fupGHDhjRq1IhTTjmFBx98kCFDhvDLL78wdOjQMn3ExMTQpk2bDGAj4ReOYgmjG4FS2IHW6ClSN3RqhPeMx+XrVHTawx8B+7fL0+JSo7CLT0VPOelCxO0ie+kMdr9zK/FtexDX6ghUXCLuvCwKtq+j+ODOipa989DelRGHEShDJFGl22pbXALxbbrR4ty72PHiFeT88CEN+40mplHzAGeJFSVV5yimcl/ifo8deeSRrZOSkq5H708A0LhxY5599llGjx5NbGwst912G82bN2fq1KmccMIJDB8+nPvuu49///vfzJw5kyZNmvDss8+WG+TkyZNnT548+bbw/RrCg1LEon0JO1qP9lQ91K4jFQgUsBib7EELn09SB40jufsgDv2szWJz136NFDuxJzYktmVnGp1wPg2OGh7oGpnA95UffvgxAmWIJKoVRWRLaEBsk7YU7fmLot1/VSBQyoqSChlOamwZyf9PEfF1t15VjkF75yV6Pjlq1Cg2bNhQpuGUKVNK/92iRQteeeWVivquSjXYWkEpuqCdSDqixSlUuVV+DWBL2Dp1jKRNSb9ehI+U8m/3FdusPU1GXl+VMTjQPpURmXJhBMoQSWyhml9c7oJc/Q9xB2yn4hLcLkfW80A2oRGGwBeMTg4RnkRXN7AnDP2Gi5FoV41QIcBuYJtSKBECikPGI2P+bHX5kk1xLXJ6KruE2hNxbiQnrRuBMkQMGY+MzW7/fx/m22Ljq7Rk4vhjOcXZe8AWQ3y7HgHbKpt9X3HO/n9U5Tr1iC3o/K5QUwB8HoZ+w0UG1RMoFzqmPMN6bBehIJgTlaIdcOG+D4/9ptWEpUfYkwsSlD0kuW5OtIHypBD0FTaMQBkiAqXU8cCz+ZtXNUrqeiKqvHdbGTzDzMVZiDNzG/l//QRA6tAJ2JMDBgNGbNRShOFEL/G9htcyXxAIej/MhZ4xua1/x6HLbawK3TCDw4qES+Xw/lGcCLODODUDOK4SlypGf/mXCNIOESq9nqwUTdG+kXGu3IQjdr930i+tLlvWz55UGK9iJPAHJDAFHHb3j2jjZCNQhlpFKdUWmApcBnBo1TwSO/VHxQX+PiwTZq5s2JNSSDziOBoeM5bETv0qumzERi1FIB8AWcBd6FlEMYeXN/Otf+ej9zLy0b/bPOvfJUugJe0PAL8AO2pi4JYgNeOwIHUEGnk0cStFnEi5OlfebKvgeKHVpkSQdor4jNwLGqVoiP5MJKGjWxu6DiU6d701+KemY1cfmdAxs5ktxl0VkXIAnwA3RLo4gakHZagllFKJwD/RX3xl9p3aTnqTmJQW4bx8xFYQNVQdpbChRbREjDpQsSvJeyJUmGugFJM5bMHl4LAYZQB7RAjZHqQVDHEF0Ao94+tD2SQqd4N+W7NTh24cZosvjlWKJAInt0dthWkjUIYaxfKtOw94HB0ZVY6EzsfQ/Owp2OJC7zYgbleRstnPjqYPqcE3XjlIJYJU2TfNYhEWBnGtY9AikAHsryiwoaooRQxwKXrWFAMMoGz6hQC/AfuxuT/vePsXCWh7sEFAE/TMVayxJqJnrUvQKwbfR9tNmREoQ42hlOoLPAMMDdwSmp5xG0ndThRbTFzIzE/dxUXk/7Fc9n/6+AUi8mGo+jXUDCHOQSohQ4S3qju2UKEUY9CiBLqoo3euxB/ogIs/gBmeQpk2JT0VXTIjHr3suDkS7YsqgxEoQ9hRSrUAHgSuITi37QMxKS0fbnP969crpdIIQc0BcTkpzt7LrnduRQodTuAcEUmvbr+GmkEpzkV/YYc67L0YeESE4hD3WyWUojF676kn0N3r8FbrkQtME6HO+3RVJxLEYAiIUipOKfVPYBNwLRWLkwu9FNHVmbX7KaXUIPQHMqiQXH+4i4sozt7L7vdvRwodoAXvQ6XUKdXp11CjOAmNOAmwC1gBzASejhRxAhDhIPAm5SsO70J/FgA+qQ/iBGYGZQgD1j7T6cBTwJFBnrYA+KeI/O75ZNqU9BR0AcOzqEISr9tZgOOPFRz48qUScfIkHxglIiGpC2WoHFakWnORiss8KEVv4NwqXKbKOUi1iVKkAg+jZ42ZwDq0uK4QYX4tDq1GMQJlCClKqR5oYRoV5CmbgFuBzwOZelpl4KehI6mSCTJqKXfdwucyP3vqAfxHc+UCp4rIiiDHa6gCPnKQOqI39V3oJbaAeUJKkYJ+n1RESHKQIgErCOQqdDSfDe2+8VokzfjCjREoQ0hQSjUG7gVuIrilmBzgfuAFEakoDwUoLQc/GLhZitUIbJKC2+YsCVnC5o5FVLayy0I8opaUUsOAL/Af4ZUNjBCRn4MZh6FigshB8uRtkdLlq0B9/oPyS18hz0GKNJSiOXoF4VMR9tb2eGoSI1CGaqGUikHvLz1AcKXaBXgduEdEqvxhU4qh9kb5o+NbZzVWMS67FNtdhbtSD7pyEr8S4Rsf4xwJzMO/dc8BYJiIrK3qmOozVcxBKuFbEb4L4hrnAF0JYw5SpBKMZ19dxDhJGKqMUupkdNh4ryBPWQz8Q0R+CcHli1w5iQWOnMRdXs/7jPgTkS+VUucDH+H7fd8E+FopNVRENvg4bvDAWn5qQ1lBqmoVwgpdvS0+B4qi/YtaKeyVneVF+2uuKkagDJVGKdUFeAI4O8hTMoB/AR+GsHicv30Fv+amIjJPKXUJ8D98R7C2AL5RSg0REVPJ0AvLuPQItKC0IwTh/xbtg/nSFqleOZZIQCnigcuVYp0Iy2p7PJGOEShD0CilGgJ3ozerg3G5dgD/BZ4SkYClrauAv32rgF+aIjJbKRUPvIvvQIs2wEKl1GARqciDrb5xLNA3xH3uQ9/AxEHg8ufRjjXrvBD9HmtjRTF+WV9nR8FgBMpQIUopGzABberaKsjT3gOmiMjfYRpWpWdQJYjI+0qpBLRLty86oEVqiIjsrOoA6yAZVE+gSuoglewfbas3+Tw6aOQsoIvH0ycBfZTiqfoUmVcZjEAZAqKUGojeZzo2yFNWArfUQNh2lWZQJYjI65ZIPe+nSRf0ct/Q6gRzRCrW3XvJ/lFb4I0g9kUyKnmZqMxBChOnAr29nktD//67KcWdIuyr8VFFOEagDD5RSrUHHgXGBXnKLuAO4IMaqi5b5RlUCSLygiVSj/tp0h0dODFcRDIrO8BIIUAOkidt0PlDgTiIrrLb0M/xOpODFEqU4kT0bMmTVPRMHXS049VK8Tmwxiz5HcYIlKEMSqkk4Da02ARTpK4QHTDxiIjkhnNsXlRrBlWCiDxhlf54wE+To4EFSqlTRCSrMn3XFpXMQSqhIxUIlAiiFBkcjtqs8zlI1UUpjgZO83o6Bn3zo9AuERvREZDnAF2VYq4Rdo0RqLpPP7Qj8hJ00IJPLHuii4DH0C7RwTAHuF1EtlR3kFXgAPA2WqicXj8ry0NoMZ7i5/gxwOdKqdNE5FAV+gc96zgRGA4MRN81/weotmGtjxykjlTeFqojwVUZXo0uOFhvcpCqilJ0xnek65HopPEc4HcoM2NKALMfVYJJ1K27DABeAbqhv7gT0QmyN3o3VEodAzyL/uIMhl/R+0wVJldGC5ZAPwX8I0Cz74DTRcSv0HvQBF2j52T0HXRn9A1CAw47bTiAueiy3pUYa0hzkEooBB41ghMalKI1cCXll5xboWdPDnR1Yc+ZUh7apbwmVyIiGiNQdY/WwNPAmei7Mc9Qagfa9eEDAKVUK7Qh5ZX4Drn2Zj86zPwNEalzSzmWSL0EXB+g2ZfAWSLivdnfChiC3gw/Ff13KEALUqCqAflAfyCo5GClGIK2ewpVDpInL4uwOwz91iuskhlXof/2niSig41cwM+Ud+mfLsIf4R9h9GDKbdQdEtDi8Sd6LTuR8qKTBLwwc+bMpkqp29FFzyb6aOdNMXp20VVEXq2L4gRgJRHfiF469MdIYNYff/xxBLpuz/voSLWt6LD1q9Czmjj0vk9FnzEbMLYSw3QQOnHaB6wCPgSeMuJUfZQiGf2+8BYnBfRAL+etobw4rTTiVB6zBxX9KPQ690vofY6Aew/FxcVJmZmZf1LedNMf6cD/icjGaowxahARt1LqarTgX1zy/JFHHsmQIUMYPXo0Q4YMOaNhw4ani4jDSl4uoSrLbPHFxZwZG8tnwB9BLLFVNtS7hHqbg1RTKEUcernWO0ISdEh5A2AtlFvC2wt8FdbBRSlmiS+6ORp9196LAMacIoJevdI4HA569erFli0BYxs2ALeKSL2pPePJV1991XjRokVfDx48uP9JJ52E3W5HRGjQwPvGODicTicHDx6kRYsW5Y4VFeFMSmKqy8UrIgRMDLYi9G6j4iAIN/A3JgepRrD2BcehraC8SUEnOK+Hcm7kxegSGnvCOsAoxQhUdNKMwzlK8fhYRsrJyWHp0qUMGzaMxMSy0eJOp5OFCxcyapTPkk1Z6LIZ00Skvoa6NgX+EJE4pVTVFAlwuVx88MEHvPbaa+zfv5++ffty1llncfHFF+N2u7HZ9J+tsJCiCy7gzXnzeFeE5RX1qxQXoZeLPDE5SLWEddNwNtDHx+EY9L7TDuvhzXwRTC0yP5g9qOgiFu2DtxW4FL3PVO5v+OKLL9K/f3+ef/55rr32WnbuLHtTHhsby6BBgxg+fLjn0250QcCuIvJcPRYn0EEmyVURpzfeeINLLrkEp9OJiNCmTRsefPBB1q9fz7XXXsukSZPKnaMU6vTTSyPygiEDHXW3CfgaeAOYKsI7IiwSYYsRpxrlZHyLE+iQ8r34Fqc/gR/CNai6gNmDih5GAa+i17f9Luf9/vvvfPrpp/zwww80bdqUY445hj179tCmTZsy7ZKTk3nttdfo1q0bLpdrIboMhqmFpPf0LiG4ootlyM3N5fXXXwdgx44ddOrUiVNOOaXM8fPOO6/M7AkgLo7YwYPpCnQIsu7PT+hNdRMSXssoxfHodAJftETv/fkqaZ8HfGJcIwJjZlCRz5HAQnSkVXt8iNP+/ftL/x0XF8f+/fvJzMxk9+7dNG7cmN9++43s7OxyHbdq1Uo++uijV4FTolGclCJFKVorRUelOEIpeipFH6UqdE0IhFB+E5v8/HxeeeUVbrrpJtLTfefWbt++nQEDBtCgQQO2bTtshJ6RkcGYMWO46KKLSE1NZfPm8t9XXbrQXimS0EnVgQcoOI041T5KcRT6xtEXieilYn/pA3NNvlPFGIGKXFLRRqar0fk15TbFCwsLue222xg9ejSPPPIICxcu5IgjjuCcc87hzjvvpFu3bpx00knMnj2bBx98kF9+KVsnMDk5WZ155pnjRCQ1/C8nLJwBXIfO4xqPLmVwDsE7rvujXOn3OXPmMG/ePI466ij+85//sHjx4tJjJfu4zz33HOPHj6d79+5lftctW7Zk5syZOBwOGjRowAsvvEBhYdnSRjYbavx49oK5o44GLPeOIfhO0VDoPLiN+P57mpDyIDECFXnY0V+6GeicmkT8LDc98sgj7N69my+++IK+ffsyfvx4HA4H99xzDwMHDmTq1Kk88MADPPXUU2RnZ/P33z4rX8QC/xeuFxNm/NkaBW0Y64d0vPJUnnvuOW6//XYmTZrEddddx9y5c1m3bh2gBSo/P5/U1FSOO+44Tj31VP744w82b96My+UiISGBpCR9fzFq1CgWL15MfHzZiPS4OJzvvovTOFpHB9YM9h187y01sZ73lS+4DxNSHjRGoCKLzsA64El0kqdfs9aCggJ27tzJDTfcQLNmzRg1ahSNGjXi5ptvBqBp06YcPHgQgC5dupCVlYXb7XNVKAG4PMSvo6bwFwgQdCKrUiQpRXelOE0prlWKlujS9KXiV1xczLHHHlsq8EOHDsVut5fOkmw2GytWrODbb79lxowZPP3007z22mtcfvnl7N27l4KCAmw2G3l5eaxZs4b+/fvjI3o2mfKmooYIRgQHuvDlJo+nY9H7S75unlzAHBPAEjwmSCJyUMBy9Lp1uRlTcXExMTGH/1wJCQkUFRUxe/ZsevTogcPh4KSTTmL+/PlkZGTQtm1b5s2bxx133MHu3bvZu3cv3bt393ftTf4ORDiVnkF51UHqiC7z7klH9BJf6RRHRGjfvn3pvlKrVq1o2bJlmb2kLl26kJuby4IFC5g4cSIul4tnn32W1q1b8/XXXzN16lQyMzNp3749jz/+eJm8NA+6VvySDZGECEVK8T+0G0gf9HvSnzPLVybfqXIYgYocjkbvM5URJ5fLxZ133klRURFnnHFGmaiwhx9+mDvuuINJkyaxZs0ann/+eZo2bcrrr7/Ogw8+SLNmzXj//fc5+uijeeedd/xdNx/tSh6NBJxBBVkHyZuO6KKLvwAngA7L79ChQ+m+U8OGDbHZbDRqdDgWo0OHDqVLfgCfffYZP/30E/369aN379689NJLdOvWLdB1i9DO1oYoQwSXUnwK7EFbYfkSqL8wIeWVxghU5NDf+wkR4eabb+bgwYOMGTOGRx99lI0bN3L11VcTHx9PmzZtePPNN9m8eTOpqam0bNmSNWvW0LRpUwD69etHnz59yoQ0e+EAZqBdzqMRXzOoJKCXtVQXTB0kbzpaod6foctsxAIMGjSIV199lXXr1tGrVy+WLVvGueeey19//cVPP/3EiBEjaNasGYWFhcTHx/PKK6+QmpoKQIsWLXw6SHjhROe4GaIQK1x8hVL8DZwLNPY47MCElFcJI1CRw1d4/T0OHTrEL7/8wvz582nUqBHNmjXj888/Z/bs2YwfPx63201sbCxHHnkkSilWrVrFRx99xF133VXahx9xykPb4FyD3m+JVorR/mapaDuZVLSgtEbfsVaFBuhZ1nfo2WUsQFpaGkOGDOHhhx8mJyeHgoICBg0aRGFhIW3btiUlJQWgNPihRJwC4LT6TwBWAPeAMWuNAuzWw+fysgjbleJl4HQOJ+/OFaGqdcTqNcbqKLJ4B100sHT/49JLL+X444/n5ptvJjc3lzlz5rBy5Uruueee0uTbgoICPvjgAx599FHuuusurrjiCn/956O/GP8FvInvKKOoQCl6AjdwuLqrJzuhWmG8H4uwAV1QrnQ/y+l0smTJEjZv3szYsWNp2bJl0B0WFhZis9kcsbGxdrSbdTqwCL2cmF+NsRrCiFJ0RBvrnou2F+uAFqeH0JWk/RYXtKrpthTh65oYa13EzKAiiEWLFt1z0kknjYuLO7zHf/bZZzN//nx27dpF69at6d27N2vWrCEzM5M2bdrw22+/cdRRRzFq1CgmTpzob/Pdhf5QvYau4lo+azf6yMH/+7eyLhAHOexhtxXtRyjoPaG+JY1iY2MZPny4t0WUT/Lz83E6ncTExLBq1SrS09P5/vvv7X/99dfpe/bsWVjJ8RlqAaUYAIx58kkKRPi3lUgNelZ9D9rt/hzAp+uyCGvR7uWGKmJmUBGAVSjvXOCJf/3rX2n33ntvqWv2rl27eOaZZ0hNTWXKFF2RfNCgQTz99NPs2rWL/fv3M2HChDIRfl440NGBk4jwaD2liEe7ZeypaEnEco9+Eg8B8WA/OlzfH/s4LEgZIuT4aXc/cAfBldFwAFJcXGxbtGhR3KeffmpfvHgxa9eu9Q7vzwVOFRFjEBrBKEUP4MIjjiD511+ZnJREgo9mLnS+3HVYRUANocUIVC2jlOqDLrc+FPRd+ubNm2nXrl1pm+XLl3PHHXcwefJkBgwYwMSJE3n66afp3bu3vxkT6H2m/egKul+G91VUDeuOtAOHI+xaoyOg5onwUxDn3waM8XHoAHoZDapXB6kjukSCd2Vi0EJjQ5u2Lga+sH5uUEoNtf7v60sN9Ax2hIiUOFYodJG7i9BLsD8AjxHFS7DRjLWsd5lSxPz2G5d17UrHmJiAs3IHsADtaFIXViciBiNQtYRSqjnwIDpQoUwkw+mnn87MmTPL1B764osvmD17NsuWLePGG29k8uTJ/rouRH/J3YMuYhgxSYFB5CCVsEaEj4Lo70L0PpQngl6m+4DQ1EEagE7G7IBeUjwIfAvMRwvSVnzY2SilTgM+xX9OViYwTETWo/cDz+Owz2Ieem9qNPrvaaghrOjPK4GE55/n2GuvZWRcXFCJ3wXoZedzgaXhHGN9wghUDaOUikOXFb8XHXnmk++++46TTjqpzNKd0+lEKeVvOc+N/jL7ALgT/QVYa1QxB6mEbBGeDuIaxwGPoO9as9F7RznAThFervyo/WJHJ9HuRc/OgkIpdSba5NfnHywmJmbvnj17Vjdp0mQQ5b0W84El6BlixNxk1GWUIgVtL9Zo5EiazZvHtUGKkyf56HIt9xIggMIQHEagahCl1OnoN++RFbXt3r07q1evdsfHxwdjR5WH3oy9hsB7L2HDEqRmlBWk6riKPyNCVgXXbA7cRPkZzAERnqvGtUOGUuoC4H94zZJjY2P5+OOPGT58uCQlJflbp3UA36Dvys2XXRhRikRgItC8YUPsf/3FDU2b0thmO7y0m52dXZpO4F2l2os8dK2nU8F4K1YHE8VXAyilegBP4d+a35s/N2zYcGtcXNwY4Ar872U40DOG69HLSTV6t2Et2R2FFqMOBKhTVQU6QmCBQkcm+nrNlb3rDRsiMlsplYBOIVCgc6XS09M54YQTCCBOoGdVJwPT0RFjpsRGGFCKWHQNsOYAs2dzakoKDT3F6fXXX2f+/Pnk5eXxxhtvlKuv5kUy0BMdnNQV41BfZYxZbBhRSjVWSj2D3rAPRpxygNuAo0TkM6XUXfjegyhCi9NDQBowl9r5EDRGv64ehE6cDqFngcFsNofLzTykiMh76EgvEhMT+frrrznxxBNJTg7qV5aETvp8E/8eb4YqYpXNOB8dPcqNN9J5+HD6ey7trVu3jscee4wHHniA9u3b8/LLL7Nr1y6cTicul984lligDdA23K+hLmNmUGFAKRUDXI0WkKZBnCLost33iIinmeRB9PJOOvrLON5q+xG6REZtOw/sRC89Ved95JmDlAEcrIQljF8vviAr09YYIvJa+/btUz788MPHe/XqVVp+w5Pi4mI2b97M1KlTady4MWeddRZDhw4FLf7no/c3bsDckYcEa1l6DNANoEsXEh97jPO9952mT5/OWWedRc+ePZk4cSI33ngj6enpHHfccVx88cUMHTrU35KfE33z5qskhyEIzAwqxCilRqDdsKcRnDgtBo4RkWu8xKmEheilh0vQkXm9gEsJsTgpRbxVlfZkpZhoZcEHRIRiKv/h2wesQgcPPCXCsyJ8IsIvIhyopKi48L3sZaMKJdvDTOr27dsv69evX7EvcXI4HEyfPp3bb7+dVq1aceKJJzJhwgRPx/RkdCj6k5iZVKgYivZbRCn49FPOi40tP/seNWoUGzdu5L777mPChAn885//ZNmyZbRt25Z7772X7Oxsf/tRsRxOdzBUATODChFKqc5o65NzgjxlG9pyaI5UHKmSi55F+a41XgUC5CCVsJ/gsuAz0MuMvqhODlKFiCBK4cR3Im0skRNYEIe2NeoWGxtb7jNXXFzMrFmzWLJkCZdffjnnnKPfQtOnT2fLli107ty5pGkyeqnQgb5ZMVQRpTgWGFby/2ee4ZguXegQG1v2xkZEGDhwIFu3bmXv3r107dqVwYMHEx8fzz333MPSpUvZsmULffv29b6EA53LZsprVAMjUNVEKdUQuAv4J8HtfTiAqcCTIlJjHmyVyEEqIS3IrjM8/u1Gm9CWCFJ1c5CCoWTp05s4Isfj7iZ05KZPR4rVq1cze/ZsJk2axNixYwH4/PPPcTqdHHvssd7Nk9Cu5w7gv+Ebct1FKbrjkeB98sk0vf56TvNc2svKyuKll16isLCQLl26MGHCBACSkpL45JNPGDhwIL/88gsHDx70JU4udLn3h8L+Yuo4RqCqiFLKBkxAi02rIE97H5giImFdk65mDlIJTZSiYRAuzDvQzt9bgR21UC202lV1a4B2+KmO7HK5ePLJJxkxYkSpOP3444+sWLGCoUOHkpSUhNvt9nalTwLuRovUM+Edet1CKTqg9/MUQHIy9unTGRcTU/a78PLLL6dp06acdtppPPTQQ0yfPp0XXniB0aNHc9999/Hdd9+RmZnJ+++/7+sy+ei9Y+MEUk2MQFUBpdRJaHuicre3flgJ3BIu/7Uw5CCV0JEK8qpEKEI7K9QWJaHmRWixKvkZSSj0GH1uVMTHx9OpUycAPv30U1auXInD4eCaa64hNrasznqIVRLwMNrBIJRJyXUWpWiB3sst/d6bPZtTUlNp5BlSvnHjRvbt28eMGTNISkrioosu4q677uK0007j7bff5vXXX8fhcFBUVETz5s29L+NA24ttrYGXVOcxibqVQCnVDm25f0mQp+xCuzq8LyJhyWGxrFkup7wTQShYLsKCMPQbMizTWHckRez5oAc6MMTn32jRokVce+21dOzYkcaNGzNgwAAuvfTSMrk2Dz30EAUFBWRkZPD000/TrFmzkkMOtBHwu2F+DVGNUjRCR9aW3rhdfz2dnn2Wcb7cIu644w769u3LhRdeiN2ut6U++eQTpk6dyvvvv0/Xrl19XaYA+Ay4ICwvoh5iZlBBoJRKQgc03EFwQlCIjraaKiK54RwbOkzbXyJvZTlE2bIT+0PUb9gQiYpllPXoCLz38bHUN2zYML744gtyc3Pp3LkzDRs2BA7Plu677z6WLl3KlClT+O677zj77LNZsGBBSR5VEnoGVQDMqrFXFEVYLhHj8RCnLl1IfPJJLvAWp5Jw8WOOOYZp06YRFxfHeeedh9Pp5Oyzz2bRokWsX7/en0Blo62SDCHCCFQArDIYF6KjcToEedqHwG0i4rNGTHDXJRad4JcvEjgKSIQipdiJ3ueoLNXJQTJUjo/QVlSv4UOkunTpgojw7LPPcuaZZ9K5c+fSfadDhw5xxx13MGLECEaMGMGuXbv4+++/OfLIUsesROAt9I3R3Bp5NdFHacK7UjB3LucG8tm78MILSUhI4D//+Q/Lly9nzJgxNGzYkM8//5wzzjjD1yn56Ahef6VbDFXACJQflFLHoDegBwV5yhr0PtOiyl+rtA5Syf5RW3Qezyr0kkFFZBCcQAVbB8kQHj5AR/I9j4+ZuFKKU045hZwc/WdxOp243W62b99ORoYOliwRJx+GwUloS6RzIbKXZWsaEfKV4l10cES3p56iv1VCo/SX6HK5sNvtZfKZzjzzTAYOHMhDDz3ERx99xM6dO7niiis4+eSTvS+Rh/6uWB7+V1O/qDN7UGlT0hsDndBfAIXAlq1TxxysbD9KqZbozeeJBJcQuR+dk/K6iAS13BREDlIJ+0R4MYj+jqT8vlhYc5AM1eIm9F6m3+XijRs3UlBQQJ8+fVi2bBkvvvgizZo145VXXmHKlCnce++95OTk8Oeff9K/f3/PUx3AWGo3cCUiUQrb6NFM+OQTXvWcPeXk5PDEE0+wf/9+hgwZQrt27Rg0qOx9aX5+PomJvgMx0YFExxI5eXd1hqgVqLQp6QoYAtyMnuU0Rn84S0hCL2EtAZ4DFm+dOsbvi1VKxVt9/RtoGMQQitF3wg+ISFaghlXIQfLk8YqERSkS0Htku6jZHCRD1bkNuA8/IrVkyRKuvPJKvvrqK9LS0vj222+59NJLufLKK3n44Yf5+eefmTJlCna7nZiYGD799FPP0/OA44Hfwv4qogwRpolwlc12WKCGDRvGcccdR7NmzSgsLGTHjh306tWLiRMnkpycTF5eHg6Hw1fEHujf9dH4KftuqB5RKVBpU9JHc9hKKJnAMx1Bv4kygUlbp475wvOgtc90BtptvEuQQ/gc+KeIbPQ+EKIcJE9mirC+okZKEWNZDxmih3vRQuXTNfbFF1/kvffeo2fPnmzYsIHLL7+c6667DoCpU6eyY8cOXnzxRW644QZsNhsvvPBCyamCXnLuG/6XEHW8ha4QAMDBgweZOHEiH3/8MQA7duzgxx9/ZMmSJXTv3p1rrrmG//73v4wdO5bevXt795WHKfceVqJqDyptSnoKOmLpTIIPq1ZAA+sxJ21K+ly0UGUrpXqh6zOdEmRfG4FbRcRL5IhH30WFMgephJKy4wEx4hSV3I9+H9+ID5G68cYbGTx4ME2bNsXtdtO2bVveeustFi9ezK5du2jbVhtl9+/fn6VLl3oaliqCTx6vb/yELl2SANCgQQOKi4u59tpreeGFF2jXrh3NmzdHKcXzzz/P6aefzgUXXOAraq8AvddnxCmMRM0MKm1Kegv0cl17qhdWXSBu1987X73+u+KsXZcTnKloNno55kURKZcEai2x3UFoTTwd6KW630RqpwihoUZQ6CXoK6mgZMknn3zCxx9/zJlnnsmYMWPo0aMHAwYMYNOmTdx1111ccMEFniK1neAjT+sTLYEf0QIeC7oQ4ZQpU0hISODqq6+mZ8+egL5B6N+/P1dd5TNyfC/aviqYsjCGKhIVAmXNnH5E+8NV28JGXMUUZ+9h1zu3IoWOQE3dwKvAf0QkYGVMpbgOHexQVcrlIJmQ74pRiibopdlYtP9eyc9dIvxcm2OrBAp4HbiIACJ15513Ehsby4MPPgjAv//9b1q0aMFFF11Es2bNPO2Q8tBWSM+Gd9iRg1I0A/oB34hUWNixOTpnbADW73vDhg3MmTOH9evX06tXL8aNG8eQIUN4++23GTFihPf5+ehquUtD+yoM3kTLEt809N2gX3HKXjaTrMXvAdDmmpeJbeo/6lrZY7A3ak6TkTeQOe8Jr6Mx6Cjv3r/AbR/C4B9EgirbnEHlBMrkIIWG1ngYf3rwO0SNQAk6RyoROAs/y9cdO3Zk5syZACxYsIBt27YxcuRIWrQoE3OTB7wHkVHyviawgpAuA1LQHpIfVrDkvQ8YgQ6K+i+Q2L17dzV58mRWrVrFE088wV9//cUNN9zgS5zy0MJvxKkGiPgZlBUQMYcAe04iwt8vX40rey8gNDruHBqPqDih211UwL5PHqdg834Obx8dcRAGfQmDN1grdrtEeKWivpSiB/oO2B8mBykM+AmxB9gkEnX7A3ZgNnAaft7vt912Gz/++CMJCQkMHz6cyy67rNQSKS8vj3Xr1n1+/PHHj6WeFDW0ltevRC/dlZABzAgyirUX8Kl1funv3Ol0lvNBRIeU/4auIWX2fGuAiBYoK5R8C1o5/JK/+Sf2zrqX5KNPoWDzT4jbRbub3kHZK14NLM4q5O+X1wB2Jxy1GE5fAQmebz4BHq3oza4UyeiIrJJzTA5SDaAUaXhEZXmQIcJbNTuakBCL/sIcih8H9LVr19K+fXtiY2NLy8bn5eXx9NNP85///MclIheIyMc1N+TaQSli0BZGaT4O7wXeD/JGMAGdaHsZgYOv8oDewOYAbQwhJNKX+IYQRFXa3F914nyDPqdhT2xIzsqPcfyxnOQeQyq8gC0phoROR2wq2HLqp9DSl2+eQi8v/hGoHxHylOIr9AfD5CDVHNFQbqMyONHLfPOBE/AhUkcfXbbYcV5eHg8//DBTp04FPQubqZQ6W0Q+D/9wawelsKFdM9L8NGmBTtz/NYjuCoDr0TZR09Ei5V3bLQ9tymvEqQaJ9JLvN1NBZJMr7yCOTSuJadKWhHY9SD5a25Dkrp4f1AVUrF1Sh3Ys9iNOJQScwZUgwlIRNhlxqlGK/DwfTPHISKUIOB1YDYHfS3l5edx9990l4lRCLPCRUqqcJ09dwMo1HAX0DNBskUhQ4uTJF0A3YBmUWfHIA2aizX4NNUikC9QgKgjdzl3zNbiLaWAJU1zzNOJaHUFBxlqcB3dWeAGlULFN8ioKx21TwXFD7VHXZlAlFKAjxX7Dj0i5XK7C22+/veDZZ30G68UD85RSg8M3xFpjIHBcgOM/oYtoVoW96ACKCcBC4BV0nto11JN9vUgiYgXK8tZrHKiNiJD765egbCT3Ohxto2dRQu7q4DwzVYwr0d4o3zO3yoFOjp2PfoO+V8nhG2qOujiDKiEP/WW5CQ83bguH3W6/5qWXXjoF/O5vJgKfK6WOD+MYaxSl6EvgxPoNQLp3RKxSxCvFuUqRGsRlBO0+fzJ66e8dqDB03RAGIlag0OvHAZOUCjJ+pThrFwlpfYlpWFrAjeSew8AeQ+7abxBXxcE24rIVJ3TcX4QWpBfR/nczRVghwq4g8ioMtUddnUGVkINeSfgWHUVWhE4OPR94T0SWoq26/C0FNgDmK6X61cBYw4pSdEW7yPhjO/Chn8/raHSAwySlONrHcUMEEskCFV9Rg5IZUsnyXgn2xIYkHXEcbkcWjk1BVFlXCG7VBhiJrv90nlIMVopuSpFqrXkbIpNifC+9xFgb6XWBHPQXbAp6878ler8EABH5Fl2LyN9sMhX4yrL2ikqUoh36s+nvb7oPmC5S/oZFKXpx2JcwHv35PteyKDNEMJEcxee9pFEGlyMbxyZdfmX/p4+z/9PHfbbLXb2A5O6BSzopQIrtLvSbv7n18Pww7wde8HGqoZYRQZTCie8lvVgqeB9FGX5TFURkvlLqQnTOoK/PdVPga6XUEBEJGJEaaShFU3Sum79ZcQ46pDzfx7kp6PIj3vQG2ivFKyaoKXKJZIHaQoCchNy134CrmLhWRxDXorPPNo4/f6Bg62qcWbuJTQ3gnWlzxxbuSg1UOyormAFbpaULzZJgjVOEb4GKo24JVEBEZK5S6lJgBr5nGi2BhZZIRUW4tIdLhL/vggK0OJXzxPMIRffn3bnFiFNkE7ECtXXqmINpU9IP4qd2UknuU5ORk4hv081nH/bF75G9bCa5v35J46ET/F5Liu35rpzEQG/UgGXXPTgXSFOKvdY5JY+9IoH30wzVoq7vQwWNiMxSSiUAb+M7ArYt8I1SaqiIbKvRwVUSyyXiUvAb2FCMdozY6+f4IPyniGSi95wNEUzECpTFEvTaepkPWkHGGooP/E1s8zS/4gTQoPdIspfNIm/t16QOvhRlK29cLoI4DyRX9EH19wHwpiX6S7Gt9ShFKQ5RVrT2oA1hg6rCawhIXY7kqzQi8q4lUv4sutI4LFIV52LUApZLxEX4Lxsi6ICIDD/ntwWG+TnXbZ3r731jiBAiXaCeQwcuNPB88lCpc8TIgCfHpLYkIa0vBVt/IX/TDyR1O6l8I6HQsaH1R0Cu93U8qHAGZS3vBaoD1dB6HOHxnEsp9lNeuHKNcWylMDMoL0TkVatKtD/T2CPQe1LDRCTYG7AawQpKOhsdyeuPz/0V8rSCH87Df0DFQhEiUpgNZakTXnzVJAPotHXqGLH89FpyOFKqJXpz+fGKCgIqRUe0aWUocGAtDVo/1/vaADZolGIC4Gsj8j0R/qrp8UQSSqnbgUcDNFkDDBeRAzU0pIBY4nQa2ubJH4tFWBigj7PxX014C/CuuQGMDiI6DHfr1DGC9r8Ky/6NuClAV9cV0H56Imy28p/mivAq8EiQ1WpbVtwkaJLQd4/Ho/M+qlOgsT5gZlB+EJHH0MU2/dEbWKCUSqmZEVXI0QQWp1/QOWE+UYqj8C9O+cDHRpyih4gWKICtU8d8gTZxDGm0jbhUcdHulD8zHhmzP2C74N/McfjfC6kOhQQRRagULZVigFJ0sDaX6xNmDyowDxB4FnUs2nHC3xJ3TfI7+F66Qxs2z/P3mbRCys8I0Pc8U+Ymuoj0PagSJqE/RGmEoqKuG7e7MCZn38fHfAaMUQq7CEFk9AboU1iiFEvREUctPR4t0MuEVU323RukSB6JtmYBQCmyKL+3daCOhsCbGVQARESUUlPQ1kc3+2l2EvCpUmqMiNTacrIIxUoxG52YPMDj0A5gtr/3bxAh5T+L8HtIB2sIO1EhUFunjslOm5I+CB3V155qLHmJSxW7C2Nydr9/0puuQ4klOTKjlMImwrLqjNMSkoPWY0PJ80oRi07+ben1CFR7poRgQ9y9lxhTrYdnmGOxRwh8aSh8LdWqikHP4J1U34TTzKAqwBKpf6A/O9f6aTYc+FgpdZaI1Fr+mAhupfgcOIT2ItyPH5cIDwZiQsrrHBEdJOFN2pT0FHT5d79lsQMhLuUs2p2yYd8n/dM9xMmThSIsru44g8HaDG5AedFqhq7pU0K6CD8G0d+NaBGsCrn4DoEPV9XQ0cAHaOueQrQZ6krgR/Sm/TprTEGhFCcDvly7a+zvGS0opWzAm8DlAZp9CpwvIoEEoTIcjRaQRmin8W+CPVEpeqMLfmYFaNMWuArfWxZu4A0R/q7MgA2RQVQJVAlWGfhp6KWzZAIvnwnaIiazaE/DKbveGtKOwDWmvkPXkqmVX4xS2NEiVRJJuCZAImLJOTHAXYR2T9GNvnNdIsKaEPZ7FTr02dcNhgM9o0oCDqCFahm6LtIadLG4cks8SjEYj+VND5aI8HVIRl2HUErZ0TcIFwVoNhu4RESqe5NyHvAu+jMai57tTke7hFc7B1Ap4qy+mvhp8rUIS6p7HUPtEJUCBaUh6IPRa+qD0G/QfLQgKfR6+wH0suBzwPdWKHlz9N1joA3hmf5yLCIRpWgNXBem7j+uQuG3QPwO9KhEexf6BsOGXrLbAqwCfkCL1lrrDvo4tLgVefz8M5r+jjWJUioWmIXON/LHe8AVIlLVfctT0QFO3lWBHejP5dlQvfQJpTgL8OfUvhUdUl4X913rBVErUN6kTUlPRefCxKOXjTZvnTomy1dbpWiGFqmGPg6vBuZGUyiqVeOmH4eXCQPW0aokL4uwu4LrJ6ILunkvEx70+j3Go5fuQrH3mY8WoUSrz9/Ry4Qb0DOtXeiZWDJ63+UgWtCM95qFlcj7Cbo6rT9eA66Tyn9R9AcW42e1QoR8pdiAru1UpRwspeiJdjj3RT76vVvOo88QPdQZgaosStEELVKe+R9r0TOGqL7jsjLpW1A24bgllQ8ucQP/rUaSchEewRgzZtD0oouYaRmAluJyuVi4cCG7d+/mlFNOoXXr1pUcpk8E/SVVso9iA/4GhqBLMxgApVQi8Bk6GMEfzwO3VEKkuqBnuamBGolQqBS70H+T7UH2DZSGlE/C/3t6tgi/VaZPQ+QRFVF84UCEA0rxFnAF+oP0O3VAnEB/8NEf+NIPvRWU0YjyQRlN8b93FWyghL8k5TignfUgI4M+RUXEx3tV4XnggQfYtGkT2dnZzJo1i6effpojjjiifG+VQ1F+n6vki/MI/Iem1ytEJF8pdSawAB3I4IvJQL5S6s4gRKoFeuYUyPYL0DdSIrRTip/REYTrKjH0JPRs2JdA/WLEqW4Q8Ym64cSKDHoLvZ/hrxJnnUAEESFbhD9E+F6EOSK8CPwXbSr6CbAcvTxWEkFX1RB3n/TtS+u4uLKh33v37mXBggXcfffdpKenExMTw6pVq8qc53brP8uuXbv46quvcLmqvLcei17+7F7VDuoiIpIHnA4Bo0VvB+6toKsGaJeH5vj4bvGlbUoRI0JTdDDMkCCHjAi7gJfRqx6eHMCjmKMhuqnXAgVgfWl/UV9dxUUotsrarxZhgQjvivAE8ATBhwMHJVDdu9NOecVbLlu2jO7du9Ojh46bOPXUU/nrr78oLj48cbPZbLhcLm6//Xbi4uKw28u70lcCN9ChOh3URUQkB+2BtzpAs3uVUnf6ORYLfI7eBy6TIL1161YyMjJQ3n98C2t23xCdq3Re8GOmAPjIehRiXMrrHPVeoAy+ESE3UO5JCdaXi8+aXd60bEkzz/+73W527NhBmzZtsNlsFBQUoJTi0KFDxMQcXn0uKCjglltu4dRTT6VRowpXjkr7drlcvu7aFYHTDOotInIQXT0gkOPCVCvh1xOFDls/Bq8ltxtuuIFbbrmF448/nlde8Vf9o5REdOSgrwq4fsaMWGkQLwOfmHynukW93YOqClbOxUC0m3K9nHF5Y5Vcf5ryARkt8HBy6NmT5JiYsu83ESEnJ4cGDXTE/8GDB9m/fz8tWx6ekB08eJA5c+aQlJTEaaedxqZNmwKOZ8+ePaxYsYIdO3ZQVFREQkICrVq14uSTTy4RNztGoPwiIvuUUiej95G6+mn2tFIqX0RKFOdJ9BJhmT2///73v2zcuJFvvvmGFStWcPnll9OhQwdGjx4daAiJ6Lwpf3lNfsZd6uBiqEMYgQoSy67oErQfYCulmB1Gp4WowioFkmE9gNKZVWMswRo3jpOLiymKiTmcE2O328nIyCiN2lu3bh379+/nxBNPLO373XffZf78+XzxxRcsX7484Djmz5/PrbfeSvfu3enWrRsJCQls2bKFWbNmsXLlSh544AHs9tiYvXtp3ro1trq851gdRGS3h0il+Wn2slKqQESaoK2Tyoh+bm4uDoeDJ598EoATTjiBSZMmsWdPUNuaKehZWf0MMTaUYgQqCCxxGsfhD2s34GKlmFmBP1i9xcp/OmA91qPzYsotKZ999tm88cYbvP3223z00UeccsopnHLKKaXHk5OT2bBhAz179iQnJ4f77ruPQYMG+bzmrbfeyrx583xGALZr146bb76ZFi1ax375JacCsUqxD6/cLZHwlHaJNkRku1JqBFqk2vlqc+GFF75VXFzsjImJKed52KBBA2688cbS2THov+WyZcu44oorcDqd7Nmzh3btfHYt6JmuuQGs5xiBqgCP0tPeBfGOAMYpxQwjUkFxHDpRtwwjR45kz549LFiwgJNPPpmbb76ZV155hZYtW3L66adz9dVXc/XVVwNw4YUX0rRpU78XaNCgAfv376d169a4XC7cbjeHDh0iIyODLl26AGC3Y0tOJg69kd/GepSiFIc4LFjfWxvx9RIR2WLNpL7Dq/T6sGHDePPNN5W3OK1bt45t27ZRUFDAueeeC+j9QJvNRqtWrSgs1BaYZ5xxBueff37p39YDB/AiRpwM1ONE3WBRiuHA0ABNtqKdlk3kUGC+BYYF0/DAgQM4HA7atWtXGlJut9sZPHgw7777Lp06+a4E/uGHH/L4449zzDHH0KVLF4qLi8nNzeWHH35g/PjxjBs3jpiYGL78kpWnnVZhKHJQScr1AaXUUcAitEckvXv3ZunSpWVmRwAbNmxg3LhxDBo0iIyMDNq0acNLL72EzaYnzrt27eLuu+/G5XIRExPDG2+84X2pAnTKx6lYeWpKMRAdvv6N+VvUP4xAVYDH8p6vkuIlbAfetxJkDb75L3ALVXChL+HPP/+kdevWJCf7j3HIzs7myy+/ZOPGjdhsNlq0aMEZZ5xRJvBi2TJ+HTiQTyq43F4RXqpoTErRHV39uGTWtbcuzqiVUn2Bbzt27Jj6888/k5qaWio8oPecRo0axfXXX8/48eP5888/efDBB3niiSdo3lyb7O/evZs2bdowevRo0tPTvS9RhHa1Px7tvYhStAGuRi8N7wHmiBgXkPqEEagg8Fjm8xfVBNpG530rYMBQnni0U0AaYVpafuaZZxg5ciQ9e/Ysd6xkmQlg9Wo29OvHzAq6WyvChxVdUynOQIdXl1Cy9+btS5gVTf6OvhgzZswpr7zyypetWrVSnmkAAJs3b+aTTz7hn//8JwBOp5MxY8Zw1VVXcdFFh03Tp02bxjXXXIPX+S5gN9pPch+URsxeh3Y6KaEY7XixKtp/l4bgMHtQQWBV+ZwJXEDZAoCetAUmKMV7ZqPdJ4Vov7fJwDno4nKF6C/0ZEKQk3f88cf73aPyvNtPTg6q0m5VXTQU+ku1KeCplEVKlROtvVG0x5Wcnp7+lNvtdtlstnLfG507d2b8+PGICG63m9jYWI4//vjSJdp58+bRp08fJk2a5H2qADmU90g8jbLiBPr7agw6qXdhaF6WIZIxM6hKYNVqOp/A5SL2oC3+a6NKbTQRgw406Q30RZccPwodYuxA51B5l2kIiszMTA4cOEBubi6FhYU0btyYbt0O31eI8IPNxih8+xKW2B18IELApCsrlH4K1avcm0352VZmhIXAxwBfAifiw/vOc3YKOr9NKcUzzzxDamoqnTt3ZvLkyXzyySd06tSJ7Hwn2w84KCx2E2tX+dn5zlGDuzYvLSypFD3wX6uqAO1SnhXKF2iITIxAVRJLpM4BegVotg94RyT4qrCGUlLRv9ve6Mi/Y9FC5kIv8SRTtuJwGT7//HM++OADMjIyEBFSUlJo0aIFw4YN48wzz6RJkyaglxqP9j7X2m9sjharjRXNhC1H/Jur8BorohjtilAZ89Rw8hy60GSZ/cPNmzfTpEkTUlNTS58rESeAGTNmcPfdd9OwUSMmP/Asq3KS+XHrQXIKnMTH2FGAw+nKc7klDp1ku8R5IOn1na8N640ofzcncyLo92IIM0agqoBS2NBl5/sEaJaJFqmcmhlVncaGXhLsjXYsmIiP5ek333yTWbNmcdNNNzFs2DAaNGhAfn4+q1ev5qGHHuLss8/mmmuuAW2I26W6g6rgTr+6vCVyOPHZz/VtgC3M0W2dgN/wms3++9//ZsWKFWRlZXHeeefRt29fRo3SZaVcLhd2u52ffvqJUVf+k/bn3k6B2MkvcgXcOBJBcKsit9NekL2062eHfuz8p1eTX0X4OKSvzhDRmD2oKiCCWynmokOR/VXzbApcoRTvmKJp1caNrqS7BfgTuBgf5Ry2bNnCcccdx9ixh63cEhMTOfHEEznqqKPYvHlzydNVjiT0YgfwIWWXCYMzC6yYYPbAWgLXKEUm5ZcJc0IUSJCAV07SmjVrmD17Nr///jvLly9n9erVfPTRR2RmZnLppZdit9vJKXDyxgah8Rm3c6AQgqnurhQKu8Tb7cXxjYdtuDC5x86N+z4+5jPXocRC9Azr8xC8HkMUYQSqilgi9Sn6k3esn2ZNgCstkTI+YaEhDz9LfF27duXjjz9m4cKFtG3blpycHHJzc1m/fj379+/nwgtLi69WtnCjT0Q4hC73UFrywaou7L231QKCCswoITvI4ImW6Nllc+vhuexc4Ccoo7L5ell47bHFxsbSrl07bDYbAwcOpEuXLqSmprJixQo6d+7MEb36MfKx+eRJHEVV3ElTdomNa5nTvfWV37fZ/cGJrxdnNnzdpHHUP8wSXzWxNspHofM3/JGNXu6rUmlrQxmaA9vwIzJffPEFDz30EHv27KFly5bY7Xbatm3L9ddfz9ChpfnW+YRuFlUh1lJciS+hp6muP0PUP0SYHkS/I9HBJZXBOwR+L3AgwGxLAdPRS9qly3wXX3wxrVq14pFHHiEhIYE9e/bw/PPP4ySG5Y0Gs+NgPu7SeJOqI27c4rLvtMW6em2dOsasRNQzjECFAEukTiXwl8UhtEjtr5lR1VmS0IIf1Ow/Ly8PEfF2PXATINCiplCKeA4HZXg+fhSpuBaXUlxGCPbS0K4NezksWttF2OlxPA6YhX6PJwH88ccfTJs2jaSkJP7xj3/QvHlzdu7cyWn3/4/C5t0pcoXue0WEAqX4eOvUMZeErFNDVGAEKkRYIjUCGBygWS7wqgmcqBYKvaxa7vY8NzeX7du307VrV+x2O0uXLuWdd97B6XRyxBFHMH78eNLS0kALVCJEnj2V9T6yBxP4oBT/QtsAhRpfSco2dDTflVgi9csvvzBr1izWrVvHP//5T15LX8Zy+1GIrfxqZsYjXiWe7DHY4pKIadScuJZdSOp2Egmd+qFsfu8bHMD5W6eOMdVy6xFGoEKI9eUyFP+ec6uAdJMFX20K8ZF79PXXX3P//ffz/fffs3btWm677TbGjh1Lnz59mDVrFsXFxUybNg20MLUiiusHWbOvyYRHoL4WYYmfY7cD92Et9x04cICZM2ey+Pvv+aXdORTYfa+clghUysBx+glx4y7Io2j/Ngr//h1cxcS16kqzM/9FbJO2/saVAXTaOnWM+fzUE0yQRAixhGeRUriAk70O/4IRp1DhU6AaNWpUaqFTXFxMXFwcN910E6DLbYwbN66kaUk+VdQKlBUw8IRSNKB8scjmVO+zHSiC8DFgpwivKkVikyZNmDRpEr1POY9r3v8ZigJH66UOvrTcc668gxz46hUcG5aw53/30Pryp7Enp/o6vSl6hWKxr4OGuocRqDAgwveWSI20nloDzDPiFDIK0HY3ZWjRogUNGjTgu+++Iz4+HpvNxrx582jSpAmLFy9m+PDhJU1d1JGqulYyeC46twsoDcpoSvm9rZQguy0VKKXoCOzyiv57f8MG8rt0YUZcnI5OfPeHbeRXIE7+sCc3ptlZt7PHkU3htrVkL59Fk1Ou9dU0GZ0YbQSqnmAEKkyIsMwSqfZoV4BIsq6Jdnwa8qalpfHYY49x3XXXkZiYSHJyMjfddBMDBgygY8eO3H///SVN3dQRgfKF9V7bZz1KXReUIoHykYQtKTsbzUcH9JSI02XAbqWY7ums0bMnv59/Pm+9+y4T4uOJX7X1oKrO3ZdSNlJOuoi929aS9/tiGp98TakjhWczwHe1SkOdxAhUGBHhB6VYaWZOIcevBVGPHj1YvHgxu3btIiMjg5SUFNq0aUNKSpnJQ/Xjn6MQK7cqw3oApfumqRwWK5sIohQt0WVmYtAVda+yjJCzrFNbzpnDrj17ePX9WUVXZuc7y81oK0tCu6PAZsftyKI4ew+xqa18NWuSNiU9devUMVm+DhrqFtV2kDYExohTWPgr0EERoXXr1pxwwgn06NGDlJQU3O4yE9g4CGwjVF8QQUQ4KMIGEb4T4VulSAEupWyuWVPgaqVKK+u2BPj+ew6eeNaB2XZsVVvf80DFxGJL1DrndvhNeconcG02Qx3CCFSEoBQJ1pKKoWLeJMAsysfSkKfbtgNdZygzHAOLdiwnjPH4tmxqAFxq1UfbhN5b3ZPtKHYVFElo/ABLbud8/A09WsSH5FqGiMcs8UUAVsjweKC1UswU4Y/aHlOE8xHaOuhtdEQfbjf24mLilULZbNhsNmwuFy63m+K4OA6h77w3A68D79fWwCMZy839EnQUoC/cwFwrR+tn60GH/1t3LHb3xVRTOKS4CHfBIQDsiX7jORQYy6P6ghGoWsaqHHopep0f4CKlmCPC+locVjQwE/gGXaHX9uOP8I9/cElWFs7MTIoOHKDYpc2zc0V4olZHGgVYkX/no4N6/DFXBG+HcWxxrr8Igb9hwY7fwe3ClpxKTKp3HchSEvGIWDTUbYxA1SLWHes4oIPH03bgAqX4UITfamdkUcN+68EJJxCH9rZzohNxS35GS8XaWsMKlBiD/2rRoJN3f/V1YOvUMQfTpqQfREcHVgkRN9nLZgKQ3HNooKYHTIBE/cEIVO1yHrrejjc24HylsIuwpobHFJVYeTrza3scUcpQ4JgAx38AllbQxxJ0Ic9KR0i68rI48NXLFG5bi71Rc1JOvNBfU7GuY6gnGIGqXVaizT59lWJQwDlKYRNhdY2OylBvUIpj8W/NBbpY4fwgolGfQyemB7Reyvr+A/0PEdyFeTj3Z+ilPVcxca2PpNmZ/8Ke5Hf/Kc+6jqGeYASqFhFhs1J8gN6YLmfdgxaps62Z1E81OzpDXUcpuqOX9vyxFfg4yFSJxVJclKti4gIKVPbSGfofHmaxDXqNIKnbQG0WqwIGFmcC3wcxFkMdwZjFRgBK0QEdKBEoCupzEVbW0JAMdRzrPTcB/zepe9Bl54Paw1NKnZXY5diZzc66M94WF5J6kN4YN/N6iMmDigBE2Aa8R+AN/dOV4kQfzyv05nQbwuNsbahjKEUL9KzdnzhlAe9XQpxuBj7O/2tVvGPTD7iLQ17FpACYa8Sp/mEEKkIQYQfwLn585ixOU6qMF1lTtHHmdnTi5E505VODwSdK0Qidc+dvmuNAi9OhivtSdqXU08CzWMERB758EVfOPsTlDNWQnej396RQdWiIHswSX4RhWclMIHBJ8m9FWIeOrmpPebPPm9BuCwZDKZZLxJX4Dwd3oqs+76i4L5WETng+x/uYLSmFVuMfw96oObYYX1urQVOAFqdBW6eO2VudjgzRiZlBRRgi7EY7JOT5a9OxI6Nzc/lFpJw4gU5kfAG4PmyDNEQdVs7dxfgXJzcwO0hxagEsxIc4gfbR2/XOrQWunH3fE8CSqgIcaMeQAUac6i9mBhWhKEUz4HK86h61a0f8jz9yVbNmNImJwW99bPQHfAr1OCzXSkCNA1zBlFCvq1guERcAPQI0myvCLxX3pboBnxPYsHUfcKaIrEibkj4amIZejk4mcJ6UoG/MMoFJZs/JYAQqglGKpmiRagTQqhVxq1YxsUULmsXGBhSnEhzAvVA/rH6U4hTgaHReWRyHgwA+qc+5ZJZQjwSfQTYAC0UqLgKolBoMzAUaB2j2B3C6iJQ6zqdNSVfoSrg3o+s5NUEvRQtasBKBA+gk3OeA701ZdwMYgYp4lKIxcHnz5jRftYorWremRWys7+irvDy9KpicXKYWnwOYCjwU7rHWNkoxBhjg45AJ0QeU4iQOV3ku4Uf07yfgF4FS6mLgHXzn65WwBDhbRAI6xadNSU9Fz8Di0cavm419kcEXRqCigGuuofW997K8ZUva+hKn3NxcPvnkE1599VUSExMZP348l112mWcTB/Ak8J+aGnNtoBQjgZN8HPpKpEKrnnqBUvRGR3ragfXofSe/1Z6Vrl1yB/omJxD/A64UEeN9aAgZxkki8kl47TVmidDSqsNThj179vDuu++ydOlS7r77btq2bcull17KcccdR7dupd6fScD/oe9+p0CdLaLoLwGnWqFkdQkR1ihFHnAC8GEF4hQDvAhcW0G3jwB3i4jfvgyGqmAEKvJ5EzhGqfJ5KwUFBUyfPp2//vqL++67j759+wLQoEEDHI5ywVNJ6PDzeOCf1E2R8pd848vrsN4iwl9UUJVYKdUQmAWMCtDMBdwgIq+GcHgGQykmzDyy6Qacjd5ELsfChQtZsGABV199NX379qWwsJAZM2bQu3fvUrHyIhm4Bh2GXmnX6SjAzKBCgFKqLToBPJA45QJnGHEyhBMjUJFNPH5mBS6Xi3fffZdLL72UY489lsLCQr7//nt++OEH+vTp47PsuUUyOjLwNere39/MoKqJUqo3sALoG6DZTmCwiJgwcENYMUt8kU0efu7+lVIkJiaSm5sLwAcffMBvv/1G06ZNmThxIgButxubzYaIeAtWMjppMw64AvzvQ0QZ9XYGZS0Bj0IHhPhN8g7chxoJzMEr986Ldegw8u1VuYbBUBlMFF/k8wZaTMpZH23YsIFx48bRqFEj0tLSGDFiBBdccAFJSeVdkoqKinA4HKSmpno+nQd8hnZSd4Vl9DWIUnRFvxZv/hTh/ZoeT01hBc+MB9LQ+UTviXCwcn2oicArBL5p/Qq4QESyqzhUg6FSGIGKfGKBD4GT8SFSWVlZ5Obm0qxZMxISDsdRuFwusrKy+PTTT7HZbMydO5dFixbx6KOPcs0113h24UBXor0IotttQSk6or3mvNkmUje9CS2XiPOBnh5P5wIfiLCr4vOVAh4A7qmg6VvAdSISMhdYg6Ei6toeRF3ECZwLfIEPX7PU1FTatWvHq6++yrfffgvAwYMHGThwIC+++CIrV66kT58+OBwORo0axYAB5fJYk9BLQx8T/Xs19WoPynKIGEVZcQJdduVKpQLaEaGUikc76FckTv8GrjLiZKhpjEBFB8XoGc6n+DHfnDhxIs2bNwcgOzub9evX0717d6ZNm8aiRYtITEzkjjvuoE+fPhQWFrJu3TrP05OAEcA8onu/pr7tQQ0EjvNzLI4A3ntKqcbomfP4AP07gctE5CExSy2GWsAs8UUXNvRSy3noQAe/bNy4kXPOOYcuXbqQmJjInXfeSZ8+fbDb7WRlZXHFFVfQq1cvHnqojANSPvAl2qU66t4YSpEC3OrjUI4IT9X0eMKJUvRFpyD4YwMwy1cirlIqDW34Gsg8Ngs4R0QWVXWMBkN1MTOo6MKN3mOZQYByHADdunXj9NNPZ/78+XTo0IH+/ftjt2t/2dTUVN566y2WL1/Ovffe63laInAKcHp4hh926sUMygoGOTNAk+34cYlQSh2LDiMPJE4ZwEAjTobaxghU9OFGW8+8TQCReuedd1i1ahUzZsygW7duuFwuRIS//vqLpUuX0rhxY9566y327t2L01lua8FfzaBIp87vQSlFW+BC/H929wHTRcr/LpRSZwLfAS0DXGIVcIKI/F7dsRoM1cXkQUUnAkxGzxiuw0d03yWXXELv3r3p168fDoejdPZks9k444wzuP322/nf//7HUUcdRWxsnfn+dqEF3PvL264UdpHoDqW3yq9cin/BzUGXa88vf666CV2aPdBN6TxgnIhUKY/KYAg1ZgYVvQjaAPY5fAROxMbG0q9fP7Kzs5k3b17p8506deKll14iMzOTiRMn8sEHH3j36QA+CevIw4RVMqJOzqKUogFwGT5uRiwK0OJUJkdJKWVTSj0JPE/gz/sL6D0nI06GiMEIVHQjaHfyJ/AT3ZeSksJnn33G9dcfrgD/7bffkpKSws0336w7ORwok4uuF1SpJM8Io87tQylFPDraLtXz+eHDabppE1fs28dNDgcXijC47HkqEZiNNgf2h1jHbxaRqJ5hGuoeJoqv7nA3cBd+7rDPO+88GjRoQGZmJjExMUybNo3WrVuXHi8qonD1akYddxyLama44UEpbkZXbPXmeRECFtKLRCyXiEvwKrH+f/9H1//+lwtiYoix2UqNfx1osXlFKdUcnZZwQoDuC4DxIvJhGIZuMFQbswdVd3gYXZ30fnyI1IcffsiMGTPIy8vj/PPPL2N5VFiIc+JE3ps+nf7A2mj8IvegzsygrETcs/ESp5deYsDVVzPSR/HKJODpr776yg3c6X2eF/uBM0VkeQiHbDCEFDODqntMRheQ87dXUYaiIpzXXcf7b7/NNuupQ8C7IuwL1wDDiVJcBbT3cehNkdLXGPFY4nQaHjMgux312WeMHDGCY+Li/O+p5eXlMWjQIFavXu2vySZgtIgErAllMNQ2Zg+q7vE8epmnXCSXN0VFOCdPZrqHOIF2sr5CqagNNa8rM6gT8RCnlBRifv6ZiysSJ4DExES+/vpr2rf3pdMsBU404mSIBswMqu5yJbpct89ihy4XhbfcwqwXX2Szn/Md6JnU7nANMBxYhrEJ6Gi+Io+fh0SiwwxXKXqj/RcB6NKFxO++Y0KLFjTzXtZzuVylKQSeFBcXs23bNvr160dOTk7J07OAy0WkIIzDNxhChhGous144FXKi1Q+cJ5SHEI7R/gjH126YWeYxmfwwsp1uhFrdWPIEJp8+ilXJCeTFBNDGSV688032bp1K3fddVcZJ/sSCgoK+Pnnnxk2bBhOp/MxYIqI1JXaX4Z6gFniq9u8j44AywWyrcchtNfeFyIsARYEOD8RuFwp2oV7oAaNFaDyJcD48bRbsIBrGzakgbc4Pfvsszz77LOMHDnSpzgBJCQk0KdPH1atWrVKRO404mSINswMqn7QFl2SwQWsh7J1gpTiOAL77y0Sie7w82hj4UL+MWgQj/nbb7rpppu46qqr6NevH7t27cLtdpOamkpysk8P4TzgaXTZDIMhajBh5vWDv62HT0RYqRQuYCygvA4vQ/u3GWoGBfxrxAjux4/7xaFDh9iwYQNut5u1a9dy7bXX0rRpUzp37szo0aMZPXq09ynJaNeRr4DFYR29wRBCjEAZABDhJ0ukzuKwSP0AfGVZCBnCjx2Yhl6WLbNvuH//fsaNG8fLL79Mly5duPDCC3nggQew2+28/PLLNGrUiAULFrBgwQIGDhxIw4YN0cVyS0lER3cagTJEDWYPylCKCKvRlXUF7Wo934hTjZGM3g+8BK9aX7/99hvnnnsuW7duZc6cOYgIp512Gj169GDNmjV07tyZTp06cfzxx7Nz506SkpK8xamEbF9PGgyRihEoQxlEWAO8DqQbcaoxWgErgZPwEqfvvvuOSy65hOuuu44nnniCTZs2oZSiY8eOnH/++fTq1avUU/Gzzz7D5XJRWFjo6xr5wMYwvw6DIaSYIAmDoRZQijgRitCFAxcBjfHac3K73TzxxBOceOKJDB48mIMHD9K3b18effRRLr74YtxuN5mZmdxxxx3YbDa2b9/ORx995CtQIh/4DRiEtsMyGKICI1CGaqEUNl+VW2sLpbCj91vi0F/4JT+LRNhRm2MrQSlOBI5bvZotffrwPtAAr+AUEUEphdvtxmaz4XQ6iY2N5fXXX2fDhg3cf//9ZYQoPz+fxESfOdl5wBJ04q9Px3uDIVIxQRKGKqMUScAEpVhh7V9FAh2Ay308nwG8VcNjKYdS9AJOe+ABevXsyeP4idTz3kMqKSrZtm1bZs2ahdut7wmKi4uJiYnxJ04O4F20P6MppWGIOswelKFKKEUiMAG9f3KWUvSv5SGVELEFC5WiM3DO5ZfT/o47ODM2tuyYnE4nH3zwAd988w1Op34ZNpvNs14Xo0ePJj4+nv/85z8AxMT4vcd0APcAN2DEyRClGIEyVBqlSEBXd21V8hRwplIMqL1RlRKRZrFK0Rq4CLCffz7dvRNwMzMz6dmzJ9988w1Tpkzh3nvv5bPPPrPO1bMpl0vrzNVXX01BQQHZ2X6D8hzo0vBPh+XFGAw1hFniM1QKj+qubXwcHqMUdhFW1PCwPIm4GZRSNEYLRjxAs2Y09G6zfv16Tj75ZF5++WW2bt3K119/zdy5c0lKSmLEiBEApaawvXr1on///qSkpHh340JbWZ2Gjgo0GKIaM4MyVJZUoFmA46OUYmANjcUXETWDUopktKA3KHlu0yb2FBWVFVKn08nChQvJzMwkLS2NUaNG0adPH+bNm8e+fbo01/r163G5XHTp0sVXKY1CYCdwDEacDHUEI1CGSiHCHuAdAtebOlUphtTQkLyJmBmUUsShE2+bej4/aRIrDh4kq7j48N7Q8OHDOeuss3jyySdxu920a9eOQYMGsW/fPvbt28fnn3/O4sWLfZbWQP8t1gJ9wW/5FIMh6jACZag0IuwC3iZw2PIIpRhuVYatSfwJVIxSNfd+t8LdL0Ab9ZYhLw/XqafyrtvtPlQSjQcwbtw4Dh06xJNPPglA3759SU5OZuXKlZx++ulcd911vi6VB3yBznE6EIaXYjDUGkagDFXCmkm9jS7l4Y+haKGqMZGy3C9qdRZlvd4zga7+2qxdO3PF8ccfcygvL6/0uX79+nH++eezevVqJk2axK+//srSpUtp0qSJv24cwEtoITQJuIY6h0nUNVQLpWiGzjsqt/HvwXLgy5qyTlKK2/CyDLJ4UoRDNXD9U9AzGj+8kQVX3wM0HTZsGJ999llp0q2IkJmZyb/+9S9sNhv9+vVj8uTJvjrJB/6BLkhpMNRJjEAZqo1SNEGLVLmwMg9Wooskhv0NpxT/QAdzePOcSHiXwZTieKBcvYvDvNgAJt8EEl/yzGWXXca0adPKOEOICC6Xy1eek6BnTucD80M4dIMh4jBLfIZqY33pvwVkBWh2HDoMvSaW+2olkk8pjgJG+T4qwFPd4eb/8xQngPfee48nnniC3NzDq6VKKV/iVIzeZxqIESdDPcAIlCEkiJCFFqlAM5Rj0Qm94X7f1fgelFJ0Qvvd+RDgYhs8MgymXIQf28L77ruPX3/9db2I+As8KUTbNfUFfg3BkA2GiMcIlCFkiJCNDpzYH6BZP+DsMItUjc6glKIVcDG64KAXuXHwyHnw8FD/w8IN3DBw4MA+SqnVlA94cKDrc/WHyDC8NRhqAiNQhpAiQg5apPYFaNYbONcKxQ4HNTaDUopUPFwiyrK7ITw9AZ7qqaPBfZIHnCki09DjHoXer3OgCwyWGL4OB3JCPHyDIaIxVkeGkCNCrlK8jTaTbemnWS/ArhRzREJuZlojMyjLzX08PiMY/2gB/7sEXk6Bg/662A2MEZGfPZ47hBajI4AmaJH6PYTDNhiiBjODMoQFEfLQjhO7AjTrYT1CTU3NoNLwconQ/NgJ/jcR3kkJ8PJ/A473EqcSXOjqt8sx4mSoxxiBMoQNERxokfrbT5PvRFgXhkvXyAxKhN+B2ZQpZ/FlH/h8PMyJD+A6tBAYJCLbQjkeg6GuYQTKEFZEKEDvoWz3OrQEXeo8HNTYHpQlUu+BuwBmD4NlZ8MCm7bG88k7wGgRyQr1WAyGuoYRKEPYEaEQeA8dJg166eqbMCbt1nAelNoJXfvBiqH6pS331/A+4EoR8RvOZzAYDmOCJAw1gghFSvE+Oo9nVZgdJWowik+lAh/C5hHwCn6i9YqBq0XknVBf32CoyxiBMtQYIjiBH2vgUjUUxac6AunAUfoZn+KUDZwnIt+E8toGQ33ACJShLuI9g3J6PEKCUuoY4DMOl733xTbgdBH5LVTXNRjqE8Ys1hCRKEUD4HQg3QpZr8y5sejE2SLAWZ3lRKVoAbhFDrtjKKXGAjOBpACn/gyMFZFAYfYGgyEAJkjCEHFYZdIvB3oCV1hiFTQiOEXIFaGomuKUgk7EnagU7fRzahIwl8DilA4MNeJkMFQPM4MyRBSWO8PllHWg2A+8a9ko1dQ4EoGJQHP9jKsYhqfB91dVcOo04GYRKQ7rAA2GeoARKEPEYInCBKC1j8MHgHcsQ9pwjyPWGkd7/Ux+DLx7DvzdE+YBq/2dejvwhJgPlcEQEswSnyGS6IdvcQLtS3elZc4aNiyX9fMpFafMJHh1AuzqqT8uZwEtvE8rBC4UkceNOBkMocMIlCGSWA78EOB4KlqkmoTj4lYxxTFAN/3M1ibwxlVwsP3hVl8Dez1PywROFpHZ4RiTwVCfMUt8hojCEolTgZMCNDuEXu4LVHeqKtceBgzT/1vdHtLHgTPxcIsf8Cpk+xfatmhTKMdhMBg0RqAMEYclUiOAwQGa5aJFKlDdqcpc8xjgDP2/RT1h8bng9qhX9RvwIRwOClwOnCUiIbm+wWAojxEoQ0RiidRQSmc0PslDR/ft8dcgbUp6Y6ATOi+qENiydeqYMgWalKI7cBGIgrknwepTy/ayFXgfD9PyOcAEEckP/hUZDIbKYgTKENEoxWDg5ABN8tEitQsgbUq6Aoa4CmKm2GJdA7BJI4RCUHryY5M4pTiIdlN/btvjo7aIyz4BiuNg+ijYPKBs93uAt/Cowv4EcIeIuEP5Og0GQ3mMQBkiHqU4CRgZoEkB8F7HO9N7o/OQmoqQbM3C/CEiOKQwxpn1fccFh35a2B/2di3bJBt4Hb2aiBuYLCIvVeOlGAyGSmAEyhAVKMXxwGhfx+wN8+Obn7vq9LhWOV2tXKpK4Xa6xfHHQXXgywyksGQZLx94Ax2khwO4SEQ+q+LwDQZDFTACZYgalGIAOgy8lJgmucmtxi+baIsvbqTsUmXzY3exG1dOEbvfX4/bkY+usbgDYDfaU++n6ozdYDBUHpMHZYgaRPgR+BQrlM7eMD++1fhlE20JztTqiBOALcZGTEocrcZ3R8XPwxKn34ETjDgZDLWDKbdhiCpE+FkpXMDZzc/5aawtoThF2creaGU8MrbCflqO+y8JHXuXeU7ZbdgbxdBk5Glkzlv7LXCuKc1uMNQeRqAMUYcIv7Y4/8eT4lrmdFM2sftrlzJwnN8+YlJa+nzeFhNLcreBrqQjjnt621MXZFV7sAaDocqYPShD1GGFkm8BOvo6XjKD6nhntWIaMoBOW6eOMR8Qg6GWMHtQhmhkCNA0zNdoSmAnC4PBEGaMQBmikZuB5DBfI9m6jsFgqCXMHpQhGhkEAZNwAcj6/gOfz6uYOFJOvKCi05V1HYPBUEuYPShDVGF56+0BYv21qSiKT8Un0+HWmcFczgm02Dp1TFYlhmgwGEKEmUEZoo1OaGeHlIoaVjNIArSdRGfg5+p2ZDAYKo/ZgzJEG/E1eC2p4esZDAYPjEAZoo3CipuEDFXD1zMYDB4YgTJEG1uApBq6ViKwuYauZTAYvDACZYgqrGKDBytsGBoOmAAJg6H2MEEShmhkCXAOFYSa+wszB0g68kTiWnYOdLpY1zEYDLWEEShDNPIcuoBhg0CNspfO8HssJqVlRQKVZ13HYDDUEkagDNHIYnQlQZ8CFYLwcqz+vw9FRwaDoWqYPShD1GEZuE5C50OFAwcwyRjFGgy1ixEoQ1SydeqYL4C5QEGIuy4A5lr9GwyGWsQIlCGamQRsR1sShQKn1d+kEPVnMBiqgREoQ9SydeqYbLSh61aqP5MqsPoZZPVrMBhqGSNQhqhm69Qxe4EBwMdUfU/KAXwEDLD6MxgMEYBxMzfUGdKmpI8GpqGLDSYTOE9K0KHkmeiACLPnZDBEGEagDHUKqxz8YHSxwUFAE7QruaAFKxE4gE7CfQ743kTrGQyRiREoQ50mbUp6KrpkRjza+HWzsS8yGKIDI1AGg8FgiEhMkITBYDAYIhIjUAaDwWCISIxAGQwGgyEiMQJlMBgMhojECJTBYDAYIhIjUAaDwWCISIxAGQwGgyEiMQJlMBgMhojECJTBYDAYIhIjUAaDwWCISIxAGQwGgyEiMQJlMBgMhojECJTBYDAYIhIjUAaDwWCISIxAGQwGgyEiMQJlMBgMhojECJTBYDAYIhIjUAaDwWCISIxAGQwGgyEiMQJlMBgMhojk/wHRMiA+nUizxwAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "execution_count": 47, "metadata": { }, "output_type": "execute_result" } ], "source": [ "import networkx as nx\n", "def plot_network(mr):\n", " G = nx.from_numpy_matrix(mr, parallel_edges=True)\n", " G = nx.relabel_nodes(G, {0:'A',1:'B',2:'C',3:'D',4:'E'})\n", "\n", " elarge = [(u, v) for (u, v, d) in G.edges(data=True) if d[\"weight\"] > 0.5]\n", " esmall = [(u, v) for (u, v, d) in G.edges(data=True) if d[\"weight\"] <= 0.5]\n", "\n", " pos = nx.spring_layout(G)\n", " nx.draw_networkx_nodes(G, pos, node_size=700)\n", " nx.draw_networkx_edges(G, pos, edgelist=elarge, width=6)\n", " nx.draw_networkx_edges(\n", " G, pos, edgelist=esmall, width=6, alpha=0.5, edge_color=\"b\", style=\"dashed\"\n", " )\n", "\n", " # node labels\n", " nx.draw_networkx_labels(G, pos, font_size=20, font_family=\"sans-serif\")\n", " # edge weight labels\n", " edge_labels = nx.get_edge_attributes(G, \"weight\")\n", " nx.draw_networkx_edge_labels(G, pos, edge_labels)\n", "\n", " ax = plt.gca()\n", " ax.margins(0.08)\n", " plt.axis(\"off\")\n", " plt.tight_layout()\n", "plot_network(m)" ] }, { "cell_type": "code", "execution_count": 48, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAFtCAYAAAAK6G3eAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/av/WaAAAACXBIWXMAAAsTAAALEwEAmpwYAAB/vElEQVR4nOzdd3xdZf3A8c9z7l7ZeydN070Xm5ZZNggiQxAVcYALB6Ii4sb1c6GIiIqyFFmyZynQQmlLd5s0TbP3zs3Nnef5/XFv0qRNaTpuk5v7vHkdzj37e5+ePN975iOklCiKoijxSxvvABRFUZTxpRKBoihKnFOJQFEUJc6pRKAoihLnVCJQFEWJc8bxDmCy2rBhQ4bRaLwfmI1KuIoyGenAtmAweOOiRYtaxzuYo6ESQZQYjcb7s7KyZqSnp3dpmqbu0VWUSUbXddHW1jazubn5fuDi8Y7naKhfqtEzOz09vVclAUWZnDRNk+np6T2Ej/pjmkoE0aOpJKAok1vkbzzm69GY/wKKoijK0VGJQFEUJc6pRKAoihLnVCKY5JYtW1b25JNPJgB86UtfyvnEJz6RP94xjQdVDmGqHFQZjEbdPnocfOPxzfkVzX32Y7nOsiyX5xdXzKs71Hzf//73G++8887c5uZm45YtW+yvvvpq5bGM43C89Kff5LfX1RzTckjLL/Sc+/mvxFQ5qP1B7QsTjToimOTOO+88t5SSP/zhD5lPPvlkldFoZMeOHeYrr7yycOXKlSXD573gggtKvve972WOV6zRNFo5/POf/0y66qqrCi+44IKSJ554ImFw3ngrh40bN1qvueaagpUrV5bcfffd6YPzTtZyGK0MAHp7e7XZs2fPeOSRRxIH552sZbA/dURwHIzll1q0rFu3ztbW1mZKSkoKJicn6wAzZ870//vf/64Zngj+9a9/JZ1//vndr732WsLB13Z0xvJrLVpGK4frrruu+7rrrutua2sz3HzzzXkf+chHeo9HOUy0/WHhwoXehx9+uDYUCvGRj3ykGGiLdjlMtH0B4Hvf+17WpZde2jk4fDz2hYlCHRFMYjU1NaaPf/zjxU888USl3W4PPf7446Pu0B6PR/znP/9Jvvnmmzt7e3sNxzvOaDtUOXz729/O/tKXvtQWz+Xw0EMPJa5YsWLqeeed1zOZy+FgZfDkk08mzJw5cyAjIyMIk/9vYn8qEUxSfX192iWXXDLl7rvvrl+4cKH3jjvuaPrRj36UM9q8d955Z1Z/f792zTXXFOzevdvmdrvF8Y43Wj6sHHRd5/Of/3zuBRdc0HPKKad44rUcAK699tqe1atX73700UdTJms5fFgZvP766653333X+dhjj6Xcf//9aXfcccekLIODEaqpyujYvHlz9bx589rHO47RNDc3G2699dbct956K+HKK6/sqKurM//73/+uAfja176WfcEFF/SeccYZ/eMdZ7T96Ec/ynjkkUdS582b15+TkxOora2Ny3J49tlnXY8//niS3+/XsrOz/fG6PwD87ne/S9U0jdWrVzvHWgabN29OmzdvXtFxDfQYU4kgSiZyIlAU5diZDIlAnRpSFEWJcyoRKIqixDmVCBRFUeKcSgSKoihxTiUCRVGUOKcSgaIoSpxTiUBRFCXOqUSgKIoS51QiUBRFiXMqEUxyqhGOMFUOYaocVBmMRr2G+nh46uZ8Wncc00Y4yJjp4dJ7YqoRjs7HK/IDzf3HtBxMWQ5PyhVlMVUOan9Q+8JEo44IJrmxNkwze/bsGddff33BiSeeWLZ+/XrreMYcDWNtmCYey2G0hmkmczmMtWGayVwG+1NHBMfDGH6pRctYGqaprKw0zZ8/v//BBx+sveuuuzJqamrMixcv9h7rWMbyay1axtIwzdy5cweORzlMtP1h/4ZpKisru6NdDhNtX4CRDdMcr7+JiUIdEUxiY22YZu3atY6qqirrNddcU7B69WrX5Zdf3nu8Y42msTZME8/lMLxhmslcDmNtmGYyl8FoVCKYpA6nYZr169fbf/WrX9U9/PDDtUajkd7e3kmzXxxOwzTxWg4wsmGayVoOh9MwzWQtg4NR7RFEyURuj2B4wzQf//jH29esWePMy8vza5pGcnJy8He/+13jeMd4PAxvmGb+/PkDzz33XGI8lsPwhmnmzJnjefHFF+OyHCDcME16enrw3nvvTR9rGUyG9ghUIoiSiZwIFEU5diZDIpjUhzuKoijKoalEoCiKEudUIlAURYlzKhEoiqLEOZUIFEVR4pxKBIqiKHFOJQJFUZQ4pxKBoihKnFOJQFEUJc6pRDDJqUY4wlQ5hKlyUGUwGvUa6uPgjnfuyK/sqjymjXCUJpd6fnjyD2OqEY6nnnoqv7W19ZiWQ0ZGhufSSy+NqXJQ+4PaFyYadUQwyY21YZrf/va3qZ/4xCfyr7/++oLPf/7zueMZczSMtWGaeCyH0RqmmczlMNaGaSZzGexPHREcB2P5pRYtY2mY5tVXX3Vs3brV9o9//KMOwOv1imjEMpZfa9EyloZpEhISQsejHCba/rB/wzSvvvqqJ9rlMNH2BRjZMM3x+puYKNQRwSQ21oZpHnjggbTbbrutdXDYarVOqlfSjrVhmnguh+EN00zmchhrwzSTuQxGoxLBJHU4DdN4vV5hMpmGdvRgMHj8Ao2yw2mYJl7LAUY2TDNZy+FwGqbxeDzaZCyDg1HtEUTJRG6PYP+GaS6//PLuu+66Kzs1NTXodrsN9957b11aWlpovOOMtv0bpjnjjDP64rEc9m+Y5uyzz47LcoB9DdNMnTrVN9YymAztEahEECUTOREoinLsTIZEoE4NKYqixDmVCBRFUeKcSgSKoihxTiUCRVGUOKcSgaIoSpxTiUBRFCXOqUSgKIoS51QiUBRFiXMqEUxy6t3rYaocwlQ5qDIYjXr76HHQ+O3v5Pt27z6m7163TJ3qyfnJj2Pq3es7dt6W3++uOKbl4HCWeWbOuDumykHtD2pfmGjUEcEkN9b2CGbPnj3j+uuvL7jsssuKfvOb36SOZ8zRMNb2COKxHEZrj2Ayl8NY2yOYzGWwP3VEcByM5ZdatIylPYLKykrTwoUL+x988MFaCB8633LLLR2DfyDHylh+rUXLWNojmDt37sDxKIeJtj/s3x5BZWVld7TLYaLtCzCyPYLj9TcxUagjgklsrO0RrF271rFo0aL+wWG73a5r2uTZNcbaHkE8l8Pw9ggmczmMtT2CyVwGo5m83yzOHU57BOvXr7efcMIJHoC1a9facnNz/ZNlpz+c9gjitRxgZHsEk7UcDqc9gslaBgejXkMdJRP5NdT7t0ewZs0aZ35+vt9isUiDwSDvvvvuxuGHzJPV/u0RPPfcc4nxWA77t0fw4osvxmU5wL72CO699970sZbBZHgNtUoEUTKRE4GiKMfOZEgEk/dYR1EURRkTlQgURVHinEoEiqIocU4lAkVRlDinEoGiKEqcU4lAURQlzqlEoCiKEudUIlAURYlzKhEoiqLEOZUIJjnVCEeYKocwVQ6qDEYzOd+pOsG89uDO/M4G9zFthCMl1+k58/oZMdUIx1d21ubv6vce03KY7rB6fjOjIKbKQe0Pal+YaNQRwSQ31oZpZs2aNeOaa64puOaaawrefPPNY/oHOhGMtWGaeCyH0RqmmczlMNaGaSZzGexPHREcB2P5pRYtY22YZt68ef0PP/xwbTRjGcuvtWgZa8M0x6McJtr+MFrDNNEuh4m2L8CBDdMcj31holBHBJPY4TRMs3v3bts111xTcOutt47aZkEsO5yGaeK1HPZvmGaylsPhNEwzWctgNOo11FEy3q+h7uvr00499dSyO++8s/Gyyy7rfeGFF5y333573qZNm3YNzrNy5cqSF198serLX/5yzpVXXtl18sknD4xXvNHyYeWg6zo333xz7rnnntt76aWX9sVrOQyfb/ny5aXz5s3zTMZy+LAy+OIXv5jb39+vlZeXW61Wqz579uyBsZbBZHgNtUoEUTLeieDDjNYwTVZWVsBgMMiVK1f2Xn/99d3jHePxMFrDNPFYDqM1TBOP5QAjG6YZaxmoRKAc1EROBIqiHDuTIRGoawSKoihxTiUCRVGUOKcSgaIoSpxTiSB6dF3XxXgHoShK9ET+xvVDzjjBqUQQPdva2toSVTJQlMlJ13XR1taWCGwb71iOlnqyOEqCweCNzc3N9zc3N89GJVxFmYx0YFswGLxxvAM5Wur2UUVRlDinfqkqiqLEOZUIFEVR4pxKBIqiKHFOJQJFUZQ4pxKBoihKnFOJQFEUJc6pRKAoihLnVCJQFEWJcyoRKIqixDmVCBRFUeKcSgSKoihxTiUCRVGUOKcSgaIoSpxTiUBRFCXOxVx7BGlpabKoqOiIlu3v78fhcBzbgKIoluKNpVghtuKNpVghtuKNpVjh6OLdsGFDu5QyfdSJUsqY6hYtWiSP1BtvvHHEy46HWIo3lmKVMrbijaVYpYyteGMpVimPLl5gvTxIvapODSmKosQ5lQgURVHinEoEiqIocU4lAkVRlDinEoGiKEqcU4lAURQlzqlEoCiKEudi7oEyRVHih5Q6UoZGdBBClyEYGhc8YB6JDlJHIiN9HaTcb1hHyh10dBgAHSl1YP955EH6kXmlHFpWokeC1vfbvty33hHDDH2OLLjf/Ow3v0RKM7D8mJezSgSKEqd0PYiue9F1L6GQb9jnAXTdR0j3ooci43QfUg+gy8Cwvh9dBtD1AHKwrwfQ9QY2b3lk5Lih5fyR+YdV5Ayv0PUR/X2VZPRs2hz1TRwzgvOjsl6VCBQlRoRCPoLBHgLBHkLBPoLBfkKhfoIhN6Ghzx5CkeFgqH9fP+SOVPBeQpHKPVzRHjkhjAhhQtNMQ31NmJEE8Hr7RkwzGlxowoTQTAhhDH8WhmGdEYSGFumH1z3YNyAwILRIf7/lBucbsRwaQhgAgRAaIMLTh/rhcZs2bWbBwoXh4cHpQ8tokT5Dywyta8Q8w9Y9OHzAvGKw1IYNi2HDjJi+/zyDy69ateqo/s0ORiUCRTnOpJSEQm78/nZ8/nYC/g78/nYCgS4CwV6CgXBlP9gP6a28sSpceR+KECYMBgdGgx2D0YHB4MRocGCxZGAw2DFoVjSDFU2zYtAswz5b0TQLBoMVTbOhGSLTNWtknAVNM4+o3Icqyf2sWrWKZUuXH+NSiw4hPCQlLhrvMMadSgSKcoxIKQkEuvD5mvB6m/D5mvH6mvH72/D72yNdB4FAB7ruH3UdBoMDkzERoykRozEBu70YT38meXnTIuOTMBkTMBpdGAwODMZwRW8wODAaHWia5Th/a2UyUIlAUcZISonf34ZnoIYBTw0D3lq83kZ83ia8vnDFr+u+EcsIYcRsTsNsTsVsTsPpKIsMhztTZLzZlIrJlISmmQ7Y7qpVq5hauvw4fUslHqlEoCj7CQS6cffvxtNfGa70IxW/Z6AWXR8Ymk8IA2ZzBlZrNi7XLNLTzsJizcZqycFqzcZiycJsTjvoKRRFmShUIlDiVijkwe0ux91fQX//bvrdu3H3V+D3tw7NI4QZm60Au62A5JSTsNkKsdsKsdkKsVpzRv0FryixRiUCJS6EQl7c7p309m2lr3crIf09Vr3ZBJF7vzXNhsNRSmrKKTgcU4c6qzU7cveJokxeKhEok5LX20h393q6u9fR07uJ/v6KyMNIYDKlArkUF12OyzULp7MMqzVPncJR4pZKBMqk4PHU0NW1lu6e9+nufh+vtwEAg8FJYuIC0lJXkJAwB5drDhZLFm+++SYlJcvHN2hFmSCilgiEEA8AFwKtUsrZo0y/FriN8BMTfcDnpZQx9IyfMp6CwX66ut+ls+MtOjrfZGCgFgj/2k9KWkJ+/idJTlqK0zldndpRlEOI5hHB34E/AA8eZPpe4HQpZZcQ4jzgPmBZFONRYtzAQB1tbS/T3v463T0bkDKAptlITj6B/LwbSEk5Bbu9ZNhTnIqijEXUEoGUcrUQouhDpq8ZNvgukBetWJTY1d9fSWvri7S1vUyfezsADkcZ+fk3kJpyKklJi9VDVIpylCbKNYJPAy+MdxDKxOD1NtHc8gzNzU/S378bgMSEBZSW3k5G+jnYbAXjHKGiTC4i/BrVKK08fETw7GjXCIbNswL4I3CKlLLjIPPcBNwEkJmZuejRRx89onjcbjdOp/OIlh0PsRTv0cYqpRfJRsIHirsIv3WyFCGWIViAEMnHKlQgvsr2eIuleGMpVji6eFesWLFBSrl41IlSyqh1QBGw7UOmzwX2AGVjXeeiRYvkkXrjjTeOeNnxEEvxHmmsfX3lcueu78k3Vs2Rr75WIt9+53S5Z89vZH//3mMa3/7ioWzHSyzFG0uxSnl08QLr5UHq1XE7NSSEKACeAK6TUlaMVxzK8afrftraXqG+4SG6u99D08xkZFxAbs5VJCYuUhd7FeU4i+bto48QbkonTQhRD9wJmACklPcC3wNSgT9G/vCD8mCHLcqkEAy6aWh8lLraB/D5W7Ba8ymdchvZ2VdgNqeMd3iKEreiedfQ1YeYfiNwY7S2r0wcfn8HdfX/oL7+nwSDvSQnn8T06T8mNfV09TSvokwAE+WuIWUS8vnaqK75E42Nj6HrPtLTz6Gw8LMkJswb79AURRlGJQLlmAsEeqip/Qt1dX9HSj9ZmZdSWHgTDkfpeIemKMooVCJQjplQyEdd/d+pqbmXYLCXzMyLKCn+CnZ70XiHpijKh1CJQDlq4VvQNvLue9/H660jNXUFU0q+hss1Y7xDUxRlDFQiUI6K211Bxe4foss1GAxTWTD/QVJSTh7vsBQlpkgpkYQfo5QSdOTQ56HxSIJRev5XJQLliIRCXvZW/57a2r9gMDgR4hqWLrkTTVO71GQT1CU+qePTJX5d4tN1/LrEL8OfQxICuiQkJQEp+UAa8bR1E5QMjQtKSVCP9KU8yDSGTZfoEkJE+pGKUpeSUKSvRyrM8PR9n/X9lglJiQ7haZHPg+Pc0oX13Z3oyKF1jKh8R62Uw7Xx/pX0yOHB6Qcuu38FfzguxspZR/nvORr1V6sctq6u99i569sMDFSTnX0FpVNuY82aLSoJHEdBXeIOhegNhugP6QzoOp6QzkBIx6NH+iGdAV0O+7yvPxCZd7AyD1fwkc8yXOH79XDlrx92dE7YVn1YSxgFmITAIAQmIdCEwChAEwKNcN8gQEOgDesPjhMCDEPjwssIAUYEFi28XjHKejr6e8l02YaW0QQIwvMKQXg9kfUPPuYYnhaZZ2h48HN4PKMsK0ZZllGW1YYP77ddrSo6z96qv1xlzILBPnZX/pTGxsewWQvUaaCjIKWkNxiiKxiiMxCkMxCiKxCkKxCkJxjCHdTpjVT0tdLB3esrhir+vkjFfjismsCmadgNGjaDhk0L981C4DQZsWgaZi1caVq08Hjz4GdNYBbhzxbtwPEmTWAU+7rNGzeybPGicKW+3zRjpIIfPmwQjNvT5KtW1bB8VtG4bPtIrNq7MyrrVYlAGZPung1s3/41vN4GCgpupKT4KxgMtvEOa0LxhnRa/QFa/UFa/QFa/EFafQFa/QE6AkG6Avsq/e5gkNCHnBdwGDQSjAachvDvylyjgTyrGZdRw2U04DIYSDBqOI0GHAYN+7BK3h6p6O2R8VaDhuE4VrRuEWK2y37ctqccPZUIlA+l6wH2Vv+B6uo/YrXmsmjhIyQlxd+bQPy6TpMvQJ3XT73XT703QL3XT4PPT4svXPF3B0MHLKcBqWYjaSYjySYjZQ4rKSYjKSYjyUYDySYjKSZDeNhkJNlkIMFoGFFxr1q1iuXz46/MleNHJQLloAYG6tm2/cv09m4iO+sjlJV9D6PRNd5hRY0npFM94GOPx0eVx8eeAS/VA37qvH6afYEDLuxlmo3kWs1MdVg4KdlJptlIptlEhsVEptlIhtlEmtl4XH+NK8qRUIlAGVV7xyq2b78V0Jk9+/dkZpw/3iEdMz2BIDulgar6Nnb1e6ny+Kga8NHoC4yYL9tiothm4bRkF3lWE7lWM/kWM3lWMzlWExZNvSdJmRxUIlBGkDLE3r1/YG/173E6pzNn9j3Y7YXjHdYRkVKyd8DPlj4PO9wD7Oj3stM9QIMvALhgdwOJRgNT7BZOSnIyxW6hxG5his1Csc2Cw6gavVfig0oEypBAoIdt279MZ+dbZGddzrRpP8BgsI53WGPW5g/wQa9nqNvU5xk6b28UUGq3sizJyQyHlWDVbq46aSlZZpNq/0CJeyoRKAB4PNVs3nIjAwMNTJ/2I3JyrprwFWSd18/abvdQVz3gB8IXaGc4rVyYnsSCBDvzE+xMtVswDzuVs2rvTrIt5nGKXFEmFpUIFLq632fLls8hhGDBggdJTloy3iGNqt0f5I3OXt7q6mNtdz913nDFn2Q0cEKSg+tz0liYYGeOy4bDoE7rKMpYqUQQ55qan2Lnztux2XKZN/f+CfWmUF1KNvcN8FpHL6919LKpz4MEUkwGTkxy8tn8dE6MnOrRJvjRi6JMZCoRxCkpJTU197Kn6pckJS1j7pw/YjIljXdY+HWd1V1unm3t5pWOXjoCQQSwMMHON4qzODM1gTlOm6r4FeUYUokgDkkpqdxzN7W1fyEz82JmzrgbTRu/8+V+XefNzj7+19bNS+299ARDuAwaZ6UmcFZqAstTEkg1q11VUaJF/XXFGSlD7Nr1XRqb/k1u7seZVnbnuLQbLCOnfR5r7uTJli66gyESjBor0xK5KD2J01Jc6j59RTlOVCKII7ruZ/v2W2lte4GiopspKf7qcb8zqM0f4PHmLh5r7mRXvxerJjgvLZHLs1I4Ldk54s4eRVGOD5UI4oSu+9m67Yu0t7/K1NJvU1Dw6eO6/Y09/dzf0M7/WrsJSMmiBDu/mJbHxelJJJrUbqgo4ylqf4FCiAeAC4FWKeXsUaYL4LfA+YAHuEFKuTFa8cQzXQ+wbftXaG9/lWlld5GX9/Hjsl2frvNMazd/rW9nU58Hl0HjhtxUrstJo8wROw+qKcpkF82fYn8H/gA8eJDp5wFTI90y4E+RvnIM6XqQ7Ttupa3tJcqm3nFckoA7GOIfjR3cW9dKmz/IVLuFn5bl8dHMZJzqtQ2KMuFELRFIKVcLIYo+ZJZLgAellBJ4VwiRJITIllI2RSumeCNliB07v0Fr6/NMLf02+fk3RHV7PYEgf21o5y91bXQFQ5ye7OILMzI4Ldk54Z9SVpR4JsL1cJRWHk4Ezx7k1NCzwM+klG9Hhl8DbpNSrh9l3puAmwAyMzMXPfroo0cUj9vtxul0HtGy4+Fo4pVSIuU/kLyFEJejiei9PbRfCp7wC143uxhAsIgAl+Jlqjjw/fwTRSztC7EUK8RWvLEUKxxdvCtWrNggpRy9YYtwhRGdDigCth1k2rPAKcOGXwMWH2qdixYtkkfqjTfeOOJlx8PRxFu551fy1ddKZOWeXx27gPYzEAzJP9a0yGmrt8is1zfKT2+tklt7+6O2vWMplvaFWIpVytiKN5ZilfLo4gXWy4PUq+N5u0YDkD9sOC8yTjlK9fX/orr6HnKyr6Sk+KvHfP26lDze0sXdVU00+AKsSHFxbmcjN8xecMy3pShK9I3nTdvPANeLsBOAHqmuDxy1ltYXKK/4PmlpZzFt2g+P+bn5Tb0eLtiwmy/trCXVbOTx+VN4ZN4UiibwaSBFUT5cNG8ffQRYDqQJIeqBOwETgJTyXuB5wreOVhK+ffST0YolXnR1vcf27beSmLiQ2bN+g6Ydu3/eDn+Qn+1t4l+NHaSZjfxuRgFXZCard/4oyiQQzbuGrj7EdAncHK3txxuPp5otW7+AzVbAvLn3YTDYjsl6pZT8p6WLO3c30BsKcVNeOl8vzsKlbgNVlElDPdI5CQSDfWzechMA8+bed8zeItrg9fON8jpe7+xjSYKDn0/LY4bz2CQYRVEmDpUIYpyuB9m67YsMDNSwYP6Dx6R9YSkl/2zs4Ad7GglJ+NHUXD6Vm6ZOAynKJKUSQYyrrPwpnZ1vMX36T0hOPvoHs9v8Ab6ys47XOns5JcnJr6bnU2izHINIFUWZqFQiiGENjY9RV/938vM/RW7Ox456fW909PKlXbX0BkP8OHIUoJ4IVpTJTyWCGNXbu4Xy8u+TknIqU0u/dVTr8us6P65q4s91bUxzWPn3vCnqWoCixBGVCGJQINDN1m23YDGnMXvW/yHEkd/B0+Tzc+O2ajb0evhUbhp3TMnBZlBtAihKPFGJIMZIqbN9x634fG0sXvQYJlPyEa9rbbebm7ZX4wnp/GVWERdlJB27QBVFiRkqEcSY6up76Oh4k2llPyAhYe4RrUNKyf317Xx/TwNFVguPzy9lmmofQFHilkoEMaSj822q9v6WrMxLyc295ojWEdAl36yo45GmTlamJfC7GYUkqIfDFCWuqUQQI3z+drZvvxWHo5Tp04/sHUI9gSA3bq/mrS43Xy3M5BvFWerZAEVRVCKIBVLq7NzxDUIhN7Nn/RODwX7Y66gd8HHtliqqB/z8bkYBV2alRCFSRVFikUoEMaCu/h90dK5mWtldOJ3TDnv5Tb0ePr6lioCUPDqvhJOTXVGIUlGUWKUSwQTX17eDysqfk5Z2Frm51x728mu63Fy/tYokk4En5paqRuMVRTmASgQTmJQ+tm3/CiZTEjOm//Swrwu80t7DZ7ZXk2818+/5U8i2mKMUqaIosUwlgglMysfweKpYMP8fmM2Hd07/qZYubtlZw0ynjUfmTiHVrP6pFUUZnXqEdILq6FiN5E0K8j9FSsrJh7XsI00dfH5HDYsTHDw+v1QlAUVRPpSqISagQKCXnbtuB7IpKfnaYS377+ZObt1Vx+nJLh6YU4xdvS5CmUCkrqPrOjLS6UP90Ij+iGmh0NAwQw2u6yDDd9TJSH9wGvs1zP5h47uqdrPbZj5wPQC6jozEPBi7RCJ1Cfv1I0tEtjH0bQ8YlnLYfEOjh42HodhGDEfm6x3wH9t/kAiVCCag3bt/hN/fhiZux2AY+yugn2jp4is7azkl2cnf5hSrdwbFAV0PEfIHCAb8BP1+gn4fwUCAoN8XHj9sONz3EwoECAUD6MEgoVAIPRiI9IPooSChYDA8LRhED4XC84ZCkeHB6SPHhyvvELqu4/P52Pbgn0at9IdXcBNF1UtPj3cIY5a5YGlU1qsSwQTT3v46Tc3/pajw89TVFY95uadbu7hlRw3Lkhz8Y06JSgITkJSSgM+L1+3G5+nH53bj9w6Eu4EBAl4vAe8Afp83Mhwe19zYQONrzw5ND/h8kUrfjx4KHnVcBqMRzWAM9yOdwWBAM5rC/aFpBgxGIyaLNTLPvnmFZkBoGppBo6mpmbz8fISmhccN9g0GhNg3vG+aAc0Q6R8wLbLc4HihIYQY6tivP/p4jfBHDSL9wWkbNmxgyZIl4fkAoWmAQGhiqC/2Gx45Pvx3Nri+QUM3dojwfIjh48XgpKFlhuYZMTxsvsgyq1evPup/79GoRDCBBAJd7Nz1bZyOaRQXf5G6urVjWu6Ftm6+sKOGJYkO/jWnRJ0OOg5CwSADvT14It1AT/fQZ6+7D19/P97+SIXf7x6q/PVQ6JDrNphMmK02TFYbZqsVPeDHbE3BkZiM2WrFZLViMJkxms0YI/3BYZPZjME8+rR9wyYMRiMGoylSuR7bp8tXrVrF8uXLj+k6o8VeU0d64dh/cI23wcRzrKlEMIFUVPyQQKCL+fMeQNPGdkpoTZebz+2oYZ7LzkNzS3Co9wYdlaDfj7uzg77Odtwd7fR1duCOdJ7ebjy9vQz0dOPtd4+6vGYwYHE4sTqcWBwOrA4niRlZWB2OA8Zb7E4sdjsmqw2T1YrZZsNksWIwjvyzjKWKNR6Fz+kPHzH0v33jR1wnGOO8kRmHzy+O/gBwVCoRTBDtHatobnma4qIv4nLNHNMy290DfGJrFQVWM/+aW4JTJYFDCni99LQ20723kg393XS3NNHb1kpfpNL39vUesIzZZseZkoojMYn0wmLsCYnhLjERe0IStsTIcEISFodjwrbqJnUJIRk+Zx+MfA7pyJCEUHicDOn7jR82bnC6LkGXSJ3whc7IcHicJKVa0D2wd8Q49HBlKUNy5DLDxh0w72AFKxlZKQ6N39eXw+cburYa+azvdyF32HzFfo3G1WuHvsvQRvbfPqNV9sdfSrGAs479eqOaCIQQK4HfAgbgfinlz/abXgD8A0iKzPMtKeXz0YxpIgoG+ynfdQd2eylFRZ8f0zI1Az6u3rwHp9HAI/OmkGJSOX1QKBiku6WJzoY6Ourr6Gqsp7ulmZ7WZvq7u4bm20O4kk/MyMSVmkZO2XScKWk4U1JxpaThTE3FlZKK2Xb473Y6HFKXyICO9IXQfcFIP4T0h5C+EAl1gr63GpDBUHi+oB7uD++COjIQ2m843BHUIxV4VL/GkGQhcNc2hM+nC4EwRM53awKhEf5siJzDHxynDX4enBfC5+Ij08KDkWsAkXUMP6cuIvMN5uDR5htcjbZvXGdjA7m56ZFz/AfOP3R+fnhuH7YNDvy4b0AMn3/kOoaf9z9w3Qdub3CwvnX3KCV+9KJWe4hws1n3AGcD9cD7QohnpJQ7hs32XeDfUso/CSFmAs8DRdGKaaKq2vt/eH2NLFr42JhOCbX7g1y9uQq/Lnlq4RTyrPH5xLCUkt62Flr27qGtpprO+lo6GuroamoccRHVmZpGcmY2xQsWk5iRRVJmFnvqGznj/AuwOl1H/Qte6hLpDaIP7Nd5Bj8H0D1B5OB4X2hfZe8LIQOhD/2FmYFGz/aqoWFh0sKdURv6jMmAMGpodtOo04VBC1fIkb4wCDAOftbCFbNBA2PkQujg56Hlhn0eXmlrDKu8w+PCp7JOO6oyPV42r6pn9vLS8Q5jzLyrxikRCCEygZ8AOVLK8yIV9olSyr8eYtGlQKWUsiqynkeBS4DhiUACCZHPiUDjYcYf83p7t1BX9w9yc68hKWnxIecfCOl8YmsVTT4//55fynRHfLQtLHWdruZGWvbuoXXvHlr3VtKydw++/n4gfCdIYmYmqXkFTFm0lJTcfFLzCkjJzcNsPbCMmletwuZKOGD80PZCklCfH93tJ+QODOsHCLn96MPH9Qc+/FSBUaDZTGg2I5rNiMFpQqRaEWYDmsWAsBjQLEaERUNYjGjmwXHh/rsb13HyaScjTIZw5TxBTz0psUvIQ9zXK4R4Afgb8B0p5TwhhBH4QEo55xDLXQGslFLeGBm+Dlgmpbxl2DzZwMtAMuAAzpJSbhhlXTcBNwFkZmYuevTRRw/jK+7jdrtxOp1HtGw0SBlElz8E3Gjihwgx8hTE/vFKCb/HzhrM3Eo/S0XgOEd8cMe6bEMBP/0tzfQ3N+BubsDd0ojuDz9MIwwGbKnp2NMysadnYE/LxJaSimY0jW3lErxd/SRqDoxeMHoFpoFwPzwMBl/kFr796AZJyAwhMwQtEDLvG9ZNEDJJdOPg53BfHuWlm4m23x5KLMUbS7HC0cW7YsWKDVLKUX9tjuXUUJqU8t9CiNsBpJRBIcSh74Ebm6uBv0spfyWEOBH4pxBitpRyxNlMKeV9wH0Aixcvlkd6B8VEu/uiuvpe9lTVM3fOn0hPP+eA6fvH+4u9TaypbuE7Jdl8sTDzOEZ6aEdbtgGvl/qd26jZuon6ndtprd4z9ERnal4Bs089g+yp08gsKSUlN/+AO2v2J3VJqNtHsGOAYIf3gD7BkUcDwqRhSLRgSLeE+4nmcN9lRnOaMDjDfc18/C/IT7T99lBiKd5YihWiF+9YEkG/ECKVyMGvEOIEoGcMyzUA+cOG8yLjhvs0sBJASrlWCGEF0oDWMaw/pg0M1LO3+vekp58zahLY3xMtXfyquoWPZaVwS0HGcYgwunQ9RGvVHmq2bqJmywc0VuwkFAxiMJnInjqNpZdcQc60GeRMnYH1Q34BSV2GK/cWD4EWD4GWfgItHoLtAxAadrRr1DCmWjGm2rCWJbO3vY7pS2ZjSApX/JrdqE65KHFrLIngVuAZYIoQ4h0gHbhiDMu9D0wVQhQTTgBXAfs3tFsLnAn8XQgxA7ACbWOMPaZV7P4hQmiUTb3jkPOu7+nnq7tqOSHRwS+m5cVshRXweqnevJE9G95jz8b3h27VTC8qYcF5F1M4Zz6502disozeZoIMhPA39ROod+Ov7yPQ1E+gzQPBfRW+IcWKKdOObXoKxjQbhhRruO8yR54KDetZVYttZmp0v7CixIhDJgIp5UYhxOnANMJ3MZVLKQ95cjpyCukW4CXCt4Y+IKXcLoT4AbBeSvkM8DXgL0KIrxI+4rhBHuqixSTQ3v4G7e2vUjrlm1itOR86b4svwKe27SXbYuKBOcWYo/RkYbQMuPvY/d4a9qx/l5qtmwgFAlgdTooXLqF4wWIKZ8/Dnph0wHJSlwSa+sMV/mDF39I/dAuk5jRhynHiLE3ClOnAlGXHmGEfl1M3ihLrxnLX0M3AQ1LK7ZHhZCHE1VLKPx5q2cgzAc/vN+57wz7vAA7vHcsxLhTyUl5xF3Z7Kfn5n/zQeYMSPrO9mr6gzmMx9KxAwOtlz4b32PnOm1Rv2ogeCpKQnsm8s85jyuITyJ0+84Bz/DKo46/vw1fdi39vD77qXqQvfClKsxsx5TpxTc/HnOvElOfCkGiO2SMjRZloxlKzfEZKec/ggJSySwjxGeCQiUA5UE3Nn/F661iw4F9o2off//8gNtb19PPnWYXMcE7s20SllNRt38LWN16hct1aAj4vzpRUFpx3ETNOPp2M4ikjKm4pJYFmD76KTrwVXfhq+iAY/rlvzLBhn5+OpSgRc4ELQ4pVVfqKEkVjSQQGIYQYPGUTeVAsPp9gOkoeTw01tfeSmXkRKcknfui8jzV18jIWPp+fziUZyccpwsPn6elm+5uvsf3Zp9jY04XF7mDGKcuZfsrp5E2fNeIlWbongHd3N96KLrwVXeh94dtBTVl2nMuysBQnYi5KwOBUu5eiHE9jSQQvAo8JIf4cGf5sZJxyGKSUVOy+CyHMTC29/UPn3dLn4baKOmYR4DslH34NYTxIKWms2MXG55+m8v130UNBnFm5rLj2BspOOHnExd5Qr4+B7R0MbGvHt7cHdBA2I9apSVjLkrFOTcaQOPY2FxRFOfbGkghuI1z5D74E5xXg/qhFNEm1t79CR8ebTC39DhbLwZ8B6A2GuHFbNSkmI1/2dWPUJs4pET0UYve6NWx49imaKsuxOpwsWHkBc85YydbKPcw6fTkQrvw9m9oY2NaOv7YPAGO6Ddfp+Vinp2DOd424g0dRlPE1lruGdOBPkU45AqGQl4rdP8LhKCMv77qDziel5NZdtTT6/Dy1YCruDybGnbQBr5ctr73ExheepretlaTMbM741OeYffpZmKzhX/9i1x76N7bg+aAVX2U3SDDlOEg4uxDbnDRMGdF9cZuiKEduLHcNnQx8HyiMzC8AKaUsiW5ok0dt7f14vQ0sXPAQmnbw1yD8o7GDZ9t6uGNKDosTHaw6fiGOKuDzsvnl53n/f0/g6ekmd/osln/iM0xZtBRNMyClxFfdQ/97zRRv0egKVWBIseJakY99QQamdFX5K0osGMupob8CXwU2AMfq1RJxw+trprrmXtLTzyU5+YSDzretz8OdlQ2ckeLi8/npxzHCAwV8Xja/8gLvP/NfPD3dFMyZz0lXXEPu9HA7CboviPuDRvrfbSLQ7EFYDPRlS0ovnIe5MEHd4aMoMWYsiaBHSvlC1COZpPbs+QUQYmrptw46jzsY4qbtNSQbjfxuRiHaOFWkeijE1tdfZu3jD9Pf3UXB7HmceOvt5E2fBUCgfQD32w14NrYi/SFMOQ6SPlKKfX4Gu9e8xayixHGJW1GUozOWRPCGEOIXwBOAb3CklHJj1KKaJHp6NtHc/BSFhZ/HZis46HzfqqinesDHfxeUkmY+/g+NSSmp3rSBN//1AB31teRMm8mFX76NvJmzAfDX9dH3Zh0D2zvAILDPTcdxQnb4oq/69a8oMW8stc6ySH/460slcMaxD2fyCN8u+iPM5nSKCj930Pmeauni8ZYuvlGUxYlJx/91uB31dbzxj/uo2fIBSZnZXHzrtyldGn7GwVveSd+b9fiqehBWA67T83GenIPBpe7zV5TJZCx3Da04HoFMNi0tz9Db+wEzZtyN0Th6Bd/g9XNbRT2LEux8+Ti/Vjrg8/LuE4+x/n9PYrJaWH79Z5h/7vkYjCa8ld30vlyNv7YPQ4KZxPOLcSzLQrPExisuFEU5PGP6yxZCXADMIvx2UACklD+IVlCxLhTyULnn57hcc8jO+sio8+hS8qWdtQSk5J6Zhcf1eYE9G97j9b/9md62VmadfianXftJ7IlJ+Gp76Xx5F77KbgyJZpIuK8WxKBNhjK0X3SmKcnjGcvvovYAdWEH4QbIrgHVRjium1dTch8/XzOxZv0WI0SvR++raeKfbza+n51NkOz5P1np6e3jtgXupWPsWqXkFfOzOn5E3czaBVg/t/9iOd2cnmsNE4oUlOJdlh9vCVRRl0hvLEcFJUsq5QogtUsq7hBC/AtRdRAcxMNBATe19ZGZceNA2iHe4B/hJVRPnpSVydVbKcYmr4t23efWvf8LX38/JH7uOJRdfjvBLuv+3B/faRoTJQMI5hThPzkWzqFc5K0o8GUsiGIj0PUKIHKADyI5eSLGtcs/dgKC09LZRp3tDOjfvqCHJZOCX0/KjftfN8KOAzJJSVt7xY1LzCulf10zvy9XoA0EcS7NIOLtQvexNUeLUWBLBs0KIJOAXwEbCdwypdw2NoqdnI62tz1Fc9MWDNjjzs71N7Oz38tDcElKjfKto7bYtvPCHX+Lp7R06Cgi1eGn94yYC9W7MxYkkXVSCOSd2Gu9WFOXYG8tdQz+MfPyvEOJZwCqlHEubxXElfLvoTzCbMygsvGnUedb39PPnujY+kZPKmakJo85zLOihEGsff5h3n/w3yVk5XHPbnaTnFtH7ci3ut+vR7CZSrp6GbW66eg5AUZSDJwIhxOi3u4SnIaV8IjohxabW1ufCt4tO/xkGw4Hv2PHpOrfuqiPHYuKOKdF7tXRvWyvP/u7nNFXsYvaKs1lxw03o9V5a/m8DoS4fjiVZJJ5XhGY/+DuPFEWJLx92RHBRpJ8BnAS8HhleAawh/KSxAoRCPir3/AKnczrZ2aPnz9/WtFDhCZ8SchqjczG2Zssmnv3dz9GDAc7/0jeYtuQUel+sxr2mEWOajfSb5mIpUa+BUBRlpIMmAinlJwGEEC8DM6WUTZHhbODvxyW6GFFf/w+83noWzH+QcANuI+10D/D7mlauyEyOyikhKSXvP/Nf3n7kQVJy87j4a9/BGUqg9fcfEGwbwHlSDgkri1TD7oqijGosVyvzB5NARAtw8BfnxBm/v5Pqmj+SmrqClJSTD5gekpJbd9XhMmrcVZp77Lc/4OHFP/2G3e+toezEUznnpi/iW9tO62ubMLjMpN04G2vpxG3qUlGU8TeWRPCaEOIl4JHI8MeAV8eyciHESuC3gAG4X0r5s1HmuZJwewcS2CylvGYs654o9lb/jlDIc9C3i95f38YHfR7+NLPwmN8l1NveypN3/4CO+lpO//inmH/6BXQ9XIGvshvbvHSSLy1Fs6nXQiiK8uHGctfQLZELx6dGRt0npXzyUMtFGrm/BzgbqAfeF0I8I6XcMWyeqcDtwMlSyi4hRMaRfInx0t+/h4aGh8nJuQqHo/SA6TUDPn5W1czZqQlcmpF0TLfdXFnBkz//AaFAgI/cfhdZ9mJaf/cBujdE8kemYl+Sqe4IUhRlTMb0czFyh9DhXhxeClRKKasAhBCPApcAO4bN8xngHillV2Q7rYe5jXFVueduNM1GSfGXDpgmpeQb5XUYBNxdlndMK+WK997hhT/8GkdSEh+948eYKwTtr24NXxC+cQ6mLMcx25aiKJPfWN419BHgbsJ3Dwn2NVV5qKueuUDdsOF69r3SelBZZBvvED599H0p5YtjC318dXatpb39NaZM+SZmc9oB0x9t7mR1l5ufleWRYz12T+y+/78nWP2vB8gum84lX7wd70st9G7vwDY/neTLpqrXQyiKctiElPLDZxCiErhISrnzsFYsxBXASinljZHh64BlUspbhs3zLBAArgTygNXAHCll937rugm4CSAzM3PRo48+ejihDHG73TidR/8UrZQ6uvwh0I8mfowQI+/J75aCr+EiH53v4eZIXyw6PF4pJQ1r36Rl83qSp0xj6gnnkbPFjNkN7dMlPYUynKLHybEq2+MlluKNpVghtuKNpVjh6OJdsWLFBinl6C9Ak1J+aAe8c6h5DrLcicBLw4ZvB27fb557gU8OG34NWPJh6120aJE8Um+88cYRLztcY+Pj8tXXSmRT09OjTv/01ipZsGqTrOwfOKrtDMYbDATk83/4lfzllRfIV//6J+nZ1SHr71wjG+5aIwcqOo9qG8fKsSrb4yWW4o2lWKWMrXhjKVYpjy5eYL08SL06lmsE64UQjwFPMbKpykNdM3gfmCqEKAYagKuA/e8Iegq4GvibECKN8KmiqjHENG5CIQ979vyKhIR5ZGZedMD059u6ebath2+XZDPFbh1lDYcn4PPy7G/upmrj+5x05bXMyV9Oxz+2Y0y3k3b9TIyptqPehjJ5SSnRdZ1QKDTUP9hnXdeHusHljqRfV1fHmjVrDjnfYHxy3w/BUT9Hc1pnZyc1NTUHzDc4vH9Zjve0xMToPBA6lkSQAHiAc4aNkxzi4rGUMiiEuAV4ifD5/weklNuFED8gnJmeiUw7RwixAwgB35BSdhzB9zhuamr/is/fwuzZvzvgAnBPIMjtFfXMdtr4fP7R3wAV8vv570/upKF8B2d9+gsUyZl0P1GJpSyZ1Gumo1nVraGxIBQKEQgE8Pv9+Hw+/H4/fr+fYDBIIBAgGAyO6PYfV1dXR3t7+wHzjFapj1bJj4c9e/aMOl4IgaZpCCGG/n4GP+8/fCymHWq+UChEMBgcMW3/eYaP3//z4Uw7FuvZP1kcK2O5ffSTR7pyKeXzwPP7jfvesM8SuDXSTXg+Xyu1tfeRkX7eqG0N/GBPI+2BIP+cW4LpKFsc83k87H72cTxtzVxw8zdIr8+gb3MdjiVZJF06BWFQjcZEm67r+P1+vF4vAwMDeL3eA7rB8YOV+2hdMBg87G0LITAajRiNxqE4BodNJhNWqxWDwTDUaZp2xJ8Hh4d3wyvsw+2vWbOGU089ddTpE82qVatYvnz5eIcxZqtWrYrKesdy19DfCB8BjCCl/FRUIprAqqr+D10PMGXKNw6Y9nZXHw81dXJzQQZzXQe+dO5wePvd/Pcn36O/rZmLbv4myTsSGKhqI2FlEa7Tj+2tqPFE13X6+/vp7++nq6uLrVu3Dg17PJ4DPvt8vkP+ArNarVgsFiwWC2azGbPZjNPpxGw2jxi3f2cymYa6wQp+sJI3Go1DlSfEXmU1mKiU2DGm9giGfbYClwGN0Qln4upz76Kx6T8U5H8Ku71wxDRPSOdru+ootpn5elHWUW1nwN3Hf398B2011Uw74zISNznwNfaQ8rFp2BfE1PN2x1UwGKS3t3dE19fXN+JzX1/fiIp98+bNQPjXt91ux26343A4yMrKwm63Y7PZsFqtWK3WEZ8HO4vFMiF/5SrK4RrLqaH/Dh8WQjwCvB21iCYgKSWVu3+C0ZhIUdHNB0z/xd4marx+/jt/CrajOGXj6e3h8R/fQWd9LZd84dvwhpeAz0PqdTOxzUg9mq8Q86SUQ7/kR+t6e3sPWMZsNpOQkEBCQgJpaWkkJCTgcrlwOBzs3r2bk046CYfDgc1mUxW6EteO5GrjVMIPl8WNjo5VdHa9w9Sp38VkGnnV/oNeD3+ua+O6nFROTnYd8Ta8bjeP/+i7dDU2cOnnvoPlHUnAC+mfnoWlJOkov0Hs0HWd7u5u2tvbaWtrG+ra29vx+Xwj5nW5XCQnJ1NcXExycjKJiYlDFb/L5frQ0xNtbW1kZMTVbqwoBzWWawR9hK8RiEi/GRi9Qd5JSNeD7K78GTZbEXm5146YFtAlX9tVS4b56Bqb8Q94eOJnd9LZUMelN30X85tBpC5pWKpTMImTgN/vp6WlhaamJpqammhubqatrW3ExVWHw0F6ejpz584lNTWV5ORkkpOTSUpKwmxWbSwryrEwllNDR/4zdxJobPo3Hk8lc+f8CU0bWfHcU9vCjn4v/5hTTMIRNjYT8Pt46uc/pHnPbi7+1LcwvxlAGDXSbppLxY51x+IrTAjBYJCmpibq6uqGKv6Ojo6hc/Y2m43s7GyWLFlCWloa6enppKenY7Op5yQUJdrGdGoo8r6hUwgfEbwlpXwqmkFNFMFgH1VV/0dS0lLS0s4eMa2i38uvq1u4OCOJc9OO7CGPUDDA/379U+p2buPC627FvlYgTBrpN80NPyi249DrmKjcbjd1dXVDXWNjI6FQCICEhASys7OZNWsW2dnZZGdnk5CQoO6GUpRxMpZTQ38EStnXHsHnhBBnSykPvGo6yVTX/JlAoJOppbePqKR0Kfl6eR12g8aPpx5ZYzN6KMTzv/8Vez9Yz8qrvohrg3VkEogxwWCQ8vJy9u7dS1VVFa2t4RfJGgwGsrOzWbp0Kfn5+eTn5+NyxfVBpqJMOGM5IjgDmBF5+AshxD+A7VGNagIYGGigru6vZGVeSkLC3BHT/t7Qzrqefn47vYB08+E3Ai+l5OX7fk/Fu29z5mWfIXlbQswlAV3XaWpqoqKigqqqKurr68PvLDEaKSgoYM6cORQWFpKdnY3JdPhlpCjK8TOWRFBJuGnKmshwfmTcpLan6peAYMqUr40YX+/18+OqJpYnu7gy68iagHznsX+yfdWrnHr+dWTszkAYNNI/M/GTQCAQYO/evZSXl1NRUUFfXx9CCHJycsjPz2f58uXk5+eril9RYsxBE4EQ4n+Erwm4gJ1CiMErl0sIv1Bu0urp3UxLyzMUFX4Bq3Xf3UBSSr5ZXocEfj7tyJ7w3fTSc7z35L9ZdOrF5NUVIgwifCSQNjGTgN/vp6Kigu3bt1NZWUkgEMBsNlNaWsq0adOYOnUqdrudVatWUVJSMt7hKopyBD7siOCXo4wThJusvCo64Yw/KSW7d/8EszmNwsLPjpj2REsXr3f28cPSXApslsNe9+51a3jtb/cyff6pTOubh5SStM/MmXBJIBgMsmfPHrZu3Up5eTmBQACn08m8efOYNm0axcXFGI3qhXeKMlkc9K9ZSvnm4GchxALCr5D+KLCXcDsCk1Jb28v09Kxn+rQfYTTuawCi3R/kjsoGFiXY+VTegS2SHUr9ru0897tfUDBlDguNK9D7g6TfNBdT+tG9l+hYam5uZuPGjWzZsgWv14vNZmPu3LnMnj2bwsJC9fStokxSH3ZqqIxwWwFXA+3AY4RbNFtxnGI77nTdT+Wen+FwTCU7+6Mjpt2xu56+oM6vpudjOMxTQh31tTz18x+QkpbHKamXoXf5SbtxDuac8W8Zyev1sm3bNjZu3EhjYyMGg4EZM2Ywd+5cpkyZgsGgmr5UlMnuw47vdwFvARdKKSsBhBBfPS5RjZP6+n8xMFDL/Hl/Q9P2Fc3L7T082drN14uymO44vNM4fR3t/Pcnd2Ix2Tmr+HpCrT7SbpiFpfBQTT5HV2trK++99x5btmwhEAiQkZHBeeedx5w5c7DbJ85RiqIo0fdhieAjhK8FvCGEeBF4lHFtFTe6AoFu9lb/gZSUU0lNPW1ofF8wxLcq6pnmsPKlwsN7N43P4+GJn30fv8fDZYu/ht7gJfXaGVinHtndRkdL13UqKyt59913qaqqwmAwMHfuXBYvXkxOTo56oEtR4tSHXSN4CnhKCOEALgG+AmQIIf4EPCmlfPm4RHic1NTeTzDYy9TS20eM/8GeRpp9Ae6fVYT5MM6R66EQz/32bjrqa7ni9NuR1T6Sr5iKbfbhX184WsFgkM2bN7NmzRo6OjpwuVycccYZLFq0CIfDcdzjURRlYhnLu4b6gYeBh4UQyYQvGN8GTJpE4Pe3U1f3dzIzL8TpnDY0/q3OPv7Z2MHn89NZmHh4FeaqB+9n76YNXLTiq2jVQRLOLsSx+OjaKjhcgUCAjRs38s4779Db20t2djYf+chHmDlzprrrR1GUIYdVG0gpu4D7Il1MCfX0YN6xA3nSSYj93lpZXfNndN1HSfGXh8b1B0PcWl5Hic3CN4uzD2tbH7z0LB+8+D9WnHQ99mozjiVZuM7IPybfYyz8fj/r169nzZo1uN1u8vPzufjii5kyZYo6/aMoygHi5meh+623Sf7d7/GdeSbWaft+9ft8LTQ0/Ivs7I9gtxcPjf9JVRP1Xj9PLSg9rMZm9m7awBt/u48Fc84jozkb67Rkki4tPS4VcCgU4oMPPmDVqlW43W6Ki4u5/PLLKSoqUglAUZSDiptEYCmbCoCvomJEIqir/ye6HqR4WMtj73a7+WtDO5/OTWNZ0thv8WyvrebZ3/yM4sIFlPnmY8qxk3LNDIQhupWwlJIdO3bw2muv0dnZSUFBAVdeeSUFBQVR3a6iKJND/CSC4mKkwYCvvBwuugiAUGiAhoZHSE8/G5stXGl6Qjpf3VVLgdXMt6eM/ZRQf3cXT/78ByTaM1jqOBfNYiTthllolujeh19TU8PLL79MQ0MD6enpXH311ZSVlakjAEVRxiyqiUAIsRL4LWAA7pdS/uwg810OPA4skVKuj0osJhPBrCy8FRVD41pbXyQY7CY/7xND4+6qbGDvgJ/H50/BMcaHqQJ+H0//8kcE+3ycN/1G8AvSPjkbgyt6LWj19vby8ssvs23bNhISErjkkkuYN2+eevpXUZTDFrVEIIQwAPcAZwP1wPtCiGeklDv2m88FfBl4L1qxDArm5eIr35cImluexmrNIylpKRB+cOwfjR18Lj+dU8bY/rCUkpf+9FtaK/dw2cKvQa9O2o1zMGVE56GsYDDI2rVrWb16Nbquc9ppp3HKKaeoZhsVRTli0TwiWApUSimrAIQQjxJ+HmH/drd+CNwNfCOKsQAQzMkl+N46Qt3dBO1BOjvfoajwcwghaPEF+MquWmY7bdxeMvZTQmv+8zAVa97iokVfwtAFKddMw1J0ZC2WHUplZSXPP/88nZ2dTJs2jXPPPZeUlJSobEtRlPghBtuMPeYrFuIKYKWU8sbI8HXAMinlLcPmWQh8R0p5uRBiFfD10U4NCSFuAm4CyMzMXPToo48eUUyh9evJuf+vtH/zBvxFLwPVwPkgLuFuEtmJkZ/RR67Qx7S+joodVL/2PCeWXklBqJi26To9RceuPN1uN06nk0AgQGVlJS0tLdhsNkpLS0lNTT1m2zkWBmONFbEUbyzFCrEVbyzFCkcX74oVKzZIKRePNm3cLhYLITTg18ANh5pXSjn07MLixYvl8uXLj2ibq7u6CbkkMu8xhOgHDEj5PC+YCtgSOJufleVxbe7Ynvyt37WdTX95maVlF1EQKMZ5cg55F005orgO5o033iA1NZUXXngBr9fLaaedxqmnnjohG35ZtWoVR/rvMh5iKd5YihViK95YihWiF280E0ED4dbMBuVFxg1yAbOBVZE7XLKAZ4QQF0fjgvHq+tV8r/83fGOlCauxB6MhgbS0M6iwXsRDNQmczNucI3KByw65ru6WZp755Y+ZmrmY4sBMbHPSSLzg2DbK0tPTw7Zt2+jo6CAnJ4eLL76YrKzj+2SyoijxIZqJ4H1gqhCimHACuIpwmwYASCl7gKGf3x92auiYcIco6smkZ3YPCa12/Bm9eBync3ttKtPtGrca17Fz1zt4BvYwpeRWwgcsB/L2u3ny7rtINmYxz7Ycc76LlCunIbRjc7umlJKtW7fy3HPPEQgEOPfcc1m2bJm6G0hRlKiJWu0ipQwCtwAvATuBf0sptwshfiCEuDha2z0Yx1adH9XdgtQSMdQG8GPmW83FSODvc0s5ccH95OR8jJqaP7F1282EQp4D1hEKBvnfr39KqN3HKRkfwZhiJe36mQjTsSlGj8fD448/zhNPPEFGRgZLlizhxBNPVElAUZSoiuo1Ainl88Dz+4373kHmXR7NWMwuK6Aj/C4INvAPcQvbPSH+OaeYwkizk9On/RiHYyq7d/+EDRuuYu68+7Basgbj47W//pHWnZVcNO1mDAYjaZ+cjWY/NufrKysreeqpp/B4PJx55pmcfPLJrF69+pisW1EU5cPEzZPFtlQX0IPmT2VN0WJWcTpfLczk7LR9t3oKISjI/yR2WxHbtn+Z9euvYP68B3A6y1j/7JPsfOMNLpxxM4aQgbRPz8KYYj3quAKBAK+88grr1q0jPT2da6+9luzsw3vJnaIoytGIm0Rgz0hCpwezP4tHCk4mP9DK14rmjTpvWtoKFi18lE2bP82GjVeSEPwybz30FOeWfRqLz0LK9TMw543tgbMP09HRwX/+8x+am5tZtmwZZ5111oS8I0hRlMktbk4+OzOS0NHZ5sinQcsnTdZz3dYqVnX2MtqzFC7XTJYs/i+B7ize+Ot/OLnwYhL9ySRfOhXb9KN/iGvr1q38+c9/pqenh6uvvprzzjtPJQFFUcZF3BwRGGQPfQYPG13Z2KWbwj31vO8c4KrNVZyc5OQnZXlMc4w81eN3m6h8Lo25abPJFWXoi9uwLznlqOLw+/28+OKLbNy4kfz8fK644goSE6PzJLKiKMpYxM0RwYBnN92GPrY7spnHB3zxLw/zcrCTH0/NZWf/AOesL+cPNS2EIkcH/gEPT919F7naVKbbTsY3ZS+7k79BecWdSBk6ohja2tq4//772bhxI6eccgo33HCDSgKKooy7uDkisFiy2JTcTp/BxiKxHYNuJfD223z69NO4OCOJb1XU86OqJp5v7+F3ZXls/eMvMHWYmJ95BtayZHKuOxFZ3URN7X34fC3MnvUbDAbbmLe/Y8cOnnzySUwmEx//+McpLS2N4rdVFEUZu7g5InC5ZrIrIXwOfoG1D/uypfRHbs9MN5u4f1YR984spKrfy09+/Qt6tjdwStZlmHOdpFw7A81opLT0NsrKvk97+2ts/OBa/P6OQ25X13VeffVV/v3vf5OZmcnnPvc5lQQURZlQ4iYRlL/7DrV2ExbpI9fqwnnqafhravDX1ADhW0cvzUzmtw2bWVZRyQm5V9HtMGH7+HQ08752CfLzrmPunD/idu9i/YaP4vHUHHSbAwMDPPzww7z99tssWrSIG264gYSEhKh/V0VRlMMRN4kglFtAs91KFk1YDdm4zlgBQO+LLw3Ns+XVF9n9zP84u/gTGC1WPj3Pwspd1WzuG/mUcXr6OSxc8C+CwV7Wb7iC3t4tB2yvpaWF++67j6qqKi688EIuuugijMa4OROnKEoMiZtEsMvooNNqJZt6zN40TLm52BYsoPe55wDY/f5aVv31L5xd8gkswkbBp+fw21PK8Oo6F27Yzf31bSNuM01MXMjiRf/BYHCwYeM1tLe/MTRt+/bt3H///QQCAT75yU+yePGob35VFEWZEOImEXQ09OE22cmmiUB1+Jd5wvnn46uooPLF53jht79ieeHVOPQEUq6dgTnfxYlJTl5dMo0VKS6+u7uBz+2owR3cd8eQ3V7M4kX/weEoYcvWz1Jf/xivvvoq//nPf8jMzOSmm24iPz//YCEpiqJMCHGTCNISrEihkU0DwT3hr51wwfm0J7l49u9/5qScS0khk+TLRj4wlmIy8vc5xXynJJv/tXZz3oYKdvUPDE23WNJZuOBhHI5T+O9/X1XXAxRFiTlxkwgGLOHXRGfRiNWdQWvNXnbv3Mr6wkyWJp5BllZIwsoiHEsPfOe/JgRfLMzkP/On0B0Mcd763TzR0jU0vaOjn3fXzqWnJ4fSqWspLV2DemGooiixIm6uXuZ2N3NG0//IyWpABKz893vfxePt4cSCiygwzMSY2k3C8g8/jXNysotXF0/js9ur+cKOGtb19HOVv5vnnn4ai8XCDTd8mkAwkerqe/D525gz+3cYDNFpxF5RFOVYiZtEUGISnNX7PJZkHbOwk2bIJTd/OQWG6eju7fStfZSUaxbir67GX1ODHPAiQyEIBdHsdjRXAprLidOVwMMZGfzKksD7b73JU3W7ycjN47qrPobL5QJuxWLJorz8TjZuvJZ58/6C2Ty25i8VRVHGQ9wkgvyZc3CX2+gQHQRtkhMzLgKgOriDnDI/wVdb2HPmWWNal89sJv3EE1iYnU1Kazv+hna2oTFv0QKs08rIy70Gizkj/CrrDR9l/ry/YbcXRfHbKYqiHLm4SQQACSY/vbrEfU0SSU94kV0+qNuK99nXCJmMGK028v/v11imlqE57AiTCaFp6B4PoT43el8vzfUNvPj+Ovr8fk4LBMnp7sZdXo7jrdeoARACy9Sp2JcsYcasW9htuZ/1Gz7KvHn3k5gw+muvFUVRxlPcJIK36t/CYuigz2/ltV83YpRmMixBivfsQqQWsakoh56BBpa+tYolJ52EGHa115CYiCExkW3bunh63XtYLBY+ed11Q7eGeoIhfrL2A7Zv3MxZ7Y2cW7+X7iefRD7kIQ0I5hmpmH01eZfcStZpn0QYDKMHqSiKMg7iJhHoIUmi0HD3FuATRmxOE+3dPtoW3UaxXXCCyUyvt5IN21dRcftXWPGpz5M7bQYAoVCI1157jTVr1pCfn8+VV14ZuR4QZjca+NGpi3m4tJhvV9TzF6OR+6blMqehhv516+hb/QaGlz+g+8Vf0ZPwRxLOXEnC+efjOGEZQrVBoCjKOIubRBDYacNs8uHoKaY1v5NPXzaPio/eyqa5X6BKFFI1ECTRWExeShHevmae+vG/yJiWx+yzV7Bmyxrq6mtZvHgxK1euPOirIq7JTmWey86N2/Zy2bZq7ijJ4bM33kjaZz6Dr7OZ8n/fRGDNLuRLz9Hz5JMYkpJwnXsuCeefj33xInWkoCjKuIibRDBzmYudH0ikN4mtCbtwv9KKOeSh4JTp7Ng8QMmCdHxuP7W7uwiRBfYsGpr62Pm/J9G1AK6eMuqed/Dgm2uxOkxY7CacKRYSUm0kpFlJSLORku1gZqKVlxZP46u7avn+nkbW9fTzmxkFJKRkMfum/7Lr9O/SWPs42U0nYd+SRs8zz9D92GMYs7JIvOxSkj7yEczqaWRFUY6juEkEDmP4PUEhr4tyy076Vm3BNn8+p31mGT2/30T1lnYu/eoC0lJ0dl/xaaoWnMlGZwdmaSajLQFvsI7cmem40rLwD4QYcAdoq+mjamMbur7vHURWp4nUXCefzHOwyOXgt3VdnNs/wP2zi5nltDFj+s+wWnPZa/otnjkLmf29p/Gv2Ur3U0/R8ef76PjTvdiXLiXpisvBaj3Y11EURTlmopoIhBArgd8CBuB+KeXP9pt+K3AjEATagE9JKQ/+XuejEGpYC8BA0ITN24xvew3pt96Kwaix8qY5PP6z9bz4l21c9s357LjsDLZ1tJNvdnB28nI0vxcPbjZsfo7utCCnX38jhXMWAKDrEneXl962ATqb+mmvd9NR72bb6kZkQOdLQJ9D455VnSyZncFlS/MoKf4SDkcpO3Z8gw3br2HuqfdRcP59BJqb6XnqabqfeILGb95GutVK09q1JF99NdaysmgUi6IoSvQSgRDCANwDnA3UA+8LIZ6RUu4YNtsHwGIppUcI8Xng58DHohFPwNcKgAeduS3dADhOOgkAq8PEeZ+bwyO/XMW99/wZr97H3L4+pr/xPzKfvhq9z4bxhb2cyuV0BJt47ee/J3lOIad//NOk5OSGTw+l2sgb9o4iPaTTXu+mqbKH6opOLBVdeKrreejZelxpVopml5JV8iBtA19j/YaPMmvWr8jIOpe0z32W1Js+g2f9esrv+SM9TzxJ9yOPYlu8iJRrrsF11lkIszkaRaQoSpyK5hHBUqBSSlkFIIR4FLgEGEoEUso3hs3/LvDxaAWTnH8pHU8G8AR1Cjs9BIXGy14HlxBuRayyfjudyR8gAxrzi5ZzwcqpVK16k6ZvfYvCB/+BddpC+tc1I14xco7xBmobdvHv225j2pmnccIVV2NzukZsTzNoZBQmkFGYwLwz8wnqOr/eWMvqjc3Mbw3R/04j+iodg+lOnFl76ax8iFkn7WHajM8hNA3H0qX0ejzMnTePnieepOvRR2m49WsY0tNI/uhHSbrySkxZB74XSVEU5XCJ4e/YP6YrFuIKYKWU8sbI8HXAMinlLQeZ/w9As5TyR6NMuwm4CSAzM3PRo48+etjxuPvc7HizCa+9hqzKd5he5eNzZ32Lb8wK0NGwF7fbTUpKCimBMnp2m8k9QZDVso7Ev/0d90UX0X/B+QBoAUjaK0iqFkhdp7z7fSoGPiBjyRLSZ84b8fzBaDZLI3/Eji8En2ztp6wpQF+DJOARIELYMxpIzMsiId+CN9iP0+kML6jrmHfswP7maszbtoEQ+ObOZWD56finTQMhDrtMjiW3270v1hgQS/HGUqwQW/HGUqxwdPGuWLFig5Ry1MZRJkQiEEJ8HLgFOF1K6fuw9S5evFiuX7/+sOPZvGErTz7zXwQhpDBg9XuQJg2fsJKQmMiZZ5zB3Llz0XXJM7/ZRGt1L5fftgjfb35A7wsvUPivf2JfsGBofcFuL70v1eD5oJUAfra0r6I7sZPTr/sURfMXfWgsbf4AX9pZyxudfZyblsCvyvIJNXnYvPpNqrd4CLgzQIA9HZadM50pC9Ox2Pc9b+Cvr6f7scfo/s/jhLq7MRcXk3z1VSReeimGcXr19apVq1i+fPm4bPtIxFK8sRQrxFa8sRQrHF28QoiDJoJoviy5ARh+H2ReZNwIQoizgO8AFx8qCRyNFEc2aa0nc1JgHUvffQ+r302h1sQJ+jr8BUuYN28eQggMBo1zPzMbi93IC/duJekb38aUlUXj179BqK9vaH3GJCspH5tGxi3zcZSksyjtHE4wnMfa3/6TJ376fToa6g4aS7rZxENzS/hBaQ5vdPRx1voKdidqnHP9xVz53WmUXXgPaTOfJ+jx8Ma/dvHAN9/mhT9vpeqDNkIBHXNeHhlf+xqlb64i+2c/RUtw0fKTn7L79OU0fe9OvOXl0SpGRVEmoWgmgveBqUKIYiGEGbgKeGb4DEKIBcCfCSeB1ijGQlNtB0IamNrTRHF1NZ2hvVTIYrq1VOZv/iH/fLucuk4PO5t62dnlJvmCPLb2evjTPyuo+/oP+cBvZdW3f0pNu5uufj/BkA6AOc9F+mfmkHr9TBIyMjk183JK2qbzv+/8iNf//mcG3H2jxqMJwU35GTy/aCouo8ZHN+3hB5WNWF3zOP2cv1B2ajvF532ZhR97h5mnZNBU2c0Lf97K3257mzf+tYuGii6EyUzSpZdS/NhjFD3+OAnnn0fP00+z95JLqb7mWnqefQ7p90ezWBVFmQSidrFYShkUQtwCvET49tEHpJTbhRA/ANZLKZ8BfgE4gf+I8DnuWinlxdGIJ32qmWDBamwV7QyQiH9BJjZxIveXl9AUdOJ9dhcFz71FKj0kCzeJoh/hkrR0wfNvg37qLDwY6fv1H3BLG33YCJhcGJ1ppCc6yU60kjXTSWq3hZRynVnma/BtKOfRd77O/MsvZO5Z52EY5Ynk2S47Ly4u467KRv5Y18orHT38bkYh8+f9lTfe+Doe+SC2gte57Jxf0dtUSMW6Zireb2HH2404ky2ULc2kbGkWqbNnYfvxj8n8xjfojlxcbvz612lJSyPpo1eQfOWVmLKzo1G0iqLEuKg+RyClfB54fr9x3xv2eWzvfT4GiouLmX3CqezY+SBZwLvubZzQ8DI3meuZbyqnUK/DIoKHvV7dI+j2JtPSkkJ9MIlGPZkdMoNntRzqndlojot44eU+El79A+ecuZSzT1tMsmPk7Z8Og4GfT8vnvLREvlZexwUbKrilIIMl4hKWLLiW7Tu+xgebr6Ko6BbO/MQXWH7NdPZubqP8vRY+eKWOjS/VkpbvpGxpFmVLMkn91CdJueET9L/zDl0PPUzHvX+m476/4DpjBcnXXIP9hBMQ43xxWVGUiSNuniwGsIQ8NPQnkGz20eHw8R3zX3GanIjsOfj7zPR21vM0ZxGyptHj1TFb7Ewrm0Z3uZmBZsFZH8sh9JcfEWquJfOrX8CcYkdzt5LS20BKbxPT+5qQPRVovt6hbQYwUm/OZLvMZ9vbO/niW6vpsJeGE1NuIgvyk1lQkITVZGBFagKrlk7nzsoGflfbSh4u/qLNYNnSZykv/z579/6Wzo7VzJr1a8qWFlC2NAtPr5/d61uoeK+ZNf+tZO0TleRNT6ZsWRYlS04k/9RT8dc30P3Yo3T/53H6Xnk1cnH5ahIvuxSDy/UhJaYoSjyIq0Tg6NuL3p1Nk70bcPPBx+7ntLKPgKZhDgXwP3AZVzU8y3f02/noZedTXbGDd3bsIhQKYUlJ4OGnQ1x28z0Evv4Zan70OIWPPIw5L29o/SLS0d8BHbuhoxJT+24KG3aSV7eFC0PvhmcMQlNFKlt2FfOWXsLvxXRkzkLmFWezrDiFu4qyuSA9iS9uqeSCDRV8MjeN26b9nNTU0ymv+B7vrbuAKSVfJy/vOuwJZuadkc+8M/Lpau6nYl0L5e8189rfd/KmuZzieelMW5ZF/le+Stott9D34ot0PvwwLT/5Ca3/938krFxJ0hWXY1u4UB0lKEqciqtEsLMzgMMTpMsFZs3Mc63vc9r0K8ITDSac1z2E5y/n88P2u7nrJQNfv/kWLrjgArZt28b69zfQ2rabvz9bydRrryP36SeQn/o0JQ8/hDFtv6YoHanhruAEIHxFXgO82/bS//wrGHp2YpXbOEGr4FxT+FbYQKuRbc3FvP9OGY/IMnrTFnOCxUxoWh4P1LTybFs3d5WeyrlLnqO8/LtU7P4BLa3PMmP6T3E4SgFIznKw7OISll5UTPOeHsrXtVC5voXd77dgc5mYuiSTactWUHTxxXi376D7scfofS78JlRzcTFJl3+ExEsuwZiefpz+RRRFmQjiKhFs7LZz2UAPhpmpBPQGXqt9DU/Ag90UaWDemoj90//D/ZcLubPzJ3z/HvjyZ7/A0qVLWbp0Kbu2VPHsI6+zp6GJ3UuWYPd4KP3uHZz23e+QNuzI4GCss4uxzPwMng9a6X5xD7IvRG3/LrrlK0ydYWKO1sq8ple4SX8OemCPzOGdxlk45Gw+SJrPzTs6mV+cwm+W3UNm1stUVPyQ99ZdRFHRFygsuAmDwQKAEILs0iSyS5M49aNTqdneQcV7zWxb3cCW1+tJzrKHryd8+XYyb/8WvS++RPd//0vrL39F6//9Bufpp5N0xeU4TzsNcZBXbiuKMnnE1V/5tl4nn/T1kTfzUlKtr9DubeeH7/6Qn576030z2VNwfuZZ+v9yAd/v/Anfv8fH5z/7ZQpS7UyfW0JqcjpP/2YDPksHlqwGttg62HL//RTl5bFgyRJmzJiB+UPeBSQ0gWNRJvZ56fRvaIGXNJI906nfUc1bchszL/kqM6YkotWvI2HDM1zb9w7XB19B7xfs6CvknZpZ/OjN2TRlLebiWX+n1PwwVVW/obnpCcrK7iAt7YwR2zOYNErmp1MyPx1vf4A9G1spf6+Z956p4r1nqsguTaRs6VKm/PkCtNYGep58gu6nnsL9+usY0tNIPP98Ei68COvsWerUkaJMUnGTCHoGAnR1h0/RJBSU8sC5V/OxZz/Gs1XPIqXkrpPvwhL5RY09BcdnnsPzwKX8qP1ufnVPOxd95k5mZCeQnu/i0q8u5unffIDWksnHFvay4/G/sdfn48n6ep5//nlmz57NggULyM3NPWjlKYwazmXZOBZm4n63EfmaIMtbRNPTVTxlXMu0S5bTOnsuK5afBg0b0PauZlrlKmY0vMxn9ecItBnY9PoU1uiz2GX5PP4UP1P2/omTS//DkjnfxG4vPmCbVoeJWafmMuvUXHrbB6h4P3yR+c2Hy1n9aAW5ZUlMWXIVxZ/8HPqm9+h+4gm6Hn6Ezn88iLmoiIQLLyTxwgswFxVF7d9JUZTjL24Swbq9naQMhO/mMWVlUZxYzBMXP8FVz17Fc3uf46XqlyhJKiHLkYXD6MBsMGNZsBx2g7PnP/zqsc2kzr6A2TlpWAwWLFdJtr7YSMVOI/OvOJOif/4JPasY74oL2bR5Exs2bCAtLY358+czd+5cEg7y6gdh0nCdmodjaTbuNQ3INyDbX0Lbk3W0+bawOeRj1oqzMRacgOn0b0JgAGrfpbvidRJ2vsYtvU9jCD3JQKuZ95unsWZHCg+9/FeSigs4b+GZnDi1CKvpwJbPEtJsLD6viEUrC2mvd7NnQyuVG1vDSUFAztQkplzxLQpv+z6hd9+g99nnaL/nHtr/8Aesc+aQeOEFuM49V734TlEmgai9ayhajvRdQ1VtblZ998uc+OYaip98AuuMcHvEvqCPb6z+Bm/UvYHL5CLdno4udfwhP76QD3/Ihy/gwY8+5m0Z0HAZEzAHzIgBgSPkID8hn3lF81hatpTCpEISLYmjHi3o/hD97zfT/WoVYgC6fC3sDW4j99wFzD3rXEzDGquRUvJmUwMvrn+Bkpa1nNm5kSkDewHokXbe1Weyjpn0ZZ9EycwlnFqWzoysBDRt9KMUKSUdDf3s2djKno2tdDV7QED2lESK56WTly3R1r1Oz7P/w7djJwDWeXNJOPtsdiQkcOqVV465jMZbLL1jJpZihdiKN5Zihei9ayhujghK0p2Y2UofYEhOHhpvMVr43Rm/47Wa17hjzR20elq588Q7WVm8ct/CUuJb+2eCr3yHapnCowXf4tqLz8JsCuHu9/Dm4ztpqu8kc1oQuf5fdAkP/uVl9GZZaeptotHdSJVexZt734RwPY3dYCc/IZ+ChAJKEksoSSxhStIUChMKcZ2ci3NZNhseXU1mbTYLezNxv9HFqy/+HwnLcpl7/vm4UtMQQrA8J49TLryRp1ov5/rqFnp6mrhyYCtXu9dzYt2bnOtbD60P0tqSxJrXZvK4cR568enMmjmbU6emkZ1oG/qaQgjS8pyk5TlZdnEJHY1u9mxso+qDVtb8txKAxIzpFF19CrlpARzlb9P/ysu0/vJXpAFVDz2M6+yzcZ19NpayqeqagqLEiLg5IgCov34efev8TNu4Ac1uP2B6o7uRb67+JpvbNnNByQV8a8m3SLImDU2Xdevo/9fHMXo7+YPpk5xx3e0sLEwhFNRZ/VgFO95qJLvIwcxN9yI3vUvazTeT9oXPIwwG+nx9rCtfx9qdaylvKqdP6yPgCOC1eWkPtKNHjjgEglxnLlOSpmDsNbJ89nKKu7Jwve3H3CoI6H5q+rcTKIFZF51DTtn0ofhCUvJ0aze/rm6m0uOjzG7lC9a9LN37S4y1O0npCeIMhd/rV61nskafyR7nIixTV7BwxlSWFKeQaDMxmt72AWq2dVC9tZ368i70oMRsNVAwK5X8PI2utx8nt6aCgQ0bQUpMOTk4Tj8N56mn4Thh2ajlPZ5i6ZdgLMUKsRVvLMUK0TsiiKtE0PDRGfTuEEzftv2gv1YDeoD7t9zPfVvuI8GSwHeWfYdzis7ZN0N/Oz0Pf4rEhjd5S5/DrqU/4RMrT8Fs1Nj1bhNvPlSOyWpgYWgtluf/hn3JEnJ+8fMR59IHBgbYvn07mzZtor6+npAI4Spw4ch3EEwMUuepo6qnir1dewkSfu2FURgptBdQ1J9FQWs6U7z5OLpDhMz9ZJ8xl+mnnIrZFq5sQ1LyTGs3f6htYbvbS4rJwFWpOqf6/oGt7r+k9kBCTzpJnQ1YQh4Adur5vC+n05gwH3PxycyYNp2lxSmkOi0HlJHfG6R+VxfVW9qp3tbBQG/4xXapeU5yC62k9O3Bvm0VvnffRno8CJMJ+5IlOE47Fedpp2MuLhr3o4VYqgBiKVaIrXhjKVZQiWDIUR0RXDQVT5OdsvWbDzlveWc531vzPXZ07OCU3FP45pJvUpwYuRNHSgbevR/t5e/i0zUesH+S5Vd/nfkFKXQ0uHnpL9voavZQmucn96kfYDboZP3gBySce84B22lvb2fr1q1s27aNjo4OhBCUlJQwe/ZsWtpamLZsGuVd5VR0VlDRVUF5VznN/c1DyycGnRR5s0lxW8lNyOSkE85h6ZwzsBgtSClZ0+3mvvo2Xm7vxSgE56doLNdfJKPrAYTuI0/MIcubi15bjb11I2Z9IFxWMo11+nRqHHMRhSdQPH0hy0rSyUq0johf6pLW2j5WP7cBkz+Rpj096EGJZhRkFSeQ5ewnoXkr5vdfIli5GwBjdjaOpUuxL1uGY9lSTLm5R/TveTRiqQKIpVghtuKNpVhBJYIhR5MI6s4uwe9PY8qb68Y0f1AP8tDOh7h38714g16unXEtN827iQRz5A6gziq6HvksyW3r2KRP4d0Zt/Oxiy/BZTKw7n972fRqLQ6Xgem1T5Ow6QWcZ51J1h13YMrMPGBbUkpaWlrYtm0b27Zto7u7GyEEU6dOZdq0aUybNm2oZaIeX084KXSUs7NmGztbtlNDA34tAIAmBdlaGjMyZzMrZy5lyWXYrEU81Sl4rLkLd0inyGrkXGslC9z34grsxmLJJjvzUnK0WZia99JX8RaWhvewBzoA6JJONuhT2WOegS9zPolTljFzSgFzchOxmgxDO2jAH6Jpdzd1u7qo39VJe50bCD/PkJFjIYV2Ehq3Yt34MqKjBQBTXh72ZUtxLFuGffFijNnZUT9iiKUKIJZihdiKN5ZiBZUIhhxxIggFqFkxDenKp+i5Nw9r0faBdn7/we95cveTOM1OPjHzE3x85sdxmBwgJZ6NjxB64bs4Ap08I07Hd/I3uWzFiXTW9vH6gzvpavaQneih8M3f4wx0kP6lL5J89dUHbYReSklDQwMvvvgifX199PT0AJCXl8f06dOZNm0a6fu9BsLv9rJzzfts2fY+NbKRaksje8w1tFv2vQDPZXZRmlSG0VpAdSiLimAGuimPE11wglzNtP5/4aCfpKRl5GRfTnra2Rh729Br1tJd/hZa3bskeaqH1lelZ7FVTqHZNYs2SwGzTjyfecWZFKU6hu5MGujz01jZTdPuHhoru2mv60NKEBqkpplI0bpwtO7Euu1trK17EEiMGRnY5s8f6qyzZqJZDjxFdTRiqQKIpVghtuKNpVhBJYIhR5wI+juoOmsppoJp5D/2/KHnH0V5Zzl/2PQHVtWtIsmSxLUzruXKaVeSYk0Bby+dL/wI5+YHkFLyjPFcTMu/zrlL5lL+diPvP7uXgD9EYbCC3Hf/jjMrkYxbv4rr3HMP+ut31apVnH766bS0tLBr1y7Ky8tpamoCICUlhSlTpjBlyhSKioqwDrutNNDmofvdGtwbm/H7Auw117PNuIM9zjqak/toNLTjCQ3s25ApE58xD2nOp8iRzELDHk6Qq0nUgqSmnEJ6xrmkp52FyZQEA93QtAl31To8e9/D1roZV6AtvF1pYJfMp0KU0JtQBllzSCqeT1lRPlMzXJiNGn5vkOaqHpoqe2jc3U1rbR9BXygchlmQYvOS4KnHUbMJe81GLL5uNJMJ68yZWOfODfdnzsBSUoIwjX5heyxiqQKIpVghtuKNpVhBJYIhR5wIOvZQefZK7IsWkvOXx44qhq1tW/nj5j/ydsPbWAwWLppyEdfNuI6SpBJkdx1N//sBmXseJygNvGBcgTzhZk5dsJRtL9Wy850mQJLfu4ncrf8haUoOqTd9BtfZZyMMIx/8Gu0fvaenh/Lycnbv3k11dTWBQABN08jLy2PKlCmUlJSQnZ2N0WhESkmwxUPHmio821oxe8JHIL2BTiqM5TQX9tNVGKLV6mZ7dyVN7jqI3L0kMWIyp5FqNFBk7GSqqY/SlNnMyT6bvMyzsdmGXfDtbWTLS/8g1+wmUPs+rp5yHMHuoZjrZRoVsoBW+1QCaTOw5s0jo2gmU7MSyXRZ6G720FLdS2tNH63VvXTUu9H18H5psUgStT4cPbVY67bh7KrG7mnCaNSwlJUNJQbrjBlYysrQbDbGIpYqgFiKFWIr3liKFVQiGHLEiaDxA8rPvYrEc08n65f3HZNY9nTv4Z87/sn/9vwPv+5nWdYyLim9hDMLzsTW20TjC78gvfJxjDLI28ynufRjzFryEVrf72TX2mak1MnsKye34lnSUiSpN9xAwoUXYnA6gEP/oweDQerq6tizZw979uwZOlowGo3k5eVRWFhIQUEBeXl5WCwWgt0+ejbW0bOhHmMHaBjQpU6nr4l+uxsxxc7ALDubrR7WduyhqnsvHm89hmALYtgDdQkGyDGbKErIoyx1AWUZJ9Na0cOlZ16KSTOBlNDXTKh5G917NzJQtxlLx05SBqoxRNbjk0aqZRY1IpduRxGh5KlYsqeTVjiLgowMLO4gbbVu2uv76Kh309nYTzAQucVWSJxGLy5vC9aW3di6anF4mrF527HmZGIpLcVSOgXzlClYppRiKSlGczhGlF0sVQCxFCvEVryxFCuoRDDkSBOB3P0auy66hbRrLyL9jp8f05g6Bjp4vOJxnqp8inp3PXajnXOLzuW84vNY5Cig/fV7sW97mKRQO20ykbcc5yCmXUPiQC573mvGPxDCFWgns3Y12b3bST/vdJIuv4L32ttYvmLFmOPo7++nurqa2tpaamtraW5uRkqJEIKsrCxyc3PJyckhJyeH9ORUvHu76Vy/F19VDxaPGYFGUA/QHWhlwOrBmG3DNyuf7bmZvNnXxIb23fh9DRgDjViD9WiBRnQZGNq+ADJtyeQnFFGYOIV8Vz55zjzyXfnku/JxCiO07cJds4ne+u2E2iqw9ewh2dcwlCAAmmQKe2UOHZY8vM58SCnGml6C3VGMPWSHngBdjf2017vp6/AOKwGJTQzgGGjD2lmD3d2EfaAVq7cDV7IF65RizPn5mAryKe/tZdH552PKy0OzjrwTaqKJp8rqeIulWEElgiFHmgiC7z7C7ht+QOYXbyDl5tuiEFn4Iu/G1o08Xfk0L1W/hCfowWV2cVreaZyRezozO3sJvPsQBe1vYUBnj8xhR8KZDDguItiRRFddPyBJ6a4go/l9kmgi/6Kzca08F+vMmYd9J43X66W+vp6amhrq6upoamrC5ws/UGY0GsnKyiInJ4eMjAzSklKxd0o825oI1vdj8pgxED5VNRDso0d24E+AuqJkdheks9OVwAaPn15fJ4ZgC9ZgAynBCqzBOgh24A168egjX8vhMrvIdmST5cgiy55FtjObTHsmWdZUsgMBnJ0tuOsr8LeUY+2uJHGgDofeN2Id3dJBncygzZRNvy0Pn70AaS5EM2QjZBJ4Tfi7/PS2DhCIXHuI/OtgDbmxetqw9oeTg22gA6u3A6dT4MpOwlKYH04U2dkYs7IxZWdhzMpC+5C3yR4P8VRZHW+xFCuoV0wctZApAwBDen7UtiGEYFHmIhZlLuL2ZbeztnEtr9e+zpv1b/Jc1XMYhIFZJbNYtvgblLb3Uly+jvP6HsLY90+aZTKbclbQaTiLfud0diVNA2D7xlpSX/kbGVorOQuLcZ1yIo4TTsCQlHTIeKxWK6WlpZSWhhuu0XWdrq4uGhsbh7pNmzbh9/uHlnG5XGQUZ5Celk6y0YWxdQBTkySpKwOLx0bBTsHJO/0E9Aa6g53sSRXszk5kqz2FpszlVKIRCLfThlXvIiuwleTgbqzBOsx6HwTbqe1s5YOWdfQFvCPiFQhSbalkZGeQXrKINNu5pBqdpEmJrX8Aa3cPzp420rqbmOWpI8W9DpN7ZDvTPmmiWSbTlpRKp6EYt6EEv5ZLUKSjh1wE/KV43KX4+0NE2pMLkxJzsxtLdRcW314svk2Y/d1YfN3YLTqORDOOdCf2rDTMOVkY09MxpqVhSE3DmJ6GISkJoWlHsNcoyviLm0SgG8OtiGkZOcdlezajjTMKzuCMgjMI6kE+aP2AtY1rea/5PR6o+DchGcKQbmBq6XKKgxaK2ttZ2vkGZ/qewOiAbeaF7AycTp99DtWu86hGoLUGSPh7NYm/fYO0xADpZemkzJ+Bff48zEWFh6yINE0jNTWV1NRU5syZA4STQ29vL62trbS2ttLW1kZrayvrN6wnGNxXyRptRpKTkkm0OrEGNGz9Aptbo8jtYF6FmY9jQuDBR4it9gG2Jwn2ppipTVzIXtsJdBv2/ao24yNdtpKr15EUqsYh27ELH0YCBPUg/SEfDb1VbGvfQpevB10OO7LQgBRwZDpJt51IislFgjBhC0rsviB2rxent5/kgT4yfJuYFXid9ICfBF0P7+wCcIHHYaMuVESrXkSPzMNDBn5SCIQS6POn0BWwENL3uyspBKIuiKnKjTnQjclfhzngxhTowxzsx2qWWO0GbAlm7El2bGkJWFMTMCUnYUhKxJCUtK9LSFCN/igTRtzsiaHe8P30hoTE475to2ZkSdYSlmQtAcDtd7OhZQOb2zazrX0bb3ds4wVrH3/KcWEUSWQJB7k+D6X9jzDb/2dyvVbwTKPJu4BWw0xqEs+mRmjQCsbn+nE+9hIufyuJzhBJ2S6SizNImlmEdcoUTFmZH1rhaJpGUlISSUlJlJWVDY3XdZ3u7m46Ozvp6uqis7Nz6HN1Z+e+JBG5vV8gcJht2DBiCxqZ0mZmXosdOxas0ozPaKHZaaPObqTWptPgctDkmMUu6xK82sgK1yb7SaaLZDoplR0kiR4c9GMVPkwE0KVOUAYZCPnpDfioD/TT4++j29dDQAbCMVkgvHtnDK3XabDhEmZsuobJF8AlNJzBKlzBXSQFB0gNDZAogzh1iUvXsYUMmAKJiEAKeigJdyATt56CO5SCN5iIX3fRL7MJShu62O85h95IVwWG4ADGYCemYD3GoAdjcABjwINJC2E2SSwWgdlqwGIzYrKbMdvNmJwWLC4b5gQ75oZq+s1mNJcLzeFEczowOJ0Im23cX9WhTA5RTQRCiJXAbwEDcL+U8mf7TbcADwKLgA7gY1LK6mjEEuoOP5RlSDr+iWB/TrOT0/NP5/T804HwtYW6vjq2tm9ld9fu8HuGevayzjSAZLA95AZMspGM4LMU+iT57mxS+nMw9xTgM+bTK5dQp1nClc9mEB8MYPOuxurtwCK8WK3gTLLiTHPgTHXgyEzAkZWCNSsVU1oaWkLCiEpF0zRSUlJISUk5IH4pJW63m97e3qFu27ZtJCcnDw039raPOKIAwBPuMhDkSyNWacYiTQRMdvqsDrptdvqsZvrsgh6bgW5zNuXmIrqMDkJi2K21ItwJoeMw9uO09eGilwzpxkU/VtmHSfowiQCaCCIIIPUAId1HIOQjIAN093URMBupCw3gCRjoD0p0ebBrAYO1ei0GKbBKwp0useshHDKIM6jhDDixB5zYAk4sQSfmoB2z347Bb0cLOCDkwB9yMCBT0SkghA1dG3ah2h/puvff/jw27+rDEGrHEPRiCPnCne7DIIMYCGHQdDRNYjSAwSgwmAyYzAZMFiMmswGjxYjRasRoNWOymTDYzJhsFoxWEyabGYPVhNFqwWAzo1nMaBYLwmxGmM1okT5Go0o8k1TUEoEQwgDcA5wN1APvCyGekVLuGDbbp4EuKWWpEOIq4G7gY9GIJ+GC89lk0JheUBCN1R8VIQQFCQUUJIyM7ZU3XqFofhF1fXU09TfR5G6isb+RZncTr/TV0enfCGwMzywFLl8yCb40MnpSKOhOJcWQjsWaiiZz0LUE8GnQQLgDwA2yF2Nwe/hXasiDQfdiwIcQQTQthGaQCBMYTKAZNQxmgcFswGgxYTCbMNlMGC0W0jv6KHFlYk7Nx5RrwWA2I42CgICAlPhlCF8wgDfgx+P14vF5GRgYwOv3ExrwYhloIrHHh70jwP4v4JCA12Sm32xjwGzGazIRtAoCFgMBiwmfycSA0Uq7IZEagw2PZsNn2O8X+v5nzSIvQxUyhBk/FunDJPsw6eHOqHsw6P1ougdNH0DTfWjSj5B+pAwi9QB+GcQnA3TIILoMIKWfkO4lKN0ECRIkxKFouoY5ZMMStGMKWcKdbsESsGD3W3D6LDj8ZuwBM9agdWgegx7uBIkIzIAJKcIdInKEFYx0/YeKYjADRWaUOpoMIfTgfv0QQoYQBMN9GUIQQqBHuvCw1EM8+pd1CAEIua8PCE2CEAhtcHy4+VahEekEaCL8ZLpBQ9P2TdcMGpomEJpAMxgifYEwGDBoBjAIhGZA0zSEwYDQNITRgMFgQDMYwGBAMxgxGDSIjO/fU0VlwBJen2FwvRqa0RAZp4XXJwSaUUMIDc0U2UYkjvD2tPB6RTh+IbTI9xTh5KmJ8KlbIcLjI/3h3Xgm2WgeESwFKqWUVQBCiEeBS4DhieAS4PuRz48DfxBCCBmFW5mEpiHt9pg6L2sSJqYmT2Vq8tRRpwf0AD2+Hjq9nXR6O+nydoU/97fS3tdMtXsnfd5u+gNuBgID6D4DmteGs9+Jy+MiwWvD4bNh0+xYNDtGacMg7QgSQVgAE1IzIzVz+DmzwfpiVDnU1Q9+/rAZNcK18LDXUksdm5TY0EHqSBFEasFwXwTRRRBECKnpQAhdhEDoCM2LEMHweCGRmo4U4S4kJH6jEZ/JgM+g4TdrBM3gNwsCJo2ASSNoNBAyaPs6TSOoGQhqLgJaMkGDCa8wEhRGAsJEQJjwCxN6+Dc4IXGIfUnqIP0IGUDIAMhApBINRMYH940f6kfGEUDoAYz6AEbdizHkxqD7Mep+NBnAoAfRZAAhQxCpmCEEkTLUQhqabsCgG9AwInQDpqABc8CEJWDEFDRiDpgx6UZMIRPGkAFjZH6jbsQQ+WzQjRikAU2GxwkMaDLcCYwIaYgkIgNgBBHpIyKf9/Xl8L7YLzPrcBhtPx0jOlBE9Y6BQ855zMnB187L8C8dBqu8cF/IkcOD4zT97ajc5RTNWjEXqBs2XA8sO9g8UsqgEKIHSAXah88khLgJuAkgMzOTVatWHVFAbrf7iJcdD4cbrw0buZH/Bi+K4hp93sHz7AEZGOoHpJ+Q7sMX7MIX8OILefEHvYRCAfRAEPw6Bl8I4ZcIv0QLSEQQCEpEUGLQBQRBkxpIgdAjfamFKwwJQg9fTxBSgAz/AhJycBjCiSLyUzE8dd84SfjXk9w3fcRvKMnQsgLBUO0ytNzwWeXg5pFCRFYZHrdvfOTHGkR+scmhX21SSEAQ0sKdboCQAaRRoBsluiYIGUE3hDtpgKBBQ9fCXydcRBJdAIbw5/A4gS4MSIxIzYYUiegiEiMi/FkTSMPgvOHYJZG+EPu+5rDvMvh9Br9vCBgA+qVERwISXYT74d9hEjl83ND4EDqByLokuq6joYMMIXQ9kvwiRxRSho8SpA56CI1Q+LkWdNAl4VULNF2iSdAkiBAYQmAMSQwhwtN0Ee4i+49BgqaL8P4kQYRkeA+R4X1Ji6x7qCM8rzZ8XGReIv3BecL76769a99+KvaNG1awYnBfHVzXiL1y5Lr27aP7lty3HkYsM3L/JzIO3Ak9UanDYuLnsZTyPuA+CD9HcKQZMZ7uGT7eYilWiK14YylWiK14YylWiF680bzxuQEYftN+HsPOTu8/jxDCCCQSvmisKIqiHCfRTATvA1OFEMVCCDNwFfDMfvM8A3wi8vkK4PVoXB9QFEVRDi5qp4Yi5/xvAV4ifCXpASnldiHED4D1UspngL8C/xRCVAKdhJOFoiiKchxF9RqBlPJ54Pn9xn1v2Gcv8NFoxqAoiqJ8OPVyFEVRlDinEoGiKEqcU4lAURQlzqlEoCiKEudirmEaIUQbUHOEi6ex31PLE1wsxRtLsUJsxRtLsUJsxRtLscLRxVsopUwfbULMJYKjIYRYf7AWeiaiWIo3lmKF2Io3lmKF2Io3lmKF6MWrTg0piqLEOZUIFEVR4ly8JYL7xjuAwxRL8cZSrBBb8cZSrBBb8cZSrBCleOPqGoGiKIpyoHg7IlAURVH2oxKBoihKnIubRCCEWCmEKBdCVAohvjXe8exPCPGAEKJVCLFt2LgUIcQrQojdkX7yeMY4SAiRL4R4QwixQwixXQjx5cj4CRevEMIqhFgnhNgcifWuyPhiIcR7kf3hscir0icMIYRBCPGBEOLZyPCEjFcIUS2E2CqE2CSEWB8ZN+H2g0FCiCQhxONCiF1CiJ1CiBMnYrxCiGmRMh3seoUQX4lWrHGRCIQQBuAe4DxgJnC1EGLm+EZ1gL8DK/cb9y3gNSnlVOC1yPBEEAS+JqWcCZwA3Bwpz4kYrw84Q0o5D5gPrBRCnADcDfyflLIU6AI+PX4hjurLwM5hwxM53hVSyvnD7m+fiPvBoN8CL0oppwPzCJfxhItXSlkeKdP5wCLAAzxJtGKVUk76DjgReGnY8O3A7eMd1yhxFgHbhg2XA9mRz9lA+XjHeJC4nwbOnujxAnZgI+G2s9sB42j7x3h3hFvzew04A3iWcMO1EzJeoBpI22/chNwPCLeAuJfITTITPd5h8Z0DvBPNWOPiiADIBeqGDddHxk10mVLKpsjnZiBzPIMZjRCiCFgAvMcEjTdymmUT0Aq8AuwBuqWUwcgsE21/+A3w/+3dX4imYxzG8e8ldtPYLNoDGoVtQ01rrewBm5SzSUuiaGPVxoFNOUCknClFygkl5Uirdklrj4Q9kvxZf4eVJcpu7GiFiZLM5eC+h7dpp3famXeeW8/1qbf3+TfTNXW//d7n9zxzPw8Cs3X9HNrNa+B1SQcl3V23NTkOgAuBn4AXatvteUljtJt3zq3A7ro8kqx9KQT/ey5fAZq611fSGcDLwH22fxvc11Je23+7nGKPA1uAS7pNtDBJ1wPTtg92nWWRttreTGm77pJ0zeDOlsYB5UFcm4FnbV8O/M681kpjeanXgrYBe+bvW86sfSkER4HzB9bH67bWHZN0LkB9n+44z78knUYpAi/afqVubjYvgO1fgAOU1spaSXNP6GtpPFwNbJP0HfASpT30NI3mtX20vk9TethbaHccHAGO2H63ru+lFIZW80IpsB/aPlbXR5K1L4XgfWBDvfNiFeVUa1/HmRZjH7CjLu+g9OI7J0mU500fsv3UwK7m8kpaJ2ltXT6dci3jEKUg3FwPayIrgO2HbY/bvoAyTt+yvZ0G80oak7RmbpnSy56iwXEAYPtH4HtJF9dN1wFf0Gje6jb+awvBqLJ2fSFkBS+4TAJfUfrDj3Sd5wT5dgM/AH9RvrnspPSG3wQOA28AZ3eds2bdSjkl/RT4uL4mW8wLbAQ+qlmngEfr9ouA94CvKafdq7vOeoLs1wL7W81bM31SX5/Pfa5aHAcDmTcBH9Tx8CpwVqt5gTHgOHDmwLaRZM0UExERPdeX1lBERCwghSAioudSCCIiei6FICKi51IIIiJ6LoUgYog6Y+U9dfk8SXu7zhSxnHL7aMQQdT6l/bYnus4SMQqnDj8kovceB9bXiesOA5fanpB0J3Aj5R9/NgBPAquA2ynTX0/a/lnSeso06Oso0wnfZfvLlf4jIhaS1lDEcA8B37hMXPfAvH0TwE3AlcBjwB8uE5q9A9xRj3kOuNf2FcD9wDMrETpisXJGELE0B2zPADOSfgVeq9s/AzbWGVqvAvaUKZoAWL3yMSMWlkIQsTR/DizPDqzPUj5fp1CeJbBphXNFLFpaQxHDzQBrTuYHXZ7T8K2kW6DM3CrpsuUMF7FUKQQRQ9g+DrwtaQp44iR+xXZgp6S5WTpvWM58EUuV20cjInouZwQRET2XQhAR0XMpBBERPZdCEBHRcykEERE9l0IQEdFzKQQRET33D/Sy/qX02qnxAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "execution_count": 48, "metadata": { "needs_background": "light" }, "output_type": "execute_result" } ], "source": [ "def plot_matrix_pop(sol):\n", " series = []\n", " ts = []\n", " for t, s in sol:\n", " ts.append(t)\n", " series.append(s)\n", " plt.plot(ts, np.array(series))\n", " plt.grid()\n", " plt.xlabel('time')\n", " plt.ylabel('Abundance')\n", " plt.legend(['$x_{1A}$','$x_{1B}$','$x_{1C}$','$x_{1D}$','$x_{1E}$',\n", " '$x_{2A}$','$x_{2B}$','$x_{2C}$','$x_{2D}$','$x_{2E}$',\n", " '$x_{3A}$','$x_{3B}$','$x_{3C}$','$x_{3D}$','$x_{3E}$',\n", " '$x_{4A}$','$x_{4B}$','$x_{4C}$','$x_{4D}$','$x_{4E}$'\n", " ], fancybox=True, loc='lower center', bbox_to_anchor=(0.5, 1.05),\n", " ncol=4,)\n", "\n", "plot_matrix_pop(T2.solution)" ] }, { "cell_type": "code", "execution_count": 33, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "inits = sel_model(0,inits,[f_i,m])\n", "# inits" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Epigrass: Modelos metapopulacionais geo-referenciados\n", "Para criar e simular modelos metapopulacionais mais complexos a biblioteca [Epigrass](https://github.com/fccoelho/epigrass) pode ajudar.\n", "\n", "\n" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ ] } ], "metadata": { "kernelspec": { "display_name": "SageMath 9.7", "language": "sagemath", "metadata": { "cocalc": { "description": "Open-source mathematical software system", "priority": 10, "url": "https://www.sagemath.org/" } }, "name": "sage-9.7", "resource_dir": "/ext/jupyter/kernels/sage-9.7" }, "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.5" }, "latex_envs": { "LaTeX_envs_menu_present": true, "autoclose": false, "autocomplete": true, "bibliofile": "biblio.bib", "cite_by": "apalike", "current_citInitial": 1, "eqLabelWithNumbers": true, "eqNumInitial": 1, "hotkeys": { "equation": "Ctrl-E", "itemize": "Ctrl-I" }, "labels_anchors": false, "latex_user_defs": false, "report_style_numbering": false, "user_envs_cfg": false }, "toc": { "base_numbering": 1, "nav_menu": { }, "number_sections": true, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": { }, "toc_section_display": true, "toc_window_display": false } }, "nbformat": 4, "nbformat_minor": 4 }