{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div style=\"width:100%\"><a href=\"https://www.dunod.com/sciences-techniques/programmation-python-par-pratique-problemes-et-exercices-corriges\"><img src=\"https://python.guillod.org/livre/banner.png\"/></a></div>\n",
    "\n",
    "Ce fichier reflète les énoncés des exercices d'un chapitre du livre *[Programmation Python par la pratique](https://www.dunod.com/sciences-techniques/programmation-python-par-pratique-problemes-et-exercices-corriges)*. Tous les énoncés peuvent être téléchargés au format [Jupyter Notebook](https://python.guillod.org/livre/) ou exécutés directement en ligne sur [GESIS](https://notebooks.gesis.org/binder/v2/gh/guillod/livre-python/v1errata). Les corrigés sont disponibles dans le livre en version papier (ISBN [9782100815142](https://www.dunod.com/sciences-techniques/programmation-python-par-pratique-problemes-et-exercices-corriges)) ou numérique (ISBN [9782100819089](https://www.dunod.com/sciences-techniques/programmation-python-par-pratique-problemes-et-exercices-corriges)) publiés chez Dunod. Ce fichier reflète les exercices tels que publiés dans cet ouvrage et diffère quelque peu des exercices présentés à la page [python.guillod.org](https://python.guillod.org/)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 10 Probabilités et statistiques"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div id=\"ch:proba-stats\"></div>\n",
    "\n",
    "Dans un premier temps, la statistique de la proportion de nombres commençant par un certain chiffre sera étudiée. Puis dans un second, des modèles probabilistes importants seront introduits et simulés, comme les marches aléatoires ou la percolation.\n",
    "\n",
    "**Concepts abordés:**\n",
    "\n",
    "* statistiques et probabilités\n",
    "\n",
    "* série harmonique aléatoire\n",
    "\n",
    "* marche aléatoire\n",
    "\n",
    "* percolation\n",
    "\n",
    "* transition de phase\n",
    "\n",
    "* histogrammes\n",
    "\n",
    "* importation de données\n",
    "\n",
    "* optimisation par compilation\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.1: Loi de Benford"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La loi de Benford prédit que statistiquement dans une liste de nombres donnés, la probabilité qu'un de ces nombres commence par le chiffre 1 est plus importante que celle qu'il commence par le chiffre 9. Plus précisément la loi de Benford prédit que la probabilité qu'un nombre commence par le chiffre $d$ est:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "p_d = \\log_{10}\\bigg(1+\\frac{1}{d}\\bigg) \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "où $\\log_{10}$ désigne le logarithme en base 10.\n",
    "Il est possible de vérifier que la loi de Benford est la seule qui reste invariante par changement d'unités, *i.e.* en multipliant les nombres de la liste par une constante les probabilités précédentes restent inchangées.\n",
    "\n",
    "**a)**\n",
    "Écrire une fonction `firstdigit(n)` qui pour un nombre `n` donné retourne son premier chiffre et une fonction `occurrences(liste)` qui retourne le nombre d'occurrences des premiers chiffres de `liste`.\n",
    "\n",
    "**Indication.**\n",
    "Faire en sorte que la fonction `occurrences` fonctionne même si la liste contient des zéros en les ignorant.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "**b)**\n",
    "Vérifier si la loi de Benford semble satisfaite pour la suite des nombres $(2^n)_ {n\\in\\mathbb{N}}$ en comparant l'histogramme empirique avec la loi de Benford.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**c)**\n",
    "Vérifier si la loi de Benford semble satisfaite pour la suite des nombres $(3n+1)_ {n\\in\\mathbb{N}}$.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**d)**\n",
    "En allant sur le site de l'INSEE à l'adresse: <https://insee.fr/fr/statistiques/4171341?sommaire=4171351>, télécharger le fichier au format TXT contenant les données de la population par sexe et âge regroupé (POP1A). Importer ces données pour avoir la population par code postal, sexe et tranche d'âge.\n",
    "\n",
    "**Indication.**\n",
    "La documentation sur comment lire un fichier est disponible [ici](https://docs.python.org/fr/3/tutorial/inputoutput.html#reading-and-writing-files). À noter que les fichiers TXT de l'INSEE sont encodés au format ISO-8859-1.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**e)**\n",
    "Déterminer si la liste de toutes les populations par commune, sexe et âge suit la loi de Benford.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**f)**\n",
    "Sommer les données précédentes pour obtenir la liste des populations par commune et déterminer si elle suit la loi de Benford.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**g)**\n",
    "<font color=\"red\">!!</font> Lire la documentation du module Pandas disponible [ici](https://pandas.pydata.org/docs/) et l'utiliser pour refaire les deux questions précédentes mais avec le fichier de population POP1B avec les âges non regroupés.\n",
    "\n",
    "**Indication.**\n",
    "Utiliser la fonction `read_csv` de Pandas.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "**h)**\n",
    "<font color=\"red\">!!</font> En allant sur le site de l'INSEE ou autre télécharger son jeu de données préféré et tester s'il suit la loi de Benford.\n",
    "\n",
    "**Indication.**\n",
    "Utiliser par exemple les comptes détaillés de l'État disponibles [ici](https://www.data.gouv.fr/fr/datasets/donnees-de-comptabilite-generale-de-letat/).\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.2: Série harmonique de signe aléatoire"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "Le but de cet exercice est de simuler la convergence d'une série harmonique donc le signe est tiré aléatoirement. Plus précisément si $(X_i)_{i\\in\\mathbb{N}}$ est une suite de variables aléatoires indépendantes valant $-1$ ou $1$ avec probabilité $\\frac{1}{2}$, alors on définit la somme partielle:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "W_0 = 0 \\,, \\qquad\\qquad W_n = \\sum_{i=1}^n \\frac{X_i}{i} \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "et la question est de déterminer si la suite $(W_n)_{n\\in\\mathbb{N}}$ converge et si oui vers quoi.\n",
    "\n",
    "**a)**\n",
    "Écrire une fonction `sign()` qui simule la variable aléatoire $X_i$.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**b)**\n",
    "Écrire une fonction `simulate(n)` qui retourne une liste avec une réalisation de $(W_0,W_1,W_2,\\dots,W_n)$.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**c)**\n",
    "Faire la représentation graphique de la fonction $n \\mapsto W_n$ pour différentes réalisations par exemple pour $0\\leq n\\leq 1 000$ et émettre une conjecture quant à la convergence de la suite $(W_n)_{n\\in\\mathbb{N}}$.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**d)**\n",
    "<font color=\"red\">!</font>\n",
    "Déterminer l'histogramme de $W_{1 000}$ pour $10^4$ ou $10^5$ réalisations pour avoir une idée de la loi de la variable aléatoire limite.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.3: Ruine du joueur"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le but est de simuler l'évolution de la somme d'argent d'un joueur jouant à pile ou face. À chaque lancer le joueur gagne un euro si c'est pile et en perd un si c'est face. La probabilité d'obtenir pile est notée $p$, celle d'obtenir face $q$. En particulier $p=q=\\frac{1}{2}$ si la pièce est équilibrée.\n",
    "\n",
    "Mathématiquement, la somme $S_i$ possédée par le joueur au temps $i$ est donnée par une marche aléatoire:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "S_{i}=\\begin{cases}\n",
    "0\\,, & \\text{si}\\:S_{i-1}=0\\,,\\\\ \n",
    "S_{i-1}+X_{i}\\,, & \\text{si}\\:S_{i-1}\\geq1\\,,\n",
    "\\end{cases}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "où les $(X_i)_ {i\\geq1}$ sont des variables aléatoires indépendantes de loi $\\mathbb{P}(X_i=1) = p$ et $\\mathbb{P}(X_i=-1) = q$.\n",
    "\n",
    "**a)**\n",
    "Écrire une fonction `simulate(p,k,N)` qui génère une réalisation de longueur $N$ du processus à partir de $S_0=k$, c'est-à-dire qui retourne $(S_0,S_1,S_2,\\dots,S_N)$. Représenter graphiquement plusieurs réalisations.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**b)**\n",
    "Simuler un joueur qui, commençant avec une somme $k$, joue jusqu'à tout perdre ou avoir la somme $n \\geq k$.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**c)**\n",
    "Si $T$ désigne le temps auquel le jeu s'arrête, *i.e.* lorsque $S_T = 0$ ou $S_T = n$, retrouver par simulation les résultats théoriques sur le temps moyen:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\mathbb{E}(T)=\\begin{cases}\n",
    "k(n-k)\\,, & \\text{si}\\:p=q\\,,\\\\ \n",
    "\\dfrac{n}{p-q}\\dfrac{1-\\rho^{k}}{1-\\rho^{n}}-\\dfrac{k}{p-q}\\,, & \\text{si}\\:p\\neq q\\,,\n",
    "\\end{cases}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "et le lieu de sortie:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\mathbb{P}(S_{T}=0)=\\begin{cases}\n",
    "\\dfrac{n-k}{n}\\,, & \\text{si}\\:p=q\\,,\\4pt]\n",
    "\\dfrac{\\rho^{k}-\\rho^{n}}{1-\\rho^{n}}\\,, & \\text{si}\\:p\\neq q\\,,\n",
    "\\end{cases}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "où $\\rho = q/p$. Pour cela on pourra faire un graphique de ces quantités en fonction de $p$ ou se contenter de considérer le cas $p=q=\\frac{1}{2}$.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.4: Urnes de Polya"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Une urne contient initialement (à $t=0$) $r_0$ boules rouges et $b_0$ boules blanches. À chaque instant, on tire une boule uniformément au hasard dans l'urne. On remet ensuite cette boule dans l'urne et on y ajoute une boule de la même couleur. Un tel système s'appelle une *urne de Polya*.\n",
    "Le but de cet exercice est d'étudier le comportement de la fraction de boules rouges dans l'urne, c'est-à-dire le nombre de boules rouges sur le nombre total.\n",
    "On appellera respectivement $r_n$ et $b_n$ le nombre de boules rouges et blanches présentes dans l'urne à l'instant $n$.\n",
    "\n",
    "**a)**\n",
    "Écrire une fonction `densite` prenant en argument un tuple représentant le nombre de boules rouges et blanches dans une urne, et qui renvoie la densité de boules rouges.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "On veut construire de manière récursive la distribution du nombre de boules rouges au temps $n$, c'est-à-dire la liste des probabilités que le nombre de boules rouges soit égal à un entier donné $k$ (qui sera l'indice de la liste). Cela se fait en écrivant deux fonctions : `next_dist_rouge`, prenant en argument la distribution à un instant $n$ et qui renvoie celle à l'instant $n+1$, qui est donc la fonction qui fait tout le travail, et `dist_rouge` qui est la fonction d'enrobage, prenant en argument $r_0$, $b_0$ et le temps $n$ et qui nous renvoie la distribution au temps $n$ par un appel récursif.\n",
    "On utilisera les faits utiles suivants (faire un petit dessin):\n",
    "* La distribution passée en argument à `next_dist_rouge` est une liste `r` et `r[k]` représente la probabilité d'avoir `k` boules rouges dans l'urne à l'instant $n$. Les indices pour `r` varient de 0 au nombre total `s` de boules à l'instant $n$.\n",
    "\n",
    "* Au temps $n+1$, pour avoir $k$ boules rouges, il faut:\n",
    "\n",
    "  * soit avoir eu $k$ boules rouges à l'instant précédent et ne pas avoir tiré une boule rouge;\n",
    "\n",
    "  * soit avoir eu $k-1$ boules rouges à l'instant précédent et avoir tiré une boule rouge.\n",
    "\n",
    "* Si $n=0$, le résultat de `dist_rouge` est complètement déterministe et les coefficients de la liste ne sont que des 0 et 1, dépendant de $r_0$ et $b_0$.\n",
    "\n",
    "**b)**\n",
    "Écrire les fonctions `next_dist_rouge` et `dist_rouge` en utilisant les indications fournies. Regarder le résultat de `dist_rouge(0,1,n)` et `dist_rouge(1,1,n)` pour différentes valeurs de `n` $(1,2,5,10,20,\\dots)$ et commenter.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "Plutôt que de calculer théoriquement pour chaque $n$ la suite des probabilités théoriques, nous allons faire des statistiques sur un grand nombre d'évolutions d'urnes de Polya, au bout d'un grand nombre d'étapes. Pour cela, il nous faut une fonction pour faire évoluer une urne de Polya.\n",
    "\n",
    "**c)**\n",
    "Définir une fonction `polya_step(r,b)` qui, à partir de la composition d'une urne passée sous forme de deux paramètres `r` et `b`, renvoie l'évolution (aléatoire) après une étape de la composition de l'urne sous forme d'un tuple. Définir également une fonction `polya(r0,b0,N)` prenant en argument $r_0$, $b_0$ et $N$ en paramètres et renvoyant la composition (aléatoire) d'une urne de Polya au bout de $N$ étapes, aussi sous forme de tuple.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**d)**\n",
    "Écrire une fonction `data_rdens_polya(r0,b0,N,nbexp)` qui renvoie une liste de longueur `nbexp` contenant les densités de `nbexp` réalisations d'urnes de Polya au temps `N` initialisées avec `r0` boules rouges et `b0` boules blanches.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**e)**\n",
    "Stocker dans une variable le résultat de `data_rdens_polya(2,3,1000,10000)` et dessiner un histogramme pour voir la répartition des densités. Attention, on veut que les hauteurs des barres soient normalisées pour que leur surface représente la proportion de points, et non pour qu'elles donnent le nombre de points par *bin*.\n",
    "\n",
    "**Indication.**\n",
    "Une bonne règle de départ est de choisir le nombre de boîtes (*bins*) pour un histogramme de l'ordre de la racine carrée du nombre de points. Consulter la documentation de la fonction `hist` de Matplotlib.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.5: <font color=\"red\">!!</font> Percolation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le but est d'étudier un modèle de percolation dans un milieu poreux. Le milieu est modélisé par une matrice aléatoire de booléens qui détermine les sites qui peuvent être envahis par l'eau et ceux qui sont imperméables. Une matrice percole s'il existe un chemin d'eau allant de la ligne supérieure vers la ligne inférieure.\n",
    "Dans les exemples suivants, les entrées d'une matrice pouvant être envahies par l'eau sont colorées et les entrées effectivement remplies d'eau sont en bleu. La première matrice ne percole pas alors que la seconde oui:\n",
    "\n",
    "<center><img src=\"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAA+gAAAH0CAIAAAAhSpB6AABPxklEQVR4XuzW0YmEMBDH4bBV+GIxVmIr+mgFVuC7BWkRiqJM4swFT9jjDvbhRHbZ/X0PEQaSwAj/jLOXBwAAAMAdXwAAAAAM7oCq2tuhpQBAUunOLgIwuMN7LyIhhF/1bdtkp6oXhDJ4dQAQ/njtxoLBHRjHcV1XOwv3V3YYBhGxpwIAwn+apnme7XrAzX0SmFlcVbWu67Is27b9WYxr13WxXhRF3/dH8SHvfQxlVX1847IsSZJUVfW9xf0f7n8qTdOmaWIlhODOmXfuSgAI/7fsapRlWZ7nR9/OEZHY2OPkP4Cb+yT4Yu+uY249jjOAx8xMsd2amRlriJmZaopJsmWZ3ZohMdfMNXPNUDOzZGaGxMycOnXrpPTTfaTR0TnfPffaSfUl/nb+eLVn393Z2X2vnpn7zO5+JHHekUceCaB32mmn7777bpRRRgkB4BVMVw++X3nllU5Y1wXmehYqpXzppZcKH1977TVvv//++7RXxgcreyqnSxUiw1NbaNjzamAxhDYppH1s6NU24Fs1OvbvXlKziMKaZm9gPeBwyl3dFfovRb3qHctP4G7NlUdGSYarNp0N1Ky77rrrrbee+nTv8yH+iqRJkyYN/L0tbX3A37PXAJUj6Uqg8bfffts7bi3ICP1ImXTmmWda2I8++sir5FS7dDY0bhTskJO2cw4czDnnnLPMMssYY4xx3XXXVVQHx2eaaSb144wzzh133BFY6aPE8/LLLwcAn3322Y/Y6THyr9SMeF69P/90bf3l/384r3provDVV1/ldM8555z6Fv379m/wd8NkpD/E4EuTJk0a+PfXoNxr/5/XtlI499xzr7nmmqn80W4lCytwTxCvPIK5DBFp0va4N+wGAdNOO+2GG264zDLLrLHGGmpyIOn6668fb7zxDj30UBFhYTcoeeKJJ0488cR999336KOPfvbZZ1NPyZ133rnaaquNPvro22+//T8Mk/fee89bnM0999yjcOGFF+qlu/Lhhx9+//33V1+Ft95667TTTtPgmGOOefTRR1OZJ2cgJN1vv/2OP/74N998sw9sXXbZZYgihVtuuWX//fc/6aST3n///bytXh9++OEZZ5xBm6dyp7Zjjz32pptuUsAeseT888/HmlT3XndlcR555BFTOP300w137bXXUtU13LvvvnvqqafSpg3WpFzg559/rruA2970k08+mT1S0tH/3HPPsUQX9eiuKIzOJ598EkN24IEH8rIMKG30APcK3FX6NKeccgol2mfNK8p/6qmnDG2bKXrssMMOU7bg6qPtuOOO+5th8o/DxEfRsZaOF6HTwprXX1rs3qRJkwb+tv1oT+HHH3+ct/W0UefKK6884IAD2ABLq55of/XVVytcddVVLDGjvPrkk0/OPfdcCnW577772DwyrmSuueaynlVjNS6++GK4fdBBBzGv0yQuxjo/88wzVp6PoI0fKbi+8cYb7boZbbTRdt55Z6vKMPZUWM8L0GmRuaGhFbs3aYF7w26bDjfffPOLLrpo1FFH/fTTT1MPd1ZYYQXAISK866670kUMOvbYY0833XSzzTbbxBNP7NURRxyRV0BzsskmAzEzzzwz2JpjjjleeOGFpPBWWWWVZZddVq9pppkm2Kpyzz33VMgpJc4A6I877rh6TT311JQ88MADUcsAA00wwQSYIfrHHHNMPqAXvvNz6aWXNoTNhVNMMcXss89uFNpeeumlYi9uuOGGCSec0BTY76kMRitmRS8tv/zyJi5sxTbpvvDCC3/55ZcVN0eUc/xIg/XXX3/JJZeUypx++ukz02+++aYi3bPPPtsok0wyCeONNdFEE91+++1RYnG032KLLbxCbrHz3nvvVQ+g1afLpJNOqhfXlUH32msvr0zNCsc20bP6rsA9btK4ZmERLJpXED9e2VPYrWbvvff23S24Ifxkaj7H2muvPdZYY3Hb88wzj1XabLPNMhcb6Mcff3xTYBjzfKxrrrlm0GP3Jk2aNPBPntBcNt10Uy0D/n6Kzgv8mQRpgRv7DaQBtqWWgnkLLrgg7gmNRQTQKi+55BKjS0foYgpsww2p7+NKoq0C90Ttk08+Oau0hMnGXWqppbiVmI1aUuMTLLHEEjPMMEP8yMYbbxyb0SjA1qeZddZZLQJMfv3119V7GsK6qWetLrvuuusQuhCsSQvcG3bDrHXWWce2PEAAndXjKoCU7Ce6pbBbPewTI5YGVIS32IJUnHfeeTnVlJ+VclWJ94VE2YaI9gCFhxxySNrcfPPNGmy77bZfffVV7EHM4GCUUcJabr311l9//bWftk7usssuGqsfMA0K/pgtGgaLKoXsYFroGePfeecd8egiiyyC+Q5fstBCCwlGhb9RPv/882OY7OxMgIt9N5ZzVIH+rnWzXLA4DULD4HISEKfZgw8+6CdGJP4Js26R+SFcu58vv/wyXwXxETzJUHtSFYUxwBM1Ze5xIV7hiowV5VB75ZVX7g3cLQWe3tzLYEySt+UOOT8ukJtE7Yfy4Vk5HoU0WHqYKJQ8/fTT8Q1uS1Dh6OpWW23F+I6ExmBKkyZNGvibgnAZ05HuIlr8ixEZD9zExCE70uuss86iDTEUeIeBfhqC2ewkgn41MPaDDz5Il+effx5uK/R3JRW4Z93Y/9BDDxVIyqyK48NbEcrRUrBUmjcTOeqoo4xbFM8JJ5zg5xdffJGfWQpBvCGSgki2QRsaylX9xKVJC9wbdqMB8CLK2NZ5551XQe4PASM+k7Yr7C5EAF7wTjrVE8TD5U4oD7gDvnAqCAOgptwJ6IloMzrE5Dw6s4clu+22G4q3y2zxrqxl7OnCbkcqReqCyzJgo402ojxxsPSucYFv9ZXtVSM/G1MRHigN5RIEUm1V7GXczR1N0lnpWCfzYoCOsLV6JetqVZFbWUNDb7fddqUcHIueV111VeVeWWyxxWaccUaFmtrBBx9MgwA92ipwz2Tjuuz88ZluvfVWy2LfS+ovuOACjdUX6XXFFVeoefzxx9MAIUQ6xxKm84JdHzGeftBdRZMmTRr42/BDszi4MNDuFDW2TSrb91hvS1DsUqbRwM4pp5ySDWW/zK14WtjdaxuX0ceVaJzAvVQR2QwraUltwrH7Bfses3FkOm655ZYxI57F+uNolEO6a6BZjZUrgF588UXlEjrJkKNRmrTrINvRZDsUwQE6AXmMaYCSwLcagGmojWZw8kYu7++HSRiF0lCFlAGQeBTPoWLAw+/oEMSGtCyqW4M6LJ+WhtOAO0EP4xiIMujPRQe69N6chVbBhVTSULaUeQrBWa8o8ZZhOUJkjuozO6bON998ZSqZaqqpbH3JRDxLUsNOAW7ap4utosxDwGggm4mQltws41dccUVtYnxM4i3qQBieCe9l2VNTtxZkXuguygvcM5YViPFqOj8T3kvi2Lj+3+IbedKA8k+D6JRNpiTrwL9Sle9YklFoC33lrbmYRWTRRRelJ3MZXGnSpEkDf2NhoMIvGJqSxRdfXEFe0Vs5Q13shIk2TxMRbaOZYqrJ0q9NgR7bRPZypLGnruEaoSuhpFZDObkCrkQoLy2A3Hn44YcNR8p4Q0e5jtKhUg18QefCFj2UuRhUJpkBmQtkllCF0lavRh8i0mT0ITvzJoAAxKB7Yd8ee+wBgKQ+C4IjQMEuFLQxHheCywyKL0XGBUAlVJGUg0T5WZWdOnWnCiRVL8akIL8pqwhtO+3cYIMNOIOUeycC8owYVRp0TsFAGJS8qjb2tcf+srAkNf1BkELP0uanQcMqQV6LGeNjCX8ghkYylfFWJpMleimzp9MYXWqVOIYYU2NlUr23CO+www58sKNLEgLysOqljOnv+uIkQD8g1nd+QXPheMxFZXU3F/v7ay6DJU2aNGngbyyqOsEfQtKMs/A26GeOncbQlj3lEWnJmnIA2eGlsi2V/V2JGL13be3JkcRwHmCTTTZBkVgQ2yntRAq2F8NVjqD8xfBWlWGsWmmllTi7skq2We6iA7p/6tKkBe5NAgrAyI5w5/cTcSZwLHoYfwzTnXyXVMUN2ClhYyL41qb3rlwyMvFczlCiCjSmqlAsJgFuGUa5wj4298J3V31IlBxURZ9DPWQJ85gtHYwxUp85/lCOKrlX+jNlC+JnaP4cMFWwmAN2710cXJG+ePcsRdQGzWE9I22OzCtj8TFJ++rVpcdHsQOen+CWMk1uxjQVftC/B32zmDo6Jov96jeXwZMmTZo08C/4TVJRUGsDujLGoZ5BsD4S1CUKovY33ngj7EbZptDHlQD8Xh9kT5HTU3YH1Svz+kFo3JWLyCguBLO1sgFAk6HrfRtwkyJr/TEOe/UwtfWKJFuauDAo6enWwtDGnsUfe6VlJ9FChjeiZhKv0nzOJIFafQklyUi6IsBxHydvungauwD7TmTgmtVXX11fRzNjnqeNj6BQBnOE3XuFKhS4uwXsPqSNExIfu3pFdtVeydwY42ipGLprMw+I7/I98QdcCzLeDniHq8TlzKulULZj/u6773YaySqF7MGK+UbSwfk0pS0FXqQ+mWtkKNGrj9vrrDG6MD0GRIO0uIxzTpKV8LWZy+BKkyZNGvjDN6+c9WcABKPHYZ6cLvWUHoSx++yzT1cvucQuqzKKJ9vefvttu/kpLNso+aGuxNQ6Md/+HIeObN83zT5oTAqNa2GD4Q4Te7odcoC5DBFp0gL3JhCkthja+uzSLkFnwKszYYodcfzFhQCOB8n0AVZHbaCSIC8NxKzAy6kamzTgHU4C4ghVo6p3/3Qg0k2FdunZLO6OW/efuCl8rbXWcptv4kVHjpyM3HHHHR2ZFxa78tZuQjcnpu/wJlKiTdFCUBXkyQWDPPs45SuJHYfopXgjLXXvo7BXUNHUIqIE3Dad46X4s+yw3H333e0kcb+k4Zz+dC2MgpnaoxnDahGIgp9uXjccesZNAtYW28TruIUm98kgWqySsTgznsOlwhrkXjZ9IzFJ2pQGp1eNu8022zi0ihnCD9VY1Zh0dSf2r7sswsFTp6C4Ve15dJZYPSfGfAUGuEnNASwetz7EYEmTJk0a+BsLSc8AOsG7LfuoDVfoujDRBG3vAVkAUx5SkK2By7vsurGJJZBY61BMCg3LLbccGyC59rQpm50GXMbwXEk49U5t8hgOHdlYaArO8irYn2PdNOtakJL60611DspASCKAbGEdbHXzjNPAXI+ne8lkEuyKzJHZjPsTlyZtq0wTh/dDsRSIwJTiUeTjAHR+OqcPxWA3shnX68C+eM7b4IWTnS7lxe9iFNS4AwtcUq5x14h02hEY/MVPi0El/lwljplQTzNoC1cBlfwdIvAE9VjlFnN3IAxvm6NZdGUPDZ2Tl2kvkLWHBFPCQwh5obwMZulRY5dkl51V7s3MCoWBu/t9hchocm7PBlBbD4Pd+rpAgHvIXwbBl9iDyG3Ybp6paYBNqbNfei2wwAIWzR/UsBqw23QE6KzyypUy/pSJgdhs0SixII46JVFLSltuOrOMvJQaayVX69aC+ogqNa5cbRkT8iwXxjujxpWywbfwdWgzNe7B0nF7+soPsI0fHfStMk2aNGngLxQGxQwAX9gNVAVToWUOm3oKduGVBgJuM4WrTvCHtQlcZx0qcGe5HKO7yHAubMB2o1QcBg3493clpc24LrvkHSzObbfdZmoMcOwVRFNSQ4dKL8kfysjCMtJEjILVAtGPPfaYVcWbCN8BsskaggFIIruYYvz/sXcmpxHEQBS1Dw6iQ/DFCXQGBgfju7GTcAhOwzl537fPFDyENFTPqA+N6f8Og5BKtQmqhUCag1VhfleJ0TUgFc2tQ/EaSTx7wknty4Zo6wBAFVMNUK152KBG3GtEAKQzHiIs1cqKPNH9HsQYikeywskECUis7JFpOYA2PJQh/dJPHogLo6VC5HnXTLvbePhMcxnd6jyaecJMSVB6q05e6sXDSpts4RKdaItOYidjGo1VYE3DFpmvVpnUSUOyEAtijHHxp2zq+DmeeQm74TZUIUhhzEpMlNVVvhHI5Kek1UaAOFB+lZSQMnDyVuUnEitn2lhwoMb4HXdjfgTlbENbQ/vUUnq4hdMOUSsZ7TM0UxUb9/gnESlBJ40+59usVm1Rqko8RJJ23ypPx7IQxhgXf+zqvFmn/nl5b7Xh86Rv1NKOTwmStPcCo0ksrMiKMN64uyIno4k87Wq0Fch1IpkPiR0Dya0H3dPjdJ//7tZ5CZK7O498XyrmLNPu00USy4IYY1z82deO4zgMA9vrjmos5lfjycy0o/tOn+4H4427MYY6q7tQuiylm1IcfiyCMca4IOsqkS7uFz3/DmO8cTfGGGOMMcYbd2NMPHa2lHVjjDFcueFC6now5lCNP/bu2AgAAARhoO4/tJbugO8G5nKUEFOP83X4996PJUAGPiBAWgRkOB8QcMveecZYVQVx/ClLE5COIgGRIkiJZinyZemCLiUkSw2GTmCBFSK9uHSNCygk4BoEgqBL7SWhKSUUIZDdpfcILE2aDRX7LzvJ1fgSP7/75v//8LL37Kf7vzNz5sycmfkr/LDClKj6dEeIroJ3fv1RgaiQsiGNkND6LCD5//5RsuE+5cEzIdE9fIK/hYiD3rReQGvYqM7iTi0+vW9p3eiZAeu4DAl+GQg0QgiJ0PL78OFDzxaMTt4wIB8lGDTh7cgd3RseM+5ZDGj/H28aIccd+QaMKOvatSvD4ZlSmZqaivtu607UG9PGcM369etXqlQJEphsl5uby7rbzAMD5xhBx6B+h8d02+eY+WfywNA+eHDlCdlrMpqRcYzP5KNt27bnzp2LcQbktTMpk2E6WDBmCe/atcuE2Y/OHjlyhLHEDEuGgerVqzN92ecBxnLmRF4QBsaLnj9/3ha9qQNzrDt16lS5cmXMONNbGZ7qhwcYADhyzIjl9U0jGEMbJxohx93kuGXLltWqVTt8+DB9+sqWLduvXz/7l5/oMq4JQ5VPnDgBCU2aNEHQ7969ayLuzeJfvHiRzQ8xYOC/N8fdXjYjI4NzC+O1z5w5s3PnToaWmy74Ocpeu3YNBhhsfvz48X379jEnvGbNmlQVq/lxzAotBy0KABhQj6M2dOhQRrifPXvWiRm38soFCxZ07tx5x44dWDCUFzYYvG/8OJQHwg3ly5cvXrw4m5orx91sFNt31apVsWCYLxggEJOTk2M8+JGBUaNGoQWc5wm7YBn4e/ny5aHXCDnuaLgFKviiR48etcUVK1bwmJeX58RTiVZmdB4GaHBrW4I3i88pbsqUKcSu2AVt0VWc5saNGwkJCRs3bvTsBa5evRoVuHPnji1u376dx8uXL4du5/Njvlq3bt2iRYtgkUDMoEGDHOpvgAEDBpAuC/7lSn8JwHHYZgsrUqTIqVOnTEhcMTBs2DBiDc7LG2DAgrCGevXqdevWTY674cmwlxvv3buXcznZVStBa9q0aYECBQ4dOsS/eIw4ANEpfpFmA3fCeISEiCfwrXllostEWydNmvTgwQNvRegm7Zs3by5WrBg2jnBFnz59CFf8oywOYB+d1zfLYIt79uwpWbIkd2bUmiAGbTjmC6tF/IW0ITacS94sct9v//79Zsf8iK7lDEkNGTOcwB32DyHYNHjw4KlTp+K7s5dBiysGEHhkYM2aNVjvTZs29ezZc/To0aShXNlw28saNGiAIwcbLBKOhISGDRtGDELYM4zDhw+3yERQ3FO4cOG5c+d6izcHhHDLHzcFz9VVzoE3vXfvXrly5bZt28ZjYmIiZQ+s4xO4+vpjx47Fca9bty7Bib59++IV9e/f348kmDBY0J2Klw4dOrRq1ap27drcmVG4PWaja2gugrp06VIz4PySNCtVqhTKa6LrLWV0/fp13JeZM2faiqt3nzhxIv4Zj1u2bIGEkydPWkjOjzrcv3+fy/01atRo1KhRWloae1np0qWxYEaREx4Apcnt2rWzQHuVKlWgIrjrKITecR8yZAiFC8HnfPz4MYWJs2bN8ua4GxukFzmMrV271pTcla/GpfZevXrZChf9iVXYv1wJwLhx4xCA2bNn2+K6det4pLDJ5MGJxYcKDjBE7DAOBK4ogMnMzIxZF1CO+61bt3DcyQ4FjvuMGTNKlChB7NnbmdNyDtRl4ryyl9mKHxt++vRpCrSsIJUCZRz3K1eusO7Kcb99+3bBggUpygz6ydBpgByUt3JtziqNGzdu3779iBEjkpKSqFm6evVq6E9xctzNrI8fP75ChQr/bgCH0C9cuNCV425vunXrVry0efPmOfTaDxw4UKhQofXr1x88eBA/tU6dOoRa0Xx2QT8jXS1UyW6HgbOZgixyT3TOnDlO1MHEftGiRSgCRX62iFTwSDFMbO58ctzpeUdhBocrc1vRaM5dFStW5NFVNyReHDRv3pyjpqvuAoE31qZNG25MUYiJttJFhOMc9YjUmrtqE0c3FWz4wIEDrfMEi9OnTyeL7iQBBQm2nWEBrMOEgTQyvnuc2HAVpxJjRr3JLllwIjs7m00a782P82pyvHv3bl6cq4GWdvBWzbNq1Spy6zQiwMDxB34A5zc0n2JNswVOeFi5ciVicOnSJXMCKHjg5hhnOSeOu70jex6JZusBzAYAM+YXxiwJmrNDJygaiQQRd1Lk1Kq68lyNiuTkZIwYB2/zXbx18W/WrBmvz80QbDi1OqgtBWwTJkwwzfWjDpgvyxjbaJrJkyfDCfbcj+NOawE2sg0bNvDWdu83PT0d2QiIkuMe+vuRuCZIti327t2bHBNK7irefOzYMa53c0HIc1M5ki1Mavg2H/Tw7tKlCytIgqvNDwbY8EhD2eLixYsxf9wTddVZb/78+by19aMAWVlZPJKTidnDvNpBUlCO3GLMWbFSPETXPqifKh3KUajI9DwyjIFxWG/CcPxyzY8kKv3iLAnjSh2mTZtGphTn1aIPVPH5aZIWRNzLlCljpxd7JOIeV4d59XFfsmQJ4VXSKCQZqeqgd7WrPu5Wj8tWxw1vKjJ79OjBL9dmPLsp7H8E7bwxgMxbR9SiRYvSX4/STFKu1Pl5m9xBfpnuUpiClJQUzAIZOdqKqbAplj8ZjholiTgo3bt3ZwhDx44dcdntX04uuVGXgg2nHhH3naADBXm0GSB36tZToTgVQrj17rCPO4c3LBghSDzXWrVqEYBnHpO3O+4EXHBsuOYOCWzohNs5xSn4YojETciZsjySKURrAkvnx3Gnlot3T8sHlRwknRlo4nDgXACKHOio5bCQxd6XEDshzDFjxhBmdigGZtxpLDNy5EiqVG0MpxDjnwwnlUEzdAnj5Bms+wmyIqXEWRFX5k9hwAEm3Vt5rllvyz9cuHCBoANRZ1v0pg6QsGzZMsSA/CGpY582gRmC+DaQQJ/AmzdvButCJM4mEKl2QdBEmygL6NHoh4wE+e6CIMiCRW1k0StCJG6+McEJsqs+P+2vUXAu4kiC54QaWgADnsWAfQ4BMF0IhSQI1sfTvpdPN+U/Nhw2nMsDJPDr3ILZnTHPG1m8kiDHXRAEQRAEQRDkuAuCIAiCIAiCIMddtyRFgkgQCSJBiKUyRJEgEkSCSJDjLgiCIAiCIAiCilNp3hzSUjyRIBJEgkgQfspHPBMgEkSCSBAJctwtjZKXl8fQbH495sdFgkgQCSIh/CNvU/MRPEZDJIgEkSAShIRIXIC4GsPe+Y2EECJBJIgEkSA8evTo/wkQCSJBJIgE4W927mAFYBCGAejSrv//w7XbWMAPEA8WkpMKOfiuhd6HAHCivVD8WgCqKiK4+pQvV/MIQQhCEAIAd2/BPscaC93M5Gd55rV7hCAEIQgBf8xsp87Lrh0bAQyDMABMkRVotCAF01OoVrgM4N62/jgNwIkOHYEkgEldrJsRmNRuvASfg5tgmTWju9VPC26Cm+Al+FXmY9eKdRMHouD4vEqILctEkZDoIf+B/+COikSiQvxAFIkCKmq6tGmQzFUgUaYhFf4IoHZBA5aJgQgn93L7BffupF3sm4ICaZ/WszOz+9ZWfje2Xq/7/T6dM+Q/f1pBCLHdbsMwbDabxWKROjOqeYbXhLBtXF0xR5vmFwm7Xfj4+EVCmp4AgzGHzYZ+M06C/tjvQW6QKlZiBxpxfQ0NbJRfO8hVs2271+uVSiUZldretdN90mj0Mwhmtv2Rph+MCpZlzedTAN1umiQJFYQW4GuGR8J0OpWvoCUJ7AjVEXwl8CNUNZQrIc92UI9vv+F5XqPRkDmJfwGD9gPlib9ararVar1eZ58zqIjrulEUUcHzNCosC7MZFou/EgngAhHAJ+HuDpeXIBYNI7ckqMftLTwPSQKOzfl2kFFLw3E8YjSCBsiGHfgN2Hg8Xi6XlUpFRqW2L7WFEK3W/fPzK/ADiAHGVD8BBwTEOvTbSjRD6+44DoA4juU+yI5QXcFXAj9C1UOhErJsB/0bMJr/ZDKp1Wq+78uczMI37lKRhUKhXC4Ph0PLspBjDAbodHBxgfd39g4KIXgLQSaB68L3s0yC/pAP3m7j4QFq8fKCKJLCgGrk1A5JkgRBQPEIQP996+amIMR3IZ5OJybnaQqCaSqXXIYjVD34SuBH6P8TRZbtoH8PdjgcAPxi525CpLrSMACfKqvwr7ttoxJNsspCDWQRHDSJQdIKLowLBUE0K00Wogbys0hQQmxRMGSEgK4MgmA0NG4HBsWNMZvgQtyIwVUyGtQZbDRqsK3+mY8+pBA3Ycy0p6rreSgut66L23X71Nev3z33TLaHU/NaPzFh5qWXXurAZ+lCo5Hq9chJETVixMf2r4SVp08q8WMMDqaentiPI+19EdpU/uD5IuQL8ozlX/1vv8XZ4zK2THAv/3UoM4cwL4LWJsv4xHPAD+KSh7/SEmqBm6ZKaIs0LAuUUImiHb8OrXwrMiaBRzGfhME9i8UT4kN2ZnCPYV2rtcRNvVotXjmpFLsIVKv5txDbMsE9zkvBr0Mug1ESU5up5Miety1DCVVCy5MoOlBk2okItFJSAQAA8PTBHQAAENwBAADBHQAABHcAAEBwBwAABHcAABDcAQAAwR0AABDcAQBAcB8dHW00GsPDw7GT2goAY2NjUcCjjI+MjMR+Kg2A2sSl9uq4x9+mNgEgtVcqlVqtlloGANWJS+0XL17ctWvXnj17rl69Gm/bqGEDILU/evTo6NGjH3300cDAQPN4KgeA6gSl9mPHjr3xxhuXLl06f/780qVLz549G82bFp8zA0BO53fv3n3rrbf2799/69atDz74YN26dTFtpmx2B6D2f6/4kdpv3769ffv2aLfv3bs3Dm7ZsmXr1q0///xz865rrv5PbIsDIDosU6ZMOXjw4JUrV3755Zc5c+bEXdNFixYdP378vffei/gelfyJMp6lVEltBEDHPffUT58+3Wg0Pvzww3gb1TzutP76668XLlxIKcVDTmlcZVyk/NhOnTo1tQAA7fZI7bFz8uTJ999/P1L7w4cPFy5cuHbt2hMnTsTx/K9ZZVy9Xm/j2fAAgvvly5fnzp07e/bsnMtffvnlqOxxMPdmmvdh/zUuMn004/PB4gB48ODBtWvXXnnllajMubC/+uqrP/30U5773jwY1TtqeBTwW7f+c+fO3TZbXxhAcI+ankN5d3d38+C0adMiuMfB/DZPlDxy5MiSJUtef/312K5ZsyaaOrmRUxYA9+7di1ze09OTS3q0YKKk379/vznNPYf7FStWRAGPp5hee+1v3377j5R649/bDIDlICOmR33PFT83ZmKGTBzMb/Md1W3btm3atCnfeL1+/Xrch82zaMoCIJfrRqMR22bDJWp1JPhmg2bmzJk//PBD1O2o8NOnz/ziix3ffHMnqnubZXcAwf2FF16I51Oj6Eehj97M4OBgFP04+HhXfta4ZpTPB4sDIIpz5PKbN2/m55TCjRs3FixYkO+LNtvwL774YvpDb2/U85HURgAE91zQ+/r6Pvvss1gLMm6hxttYETK6MsuXL8/99SeWI4jIPjQ0lEoDIK/bG82UmMd45syZqOTTp0+P47Gk7+rVq/MCA0+U8bzOTGzbbFUZAME9CnoU/WXLlq1atSpmwnz33XcxDzLWAI5FxKI3k5d4b0b85k5ILQPAOu67d++OMn7gwIENGzYcPnw4nkCN9cEeL915P2+z1HYATJXJaw4MDAzs3Llz/fr1kdQ3btwYSwLHwRav7ADkKY4rV66MFSH37dt36NCh+fPnR/d98eLFzebL5AEguMd23rx5p06dilUIosrPmDEjtQkAcv/l3Xff3bx5cywI1tvbmzvxZVM7ANWJu9kavZmurq5I7fnxptQ+AGT3mM4e25za834qCoDaxBX9kPN6mzZpAMyZae6nyQpAcM9atEMDgBoOYKoMAAAguAMAgOAOAAAI7gAAgOAOAACCOwAAILgDAACCOwAACO4AAIDgDgAA1FLLGB5OjUYaG0uVSiqiXk9lVavxM5T5MeKax5Wv1Ur+FuK8YXQ0FVevFz/1aKMxEhfk2cu/+hgGtVo9rkMeGEWMjPxvgwGihE6ZkkKHl1AlNEZCK6DRKJjrymfa+OARaarVSRrco8Y9/3xOLZ3r/v08ylMpg4PpuecKF9yurlaoNcVPXa3X45VKiWEQg6G9QGTWeIUOL6FKaPwxpRVy3Zw5qWPlSDNjRvr998kY3CO1P3qUPv88zZoVjYoynZKpU9Mnn8S2TKcktwfeeSd+x9Exil5jKqW/P5WS/2N67lwKBa9AjIGPP45tKiQ++ZSUzvX3nxvfefYXopLSWFyGPXs+ie34fqXUSPj++9j58ygGeZy8/Xbq68vjp0NLqBKa/4D29ZXsu5OD3MOHadeuGAmd20eo1dKPP6Y335yMwT1qTQT3r79OZe3YUbLWhKg18SqrUmmREV/sZve0aenAgZIfPYZDf/8/9+79e+nKu2M8uBfz1VcR3OOr8eeDAfI4id7Hp5+2TglVQsuPioIYGkpffplYvXryTpWp1SLBF+u49/TETvFSm9skxWY237sX93lj29Ezm+OD374d46HULNV8f62rPq4R75917KiMjY11d/cMDla6uwtehKe/2Y3ZhvV6Z5ZQJbQZ2Vuh3U6nz3+OWDs0FJFmUj+cWqkUC+5R64urVuNVMrg3n6zK16RT5UeE8zUpVelGc2QvFdyHhxu1WtmL4O8uT/l8f4hti5RQJbTDeTi1g+WHU9tvOUgAAEBwBwAAwR0AABDcAQAAwR0AAAR3AABAcAcAAAR3AAAQ3AEAAMEdAAAEdwAAQHAHAAAEdwAAENwBAADBHQAAENwBAEBwBwAABHcAAEBwBwAAwR0AABDcAQBAcAcAAAR3AABAcAcAAMEdAAAQ3AEAAMEdAAAEdwAAQHAHAADBHQAAENwBAADBHQAABHcAAEBwBwAABHcAABDcAQAAwR0AAEi15h5jY/nVoZ89b8syEprnrfwhPVv5jE93XiPB97FS+S979xci11n3AfyZ3TPZ/JnstgTbi9BL0dgLazQUpLVB1NemQqliFMQLa7XWIpgiNSnWVqkk9SIpRhsqWi9siyA0F61Vi4qvEVSipaIY/4AWbS2tpuluJtlsZmbn/e08vENI3lcaMfvM7Pl8OBxOhg1n5swzv/M9v3PmTEwF1+styyj/XoxCEW008rRC3wXBnWaz7BAvX2JiC5RlJORVp7TY7/c7nU7M0/LKa4xVp3KMhLzGqkpjp9Mptd7hEW+u5HUvoWg+LC5GPR9+NFYiwV1r7dix1O3GQh2ze37Vx48X/oQbCd1ujmtrpqenq6rqdrvL33GP7L5+/fqCTXcjIa9xdjYWxmlzRWxttWKLFSihVbW03jA3lz9EtS6hLCzEFOOwWGaNVbdahbPEmjVpejp/NIrVhHZbqLhgwd0JtRhemzbV/Rg9f8zyQhFGQlU1o8ju2HHb3/52a7fbqKp+ob1Oo9Vq5QUjodSOv9dbeiaTk2ksRONj+/Z04ECkhALPOYf1vXvTZZflpFLTEkqnk5rNpZGwZ0+ZkZCPFtavT0eOxLxMAyif+bnttnTrrfn5lGoFRgmN+QV7DoK7+B7DqyiMhGGBm5qensq1r87UhLVrx2lzrVoVSaXkOcPccZcSmJ8vPBLymbqypqaWpoKqyi5McL/wmYmytcZIaDZTp5ODSD+3jooo3Ws3EnLmGMfv8ubmd5E+a95u+UNU6xLKxESxkZA/uc1mnXfoueNeYOML7jIrRkLZu8poV4zZSFCRGo08lb+rTFm4OZjWQ0w15j7uAAAguAMAAII7AAAguAMAgOAOAAAI7gAAgOAOAACCOwAAILgDAACCOwAACO79fr/X63UGYiGNFQAWFxe73W7U8JhHSU8jAIAqXQCNRmNyYJjj45E0JgCYGDgzx+d/rjQAOu7/+Mc/Hn744dtvv33Xrl0/+clPIrWPS8MGgF6vd/jw4d27d3/yk5/8yle+8tJLL6201A4guOd0furUqTe/+c133nnnX/7yl6effnrr1q1R+kc/uwPQ7XZjfuDAgauvvvrJJ5988cUX9+/f//rXv/43v/lN7runcgCoLsQF7vfff3/k9WazGQ888sgjH/jAB7Zv3x5pPro4Z14/c9a8LAByZz3K9a9+9avLL788lk+fPn3FFVfs2LHjBz/4wbnVPs/DmF0SCSC456q9Zs2at7/97blzE4+8+93vrqoqujUR3HOVH/7xcD41NZVGA4Dgvnnz5pzII7VHfb7hhhu+/vWv51h/ZkDPC7lHE3V+LG9FAOAa92H5jv76j370o0jweTeQq3yO77Ozs38deO6555555pn8YHEALA7kOB7z733ve9F0z4/nMp6Xo3pHDY8CHt9ripI+ltfBAwjukdejpkfFn5ubu+mmm6Jbs2XLlngkHh9eQ/nAAw9Emr/yyivf8IbN//Vf18aV8fH/iqd3ACYGolZH0b7vvvueeuqpe+6556yLG0+cOBHXwUcZj/L+xje+8bHHHrvooovS2AFwO8h8LrXT6bzjHe+IUv7ggw8O+zTDFs7NN9/8/ve/P6f8559/9i1vuW5+vhd/Ujy7AxCpvdlsPvroo3F1+1e/+tU3velNw+ZLLubr1q07dOhQr9eLx2P54x//+Msvv5zGDoAvp4ZYuO6665599tm4p1hk9zODe16YGUgDU1NVfrA4ACKOV1X1/e9//z3vec+99977kY98JBoxkePP6spv3LgxL+eSPmbXuAMI7v2BKOjve9/7fv/738cJ1ksuuSSX+HP/ctibX1hYSCMAgNxZ/+lPf3r99dd/4QtfiF/kGH4D9dwynnvzkfJjPmb9FwD3cc9BPK5r//a3vx13Ap6fn4/4/sc//jGfQj3rrjJnSiMAQK89+izRc3nnO9+5bdu2G2+8MWr4H/7wh/hdjuFN3FdCGQcQ3HNnvd1uHzx4cMOGDbfccktcFvnWt771qquu+ta3vpV3CWlUAZDT+eOPPx7zX//61/HTS9dcc018CTXuMRC3hlxpP7sB4D7urVYrbhDW7XZz5yaFwc3dh19LHU0A5Etidu7cGd9JjTI+LOxRzPMPbqyQzjqA4D6M73GHgTSeAFg1kGoCwO0gz03zqSwAlHEAwX3M6jsAyjiAX04FAAAEdwAAENwBAADBHQAAENwBAEBwBwAABHcAAEBwBwAAwR0AABDcAQCAKo0MJibS5GQqq9NJNJuFVx0jobhOp/CqFxdrPRIajX9nCyihzWZMSigjUUJj6veXPst1k191p1P3WFVVUcZXdHAnkkpMyKzttqOX1GqN10go/y7Qbue0lMqCiK0bNqSaiy3QaNQ6VuVaNDGxQoO7RlGMrWuuSVu3xkKxhuvCQtq3b2leZ1NTaceOpXkR0Rvo9WIY5FFRrFMSY2Dv3pjH8tIjyyx/BH784xR6vZqOhNjys7PpwQdf6V5HCQ3btqW1a/OHSAktg1wwT51Ku3bl6lFfCwuxHfI2qWOsipVWVTp4MLXbKzG4E3EtBlbsdW6/PZV14ECpuFZeftWrV6fdu0dlVJQS1faOO0bkNFR9R8LJk+mhh17pJ1EJDVu3xqSEMhKZdc+ehFj1/PPpxIkaBHfneZvNMn3WubmlBf2So0fT9HTeJsUqzsRE+fO8c3Ml80evF2W3piMhr/GFF4S/821xDXvtSiiuc3O1YclY1e2mqor+Sw0ulfHNqlAouOf1kr/flrdJzb9ZVSa4Gwl5jVWVzpcSGlNWsISCrymLVbmGT0yM3+0gAQAAwR0AAAR3AABAcAcAAAR3AAAQ3AEAAMEdAAAQ3AEAQHAHAAAEdwAAENwBAADBHQAAENwBAEBwBwAABHcAAEBwBwAAwR0AABDcAQAAwR0AAAR3AABAcAcAAMEdAAAQ3AEAAMEdAAAEdwAAQHAHAAAEdwAAENwBAADBHQAABHcAAEBwBwAABHcAABDcAQAAwR0AABDcAQBAcAcAAAR3AAAgVWlkNBpLU1HlX3t/SYppmcUa87wRir4R8RyKi+eQpxoajgT+jZHgXej/r0KFtDEaxTymmu7IhnNvROkqWvq9EKsGVnZwr29UCp1Ofu1LubnZjHmZgltVqdut9RuRNZt13/XGFhgyEgp+HsdLY6BZ79HT6dT9sD+2QHGLi/Fk8ntR99ZDKWJVVMK89v94fK9GJy60WmlysqbjLCfmlBbm5ha63UZV9Ze/4MawbrcbMzOtyclGLMcjRZ7DwkJMhSvdsWP5AKaO2T2/6uPHNd1LjoS8xtnZWBin1H769Onjx493u93JKOXLK+8dpwbK7sguvnipnucSWjf5Va9fX75yrlmTpqfzjrVY9Wi3C4+BqamYiu3IxKper1dVVVTFVatWrcDg3uul1avToUNp48aaDrJut1NVzb1791522Z6qqmLPV2S312qtP3z4SMxjOR5Z9o2wVGd370579ix1OjudMs2JqLabNo3OICy518kL9VN+JOQMFIUxnklk4HEQgfXi/fv3P/HEE1G+lr965LK5c+fOXbt2dTqd6HUVOUnSaqUjR2re5c1tuLxQ5ixZuO22dOut+UCiVO8jqkfMyzyHvAPdsSPt2pV3rGn5iVWRo2K97Xb7E5/4xIrtuM/MLB0i19KwpTc/NzeXM3ShI6juzEx/3bqSG2HNmvKhLartKMBIWLt2nDrukZiPxUmKEnLZnJ+fL55Zo9+cUbzZXFBVjchph8LNULEqV6eVeo17PiVd00HW6eQ+wUS+QjT2f0V2e1XVLPgu5I2wuDgq13mLrZQbCeN6rUWUr+XvbIVcNicmJnxwMiW03y+YF0fqQv8y31kSq/KZwFhvHe4qE/Naf/05K3U7iILvQl71+GZWHL0Yh2XL14gUcxx1u7VOJla5jzsAAPgBJgAAQHAHAAAEdwAAENwBAADBHQAAENwBAEBwBwAABHcAAEBwBwAAwT2lxcXFTqfT6/XSGAKg2+1GGc+/3V0WAFW6kCYG0hgCIPJ6VVVpZQPQcY9ee8yfeeaZu++++7HHHhs+Mi4ApPZGo/GNb3zjrrvu+uc//5kfSSsMgODeH4iFW2655XOf+9wjjzwyRhUfgF6vF6n9hz/84Y033vj5z3/+73//e/EyDkB1gdrtk5OTX/rSl15++eWrrroqXy1zbrg/az46APTaT5w48alPfepDH/rQ448//q/LeJbGDoDgHqk9Svxzzz13zz33PPnkk1H3FxYW0jlirzCcT01NpdEAQG6+fPrTn37ta1/70Y9+NK6WGQb3c8t4s9mM+TJcCg/AxAVq1XzsYx+74YYbrrjiimPHjkVBP7dJMzs7+9eBiPhxKfzodGsApPZf/vKXDz/88H333Tc/P////VlU76jhUcBfeOGFOL+axguA4N7r9aLiP/roo4cPH7733nvjn9GniVAeJf7Mm4vF/IEHHti8efOVV14Z82uvvfbUqVPxH1M5AAyvfvnwhz+8Y8eOSy+99Pjx49GL6Q+c+TdxIc3VV18dBXzLli3Ro/nmN785LO9jAEBwj2oe4full17auXPn1772tYsuuij+2Ww2V61aFfE9n1QdnlG9+eabn3rqqV/84hcx/+53v7t69eqyt3sHIF/ruG/fvnXr1n32s5+NR6anp2PearWGNTwvxB8cOnQoCni0aZ5++ukPfvCDF/qCGQCq//gJ1u985zt/+tOfvvzlL8c51sjrv/vd7/785z9ff/31EeVf9apXRbjPRX9mYBjl84OlAJCbLzHfv39/9M7f9a53dTqdF198MSr5e9/73ptuuimud8+nVfPPdGzcuHG4gujUpDECILjnby/FmdO46e/c3Fw+nfrb3/42Gjavec1r/s8r3XOOH357tTgA7rzzzui/xNXtUZ+jsx5l/HWve90ll1ySK/ZZZTwifpT3mKcxAiC454J++cDwwZ///OevfvWrv/jFL+YSn//mrIWQigJgeLOvuHf78MGf/exnBw8ejBoewT1fSHPu32dp7AC4q0xU9s4Z4uL1KOixMBa/nApAtM+jaMc9A2IedwCLyxrjgpkV+DUkAPdxnxgYPpB/uSPf6Hf0AZCvbMx1+21ve1vc87HVauXCXhAAVbrA4vrINLYAhPh8Y5niAJhIF9hY/xQ2AHWp4QA67qP4jSUAlHEAHXcAAEBwBwAAwR0AABDcAQAAwR0AAAR3AABAcAcAAAR3AAAQ3AEAAMEdAACo0sjodrudTqff79fw57XjheeF5kCRnzSPLZ9XXdbERGyEmMq+HYmyY4GogouL6bywuLjYGUg1Njk5ORFltBAltN9PjcaobAHis1AqVoWqqqIordjgHsHx0ksvzRu3hoYvvOxe5+jRoxHfU1HtdmyEpQm7Xs6rItJqtZoDqSiM2w0bluJ7cbTb7YKxKq/3P34gXY1Iaj99+vRnPvOZmZmZHBxr66677kqhUMd9ampq9erVBXvtYdu2tHZtNI1Sr5eKWFhI+/Ytzetsairt2LE0L0jHfXY2PfjgK+q70xvUiyeeeOLkyZPR4oo9ZT177bEdtg7EQvyzSLM5iufevTFfWq7t/jxe/qlTKRTcAnrtg0Sxbe3atfmjsfwnAKuqOnjwYBw8rMDgHts3gvu+iEuh3qn97rvvHpFDqRJ7nRS2bl2ayjpwoL57nfyq4/Bt9+5UFidPpoceekWDkHw++r8Haj9sTkZwjw0yCO4FRGC9445UFuTxv3UglfP888+fOHFixV4qU1VVJPh6dtybzWY+pdIZiH8WfCald8DFeu25XTQ3JyotbYGjR9P0dN4miSJD8YUXznso6rEN02pt9yOtVqv4wf+GDVFIR6L34WpDh/TRay/1vc2ItXEgHXVpJX85NV+wUdcPeSdX3jyv8a43ppJpKW97ms2YSgZ3wb2q0nlhcSDVWAT3EdgC+XtK5YM7TAyU7UeP3e0gAQAAwR0AAAR3AABAcAcAAAR3AAAQ3AEAAMEdAAAQ3AEAQHAHAAAEdwAAENwBAADBHQAAENwBAEBwBwAABHcAAEBwBwAAwR0AABDcAQAAwR0AAAR3AABAcAcAAMEdAAAQ3AEAAMEdAAAEdwAAQHAHAAAEdwAAENwBAADBHQAABHcAAEBwBwAABHcAABDcAQAAwR0AABDcAQBAcAcAAAR3AAAgVekM0O/3U1GNRiOVFk8hTzUUrzrP+0tSTGWN12ikMVC8fCnj2bKvN89TIxQtoQYCgjuSSo10Oqm2mTW/6m43VdV573qh0+mkopiYmIgy3mw2Y17ksL+qooAsSwkFwR19mna7HfNSvbpY9dRA2ZbzxRenqoqFOmb3/Kpbrf7sbLvX6w+aZ/VtHrdarVeeflSPCIuxxbrdbmy05V97rHRhINXb/Pz83NxcvAtVVaXllctFu92YmWlNTjZyMSnyHGIUxLRiIbhD3u1Fat+0adPx48dzhk7Lq6qq2Nns3Llz165d0beLEFCkXdRqpSNHYoPUfSRs2bKp3R6OhHrJr3r9+vVHjhyJed4m6V/i2LFj27dvP3DgQETGycnJtLxyTt29e/eePXuidNSz8Z9f9b59++6///4igzaXi1Zr/eHDR2Je4jnks4Vp9+60Z09qNtMKHwgI7ojvsfc9efJkwbQU7aLiLef161PNTU5Gx/3YiRMnU41FFjyvgxbVY9WqVXGcU7D1sGbNGmW8+GmHXq87M9Nfty5lRTruozAQ4H/Yu5fQKLIojsO3Kt0mtkxGIxriAxR1JQoqIhEVN0ERlWzEnbjwsXIRdaMLW1GMGHWhK3EhvnA5CyEigk5QfOALd7pIlCGCIolGYsyzM3/qQDHjzGZg2mMlv49QFBWl4aY5/c+pcys/IrgDNhnp0me1Plkcx8ENG5vsY8+6VvlIxnHH/b/e88FowmVIw+7RlUql4IetSlYuVDrSGfcoctihpEb7+HojgOAOPnrtxOml/fmPRfg/VKcM7wQHDu9GIqPxfZoNZdz3E8RqiOjoUr68ATzHHQAAACC4AwAAACC4AwAAACC4AwAAAAR3AAAAAAR3AAAAAAR3AAAAgOAOAAAAgOAOAAAAgOAOAAAAENzTP308MjIyNDQ0PDxcKpVCpgAAVLpVwFXGVcxV0oMrAECufOU+juOKRMgaAICVcQnGGwAgLlOvXbV+YGDg4sWLe/fuPXbs2JMnT0KmAACp/fnz58Vicf/+/ZcuXfry5YuV9+AEABCXo9xHUfTmzZtly5adPn36w4cPd+/ebWlpUbnPStEHAFL70aNHV61a9fjx487OThXz+/fv27eCEwBArkxFf8uWLXPnzr1x40ZI9PX12YkyfdqV/+7oDgCgcXaNOKp6Hzp06MGDB/X19VbVNemuE5t+/K6MGyvvmQEAdNxV8eM4bmtr0w3Wq1ev6tZqR0eHrhcKhfB3UUL/WMfKysrgDQBgg446Ubt9586dSu2q4V1dXbr4z0IdJfL5vI65XC5kCwAQ3K0Bo+BeU1Nz6tSpRYsWqe4vXrz44cOHquzq2aT/pqen54/Eu3fv3r59axd9AQBUq799+/b69WuF9YaGhpUrV86fP3/fvn3quI8m0nEaVW/VcBXwjx8/qqRnchsrALA5tb+/v7u7+86dO62trU+fPp0xY8bGjRt1RWVdRV8PFwshnD9/funSpStWrFiyZMm6dev0X9xnZgAA1lhRLr9w4cLmzZtfvXqlxwycOXNGvRjrv1ih/vr16+rVq1XGly9frh1NGq2ZPHlyyBwAILhbQD979uzChQtnz5595coVpXaFeJulsTuqu3fv1jiNtj29ePHi1q1bVVVV6QS8IwCApmJ6e3vXr1+/Z88exfHGxkaNzajbosKezrhPmjTp3r17KuN6aNizZ882bdr0+fPnkDkAwObUadOmKbvX1taqua4srrqvEq+RmL/uT/01kX5I2EV3AIBfEnPmzNF4zODg4IQJExYsWHD9+nXrvKQNmpkzZ9q5lXR9N2QIABDcLX+vWbNG44/t7e2zZs3Sxffv3+um6rx58/71qTK6oie+B28AABuGUTrXgyDVTc8ndP3ly5d1dXW6bqXb2LkaNLquY8b6LwBAcNddVBV9zTtu2LBh+/btly9fLhQKu3bt0t4mXbHnFaQRPz2R4A0AkMbxgwcPaoS9WCxu27bt9u3b165d08i77UlVnU+rtx1NyBwAYMZd5Vt1X1VeffetW7dqW+r06dNv3rypBP+TT7EDAJTLVcPVcf8toRP99aVz587t2LHDvhvGDABgxt2iueba9SeyNSGj9oxmJdOpmPBzAwBY/6Ux8enTp+rqauX1MVjDAYDgbkYT2pOanmfoEb8AQHa3P6E6ZcoUXbDzMCYBAME9SiivZ3eKHQCYmUnPw1gFAAR385PmdQAANRwA2JwKAAAAgOAOAAAAENwBAAAAENwBAAAAENwBAAAAgjsAAAAAgjsAAAAAgjsAAABAcAcAAABAcAcAAACQS8+Qz+d9XzqOY/9FSERRNDo66rYI416pVBoZGQke9HOPomh4eDiXy3m9E4xWQOsQXA0lbE1cfgohqwAShf9LDw3pK7iA6rdiVKmkVENwL+eHtO9L9/b2BleKC11dXf6LwARbIvipqanp7u4O45ui89SpU4Of2tpaxfeQSQCJwv+lleD15QK28oVC6OsrT3BHZWVlU1OTjsFDRUWF+otr1671bTlXVVUdP358YGDA+qzOizAuaQW0Dr8nbEGCk2Kx6HjDIY5jrUBbW5tOfnzf3d78/f39Bw4c8KoJiuw9PT2Dg4OZy+4AVDaamnQMvg4fDl6gXnsuFx49CvX1BPf/m4VUZdbm5ubgzfKryyLYby9KKtlehOwPySivt7a2trS0uN9+Ca5Onjyp4K7V8BqY0W+wJ06cCK4mTpyYoeExAFGk4qkuWGhu9k/tR44Ed2hoILiXc0qkurrapku90mocx9mdF/JfhOz7k507uGkcCAMwaiBmQbIjcaYCCqAqTnCCCsiJC73QAEWkAM67CGEpCgay/2aEK9jwY/k9RSNrL7OxPdZny6FomqbeisORch/btm28KhNjynKIbx3fveu6Kf/uJXZ73LGM8XE7EO0er53O57ERKZ/wqkxdD6/KpL20w2xWrdfxPG5n4U69VUrFL3TJfe5ekj0r3MuPUxOXQ0xd7t8mex9bDsTh4WE1PkCJ5hLuObOXgi9jBsqPU8f25yABAADhDgAAwh0AABDuAACAcAcAAOEOAAAIdwAAQLgDAIBwBwAAhDsAAAh3AABAuAMAAMIdAACEOwAAINwBAADhDgAAwh0AABDuAACAcA8AACDcAQAA4Q4AAMIdAAAQ7gAAgHAHAADhDgAACHcAAEC4AwCAcAcAAIQ7AAAIdwAAQLgDAADCHQAAhDsAACDcAQAA4Q4AAMIdAAAQ7gAAwKz6Mfa2UuYt4+ZL4h6IkeEQTPlA7H2Z4HIY5h0dNl9S5k1cOJSdP4xsNuWTNvVPkB41u60J4Z5YCTG+v7/PZrOfcMXHIfj8/IwK6fs+xgkuh7quY979/f1qVNjbisOXdc7ERtbCYbh6xP6PkbqObI1P2tTlCqqbhfsOW61pmoODg+8vlfJwMWZ/eXn5+PiI7fiXlP/Dr61q2tZbZYcknoq5Nw/Hx8fz+TzSOQJ6gsuh3DasVqtqVFT729vb6+trHL64kqecM7GRtXAYrh5t28aGx+3Pz3FaxkZCu8e8sRpWq/zH7U3zb0xR9vx6HR/hvgPRB0dHR4+Pj6enp98fCmXGruvOz89jTEnGkmjX19c3Nzd935dHVlNTvvjd3d1isUhs1rZtl8tljCnNWg791dXVxcVFxqmYvxyGJRn3byN6esfJycn9/f3Dw0NEc3RbyiX08vLy6elpiHgSn32UjWkme+i66uwsv1lD3+cke2RUVPtyWbVt5t3L7W21WFR1/f/2g3D/y969hDa1rnEc/tqukpzai1qxgk7EgQMdqCA68FIFBTsQ4aCCIDjyhoj1hg1e0oHWa3sGgojYiTpRoSCoM6mo4A1xpOJE3SDaQVqbvePuNTl/1gcL2Z3u+mYlv4fyEYMQmma9+ef9Lvk1MzU0NKhZ4oyoRaQWYy6XM2xU0GIUPQnZbNaw467YYT7PHs29xOhyAHFNTYeBgQHDEuo77lH3BzCM73/+SU3Q53lXU2PYcdf0NUtlJj8wWXXco0W9JpFRfVZ1FlnUK3oSlAD8E2LScddDx3GHrv3lYP8M2INKmZ8rMymhvp6bT1qCD07sy/R0URYKPkMbPLTKQD7vfgdOldFo8tDR7ljDIxEEnGcy4SqIx+UARKWMU2VA0x3R9lyNFg/NOe4AAAAAzIM7AAAAAII7AAAAAII7AAAAQHAHAAAAQHAHAAAAQHAHAAAACO4AAAAACO4AAAAACO4AAAAAwX18fHw0pBsuVgAA+Xx+bGxMNVxjoVBw1gAAgZs0VSEXUtGvqKhwMQEApPbKkJOSL+MAQMf93r17hw8fTqVST548UbmPRcMGAKByrcj+6dOnc+fOHThw4MqVK/39/SWY2gGA4K6Kr3nVLVu2bN68+ePHjy9evFi1atXJkydV9NXCcUUMAKDFjSrXvb29ixcvvnXr1vfv38+cObN06dLPnz/7TrwzAgAI/vWKr+UxT58+vXPnTk9Pz6ZNm3RnOp1ub2/ft2/fzJkz/fSrz/f/GIsBANBr140TJ07Mmzfv9evXCvGZTGbu3Ln/C0U1/NcC7rk4AgDWuOdyOY3q0LjQkiVLNI6MjLhfqKMTjYlEwhUNAMDg4ODChQt9iW5sbFTbJZvNOplQxqurqzUGQeDiBQBYKqNOjPouzc3Nq1evPnLkyNu3b589e3bp0qUdO3bMmTPH93I0+neFP0Jfv37VDKx5twYAEK1pPHbsmOZOb9++rRWPXV1d6rzs2bPHF/lo96qqt2q4CnhfX9+PHz9cvAAAwd1vQp0yZcrFixfv3r27du3adevWffny5cKFC9G8qlbAa7x69ao68cuWLdO4YcOGoaEhfwSNIQCAb69s27ZtzZo1W7du1SalgwcPaq2jJlGj5oufWV25cqUKuO5ftGjRjRs3ovIeDwBAcPfLH9+9e7d+/fq2tja1at6/f6+yvmLFioGBAR/r/Yzqrl273rx5o62rGh8+fJhMJm2PewcA+FCuWr179+7nz5+/fPnyw4cP9+/f7+jouHz5su73W1edc2rQ6MQwFfBXr15pcnX79u2TvWAGABBMRnC/fv262ufakKp7ZsyYoea6lkgqnauFo36Mr+wNIRfSPeYHjQEAfGrXuhfV7WvXrqntoqLd0tKi1rsmTtV3V233B7qr1M+ePdtJaOrUqS6OAIDNqdE6SE9132f6iX0d/wYwPDzsigkAQEU7aqJHh8lMLOO+HaPRxQsAENx941wnuHd2du7cuXP//v1aB3no0KFp06Zpsbv/OtXof0Y3xBUHACCvq32+ceNG7U+dNWvWggULHj16pO57Op2Ozvz9x+FgnosXACC4+1nU5cuXa2fq6dOnfVifP3/+gwcPmpqaivkbswEAvkSrVnd3dx89enTv3r1K6olEQl+il0qlolPeSwEAENyjg2X+G9KZj6rydXV1ur/4UzsAIDq7XbuVhkK1tbWluesUAAjuvu776VS//bQg8enTAAAKoWTIr5BRDS/B5gsAENyjNTNxXMUOAKgIRWW8lL9nAwAI7lKkYR0AQBkHAL6ASQAAAAAQ3AEAAACCOwAAAACCOwAAAACCOwAAAEBwBwAAAEBwBwAAAEBwBwAAAAjuAAAAAAjuAAAAAAJXNMbGxkZHRwuFwm/+km3/iHr0IAiqq6t1W/e430uP60KjIVeW/C+ez+edtVExeCnaiy4H/7eAvzBNqBD5ayFeYF7Aq6qqKisrnTWoeNgaH3eUEF0KquL6MREEijSlG9yVFZqamgzfJqdPn97f329e66tD5RySamtrzV+KjY2NrrzpGfAfWghh9hcFYljHABUPc/jrL5Vx/Vi+BiorSzG4KyKMjIwcP368oaHBd7tNnDp1yllLp9OuXKm/qEZRb29v2CoYN2k2axwaGmpra0skEq6MDQ8P63mInpPypNdAa2urRquqODg42N3dTd89RlNVunA6Ozs1mszcqteuytkc0g390xmBykZrq0arN1OFRaf30sePdaNM++4+Lre0uJoaXRoKFQZ/hSBwPT368FCKwV1xTcG9q6vLvPKap/b29naHMMQbZtazZ8+6MgafupLJZEdHh7Pz8+fPmzdvxuizE/RxN5VKOTF92Si4q4QaBXeqh7KESyadioet8+cV3JVZyzS4+5d/c7N+nKFv31wuV7pLZYIgUIJXejZ5k66rq9NSGY0mK5s1I+8nWDVKOa8wVqPIp3bzyW5WidBDzWQy9fX1Vhtv+vr64pXa4RfaZbNZq71SumzNVxsKFGQyGVdfrxvKGAYrNBQo1OhFPq9QYbVv0wWB+i/qTf+ffTvWdRCGAShqo/T/P7jFKSIS+xtoeHCOGNiMYgndJbe+nJqZs8L9uJw6wn1KLI5aukgzaVZ47WaFe2st/hveuxHuU6ZXVczHuBM5wn3O9GUJlmXaOWRGa6dMP3uxAACAcAcAAOEOAAAIdwAAQLgDAIBwBwAAhDsAACDcAQBAuAMAAMIdAACEOwAAINwBAADhDgAAwh0AABDuAACAcAcAAOEOAAAIdwAAQLgDAIBwBwAAhDsAAAh3AABAuAMAAMIdAACEOwAAINwBAADhDgAAwh0AABDuAAAg3AEAAOEOAAAIdwAAEO4AAIBwBwAAhDsAAAh3AABAuAMAANGONz673ntm/n50ZlaVLYB/wpi4rmv8BVX12Y19/X5l2+hrHcLzbBvIjKpLfMn29B5jDw7hgUl5/BluGO5f9swfR0EoCOMsLiHwQkALTmC8h17AGhM7r0HHNayMdph4AUoOAfECNpAQQBP3z7dOseUmL4Zh5X0FGSmYeT/D9+YNWNX7QwiwSGlGeIJcdiQ1DGMymXARoMpd10WASrqE0H8pCArCAD0B8n0fNfyX9ms0GoG2tP/QQYWeI2GAZJtkoYweiE0ElfzuZbwQBilat+siABC0Td33M0j6k5rMQ0EYrMiIhBDX6/XJeOEOvLMBXdfP5/N0OmWnvFqtTNNESXDP7jnYtp0kSZZldEdJSWmwngBbhumXZRnHcZ7nsEeyyj4fkDabzXa71Vg1m80Wi0XTNN2zwl+GpLfb7XA4vBiE7kk6joOgqiqJ9w5jZtvWkkTLMnb30EwT9aCLZYCgILALpo36j8fjfD7f7/fkk6/QuNNHhMvlEkWRYRgwGtyRmPRgh9vtduv12vM8TCzkKimKAldGFEIIy7LoJwOEvkpBUBAG6Am0UQkhwjD0ff/P763s8xc0rGmaomDQlptcnE4nxMvlUrrpbNu2rmtGUEg9Ho+pAH4I3OKCgFdWCM2yeI8fWlHgyghBQeCX/hAO0kEQPHP48s0OHdoAAAJRDD3BCn9GBNN3giNBswCktXWvvwhIAvSDiSCCCCLYPIkgwjpdlwgi2Kiqzb4d4wAMglAAbUsn7n9OJgaGohId3I2D2v8HQ5zIG0FXiLsHwPDMXkSijpOZ2821eYAABCAAIdomol1mpblm+LGNmbX6lz80OgIQVPUMBCAA4a6ZvC4o7dqxCsIwFIbRgm5Zkhf03fMKGVxCWkOELoKjePWc4YKO39J/yBHc82tRa805z3v+E4IIIoggAr33eW/L+fOVCCKIIALXn3nSWkqZdwtIBBFEEIGU0vYPRBBBBBG8cR9jtNbmPQISQQQRROC+vAkggggiiMCnhjsAAGC4h37GKoIIIojAvhxfQgQRRBDBcAcAAAx3AAAw3AEAAMMdAAAw3AEAwHAHAAAMdwAA4AHjDU0WOGjZzgAAAABJRU5ErkJggg==\" style=\"width:90%;max-width:800px;\"></center>\n",
    "\n",
    "**a)**\n",
    "Écrire une fonction `generate(n,p)` qui génère une matrice de booléens de taille $n \\times n$ telle que chaque entrée ait probabilité $p$ d'être juste et $1-p$ d'être fausse.\n",
    "\n",
    "**Indication.**\n",
    "La fonction `random.binomial` de Numpy peut être utile.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**b)**\n",
    "Définir une fonction `fill(isopen)` qui pour une matrice de booléens donnée renvoie une autre matrice de booléens avec les entrées envahies par l'eau.\n",
    "\n",
    "**Indication.**\n",
    "Définir une matrice de booléens `isfull` pour stocker si une entrée est remplie par l'eau ou pas, puis définir une fonction récursive `flow(isopen, isfull, i, j)` permettant d'envahir toutes les entrées possibles à partir de $(i,j)$.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**c)**\n",
    "À l'aide de Matplotlib représenter le remplissage de différentes matrices générées aléatoirement.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**d)**\n",
    "Définir une fonction `percolate(isopen)` permettant de déterminer si une matrice de booléens percole ou non.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**e)**\n",
    "<font color=\"red\">!!</font> Calculer le temps nécessaire pour déterminer si une matrice de taille $50 \\times 50$ avec $p=0.9$ percole ou non. Lire la documentation du module Numba pour réduire le temps de calcul en compilant une des fonctions: <https://numba.pydata.org/>.\n",
    "\n",
    "**Indication.**\n",
    "La fonction qui est la plus utilisée est la fonction récursive, donc c'est celle qu'il faut optimiser en la compilant.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**f)**\n",
    "En faisant des statistiques, déterminer la probabilité qu'une matrice aléatoire booléenne de taille $n \\times n$ avec probabilité $p$ percole. Étudier cette probabilité en fonction de $p$ et de $n$.\n",
    "\n",
    "**Indication.**\n",
    "Faire le graphique de cette probabilité de percolation en fonction de $p$ pour différentes valeurs de $n$.\n",
    "\n",
    "**Réponse.**\n",
    "Dans la limite des $n$ très grands, une matrice percole presque sûrement si $p>0.592746$ et presque jamais sinon.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**g)**\n",
    "<font color=\"red\">!!!</font> Les statistiques effectuées au point précédent sont un exemple typique de calculs pouvant être facilement exécutés en parallèle, car chaque cas est indépendant des autres. Paralléliser l'algorithme précédent de manière à utiliser tous les cœurs de son processeur, par exemple à l'aide du module [`mpi4py`](https://mpi4py.readthedocs.io/).\n",
    "\n",
    "**Indication.**\n",
    "L'utilisation de Jupyter Lab pour faire du calcul parallèle est assez complexe à mettre en œuvre, il vaut mieux utiliser la ligne de commande pour exécuter un script en parallèle, par exemple pour quatre cœurs: `mpirun -n 4 script.py`. À noter que [Open MPI](https://www.open-mpi.org/) ou [MPICH](https://www.mpich.org/) doit être installé sur l'ordinateur.\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}