{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Aiyagari (1994) in Continuous Time\n",
    "\n",
    "#### By [SeHyoun Ahn](http://www.princeton.edu/~sehyouna/) and [Benjamin Moll](http://www.princeton.edu/~moll/)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The material in this notebook is based on [Achdou et al. (2015) \"Heterogeneous Agent Models in Continuous Time\"](http://www.princeton.edu/~moll/HACT.pdf) and follows closely the material in the paper's [online Appendix](http://www.princeton.edu/~moll/HACTproject/HACT_Numerical_Appendix.pdf). Additional codes (mainly in MATLAB) can be found [here](http://www.princeton.edu/~moll/HACTproject.htm). The code and its performance should be compared to [QuantEcon's discrete-time version of the Aiyagari model](http://quant-econ.net/py/aiyagari.html). The structure of the code in the present notebook is purposely kept as close as possible to that of the discrete-time code. The continuous-time code is considerably faster."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We begin by importing some packages"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy import sparse\n",
    "from scipy.sparse.linalg import spsolve"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Aiyagari Economy\n",
    "The economy can be represented by the following system of equations which we aim to solve numerically:\n",
    "$$\n",
    "\\begin{align*}\n",
    "\\rho v_1(a) &=  \\max_c \\ u(c) + v_1'(a)(wz_1 + ra - c) + \\lambda_1(v_2(a) - v_1(a))\\\\\n",
    "\\rho v_2(a) &=  \\max_c \\ u(c) + v_2'(a)(wz_2 + ra - c) + \\lambda_2(v_1(a) - v_2(a))\\\\\n",
    "0 &= - \\frac{d}{da}[s_1(a)g_1(a)] - \\lambda_1 g_1(a) + \\lambda_2 g_2(a)\\\\\n",
    "0 &= - \\frac{d}{da}[s_2(a)g_2(a)] - \\lambda_2 g_2(a) + \\lambda_1 g_1(a)\\\\\n",
    "1 &= \\int_{\\underline{a}}^\\infty g_1(a)da + \\int_{\\underline{a}}^\\infty g_2(a)da\\\\\n",
    "K &= \\int_{\\underline{a}}^\\infty a g_1(a)da + \\int_{\\underline{a}}^\\infty a g_2(a)da\\\\\n",
    "r &= \\alpha K^{\\alpha-1} - \\delta, \\quad w=(1-\\alpha)K^\\alpha\n",
    "\\end{align*}\n",
    "$$\n",
    "\n",
    "where $z_1$ < $z_2$ and $s_j(a)=wz_j + ra -c_j(a)$ and $c_j(a) = (u')^{-1}(v_j(a))$ are optimal savings and consumption. Finally, there is a state constraint $a\\geq \\underline{a}$. The first order condition $u'(c_j(\\underline{a}))=v'_j(\\underline{a})$ still holds at the borrowing constraint. However, in order to respect the constraint we need $s_j(\\underline{a}) = z_j + ra - c_j(\\underline{a}) \\geq 0$. Combining this with the FOC, the state constraint motivates a boundary condition\n",
    "$$\n",
    "\\begin{align} v_j'(\\underline{a}) \\geq u'(z_j + r \\underline{a}), \\quad j=1,2 \\end{align}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We use a finite difference method, an in particular an \"implicit upwind scheme.\" The details are explained in the paper's [online Appendix](http://www.princeton.edu/~moll/HACTproject/HACT_Numerical_Appendix.pdf). We here provide a brief summary. We approximate the functions $(v_1,v_2,g_1,g_2)$ at $I$ discrete points in the space dimension, $a_i,i=1,...,I$. We use equispaced grids, denote by $\\Delta a$ the distance between grid points, and use the short-hand notation $v_{i,j} \\equiv v_j(a_i)$ and so on. The derivative $v_{i,j}'=v_j'(a_i)$ is approximated with either a forward or a backward difference approximation\n",
    "$$\n",
    "\\begin{align*}\n",
    "v_j'(a_i) \\approx \\frac{v_{i+1,j} - v_{i,j}}{\\Delta a} \\equiv v_{i,j,F}' \\\\\n",
    "v_j'(a_i) \\approx \\frac{v_{i-1,j} - v_{i,j}}{\\Delta a} \\equiv v_{i,j,B}'\n",
    "\\end{align*}\n",
    "$$\n",
    "An upwind scheme means that we approximate the derivative $v_j'(a_i)$ with a forward difference approximation whenever the drift of the state variable is positive and the backward difference approximation whenever it is negative. An implicit scheme is a particular way of iterating on the value function.\n",
    "\n",
    "In a nutshell, the discretized HJB equation can be written as \n",
    "$$\\rho v = u + \\mathbf{A} v$$\n",
    "where $\\mathbf{A}$ is $N \\times N$ transition matrix with $N = 2 \\times I$ and where $I$ is the number of wealth grid points. The matrix $\\mathbf{A}$ depends on $v$, i.e. this is a nonlinear problem and we therefore need to iterate (this is where the implicit scheme comes in). The matrix $\\mathbf{A}$ has the interpretation of a Poisson transition matrix (or \"intensity matrix\") on the discretized state space $(a_i,z_j)$.\n",
    "\n",
    "Similarly, one can show that the discretized Kolmogorov Forward equation is\n",
    "$$0 = \\mathbf{A}^T g$$\n",
    "which is an eigenvalue problem. That is, the discretized stationary distribution $g$ is the eigenvector corresponding to a zero eigenvalue of the transpose of the Poisson transition matrix $\\mathbf{A}$. The transpose comes from the fact that the differential operator in the KF equation is the \"adjoint\" of the operator in the HJB equation. And an \"adjoint\" is the infinite-dimensional analogue of matrix transpose.\n",
    "\n",
    "The matrix $\\mathbf{A}$ is found from the discretized HJB equation. Skipping a number of steps it can be written as\n",
    "$$\n",
    "\\begin{align*}\n",
    "&\\frac{v^{n+1}_{i,j} - v^{n}_{i,j}}{\\Delta} + \\rho v_{i,j}^{n+1} = u(c_{i,j}^n) + v_{i-1,j}^{n+1}x_{i,j}   +  v^{n+1}_{i,j} y_{i,j} + v_{i+1,j}^{n+1} z_{i,j} + v_{i,-j}^{n+1}\\lambda_j \\quad \\mbox{where}\\\\\n",
    "&x_{i,j} = -\\frac{(s^n_{i,j,B})^-}{\\Delta a},\\\\\n",
    "&y_{i,j} = - \\frac{(s^n_{i,j,F})^+}{\\Delta a}  +  \\frac{ (s^n_{i,j,B})^-}{\\Delta a} - \\lambda_j,\\\\\n",
    "&z_{i,j} = \\frac{(s^n_{i,j,F})^+}{\\Delta a}\n",
    "\\end{align*}\n",
    "$$\n",
    "where $s^n_{i,j,F}$ and $s^n_{i,j,B}$ are the discretized optimal household savings at grid points $(a_i,z_j)$ using forward and backward approximations, and where for any number $x$, the notation $x^+$ means \"the positive part of $x$\", i.e. $x^+ = \\max\\{x,0\\}$ and analogously $x^{-} = \\min\\{x,0\\}$. This part is what makes it an upwind scheme.\n",
    "\n",
    "\n",
    "This is a system of $2 \\times I$ linear equations which can be written in matrix notation as:\n",
    "\\begin{equation}\\frac{1}{\\Delta}(v^{n+1} - v^n) + \\rho v^{n+1} = u^n + \\mathbf{A}^n v^{n+1} \\end{equation}\n",
    "where\n",
    "$$\\mathbf{A}^n = \\left[\\begin{matrix}y_{1,1} & z_{1,1} & 0 & \\cdots & 0 & \\lambda_1 & 0 & 0 & \\cdots & 0 \\\\ x_{2,1} & y_{2,1} & z_{2,1} & 0 & \\cdots & 0 & \\lambda_1 & 0 & 0 & \\cdots \\\\ 0 & x_{3,1} & y_{3,1} & z_{3,1} & 0 & \\cdots & 0 & \\lambda_1 & 0 & 0 \\\\ \\vdots & \\ddots & \\ddots & \\ddots & \\ddots & \\ddots & \\ddots & \\ddots & \\ddots & \\vdots \\\\ 0 & \\ddots & \\ddots & x_{I,1} & y_{I,1} & 0 & 0 & 0 & 0 & \\lambda_1\\\\  \\lambda_2 & 0 & 0 & 0 & 0 & y_{1,2} & z_{1,2} & 0 & 0 & 0 \\\\ 0 & \\lambda_2 & 0 & 0 & 0 & x_{2,2} & y_{2,2} & z_{2,2} & 0 & 0  \\\\ 0 & 0 & \\lambda_2 & 0 & 0 & 0 & x_{3,2}  & y_{3,2} & z_{3,2} & 0 \\\\ 0 & 0 & \\ddots & \\ddots & \\ddots & \\ddots & \\ddots  & \\ddots & \\ddots & \\ddots \\\\ 0 & \\cdots & \\cdots & 0 & \\lambda_2 & 0 & \\cdots  & 0 & x_{I,2} & y_{I,2}  \\end{matrix}\\right], \\quad u^n = \\left[\\begin{matrix} u(c_{1,1}^n)\\\\ \\vdots \\\\ \\vdots \\\\ u(c_{I,1}^n)\\\\  u(c_{1,2}^n) \\\\ \\vdots \\\\ \\vdots \\\\ u(c_{I,2}^n)\\end{matrix}\\right]$$\n",
    "This system can in turn be written as\n",
    "$$\n",
    "\\begin{align*}\\mathbf{B}^n v^{n+1} = b^n, \\quad \\quad \\mathbf{B}^n = \\left(\\frac{1}{\\Delta} + \\rho\\right)\\mathbf{I} - \\mathbf{A}^n, \\quad  b^n = u^n + \\frac{1}{\\Delta}v^n \\end{align*}\n",
    "$$\n",
    "This system of linear equations can be solved very efficiently using sparse matrix routines. In Python this is implemented with the function \"spsolve().\"\n",
    "\n",
    "#### Summary of Algorithm\n",
    "First consider the algorithm for solving the HJB equations. Guess $v^0_{i,j},i=1,...,I,j=1,2$ and for $n=0,1,2,...$ follow\n",
    "1. Compute $(v^n_{i,j})'$ using the current guess of the value function and the upwind scheme (forward difference when drift is positive, backward difference when drift is negative)\n",
    "2. Compute $c^n$ from $c_{i,j}^n = (u')^{-1}[(v_{i,j}^n)']$\n",
    "3. Find $v^{n+1}$ by solving the system of linear equations involving the matrix $\\mathbf{A}$ described above (implicit scheme)\n",
    "4. If $v^{n+1}$ is close enough to $v^n$: stop. Otherwise, go to step 1.\n",
    "\n",
    "After solving the HJB equations, solving the Kolmogorov Forward equation: simply solve the eigenvalue problem $0=\\mathbf{A}^T g$ described above. That is, once the HJB equation is solved, we basically get the Kolmogorov Forward equation \"for free.\"\n",
    "\n",
    "#### Overview of Code\n",
    "Given this overview of the algorithm, we now briefly describe the code. As in the discrete-time version, we define a \"household\" class. Household objects contain all the data relevant to solving a household's decision problem. The household class  contain:\n",
    "* economic parameters (e.g., w, r)\n",
    "* utility paramters (e.g., discount factor $\\rho$)\n",
    "* asset and skill level represented on a grid\n",
    "\n",
    "The household's decision problem is solved by invoking the function solve_bellman() which solves the HJB equation given the relevant parameters saving the value-function as v.\n",
    "\n",
    "We also include the stationary wealth distribution of households in the household object. This is a natural thing to do since this stationary distribution is computed using the household's decision rule. Hence, after computing the decision problem of households, the stationary distribution can be found by invoking compute_stationary_distribution()\n",
    "\n",
    "The definition of the household class is given below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "class Household(object):\n",
    "    def __init__(self,\n",
    "                 dep=0.05,   ###ADDED DEPRECIATION\n",
    "                 r=0.03,         # interest rate\n",
    "                 w=1,            # wages\n",
    "                 rho=0.04,       # discount factor\n",
    "                 a_min=1e-10,    # minimum asset amount\n",
    "                 pi=[[-0.33, 0.33], [0.33, -0.33]],  # poisson Jumps\n",
    "                 z_vals=[1.0, 2.0],                # exogenous income states\n",
    "                 a_max=40,\n",
    "                 a_size=1000,    # number of asset grid points\n",
    "                 delta=1000.0):\n",
    "        # Initialize values, and set up grids over a and z\n",
    "        self.r, self.w, self.rho, self.dep = r, w, rho, dep\n",
    "        self.a_min, self.a_max, self.a_size = a_min, a_max, a_size\n",
    "        self.da = (self.a_max-self.a_min)/(self.a_size-1)\n",
    "        self.k = 10\n",
    "        self.pi = np.asarray(pi)\n",
    "        self.z_vals = np.asarray(z_vals)\n",
    "        self.z_size = len(z_vals)\n",
    "        self.a_vals = np.linspace(self.a_min, self.a_max, self.a_size)\n",
    "        self.n = self.a_size * self.z_size\n",
    "        self.delta = delta\n",
    "        ###### ADDED TO MATCH LABOR SUPPLY IN .m\n",
    "        self.z_ave = (self.z_vals[0]*self.pi[0, 1] +\n",
    "                      self.z_vals[1]*self.pi[1, 0]) / \\\n",
    "            (self.pi[0, 1] + self.pi[1, 0])\n",
    "        \n",
    "        # Initial Guess of Value Function\n",
    "        self.v = np.log(np.tile(self.a_vals,(self.z_size,1))*self.r\n",
    "                        +self.w*np.tile(self.z_vals,(self.a_size,1)).transpose())/self.rho\n",
    "        \n",
    "        # Build skill_transition, the matrix summarizing transitions due to the Poisson income shocks\n",
    "        # This is analogous to the Q matrix in the discrete time version of the QuantEcon Aiyagari model\n",
    "        self.z_transition = sparse.kron(self.pi,sparse.eye(self.a_size), format=\"csr\")\n",
    "        \n",
    "        # Preallocation\n",
    "        self.v_old = np.zeros((self.z_size,self.a_size))\n",
    "        self.g = np.zeros((self.z_size,self.a_size))\n",
    "        self.dv = np.zeros((self.z_size,self.a_size-1))\n",
    "        self.cf = np.zeros((self.z_size,self.a_size-1))\n",
    "        self.c0 = np.zeros((self.z_size,self.a_size))\n",
    "        self.ssf = np.zeros((self.z_size,self.a_size))\n",
    "        self.ssb = np.zeros((self.z_size,self.a_size))\n",
    "        self.is_forward = np.zeros((self.z_size,self.a_size),'bool')\n",
    "        self.is_backward = np.zeros((self.z_size,self.a_size),'bool')\n",
    "        self.diag_helper = np.zeros((self.z_size,self.a_size))        \n",
    "        self.A = self.z_transition.copy()\n",
    "        self.B = self.z_transition.copy()\n",
    "        self.AT = self.z_transition.copy()\n",
    "\n",
    "    def set_prices(self, r, w):\n",
    "        \"\"\"\n",
    "        Resets prices\n",
    "        Calling the method will resolves the Bellman Equation.\n",
    "        \n",
    "        Parameters:\n",
    "        -----------------\n",
    "        r : Interest rate\n",
    "        w : wage\n",
    "        \"\"\"\n",
    "        self.r, self.w = r, w\n",
    "        self.solve_bellman()\n",
    "\n",
    "    def reinitialize_v(self):\n",
    "        \"\"\"\n",
    "        Reinitializes the value function if the value function\n",
    "        became NaN\n",
    "        \"\"\"\n",
    "        self.v = np.log(np.tile(self.a_vals,(self.z_size,1))*self.r\n",
    "                        +self.w*np.tile(self.z_vals,(self.a_size,1)).transpose())/self.rho\n",
    "\n",
    "    def solve_bellman(self,maxiter=100,crit=1e-6):\n",
    "        \"\"\"\n",
    "        This function solves the decision problem with the given parameters\n",
    "        \n",
    "        Parameters:\n",
    "        -----------------\n",
    "        maxiter : maximum number of iteration before haulting value function iteration\n",
    "        \n",
    "        crit : convergence metric, stops if value function does not change more than crit\n",
    "        \"\"\"\n",
    "        dist=100.0\n",
    "        for i in range(maxiter):\n",
    "            # compute saving and consumption implied by current guess for value function, using upwind method\n",
    "            self.dv = (self.v[:,1:]-self.v[:,:-1])/self.da\n",
    "            self.cf = 1.0/self.dv\n",
    "            self.c0 = np.tile(self.a_vals,(self.z_size,1))*self.r \\\n",
    "                            +self.w*np.tile(self.z_vals,(self.a_size,1)).transpose()\n",
    "\n",
    "            # computes savings with forward forward difference and backward difference\n",
    "            self.ssf[:,:-1] = self.c0[:,:-1]-self.cf\n",
    "            self.ssb[:,1:] = self.c0[:,1:]-self.cf\n",
    "            # Note that the boundary conditions are handled implicitly as ssf will be zero at a_max and ssb at a_min \n",
    "            self.is_forward = self.ssf>0\n",
    "            self.is_backward = self.ssb<0\n",
    "            # Update consumption based on forward or backward difference based on direction of drift\n",
    "            self.c0[:,:-1] += (self.cf-self.c0[:,:-1])*self.is_forward[:,:-1]\n",
    "            self.c0[:,1:] += (self.cf-self.c0[:,1:])*self.is_backward[:,1:]\n",
    "            ######\n",
    "            # UNCOMMENT FOR DEBUGGING\n",
    "            #plt.plot(self.a_vals, self.c0.transpose())\n",
    "            #plt.show()\n",
    "            self.c0 = np.log(self.c0)\n",
    "            \n",
    "            # Build the matrix A that summarizes the evolution of the process for (a,z)\n",
    "            # This is a Poisson transition matrix (aka intensity matrix) with rows adding up to zero\n",
    "            self.A = self.z_transition.copy()\n",
    "            self.diag_helper = (-self.ssf*self.is_forward/self.da \\\n",
    "                               + self.ssb*self.is_backward/self.da).reshape(self.n)\n",
    "            self.A += sparse.spdiags(self.diag_helper,0,self.n,self.n)\n",
    "            self.diag_helper = (-self.ssb*self.is_backward/self.da).reshape(self.n)\n",
    "            self.A += sparse.spdiags(self.diag_helper[1:],-1,self.n,self.n)\n",
    "            self.diag_helper = (self.ssf*self.is_forward/self.da).reshape(self.n)\n",
    "            self.A += sparse.spdiags(np.hstack((0,self.diag_helper)),1,self.n,self.n)\n",
    "            # Solve the system of linear equations corresponding to implicit finite difference scheme\n",
    "            self.B = sparse.eye(self.n)*(1/self.delta + self.rho) - self.A\n",
    "            self.b = self.c0.reshape(self.n,1) + self.v.reshape(self.n,1)/self.delta\n",
    "            self.v_old = self.v.copy()\n",
    "            self.v = spsolve(self.B,self.b).reshape(self.z_size,self.a_size)\n",
    "\n",
    "            # Compute convergence metric and stop if it satisfies the convergence criterion\n",
    "            dist = np.amax(np.absolute(self.v_old-self.v).reshape(self.n))\n",
    "            if dist < crit:\n",
    "                break\n",
    "\n",
    "    def compute_stationary_distribution(self):\n",
    "        \"\"\"\n",
    "        Solves for the stationary distribution given household decision rules\n",
    "        \n",
    "        Output:\n",
    "        Capital level from the stationary distribution\n",
    "        \"\"\"        \n",
    "        self.AT = self.A.transpose().tocsr()\n",
    "\n",
    "        # The discretized Kolmogorov Forward equation AT*g=0 is an eigenvalue problem\n",
    "        # AT is singular because one of the equation is the distribution adding\n",
    "        # up to 1. Here we solve the eigenvalue problem by setting g(1,1)=0.1\n",
    "        # and the equation is solved relative to that value.\n",
    "        # Alternatively, one could use a routine for solving eigenvalue problems.\n",
    "        b = np.zeros((self.n,1))\n",
    "        b[0] = 0.1\n",
    "        self.AT.data[1:self.AT.indptr[1]] = 0\n",
    "        self.AT.data[0] = 1.0\n",
    "        self.AT.indices[0] = 0\n",
    "        self.AT.eliminate_zeros()\n",
    "        self.g = spsolve(self.AT,b).reshape(self.z_size,self.a_size)\n",
    "        \n",
    "        # Since g was solved taking one of g(1,1) as given, g needs to be\n",
    "        # renormalized to add up to 1\n",
    "        self.g = self.g/np.sum(self.g)\n",
    "        return np.sum(self.g*(np.tile(self.a_vals,(self.z_size,1))))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "For example, if interest rate is 0.05 and wage is 1, we can initialize a household, solve it decision problem, and find the stationary distribution by running"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.69274641340853271"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lam = 0.11\n",
    "PI = [[-lam, lam], [lam, -lam]]\n",
    "am=Household(rho=0.05, r=0.02,w=1, pi=PI)\n",
    "am.solve_bellman()\n",
    "am.compute_stationary_distribution()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once, a household object is created, the object can be reused to solve a different problem by changing parameters. For example, if the interest rate were to change to 0.03 and wage to 0.9, the new problem can be solved by setting parameters directly by"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "am.r=0.02\n",
    "am.w=0.9"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "and running"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.62323720534758664"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "am.solve_bellman()\n",
    "am.compute_stationary_distribution()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Alternately, a helper function can be created to automate this process. For example, set_prices(r,w) function resets parameters, and solves the decision problem again. For example, you can run the following to the same effect."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.129833308836365"
      ]
     },
     "execution_count": 44,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "am.set_prices(r=0.03,w=0.9)\n",
    "am.compute_stationary_distribution()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Given the household class, solving for the steady state becomes really simple. To find the equilibrium, we solve for the capital level that is consistent with household's decisions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Before solving for the actual steady state, we can visualize the capital demand and capital supply."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "A = 0.1\n",
    "##### NOT NECESSARY\n",
    "#N = 0.05\n",
    "alpha = 0.33\n",
    "\n",
    "def r_to_w(am, r):\n",
    "    return A*(1 - alpha) * \\\n",
    "        (alpha*A / (am.dep + r))**(alpha / (1 - alpha))\n",
    "\n",
    "def rd(am, K):\n",
    "    return A*alpha*(am.z_ave / K)**(1 - alpha) - am.dep\n",
    "\n",
    "def prices_to_capital_stock(am, r):\n",
    "    \"\"\"\n",
    "    Map prices to the induced level of capital stock.\n",
    "    \n",
    "    Parameters:\n",
    "    ----------\n",
    "    am : Household\n",
    "        An instance of the Household class\n",
    "    r : float\n",
    "        The interest rate\n",
    "    \"\"\"\n",
    "    w = r_to_w(am, r)\n",
    "\n",
    "    # Set new prices and solve the Bellman equation\n",
    "    am.set_prices(r, w)\n",
    "\n",
    "    # Compute the stationary distribution and capital\n",
    "    return am.compute_stationary_distribution()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We make a grid of interest rate points, and plot the resulting capital."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAqUAAAHuCAYAAACmr3UAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xl4VOXZx/HvM0lICBASEtawhEVAZRVFFCtRqqIoorxV\n1GrBpXRBq91AUNBi6/pal2pVXCoq8lZrVcSdTVEQZFM2WcMaZA3ZJOvz/pHFAAlMlpnzzOT3ua65\nmmfmnJk7/UF7c859zhhrLSIiIiIiXvJ5XYCIiIiIiJpSEREREfGcmlIRERER8ZyaUhERERHxnJpS\nEREREfGcmlIRERER8VzAm1JjzBBjzDpjzHpjzLjjbHeGMabAGHNlhefSjDErjTHLjTGLA12riIiI\niHgjMpBvbozxAf8ABgO7gCXGmHestesq2e4B4KOj3qIYSLXWHgxknSIiIiLirUAfKe0PbLDWbrXW\nFgAzgMsr2e5W4E1gz1HPGzRiICIiIhL2At3wJQPbK6x3lD5XzhjTBhhurf0nJU1oRRb4xBizxBhz\nS0ArFRERERHPBPT0vZ8eAyrOmlZsTAdaa9ONMc0paU7XWmsXHP0Gxhh9V6qIiIiIA6y1Rx9k9Eug\nj5TuBNpXWLctfa6i04EZxpgtwP8ATxljhgFYa9NL/3Mv8F9KxgEqZa3Vw7HH5MmTPa9BD+USKg/l\n4u5D2bj5UC5uPmoj0E3pEqCLMaaDMaYBMBJ4t+IG1tpOpY+OlMyV/sZa+64xJtYY0xjAGNMIuBBY\nFeB6pQ6lpaV5XYJUQrm4Sbm4S9m4SbmEn4CevrfWFhljxgIfU9IAv2CtXWuMGVPysn3u6F0q/NwS\n+G/pqflI4DVr7ceBrFdEREREvBHwmVJr7YdAt6Oee7aKbW+s8PMWoE9gq5NAGjVqlNclSCWUi5uU\ni7uUjZuUS/gxtT3/7wJjjA2H30NEREQklBljsI5e6CT12Lx587wuQSqhXNykXNwVCtmkpKRgjNFD\nj6A9UlJS6vzPsQu3hBIREZFa2Lp1a62vfBapDmNqdDD0+O8ZDn+IjU7fi4hIPWaMUVMqQVXVn7nS\n53X6XkRERERCk5pSCZhQmMOqj5SLm5SLu5SNSHCoKRURERERz6kplYBJTU31ugSphHJxk3Jxl7IJ\nPaNHj2bSpEl1/r5ffvklXbt2JS4ujnfffffEO9SBBQsWcPLJJ9d4/0D9dxEIakpFRERE/DBp0iRu\nu+02MjMzGTZsWFA+85xzzmHt2rXl644dOzJnzpygfHawqSmVgNEclpuUi5uUi7uUjZTZunUrp5xy\nitdlhC01pdUwa/0sJs6eyFc7vvK6FBERkZDw4IMP0rZtW+Li4jj55JOZO3cucOxp5fnz59OuXbvy\ndceOHXnggQc49dRTSUxM5KabbiI/P/+Ibe+//36aN29Op06dmD59eqWf37NnT2bNmlW+LiwspHnz\n5qxcubLS7adOncpJJ51EUlISw4cPZ/fu3QB06dKFLVu2cOmllxIXF0dBQcEx++7YsYMRI0bQokUL\nmjdvzm233QbA5s2bGTx4MElJSbRo0YKf//znZGZmVut3BbjhhhvYtm0bl112GXFxcTzyyCMAXHXV\nVbRu3ZqEhARSU1NZs2bN8SJxlm6eXw2HCw+zL3cfB3444HUpIUFzWG5SLm5SLu4K9WzGjKm793r2\n2eptv379ep566imWLl1Ky5Yt2bZtG0VFRVVuf/QN2adPn84nn3xCbGwsl156Kffddx9/+ctfANi9\nezcHDhxg165dLFy4kEsuuYQzzjiDk0466Yj3uOGGG3jllVcYOnQoALNmzaJNmzb07t37mM+fM2cO\nEyZM4NNPP+WUU07hD3/4A1dffTXz589n48aNdOzYkRdffJHzzjvvmH2Li4u59NJL+elPf8prr72G\nz+fj66+/BsBay4QJExg0aBCHDh1ixIgR3HPPPTz66KN+/a5l/71MmzaNzz///JgaLrnkEv71r38R\nFRXFuHHjuO6661i+fHnVwThKR0qrIS46DoCs/CyPKxEREXFfREQE+fn5rFq1isLCQtq3b0/Hjh39\n3v/WW2+lTZs2xMfHM3HiRF5//fXy14wxTJkyhaioKM4991yGDh3Kv//972Pe4+c//zkffPAB2dnZ\nALz66qtcf/31lX7e9OnTuemmm+jduzdRUVHcf//9LFy4kG3btpVvU9WXFCxevJj09HQeeughYmJi\naNCgAWeffTYAnTt3ZvDgwURGRpKYmMgdd9zB/Pnz/f5dj3Z0DaNGjSI2NpaoqCgmTZrEypUrycoK\nvV5FR0qroawpPXT4kMeVhIZ58+aF/BGGcKRc3KRc3BXq2VT36GZd6ty5M4899hj33HMPa9as4aKL\nLuLRRx+lVatWfu3ftm3b8p87dOjArl27ytcJCQnExMRU+XqZ1q1bM3DgQP7zn/8wfPhwPvjgA554\n4olKP2/Xrl3069evfN2oUSMSExPZuXMn7du3P26t27dvp0OHDvh8xx7v27NnD7/73e/4/PPPyc7O\npqioiGbNmvn9ux5PcXExEyZM4M0332Tfvn3l302/b98+mjRp4td7uEJHSquhrCnNzMs8wZYiIiIC\nMHLkSD7//HO2bt0KwLhx44CShi83N7d8u/T09GP23b59e/nPW7dupU2bNuXrgwcP8sMPP5Svt23b\ndsTrFZWdwn/jjTc4++yzad26daXbtWnTprxOgJycHPbv339Ew1iVdu3asW3bNoqLi495bcKECfh8\nPlavXk1GRgavvvrqMUc7j/e7VlTZiMPMmTOZM2cOGRkZpKWlYa0Nya+dVVNaDTp9Xz2hfGQhnCkX\nNykXdymbmlu/fj1z584lPz+fBg0a0LBhw/IjiX369OH999/n4MGD7N69m8cff/yY/Z966il27tzJ\ngQMH+Nvf/sbIkSPLX7PWMnnyZAoKCvj888+ZNWsWV111VaV1DB8+nGXLlvHEE09www03VFnvNddc\nw0svvcQ333xDXl4eEyZMYMCAAUdcgFWV/v3707p1a8aPH09ubi55eXl8+eWXAGRlZdG4cWOaNGnC\nzp07efjhh6v1u1bUqlUrNm/eXL7OysoiOjqahIQEcnJyuPPOO49pXEOFmtJq0JFSERER/+Xl5TF+\n/HiaN29OmzZt2Lt3L/fffz8A119/Pb169SIlJYUhQ4ZU2oRde+21XHjhhXTp0oWTTjqJiRMnlr9W\ndrV5mzZtuP7663n22WfLL3I6uimLiYlhxIgRbNmyhSuvvLLKegcPHsyUKVO48sorSU5OZsuWLcyY\nMaP89eM1ez6fj5kzZ7Jhwwbat29Pu3btymdcJ0+ezNKlS4mPj+eyyy5jxIgR1fpdKxo/fjxTpkyh\nWbNmPProo/ziF7+gffv2JCcn06NHj/I51lBkQvHw7tGMMTYYv0exLeY3s36DtZanhz5NhC8i4J8Z\nykJ9DitcKRc3KRd3hUI2xpiQPF17PB07duSFF17g/PPPP+a1+fPnc/311x9xAdKJTJkyhQ0bNjBt\n2rS6LLNOHO93dVVVf+ZKn6/RoVodKa0Gn/HRpEHJ0LBO4YuIiISGAwcO8MILLzCmLu+PJXVOTWk1\n6Qp8/7l+ZKG+Ui5uUi7uUjbeqKu5yOeff5727dszdOhQBg4cWCfvWddCdQa0run0fTU9vuhx1uxd\nw61n3kqPFj2C8pkiIiLHE46n78VtOn3vAF3s5D99X7SblIublIu7lI1IcKgprSY1pSIiIiJ1T6fv\nq+njTR/znzX/YXCnwVx1auX3QxMREQkmnb6XYNPpewfoSKmIiIhI3VNTWk1qSv2nOSw3KRc3KRd3\nKRuR4FBTWk1qSkVERGpu9OjRTJo0yesyTui8887jxRdfrPL10aNH06xZMwYMGBDEqqrn/vvv55e/\n/GWN9+/YsSNz5sypw4qOLzJonxQm1JT6T/f2c5NycZNycZeykaMtWLCA2bNns2vXLmJiYrwup0p3\n3nln+c9bt26lY8eOFBYW4vO5eUzSzaoc1rhBY3zGR05+DoXFhV6XIyIiIkGWlpZGSkpKjRrSoqKi\nAFR0YtZa5y+IU1NaTT7j07c6+UlzWG5SLm5SLu5SNrWzfPly+vXrR9OmTRk5ciSHDx8+4vX33nuP\nvn37kpCQwDnnnMO3335b/lrHjh155JFH6N27N02aNOGWW25hz549XHLJJcTFxXHhhRdy6NCP/198\n1VVX0bp1axISEkhNTWXNmjXlr40ePZqxY8dy6aWXEhcXx1lnncWWLVvKX//kk084+eSTSUhI4NZb\nb62yeXvxxRe55ZZbWLhwIXFxcdx7770ATJ06lZNOOomkpCSGDx9Oenp6+T4+n4+nn36arl270rVr\n10rfd8GCBQwcOJCEhAQ6dOjAtGnTAHj//fc57bTTaNq0KR06dCj/PCg5+unz+Zg6dSrJyckkJyfz\nv//7v+Wv33vvvdxwww0ADBo0CID4+Hji4uL46quv2Lx5M4MHDyYpKYkWLVrw85//nMxM784E6/R9\nDcTHxJNxOIOMwxkkxiZ6XY6IiEiVxsysu+97f/ayZ6u1fUFBAVdccQW///3v+e1vf8vbb7/NNddc\nw/jx44GShvWmm25i1qxZ9OvXj1dffZVhw4axfv16oqKiAHjrrbeYPXs2BQUF9OnTh+XLl/Piiy/S\nvXt3Lr74Yp544gnuvvtuAC655BL+9a9/ERUVxbhx47juuutYvnx5eT3/93//x4cffkjfvn254YYb\nmDhxItOnT2f//v2MGDGCl19+mWHDhvHkk0/yzDPPlDd0Fd14441ERETwwgsv8NlnnwEwZ84cJkyY\nwKeffsopp5zCH/7wB0aOHMn8+fPL93vnnXdYsmRJpUdXt23bxiWXXMLzzz/PiBEjyMzMZPv27QA0\nbtyYV155hVNPPZVVq1ZxwQUX0LdvX4YNG1a+/7x589i0aRMbN27k/PPPp2/fvpx//vlHfMZnn31G\np06dyMzMLP9a002bNjFhwgQGDRrEoUOHGDFiBPfccw+PPvpotXKuKzpSWgMJDRMAOHj4oMeVuE1z\nWG5SLm5SLu5SNjW3aNEiCgsLue2224iIiGDEiBGcccYZ5a9PnTqVX/3qV5x++ukYY7j++uuJjo5m\n0aJF5dvceuutJCUl0bp1a37yk59w5pln0qtXLxo0aMAVV1xxRNM5atQoYmNjiYqKYtKkSaxcuZKs\nrKzy16+44gr69euHz+fjuuuuY8WKFUDJ0cgePXpwxRVXEBERwe23306rVq38/j2nT5/OTTfdRO/e\nvYmKiuL+++9n4cKFbNu2rXybCRMm0LRpU6Kjoyvd/4ILLuCqq64iIiKChIQEevXqBcC5557Lqaee\nCkCPHj2OaXYB7rnnHmJiYujRowejR4/m9ddfr7LWikeAO3fuzODBg4mMjCQxMZE77rjjmPcOJh0p\nrYH4mHgAMg5neFyJiIjI8VX36GZd2rVrF8nJyUc816FDh/Kft27dyrRp03jyySeBkoapoKCAXbt2\nlW/TsmXL8p8bNmx4zDo7OxuA4uJiJkyYwJtvvsm+ffswxmCMYd++fTRp0gTgiEYzNja2fN9du3bR\nrl27I+o8en2i37Nfv37l60aNGpGYmMjOnTtp3749AG3btq1y/+3bt9O5c+dKX1u8eDHjx49n1apV\n5Ofnk5+fz89+9rPy140xR7x3hw4dWLVqlV9179mzh9/97nd8/vnnZGdnU1RURLNmzfzaNxB0pLQG\nEmJKj5T+oCOlx6M5LDcpFzcpF3cpm5pr3bo1O3fuPOK5ikcP27Vrx8SJEzlw4AAHDhzg4MGDZGdn\nc/XVV1f7s1577TVmzpzJnDlzyMjIIC0tDWutXxf2tG7d+oi6gPLT5/5o06YNW7duLV/n5OSwf//+\nI5rFslPmlWnXrh0bN26s9LVrr72W4cOHs3PnTjIyMhgzZswRv5O19ohat23bRps2bY55n8o+f8KE\nCfh8PlavXk1GRgavvvqqpxdCqSmtAR0pFRERObGzzjqLyMhInnzySQoLC3nrrbdYvHhx+eu33HIL\nzzzzTPlzOTk5vP/+++Tk5FT7s7Kzs4mOjiYhIYGcnBzuvPPO4zaCFQ0dOpQ1a9bw9ttvU1RUxOOP\nP87333/v92dfc801vPTSS3zzzTfk5eUxYcIEBgwY4PfR1uuuu47Zs2fz5ptvUlRUxIEDB1i5cmX5\n75WQkEBUVBSLFy9m+vTpx+w/ZcoUfvjhB1avXs1LL73EyJEjj9mmefPm+Hw+Nm3aVP5cVlYWjRs3\npkmTJuzcuZOHH37Y7985ENSU1oBmSv2jOSw3KRc3KRd3KZuai4qK4q233uKll14iMTGRN954gxEj\nRpS/3q9fP6ZOncrYsWNp1qwZXbt25eWXXy5//eim8nhN5g033ED79u1JTk6mR48enH322X7XWVbb\nuHHjSEpKYtOmTQwcONDv/QcPHsyUKVO48sorSU5OZsuWLcyYMcOvuqHkSOn777/PI488QrNmzejb\nty/ffPMNAE899RR33303TZs25b777qv0KPKgQYPo0qULF1xwAX/+858ZPHjwMds0bNiQiRMnMnDg\nQJo1a8bixYuZPHkyS5cuJT4+nssuu+yIbPypu64Zl+9X5S9jjA3m77EnZw93z7mbxNhE/jb4b0H7\nXBERkcq4fv9JCYytW7fSqVMnCgoKgn5D/Kr+zJU+X6NuVkdKa6Ds9P2hw4f0PwLHoTksNykXNykX\ndykbcVk49SFqSmugQUQDYqNiKSwuJDs/2+tyREREpJ4K9in2QNLp+xq6d9697MraxV3n3kW7pv7f\nNkJERKSu6fS9BJtO3ztEFzuJiIiI1J2AN6XGmCHGmHXGmPXGmHHH2e4MY0yBMebK6u7rBd0W6sQ0\nh+Um5eIm5eIuZSMSHAFtSo0xPuAfwEXAqcA1xpjuVWz3APBRdff1ippSERERkboT6K8Z7Q9ssNZu\nBTDGzAAuB9Ydtd2twJvAGTXY1xNqSk9M9/Zzk3Jxk3JxVyhk06FDh7C64EXcV/HrYutKoJvSZKDi\n93TtoKTZLGeMaQMMt9aeZ4zpX519vaSvGhUREVekpaV5XYJIrQW6KfXHY0Ct50VHjRpFSkoKAPHx\n8fTp06f8X7dl80B1ud6TvQcoOVIaiPcPh3XZc67Uo3XJ+rHHHgv43w+t9fclnNYrVqzg9ttvd6Ye\nrY/8u+JKPfV1vWLFCjIySs4a1/YfRwG9JZQxZgBwj7V2SOl6PGCttQ9W2GZz2Y9AEpAD/BLYc6J9\nK7xH0G8JlZWXxR8//iOxUbH8fcjfg/rZoWLevHnlf3DFHcrFTcrFXcrGTcrFTbW5JVSgm9II4Dtg\nMJAOLAausdaurWL7l4CZ1tq3qrOvF02ptZax74+lsLiQJy5+gujI6KB+voiIiIhrnL1PqbW2CBgL\nfAysBmZYa9caY8YYY35Z2S4n2jeQ9VaHMYamMU0BXewkIiIiUlsBbUoBrLUfWmu7WWtPstY+UPrc\ns9ba5yrZ9kZr7VvH29clZRc7qSmtXMV5H3GHcnGTcnGXsnGTcgk/AW9Kw5luCyUiIiJSNwI6Uxos\nXsyUAryx+g0+3fwpV558JRd1uSjony8iIiLiEmdnSsNdQsPSe5Ue1r1KRURERGpDTWkt6PT98Wne\nx03KxU3KxV3Kxk3KJfyoKa0FfauTiIiISN3QTGkt7M/dz4TZE4iPiefBC465p7+IiIhIvaKZUo+U\n3ac0My+TYlvscTUiIiIioUtNaS1E+iJpEt2EYltMZl6m1+U4R/M+blIublIu7lI2blIu4UdNaS3p\nYicRERGR2tNMaS09tfgpvvn+G351+q/o27qvJzWIiIiIuEAzpR7SkVIRERGR2lNTWktlTaluoH8s\nzfu4Sbm4Sbm4S9m4SbmEHzWltVT2rU46UioiIiJSc5opraU1e9fw+KLH6ZbUjd+f9XtPahARERFx\ngWZKPaSZUhEREZHaU1NaSxWb0nA46lyXNO/jJuXiJuXiLmXjJuUSftSU1lLDyIZER0aTV5jHD4U/\neF2OiIiISEjSTGkdmDR3Et9nf8/k1Mm0adLGszpEREREvKSZUo9prlRERESkdtSU1oGEmJLbQu3P\n3e9xJW7RvI+blIublIu7lI2blEv4UVNaB5o3ag7A/h/UlIqIiIjUhGZK68CiHYt4aflLnJF8Bjef\ndrNndYiIiIh4STOlHktsmAjA3py9HlciIiIiEprUlNYBnb6vnOZ93KRc3KRc3KVs3KRcwo+a0jrQ\nNLopkb5IsvKyyCvM87ocERERkZCjmdI6Unav0kmDJpEcl+xpLSIiIiJe0EypA5JikwDYl7vP40pE\nREREQo+a0jqipvRYmvdxk3Jxk3Jxl7Jxk3IJP2pK64iaUhEREZGa00xpHVmWvoxnv36WXi178dv+\nv/W0FhEREREvaKbUATpSKiIiIlJzakrrSMWm1Oujtq7QvI+blIublIu7lI2blEv4UVNaR2KjYomN\niiW/KJ/s/GyvyxEREREJKZoprUP3fXYf2w9tZ/w54+mY0NHrckRERESCSjOljmgeW/J1o5orFRER\nEakeNaV1SBc7HUnzPm5SLm5SLu5SNm5SLuFHTWkdSoxNBNSUioiIiFSXZkrr0Oo9q3niqyfontSd\nO866w+tyRERERIJKM6WO0Ol7ERERkZpRU1qHmjVsBsCBHw5QbIs9rsZ7mvdxk3Jxk3Jxl7Jxk3IJ\nPwFvSo0xQ4wx64wx640x4yp5fZgxZqUxZrkxZrExZmCF19IqvhboWmsrKiKK+Jh4im0xB3444HU5\nIiIiIiEjoDOlxhgfsB4YDOwClgAjrbXrKmwTa63NLf25J/Bva+3JpevNQD9r7cETfI4TM6UAD3/x\nMBsPbOSOs+6ge1J3r8sRERERCRqXZ0r7AxustVuttQXADODyihuUNaSlGgMVz3ubINRYp8rmSvfn\n7ve4EhEREZHQEeiGLxnYXmG9o/S5Ixhjhhtj1gIzgRsrvGSBT4wxS4wxtwS00jpS1pTuzd3rcSXe\n07yPm5SLm5SLu5SNm5RL+HHiKKS19u3SU/bDgfsqvDTQWnsacAnwW2PMOZ4UWA26Al9ERESk+iID\n/P47gfYV1m1Ln6uUtXaBMaaTMaaZtfaAtTa99Pm9xpj/UjIOsKCyfUeNGkVKSgoA8fHx9OnTh9TU\nVODHf00FY50Um8Sub3exOG0xN592c9A/X2utT7Que86VerTWOhTWZVypR+tUUlNTnaqnvq5XrFhB\nRkYGAGlpadRGoC90igC+o+RCp3RgMXCNtXZthW06W2s3lf58GvCOtbadMSYW8Flrs40xjYCPgXut\ntR9X8jnOXOh08IeDjP90PE2im/DIhY94XY6IiIhI0Dh7oZO1tggYS0lDuRqYYa1da4wZY4z5Zelm\nI4wxq4wxy4AngatKn28JLDDGLAcWATMra0hdEx8TT6Qvkqy8LPIK87wux1NHH2EQNygXNykXdykb\nNymX8BPo0/dYaz8Euh313LMVfn4IeKiS/bYAfQJdX10zxtC8UXPSs9LZk7OHdk3beV2SiIiIiPMC\nevo+WFw6fQ/w9JKnWbl7Jbf0u4XT25zudTkiIiIiQeHs6fv6qlXjVgDszt7tcSUiIiIioUFNaQC0\nbNQSgO+zv/e4Em9p3sdNysVNysVdysZNyiX8qCkNgLIjpd/n1O+mVERERMRfmikNgOz8bP7w0R+I\njozm8SGPY0yNRitEREREQopmSh3TuEFjGjVoRF5hHpl5mV6XIyIiIuI8NaUBUj5XWo9P4Wvex03K\nxU3KxV3Kxk3KJfyoKQ0QXYEvIiIi4j/NlAbIBxs+4O11b/PTTj/lZ6f+zOtyRERERAJOM6UO0pFS\nEREREf+pKQ2Qlo01U6p5HzcpFzcpF3cpGzcpl/CjpjRAWjRqgTGGfbn7KCwu9LocEREREadppjSA\nJs6eyL7cfdyTeg+tm7T2uhwRERGRgNJMqaP0zU4iIiIi/lFTGkDlc6XZ9bMp1byPm5SLm5SLu5SN\nm5RL+FFTGkC6Al9ERETEP5opDaB1+9bx94V/p3Ozzvx54J+9LkdEREQkoDRT6igdKRURERHxj5rS\nAGoa3ZToyGhy8nPIyc/xupyg07yPm5SLm5SLu5SNm5RL+FFTGkDGGFo3LrkV1K6sXR5XIyIiIuIu\nzZQG2LSV0/hi2xeM7DGS8zqe53U5IiIiIgGjmVKHJTdJBmBn1k6PKxERERFxl5rSAEuOK21KM+tf\nU6p5HzcpFzcpF3cpGzcpl/CjpjTAKh4pdXXEQERERMRrmikNgj99/Ccy8zL56+C/khSb5HU5IiIi\nIgGhmVLH1edT+CIiIiL+UFMaBG3j2gL172Inzfu4Sbm4Sbm4S9m4SbmEHzWlQVA+V6ojpSIiIiKV\n0kxpEGw/tJ37PruPVo1bce9593pdjoiIiEhAaKbUca0at8JnfHyf8z0FRQVelyMiIiLiHDWlQRAV\nEUWLRi2w1pKene51OUGjeR83KRc3KRd3KRs3KZfwo6Y0SMovdtJcqYiIiMgxNFMaJO9veJ931r3D\nBZ0v4H9O+R+vyxERERGpc5opDQFlV+DvyNzhcSUiIiIi7lFTGiT18Qb6mvdxk3Jxk3Jxl7Jxk3IJ\nP2pKgySxYSLRkdFk5mWSlZfldTkiIiIiTtFMaRA9uOBBNh/czB1n3UH3pO5elyMiIiJSpzRTGiLq\n4yl8EREREX+oKQ2i8ttCZdWPplTzPm5SLm5SLu5SNm5SLuFHTWkQlV2BryOlIiIiIkfSTGkQ5Rbk\ncseHdxAVEcUTFz+Bz+jfBCIiIhI+nJ4pNcYMMcasM8asN8aMq+T1YcaYlcaY5caYxcaYgf7uG2pi\no2JJaJhAQVEB+3L3eV2OiIiIiDMC2pQaY3zAP4CLgFOBa4wxR192/qm1tre1ti9wE/B8NfYNOfXp\nFL7mfdykXNykXNylbNykXMJPoI+U9gc2WGu3WmsLgBnA5RU3sNbmVlg2Bor93TcUlV2Br292EhER\nEflRoJvSZGB7hfWO0ueOYIwZboxZC8wEbqzOvqGm/EhpPbgCPzU11esSpBLKxU3KxV3Kxk3KJfxE\nel0AgLUT7BFmAAAgAElEQVT2beBtY8w5wH3ABdV9j1GjRpGSkgJAfHw8ffr0Kf8DW3aI34V1+6bt\n2fXtLrLXZ/Or03/leT1aa6211lprrbXWNV2vWLGCjIwMANLS0qiNgF59b4wZANxjrR1Suh4PWGvt\ng8fZZxNwBtDV331D5ep7gGJbzB0f3sHhwsM8fOHDxEXHeV1SwMybN6/8D664Q7m4Sbm4S9m4Sbm4\nyeWr75cAXYwxHYwxDYCRwLsVNzDGdK7w82lAA2vtAX/2DUU+46NDfAcA0jLSvC1GRERExBEBv0+p\nMWYI8DglDfAL1toHjDFjKDnq+Zwx5s/ADUA+8APwR2vtwqr2reIzQuZIKcBba9/io40fMbTrUIZ1\nG+Z1OSIiIiJ1ojZHSgM+U2qt/RDodtRzz1b4+SHgIX/3DQcd4zsCsOXgFo8rEREREXFDoE/fSyVS\n4lOAktP3oXSEt7rKBqLFLcrFTcrFXcrGTcol/Kgp9UB8TDxNY5qSW5DL3ty9XpcjIiIi4rmAz5QG\nQ6jNlAI8veRpVu5eyY19b+TMtmd6XY6IiIhIrbl89b1UoWyuVFfgi4iIiKgp9UzFudJwpXkfNykX\nNykXdykbNymX8KOm1CNl9yrddmgbhcWFHlcjIiIi4i3NlHpo0txJfJ/9PRN+MqG8SRUREREJVZop\nDVGaKxUREREpoabUQ+E+V6p5HzcpFzcpF3cpGzcpl/CjptRDHRNKv9kpQ9/sJCIiIvWbZko9VFhc\nyG0f3EaxLeaxIY8RExnjdUkiIiIiNaaZ0hAV6YukXVw7rLVszdjqdTkiIiIinlFT6rFwnivVvI+b\nlIublIu7lI2blEv4UVPqMc2VioiIiGim1HO7s3czee5k4mPiefCCB70uR0RERKTGNFMawlo2aklM\nZAwZhzPIOJzhdTkiIiIinlBT6jFjTNjOlWrex03KxU3KxV3Kxk3KJfyoKXVA2Vzp5oObPa5ERERE\nxBuaKXXAqj2rePKrJ+mU0Ilx54zzuhwRERGRGtFMaYjr0qwLPuMjLSONvMI8r8sRERERCTo1pQ6I\niYyhfdP2FNtiNh7Y6HU5dUbzPm5SLm5SLu5SNm5SLuFHTakjuiV1A2D9/vUeVyIiIiISfJopdYTm\nSkVERCTUaaY0DGiuVEREROozNaWOCMe5Us37uEm5uEm5uEvZuEm5hB81pQ7RXKmIiIjUV5opdYjm\nSkVERCSUaaY0TGiuVEREROorNaUOCbe5Us37uEm5uEm5uEvZuEm5hB81pY7RXKmIiIjUR5opdczq\nPat54qsnNFcqIiIiIUczpWGkc7POmisVERGRekdNqWNiImPoEN8hLOZKNe/jJuXiJuXiLmXjJuUS\nftSUOqhrYlcAvtv/nceViIiIiASHZkodVDZX2jGhI+PPGe91OSIiIiJ+0UxpmCmbK92asZXDhYe9\nLkdEREQk4NSUOqjiXOmmA5u8LqfGNO/jJuXiJuXiLmXjJuUSftSUOkpzpSIiIlKfaKbUUWVzpSnx\nKdz5kzu9LkdERETkhDRTGoY6N+tMhC+CrYe2kpOf43U5IiIiIgHlV1NqjDnHGDO69OfmxpiOgS1L\nYiJjOKnZSVhrWbVnldfl1IjmfdykXNykXNylbNykXMLPCZtSY8xkYBxQdg45CnjV3w8wxgwxxqwz\nxqw3xhzzvZnGmGuNMStLHwuMMb0qvJZW+vxyY8xifz8zXPRqWfJfxbd7vvW4EhEREZHAOuFMqTFm\nBdAXWGat7Vv63DfW2l7H3bFkOx+wHhgM7AKWACOttesqbDMAWGutPWSMGQLcY60dUPraZqCftfbg\nCT4n7GZKAfbm7OWuOXcRGxXLIxc+QoQvwuuSRERERKoU6JnS/NKOz5Z+WKNqvH9/YIO1dqu1tgCY\nAVxecQNr7SJr7aHS5SIgucLLxs8aw1LzRs1p1bgVuQW5bDoYureGEhERETkRfxq+fxtjngXijTG3\nAJ8Cz/v5/snA9grrHRzZdB7tZuCDCmsLfGKMWVL62fVO+Sn870PvFL7mfdykXNykXNylbNykXMJP\n5Ik2sNY+Yoy5AMgEugGTrLWf1HUhxpjzgNHAORWeHmitTTfGNKekOV1rrV1Q2f6jRo0iJSUFgPj4\nePr06UNqairw4x/cUFz3bNmTf739L97d+C4jThnheT3VWZdxpR6tS9YrVqxwqh6tS9ZlXKlH6x/X\nK1ascKoerbV2ab1ixQoyMjIASEtLozb8mSl90Fo77kTPVbHvAEpmRIeUrscD1lr74FHb9QL+Awyx\n1lZ6nrr0gqssa+2jlbwWljOlAEXFRfzx4z+SW5DLlPOn0KJRC69LEhEREalUoGdKL6jkuYv9fP8l\nQBdjTAdjTANgJPBuxQ2MMe0paUivr9iQGmNijTGNS39uBFwIhOa9kWohwhfBqS1OBULzFL6IiIiI\nP6psSo0xvzbGfAt0M8Z8U+GxBfjGnze31hYBY4GPgdXADGvtWmPMGGPML0s3uxtoBjx91K2fWgIL\njDHLKbkAaqa19uMa/ZYhrmeLnkDo3Rqq7DC/uEW5uEm5uEvZuEm5hJ/jzZROp+Sio/uB8RWez7LW\nHvD3A6y1H1Iyi1rxuWcr/HwLcMxFTNbaLUAffz8nnPVo0QNjDOv3r+dw4WFiImO8LklERESkTp1w\nprR8Q2NaAOXdkLV2W6CKqq5wnikt8/AXD7PxwEbGnD6G01qf5nU5IiIiIscI6EypMeYyY8wGYAsw\nH0jjyNs2SRCU3Rrqm+/9mpwQERERCSn+XOh0HzAAWG+t7UjJtzMtCmhVcoyeLUvmSlftWUWxLfa4\nGv9o3sdNysVNysVdysZNyiX8+NOUFlhr9wM+Y4zPWjsXOD3AdclRWjduTVJsEll5WaRlpHldjoiI\niEid8uc+pZ8Cwym54CkJ2AOcYa09O/Dl+ac+zJQCzFg1g7lb5nLJSZdweffLT7yDiIiISBAF+j6l\nlwO5wB3Ah8Am4LKafJjUjuZKRUREJFwdtyk1xkQA71lri621hdbal621T5Sezpcg65rYlejIaHZk\n7uDgDwe9LueENO/jJuXiJuXiLmXjJuUSfo7blJbe/L7YGNM0SPXIcUT6Ijml+SlA6N1IX0REROR4\n/JkpfQfoC3wC5JQ9b629LbCl+a++zJQCfLHtC6atnEbPlj0Z23+s1+WIiIiIlKvNTOnxvtGpzFul\nD3FA2a2h1u5dq293EhERkbBxwgudSudIj3kEozg5Vlx0HF0Tu1JYXMiy9GVel3Ncmvdxk3Jxk3Jx\nl7Jxk3IJP/5cfS+O6Z/cH4DFOxd7XImIiIhI3TjhTGkoqE8zpQC5Bbn86eM/UWSLePCnD9I0Rteh\niYiIiPcCep9SY8zP/HlOgic2KpYeLXpgrWXJriVelyMiIiJSa/6cvr/Tz+ckiM5seyYAX+34yuNK\nqqZ5HzcpFzcpF3cpGzcpl/BT5dX3xpiLgUuAZGPMExVeigMKA12YHF/PFj2JiYxh26Ft7M7eTavG\nrbwuSURERKTGqpwpNcb0BvoAfwEmVXgpC5hrrXXmK4Xq20xpmWkrp/HFti+45KRLuLz75V6XIyIi\nIvVcQGZKrbUrS2/91KXCbaDeBTa61JDWZxWvwq+PTbmIiIiED39mSj8xxsQZY5oBy4Cpxpi/B7gu\n8UPXxK7Ex8SzL3cfWzK2eF3OMTTv4ybl4ibl4i5l4yblEn78aUqbWmszgSuBadbaM4HBgS1L/OEz\nvvKjpS5f8CQiIiJyIie8T6kx5lvgQuBlYKK1dokx5htrba9gFOiP+jpTCrD90Hbu++w+GjdozEMX\nPESEL8LrkkRERKSeCuh9Sim50OkjYFNpQ9oJ2FCTD5O61zauLa2btCY7P5s1e9d4XY6IiIhIjZyw\nKbXWvmGt7WWt/XXperO1dkTgSxN/GGM4M7n0nqU73TqFr3kfNykXNykXdykbNymX8OPPNzp1NcbM\nNsasKl33MsbcFfjSxF9lc6Urd68krzDP42pEREREqs+fmdL5wJ+AZ621fUufW2Wt7RGE+vxSn2dK\nyzz8xcNsPLCRG/veWP5tTyIiIiLBFOiZ0lhr7eKjntM3Ojmm/GtHHTuFLyIiIuIPf5rSfcaYzoAF\nMMb8D5Ae0Kqk2vq17ofP+Fizdw2ZeZlelwNo3sdVysVNysVdysZNyiX8+NOU/hZ4FuhujNkJ3A78\nKqBVSbU1atCIHi16YK3l611fe12OiIiISLUcd6bUGOMD/sda+29jTCPAZ63NClp1ftJMaYmvd33N\n1KVTSYlP4c6f3Ol1OSIiIlLPBGym1FpbDPy59OccFxtS+VHvlr2JiYwhLSONPTl7vC5HRERExG/+\nnL7/1BjzR2NMO2NMs7JHwCuTaouKiOK01qcB8PnWzz2uRvM+rlIublIu7lI2blIu4cefpvRqSuZK\nPwOWlj40tOio1JRUABZsW6B7loqIiEjIOOF9SkOBZkqP9NAXD7HpwCau63Ud53Y41+tyREREpJ4I\n6H1KjTGxxpi7jDHPla5PMsZcWpMPk+A4v+P5AMzZMgc16yIiIhIK/Dl9/xKQD5xdut4J3BewiqTW\n+rbqS3xMPOlZ6azbt86zOjTv4ybl4ibl4i5l4yblEn78aUo7W2sfAgoArLW5QI0Oy0pwRPgiymdL\nZ2+Z7W0xIiIiIn444UypMeZLYDDwhbX2tNJvd3rdWts/GAX6QzOlx8rOz2bcJ+MoLC5kyvlTaNGo\nhdcliYiISJgL6EwpcA/wIdDOGPMaMBsYV5MPk+Bp3KAxZ7Y9E4C5W+Z6XI2IiIjI8Z2wKbXWfgxc\nCYwCXgdOt9aqywkBZRc8fbn9Sw4XHg7652vex03KxU3KxV3Kxk3KJfz4c/X9bGvtfmvtLGvte9ba\nfcYYDSqGgLZxbema2JXDhYdZuH2h1+WIiIiIVKnKmVJjTAwQC8wFUvnx4qY44ENrbXe/PsCYIcBj\nlDTAL1hrHzzq9Wv5cRwgC/iNtfYbf/at8B6aKa3C8vTlPPP1M7Ro1IK/nPcXjNE1aiIiIhIYgZop\nHUPJtzd158dvcloKvAP8w8/CfKXbXgScClxjjDm6md0MnGut7U3Jraaeq8a+cgK9W/UmMTaRPTl7\nWL13tdfliIiIiFSqyqbUWvu4tbYj8EdrbSdrbcfSR29rrV9NKdAf2GCt3WqtLQBmAJcf9TmLrLWH\nSpeLgGR/95UT8xnfj7eH2hzcqQvN+7hJubhJubhL2bhJuYQffy50etIYc7Yx5lpjzA1lDz/fPxnY\nXmG9gx+bzsrcDHxQw32lCue0P4cGEQ1Ys3cN6VnpXpcjIiIicozIE21gjHkF6AysAIpKn7bAtLos\nxBhzHjAaOKcm+48aNYqUlBQA4uPj6dOnD6mpqcCP/5qqr+vFXywmLj2OfS32MTdtLm32t3GqPq2D\nuy57zpV6tNY6FNZlXKlH61RSU1Odqqe+rlesWEFGRgYAaWlp1IY/N89fC5xSkyuJjDEDgHustUNK\n1+MBW8nFTr2A/wBDrLWbqrNv6Wu60OkE0rPSuWfePTSIaMCDFzxIbFSs1yWJiIhImAn0zfNXAa1q\n8ubAEqCLMaaDMaYBMBJ4t+IGxpj2lDSk15c1pP7uK/5r3aQ1Jzc/mfyifL7Y9kVQPvPoIwziBuXi\nJuXiLmXjJuUSfvxpSpOANcaYj4wx75Y9/Hlza20RMBb4GFgNzLDWrjXGjDHG/LJ0s7uBZsDTxpjl\nxpjFx9u3Wr+dHKHsZvpz0+ZSbIs9rkZERETkR/6cvh9U2fPW2vkBqagGdPreP9Za7p57N3tz9vLr\nM35Nn1Z9vC5JREREwkhAT99ba+dX9qjJh4m3jDGcl3IeAJ9s+gQ18iIiIuKKKptSY8yC0v/MMsZk\nVnhkGWMyg1ei1KWB7QfSqEEjNh7YyJq9awL6WZr3cZNycZNycZeycZNyCT9VNqXW2nNK/7OJtTau\nwqOJtTYueCVKXYqJjGFIlyEA/Hfdf3W0VERERJxwwpnSUKCZ0uopKCrgrjl3kXE4g1/2+yX92vTz\nuiQREREJA4G+JZSEmaiIKIZ2HQrAO9+9oyvxRURExHNqSuupge0G0rxRc77P/p6F2xcG5DM07+Mm\n5eIm5eIuZeMm5RJ+1JTWUxG+CIZ1GwbAzPUzKSgq8LgiERERqc80U1qPWWu577P72JG5g6tOvYrB\nnQZ7XZKIiIiEMM2USo0YY7i8++UAvL/hfQ4XHva4IhEREamv1JTWcz1b9KRzs85k52fz6eZP6/S9\nNe/jJuXiJuXiLmXjJuUSftSU1nPGGIZ3Hw6UfMtTTn6OxxWJiIhIfaSZUgHgia+eYPWe1VzY+UJG\nnDLC63JEREQkBGmmVGqt7Gjp3LS5ZBzO8LgaERERqW/UlAoA7Zu2p1+bfhQUFfDe+vfq5D017+Mm\n5eIm5eIuZeMm5RJ+1JRKuWHdhmGM4YttX7AnZ4/X5YiIiEg9oplSOcK0ldP4YtsXnJF8BjefdrPX\n5YiIiEgI0Uyp1JlLu15KpC+SJTuXsCNzh9fliIiISD2hplSO0KxhMwalDALg7XVv1+q9NO/jJuXi\nJuXiLmXjJuUSftSUyjEu7nIx0ZHRfPv9t2w6sMnrckRERKQe0EypVOrd795l1vpZdG7WmT+d/SeM\nqdF4iIiIiNQjmimVOndBpwuIi45j04FNzN863+tyREREJMypKZVKNYxqyLU9rwXgrbVvsT93f7Xf\nQ/M+blIublIu7lI2blIu4UdNqVSpb+u+nNb6NPIK83jt29fQiISIiIgEimZK5bgy8zKZPHcyuQW5\njO47mgFtB3hdkoiIiDhKM6USMHHRcVzd42oA/m/V/5GZl+lxRSIiIhKO1JTKCZ2ZfCantjiV3IJc\nXv/2db/307yPm5SLm5SLu5SNm5RL+FFTKidkjOG6ntcRHRnNsvRlLE9f7nVJIiIiEmY0Uyp+m5c2\nj9e/fZ246DjuPe9eYqNivS5JREREHKKZUgmKQR0G0aVZFzLzMnlj9RtelyMiIiJhRE2p+M0Yww29\nbyDSF8mX279k9Z7Vx91e8z5uUi5uUi7uUjZuUi7hR02pVEvLxi0Z1m0YAK9+8yqHCw97XJGIiIiE\nA82USrUV22Lu//x+th3axnkdz2Nkj5FelyQiIiIO0EypBJXP+PhFn1/gMz7mpc1j44GNXpckIiIi\nIU5NqdRI27i2DOkyBGst01ZOo6Co4JhtNO/jJuXiJuXiLmXjJuUSftSUSo0N7TqUVo1b8X3298za\nMMvrckRERCSEaaZUamXzwc089MVDGAx3/uRO2jdt73VJIiIi4hHNlIpnOiV04vyO51Nsi5m2chpF\nxUVelyQiIiIhSE2p1Nrl3S4nMTaR7Ye28/Gmj8uf17yPm5SLm5SLu5SNm5RL+FFTKrUWHRnN9b2u\nB+C99e+x+eBmjysSERGRUKOZUqkzM1bNYO6WucRFxzHhJxNIaJjgdUkiIiISRJopFSf87JSf0TWx\nK5l5mfzz639WepsoERERkcoEvCk1xgwxxqwzxqw3xoyr5PVuxpgvjTGHjTG/P+q1NGPMSmPMcmPM\n4kDXKrUT4YtgzOljSIxNZGvGVia8MAEdwXaP5rDcpFzcpWzcpFzCT0CbUmOMD/gHcBFwKnCNMab7\nUZvtB24FHq7kLYqBVGttX2tt/0DWKnWjcYPG/PaM3xIdGc26feuOuPBJREREpCoBnSk1xgwAJltr\nLy5djwestfbBSradDGRZax+t8NwW4HRr7f4TfI5mSh2zPH05z3z9DMYYxvYfS48WPbwuSURERALM\n5ZnSZGB7hfWO0uf8ZYFPjDFLjDG31GllElB9W/flsm6XYa1l6tKp7M7e7XVJIiIi4rBIrws4gYHW\n2nRjTHNKmtO11toFlW04atQoUlJSAIiPj6dPnz6kpqYCP86daB3cdSPbiNNan8Z7H7/Hn9f9mWfG\nPkNsVKwz9dXX9WOPPaa/Hw6uy55zpR6tf1yvWLGC22+/3Zl6tD7y74or9dTX9YoVK8jIyAAgLS2N\n2gjG6ft7rLVDStfVOn3v7+s6fe+mefPmcdY5Z/HQFw+xI3MHp7Y4lbH9x+IzPq9Lq9fmzZtX/j8o\n4g7l4i5l4ybl4qbanL4PdFMaAXwHDAbSgcXANdbatZVsOxnIttb+b+k6FvBZa7ONMY2Aj4F7rbXH\nXDmjptRt+3L38bfP/0ZOfg4Xdr6QEaeM8LokERERCQBnm1IouSUU8Dgl86svWGsfMMaMoeSI6XPG\nmJbA10ATSq62zwZOAZoD/6VkrjQSeM1a+0AVn6Gm1HHr96/n7wv/TrEtZnTf0QxoO8DrkkRERKSO\nuXyhE9baD6213ay1J5U1ldbaZ621z5X+/L21tp21Nt5a28xa295am22t3WKt7VN6O6ieVTWk4q6K\n8z5dE7syssdIAF5Z+QppGWneFCVH5CLuUC7uUjZuUi7hR8N9EjSDUgZxbodzKSwu5OklT5NxOMPr\nkkRERMQRAT99Hww6fR86CosLeWzRY2zYv4GU+BT+ePYfiYqI8rosERERqQNOn74XqSjSF8mYfmNo\n1rAZaRlpvPrNq/oqUhEREVFTKoFT1bxPk+gm/OaM39AgogGLdixi9pbZwS2sntMclpuUi7uUjZuU\nS/hRUyqeaNe0HaP6jALgzTVvsnrPam8LEhEREU9pplQ89e537zJr/Sxio2K58yd30qJRC69LEhER\nkRrSTKmErMu6XkafVn3ILcjlqcVPkZ2f7XVJIiIi4gE1pRIw/sz7GGMY3Xc0bZq0YXf2bh758hEO\nHT4U+OLqMc1huUm5uEvZuEm5hB81peK5mMgYbh9wO22atCE9K52Hv3yY/bn7vS5LREREgkgzpeKM\n7PxsHl/0ONsObSOhYQK/P+v3mjEVEREJIbWZKVVTKk75oeAHnlz8JJsObCIuOo47zrqDNk3aeF2W\niIiI+EEXOomTajLv0zCqIb8783d0T+pOZl4mj3z5CFszttZ9cfWY5rDcpFzcpWzcpFzCj5pScU50\nZDRj+4+lZ8ue5OTn8OjCR9l0YJPXZYmIiEgA6fS9OKuwuJAXl7/I0l1LaRDRgN/2/y3dk7p7XZaI\niIhUQTOlakrDVrEtZtrKaSzcvpBIXyRjTh9Dr5a9vC5LREREKqGZUnFSXcz7+IyPX/T+BakpqRQW\nF/LPJf9k6a6ltS+uHtMclpuUi7uUjZuUS/hRUyrOM8YwssdILux8IcW2mKnLpvLl9i+9LktERETq\nkE7fS8iw1vL+hvd597t3Abim5zWkpqR6W5SIiIiU00ypmtJ65dPNn/LG6jcAGHHKCC7sfKHHFYmI\niAhoplQcFah5n592+inX9boOYwz/WfMfZn43E/2jxH+aw3KTcnGXsnGTcgk/akolJJ3b4VxG9RmF\nMYb31r/Hm2veVGMqIiISwnT6XkLasvRlPL/seYqKizi3w7lc2/NajKnRWQMRERGpJc2Uqimt11bt\nWcUzXz9DQVEBA9oO4Bd9foHP6CSAiIhIsGmmVJwUrHmfHi16cGv/W4mOjGbRjkVMXTqVwuLCoHx2\nKNIclpuUi7uUjZuUS/hRUyphoVtSN24fcDsNoxqyLH0Z/1zyTwqKCrwuS0RERPyk0/cSVrYf2s5j\nix4jOz+bdk3bcfNpN9OqcSuvyxIREakXNFOqplQqSM9K5x+L/8G+3H00iGjA1T2uZmC7gboASkRE\nJMA0UypO8mrep3WT1tw96G7ObHsm+UX5vLLyFZ5b+hy5Bbme1OMazWG5Sbm4S9m4SbmEHzWlEpZi\nImO4se+N3Nj3RqIjo1mWvowp86ew8cBGr0sTERGRSuj0vYS9vTl7eX7Z86RlpGGMYehJQxnadahu\nGyUiIlLHNFOqplROoKi4iJnrZ/Lhxg+x1tK5WWdu6nsTibGJXpcmIiISNjRTKk5yad4nwhfB8O7D\nuX3A7TSNacqmA5uY8tkUlu5a6nVpQedSLvIj5eIuZeMm5RJ+1JRKvdI9qTuTBk2iV8te/FDwA88t\nfY5pK6eRV5jndWkiIiL1mk7fS71krWX+1vm8sfoNCosLadm4JbecdgvtmrbzujQREZGQpZlSNaVS\nQzszd/L8sufZlbWLSF8kV558Jed3PF/3NBUREakBzZSKk0Jh3ic5LpkJP5nAoJRBFBYX8u/V/+Yf\ni/9BVl6W16UFTCjkUh8pF3cpGzcpl/CjplTqvaiIKK7teS2/PuPXNGrQiFV7VvGX+X9hzd41Xpcm\nIiJSb+j0vUgFB384yIvLX2T9/vUAXND5AoZ3H06kL9LjykRERNynmVI1pVKHim0xH238iHe/e5di\nW0z7pu25+bSbadm4pdeliYiIOE0zpeKkUJ338RkfF590MX8a+CeSYpPYdmgbf/38r3y5/UvC4R8/\noZpLuFMu7lI2blIu4SfgTakxZogxZp0xZr0xZlwlr3czxnxpjDlsjPl9dfYVCaROCZ2469y76J/c\nn7zCPF5e8TLPL3ue3IJcr0sTEREJOwE9fW+M8QHrgcHALmAJMNJau67CNklAB2A4cNBa+6i/+1Z4\nD52+l4Cx1vLVzq+Y/u108grzSIxN5Ka+N9G5WWevSxMREXGKy6fv+wMbrLVbrbUFwAzg8oobWGv3\nWWuXAoXV3VckGIwxDGg7gLvOvYsO8R3Yn7ufR758hFnrZ1Fsi70uT0REJCwEuilNBrZXWO8ofS7Q\n+4oDwm3ep0WjFvx54J+5qMtFFNti3v3uXR5d+CjpWelel1Yt4ZZLuFAu7lI2blIu4UcXOolUQ9m3\nPt0+4HbiouPYsH8D986/l5dXvMz+3P1elyciIhKyAn3zxZ1A+wrrtqXP1fm+o0aNIiUlBYD4+Hj6\n9OlDamoq8OO/prTWui7XkwZNYub6mcx4bwZvfvMmX+38inM7nEuTXU1o1KCR5/VVtS57zpV6tNY6\nFNZlXKlH61RSU1Odqqe+rlesWEFGRgYAaWlp1EagL3SKAL6j5GKldGAxcI21dm0l204Gsq21/1uD\nfUvlq70AAB2rSURBVHWhk3hmb85eZq6fyeKdi7HW0iCiAYM7DebCzhcSGxXrdXkiIiJB4+yFTtba\nImAs8DGwGphhrV1rjBljjPklgDGmpTFmO3AHMNEYs80Y07iqfQNZr9Sto48whKvmjZpzY98bufvc\nu+ndqjf5Rfl8sOEDJs6eyIcbPyS/KN/rEo9QX3IJNcrFXcrGTcol/AT8uxOttR8C3Y567tkKP38P\ntPN3XxFXJccl85szfsPmg5t5e93bfLfvO/679r/M3jyboV2Hck77c/R1pSIiIlXQ14yKBIC1lnX7\n1vH2urdJy0gDICk2icu6XUb/5P74TEBPUoiIiHiiNqfv1ZSKBJC1lhW7V/DOd++U3zqqTZM2XN79\ncnq37I0xNfp7KyIi4iRnZ0qlftO8T8lfzr6t+zJp0CRG9x1NYmwiu7J28c8l/+SBBQ+wbt8xX1AW\ncMrFTcrFXcrGTcol/GjATSQIfMbHgLYDOL3N6SzYtoBZ62eRlpHG3xf+ne5J3RnefTgdEzp6XaaI\niIhndPpexAN5hXnMTZvLRxs/IrcgF4DerXozvPtw2jRp43F1IiIiNaOZUjWlEqJyC3L5eNPHzN48\nm/yifIwx9E/uz7Buw0iKTfK6PBERkWrRTKk4SfM+JxYbFcvw7sP56+C/cl7H8/AZH1/t+Iq759zN\n9G+nc+jwoTr/TOXiJuXiLmXjJuUSftSUijggLjqOkT1GMuW8KZzV7iwslvlp85k4ZyJvrX2LnPwc\nr0sUEREJKJ2+F3FQelb6/7d3r8GNnfd9x79/ECDBO8i98mLtkpRW2pVXWt2W1FiRZKlVtO5FiZ24\nTtK4dp2pWzdNM56p6unUdZy4E7fTycRuJskokxfyjDN2nYx1cerEVi3asuSlVnvTxitppeXeSO6d\nAO8ESeDpiwOQ4B2rJXgOiN9n5gzOOTgAH+i/IH96znOew/PvPM/Ri0cBqIxU8kTHEzze9jgV4Qqf\nWyciIrI0jSlVKJUN6mziLM+9/RxvXfXusFtbUctHbvsID+94WHeHEhGRwNGYUgkkjfe5eTtjO/nd\nrt/l8w9+nvaGdkaSI3z7H77NF3/0RV678Bppl77h91Rdgkl1CS7VJphUl41HoVSkCNy++Xae/tDT\nfO6Bz9Fc28zgxCDPHnuWL3d/mSMXj6AzBSIiUux0+l6kyKRdmkP9h3jhnRe4Nn4NgB2xHTx1+1Ps\n2bJHty4VERHfaEypQqmUoFQ6xasXXuV7p743O3XUrk27OHDbAe7YfAch04kQERFZXxpTKoGk8T6F\nVRYq4+EdD/OVx77Cx/Z8jOryak5dP8XXDn6NL7z0Bb7z8+9wYejColP7qkswqS7BpdoEk+qy8ejy\nXZEiV15WzhMdT/ALt/wCL599mZ9d+BlXxq7wUu9LvNT7Ek21TXS2dLK/ZT+bqjb53VwREZEl6fS9\nyAbjnONs4iw9/T0c6j/E6NTo7HO3bbqNzpZO7mu+j6pIlY+tFBGRjUhjShVKRZaUSqc4efUkPf09\nHLt0jOnUNADhUJi92/bS2dLJ3m17NeepiIisCYVShdJA6u7u5tFHH/W7GZIxOTPJ0YtH+eaL32Tm\nlpnZsaZVkSrua76PzpZObm28VVfv+0Tfl+BSbYJJdQmmmwml6h4RKRHRcJQHP/AgyT1J9nXt41D/\nIXr6e7gwdIFXzr3CK+deYVPVJva37KezpZOm2ia/mywiIiVEPaUiJW5gZICevh56+nuIT8Rn93+g\n/gN0tXbxQPMD1EfrfWyhiIgUC52+VygVuWnOOd4dfJeevh4OXzzMxPQE4P2CuWPzHXS1drFv+z6i\n4ajPLRURkaBSKFUoDSSN9wmmfOoynZrmxJUT9PT1cOLKCVLpFACRsgj7tu+jq7WL3Zt3UxYqW4cW\nlwZ9X4JLtQkm1SWYNKZURNZUpCzCvU33cm/TvYxNjXH44mF6+np4b/A9DvUf4lD/IWorarm/+X66\nWrvYUb9DF0iJiMhNUU+piOTt2vg1Xu9/nYN9B7k8enl2/9bqrXS2dtLZ0smW6i0+tlBERPyk0/cK\npSLryjnH+aHzsxP0DyeHZ59rb2ins7WT+5vvp6a8xsdWiojIelMoVSgNJI33Caa1rkvapXnr6luz\nE/QnZ5IAhCzEB7d+kM7WTu7edjeRssia/cyNSN+X4FJtgkl1CSaNKRUR34QsxJ1b7+TOrXeSnEly\n/PJxDvYd5OTVk7x5+U3evPwm0XCUe5vupbO1k12bdhGykN/NFhGRgFFPqYgUxHBymDcG3qCnr4ez\nibOz+2PRmDdBf2snrXWt/jVQRETWnE7fK5SKBNrl0cv09PfQ09fDtfFrs/tb6lrobOlkf8t+Giob\nfGyhiIisBYVShdJA0nifYPKzLs45euO99PT38MbAG4xNjQHeL7Fdm3bR2dLJPU33UBWp8qV9ftL3\nJbhUm2BSXYJJY0pFpCiYGR2NHXQ0dvDxOz/Oz6/8nJ7+Ho5fOs47197hnWvv8Fcn/oq7tt1FV2sX\nd269k3BIv6ZEREqBekpFxHcT0xMcuXiEnv4eTl0/Rfb7XF1ezX1N99HV2kV7Q7sm6BcRCTidvlco\nFdkw4hPx2Qn6B0YGZvdvrtrM/pb97N6ymx31O6gIV/jYShERWYpCqUJpIGm8TzAVU136hvvo6evh\n9f7XSUwmZveHLERLXQsdDR20N7TT3tDO5qrNRd2TWkx1KTWqTTCpLsGkMaUisiG11rXSuqeVX979\ny5y6foqjF49yOn6a/uF+Lgxd4MLQBbrPdgNQW1FLW6yN9oZ2Oho71JsqIlJk1FMqIkUnOZPk3NA5\nTg+epjfeS2+8l9Gp0XnHZHtTsz2p7Q3tbKnaUtS9qSIiQafT9wqlIiXNOce18WuzAbU33kvfcB9p\nl553XE15zWxAbWtooy3Wpt5UEZE1pFCqUBpIGu8TTKVSl2xvam5QHUmOzDvGzGipnd+burV6qy+9\nqaVSl2Kk2gST6hJMGlMqIrJARbiCXZt2sWvTLsDrTb0+cZ3Tg6c5kzhDb7yXC0MX6Bvuo2+4j5+c\n+wngTUOVG1J3xnYSDUf9/CgiIiVBPaUiUrKmUlOcS8zvTR1ODs87Jrc3ta2hjY6GDt96U0VEgi7Q\np+/N7Engj4EQ8JfOuf+xxDFfBw4AY8CnnXNHM/vPAkNAGph2zu1f5mcolIrITcv2pp6Jn+F0/PRs\nb+rCsanV5dWzV/pnw6p6U0VEAhxKzSwEnAIeBwaAQ8AnnHNv5xxzAPht59w/MbNO4GvOua7Mc73A\nfc65+Co/R6E0gDTeJ5hUlxszlZri/ND52Z7U04Onl+xNbappoqNxbt7UbdXbbqg3VXUJLtUmmFSX\nYArymNL9wLvOuXMAZvYt4Cng7ZxjngK+AeCc6zGzejPb5py7DBheD6uIiC/Ky8q5tfFWbm28FfB6\nUwcnBued8j8/dJ6BkQEGRgZ45dwrAFRFqmhryMyb2tDBzthOKiOVfn4UEZFAK3RP6ceAX3TO/ZvM\n9r8E9jvnfifnmBeBP3TOvZbZfgl42jl3JNNTmgBSwDPOub9Y5ueop1REfDOdmub80PnZU/698V6G\nJofmHZPtTc29iGp7zXaNTRWRDSXIPaU360POuYtmtgX4oZm95Zz76VIHfupTn2Lnzp0AxGIx9u3b\nN9ut393dDaBtbWtb2wXd7mjsoLu7m9vDt3P3P7qb3ngvL/79i1wcvUhZexkDIwO88dobADTvbaYq\nUsVM7wzNtc089eRT7Izt5PVXXw/M59G2trWt7dW2jx07RiLh3Qb67Nmz3IxC95R2Ab/nnHsys/0F\nwOVe7GRmfw687Jz7dmb7beCRzOn73Pf6EjDinPujJX6OekoDqLu7e/YfrgSH6uKPbG9q7mn/xGRi\n9vmBEwO03NXC9prts6f82xraaKppUm+qz4L4nXHOWwDS6cXrQVvS6bl2r7S+0mda+PmOHu3m7rsf\n9f2zFerz5S4dHfAbv7F2/34KKcg9pYeAW81sB3AR+ATwawuOeQH498C3MyE24Zy7bGZVQMg5N2pm\n1cATwJcL3F4RkYKIlEXoaOygo7Fjdl98Is7p+GnOxM/ww7M/pMzKuDhykYsjF3n1/KsAVEYqF13p\nXxWpWvFnZf8AptOQSs2t526nUt5xudv5vib7h3Ot/uAGOUQAvPcevP56sD6fwMAAnDrldyvWRyzm\ndwvWx3pNCfU15qaE+qqZfRZwzrlnMsf8CfAkc1NCHTGzNuC7gMMLz990zn11mZ+hnlKREpcbNpYL\nVfmGMr9eMzUzw2DqPNdSvVxL9XI93cu4iy8KJNXpJurS7dSl2qid7qAy1YRL27yfJxufmb9LKDTX\njpXWsx39q60HbQnS54tGYfPmtf83VAiBnRJqvSiUihRWNlTNzMD09PLrMzPe9lLra7Ws1Ku3ESUt\nzkj4DMNlpxkO9zJadp40M/OOCbsotak2amfaqUt5S7lVEQpBWZn3R3Hheu72cvuXey73D3OQAkW+\nISJ3vZgCU27bs+0TCRqFUoXSQAriOKxilUrB1JS3JJNz6wuXZHL10HjyZDft7Y/eUGhMpfz+L5Af\nsxsPWO8nlBXiNYcPd9PV9eiqr3E2w8DYeS6MnOH8SC9nh0+TSMbnhReApszY1OzSVNtEyEL+FqhI\n6XdZMKkuwRTkMaUiG146vXRYXC085vtcMrm2vYADAzAycuOvKyuDcHjlJRKZO27hejjsbS+1fiPL\nSuGvmHuPBgdh9+58jgyznXbupR3vviSQmEzMu4DqXOIcl0YvcWn0Eq9deA2AaDhKW0MbbbE2Oho7\naIu1UV1eXbDPIyJyo9RTKiUrnYaJCW+ZnJxbz91ebn/u9nr0IppBRQWUl89fltqXb9C70QBZzIGv\n1MykZ+gb7uP0oDdv6pnEGa6PX1903LaabfN6U5trm9WbKiI3RafvFUpLknNeL+LYGIyOekt2fWxs\n9TA5NbU27TBbOSQu99yNHKtQKDcrMZngTPyMd6vU+GnOD51nOjU975iKcMW8K/3bG9rVmyoiN0Sh\nVKE0kG5kvI9zMD6+fMBcbt/N9FKaQWWld1Vj7mN2WW47d3806vUkFlNg1DisYFrvumR7U3NP+y/V\nm7q1eqs3b2pjR8n2puo7E0yqSzBpTKkE1syMN35xaGjxMjw8f/39BMzycqiuhpoab8muV1dDVdXK\nobO8vLjCpMhaCofC7IztZGdsJ4+1PQbA0OTQbEA9kzjD2cRZroxd4crYFQ72HQS83tSW2hYaKxtp\nqGygIdpAQ2UDjZWNxKIx6irqSi60isjaUE+pvG/JJFy/7l2gMTjorScS84Pn2Fj+7xeNzoXL3ICZ\nXV8qfEYihft8IqUulU4t6k29Nn5txdeELERDZQOxaMwLrpnQmhtea8tr0V2qRDYmnb5XKF1zznmn\nx69fnwueCwPo+Pjq7xMKQV2dt9TXe0vueu4+BUyR4BtODnN59DKDE4MkJhMMTgwSn4wTn4gTn4wz\nklx9aoeyUNlsaF0uvCq4ihQnhVKF0vdtbAyuXJlbLl+ee5ycXPm1kQg0Ns4tmzZBQ8Nc0Dx+vJsD\nBx6dnfBZgkHjsIJpo9RlOjVNYjIxG1QXhtb4RJzRqdFV3yccChOLxuYNDVgYXmvKa9YluG6U2mw0\nqkswaUyprMg577T6wAD093uPly97y0qn16NR2LJlLnDmPjY2Qm3tymMy33sPBVKREhMpi7Clegtb\nqrcse8x0anpeUB2cGJwXWuOTccamxrg2fm3F4QLhUHjR0ICFPa7VkWr1uIoUCfWUbjATE9DXNxc+\ns4/LnWovL4etW71l27b5jzU1uhBIRPyRnEnO9rjOGyqQE17Hp1cfQxQpi8wLqgsv0GqINlAVqVJw\nFVkjOn1foqF0fBzOn/eWc+e8xytXlj62uhpaWqC52Xvcvt0LnvX1Cp4iUpySM8lFPa4Lw+vE9MSq\n71NeVr4oqC4Mr5XhSgVXkTwolJZAKHUOLl6Ed9/1ljNn4NoSZ7XCYS905i7Nzd6FROv9+1TjfYJJ\ndQkm1aUwJmcmFw0NWBheJ2dWHkA/cGKAtnvaFoXW2XGumX2Vkcp1+lQC+s4ElcaUbkDOwYULcOrU\nXBBdOP4zEoHWVrjlFtixw3tsbvbu/iMiIhANR2mqbaKptmnZYyZnJpcc15rdd7XsKsmZJJdGL3Fp\n9NKKP2upHtfc8BoNRwvxMUU2BPWUBsjgILz1Fpw8CW+/7U3JlCsWg9tug127oL0dmpoUQEVECs05\nx8TMxJLjWnPD61Rq9XsXR8PRJce15obXinDFOnwqkcLQ6fsiDaXpNJw+DcePw4kTcGnB/4A3NsId\nd8wF0U2bNP5TRCSIssF1uYuyslNjTaemV32vqkjVoqEBC8OrgqsElUJpEYXSZNLrCT12zAuiuafk\no1G4/XbYswd27/YuRCrmEKrxPsGkugST6hJca1Ub5xzj0+PLToOVfcw3uC41rjX3sbys/KbbHGT6\nzgSTxpQG3NiY1xt67JgXSKdzft9s3Qp33+0t7e06HS8islGZGdXl1VSXV9Na17rkMc45xqbHlr3x\nQLYndnx6nPHpcfqG+5b9edXl1ctelJV9jJTpVnoSHOopLaDeXvi7v/N6RNPpuf3t7XNBdPv24u4N\nFRGR9eWcY3RqdMW7ZsUn46TSqVXfq6a8ZsU5XBsqGwiH1H8l+dPp+wCF0nTa6w39wQ/gnXe8faGQ\nd1r+nnu8IBqL+dtGERHZ2JxzjEyNzOtdXSq8pl161feqrahdcQ7XWDSm4CqzFEoDEEoHB+HVV70l\nHvf2RaPw4Q/D4497t+QsNRrvE0yqSzCpLsG1UWvjnGM4ObxoaEBueE1MJvIKrnUVdSve8jUWjVEW\nWtvxaRu1LsVOY0p9MjICR4/CG29484lmc/HmzfDww95SqbmURUQkgMyM+mg99dF6dsZ2LnlM2qW9\n4LrETALZfYnJBMPJYYaTw5zj3LI/q66ibtHQgNzwGovGCFmogJ9Ygk49pe/DkSPw4x97p+ezPzYc\n9k7PP/SQd6pe40RFRKQUpF2aocmhJS/Kyq4PJYdY7e+0mVFfUb8otM5epBVtoD5ar+AacDp9v06h\n1Dn4/vfh+ee97bIyb/qm+++Hu+6CqqqCN0FERKTopNIphpPDi3pZc6fGGk4O5xVcY9HYkj2u2fBa\nV1Gn4OojhdJ1CqXf/a53Nb0ZfPSjXq+ogujyNN4nmFSXYFJdgku1WR+pdIrEZGLJmQSy4XU4OTx7\n/MCJAZr3Ni96n5CFqI/WLxrXmhteFVwLR2NK10ldnXcl/Wc+4/WOioiIyNooC5WxqWoTm6o2LXvM\nTHrGC64TcV6afIn23e2zoTV7kdZI0pt1ID4RX/Z9QhaaHce6XHitLa/FNBZvXamn9AZdueJNeC8i\nIiLBkw2uS901KxteR5Ijq75PWahsNrQuF14VXBfT6fsATAklIiIixWE6NT07c8By4XVsamzV9wmH\nwt4Y1yUuysqG15rympIKrgqlCqWBpHFYwaS6BJPqElyqTTAVui5TqanZ0LpceM03uK40h2tDZQPV\nkeoNE1w1plRERERkDZWXlbO1eitbq5cfs5ecSc5enLVceB2fHufq2FWujl1d9n0iZZEV53BtrGyk\nKrLxr6xWT6mIiIhIgSRnkstOg5V9nJieWPE9Oho7ePpDT69Ti2+OekpFREREAqgiXMH2mu1sr9m+\n7DGTM5PLToMVn4zTVNO0ji32j3pKpWA0DiuYVJdgUl2CS7UJJtUlmG6mp1Qzx4qIiIiI79RTKiIi\nIiJrQj2lIiIiIlLUFEqlYLq7u/1ugixBdQkm1SW4VJtgUl02HoVSEREREfGdxpSKiIiIyJrQmFIR\nERERKWoKpVIwGu8TTKpLMKkuwaXaBJPqsvEUPJSa2ZNm9raZnTKz/7zMMV83s3fN7JiZ7buR10pw\nHTt2zO8myBJUl2BSXYJLtQkm1WXjKWgoNbMQ8CfALwJ3Ar9mZncsOOYA0OGcuw34LPDn+b5Wgi2R\nSPjdBFmC6hJMqktwqTbBpLpsPIXuKd0PvOucO+ecmwa+BTy14JingG8AOOd6gHoz25bna0VERERk\nAyh0KG0BLuRs92X25XNMPq+VADt79qzfTZAlqC7BpLoEl2oTTKrLxhP2uwFLeF/TCJi9r5dJgT37\n7LN+N0GWoLoEk+oSXKpNMKkuG0uhQ2k/cEvOdmtm38JjPrDEMeV5vBbgfc+HJSIiIiLBUOjT94eA\nW81sh5mVA58AXlhwzAvAJwHMrAtIOOcu5/laEREREdkACtpT6pxLmdlvAz/AC8B/6Zx7y8w+6z3t\nnnHO/V8z+4iZvQeMAZ9e6bWFbK+IiIiI+GND3GZURERERIpb0dzR6WYm4ZfCWq02Zna7mb1mZpNm\n9nk/2liK8qjLr5vZ8czyUzPb60c7S00edfnnmZocNbPXzexDfrSz1OR7sxYze8DMps3so+vZvlKW\nx3fmETNLmNmRzPJf/Whnqckzlz2a+V32D2b28qpv6pwL/IIXnt8DdgAR4Bhwx4JjDgB/m1nvBA76\n3e5SWPKszWbgPuAPgM/73eZSWPKsSxdQn1l/Ut+ZwNSlKmd9L/CW3+3e6Es+dck57v8B3wM+6ne7\nS2HJ8zvzCPCC320tpSXPutQDPwdaMtubV3vfYukpvZlJ+KWwVq2Nc+6ac+4wMONHA0tUPnU56Jwb\nymweRPMAr4d86jKes1kDpNexfaUq35u1/Afgr4Er69m4EpdvbTQLz/rKpy6/DvyNc64fvCyw2psW\nSyh9P5Pw9y9xjKw93eQgmG60Lr8FfL+gLRLIsy5m9ktm9hbwIvCv16ltpWzVuphZM/BLzrk/QwFo\nPeX7u+zBzNC9vzWzPevTtJKWT112AY1m9rKZHTKz31ztTYM4eb6IrCMz+zDerBcP+d0W8TjnngOe\nM7OHgK8A/9jnJgn8MZA7bk7BNDgOA7c458bN7ADwHF4gEn+FgXuBx4Bq4Gdm9jPn3HsrvaAY3Mwk\n/FJY+dRG1l9edTGzu4BngCedc/F1alspu6Hvi3Pup2bWbmaNzrnBgreudOVTl/uBb5l3+8DNwAEz\nm3bOaf7swlq1Ns650Zz175vZn+o7U3D5fGf6gGvOuUlg0sx+AtyNNxZ1ScVy+v5mJuGXwrrRmxyo\nd2F9rFoXM7sF+BvgN51zp31oYynKpy4dOev3AuX641pwq9bFOdeeWdrwxpV+ToF0XeTzndmWs74f\nb7pLfWcKK5+//c8DD5lZmZlV4V2EvuJ880XRU+puYhJ+Kax8apP5hfEGUAukzew/Anty/+9W1lY+\ndQG+CDQCf5rp/Zl2zu33r9UbX551+ZiZfRKYAiaAj/vX4tKQZ13mvWTdG1mi8qzNr5jZvwOm8b4z\n/8K/FpeGPHPZ22b298CbQAp4xjl3cqX31eT5IiIiIuK7Yjl9LyIiIiIbmEKpiIiIiPhOoVRERERE\nfKdQKiIiIiK+UygVEREREd8plIqIiIiI7xRKRUQCwsyazOz/ZNbvztwycbXXPGJmLxa+dSIihaVQ\nKiISEM65i8657GT5+4CP5PvSAjVJRGTdKJSKiKwRM/ukmR03s6Nm9qyZ/VMzO2hmh83sB2a2JXPc\nl8zsG2b2mpm9Y2a/ldm/w8xOmFkY+H3g42Z2xMx+1cweyBx/2Mx+ama3+flZRUTWWlHcZlREJOjM\nbA/wX4AHnXNxM4vh3W6vK/P8Z4Cngf+UeclevHtB1wJHzex7mf3OOTdjZv8NuM859zuZ19cADznn\n0mb2OPCHwK+s1+cTESk0hVIRkbXxGPAd51wcwDmXMLMPZsaINgER4EzO8c8756aA62b2I2A/cHyF\n948B38j0kDr0+1tENhidvhcRKZz/DXzdOXcX8G+BaM5zueNAjdXHhf4B8CPn3F7gny14LxGRoqdQ\nKiKyNn4E/KqZNQJkHuuAgczz/2rB8U+ZWbmZbQIeAQ4teH4k8/qsOqA/s/7ptWy4iEgQKJSKiKwB\n59xJ4L8DPzazo8D/An4P+GszOwRcXfCSN4Fu4DXg951zlxY8/zKwJ3uhE/A/ga+a2WH0u1tENiBz\nTjOJiIisJzP7EjDinPsjv9siIhIU+r9tEREREfGdekpFRERExHfqKRURERER3ymUioiIiIjvFEpF\nRERExHcKpSIiIiLiO4VSEREREfHd/wd+VtNrUAYvQQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fe87ff72e48>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "num_points = 20\n",
    "r_vals = np.linspace(0.02, 0.048, num_points)\n",
    "\n",
    "# Compute supply of capital\n",
    "k_vals = np.empty(num_points)\n",
    "for i, r in enumerate(r_vals):\n",
    "    k_vals[i] = prices_to_capital_stock(am,r)\n",
    "\n",
    "# Plot supply and demand of capital\n",
    "fig, ax = plt.subplots(figsize=(11, 8))\n",
    "ax.plot(k_vals, r_vals, lw=2, alpha=0.6, label='supply of capital')\n",
    "ax.plot(k_vals, rd(am, k_vals), lw=2, alpha=0.6, label='demand for capital')\n",
    "ax.grid()\n",
    "ax.set_xlabel('capital')\n",
    "ax.set_ylabel('interest rate')\n",
    "ax.legend(loc='upper right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, the equilibrium interest rate can be found by using the bisection method, and we can see the equilibrium distribution of assets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAqUAAAHuCAYAAACmr3UAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xl8VdXV//HPCgmBICQhlCgIQUAElVFGFUFtFRyqYgdo\n1YJaEUVR0J9Wa9Wn1uEp4oQKIm2xPoq1aqtolVoFVFBEiKBlEgyzQYYwEzLs3x/3JsQYwiXJGXLz\nfb9eeYV97rn3rEuXz7PYZ+19zDmHiIiIiEiQEoIOQERERERERamIiIiIBE5FqYiIiIgETkWpiIiI\niARORamIiIiIBE5FqYiIiIgEzvOi1MwGmdkyM1thZrcd4pzHzWylmWWbWbcyx8eY2ZLoz41exyoi\nIiIiwfC0KDWzBGAicC5wEjDMzDqWO2cw0M45dzwwEpgUPX4ScBXQE+gGXGBmbb2MV0RERESC4fVM\naW9gpXNujXOuAJgOXFTunIuA5wCcc58AqWaWCXQCPnHO5TvnioA5wBCP4xURERGRAHhdlLYE1pUZ\nr48eq+ycDdFjXwD9zSzdzFKA84BWHsYqIiIiIgFJDDqAQ3HOLTOzh4B/A7uBRUBRReeamZ6VKiIi\nIhICzjmryvu8Lko3AK3LjI+NHit/TquKznHO/Rn4M4CZ/YHvzqh+h3OqS8U799xzD/fcc0/QYUgc\nU46JH5Rn4jWzKtWjgPe37z8F2ptZlpnVB4YCr5c753XgCgAz6wvkOedyo+MfRH+3Bi4BXvA4XpEK\n5eTkBB2CxDnlmPhBeSZh5ulMqXOuyMxGAzOJFMBTnXNLzWxk5GX3jHPuLTM7z8y+AvYAI8p8xCtm\n1hQoAK5zzu30Ml4RERERCYbFw21vM3Px8D0kvGbNmsXAgQODDkPimHJM/KA8E6+ZWZV7SlWUioiI\niEiNqE5RqseMisRg1qxZQYcgcU45Jn4IKs/atGmDmeknjn7atGlT43kS2i2hREREJD6sWbNGu+TE\nGavGKvtDfmY8JIlu34uIiIRX9JZu0GFIDTrU/6a6fS8iIiIitZqKUpEYqN9PvKYcEz8ozyTMVJSK\niIiISODUUyoiIiKeCntP6XHHHcfUqVM566yzPPn8Bx54gK+//ppnnnnGk88Pghc9pVp9LyIiIuKh\n3/zmN0GHUCvo9r1IDNSHJV5TjokflGcSZipKRURERKIOHDjATTfdRMuWLTn22GO5+eabKSgoAGDg\nwIG89tprAHz00UckJCTwr3/9C4D33nuP7t27V/iZ9957L5dffjkQ2bM1ISGB5557jqysLJo3b879\n999fem5xcTH3338/7du3JzU1lV69erFhwwYA5s6dS+/evUlPT6dPnz7Mmzev9H1nnnkmd911F6ed\ndhqNGzfmoosuYtu2bVx22WWkpqbSp08f1q5dW3r+smXLOOecc8jIyKBTp068/PLLNfi3WDW6fS8S\nAz0rWrymHBM/hDXPRo6smc+ZPLn6n3Hfffcxf/58Fi9eDMCPf/xj7rvvPu69914GDBjArFmzuOSS\nS5gzZw7t2rVjzpw5DB48mNmzZ1f691t+s/mPPvqIlStXsmzZMnr37s2ll17KCSecwMMPP8xLL73E\n22+/Tfv27VmyZAkpKSls376dCy64gIkTJzJ06FD+9re/cf7557Nq1SrS09MBeOmll5g5cyYZGRn0\n7duXU089laeffprnnnuOESNGcO+99zJ16lT27t3LOeecw3333cc777zD4sWL+dGPfkTnzp3p2LFj\n9f8Sq0gzpSIiIiJRL7zwAnfffTcZGRlkZGRw991389e//hWAAQMGMHv2bADmzJnDb37zm9Lx7Nmz\nGTBgQEzXMDPuuece6tevT5cuXejatSuff/45AFOnTuUPf/gD7du3B6Bz586kp6fz5ptv0qFDB37x\ni1+QkJDA0KFD6dixI2+88Ubp544YMYI2bdrQuHFjBg8eTLt27TjzzDNJSEjgpz/9KYsWLQJgxowZ\nHHfccVxxxRWYGV27dmXIkCGBz5ZqplQkBrNmzQrtDIPEB+WY+CGseVYTM5w1ZePGjbRu3bp0nJWV\nxcaNGwHo168fK1asYPPmzXz++ee88cYb3H333WzdupX58+dzxhlnxHydzMzM0j+npKSwe/duANat\nW0fbtm0rjCsrK+s7x7Kyskpv7Zf/zIYNG35vXHKNNWvW8PHHH9O0aVMAnHMUFRWVthgEJW5mSouL\ng45AREREarsWLVqwZs2a0vGaNWto0aIFECnsTjnlFB577DFOPvlkEhMT6devHxMmTKB9+/alRV51\ntGrVilWrVlUYV05OzneOrV27lpYtW1bpGgMHDmTbtm1s27aN7du3s3PnTp588smqhl0j4qYo3b8/\n6AgknoVxZkHii3JM/KA8O7xhw4Zx3333sWXLFrZs2cLvf//778wgnnHGGUycOLH0Vv3AgQO/M45F\nZXu2Xn311dx111189dVXACxZsoTt27dz3nnnsXLlSqZPn05RUREvvfQSS5cu5cILLzzi73jBBRew\nYsUKnn/+eQoLCykoKGDBggUsW7bsiD+rJqkoFRERkTqt7CKk3/72t/Ts2bO017Nnz57ceeedpa8P\nGDCA3bt3l96qLxkfSVFaftFT2fHYsWP52c9+xjnnnENqaipXX301+/bto2nTpsyYMYPx48fTrFkz\nxo8fz5tvvlm6yKn8Z1bmqKOOYubMmUyfPp0WLVrQokULbr/9dg4cOBDzZ3ghbp7otGGDIzq7LlLj\nwtqHJfFDOSZ+CCrPwv5EJzlyXjzRSTOlIiIiIhK4uJkp/eILx0knBR2JiIiIlKeZ0vijmdJK5OcH\nHYGIiIiIVFXcFKX79gUdgcQzPS9avKYcEz8ozyTM4qYoVU+piIiISO0VNz2lM2Y4zj8/6EhERESk\nPPWUxh/1lFZCM6UiIiIitVfcFKXqKRUvqQ9LvKYcEz8ozyTM4qYo1UypiIiIVMVxxx3He++9V+Fr\nH374IZ06dYrpc2bPnk2rVq1ivu55553HX//615jPj3eJQQdQU1SUipf0pB3xmnJM/KA8O3Knn346\nS5cujfn8I3nc51tvvVWVkOKWZkpFREREJHAqSkVioD4s8ZpyTPygPDu0RYsW0bVrV9LT0xk2bBgH\nDhwAvn9LfuHChfTo0YPU1FR+9rOfMXToUH73u9+Vvu6cY8KECWRmZtKyZUv+8pe/HPKaZ555Jn/6\n058AmDZtGv379+fWW2+ladOmtGvXjrfffrv03O3bt3PllVfSsmVLMjIyGDJkSOlrU6ZM4fjjj6dZ\ns2ZcfPHFbNq0qfS1hIQEnn76aTp06EBqaiq/+93vWL16NaeddhppaWkMHTqUwsLC0vNnzJhB9+7d\nSU9P5/TTT2fJkiVV/0s9QnFz+14LnURERGqnkW+MrJHPmXzh5Cq/9+WXX2bmzJkkJydz6qmn8pe/\n/IVrrrkGOHhLvqCggCFDhnDLLbcwatQoXn/9dYYOHcptt91W+jnffPMNu3btYuPGjcycOZOf/OQn\nXHLJJaSmph42hvnz5zNixAi2bt3K5MmTueqqq9iwYQMAl112GU2aNGHp0qU0atSIuXPnAvDee+9x\nxx138O6773LiiScybtw4hg4dyuzZs0s/d+bMmSxatIi1a9fSvXt35s2bxwsvvEDTpk3p27cvL774\nIpdffjmLFi3iqquu4s033+SUU07h+eef58c//jErVqwgKSmpyn+3sYqbmVI9ZlS8pD4s8ZpyTPyg\nPDu0MWPGkJmZSVpaGhdeeCHZ2dnfO2fevHkUFRUxevRo6tWrxyWXXELv3r2/c079+vW56667qFev\nHoMHD+aoo45i+fLlMcWQlZXFlVdeiZnxq1/9ik2bNrF582a++eYb3nnnHSZPnkyTJk2oV68e/fv3\nB+CFF17gqquuomvXriQlJfHAAw8wb9481q5dW/q5t912G40aNaJTp06cfPLJnHPOOWRlZdG4cWMG\nDx7MokWLgMiM67XXXkvPnj0xMy6//HKSk5P5+OOPq/rXekQ0UyoiIiKBqs4MZ03JzMws/XNKSsp3\nboGX2LRpEy1btvzOsfKr7TMyMkhIODjnl5KSwu7du2OK4eijjy79c8OGDQHYvXs3W7dupWnTpjRp\n0uR779m4cSOnnHJK6bhRo0ZkZGSwYcMGWrduDUDz5s2/87llv2vDhg3Jzc0FYM2aNTz33HM88cQT\nQKQVoaCggI0bN8YUf3XFzUxpYWHkR8QL6sMSrynHxA/Ks+o55phjSm+nl1i3bp3n123VqhXbtm1j\n586d33utRYsWrFmzpnS8Z88etm7dyrHHHlul69x5551s27aNbdu2sX37dnbv3s3Pf/7zasUfq7gp\nSkG38EVERMQ7/fr1o169ejz55JMUFRXxz3/+k/nz53t+3aOPPprBgwdz3XXXkZeXR2FhIR988AEA\nw4YN489//jOLFy8mPz+fO+64g759+x7Rfqklfv3rXzNp0qTS77Rnzx7eeust9uzZU6Pf51DiqijV\nLXzxivqwxGvKMfGD8qxise4tmpSUxKuvvsqzzz5Leno6L7zwAhdeeCHJyclV+uzDXbfs63/9619J\nTEykY8eOZGZm8thjjwFw9tln8/vf/54hQ4bQsmVLvv76a6ZPn37Ia1R2zVNOOYUpU6YwevRomjZt\nSocOHZg2bVqlMdYkc875djGvmJm75hrHXXdBFWarRURExENmRjzUGxXp27cvo0aN4le/+lXQofjq\nUP+bRo/H/gSBMuJqplR7lYpX1IclXlOOiR+UZ9U3Z84ccnNzKSoqYtq0aSxZsoRBgwYFHVZciJvV\n96CiVERERLy1fPlyfvazn7F3717atm3LK6+88p3V7FJ1cXX7/uqroVevoKMRERGRsuL59n1dpdv3\nh6GZUhEREZHayfOi1MwGmdkyM1thZrcd4pzHzWylmWWbWbcyx282sy/MbLGZ/Z+Z1a/sWipKxSvq\nwxKvKcfED8ozCTNPi1IzSwAmAucCJwHDzKxjuXMGA+2cc8cDI4FJ0eMtgBuAHs65LkT6X4dWdj0V\npSIiIiK1k9cLnXoDK51zawDMbDpwEbCszDkXAc8BOOc+MbNUMyvpGK4HNDKzYiAFqPQ5VypKxSva\n20+8phwTPwSVZ1lZWTHvBSq1Q1ZWVo1/ptdFaUug7PO31hMpVCs7ZwPQ0jm30MweBtYCe4GZzrl3\nK7uYNs8XEREJn5ycnKBDkFogtFtCmVkakVnULGAH8Hcz+4Vz7oWKzn///eFs3NiG1ashLS2Nbt26\nlf6LsKSHRmONqzrOzs7mpptuCk08GsffuORYWOLROD7Hjz76qP7/o8Y1Os7OziYvLw+o/j8+PN0S\nysz6Avc45wZFx7cDzjn3UJlzJgHvO+deio6XAQOA/sC5zrlfR49fDvRxzo2u4DrummscJ50EN97o\n2deROmzWrFml/xGKeEE5Jn5QnonXwrwl1KdAezPLiq6cHwq8Xu6c14EroLSIzXPO5RK5bd/XzBpY\npBHlbGBpZRdTT6l4Rf9HXLymHBM/KM8kzDy9fe+cKzKz0cBMIgXwVOfcUjMbGXnZPeOce8vMzjOz\nr4A9wIjoe+eb2d+BRUBB9PczlV1PRamIiIhI7RRXT3TKyID77w86GolHuuUlXlOOiR+UZ+K1MN++\n95VW34uIiIjUTnE1U5qQAE89BdoKTURERMR/mikFEhOhuBgKC4OORERERESOVNwUpQ0bRn7rFr54\noWRvNhGvKMfED8ozCbO4KUqTkyO/tQJfREREpPaJm6K0ZKZURal4QatVxWvKMfGD8kzCLG6K0gYN\nIr9VlIqIiIjUPipKRWKgPizxmnJM/KA8kzCLu6JUC51EREREap+4KUrVUypeUh+WeE05Jn5QnkmY\nxU1Rqtv3IiIiIrWXilKRGKgPS7ymHBM/KM8kzOKuKFVPqYiIiEjtE3dFaX5+sHFIfFIflnhNOSZ+\nUJ5JmMVNUarHjIqIiIjUXnFTlKqnVLykPizxmnJM/KA8kzBTUSoiIiIigVNRKhID9WGJ15Rj4gfl\nmYRZ3BSl2jxfREREpPaKm6I0OTnyWwudxAvqwxKvKcfED8ozCbO4KUrLbgnlXLCxiIiIiMiRMRcH\nFZyZOeccN9wABw7A448fnDkVEREREX+YGc45q8p742amFPRUJxEREZHaKq6KUi12Eq+oD0u8phwT\nPyjPJMziqijVtlAiIiIitZOKUpEYaG8/8ZpyTPygPJMwi8uiVD2lIiIiIrVLXBal+fnBxiHxR31Y\n4jXlmPhBeSZhFldFaclCJ82UioiIiNQucVWUqqdUvKI+LPGackz8oDyTMFNRKiIiIiKBU1EqEgP1\nYYnXlGPiB+WZhFlcFaXqKRURERGpneKqKC153r1mSqWmqQ9LvKYcEz8ozyTM4qoo1WNGRURERGqn\nuCpK1VMqXlEflnhNOSZ+UJ5JmKkoFREREZHAmXMu6Biqzcycc468PLjtNmjSBP74x6CjEhEREalb\nzAznnFXlvXE5U6rHjIqIiIjULnFVlCYng1mkKC0uDjoaiSfqwxKvKcfED8ozCTPPi1IzG2Rmy8xs\nhZnddohzHjezlWaWbWbdosc6mNkiM1sY/b3DzG6s/FraFkpERESkNvK0p9TMEoAVwNnARuBTYKhz\nblmZcwYDo51z55tZH+Ax51zfCj5nPdDHObeuguu4ku9x222QlwcPPgjp6V59MxEREREpL8w9pb2B\nlc65Nc65AmA6cFG5cy4CngNwzn0CpJpZZrlzfgisqqggLU9PdRIRERGpfbwuSlsCZQvJ9dFjlZ2z\noYJzfg68GMsFtS2UeEF9WOI15Zj4QXkmYRb6hU5mlgT8GHg5lvNVlIqIiIjUPokef/4GoHWZ8bHR\nY+XPaVXJOYOBz5xz31Z2oeHDh9OmTRsWLICtW9P48MNunHjiQODgvwxLnvmrscZVGZcISzwaa6yx\nxkc6LjkWlng0rv3j7Oxs8vLyAMjJyaE6vF7oVA9YTmSh0yZgPjDMObe0zDnnAddHFzr1BR4tu9DJ\nzF4E3nbOTavkOqULnaZNg7lz4Yor4LTTPPlaIiIiIlKB0C50cs4VAaOBmcCXwHTn3FIzG2lm10TP\neQv42sy+AiYD15W838xSiCxyejXWa2pLKPFCyb8ORbyiHBM/KM8kzLy+fY9z7m3ghHLHJpcbjz7E\ne/cCPziS65WsvldRKiIiIlJ7eDpTGgQtdBIvlO3HEvGCckz8oDyTMFNRKiIiIiKBi9uiVJvnS01S\nH5Z4TTkmflCeSZjFXVFa0lOanx9sHCIiIiISu7grSjVTKl5QH5Z4TTkmflCeSZjFbVGqnlIRERGR\n2kNFqUgM1IclXlOOiR+UZxJmKkpFREREJHCePmbUL2UfM1pQAKNHQ2IiPPlkwIGJiIiI1CGhfcxo\nEBITISEBCgsjPyIiIiISfnFXlJodvIWvbaGkpqgPS7ymHBM/KM8kzOKuKIWDe5VqWygRERGR2iEu\ni9Lk5MhvLXaSmqK9/cRryjHxg/JMwiwui9KSmVIVpSIiIiK1Q1wWpdoWSmqa+rDEa8ox8YPyTMJM\nRamIiIiIBC6ui1ItdJKaoj4s8ZpyTPygPJMwi8uiVD2lIiIiIrVLXBalun0vNU19WOI15Zj4QXkm\nYaaiVEREREQCp6JUJAbqwxKvKcfED8ozCTMVpSIiIiISuLgsSvWYUalp6sMSrynHxA/KMwmzuCxK\n9ZhRERERkdolLotSbQklNU19WOI15Zj4QXkmYRaXRal6SkVERERqFxWlIjFQH5Z4TTkmflCeSZjF\ndVG6bx84F2wsIiIiInJ45uKgajMzV/57XH89FBbCxImQlBRQYCIiIiJ1iJnhnLOqvDcuZ0pBt/BF\nREREapO4L0q1V6nUBPVhideUY+IH5ZmEWdwWpdoWSkRERKT2iNuiVLfvpSZpbz/xmnJM/KA8kzBT\nUSoiIiIigYv7olQ9pVIT1IclXlOOiR+UZxJmcV+U5ucHG4eIiIiIHF7cFqUlC500Uyo1QX1Y4jXl\nmPhBeSZhFrdFqXpKRURERGoPFaUiMVAflnhNOSZ+UJ5JmKkoFREREZHAqSgViYH6sMRryjHxg/JM\nwszzotTMBpnZMjNbYWa3HeKcx81spZllm1m3MsdTzexlM1tqZl+aWZ9Yr6stoURERERqD0+LUjNL\nACYC5wInAcPMrGO5cwYD7ZxzxwMjgUllXn4MeMs51wnoCiyN9dp6zKjUJPVhideUY+IH5ZmEmdcz\npb2Blc65Nc65AmA6cFG5cy4CngNwzn0CpJpZppk1Afo75/4cfa3QObcz1gvr9r2IiIhI7eF1UdoS\nWFdmvD56rLJzNkSPHQdsMbM/m9lCM3vGzBrGemEVpVKT1IclXlOOiR+UZxJmiUEHUIlEoAdwvXNu\ngZk9CtwO3F3RycOHD6dNmzYApKWl0bFjN2Ag+/YdvF1R8h+jxhprrLHGGmusscbVH2dnZ5OXlwdA\nTk4O1WHOuWp9QKUfbtYXuMc5Nyg6vh1wzrmHypwzCXjfOfdSdLwMGBB9eZ5zrm30+OnAbc65Cyu4\njiv/PYqLYdQoMIOnn478FqmqWbNmlf5HKOIF5Zj4QXkmXjMznHNVqroSajqYcj4F2ptZlpnVB4YC\nr5c753XgCigtYvOcc7nOuVxgnZl1iJ53NvDfWC+ckAD164NzcOBAtb+HiIiIiHjI05lSiGwJRWQV\nfQIw1Tn3oJmNJDJj+kz0nInAIGAPMMI5tzB6vCvwLJAErI6+tqOCa3xvphTg1lth50546CFIS/Pm\n+4mIiIhIRHVmSj3vKXXOvQ2cUO7Y5HLj0Yd47+dAr6peu0GDSFGan1/VTxARERERP3h9+z5QJXuV\nagN9qa6S5m4RryjHxA/KMwmzuC5KtS2UiIiISO2golQkBlqtKl5TjokflGcSZipKRURERCRwdaIo\nVU+pVJf6sMRryjHxg/JMwqxOFKWaKRUREREJt7guSktW36solepSH5Z4TTkmflCeSZjFdVGqmVIR\nERGR2kFFqUgM1IclXlOOiR+UZxJmdaIo1UInERERkXCrE0WpHjMq1aU+LPGackz8oDyTMIvrolSP\nGRURERGpHeK6KFVPqdQU9WGJ15Rj4gflmYSZilIRERERCZw554KOodrMzFX0PfbvhzFjIDkZHn88\ngMBERERE6hAzwzlnVXlvXM+UJieDWWShU3Fx0NGIiIiIyKHEdVFqFilMQSvwpXrUhyVeU46JH5Rn\nEmZxXZSC+kpFREREaoO4L0q1LZTUBO3tJ15TjokflGcSZnFflJbcvtdMqYiIiEh4xX1RWjJTqqJU\nqkN9WOI15Zj4QXkmYRb3Ral6SkVERETCT0WpSAzUhyVeU46JH5RnEmZ1pijVQicRERGR8KozRan2\nKZXqUB+WeE05Jn5QnkmYxX1Rqi2hRERERMIv7otS9ZRKTVAflnhNOSZ+UJ5JmKkoFREREZHAqSgV\niYH6sMRryjHxg/JMwizui1L1lIqIiIiEX9wXpXrMqNQE9WGJ15Rj4gflmYRZ3BelesyoiIiISPjF\nfVGqnlKpCerDEq8px8QPyjMJMxWlIiIiIhI4c84FHUO1mZk71PdwDq67DoqL4cknITHR5+BERERE\n6ggzwzlnVXlv3M+UmulRoyIiIiJhF/dFKegWvlSf+rDEa8ox8YPyTMKsThWl2qtUREREJJzqVFGq\nmVKpKu3tJ15TjokflGcSZnWiKNVepSIiIiLhVieKUs2USnWpD0u8phwTPyjPJMw8L0rNbJCZLTOz\nFWZ22yHOedzMVppZtpl1L3M8x8w+N7NFZja/qjGop1REREQk3DzdtdPMEoCJwNnARuBTM/unc25Z\nmXMGA+2cc8ebWR/gaaBv9OViYKBzbnt14tCWUFJd6sMSrynHxA/KMwkzr2dKewMrnXNrnHMFwHTg\nonLnXAQ8B+Cc+wRINbPM6GtWEzGW9JRqplREREQknGIu+MyspZmdamZnlPzE8LaWwLoy4/XRY5Wd\ns6HMOQ74t5l9ama/jjXW8pKTI7/VUypVpT4s8ZpyTPygPJMwi+n2vZk9BPwc+C9QFD3sgDkexVXi\nNOfcJjP7AZHidKlz7sOKThw+fDht2rQBIC0tjW7dupXepli2bBYbN8L+/ZFxyX+UJa9rrPHhxtnZ\n2aGKR+P4G5cISzwax+c4Ozs7VPFoXPvH2dnZ5OXlAZCTk0N12KGeGf+dk8yWA12cc0fUlWlmfYF7\nnHODouPbAeece6jMOZOA951zL0XHy4ABzrnccp91N7DLOTehguu4yr7Hp5/Cs89Cz57w6yrPt4qI\niIhIZcwM55xV5b0JMZ63Gkiqwud/CrQ3sywzqw8MBV4vd87rwBVQWsTmOedyzSzFzI6KHm8EnAN8\nUYUYtCWUiIiISMjFWpTuBbLNbHJ0+6bHzezxw73JOVcEjAZmAl8C051zS81spJldEz3nLeBrM/sK\nmAxcF317JvChmS0CPgbecM7NPKJvF6UtoaS6Sm5ZiHhFOSZ+UJ5JmMW6JdTrfH+GMybOubeBE8od\nm1xuPLqC930NdKvKNcvTE51EREREwi2mnlKA6O33DtHh8ugWT6FwuJ7SLVvgzjshIwPuv9/HwERE\nRETqkOr0lMa6+n4gMA3IIbJ3aCsz+5Vzbk5VLuo39ZSKiIiIhFusPaUPA+c45wY4584AzgUe8S6s\nmlW2pzTGiWGR71AflnhNOSZ+UJ5JmMValCY555aXDJxzK6jaavxAJCZGfoqLobAw6GhEREREpLxY\n9yn9E5Hn0D8fPfRLoJ5z7koPY4vZ4XpKAcaNg927Yfx4aNzYp8BERERE6hA/9ikdReRpTjdGf/4b\nPVZraFsoERERkfCKqSh1zuU75yY454ZEfx450qc7Ba2kKM2vVVFLWKgPS7ymHBM/KM8kzCpdfW9m\nf3PO/czMlhB51v13OOe6eBZZDSvZq1QzpSIiIiLhU2lPqZkd45zbZGZZFb3unFvjWWRHIJae0iee\ngC++gOuvhy61ppQWERERqT086yl1zm2K/vE659yasj8cfBxoraCnOomIiIiEV6wLnX5UwbHBNRmI\n17SBvlRqoXUFAAAgAElEQVSH+rDEa8ox8YPyTMLscD2lo4jMiLYzs8VlXmoMfORlYDVNRamIiIhI\neB2upzQVSAceAG4v89Iu59w2j2OLWSw9pTNmwBtvwHnnwUUX+RSYiIiISB1SnZ7SSmdKnXM7gB1m\ndlu5l44ys6Occ2urctEgqKdUREREJLxi7Sl9E5gR/f0fYDXwL6+C8oJu30t1qA9LvKYcEz8ozyTM\nKp0pLeGc61x2bGY9qGWr71WUioiIiIRXrDOl3+GcWwj0qeFYPKXHjEp1DBw4MOgQJM4px8QPyjMJ\ns5hmSs1sbJlhAtAD2OhJRB7RY0ZFREREwivWmdLGZX6SifSW1qo17HrMqFSH+rDEa8ox8YPyTMIs\n1p7SewHMrElk6HZ5GpUH1FMqIiIiEl6V7lNaepJZT+DPRGZKAXYAVznnFngYW8xi2ad07164+ebI\njOmjj/oUmIiIiEgd4tk+pWX8CbjOOfdB9IKnR491qcpFg1B2ptQ5sCr9dYmIiIiIF2LtKS0qKUgB\nnHMfAoXehOSNhARISooUpAcOBB2N1DbqwxKvKcfED8ozCbNKZ0qj+5ECzDazycCLgAN+DszyNrSa\n17AhFBREZkuTk4OORkRERERKVNpTambvV/Je55w7q+ZDOnKx9JQC3HUXbN4M//M/kJnpQ2AiIiIi\ndYhnPaXOuTOrFlI4aQN9ERERkXCqtKfUzC6L/h5b0Y8/IdYcbQslVaU+LPGackz8oDyTMDvc6vtG\n0d+NKz2rlijZQF9FqYiIiEi4HO72/WQzqwfsdM494lNMntFMqVSVnhctXlOOiR+UZxJmh90SyjlX\nBAzzIRbPqSgVERERCadY9yn9yMwmmll/M+tR8uNpZB7QQiepKvVhideUY+IH5ZmEWaxPdOoW/f0/\nZY45IBRbQsWqpCjNzw82DhERERH5rkr3KS09yaytc2714Y4FJdZ9St9/H6ZPh4EDYVhcNCSIiIiI\nhEd19imN9fb93ys49nJVLhgk9ZSKiIiIhNPh9intaGaXAqlmNqTMz3CggS8R1iAVpVJV6sMSrynH\nxA/KMwmzw/WUngBcAKQBF5Y5vgv4tVdBeUVFqYiIiEg4xdpT2s85N8+HeKok1p7SnBx44AHIyoI7\n7vA+LhEREZG6xI+e0kvMrImZJZnZf8zs25JHkNYmJTOl27ZptlREREQkTGItSs9xzu0kcis/B2gP\n3OpVUF5p3hxatIBdu2DaNIhhclUEUB+WeE85Jn5QnkmYxVqUJkV/nw+87Jzb4VE8nkpIgGuvjcyY\nLlwI77wTdEQiIiIiArH3lD4IXAzsA3oTWfg0wznXJ4b3DgIeJVIAT3XOPVTBOY8Dg4E9wHDnXHaZ\n1xKABcB659yPD3GNmHpKS3z+OTz1FJjBjTfCiSfG/FYREREROQTPe0qdc7cDpwI9nXMFRIrHi2II\nLAGYCJwLnAQMM7OO5c4ZDLRzzh0PjAQmlfuYMcB/Y4kzVl27wgUXRG7fT5kCW7bU5KeLiIiIyJE6\n3D6lZ0V/DwEGAhdF/zyISJF6OL2Blc65NdFidjrfL2YvAp4DcM59QmRP1MzodY8FzgOejfULxeqC\nC6BLF9i7F55+Gg4cqOkrSDxRH5Z4TTkmflCeSZgdbqb0jOjvC4kscir/+3BaAuvKjNdHj1V2zoYy\n5zxCZEFVjS9JMoMRIyKLn9avh+ef18InERERkaAcbvP8XWY2FviCSGFY0iPgeflmZucDuc65bDMb\nWObaFRo+fDht2rQBIC0tjW7dujFw4EDg4L8MKxqPGgWjR8/itdcgK2sgZ59d+fka191xibDEo7HG\nGmt8pOOSY2GJR+PaP87OziYvLw+AnJwcqqPShU5mdnf0jycAvYB/EikOLwTmO+cq3avUzPoC9zjn\nBkXHtwOu7GInM5sEvO+ceyk6XgYMINJLehlQCDQEGgOvOueuqOA6R7TQqbzPPoNnnomszr/5ZujQ\nocofJSIiIlJnebbQyTl3r3PuXuBYoIdz7hbn3DjgFKB1DJ//KdDezLLMrD4wFHi93DmvA1dAaRGb\n55zLdc7d4Zxr7ZxrG33fexUVpDXhlFPg3HOhuDhSnG7f7sVVpDYr+dehiFeUY+IH5ZmEWaVFaRmZ\nQNmlQAeixyrlnCsCRgMzgS+B6c65pWY20syuiZ7zFvC1mX0FTAauO4L4a8zFF0OnTpGN9SdNgoKC\nIKIQERERqZti3af0TuBnwGvRQxcDLznnHvAwtphV9/Z9iT174A9/gK1b4bTT4PLLIwuiREREROTw\nqnP7PqaiNHqRHkD/6HCOc25RVS7ohZoqSgHWrYOHHorMlP7yl3DGGYd/j4iIiIj4sHk+gHNuoXPu\nsehPaArSmtaqVWSGFGD6dFi9Oth4JBzUhyVeU46JH5RnEmYxF6V1SZ8+cNZZUFQEkyfDzp1BRyQi\nIiIS32K+fR9mNXn7vkRRETzyCKxcCe3awQ03QMOGNXoJERERkbjiy+37uqZePbjmGkhLg1Wr4P77\nYePGoKMSERERiU8qSivRpAmMGwctW8LmzfDAA/Dpp0FHJUFQH5Z4TTkmflCeSZipKD2M5s3h9tsj\nfaYHDsCzz8Lf/ha5vS8iIiIiNUM9pTFyDmbPPliQtmsHI0dCaqqnlxURERGpNXzZpzTM/ChKS6xe\nHVmRn5cXub1/zTVw/PG+XFpEREQk1LTQyUdt28Kdd0KHDpGtoiZMgHffjcykSvxSH5Z4TTkmflCe\nSZipKK2CJk3g5pvhnHOguBhefhmmTIH9+4OOTERERKR20u37alq4EKZNixSkxxwD114LRx8dSCgi\nIiIigVJPaYBFKUBuLjz9NGzaBMnJMHw49OgRWDgiIiIigVBPacAyM+E3v4GePSE/P7IQ6u9/17ZR\n8UR9WOI15Zj4QXkmYaaitIYkJ8PVV8NPfwoJCfDvf0ceU7pjR9CRiYiIiISfbt97YOXKyMKnHTu0\nbZSIiIjUHeopDVlRCpHtoqZMgRUrIjOnl14KZ58NVqX/mURERETCTz2lIVTRtlHPPKNto2or9WGJ\n15Rj4gflmYSZilIPlcyQXnstNGgQ2T7q/vth48agIxMREREJF92+90luLkyaFClI69eHK66AXr2C\njkpERESk5qintBYUpRDZLur//g8++SQyPvNM+MlPIDEx2LhEREREaoJ6SmuJ5GQYMQJ+8QuoVw/e\nfx8efhi2bw86Mjkc9WGJ15Rj4gflmYSZilKfmcGAAXDrrZCeDqtXwx/+EFmlLyIiIlJX6fZ9gHbt\ngqlTYenSyKKoyy+HU08NOioRERGRqlFPaS0tSiGyXdRrr8HMmZHxoEFw8cWV72fqnONA0QEcDudc\n6W+AYlccOafca/Xr1adR/UZefx0RERGpw1SU1uKitMQHH8ALL0SK1B494MorISnpu+es27GOj9d/\nzPwN89mZv/OIr9E6tTVdMrvQJbMLrVNbY9rJP2azZs1i4MCBQYchcUw5Jn5QnonXqlOUat13SPTv\nD82aRbaNWrgQtm2D668HknfyyfpP+Hj9x6zfub70/Pr16pNgkZbgkt9mhmGlxWbJnw1jT8Ee1u5Y\ny9oda5mxYgapDVJLC9SOzTpSv15937+ziIiISAnNlIbMxo3w2MQCVu78nN1p8zi6y39p2DByS75R\n/Ub0atGLfq36kZWadUQznQVFBSzbsozFuYtZnLuYvP15pa8l1UuiY7OOpUVqWoO0Gv9eIiIiEv90\n+z4OilLnHKu3r2be+nnMzVnAZ5/vY+cuSKqXwEX9OjOkVz86Z3YmMaH6k9vOOdbvXM/nuZ+zJHcJ\nOXk533m9VWorumZ25fTWp5PeML3a1xMREZG6QUVpLS5KC4oKeHf1u3y07iO+3fNt6fFjG2eRt6Qf\n27/sRbIdxbBhcMYZ3sSwY/8OlmxewuLcxSz9dikHig4AkRnUs447i0HtB5GSlOLNxWsJ9WGJ15Rj\n4gflmXhNPaW11PZ925m0YFLpTGVagzT6HNuHvsf2pUXjFrgB8M9/wr/+FXkSVG4uXHppZPuompTa\nIJXTW5/O6a1Pp6CogOVbl/PR2o9YuGkh73z1Dh+u/ZDzjj+PgW0G1shMrYiIiEh5mikNyPIty5my\ncAq78neRkZLBLzr/ghN/cGLpoqWy5s6F55+HoiLo2hWuuirydCiv5eTl8Mp/X2HF1sjO/hkpGVzS\n8RJ6tuiplfsiIiLyPbp9X4uKUucc765+l1eWvoJzjhN/cCJX97j6sHuIrlgBTz8Ne/dC69aRlflp\nPqxHcs6xZPMSXl36Kpt2bQIgKy2LIZ2G0LFZR+8DEBERkVpDRWktKUrzC/N57vPnWLBxAQCDjx/M\nj0/4cYWzoxXJzYUnnoBvv408ovSmm+Doo72M+KBiV8y8dfN4ffnrpSv3T25+MkM6DaFlk5b+BBEg\n9WGJ15Rj4gflmXitOkVpDXcnyqFs3rOZBz98kAUbF9AgsQHX9ryWizteHHNBCpCZCbffDu3awfbt\n8Mc/wtq1HgZdRoIlcFrr0/j9Wb/n4o4X0yCxAV9s/oLfz/k907KnsX3fdn8CERERkbikmVIfLM5d\nzNSFU9lfuJ+jjzqaUb1GcfRRVZ/iPHAgssn+l19CgwYwejQcf3wNBhyDXfm7eGvlW8xeM5ui4iKS\n6iVxbrtzOb/D+UdUaIuIiEj80O37kBalzjlmrJjBjBUzAOh+THeGdxtOg8QG1f7swkL405/gs88i\njyO99lo4+eRqf+wR+3bPt/xj2T9KWxLaNW3H1T2upmnDpv4HIyIiIoFSURrConRvwV6mLpzKF5u/\nwMy4uOPFnNvu3BpdtV5cHNkq6sMPI9tEXXUV9OxZYx9/RFZsXcHUhVPJ259HSlIKV3S9gu7HdA8m\nGA+oD0u8phwTPyjPxGvqKQ2Z9TvX84c5f+CLzV/QqH4jxvQZw6D2g2p8G6WEBLjsMjjnnEiB+uyz\n8MEHNXqJmHXI6MBdA+6iS2YX9hbsZdKCSby45EUKigqCCUhERERqFc2U1rAvN3/J0wuepqCogFap\nrRjVcxQZKRmeXtM5ePtt+Mc/IuMhQ+Dccz29ZCWxON7PeZ9X/vsKhcWFHNvkWK7ucTXHND4mmIBE\nRETEN7p9H5KiNHd3Lvd/cD/7C/fT99i+XNblMpLqJfl2/Vmz4MUXI38eNAguvhiC2uN+7Y61TPls\nCpv3bKZ+vfoM6zyMfsf206b7IiIicSzUt+/NbJCZLTOzFWZ22yHOedzMVppZtpl1ix5LNrNPzGyR\nmS0xs7u9jrU69hfu5+kFT7O/cH/pgiY/C1KAgQPhyisjt/XffhteeCEyixqE1qmtufOMO+l7bF8O\nFB1gWvY0pi6K7EBQG82aNSvoECTOKcfED8ozCTNPi1IzSwAmAucCJwHDzKxjuXMGA+2cc8cDI4FJ\nAM65fOBM51x3oBsw2Mx6exlvVTnnmJY9jU27NnH0UUczvNvwwGYE+/SJrMRPTIQ5c2Dq1MjjSYPQ\nILEBI7qPYET3ESQnJvPphk+5b859rMlbE0xAIiIiElqe3r43s77A3c65wdHx7YBzzj1U5pxJwPvO\nuZei46XAQOdcbplzUoA5wCjn3KcVXCfQ2/fvfPUOry59lQaJDbij/x1kHpUZWCwlli+HJ5+E/Hzo\n3BlGjoxsHRWU3N25TFk4hXU71lEvoR5DOg3h7OPO1u18ERGROBLm2/ctgXVlxuujxyo7Z0PJOWaW\nYGaLgG+Af1dUkAZt6bdLeW3ZawBc2f3KUBSkACecAGPHQqNGsGQJPP447A/wznnmUZncfvrtnHXc\nWRQVF/Hyly8zcf5E9hXsCy4oERERCY3EoAOojHOuGOhuZk2Af5jZic65/1Z07vDhw2nTpg0AaWlp\ndOvWrXQvtpIempoed+7dmSkLp7Bh8Qb6HNuHrkd39fR6VRnfcguMHTuLWbOgqGggN94IH38cXDw/\nP/nn7Fqxi3e+eocv+IL//eh/6Z7fnSbJTULx93WocXZ2NjfddFNo4tE4/sYlx8ISj8bxOX700Ud9\n+f+PGtedcXZ2Nnl5eQDk5ORQHX7cvr/HOTcoOo7l9v0yYEDZ2/fR43cBe5xzEyq4ju+37wuKCnjo\no4dYt2MdJzU/idG9R4f28ZrffgsPPwzbt0PbtjBmTOTxpEHauncrT8x/gk27NtEkuQk39LmB1qmt\ngw2qErNmzSr9j1DEC8ox8YPyTLwW2i2hzKwesBw4G9gEzAeGOeeWljnnPOB659z50SL2UedcXzNr\nBhQ453aYWUPgHeBB59xbFVzH16LUOcdfsv/Cx+s/pllKM+7ofweN6jfy7fpVsWVLpDDdti1SmN54\nIzRsGGxMJZvsL9+ynOTEZK7ucTVdMrsEG5SIiIhUWWh7Sp1zRcBoYCbwJTDdObfUzEaa2TXRc94C\nvjazr4DJwHXRtx8DvG9m2cAnwDsVFaRBmL1mNh+v/5ikekmM6jUq9AUpQLNmMG4cZGTA6tXw2GOw\nL+B2zpSkFG7scyN9j+1LfmE+T336FLNyZgUblIiIiARCm+cfoVXbVjF+7niKXTFX97iaXi17+XLd\nmrJ1a2TGdOtWOO64yK38oGdMnXO8ufJN3lj+BgA/avcjLu10aahW5uuWl3hNOSZ+UJ6J10I7Uxpv\n8vbnMWnBJIpdMT9s+8NaV5BCZKa0ZMb066/h0Udh795gYzIzLuhwAcO7DSfBEvj3qn/zzGfPUFBU\nEGxgIiIi4hvNlMaosLiQCfMmsGrbKjpkdOCmvjdRL6Gep9f0UtkZ0zZtIjOmKSlBRwXLtixj0oJJ\n7CvYR9v0tlzX6zoaJzcOOiwRERGJgWZKffDyly+zatsq0hqkcc0p19TqghQiM6W33BLpNc3JCceM\nKUDHZh35f6f9P5o2bMrq7at58MMHyd2de/g3ioiISK2mojQGc9fNZVbOLBITEhnVa1TczNw1bRq5\nld+sGaxZE57CtEXjFvym/2/ISstiy94tPPjhg6zcujLQmEr2ZhPxinJM/KA8kzBTUXoYa/LW8H+L\n/w+AYZ2H0SatTbAB1bCmTQ/OmK5ZA488Anv2BB0VNEluwrh+4+iS2YW9BXt59ONH+XRD6B7oJSIi\nIjVEPaWVKCwu5N5Z97J5z2b6Z/Xnsi6X1fg1wmL79kiP6bffQqtWcPPNkUeUBq3YFfO3L//G+1+/\nD8BPT/opP2z7w4CjEhERkYqop9Qjc9bMYfOezWQelcnQk4cGHY6n0tMjt/KbN4d168IzY5pgCfz8\npJ/z05N+CkR6e99c8Sbx8I8pEREROUhF6SHsLdjLjBUzALi006UkJiQGHJH3yhemYekxNTN+2PaH\n/KrbrzAzXl/+Oq8te83XwlR9WOI15Zj4QXkmYaai9BD+tfJf7Dmwhw4ZHerUoy/T0g4WpmvXRp78\ntH9/0FFFnNrqVK7ucTUJlsA7X73D9C+ma8ZUREQkTqintAJb927ld+//jsLiQu7ofwdZaVk19tm1\nxfbt8Mc/RvYxbdcuso9pcnLQUUUszl3M5AWTKSwupF+rflzR9QoSTP++EhERCZp6SmvYP5b9g8Li\nQnq37F0nC1KI3MofOzbye9UqmDgRDhwIOqqILpldGN17NPXr1Wfeunk8u/BZCosLgw5LREREqkFF\naTk5eTnM3zCfxIRELu54cdDhBKpZs0hhmpoKK1bAU09BQUie/NnpB524qe9NNEhswGcbP2PSgkme\nPpZUfVjiNeWY+EF5JmGmorQM5xx//+/fATi77dlkpGQEHFHwmjePFKaNG8PSpTB5MhSGZFKyXdN2\njO03lkb1G7EkdwkT508kvzA/6LBERESkCtRTWsbn33zOU58+RaP6jbjvrPtISQrBw+BDYuNGGD8+\nsk1Ut25wzTVQLyRPWt24ayOPzHuEnfk7aZvelhv63KD/7URERAKgntIaUFRcxCtLXwHgwg4Xqqgp\np0WLyIb6KSmQnQ1/+hMUFwcdVUSLxi249bRbSW+Yzurtq5kwbwK7D+wOOiwRERE5AipKoz5Y+wG5\nu3Np3qg5/bP6Bx1OKLVqBTfdBA0awIIF8Je/hKcwbd6oObeeeivNGzVn3Y51jJ87nrz9eTX2+erD\nEq8px8QPyjMJMxWlwL6Cfbyx/A0AhnQaUic2yq+qrKyD20N98gk8/zyEpQMkIyWDW069hWMaH8Om\nXZsYP3c8W/duDTosERERiYF6SoHXlr7G21+9Tfum7bnl1Fswq1IrRJ2yciU8/nhkm6gBA2DYMAjL\nX9vuA7t57OPHWLtjLRkpGYzrN06L1kRERHygntJq2LZvG++ufheAn5z4ExWkMTr+eLj+ekhKgtmz\n4eWXwzNjelT9oxjbbyxt09uyde9Wxs8dz5a9W4IOS0RERCpR54vSko3ye7XsxXHpxwUdTq3SsSOM\nGgWJifCf/8Crr4anMG2Y1JAxfcfQNr0t2/Zt4+G5D1erMFUflnhNOSZ+UJ5JmNXponRN3ho+Wf+J\nNsqvhpNOgpEjISEBZs6EGTOCjuigBokNvlOYasZUREQkvOpsT6lzjgnzJrBi6wrOaXcOl554qUfR\n1Q0LF8KUKZHV+BdfDIMHBx3RQfsL9/P4J4+zatsq0humM67fOH7Q6AdBhyUiIhJ31FNaBUs2L2HF\n1hU0qt+IwceHqIKqpXr0gBEjIoud/vEPePfdoCM6qEFiA27scyPtmrZj+77tPDzvYb7d823QYYmI\niEgZdbIoLSouKn2c6PnHn6+N8mtI795wxRWRP7/8MoSpdamiwnTzns0xv199WOI15Zj4QXkmYVYn\ni9IP135I7u5cftDoBwxoMyDocOLKqafCL38Z+fOLL8JHHwUbT1klhWn7pu0jhencIytMRURExDt1\nrqd0f+F+fvveb9mVv4uRPUfS45geHkdXN/3nP/C3v0Vu548YAX36BB3RQfsL9/PEJ0/w1bavSGuQ\nxrhTx9G8UfOgwxIREan11FN6BN7+6m125e+iXdN2dD+6e9DhxK2zz4ZLLolsEfXnP8NnnwUd0UEl\nM6bHZxxP3v48zZiKiIiEQJ0qSnfm79RG+T4aNAguuCBSmD77LHz+edARHZScmMwNvW+IuTBVH5Z4\nTTkmflCeSZjVqaJ03rp5FBQV0DmzM23T2wYdTp1wwQVw7rmRraImT4Yvvgg6ooMqKkxzd+cGHZaI\niEidVGd6Sp1z/Pa937Jl7xZG9x5N58zOPkUnzkVW4//nP5HHko4eHXkaVFjkF+Yzcf5EVmxdQWqD\nVMb1G0fmUZlBhyUiIlLrqKc0Bsu2LGPL3i2kN0znpOYnBR1OnWIGP/0pnHEGFBTAk0/CypVBR3VQ\ncmIyo3uPpkNGB3bs38GEeRPUYyoiIuKzOlOUfrD2AwBOb306CVZnvnZomMEvfhHZMurAAXjiCfj6\n66CjOqhsYVpRj6n6sMRryjHxg/JMwqxOVGc783eyaNMizIzTWp0WdDh1lhlcfjn06gX5+fDYY7Bm\nTdBRHVRSmGpVvoiIiP/qRE/pO1+9w6tLX6Xr0V25rtd1PkYmFSkuhilTYOFCSEmBsWOhVaugozoo\nvzCfJ+Y/wcqtK7WPqYiIyBFQT2klnHOlt+77t+4fcDQCkJAAV10FXbvC3r3wyCOwfn3QUR1UflX+\nhHkT+HbPt0GHJSIiEtfivihdtmUZ3+75VgucQiYxEa65Brp0gT17YMIE2LAh6KgOKilMSx5JOmbS\nGBWm4in1+okflGcSZnFflGqBU3glJsLIkXDyyZHC9JFHYOPGoKM6KDkxmRv63EC7pu3Ylb+Lh+c9\nzJa9W4IOS0REJC7FdU/pzvyd3Pbv23A4Hjj7AdIbpgcQnRxOQQE8/TR8+SU0bgzjxsExxwQd1UH7\nC/fz+CePs2rbKtIbpnPLqbfQLKVZ0GGJiIiEjnpKD2HeunkUu2I6N++sgjTEkpJg1Cjo1Al27Yrc\nyv/mm6CjOqhBYgNu7HMj7Zq2Y/u+7Tw8VzOmIiIiNS1ui9KyC5zOyDoj4GjkcJKS4LrrIk962rkz\nUpjmhuiJnx9/+DE39rmRtult2bZvmwpTqXHq9RM/KM8kzDwvSs1skJktM7MVZnbbIc553MxWmlm2\nmXWLHjvWzN4zsy/NbImZ3Xgk112+dbkWONUy9evD9ddDhw6wY0ekMN0com1CGyQ2YEzfMaWF6YR5\nE9i6d2vQYYmIiMQFT3tKzSwBWAGcDWwEPgWGOueWlTlnMDDaOXe+mfUBHnPO9TWzo4GjnXPZZnYU\n8BlwUdn3lvmM7/WUPvPZM3y28TMuPOFCLuhwgWffUWpefn7kiU8rV0J6OtxyCzQLUQvn/sL9PPbx\nY6zevpqMlAzG9RtHRkpG0GGJiIgELsw9pb2Blc65Nc65AmA6cFG5cy4CngNwzn0CpJpZpnPuG+dc\ndvT4bmAp0DKWi+oJTrVbcjLccAO0awfbt8PDD8PWEE1IlvSYtk1vy9a9Wxk/d7y2ixIREakmr4vS\nlsC6MuP1fL+wLH/OhvLnmFkboBvwSSwX1QKn2i85GW68Edq2hW3bgi9My/dhNUxqyJi+Y2jXtB3b\n9m1j/NzxeiSpVIt6/cQPyjMJs9AvdIreuv87MCY6Y1qp7zzBKUtPcKrNGjSAMWPguOMiBemECZEC\nNSxKZkxLnvw0fu54cneHaHWWiIhILZLo8edvAFqXGR8bPVb+nFYVnWNmiUQK0r865/5Z2YWGDx9O\nmzZt2LJ3C59t+4zjOh3HyRecDBz8l+HAgQM1rmXjBg2gW7dZfPUVbNkykPHjoU+fWaSm+h9Pie/E\nl9iALvu6sGrNKnZk7WD83PH0LexLRkpGKP7+NNZYY43LjkuOhSUejWv/ODs7m7y8PABycnKoDq8X\nOtUDlhNZ6LQJmA8Mc84tLXPOecD10YVOfYFHnXN9o689B2xxzo09zHVKFzpN+WwKCzYu0AKnOLN3\nLzz2GOTkRBY/jR0LzZsHHdVBB4oO8OT8J1m2ZRmNkxsztt9YWjRuEXRYIiIivgrtQifnXBEwGpgJ\nfIBq1TMAACAASURBVAlMd84tNbORZnZN9Jy3gK/N7CtgMjAKwMxOA34JnGVmi8xsoZkNqux6u/J3\nsegbLXCKRykpcPPNBxc/jR/v7wb7Jf86PJT69eozuvdoTvzBiezK38X4ueNZv3O9P8FJXDhcjonU\nBOWZhJmnRSmAc+5t59wJzrnjnXMPRo9Nds49U+ac0c659s65rs65RdFjHznn6jnnujnnujvnejjn\n3q7sWnPXzaWouEgLnOJUSY9pyT6m48fDxo1BR3VQUr0krut1HSc3P5k9B/YwYd4E1u5YG3RYIiIi\ntYKnt+/9YmauuLiY373/Ozbv2cz1va+nS2aXoMMSjxw4AE89BUuXQqNGkRnUVq0O/z6/FBYX8sxn\nz/D5N5+TkpTCmL5jaJPWJuiwREREPBfa2/d+Wr51OZv3bCa9YTonNz856HDk/7d35+FVlnf+x9/f\nk5OVJYFEQti3iAoEiYKIotQFtWNrsa3Vq46jdRuLYzudcRl/M23n1850sdZ1ptROpVVbW1uxtfKr\nVEcRrQrILoJEMZGwExIg+3Lu3x/3SU4IiwFyznOSfF7X9Vzn2U6eL7lunuube42jtDS/8tPEiVBT\n40fll5UFHVVMOBTmljNuYUrBFGqbanngrQfYXLk56LBERESSWo9JSl8v89NAnTviXELWY/5ZcgSp\nqXDbbTB5sh8E9eMfw+Y45n3H2g8rHApzc/HNnDHkjLYVoD7c+2F8gpMeQX39JBFUziSZ9ZjsTQOc\nep9wGG69FYqLob4eHnzQL02aLFJCKdxUfBNTh071ienShyipSKIARUREkkiP6VN6y/O3UJRfxNxp\nc4MORxIsEoH582HZsljT/imnBB1VTMRF+MXqX7C0fClpKWnMnTaXU/KSKEAREZEuoj6lUVrBqXcK\nheCGG2DGDD8I6tFHYf36oKOKCVmI60+/nhnDZ9DY0sgjSx9hzY41QYclIiKSVHpMUqoBTr1bKATX\nXQczZ0JTkx+dv3Zt1/38E+2HFbIQ102+jlmjZtEcaWbeO/NYWr60a4KTHkF9/SQRVM4kmfWYpFQD\nnMQMvvxl+NSnoLkZ5s2DlSuDjirGzLh64tVcVngZERdh/ur5vFb6WtBhiYiIJIUe06d0b+1eTZgv\nADgHzz4LL73kE9W//Vs4J8nGvy36YBELNiwAYM6pc7h03FEXKxMREekWTqRPaY9JSnvCv0O6jnOw\ncCH86U/++MorYfZsn6QmiyVlS/j1ul/jnOOScZcw55Q5WDIFKCIicow00EmkAzO4/HK45hq/v2CB\nrz093r9d4tEP67yR53HjlBsJWYhFHyxqS1Cld1JfP0kElTNJZkpKpUebNQtuvBFSUnxz/i9/CS0t\nQUcVM3XoVG6behvhUJglZUt4fNXjtESSKEAREZEEUfO99Arr1/uBT42NUFQEt9ziV4VKFpsqNvHo\nskdpaG6gKL+IW864hdSUJApQRESkE9SnVEmpdMLmzX4O05oaKCyEr34VsrKCjiqmtKqUh5c+TE1j\nDePzxvPVqV8lI5wRdFgiIiKdpj6lIp0wZgzceSfk5PjlSO+/H/bv79x3E9EPa1TOKP55xj+TnZHN\n+3ve54G3HqCmsSbuz5XkoL5+kggqZ5LMlJRKr1JQAHfdBfn5UF4OP/gB7N4ddFQxQ/oN4c4Zd5KX\nlUdpVSk/evNH7KvfF3RYIiIicafme+mVDhyARx6BsjLo3x++9jUYNizoqGKq6qt48O0H2X5gO3lZ\nedxx1h3k980POiwREZGjUp9SJaVyHOrr4Sc/gY0bITMT5s71fU2TRXVjNY8sfYTSqlL6pPVh7tS5\njB04NuiwREREjkh9SkWOQ0YG/MM/QHEx1NXBQw/B2rWHvzeIflh90/ryjbO/QVF+ETWNNfz4rR+z\ncnsSrZsqXUp9/SQRVM4kmSkplV4tHIabb4aZM6GpCf77vyGZ3tnp4XRum3ob5408j+ZIM4+teIyX\nN78cdFgiIiJdTs33IviVnl54wW8AF1wAX/wihJLkzzbnHIs+XMRzG54D4MIxF/KF075AyJIkQBER\nEdSnVEmpdJm334YnnvCrPk2aBDfd5Jv5k8Wyrcv4xepf0BJpYUrBFG6ccqMm2RcRkaShPqUiXWT6\ndPjHf4Q+fWDdOrjvPqisTJ5+WNOGTuPr079OZmomq7av4oG3H6C6sTrosKQLJEsZk55N5UySmZJS\nkQ4KC+Gee2Jzmf7nf8KOHUFHFXNy7sncdc5dDMgcwId7P+QHb/yA3TVJNNmqiIjIcVDzvcgR1NTA\nvHmwaROkpsKNN8KUKUFHFVNVX8Wjyx5ly74t9Evvx+3TbmdUzqigwxIRkV5MfUqVlEqcNDfDr34F\nb74JZnDllXDxxX4/GdQ31/PTd37Ke7vfIzUllZuLb2by4MlBhyUiIr2U+pSKxEk4DNddB2PHLsY5\nePZZeOopPxAqGWSEM7h92u3MGD6DppYmfvLOT1hcujjosOQ4qK+fJILKmSQzJaUin8AMpk2DW2/1\nzfhvvAEPPwy1tUFH5qWEUrhu8nV8ZvxncM7x9LqneXrd0zRHmoMOTUREpNPUfC9yDEpL4b/+C/bv\nh8GD4fbb4aSTgo4q5q0tb/HU2qdojjRTmFvIrWfcSr/0fkGHJSIivYT6lCoplQSqqPCJ6dat0Lcv\n3HILjB8fdFQxH1V+xLx35lFVX8WAzAHcduZtjMwZGXRYIiLSC6hPqUicte+HlZsLd90FEydCdTU8\n8AAsWuRXhUoGoweM5t6Z9zJmwBgq6yq57837WFq+NOiw5BOor58kgsqZJDMlpSLHISMD5s6FT3/a\nJ6MLFsBPfwp1dUFH5mVnZPNPM/6Jc0ecS1NLE4+vepzfv/d7Ii4SdGgiIiKHpeZ7kRO0di08/rhP\nSAcNgttugyFDgo7Kc87xWtlr/Pbd3xJxEU496VRuLr6ZPml9gg5NRER6IPUpVVIqAdu1y9eUlpdD\nWhpcey2cdVbQUcWUVJTw0xU/5UDDAfKy8pg7bS5D+iVJ5iwiIj2G+pSKxNkn9cMaNAjuvhumT4fG\nRl9z+pvf+Mn3k0FhbiH3zryXEdkj2FO7h++/8X1WbV8VdFjSjvr6SSKonEkyU1Iq0kXS0uD66+HL\nX4aUFHj1Vbj/fqisDDoyb2DmQO48506mDZ1GQ3MD896Zx/PvP49aGUREJBmo+V4kDkpLYd48n5D2\n6wc335w800Y553h588s8u+FZnHMU5RfxlSlfITM1M+jQRESkm1OfUiWlkoSqq+F//gc2bPCrQs2Z\nA7Nn+/1k8N7u9/jZip9R21RLXlYeNxXfxOgBo4MOS0REujH1KRWJs+Pph9W3L9xxx8HTRs2blzzT\nRp120mkH9TP94V9/yKIPFqk5PyDq6yeJoHImyUxJqUgchUJwxRV+TtPMTFi9Gr7zHSgpCToy76Q+\nJ3H3uXdz0ZiLiLgICzYs4KGlD7G/YX/QoYmISC+j5nuRBNm9G372Mygr8034s2fDZz8L4XDQkXnv\n7nqX+avmU91YTb/0ftxw+g1MGDQh6LBERKQbSeo+pWZ2KfAgvlb25865HxzmnoeBy4Aa4Abn3Kro\n+Z8DlwM7nXNFR3mGklLpFlpa4IUX4M9/9k36w4fDV76SPJPtV9VXMX/VfDbu2QjA7LGzueKUKwiH\nkiRzFhGRpJa0fUrNLAQ8ClwCTACuMbNTOtxzGTDWOVcI3Ar8pN3l+dHvigSqq/phpaT45vw774S8\nPNiyBf7jP+B//9cnqUHLycjha9O/xudO+RwhC/GXD//CD//6Q3bV7Ao6tB5Pff0kEVTOJJnFu0/p\nNKDEOVfmnGsCfgNc0eGeK4AnAJxzS4FsM8uPHr8BJMksjyJdZ+xY+Ld/g3PO8RPsP/MMPPQQVFUF\nHRmELMRlhZdx5zl3kpuVS1lVGd9d8l2Wli8NOjQREenB4tp8b2afBy5xzt0SPb4WmOacu6PdPX8C\nvuecezN6/DJwl3NuZfR4JPAnNd9LT7V6NTz5pJ9CKivLL1F6xhlBR+XVNtXy1NqnWLFtBQBnDz+b\nayZeQ3o4PeDIREQkGZ1I832P6Sh2/fXXM2rUKABycnI4/fTTmTVrFhBrrtCxjpPxuKpqMRdcAJs3\nz+Ldd+Hb317MaafBt789i8zMYOPLSs2icH8h9TX1bOq/ibe2vMUrr77Cpws/zdWXX50Uvz8d61jH\nOtZxcMerV6+mKtrMV1payomId03pdODbzrlLo8f3AK79YCczmwe86pz7bfR4I3C+c25n9Fg1pRK4\nxYsXt/0njBfnYMkS+N3voKkJcnPhhhugsDCuj+207Qe289iKx9h2YBshCzF77GwuP/lyUlNSgw6t\nR0hEGRNROZN4S9qBTsByYJyZjTSzNOBq4PkO9zwPXAdtSWxVa0IaZdFNpEczg/PP931NR46Eigq4\n/374/e+hoSHo6KCgXwH3zryXC8dciMPx4gcv8p0l3+GDvR8EHZqIiPQAiZoS6iFiU0J938xuxdeY\nPha951HgUmJTQrX2J/01MAvIBXYC33LOzT/MM1RTKj1Kx6mjcnPh6quh6IjtBYm1uXIzT6x5gu0H\ntgMwa9Qs5pw6h4xwRsCRiYhIkJJ6ntJEUFIqPVVpKTz1lJ86CqC4GL70JcjJCTQsAJojzSzctJAX\nP3iRiIswMHMg1xZdqwn3RUR6MSWlSkolzoLshxWJwCuvwPPP+2b8jAw/1+msWX4Z06CV7y/niTVP\nUFZVBsD0YdO5asJV9EnrE3Bk3Yv6+kkiqJxJvCVzn1IROUGhEFx0Efz7v8Ppp0N9Pfz2t/D97/sl\nS4M2rP8w7jn3Hq489UpSU1J5u/xtvrX4W6zcvjLo0EREpBtRTalIN7NmDTz9NFRW+sFRF1wAn/2s\nr0EN2s7qnTy59klKKkoAmFIwhWsmXkN2RnbAkYmISCKo+V5JqfQyDQ3wpz/55UkjEd/H9OqrfU2q\nBTxXhXOOJWVLeHbDszQ0N5CVmsUXTvsCM4bPwIIOTkRE4kpJqZJSibNk7Ye1ZYsfCNU6X3FRkU9O\nc3MDDQuAvXV7+dXaX/HurncBGJkzkqsmXMW4geMCjiw5JWsZk55F5UziTX1KRXqp4cPh7rvhmmt8\n8/3atfCtb8Ef/gB1dcHGNjBzILdPu50bi28kJyOHsqoy7vvrffxsxc+oqK0INjgREUk6qikV6SH2\n7fOrQS1f7o/79oW/+Rs47zwIB7ygcENzA3/58C8s+nARTS1NhENhLh57MZeNu4z0cHqwwYmISJdR\n872SUpE2mzfDggVQ4scakZcHc+bAGWcE39+0sq6SBRsWsGzrMgD6p/dnzqlzOHvY2epvKiLSAygp\nVVIqcdbd+mE555vyn3sOtvtFlxg1Cj7/eTj55EBDA/yKUM+sf4aPKj8CYET2CK6acBWFuYUBRxac\n7lbGpHtSOZN4U59SETmIGUyeDN/8Jlx7LWRn+8FQ998PjzwC27YFG9+YAWO4+5y7+cqUr5CTkcPH\n+z7mR2/+iMdWPMae2j3BBiciIoFQTalIL9DQ4KePWrTIT75vBjNmwGc+AwMGBBxbcwMvbX6JFz94\nsa2/6UVjLuKScZeQlZoVbHAiInJM1HyvpFSkUw4cgIUL4bXX/Pymqalw4YUwezb0CXhV0Mq6Sp7b\n+BxLy5cCkBHO4ILRF3DRmIu0ZKmISDehpFRJqcRZT+uHtWuXnzZqxQp/nJ4OM2fCxRf7ifiDtLly\nM3/c+Ec27tnoYwuntyWnfdP6BhtcHPW0MibJSeVM4k1JqZJSibOe+iIvLYXnn4f16/1xSgqcfTZc\ncgkMGhRoaHy490Ne2PQC7+1+D/DJ6axRs7h4zMX0S+8XbHBx0FPLmCQXlTOJNyWlSkpFTsjHH8OL\nL8LKlX7kvhkUF8Nll/kJ+oO0uXIzCzctbFsZKi0ljfNHnc/ssbPpn94/2OBEROQgSkqVlIp0iV27\n/GCot96ClhZ/bsIEuPRSKCwMdp7T0qpSFm5ayNqdawFITUnlvJHnccnYS8jOyA4uMBERaaOkVEmp\nxFlva/KqqoKXXoLXX/cj9wHGjPHJaVFRsMlpWVUZC0sWsmbHGgDCoTAzR87k4jEXk5uVG1xgJ6i3\nlTEJhsqZxJuSUiWlEme99UVeUwOvvgqvvOL3AYYMgYsugjPP9AOkgrJl3xYWlixk1fZVgH8RTs6f\nzKdGf4rxueO73QpRvbWMSWKpnEm8KSlVUioSVw0N8MYbvva0stKfy8yE6dPhvPN8ohqU8v3lLPpg\nESu2r6Al4vscFPQrYNaoWUwfNp2McEZwwYmI9DJKSpWUiiREczO8846f53Tz5tj5sWN9cnrGGX7u\n0yDsb9jP62Wvs6RsCVX1VYCf6/Ts4WfzqVGfIr9vfjCBiYj0IkpKlZRKnKnJ61Dl5b7P6dtv+1Wi\nALKy/EpRM2fC4MHBxNUSaWHVjlUsLl1MSUVJ2/lTTzqVC0ZfwMRBEwlZ8q2wrDImiaByJvF2Iklp\nuKuDEZHeYdgwuOYauPJKWL4cliyBsjJ4+WW/nXyyrz2dMgXCCXzTpIRSOHPImZw55EzK95ezuHQx\nb5e/zYbdG9iwewO5WbnMGjWLc4afo5WiRESSiGpKRaTLlJX52tNly2Kj9vv29RPyT50KI0YEM3K/\ntqmWN7e8yasfvcqe2j2AT14nDprItKHTmJw/mdSUgPodiIj0IGq+V1IqklTq62HpUl97Wl4eO3/S\nSX7U/tSpfnBUohPUiIuwftd6FpcuZv3u9bS+NzLCGRQXFDNt6DTG541PyuZ9EZHuQEmpklKJM/XD\nOj7OwUcf+ZrTFStg//7YtYKCWIKaH8AYpP0N+1m+dTnLti6jtKq07Xx2RjZTh0zlrGFnMbz/8IRN\nLaUyJomgcibxpqRUSanEmV7kJy4SgZIS3/905crYvKfg+6dOneqT1Ly8xMe2s3ony7YuY+nWpeyu\n2d12vqBfAdOGTmPa0GnkZcU3MJUxSQSVM4k3JaVKSkW6lZYW2LjRTy+1ahXU1cWujRrlk9OiIhg0\nKLFN/M45Pqr6iGVbl7F863KqG6vbro0ZMIbJgydTlF9EQd+Cbjc5v4hIIigpVVIq0m01N8P69T5B\nXbMmNkAKfK3pxIkwYQKMH5/YFaRaIi1s2LOBpeVLWbVjFU0tTW3XcrNyKcovYtKgSYzPG084pIlM\nRERASamSUok7NXklRlMTrFvna0/Xrz+4iT8c9tNMTZzot0TWojY0N7B+93rW7VzHul3rONBwoO1a\nejid0046jUmDJjEpfxL90/sf1zNUxiQRVM4k3jRPqYj0CKmpUFzst0jETzH17rt+KyuD997z2zPP\nJLYWNT2cTnFBMcUFxW1N/Ot2rmPtzrWU7y9n1fZVrNq+CoBROaOYlD+JSYMmMSJ7hJr5RUQ6STWl\nItItHDjgE9J33z18Lero0TBuHBQWwpgxkJmZmLgq6ypZt2sda3asYeOejTRHmtuuZaVmUZhbyMm5\nJ1M4sJDh2cM13ZSI9GhqvldSKtKrHK4Wtf0rwMyP6B83Lrbl5MQ/robmBt6veJ+1O9fy7q53qayr\nPOh6RjiDcQPH+SQ1t5CR2SNJCaXEPzARkQRRUqqkVOJM/bCSW00NfPBBbCsr8yP828vLi9Wkjhvn\n50aNZ8u6c46KugpKKkrYVLGJTRWb2laTapWWksbYgWM5Ofdk9r63ly9d/iWtLCVxpXeZxJuSUiWl\nEmd6kXcvjY1QWhpLUj/80K8y1V5WFgwf7pc+bf3Mz4dQHFvXK+sqKdkbS1J3Vu9su7Zt3TaGFQ1j\nSL8hjMgewcickYzMHsmw/sOUqEqX0btM4k1JqZJSETmKSAS2bvUJakmJ/9y379D70tJ8s/+IEbFk\ndcgQ32c1HvY37G+rSS3ZW8K2A9vo+C4LWUiJqoh0G0pKlZSKyDFwzielH38MW7bEPisqDr03JcUn\npsOG+aVRBw/2W16ev9aVGpobKN9fzsf7PqZsXxllVWVsr95+xER1ePZwCvoWMLjvYAr6FZCXlaeB\nVCISKCWlSkolztTk1TvU1MSS1NZEdefOgwdRtQqF/FyprUlq++14Rv4fqYx1NlEFCIfCDOoziIJ+\nBRT0LaCgn09Y8/vkq2ZVAL3LJP40T6mISBfo0wdOOcVvrRoaoLzcN//v2BHbKipi+x317x+rTc3N\njX3m5vpZAI6l32p6OJ2xA8cyduDYWEzRRHXrga3sqN7B9gPb2V69ncq6SrYd2Ma2A9sO+hlmRl5W\nHoP7DiYvK4+8rDxyM3Pb9jNTEzR/lojIUaimVETkODQ2wq5dByeqrVtT05G/FwrBwIGxJLV90jpg\nAGRn+0UEjkd9cz07qncclKjuqN7B7prdRFzkiN/LSs0iNyv3sAlrTkYOGeEMLQIgIp2i5nslpSKS\nJJyDvXt9s39Fhd/27IntH26AVUeZmb5GNTv70K39+c6uYtUcaWZXzS52Vu9kT+0e9tTuoaKuom2/\nqeUoWTR+6qqcjByyM7LJTs8+ZL/1OD0lXcmrSC+npFRJqcSZ+mFJV2lq8klra5LamrQuW7aY3NxZ\n7Nt36ByrR5KeDn37xrY+fQ4+Pty5jjMJOOeobqxuS1DbJ6wVtRVU1VfR2NLYuXjC6fRL60fftL6d\n2vqk9dHArATTu0ziLan7lJrZpcCDQAj4uXPuB4e552HgMqAGuN45t7qz3xVJhNWrV+tFLl0iNdXP\nh5qff/D56urVfP3rs3DOD7jat89vVVWx/fZbVZXv79rQcPhZA472/MzM9puRmdkvuo0mMxOGZMHY\nTMgcCOnpDsINNFgV9eyjzlVRF9lHdfM+9tVXsa9hH/vq91FZX0lDcwMNzQ2HLBJwJGZGZjiTzNTM\nts+s1KxDzmWGo+ej++nhdNJT0ts+w6Gwamg7Se8ySWZxTUrNLAQ8ClwIbAOWm9kfnXMb291zGTDW\nOVdoZmcB84DpnfmuSKJUVVUFHYL0cK1lzCxWqzl06JHvd84npNXVh241NUc+39Tkt/37OxuZARnA\n4OhGW5zp6ZCR4ed3zU9zhDPqiaQeIBKupiWlmuZQ7LPJ/NaI3xpcNU3UYqFaQlaLhSBktH2GQv4Z\nFvL7rdfMfEQHRWhGeko6GeEM0lLSDkla01LSSE1JJTWU2vbZ8VxaShrhUNifD6WSEkohHAoftKXY\nwee6YyKsd5kks3jXlE4DSpxzZQBm9hvgCqB9YnkF8ASAc26pmWWbWT4wuhPfFRHplcx8QpiR4QdK\ndYZz0NwMtbVQVxf77Li1P99aG9txa272q2TFVsoyIDO6DTpiDCF8ipsBOCI0Wx0tVkez1dJMHc3R\nfX+uru2z/X4k1IALNfhPa8CFmjGrx6z+4OT1cBvR69HfIYe5p/X3a9F/1kH3EvtOiBAp0QQ1ZCFS\nLIWQpbTttx6nhFJIsRChUAop+OOQhQ67GSFSQh3Oh0KECGEYZubvMyMl+um/Z4RC0c/oeYNo9wgj\nFDIMY01ZGU8teSN63Q7+NNr2/T8z+mlGKHpf6y/H38dB94Us9t1YmejwMzt8J3Zn7FkHnT9C4h9q\nd77jz+rM9w/3czovvn+MJOMfOzl9Mpkw6sj/r7tKvJPSocCWdsfl+ET1k+4Z2snviiREaWlp0CFI\nD5eIMmbmm+9bB0qdiJaWwyerjY2+Jrb1s3U7/HGIpqY+tLT0oanJJ7otLRyy39Lij5ua/OpcHTki\ntNBAi/ktYg0HH9NExBqjn01EaKLlcMfR/Yg14WghQjPOmv0nLbHz0XMQARqjW/ewbflf2XhSWtBh\nSDczaVARz/zz3Lg/JxnnKT2+EVtJ+JeF9Cy//OUvgw5BejiVMUmE/etLgg5BupmNgN15e9yfE++k\ndCswot3xsOi5jvcMP8w9aZ34LsBxj/ISERERkeQQ77k4lgPjzGykmaUBVwPPd7jneeA6ADObDlQ5\n53Z28rsiIiIi0gPEtabUOddiZrcDfyE2rdMGM7vVX3aPOef+n5l92sw+wE8JdcPRvhvPeEVEREQk\nGD1i8nwRERER6d66zVIaZnapmW00s01mdvcR7nnYzErMbLWZnZ7oGKX7+6RyZmbnm1mVma2Mbv8a\nRJzSfZnZz81sp5mtPco9epfJcfukMqb3mJwoMxtmZq+Y2XozW2dmdxzhvmN6l3WLpLTdRPqXABOA\na8zslA73tE3CD9yKn4RfpNM6U86iljjniqPbdxMapPQE8/Fl7LD0LpMucNQyFqX3mJyIZuAbzrkJ\nwNnA3K7Iy7pFUkq7Sfidc01A60T67R00CT/QOgm/SGd1ppxBvGdOlh7NOfcGUHmUW/QukxPSiTIG\neo/JCXDO7WhdEt45Vw1swM8v394xv8u6S1J6pAn2j3bP1sPcI3I0nSlnAGdHmyIWmtlpiQlNehG9\nyyQR9B6TLmFmo4DTgaUdLh3zuywZJ88XSWYrgBHOudpo08QfgJMDjklE5FjoPSZdwsz6Ar8Hvhat\nMT0h3aWm9EQm4RfprE8sZ865audcbXT/z0CqmQ1MXIjSC+hdJnGl95h0BTML4xPSJ51zfzzMLcf8\nLusuSemJTMIv0lmfWM7a94cxs2n4adX2JjZM6QGMI/fp07tMusIRy5jeY9JFHgfec849dITrx/wu\n6xbN9ycyCb9IZ3WmnAFfMLPbgCagDvhScBFLd2RmvwZmAblm9jHwLfyyynqXSZf4pDKG3mNygszs\nHODLwDozWwU44F5gJCfwLtPk+SIiIiISuO7SfC8iIiIiPZiSUhEREREJnJJSEREREQmcklIRERER\nCZySUhEREREJnJJSEREREQmcklIRkTgxs8nRZRzj+YzHzOyU6P6/dLj2RjyfLSLSlTRPqYhInJjZ\n3wFnOuf+IUHPO+Cc65eIZ4mIdDXVlIpIr2dmz5nZcjNbZ2Y3Rc+FzGy+ma01szVm9rXo+TvMbL2Z\nrY6unIOZZZnZz83sbTNbYWafMbNU4P8CV5nZSjP7Yodn/p2Z/cHMXjWz983sm+2ufSMay9p2jxIB\nHQAAAvBJREFUz80ysxfMbFX0/Bej5181s2Iz+x6QGX3Wk9FrB9r9zPuiP3ONmV0VPXd+9Pu/M7MN\nrd8TEQlCt1hmVEQkzm5wzlWZWQaw3MyeBUYDQ51zRQBm1j96793AKOdcU7tz/wf4X+fcjWaWDSwD\nXga+CZzhnLvjCM+dCkwA6qPPfSF6/u+i11KApWa2GBgLbHXOXR6N56AaUefcv5jZXOdccfvT0Xs/\nDxQ55yaZ2aDos16L3nM6cBqwA/irmc1wzr3Z6d+ciEgXUU2piAh83cxWA28Dw4BCYDMw2sweMrNL\ngNZaxzXAr83sy0BL9Nxs4J7oGtCL8euMj+jEc19yzlU55+qBZ4GZwLnAc865eudcDbAgen4dcLGZ\nfc/MznXOHTjyjz3EOcDTAM65XdEYp0avLXPObXe+L9dqYNQx/FwRkS6jpFREejUzOx+4ADjLOXc6\nPjHLcM5VAZPxCdytwP9Ev/I3wKNAMb7GMQUw4PPOuSnRbbRz7v1OPL59p34DItFzdsiNzpVEn7kO\n+K6Z/evh/jmdeGbH+xra7begFjQRCYiSUhHp7bKBSudcQ3QU+3QAM8sFUpxzzwH/BkyJ3j/COfca\ncA/QH+gDLALamujN7PTo7oHoPUdysZnlmFkm8Dngr8AbwBVmlmFmfYA5wOtmVgDUOed+DdyHT1A7\najSz9klla/L5OvClaD/Zk/A1r8s+8TcjIpJA+otYRHq7F4G/N7P1wPvAW9HzQ4H5ZhbC117eE034\nnor2JTXgIefcfjP7DvCgma2Nnv8I+CzwavR7K4HvOed+1+HZy/DN80OBJ51zKwHM7BfA8uhzH3PO\nrTGz2cB9ZhYBGoG/j/6M9rWtjwFrzWyFc+5vW685554zs+n4rgcR4E7n3C4zO7VDPJqORUQCoymh\nREQCEJ0u6miDoEREehU134uIiIhI4FRTKiIiIiKBU02piIiIiAROSamIiIiIBE5JqYiIiIgETkmp\niIiIiAROSamIiIiIBO7/A9lU7fMtEdzDAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fe87ff8bcc0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Set parameters for bisection method\n",
    "crit = 1e-6\n",
    "r_min = 0.02\n",
    "r_max = 0.05\n",
    "r = 0.03\n",
    "\n",
    "# Bisection loop\n",
    "for i in range(100):\n",
    "    am.set_prices(r,r_to_w(am, r))\n",
    "    r_new = rd(am, am.compute_stationary_distribution())\n",
    "    if np.absolute(r_new-r)<crit:\n",
    "        break\n",
    "    elif r_new > r:\n",
    "        r_min = r\n",
    "        r = (r_max+r_min)/2.\n",
    "    else:\n",
    "        r_max = r\n",
    "        r = (r_max+r_min)/2.\n",
    "\n",
    "# Plot stationary distribution at the equilibrium\n",
    "fig, ax = plt.subplots(figsize=(11, 8))\n",
    "n=50                      # Determine the max asset level to show in the plot\n",
    "ax.plot(am.a_vals[0:n], am.g[0,0:n], lw=2, alpha=0.6, label='low income')\n",
    "ax.plot(am.a_vals[0:n], am.g[1,0:n], lw=2, alpha=0.6, label='high income')\n",
    "ax.grid()\n",
    "ax.set_xlabel('asset position')\n",
    "ax.set_ylabel('distribution')\n",
    "ax.legend(loc='upper right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.04605979919433595"
      ]
     },
     "execution_count": 48,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "am.r"
   ]
  }
 ],
 "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.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}