piti118/fix numerical instability.ipynb

## fix numerical instability.ipynb
{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# function t = taufunW(m,a) \n",
    "# % This function calculates T_m(a) as defined in the \n",
    "# % last line of page 2.\n",
    "# t = pi^(m*(m-1)/4);\n",
    "#   for i = 1:m\n",
    "#       t = t*gamma(a - (i-1)/2);\n",
    "#   end\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 265,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "from matplotlib import pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from math import pi\n",
    "from scipy.special import gamma, gammaln, gammainc\n",
    "from numpy import log, exp, sqrt\n",
    "import math\n",
    "\n",
    "def org_t(m, a):\n",
    "    t = pi**(m*(m-1)/4.0)\n",
    "    for i in range(1, m+1): #matlab is inclusive\n",
    "        t = t*gamma(a - (i-1)/2.0)\n",
    "    return t\n",
    "\n",
    "def log_t(m, a):\n",
    "    lt = log(pi)*(m*(m-1)/4.0)\n",
    "    #need to take care of the sign?\n",
    "    return lt + sum(gammaln(a - (i-1)/2.0) for i in range(1, m+1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1.0, 362880.0, 76727486417.585541, 9718995188766424.0, 7.5952181508586655e+20, 3.7780746438430019e+25, 1.2367389208381925e+30, 2.7609601795779169e+34, 4.3682164847183517e+38, 5.1060479665939196e+42]\n",
      "[  1.00000000e+00   3.62880000e+05   7.67274864e+10   9.71899519e+15\n",
      "   7.59521815e+20   3.77807464e+25   1.23673892e+30   2.76096018e+34\n",
      "   4.36821648e+38   5.10604797e+42]\n"
     ]
    }
   ],
   "source": [
    "a = 10\n",
    "ms = np.arange(10)\n",
    "org_ts = [org_t(m,a) for m in ms]\n",
    "log_ts = [log_t(m,a) for m in ms]\n",
    "print org_ts\n",
    "print np.exp(log_ts)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 301,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3.54488054185e+22\n",
      "3.54488054185e+22\n"
     ]
    }
   ],
   "source": [
    "def org_k(nmin, nmax):\n",
    "    nmat = nmin if nmin%2==0 else nmin + 1\n",
    "    alp = (nmax - nmin - 1)/2.0\n",
    "    K = (pi**(nmin**2/2.0)/2**(nmin*nmax/2.0)) #that's gonna overflow real quick\n",
    "    K = K/org_t(nmin, nmax/2.0) #I can see a lot of gamma cancellation here too\n",
    "    K = K/org_t(nmin, nmin/2.0)\n",
    "    Kd = K*2**(alp*nmat+nmat*(nmat+1)/2.0)\n",
    "    for k in range(1, nmat+1):\n",
    "        Kd = Kd*gamma(alp+k)\n",
    "    return Kd\n",
    "\n",
    "def log_k(nmin, nmax):\n",
    "    nmat = nmin if nmin%2==0 else nmin + 1\n",
    "    alp = (nmax - nmin - 1)/2.0\n",
    "    logK = nmin**2/2.0 * log(pi) - nmin*nmax/2.0*log(2.0)\n",
    "    logK = logK - log_t(nmin, nmax/2.0)\n",
    "    logK = logK - log_t(nmin, nmin/2.0)\n",
    "    logKd = logK + (alp*nmat+nmat*(nmat+1)/2.0)*log(2.0)\n",
    "    for k in range(1, nmat+1):\n",
    "        logKd = logKd + gammaln(alp+k)\n",
    "    return logKd\n",
    "\n",
    "print org_k(5,25)\n",
    "print exp(log_k(5,25))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 268,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# function Kd = KfunW(nmin,nmax)\n",
    "# % This function calculates K' as defined in Theorem 1 (page 3)\n",
    "# % using K in equation 2 (page 3).\n",
    "# nmat = nmin;\n",
    "# alp = (nmax - nmin - 1)/2;\n",
    "# K = (pi^((nmin^2)/2))/2^(nmin*nmax/2);\n",
    "# K = K/taufunW(nmin, nmax/2);\n",
    "# K = K/taufunW(nmin, nmin/2);\n",
    "# Kd = K*2^(alp*nmat+nmat*(nmat+1)/2);\n",
    "# for k = 1:nmat\n",
    "#     Kd = Kd*gamma(alp+k);\n",
    "# end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 269,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])"
      ]
     },
     "execution_count": 269,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.zeros((10,10))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 316,
   "metadata": {
    "collapsed": false,
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEnCAYAAACnsIi5AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXdYVNfWxt9Dr0NHEQV7wC7RqFETjImxxMQWb6IxMSb3\npt6YcuO96Zqr6V/ijZpibMQeRY2KBQVRUIoiKCBIEelShhmYGdqU9f1xwCChz96jxPN7nvMoM2fW\n3rM5rLPPu9deSyAiSEhISEj8dTG73R2QkJCQkOCL5OglJCQk/uJIjl5CQkLiL47k6CUkJCT+4kiO\nXkJCQuIvjuToJSQkJP7itOnoBUHYKAhCsSAIl1s553tBEDIEQUgUBGEE2y5KSEhISBhDe2b0mwE8\n2tKbgiBMA9CPiAYAeAnAT4z6JiEhISHBgDYdPRFFAVC0csoTAH6tPzcWgJMgCN3YdE9CQkJCwlhY\naPTeAPIa/VxQ/5qEhISExB2AtBgrISEh8RfHgoGNAgC9Gv3cs/61PyEIgpRYR0JCQqITEJHQ2c+2\n19EL9UdzHATwGoDdgiCMBaAkouKWDElJ1ESWL1+O5cuXM7d74wZw7hwQGwtcuQKkpgL5+YC3N9Cr\nF9Czp/h/Dw/A3V08nJ0BmUw87O0BOzvAxgYwN+98PwxkQGppKuKL4nGl9ApSy1KRVpaGHGUOZNYy\n+Dj5oIdjD3Sz74a0vWmY++pcuNi4wNnGGU42TpBZy+Bg5QAHKwfYWdrB1sIWVuZWEIR2XOulpcDF\ni+KXv3oVSE8HsrKAoiLAxUUchO7dAU9PcSDc3MRBcHYGnJwAR0fAweGPwbC1FQ8rK6CZ9it0OqRq\nNEirqkJaVRUyqquRU1ODnNpaqHQ6eFlbo7uVFbpZWsLd0hJulpZwtbSEs4UFnMzN4WhhAUdzcziY\nm2PT55/j3x9/DFszM9jUHxZmfzx415XVoSqtCtXp1ahKr0JNVg1q8mpQm1sLbZkWlu6WsPS0hFU3\nK/H/bpawcLWAhbMFLJwsYCGzgLmjOcwdzGFmZwZzO3OY2ZqJh40ZzKzNAHNApyuFRnMF1dWZqKm5\nhurqbNTW5qC2thB1dTdgbu4AK6tusLT0gKWlOywt3WBh4QILC2dYWDjD3NwR5uYOjQ5bmJnZwczM\nFmZmNjAzs0F5uQpZWdeRlZWFrKws5OTkIC8vD3l5eSgoKEB1dTV69OgBT09PuLu7w9XVFW5ubnBy\ncoKTkxNkMhkcHR3h4OAAe3t72Nvbw87ODra2trCxsYGFhQ1KSmyQl2eNnBxL5OUJyM0F8vKAwkKg\nuBgwGIBu3cS/Aze3Py4FJyfxkMlavhysbQwo1+WjoDoT+ZprKKrKRYEqD3mVeShSF+GG+gZUtSq4\n2bnBzdYN7nbucLV1hbONs3idW/9xnTtaO8Le0h72VvawtbCFraUt7HRmcMgvgW1eEVznP9v5P0a0\nw9ELgrADQCAAN0EQcgF8AsAKABHReiI6IgjCdEEQMgFoADxvVI/uEq5fv87ETl0dEB4O7N8PnDgB\nKJXAuHHA2LHAkiWAvz/Qrx9gacmkuRYxkAEJRQk4knEEp3NO43zheXjYeWC092gM9hiMZ4c9C38P\nf/Rx7gNbS9tbPrt432K8OfbNTjZsAFJSgLAw4MwZ4MIFoLISCAgABg8Wj9mzxUHo2ROwtjbyexIS\n1WqcVipxXqXCBZUKhbW18Le3h5+dHfzt7PCUpyd8bWzga2MDT0vL9t2c6vm+sBC9bGwAALoKHSrO\nKVB5rhKqBBXUiWro1XrYD7KH7UBb2A2wg/tcd9j42sC6lzWsvawhmHds0kdEqK7ORIXqPFQ3zkOl\niodGkwKAYGc3CHZ2A2Fj0xfu7jNhY9MbVlbesLLqDnNzm/aPmcGAzMxMxMbG4OLFi0hKSkJycjJq\namowYMAA9O/fH/369cP999+PXr16oWfPnvD29sbSpUsRFBTUrjbKyoC4OPFvISlJPDIzRSfev7/4\n6+/dG5g+XZzw9Ogh3usdHJq9d/+JQlUhLhReQFxxEpKzk5FckozM8ky42bqhn2s/9HXpC18nX0zw\nmYBeTr3Qw7EHujt0h6utK8yENhTy2lqxwzGJ4uys4SgpAXx8xM4bSZuOnogWtOOc143uiUSHiI8H\n1qwBDh4E/PyAOXOApUvF/5uZaOWFiHAu7xyCLgXhUPohyKxlmDFgBpaOWYqxPcfCw96DT8N1dcDJ\nk8Du3cDRo+K0a/Jk4Mknga++Ev8wOuBc20Ku1eJAWRmOyOWIUCrhaWWFQGdnTHFxwfs+PvCzs7tl\n1t1ZDDoD6m7UIevdLChOKlCVUQXZaBlk42XwesELDiMcYNPbpkM3juaors6GQnECSuUpKBSnYGZm\nBUfH++DoOAq+vh/DwWEoLC09O92OwWBAUlISTp48ifDwcERHR0Mmk2HMmDEYNWoUHn30UQwZMgTe\n3t6tttHae9nZ4r391CkgJkZ09KNGiff3adOAZcvEvwVb2xZNtIjeoMel4ks4lX0KUXlROF9wHjW6\nGoz2Ho3h3YZjxoAZWHb/Mgx0Gwh7K/uON5CfL05KIiPFu1Nqqng3GjlSnJhMnizO0Hx9/3isNvJ3\nLphSShEEgSTpRiQiIgKBgYEd+ozBAPz+O/Ddd8D168DrrwOLFgFeXly62CKKagU2JWzChoQNICIs\nGbkEc/znoL9r/07Za/dYXLoErFsH7NsH3HMPMH8+8MQT4lSNMVV6PXaXlGBXSQliKivxiIsLnnB3\nx0MuLvA28qmgMaQnKMIVKP61GPIQOVLcUzDl6SlwfdQVjqMcYWbF5q6tViehrGwfSkv3oa7uBlxd\np8DZeRKcnSfB1raP0fbr6uoQFhaG4OBgHDp0CE5OTpg8eTImT56MiRMnolu3jkdcN74utFrRN+7f\nDxw5AlRVif7woYfEJ1hjJzglmhKEpIfgUPohnLp+Cl4OXpjUexIm+k7EGO8x6O3cu/M3WI1GfNQ4\nehQ4dgxQqYCJE8Vj3Dhg2DBRE2oFQRCM0uglR99FuHQJePVVoKYG+Pe/xRm8BYul9A6gqlXhf7H/\nw+qY1Zg2YBpeuvcljO813ugZZqsYDMDhw8Dq1aLe/uqr4t3Nx4dLc5lVVfixsBBBN25gnJMTFnXr\nhhlubrA3ZsGiGWryalD4UyGKfy2Gpacluj/XHR5zPWDtze4motUqUVKyA0VFG6DVyuHhMRfu7rPh\n5HQ/BMH470NEiI2NxYYNG7Bv3z74+/tj7ty5mDVrFvr27cvAPhAVBWzeLD659u0rqnAzZ4oTX2Mv\nO3mVHLtTdmNH0g4klSRhSr8pmDlwJqb0m4LuDt2NM65SiZ3etQuIiABGjxZ1o6lTO9V5Yx09iMhk\nh9icBBHRqVOn2nVeRQXR0qVEnp5E69cT6fV8+9UcOr2Ovo/5njy/9qSFwQspvSydqf0WxyIiguje\ne4kCAoi2byeqrWXabmMyNBp6OiWFPKKiaFlmJl2rquLSjjpZTVeevUKRLpGUvjSdVJdVt7zf3uui\nNTSaDEpL+zudOeNEycnzSS4PJYOB3YWjVqtpzZo1NGTIEOrfvz998cUXlJ+fz8x+SQnRZ58R9ehx\nigYNIvr6a6JevXwJwF/+8PX1bXZM6n1np32vieeEEh0hLQ14/HFgwgRxvdHd3fR9yCzPxPO/Pw8B\nAsKeDcMQzyH8G83KAv71LyAhAfjiC+Bvf2OquTemsLYWn16/jr2lpXirVy+sHzgQDhwelaqvVSNr\nWRYqz1bC+5/e6L+6Pyxd2K6Qq9XJyM39DOXlofD2fhVjxmTAyordOolSqcS6devw/fffY/z48Viz\nZg0efPBBZk90164B334LbN8uPrF+8AHwyivir/7dd3Puioi9pmOZVV2N/+XnG2/YmLtERw9IM/p2\nc/QokYcH0caNt6d9g8FAa2PXktuXbvRd9HekZzgjbKVR8bHF3Z3o88+Jqqs5NmWgXwoKyD0qit7N\nzCR5XR2XdrQVWspclkmRbpF0/bPrpKvSMW+jtraU0tJeoqgoT7p+/XPSaiuY2q+urqbPPvuMXF1d\n6bnnnqMrV64wtX/9OtGiRURubkTvvUdUWPjnc+4W39HwPaOUSpqdlETuUVH0n6wso2f0kqO/A/nu\nOyIvL6LIyNvTfp2ujhYfWEwBPwdQWmmaaRotLSWaNYtoxAgixo6kKdlVVfRwYiKNunCBLqtUbX+g\nk5SFlNHZHmcp9flUqimsYW5fr9dSXt7/KCrKg9LTl1JdXTlT+waDgX777Tfq3bs3zZkzhzIzM5na\nVyiI3n2XyNWV6KOPRJmyJe4W3wGAHklMpD7R0bQuP5/UOt3N10ly9F2PlrTYdeuI+vcXZzm3g4qa\nCnrk10fosR2PkbpWbZI2T/38M1HPnkT/+hdRDXuH2Jj9JSXkHhVFX+TkkJbTgodOraOrL1+lc77n\nSBGh6NBn26vRV1VlUnz8WEpImERqdUonetk6OTk5NHnyZBo+fDiFh4cztW0wiEsu3boR/f3vRAUF\nzZ/XeCzuFt8BgH4qKKC6Jtem5Oi7KM39Qe/ZQ9SjB1FWlun7Q0SUX5FPw34cRi8fepm0eq1pGg0N\npVNOTkTBwVyb0RsMtCI7m3qeO0dxrU0djUSVqKKYgTF0ZdEV0io7PoZtOXqDwUCFhRsoKsqd8vJW\nM11kbbC/efNmcnd3p88++4y0WrbXQV4e0YwZREOHEp0/3/q5d6ujb+V1ydF3dU6dEjX5ixdvT/vy\nKjn5r/WnTyM+JYPBYJpGt24Vw4nOnOHajEqrpblJSTQuPp6KOD4xlB4spSj3KLqx/QYX+zqdhpKT\n51Nc3DBSqZKY21coFDRr1iwaNmwYJSYmMre/Y4e4/LJiRccDqO503/Hggw+SjY0NOTo6koODA/n5\n+XXKjuTo/8JcuSI6+bCw29N+tbaaJmyaQG8de8t0jf7wA5GPD1FyMtdmKrVaGh8fT89duUI1nKQa\ng8FAud/l0lmvs1QRw+dpoaamgC5cGEVXrjxDOh37ReorV67QwIED6Z///CfVML4Z1tQQvfqqKEkm\nJHTOxp3uOwIDA2nTpk1G2+Hl6KU0xbeJiIgIAOKOv0WLgE8/FXf5mRq9QY9n9j0Db0dvfDPlG9M0\nuncvsHKluH998OCbY8GaSp0OUy9fxhB7e2zy84M1h9wQRITMpZm4sfEGAqIDIBsjM8pec2OhUiXg\n4sWxcHefDT+/XzuUZ6Y9HD58GA888AD+/e9/4/vvv4c1w52/16+L4cE3bohpiEZ0oNAor+uCF6I/\nvjORHP1t5vPPxfj4l166Pe2/E/oO5NVyBM0Kajv5EgsiIsTdrYcPi1sdOdHg5Ic5OOCHgQNhxiEO\nn4iQ+UYmVOdVGBk1Eja+bB0wACiVZ3D58qPo1+9b+Pq+z3wX8po1a/DSSy/h4MGDWLJkCVPb8fHA\n/fcDTz0l3tudnJiav+N477334OnpiYkTJ+L06dO3uzu3YszjQEcP3OGPX6YmPl6UbBhuKuwQ+1P3\nU9//9SVFdcciQzrNpUviFz55kmszNXo9Tbx4kV65epX0nNYbDAYDZbyVQRdGXejUomt7UChOU1SU\nB5WX89H0PvvsM+rXrx9lZ2czt338uPir3r+fjb073XfExcWRWq2muro6CgoKIkdHR7p27VqH7bT0\nPSFp9F2TmhqiwYPF9cjbQWFlIXX7uhudzT1rmgblclGT37GDazMGg4GeT02lOUlJXJ185rJMOj/y\nPNWV89loxdPJGwwG+uCDD8jf358KWoptNIKGNXaW+0Da8h1iZhzjD1ZMnTqV1q5d2+HP8XL0knRz\nm3jxxQgMHAgsXGj6tokISw4uwUv3voT7e91vigaBF18UM1I9/fSf3mapxX6Xn4+LKhV+9ffnItcA\nQMG6AsgPyzH8xHDmaQwiIiJQUXEOKSnzMGjQLri4sF24ISIsW7YMISEhOH36NHr06MHU/ubNwHvv\nickaJ0wwzlZHrgtWrp4V9UnI2Bk0EsnR3wYKC8V0w2vWcEvh0irrzq9DeXU5PnzgQ9M0+NNP4qrc\nl19ybeaIXI5v8vJwcOhQ5tkmG1CEKZCzMgdDDw2FpRv7ai61tQVISZkLP79fmTt5APjmm29w7Ngx\nhIeHw8ODbb2AXbuADz8U88QPHszOrt6gZ2eMAxUVFQgNDUVtbS30ej22b9+OyMhITJ069XZ37Q+M\neRzo6AFJuiEiMdTsnXduT9tXy66S+1fuzDNQtsjly2Lw9NWrXJvJrKoij6goOqtUcmujKrOKojyj\nqPwU21QDDdTVKSg21o/y83/gYn/r1q3k4+NDeXl5zG0fOCDudE1iHN6vN+hp0b5Fd7RGX1paSqNH\njyaZTEYuLi40btw4CutkrHRL3xOSRt+1yM4Wc3uUlt6e9mfumElfRX1lmsaqqogGDSLasoVrMzqD\nge6Pj6fvcnO5taGt0FLsoFjK/4HPyrleX0cJCZMpPX0pF/uhoaHk6elJyRz2LZw8KS68XrjA1q7B\nYKC3jr1FEzZNuKMdPUt4OXpJujExK1YAr70GJCdHmLztsGthSC5Jxhtj3jBNg59/LpZEe7b1wsbG\navRf5+bC2swMb/TsaZSd1kh/OR1OE5zg/Yo3F/uZmUthZmaD/PyZzG1fuXIFCxcuxN69ezGYpaYC\nsRbMggXAnj3AvfcyNY1X1r6C0KxQHHzqIFvDXYyCHwqMN2LMXaKjB+6Su3JLpKaKMx+lkk2BiY6g\n0+to+I/DaU/KHtM0mJEh5p1th0xgzFgkVFaSR1QU5XBMaXxjxw2K9YvlkmKYiKi4eDfFxPQnrbaC\n+XVRWVlJfn5+tHnzZqZ2iYjKy4kGDCDasIG5afo18VfyfM2T8irE6+du8R1Nv+eNbTforPdZSbrp\nSsyfL6ZZvx1svLiRJmyaYJo8NgYD0bRpRF9+ybWZap2OhsTFUVBREb82cqspyiOKKi9UcrFfVZVN\nUVEeVFHRRoavTmAwGGj+/Pn04osvMret1RI9/DDRm28yN00R2RHk+bUnpZT8kZXzbvEdjb9n6cFS\niuoWRepkteTouwrp6WJssdo0mX9vobKmkry+8aK4/DjTNHjgAJGfH9fSf0REH2Rl0ZykJG43L4Pe\nQAmTE+j6yutc7Ov1WoqPH0c5OV9zsf+///2PRo4cSdUcnnbeeINo6lTR4bMkryKPvL7xotDM0Fte\nv1t8R8P3VEQoKMojiipiKxq/Ljn6O51//Yto2bI/fjaldPNx+Mf0zL5nTNOYRkPUu3eHdr92Ziwy\nq6rILTKSCjhmo8xbnUfx4+JJr+WTDO3atY8oMXHKLamGWV0X0dHR5OnpSVkccl4HBxP16SMWDmFJ\njbaG7vvlPvoi8gsiunvTFKuvqCnKI4rKw8pveZ2M8L1SzVgTUFMDBAUB0dGmb1tdp8YPF35AzAsx\npmnw66+B++4DJk/m2sw7mZl4p1cv9GCYgKsxNbk1uP7pdQTEBsDMgn3MQkVFDIqKfsG99yZAYJxj\nSKPRYNGiRfjpp5/Ql3E+obw8sY7rwYOAszNT03jj6BvoJeuFZeOXsTXcxUiamYR+X/eDy0Mu7Iwa\nc5fo6IG75K7clG3biB555Pa0/V30dzTvt3mmaay8XFyA5Vw5JVQup77R0VSt47M4SkSUPC+Zspdn\nc7Gt19dRXNwQunFjJxf7b7zxBj3zDPsnOJ2OaOJEPutMmy5uIv+1/lRZ0/xayN3iOwBQ5rt/LtkI\naUZ/5/Pzz8Cbb5q+Xa1ei+9ivsOeJ/eYpsHVq4EnnuCalVJrMGBpZib+r18/2HDa/Vp+shyqeBX8\nfvXjYj8//1tYWXnD0/NvzG2fPn0ae/fuRVJSEnPbq1YBlpbAMsYT7szyTCw7uQwRz0XA0dqRrfEu\nSN/POfz9GHOX6OiBu+Su3JjkZLHQd12T3Fem0Oi3X95OD2x+gHs7RGTUbL4jY7E6L48eTkzktgCr\nr9VTrF8slf7OZ0dbVVUWRUa6UVVV85kNjbkuVCoV9e3blw4ePNhpGy0REyPufGWdA02r19K4DeNo\ndfTqP713t2r0rbzead8rbZjizPr1wAsviDMhU0JE+Prc11h2v4n0ThPM5it1OqzMycHq/v2Z52Vv\nIP9/+bDpYwO3mW7MbRMR0tNfhY/Pu7C17cPc/n/+8x9MnDgRM2ey3XRVVyfmpFu9GmCcAw1fRH0B\neyt7/HPMP9kaNiHr1q3D6NGjYWNj86ec/mFhYfD394eDgwMmT56M3Nzc29NJY+4SHT1wl9yVG9Bo\nxHQH16+bvu0TWSdo0LpBpGdcPLpZTKTNr7x+nRampLR9YiepKayhSLdI0qRruNgvLt5FcXFDSa9n\nn9r4/Pnz1L17dyovZ5+HZ+VKounTxe0RLDlfcJ48vvK4uSmqNe5k37F//376/fff6dVXX6Xnn3/+\n5utlZWXk5OREwcHBVFtbS++++y6NHTu2VVstfU9I4ZV3Llu2iH8gt4MpW6fQpovG17BsFx9/TLRk\nCdcmKrVa8oiKolSOGxGuvnaVMt7O4GJbp6umc+d8SKFgXwhdr9fT2LFjmdQsbUpamngPz8lha7dG\nW0N+a/1oZ1L7FqS7gu/48MMPb3H069evp/Hjx9/8WaPRkK2tLV1tJcEfL0cvSTcc2bMHeOaZ5t/j\nWQ8zXZ6OxBuJWDB0Abc2bqJSAWvXAh980GkT7RmLtQUFeNjFBX729p1upzVqcmpQsrMEPv/x4WK/\noGAtHBxGwtl5Yqvndea62Lp1KwwGA5577rlO9q55DAbgH/8APv4Y8GE8LF9EfYF73O7BU0OeavGc\nrlYztikpKSkYPnz4zZ/t7OzQv39/pKSkmLwvUtQNJyorgTNngB07TN92UGIQFg5dCGsLPjHmt7Bt\nGxAYyFWbV+t0+C4/HxEdqSzdQXJW5qDHyz1g5WHF3LZWW468vC8xYkQkc9uVlZV47733cODAAZgx\nLn6+aRNQWysm4WNJujwda+LWIOGlBLaG7zDUajU8PT1veU0mk0GlUpm8L5Kj58SRI8DEiYBM1vz7\ngYGBXNrVG/QIuhSEowuPcrF/C0TibH7tWqPMtDUWPxQWYpKzMwZxms1XZVahdH8pxmSM4WI/N/dz\nuLvPgb192+GaHb0uPv30U0ydOhX33XdfJ3vXPEqlWETk2DGAZRQrEeGVkFfw/sT30cupV6vndmQs\nhBVsFufpE3ZVoRwcHFBZWXnLaxUVFXB0NH0IqeToObFvHzBnjunbPXntJLo7dMfQbkP5N3bqlPgv\np5sWAGj0enybl4eTjR6BWZOzIgc93+jJvCwgAFRXX0dR0SaMHp3M3HZ6ejqCgoKQnMze9qpVwGOP\nAawforYnbUd5dTnzVNksHTQrBg8ejKCgoJs/azQaZGVlMU8V3R4kjZ4D1dXA8ePA44+3fA4v/XFz\n4mY8P+J5Lrb/xNq1wOuvG10PsbWx2FxUhPudnDDEwcGoNlpCc0WD8uPl6Pkmn1z2169/BG/v12Ft\n7dWu8ztyXXz00Ud4++230a1bt072rnmyssTarytXMjULRbUC7554Fz8/9jMszNqeY3YVjV6v16Om\npgZ6vR46ne5mScHZs2cjJSUF+/fvR21tLVasWIERI0Zg4MCBpu+kMSu5HT3QBVbOWXDgANGkSa2f\nw2PDVHlVOTl97kTyKjlz238iJ0eMHVWpjDbV0ljoDQa6JyaGTrPOntWIlKdT6PrnfOJf1eorFBXl\nQVpt+1Mct/e6SEhIoO7du5OaQxTSnDlEn33G3Cy9efRNeunQS+0+v6tsmFq+fDkJgkBmZmY3jxUr\nVhARUVhYGPn5+ZGdnR1NmjSJctoIX2rpe8IU4ZUApgJIA5AO4N/NvO8G4CiARABJABa3YKeDQ9g1\nefZZojVrTN/uD3E/0Pw9803T2HvvES3lU/augWNyOQ2Li+O2C7Y6p5oiXSNJq2Sca7eeK1eepezs\n/3KxPWPGDPr++++Z2z11isjXl4h1ZuMMeQa5felGxeriTn3+bvEdt83RQ5R3MgH4ArCsd+Z+Tc75\nBMDn9f93ByAHYNGMLcbDcudRVydOdDnUX26T0etH05H0I/wbqq4WS2VxLvg949Il2lBYyM1+xtsZ\n3OLmq6quUWSkK9XVsX8aOXfuHPn4+FAN4xTNej1RQADRrl1MzRIR0dzdc2nVmVWd/vzd4DuIbm8c\n/X0AMogoh4i0AHYBeKLJOTcANCwlOwKQE5GuHbb/ckREAAMGAG2VL2WtP6aUpKBAVYAp/aYwtdss\nwcHiKh0jrbG5scisqkKsSoUFTcLTWKGr0OHGlhvouZSPNp+X9xV69HgJlpYdy+Xb1nVBRHj//ffx\n8ccfw5pxiuZ9+wAzM2D+fKZmcTb3LOIK4vDW2Lc69LmuotF3BdoTdeMNIK/Rz/kQnX9jfgEQJghC\nIQAHAOzT8nURgoNvT7TNruRdWDBkAczN+GR0vIWtW4Hn+S74rissxAvdu8OWU4bKog1FcJ3iChsf\nG+a2a2sLUFKyG/fdd5W57bCwMBQUFDDfHKXXA598AnzzjdFr67dARHgn9B2semgVbC1t2RmW6Bht\nTfkBzAWwvtHPzwD4vsk5HwBYXf//fgCuAXBoxharJ5w7EoOByNtb3DZuagatG0TRedH8GyosJHJ2\nFhP5cEKl1ZJrZCRd51TwW1+np3O9zlHF+Qou9jMy3qKMDA7FVIlowoQJtH37duZ2d+4kGjuWfT6b\nXUm7KODnAKNzLv3VfUcDLX1PmCAffQGAxhuge9a/1pjxAFbVe/IsQRCyAfgBuNDU2OLFi9G7d28A\ngLOzM0aMGHFzY0TDo1pX/Xnr1gjU1gIDB5q2fa8hXlBUK1CVUYWIzAi+7f32GwKfeAKws+P2fa4M\nGIBAZ2dkx8Qgm4N9/yJ/2PSxwUX1RSCCrX2ttgK2tlswenQS8/FZs2YNsrKyML9eW2Flf+LEQCxf\nDrz4YgROn2bX37DwMLx94G1sfWsrzAQzo+3dTURERGDLli0AcNNfGkVbdwIA5vhjMdYK4mKsf5Nz\n/g/AJ/X/7wZR6nFtxhbTu9+dxrp1RM89175zWYZXfh75Ob1y+BVm9lpl5MgO1YNtD03HYnhcHJ3k\nkIWRiMiv8fXOAAAgAElEQVRgMND5e89T6UE++eazsz+l1NQXOv351q6LmTNn0g8//NBp2y2xdSvR\nhAnsZ/ObEzbTg5sf7PTnu0p4JUta+p7gPaMnIr0gCK8DCIUYgbORiFIFQXipvvH1AD4HsFkQhEsA\nBADLiKjc+NtQ1yIsDJg1y/Tt7k/bj1UPreLfUEoKUFLCdSdsokoFhU6HSawLktajOq+CTqGD2wz2\n+eYNhjoUFv6IYcOOM7edkpKCuLg47N69m6ldnQ5YsUKsm8BSm9fqtfjvmf9i0+Ob2BmV6DSCeLMw\nUWOCQKZsz5QYDICHB5CUxL44Q2vkV+Zj+E/DceOdG7A051zd5D//EfPbfPkltybezMiAzMICn/Zh\nX5gDANJeTINtf1v4/seXue3i4p0oKvoFI0aEM7e9ePFiDBgwAB8YkSW0OX79VUxexjrAZVPCJmy7\nvA3hz7EZC0EQ8Ff1HY1p6XvWv97pW7GU64YRiYmAp6dpnTwAHEg7gMcGPsbfyRsMwPbtwFF+ydLq\nDAbsKClBdEAAF/u6Sh3KgsswOnU0F/sFBWvQqxf7il65ubk4ePAgsrKymNolAr76Cvj2W6ZmodVr\nsfLMSgTNCmr7ZAmTIOW6YURYGPDQQ+0/n1WM8L7UfZjtN5uJrVY5fRpwcwOGDGFuumEsjsjl8LOz\nQz9bPmF4xTuK4TzZGdbd2advrqw8j9raQri7G1fGr7nr4rvvvsOSJUvg4uJilO2mHD0KWFgAjzzC\n1Cx+vfQr+rr0xUTf1nPvtwWrvxHetFRKMCcnB2ZmZpDJZHB0dIRMJsOqVSaQWJtBmtEzIixMLNJg\nSsqqyhBfFG+aTVLbtgGLFnFtYsuNG1jcvTsX20SEop+L0PdLPnnzCwrWwNv7dQgC27h/pVKJoKAg\nJCUlMbULAF9/DSxbxl6bXxm5Eltnb2Vn9A7H29sbH330EY4fP47q6upb3hMEARUVFdxqHLcXaUbP\ngLo64Ny5jq1RsggZO3T1EB7u+zDsLO2MttUqOh3w++/st0zWExgYiJK6OkQolZjn4cGlDVW8Cjql\nDi4Ps50VA0Bt7Q3I5Yfg5bWk7ZPboOl1sWXLFkydOhXe3t5G225MXBxw7Rrw5JNMzWJ3ym70du6N\nCT4TjLbVVcIqZ82ahccffxyurq5/eo+IYDAYjLJfXm784r7k6BkQGytmA2jm98yVA1cPmEa2OXMG\n6NMH6NV6oQhj2FFcjJnu7pBZ8HnILFpfBK+/e0EwYz+zKipaDw+P+bC0ZHsBGAwGrFu3Dq+//jpT\nu4A4m3/7bcCS4dIOEeGrs19h2f2M1yka6h50QQRBQO/eveHj44MlS5ZALpd36PNEemRmvml0PyRH\nz4CwMGDy5I59xlj9sU5fh4jrEZjWf5pRdtrF/v3AbH43lIiICAQVF3OTbXQqHUr3lKL78+ztE+lR\nVLQB3t6vMrHX+LoIDQ2Fo6Mjxo0bx8R2A5mZYpTNCy8wNYvQrFAQCFP7T2Vi7+ZYrFjBxJ6pcXd3\nx/nz55GTk4P4+HioVCosXLiwQzaKi7fB0tL4p1zJ0TOgowuxLIjOi8Y9bvfAzY59PPgtEAEHDnDd\nIJBdXY3SujpusfMlu0vgHOgMay/2i7Dl5SdgZdUNDg7sK2CtXbsWr7/+OnN999tvgZdfBljXcvnq\nnDibZ9rfmBggJ6ft8wSBzcEQe3t7BAQEwMzMDB4eHli7di1CQ0Oh0Wja9XmDoQ7Xry9Hnz7GL+BK\ni7FGolYDCQnAhA5Kksbqj6FZoaZZhL1wAbCzA/z9uTWR4+eHeTodzDgtWBVvK+aWpbKoaAO6d2c3\nNW64LrKyshAbG4s9e/Ywsw2ItWB37gRSU5maRXxhPNLl6XhqyFPMbAYGBooZAt95B/jnP1s/uYvE\n2AuC0G7NvqhoA+zs/ODsbFz0EiDN6I0mJkbM2MupbnWLhF4zkaNvkG04Rg3sKS3FfE7piGtya6BJ\n0sBtOvsnn7q6UigUJ9Gt29PMbf/4449YvHgxbBmHmgYFAdOmAaxVsq/PfY23xr7Fdj/H1atAVBT3\nTKnG0lIpwbi4OKSnp4OIIJfLsXTpUkyaNKldxcH1+irk5KxEnz5s6jlKjt5IYmKAzkioxmj0ZVVl\nSJenY2zPsZ220W446/NXNBqUxMVhrEzGxX7JzhJ4zPWAmTX7S724eCvc3R+HhYUTM5sRERGoqqrC\nli1b8MorrzCzC4h73tatA157jalZZCuycfLaSfw94O9M7Ua8/bbYWVPPojrIypUrYWdnhy+//BLb\nt2+HnZ0dVq1ahWvXrmHq1KmQyWQYNmwYbGxssGPHjnbZLChYC5nsfjg63sukj5J0YyQxMaafcIRd\nC8ODvg/CytyKb0NpaUBlJTCaz05SQJzNP+jszE+22V6MAWsGMLdLRCgq2oiBA39gbnvXrl0YO3Ys\n+vZlG/N/8qToM++/n6lZfB/7PV4MeBGO1m3PVNtNUZEY7fXrr+xscuKTTz7BJ5980ux7Tz3VcSlL\np1MhL+8bjBgRYWTP/kCa0RsBkRhaObYTE2tjNHqT6fP794uLsGb8LpM9JSV4a6Zxu0lbQn1ZDZ1S\nB6eJ7GbcDVRWxoKoDk5ODzC1GxgYiA0bNuAfHHbfrV0rTpBZ3lPVdWr8evlXvDqaTdTRTVavRuDz\nz4u7se8yCgt/grPzQ7C3H8TMpuTojeDaNcDGBmC8l6VViMh0+vyBA9xlG6VOh3GcZJvi7cXwXODJ\nJXb+xo2N6N59CfOImNTUVGRnZ2P69OlM7WZni5v6FixgahbbLm/Dg74PwsfJp+2T20tlJbBhgxjo\nfxeSn/8tfH3ZJq+THL0RxMR0bjYPdF6jTytLg7lgjgGu7OWIWygsBDIygAcf5NbEntJSzPPwwJnT\np5nbJgOhZEcJuj3TjbltvV6D0tK96N6dbTk/AFi+fDmeffZZWDDeOPbTT8DixWIAFSuICGvj1uKf\n97UREdNRNm8GHn4YEdevs7XbRZDJxsLBYShTm5JGbwTGOPrO0iDbcM+dcfw48PDDbLdONmFPSQl+\nvuceaAuaFiwzHuUZJSxcLeAwhHGwOICysoOQycbC2pptqlKtVosTJ04gOjqaqd2aGjEVcUwMU7OI\nuB4BAAjsHcjOqMEArFkjavN1dezsdiF8fT9kblNy9EYQEwP8rZNl0Dur0YdeC8XzI0yw+nvsGDCV\nzQ7H5khtJNuYcchpUrytGN0Wsp/NA0BJyQ54enZsh2N7CAkJweDBg3HPPfcwtbtvHxAQAPTrx9Qs\n1sStwev3Md7QdfQo4OwMjBuHwNucCOx2wSrSpjGSdNNJqqvFgkucUqc3S62uFpE5kXioD+dtuHo9\ncOIE8Oij3JrYV1aGuR4eXKJtDFoDyg6UwfNv7GPz6+rKoFRGwt39Cea2N23adEuaW1Zs2AC8+CJb\nm7kVuTidcxrPDHuGreHvvwfeeIPrvo27EcnRd5KEBGDQoM5rnp3R6OMK4nCP+z1wteWcPe38eaBn\nT66rzAfLyvC4uzsA9nnHlRFK2Pa1hY2vDVO7AFBauhdubtNgYcEwlBBAUVERIiMj0Z3xTqasLCA5\nGXj8caZm8eP5H/HssGfhYMVQGktNBS5duvmY3FXy0XcFJOmmk9wOff5Mzhk84MM2nK9ZOMs2xXV1\nuFpVhYlO7MMeAaA0uBQec/mkOy4p2c6lilRQUBDmzZvHfCfs5s3AwoWANcM0P3X6OmxK3ITI5yPZ\nGQXE+M9//INtZyUASDP6TmOso++MRh+ZG2l01Z52wdnRh8jlmOLqCqv6+HyWecdJTyg7UAb3ue7M\nbDZQU5MDjSYVrq5sJS0iwpYtW/D8888zHQudTnT0rLNU/p72OwZ5DMJAt4HsjCqVwI4dYra1erpK\nPvqugOToO4mpZ/Q6gw7R+dFMCjq0ilwuPkKPH8+tiUNyOWZy2ghTcbYCVt2tYNeffTGW4uKd8PCY\nBzMztjuSL168CJ1Oxzwd8fHjYgkB1tUff7n4C/N0B9i8WZxcmLroMgPq6urw4osvonfv3nByckJA\nQACOHTt28/2wsDD4+/vDwcEBkydPRm5ursn7KDn6TlBQAFRVGRfF0FH98dKNS+gp6wl3O/Yz1Vs4\ncUKMnef0+Fyj1yNcocD0Ro6epRbLV7bZgW7dGO84ArBt2zYsWLAAgiAwHYuNG9nP5rMV2Ui4kYA5\n/nPYGSUCfv75T0l4uopGr9Pp4OPjg8jISFRUVOC///0v5s+fj9zcXMjlcsydOxerVq1CeXk57r33\nXvyts6F6RiBp9J0gNhYYM8a0gQGRuZF/CX3+lFKJYQ4OcOMQn08GQmlwKYaHss8Nr1YnQadTwsmJ\n7ROVXq/Hrl27mDu14mIgPBzYsoWpWWxM2IiFQxfCxoLhQvfp04C5OdenSJ7Y2dnh448/vvnzjBkz\n0KdPH8THx6OsrAxDhgzBnDnijXH58uVwd3dHeno6Bg5kKH21gTSj7wQNjt4YOqo/nsk5w1+fJxKf\n9zk6+uZkG1ZabGVcJSxkFrAfxD7bYUnJLnh6PgVBYPsnEx4eDm9v75ux86zGYutWMXsFy+wSOoMO\nmxM3s5dtfv4ZeOmlP82cuqpGX1xcjIyMDAwePBgpKSkYPvyPiYednR369++PlJQUk/ZJcvSdICEB\nGDXKdO0RkbgQ68PZ0V++LJYdYpw1sQEiwmGO+nxZcBkX2YaIUFq6Bx4ejCtpA9i+fXuHy8u1h19/\nBZ5jnKHhaMZR+Dr5YrDnYHZGS0vFTVKLFrGzeRvR6XR45plnsHjxYgwcOBBqtRpOTaLLZDIZVCqV\nSfslSTcdhEh09CNGGGcnIiKi3TOWtLI0OFg5oJcTv+LcAMTZPMdNUpc1GlgKAvyabD7oyFi0BJEo\n2wzex9AJ1aPRJMNgqIWjI9u7e1VVFX7//Xd8/vnnN19jMRaXLgEVFcADjJW+Xy7+ghcDGO+82rIF\neOIJwMXlT291ZCwERtIXGTH2RIRnnnkG1tbWWLNmDQDAwcEBlZWVt5xXUVHRruIjLJEcfQcpLBSf\nML28TNdmZG4kHvA1gT4fHn5LeBtrDpWVYaabG5c8PZrLYh1Oh+Hsc9uUlgbDw2Me834fOnQIo0eP\nhhfji2nbNjF2nmV26UJVIaJyo7Bz7k52Rg0GYP16seyVkRjjoFnxwgsvoKysDEeOHIG5uTkAYPDg\nwQhq9P00Gg2ysrIweDD7CUlrSNJNB0lMFGfzxv7Nd2TWdibnDH/Zpq5OzGPLMVvlYbkcjzUj27DQ\nYssOlsHtcT43kdLSvfDwmMfcbnOyjbFjodeL4eislZCtl7Zirv9c2FsxXP84dUrM891CWGlX0uhf\nfvllpKWl4eDBg7Cy+iP8dvbs2UhJScH+/ftRW1uLFStWYMSIESZdiAUkR99hGhy9KTHJjP78eaB/\n/2YfoVkg12pxpaoKE52d+dg/JIf7TPahpxpNKnQ6BWQyI1ffmyCXy3H69GnMZpzvPzxcfNpkWcud\niBB0KQjPjWAs+rewCNvVyM3Nxfr165GYmIhu3brB0dERMpkMO3fuhLu7O4KDg/H+++/D1dUVFy5c\nwK5du0zeR8nRdxBWjr694XQ5yhzU6Gr4558/dQqYNImb+TCFAg84OcG6GT3B2NDC2qJaVGdUc6kk\nJco2c5lH2xw4cACPPPIIZE3CYowdi61bgWcY5xm7UHgBdfo6jO/FMPyxrAwIDW21s10ljt7HxwcG\ngwFVVVVQqVRQqVSorKzE00+LReMfeughpKamQqPRIDw8HD4+DIu0tBPJ0XeQxERg5EjTtdcwm+ee\nf56zoz9eXo5HXfkkYys/Ug6XKS4ws2J/OZeVBXORbX777TfMnz+fqU2NBjh4EKj3L8wIuhSEZ4c/\ny/Ya3L4deOwxMSWxBHckR98BVCpxMZaFvNZe/fFs7lm2M6nmqKkRNwewDtOoh4hadfTGarFlh8rg\nNpN9yGZVVSZqa4vg5MR2/OVyOWJiYjBjxow/vWfMWOzfLxb+7sYwDX+trha7U3bj2eHPsjNKJFZC\naSMlc1fS6O90JEffAS5fFvOG1C+om4SYghiM68k2B8qfG4kBBg9mu7umEVeqqmAhCBjAODMjAOhr\n9FCeUsJtGntHL87m50AQ2P7CDxw4gClTpsDenu3Grm3b2C/ChmSEYLDHYPR27s3OaEKCWBdWcuQm\nQ3L0HYBF/HwD7dEfNXUapMvTMaI759VfE8k2LT36G6PFKsOVcBjuAEs39ikVGvR51vz222948snm\nN191dixKSoDoaPZ554MuBeG54YwXYTdvFndztRH/2VU0+q6A5Og7gKkjbuKL4jHUcyisLTjn5+bs\n6EM56vPyQ3Iusk1tbQGqq7Pg5MRWzmpNtjGG4GBgxgyA5UNCqaYUp6+fxrxBDNcoamqAnTvFSuUS\nJqNdjl4QhKmCIKQJgpAuCMK/WzgnUBCEBEEQkgVBOMW2m3cGLB19e/TH2PxYjPFmG9b3J6qqgIsX\nuSWUqtbrca6yEpNbCdvsrBZLRJAf5uPoy8oOwdV1GszM2D4ptCXbdHYsdu3qfP3iltiZvBOPDXwM\njtYMd3EePAgMHw707t3mqZJGz442Hb0gxpWtBfAogMEAnhYEwa/JOU4A1gF4jIiGAGCfFOQ2o9UC\nV64AQ4ears2YghiM7ck56f3Zs+Ldy4H9jlIAiKyowDB7ezhZsN+ErU5Uw8zGDHb3sM89L5f/zqUu\nLI9om4ICICmJfS66bZe3YdEwxqJ/OxZhJdjTnhn9fQAyiCiHiLQAdgFo+hewAEAwERUAABGVse3m\n7efqVbGIAyt/2B79MTY/FmN6cp7Rm0Cfn9KGbNNZLbZhNs869FSnU6Gi4izzSlINss306dNbPKcz\nY7Fnj5guhmUJgQx5BnIrcjG572R2RvPzgbg4Ma1mO5A0ena0x9F7A8hr9HN+/WuNGQjAVRCEU4Ig\nnBcE4a+Riq4Rpo6fz6/MR62+Fn2c+/BtiLc+r1Dwi58/Wg7X6extl5cfh0x2Pyws2EYh8Yq24SXb\n/G3w32BhxvBJbNs2YN48wI79E5hE67BajLUAEABgGoCpAD4SBKE/I9t3BKwXYtvSH2PzYzG251i+\nG6U0GjFmlFNNxMLaWhTW1mJUG5n6OqPFasu10CRr4DSB/W5YXrJNcHAw5s1rfWGzo2ORnQ1kZQGT\nGU68iQjbk7ZjwVCG1bSIxG27HYj/7CoafWulBHNycmBmZgaZTHYzNcKqVatM3sf23K4LADTes9uz\n/rXG5AMoI6IaADWCIJwBMBxAZlNjixcvRu/6hRhnZ2eMGDHi5i+04VHtTvw5MRGYMiUCERGmaS8m\nPwaeJZ63pGpl3t769YCvLwLrZ1is7a8NCcGQigqYT5jA3L7ihALpg9Ohj9EzHR+DQQ8rqyPo0+dz\npv2trKxEREQEXmtULo+F/R07gLlzA2FpyW58HQc6QmfQoSqjChGZjK6/hAREyOWAVovATn7/O5XG\npQR79eqFkJAQzJ8/H8nJyQAAQRBQUVHRoUlbREQEttSXB+vdjoXrNiGiVg8A5hAdti8AKwCJAPyb\nnOMH4ET9uXYAkgAMasYWdUUMBiI3N6LCQnY2T5061er7EzdNpBNZJ9g12BwrVhC9+y4384tTU+mH\n/Pw2z2trLJojdXEq5a3J60SvWqe8/BSdP38vc7s7d+6k6dOnt3leR8dixAii8PBOdqoF3j72Nn0Y\n9iFbo2++SfRhx2w2Houu5juGDRtG+/bto+vXr5MgCKTT6dr1uZa+Z/3rbfrrlo42pRsi0gN4HUAo\ngBQAu4goVRCElwRB+Ef9OWkAjgO4DCAGwHoiumL8bejO4MYN8d/u3U3TnlavxcWiixjdYzTfhqKi\ngAlsa6A2QEQIUyjwEIdsmESE8mPlXHbDyuUH4e7OeNcRgP379zPPVJmeLl6bLDNX6A167ErZxVa2\n0enE2Pm/SBWptiguLkZ6ejqGDBkCQJzR9+7dGz4+PliyZAnkcrnpO2XMXaKjB7rYXbmBEyeIHnjA\ndO1dLLxIg9YN4tuIVkvk6EhUVsbFfIZGQz3OniWDwcDctipRRTH9Y5jbNRgMFB3dl1SqRKZ2q6ur\nycnJiYqLi5naXbmS6LXXmJqksGthNPKnkWyNHjlCNGaMUSa6iu/QarX08MMP0yuvvEJERGq1muLj\n40mv11NJSQnNmzePHn300RY/39L3hJEzeqnCVDtISRFz3JiKmPwY/hulLl0S40U51W8NVyox2cWF\ny2Ky/KgcrlPZR9tUVV0BkR729sOY2g0LC8OwYcPg6enJ1O7evcDq1UxNYkfSDrazeaDDi7CdIUKI\nYGInkAI7/VlqppSgvb09AgICAAAeHh5Yu3YtvLy8oNFomEdftYbk6NtBSgr71AeNF1mbElsQi/t7\n3c+2waZERgIT+VWtClcoMLWdYZWtjUVzlB8rh8+77HN6y+WH4eb2GPObU0dkm/aORVYWUFTEVnmr\n1dVif9p+XHr5EjujKhVw5Ajw/fcd/mhHrgtjHDQrmisl2ByCIMBgMJiwZ1Kum3aRkiImdzQVJpnR\nc9TnDUQIVyoxiYM+r6vUQR2vhnMg+zzmcnkI3NzY5qDR6XQ4ePAgc30+OFjcd8Qyk2poVigGewxG\nT1lPdkaDg8XylO7sq3/dSbRUSjAuLg7p6eliug65HEuXLsWkSZNMXhxccvRtQMTH0bc0U6msrURe\nZR4Ge3K8sxCJM3pOjj5Fo4HM3By+NjbtOr8js3lFuAKycTKY27NNHazVKqBWJ8LZuf19aQ9nz55F\nz5492x0i196x2LtX3HvEkt0pu/G3wYx3Xm3b1umSV3d6WGUDrZUSvHbtGqZOnQqZTIZhw4bBxsYG\nO3bsMHkfJemmDQoLxa3lppqQXCy6iOHdhrPdkdiUzEzAygrw9eViPlyp5BJtA9TvhuWgz5eXH4eT\n0wMwN2ebM59HtE1OjrhRimUd92ptNUIyQvB/U/6PndGiIiA+Xqwk9RemoZRgSzz11FNG2Y+LizPq\n84A0o2+T5GQ+sk3DZpCmxBfG416ve9k32JgG2YbTrttwhaLVbJVNaWksmkL1YZV8HD172YaIcODA\nAcyaNavdn2nPWAQHi7ltWOaJO5Z5DAFeAejmwLA81e7dYkc7WXCmvdfFX50VK1YYbUNy9G1g6oib\nC0UXMKrHKL6NcFyI1RkMOK1UIpBDLdDqjGqQnmDnzzZXCpEe5eXHmDv6pKQkmJub34ynZkWXkW12\n7AAWMI7gucvIyclBbGys0XYkR98GvBZiW9If4wvj+Tt6jguxF9Vq9LKxQbdGC1Jt0V4tVnFCAddH\nWq5U1VkqK8/Dyqo7bGzYRvIcPHgQjz/+eIf629ZY5OcDaWnAQw8Z2blGaOo0OJZ5DHP857AzmpEh\nakxGdLSraPQ82bhxIxYuXGi0HcnRt4EpI24qaipQqCqEn7tf2yd3lhs3gNJSbo8p4QoFHuIwmweA\n8tByuDzCXvsvLw+Bqyvb2Tzwh6Nnyf79YrnADtxH2yQkIwRjeo6Bux3DhaidO8WUmhzqENwt6HQ6\nbNy4EX//+9+NtiU5+lYgEouNmEqjv1h0ESO6j4C5Gcfq4+fOAePGtVmvs7NEKJWY1EFH3x4t1qA1\nQHlaCZeH2Tt6HmGVhYWFyMzMxIQOPjm1NRbBwcBcxmVsmcs2RExkm7tdoz9y5Ah69+7NRPqTHH0r\n5OWJhUY4BZD8iQuFF/gvxJ47B9zPZzOW1mDAucpKTOQwo1fFqWDbxxZWngynsgBqawtRU3MdMtk4\npnYPHTqEadOmwdKSXSnC0lIxXfYjjzAzCVWtCievncRsP4aRQQkJQF0dMIbzXpC/OOvXr8c//vEP\nJrYkR98KPGWb5vTH+CIT6PPR0dwc/UW1Gn1sbODWQefWHi22/AQf2UYuPwIXlykwYxzO2lnZprWx\nOHgQmDIFaOf2hHZxKP0QJvpMhIstw7FtmM0buZbSMBZ1Jt5FeieQm5uL6OhoPPkkm6qskqNvBV6h\nlS1xofAC7u3BcUZfWytOCe+7j4v500olHuSkzytOKDjp80eYyzZqtRqRkZGYyriI6/797a7C1272\nXNmDJwcxLPFsMIglrxhG2xwtL2dmq6uwceNGLFiwAHaMqnFJjr4VeIZWNtUfFdUKFGuKcY/bPXwa\nBMRH6oEDuRUCj+iko29Li9VV6KC5zL6alMFQB4UinHlt2BMnTmDMmDFwcup4f1sai8pK4MwZYAbD\ne5KqVoXw7HA8fg/DBeOoKMDVFRg0yGhTDWOxrbjYaFtdCZaLsA1Ijr4VTBlxY5KF2OhocSGWAzqD\nAWcrKvBAJ5xbWygjlJCNlcHclu3YVFZGw85uAKys2GaV5BFtc/SoGBErY1jGNiQjBON7jWcr2+za\nBTz9NDNzFTodQrvAjH7RokXw8vKCk5MT+vXrd0u5wLCwMPj7+8PBwQGTJ09Gbm5uq7aOHz+Onj17\nYtgwdllUJUffAgYDkJrKZGLSLE212AuFFzDKi7M+z3EhNkGtRi9ra3h0Iu6vLY2enz5/FK6ubOUV\nvV6Pw4cPY+bMmZ36fEtjwUO22XtlL1vZRqcTd3MxqlQeGBiI4NJSbuk0WPLee+8hOzsbFRUVOHr0\nKNasWYPjx49DLpdj7ty5WLVqFcrLy3Hvvffib22Mz4YNG/Diiy8y7Z/k6FsgJ0eMtuEwQW0W7gux\nRH+EVnKAuz4/hYc+fxSurtOY2oyJiYGXlxebOp/11NQAx46J8fOs0NRpcOLaCTzhx7AIeng40KcP\n0LcvM5Pbi4uxkHEefx4MGjQINvWr5EQES0tLeHh4YN++fRgyZAjmzJkDKysrLF++HJcuXUJ6enqL\ntiIiItq8GXQUydG3AG/ZpqkWy30hNi9PnHEx/CNszOmKik6nPWhNo6/JrYFOoYPDMLbrCrW1hait\nzYOjI9uF6cOHDxsl2zQ3FmFhwLBhQDeGaWhCMkIwruc4uNoyzBu0axez2TwA7AkNRYJajcc4Fcdh\nzQkmS4MAACAASURBVGuvvQZ7e3sMGTIEH3zwAQICApCSkoLhw4ffPMfOzg79+/dHSkpKi3bmzp3L\nPI2x5OhbgKds05Ty6nKUVZVhoNtAfo006PMcEpnpiRCpVOIBDjN6xQkFXCa7QDBj2+/y8mNwcXmE\neVjloUOH8BjjbI1dQraprQUOHADmz2dmMlypxCx3d9iwTLrPkXXr1kGtVuPEiRP48MMPERcXB7Va\n/adFeZlMBpVK1aId1rINIKUpbpG0NL77PRprsfGF8RjpNRJmAsf7Lkd9/pJajR7W1h3Kb9OY1jR6\nRZiCy27Y8vJjcHWdztRmdnY2SktLcZ8R4atNx0KvF+PnP/jAyM41okpbheNZx/HDjB/YGQ0NFUPU\nerIrWhLbty++budjTEQEm4lAYCAZ9XlBEBAYGIgnn3wSO3fuhIODAyorK285p6KiotUZ+xgOjkdy\n9C2QmgosXmyathJuJCCgewDfRqKjgf9jmGu8EZ0Nq2wLIoIiTIE+n/Vhatdg0EGhOIn+/f/H1G5I\nSAimT58OM4bpJc6dA7y9RembFUczjmKMN+PcNrt2AUbmXW9MqkaD4rq6dsuBxjpo1uh0Otjb28PX\n1xdBQUE3X9doNMjKysLgVnRhHnWWJemmGYjEGb0fx9xijbXYi0UXEeDF0dFXV4uLDqP4LPaeVirx\noBGr1i1p9JpkDcwdzGHbm20xkMrKGNjY+MLa2oupXRayTdOxOHAA6EA6+3ax58oezBvEMM9xVRUQ\nEsI0d/L24mKMz86GOaeaCSwpLS3F7t27odFoYDAYcPz4cezZswezZs3C7NmzkZKSgv3796O2thYr\nVqzAiBEjMHAgR5m2GSRH3wwlJWLOLw8P07SXcCMBI71G8mvgwgXxsbqTBSBaw0CEyIoKLjN6vrIN\n22gblUqFc+fOYcqUKcxsErF39NXaahzLPMY2t01ICDB6NMAoOoaIsKOkBA93gbBKQJyB//jjj+jV\nqxfc3Nzw0UcfYevWrRg1ahTc3d0RHByM999/H66urrhw4QJ27dpl8j5K0k0z8J7NA39osapaFfIr\n8/mmJuYYVpms0cDd0hJe1tadttGSRq8MU6LbIoahJvWUlx9F//6rmdo8efIkxo0bZ3S0ROOxSE4W\n93Mw3DeD0KxQBHgFwMOe4Sxm926m0TbRlZWwNjPD36exvRnzwt3dvdXIsYceegipqamm61AzSDP6\nZkhN5e/oG7hUfAlDPIfwrRHLcUfsGaWSy25Yg9YA5RklnCexfVKoqytBTc01yGRjmdo9fPgw82ib\nhtk8S/UiODUYc/0Z5jlWqYATJ4A57IqW7KiPneehVd+tSI6+GdLSAH9/vm00zAAuFl3EyO4cZRsi\nICaGn6OvqDA6rLK52ZDqvAq2fW1h5cE2LXF5eSicnR+CmRm79MEGgwEhISFMHH3jsWAt29Tp63A4\n/TBm+zOUbQ4dAsaPF/PbMEBrMOC30lI83a3bXZ+PniWSo28GU87oE24k8F2IzckRFxx69WJumoi4\nzegVYQo4T2av+4v6PNskZhcuXICbmxv6MtyMlpsrHuPHMzOJsGthGOQxCD0ce7Azyli2CVMo0NfG\nBv04rCfdzUiOvhlMMaNv0GK5z+hjYoCxY7lslMqoroaVmRl8jUyQ3pxGrzgpbpRiCZEBCkUoc0fP\nUrZpGIvffwcee4xtJb69V/ayjbZRKoFTp5g+duwoKcGC+th5qWYsOyRH3wS1Wqzk4+vLv60aXQ0y\n5BkY2m0ov0YaHD0HGmbzrLVUvUYPVbwKThPZPimo1QmwtHSHjQ3bXy5PfZ4VWr0Wv1/9nW0B8IMH\ngUmTmCWEqtLrcUgux3xThbvdRUiOvgnp6cCAAQDvXdcRERFILklGf9f+sLFgWDKoKdHR/Bw9A30e\n+LNGr4xUwjHAERYObBeoecg2BQUFyMnJwThGayARERFQKMSIWJYlA0/nnEY/137wcfJhZ5SxbHNY\nLsd9jo7oXh/BJWn07JAcfRNMEVrZQEIRZ32+pkaM0eO0UYqXPq8MUzKXbYAGR882LXFISAimTp0K\nC4YaS0iIOFFmVFwIABB8hXG0TXm5WGSkk+mYm2NHcfFN2UaCLZKjb4KpFmIDAwP56/MJCeKXYekx\n6smpqUGVwYB7GNhuqsXyWIjV6SqgVifCyekBpnZZyzaBgYHMZRu9QY/9afvZOvr9+8VHDkZZFhVa\nLU4plZjt/kdaBkmjZ4fk6JtgioXYBrjviOWoz0dy0ue1ci2qM6shu49hKSUACkUYZLLxMDdnF81R\nXV2N06dP49FH2clB1dViWDpLyf9s3ll4OXqhn2s/dkYZyzbBpaWY4uoKGcvVZ4mbSI6+Caaa0YeF\nhyGpJAkjuo/g1wjPhVhG+jxwqxarOKWA0wQnmFmxvTTLy48zl20iIiIwfPhwuDKKIQeA1asjMHIk\n4M4w39jeK3sxz59htE1pKRAXB0xnl/1ze0kJFjRJodCVNPqWSgnm5OTAzMwMMpkMjo6OkMlkt5QZ\nNBXS7bMROh2QlSXWz+ZNbkUuejj2gMya7cz1FqKjgZUruZg+o1Ti1R4M47Hr4aHPExHKy4+hZ883\nmdrlEW0TFcVWtjGQAftS9+HksyfZGQ0OBqZNA+ztmZjLr6nBZbUa07tIgZHmeO+99/DLL7/AxsYG\n6enpeOCBBzBq1Cj4+flBEARUVFTc1p2+0oy+EdnZQPfuXCTtP2HR14LvQmxBgZhVsH9/5qaL6+pQ\nrNViqAObqk+NtVge+nxVVRoAwM6O3aMaETF39Ho9cOFCIJ5gWN0vNj8WzjbObHMpMZZtdpeWYra7\nO6ybpHfuShp9S6UEG342GAy3s3vtc/SCIEwVBCFNEIR0QRD+3cp5owVB0AqCwDBY13SYNOLmRgLf\nhdjYWG4bpSKVSoyXyZinkOVVNlCUbR5lOqNKSUmBubk5/Bku6MTEiBMNlrnnmee2KSoCEhOBqexk\nsO1/kWib5koJAmJ2y969e8PHxwdLliyBXC43ed/adPSCIJgBWAvgUQCDATwtCMKf3GH9eV8AOM66\nk6YiNdV0C7Hhp8JNsyOWA2cYpyVu0GIVYQo4T3LmUjaQtT5/+PBhzJgxg+nN48ABYMSICGb2iIj9\nbti9e8WQSiN3QzfQUGCkueupK2n0wB+lBE+ePIkPP/wQ58+fh7u7O86fP4+cnBzEx8dDpVJh4cKF\nJu9be2b09wHIIKIcItIC2AWguYfLfwLYC6CEYf9Miqlm9ESEzPJMvtIN7x2xXST/vF5fjcrKs3Bx\nmczULmvZpiH3PMvcNheLLsLK3ApDPIewM8pYttlRUoKnPD2NfjoUBIHJYSyCIODBBx+8WUrQ3t4e\nAQEBMDMzg4eHB9auXYvQ0FBoNBqj2+oI7XH03gDyGv2cX//aTQRB6AFgFhH9CKDL5hY1laPPVmbD\n2c+ZbU7wxmi1wMWLYjEIxii0WmTV1CCAkT4PiFosEXFZiP3/9s47Popq7eO/SUJIQhISSIXQe0fp\nINJUbGDBhiKCF/V6LahXverVK3rFV7GL4LXQBSIQSpDeAiQQQkgCIYQkJKT37b3MPO8fk2CAlM3O\n2RWS+X7cj+zs7G/Pnp08c85znvM8avVR+PvfAi8vdhu7FAoF0tPTMWnSJGaamZlife3nnpvMTLN2\nNM9s1lFUJDaU0ZZdImp0k1RzfPRExOTBCrvdDr8GFvs4jnO7z57VYuy3AOr67m86Y19bPtAdrhuX\nlw48d0509Lpg12q8RoOxgYFow7AuKgAYM43g2nLw6ck2HYQr3Da7d+/G1KlTryy+sYB17nkiwpbM\nLWz985s3i410sgj8tSTpdGjDcUwHDX8FDZUSfOCBB5CUlITs7GwQERQKBRYuXIgpU6ZILlDTXBwJ\nrywBUDdBRlTNsbqMBBDNiUOHEAD3cBxnI6LYa8XmzZuH7t27AwCCgoIwfPjwK3fuWp/cX/G8qgqw\n2+Nw/jwwZYprPy+VT0VwebDrvs/588DYsS7R/62kBLfX+BdY6QNA7/TeyBmYA8tRC9P2XriwFU8/\nvZ1pe//44w/MmDGDaf9u3w48/ngcvv02Da+99ppkvfTKdOiydNBmaYGaKFjJ7f3lF+DZZyE+k673\neWwsxnl5gRszpt7Xv/3226vsw41KbSnBF198EUSEPn36YN26dRg1ahSio6Px3nvvoaqqCoGBgbjz\nzjuxYcOGJjXj4uKwevVqALhiLyXhwFTGE8AlAN0AeANIAzCgkfNXAXi4gdfoRiUujmj8ePd81j2/\n3UP/XfNf133AU08RrVjhEunRycl0VKViqnnkyBE698A5Kl9fzlTXaMyj+PgwEgSemabFYqGgoCAq\nKytjpllURNShA5HVKvYFCz44/AH9c98/mWgREVFeHlFICJHNxkTOxvMUHh9P2QZDg+fU7Ysb2Xaw\npKHvWXPcaZdUk/NvIuIBvAxgP4AMANFElMlx3Ascxz1f31uk3nz+CtwZWplSloK5D8x13Qe4KGOl\n3m5HhsGA0YynnZMmToLmqAZBU9ku8P4ZVsnOzRQfH4++ffsiIiKCmeb27WLKgzZt2MWOb7mwBY8O\nfJSJFgBxEfaRR5glyD+oUqG7jw/6NLJp5UYfyd9MOPSrEdFeAP2uOfZTA+c+y6Bdbsddhr5MVwab\nYEOXQPYVnwAAlZWAQuGSL3NSq8WtAQHwYZzDWXdGB+/O3mgb4XyB8fpQKvciLOwxppo7d+5kvht2\n2zbglVfY6WVUZkBv1WN059HsRKOjge++Yya3vrIST7WA2PmbBXlnbA3uMvS1pQOPHj3qmg9ITATG\njBHLBzLmqIvSEu/+aTc63MkuXwwACIIVavURBAezS+pORNi5cydmMEzNq1SKuefvukt8ziJ2fPOF\nzWyjbTIzxfw2EycykdPb7dhZXY3Hr8ltcy03Wxz9jYxs6GtwZ8SNyzdK3cCFwOtDf0bPPH5eozkB\nP7++8PZmF8KalZUFi8WCYcOGMdP84w9g6lS2aTeYu22io8XYeUaDhx0KBca3b48wRtE7Mk0jG3qI\nKWHKywEWi9tNUTuid5n/8eRJlxh6M88jRafD+EC2Sdh4A4++l/qi/e1sZwquCKusdduw3A27bRvw\n0EN/Ppd6XWRWZUJtVmNM1BhpDauFSDT0TzzBRg9iyoM5DrhtZB89O2RDDyAnB+jVi20h5oZILXNh\njhu7XfQDjGbom60hSafDoHbt4M+4kzTxGgTcEgCvANZlA/e4JO0BS7eN0QgcPsw293yt28aD1QJ0\nWpp4XTHafFdpteKERoMHWOZhlmkS2dDDff55lUmFKmMV+nTs4xr/4/nzQJcuQDD7MnwuS3twUIWs\n3llMNS2WElgsxQgIYDSqBaBUKpGWloYpU6Yw09y/X7SfddPZS70uNl/YzN5t88QTzHZy/V5ZiRkh\nIWjnwIK+7KNnh2zo4T5Dn1aehmHhw9iNtq7FhYXAj2o0LlmIVR1UwX8E252RCsUedOhwFzw82M0S\n9u7di8mTJ8PXl12FqmvdNlK5WH0RSpMS47owct0JAnO3zW8Oum1k2CIberg34qbWbeMS/6OLFmKt\ngoBErRYTGRt6a5UVpjwT7n2BXaUiAFAqd6NDh3uYau7YsYOp28ZmExdir809L+W62HJBTHnAbCCR\nmAj4+wOD2SRFyzEaUWA2Y5qDM0PZR88O2dDDfeUDU8pSXFsj1kULsad1OvT19UVQmzZMddWH1Qia\nFASPNuwuQ0GwQqU6zNQ/b7VasX//fqaG/uhRcV0oKoqZ5BX/PDM2bmTqtllXUYEnwsLg5YLQ3xuF\nnJwc+Pr6Yu7cPzdEHjp0CAMGDIC/vz+mTZuGwsJCt7er5fa4gwgCkJ0N9OvX9LlSqY24AVzgf1Qo\ngIoKl8SIxqnVmOwC/7zygBLBdwQz7QuNJqEmrLLxGO3mEBcXhwEDBiCcocshJkbcaFrfZznDxeqL\nqDZW47aut0lrWC02G7BpE/Dkk0zkBCKsq6jAM83YUXwz+uhffvlljK4TDFFdXY1Zs2Zh8eLFUCqV\nGDFiBB5nmObZUVq9oS8sFBfDXJ1Mzmgz4rLqMgaGDnTNByQmitE2jHetAqKhn8LY0BMRVAdUCL6T\n7cKxq9w2M2fOZKbH86J/fhbDxJK/n/8djw58lJ3b5uBBMQNqr15M5OI1Gvh7emL4TZ6psjGio6MR\nHByMadP+rH2wbds2DB48GA8//DC8vb2xaNEinD17FtnZ2W5tW6s39O5ciB0YOhDenuImEeb+xxMn\nXOqfv42xf96UawLZCX79/Zj2hUKxBx07svP5ExFiY2PxAMNCrgkJYsnA+myoM31BRIjOiMbjgxiO\nFNevBxhWQlpbXo5nwsObtQfhZvLRa7VafPjhh/j666+vymufkZFx1QY7Pz8/9O7dGxkZGW5tn2zo\n3WToz5SewchOI133AQkJbMsT1XBap0M/F/jnVQdUCJ4WzHTzkdlcAJutAgEB7Po5NTUVvr6+6M/w\nIomJYTuaT69Mh8lmwtgoRhFXBoO4UszIxWDkeWytrm4RdWEb4j//+Q+ee+45dOrU6arjer0e7a8Z\nJAUGBkKn07mzebKhd1fqg+SyZIyIHHHlOVP/o80GnDnjktDKIyqVa/zz+5ToMF0MIGfVF2JY5d3g\nOHbuq1q3DasbkiAAW7fW758HnOuL38//jscGPcbuphkbK84Om8hF4yg7qqsxOiAAndo2L2ldc/qC\n49g8nCEtLQ0HDx68UkegLv7+/tBqtVcd02g0bi88Ihv6ljCiT00FevZ0SUUpVyzECjYB6ji1C/zz\ne1zin2fptklKEteDWA0uiAi/Z/x+Y7ttKiowl2Fa5/ogYvNwhqNHj6KgoABdu3ZFZGQkvvzyS8TE\nxGDkyJEYPHgw0tLSrpxrMBiQm5uLQYMGMfrmDiIlmX1zH7gBiweEhxMVF7v2M/QWPfl+4ksWu8U1\nH/D110Qvvshc1szz5H/sGKmsVqa6qqMqOn3raaaaPG+mY8cCyGqtZqaZn59PISEhZLfbmWm++SbR\n++8zk6PTJaep13e9SBAENoJVVUTt2xPpdEzkSs1mCjp+nAwS+/BGtB21mEwmqqiouPJ488036dFH\nHyWFQkFVVVUUFBREW7duJbPZTG+99RaNGzeuQa2GvickFh5xQ3aXGxeVSnRHXuNWY05aeRoGhQ26\nshDLnIQEtlssazit1brEP1/XbcMKtToO7doNRps2HZlpxsbG4v7774cno0gmItE/v3UrEzkAotvm\nicFPsHPbbNoE3HuvuFGKARsqK/FQSAj8XBANdqPg4+NzVf1gf39/+Pj4oENNbouYmBi89NJLmDNn\nDsaMGYPo6Gi3t7FVu24uXBCn0AzXA+sluTQZIyOvdtsw89ETuWwh1mXx89cYehZ9UV29Ex07sguB\nBNiHVaamipl+G8ty3Jy+EEjApgub2LttGMXOExFWlZVhrpOLsDdjHD0AfPjhh1i7du2V51OnTkVm\nZiYMBgMOHz6Mrl27NvJu19DqDb07XGVnys5gRKcRTZ/oDJcvi9ajWzfm0q4w9NYqK0w5JgSOY5fu\nmIigUMQiJIRlQRAlTp8+jbtqK4IwoDbahtXAIrE4Ee3atMPgMDYpCpCbK6ZynT6didxpnQ5mQcAk\nFwwWZJpHqzb0GRnAQBftX6pLcmnydQuxzGKEa0fzjKclFkHAKZ2Oefy86oAKQZOD4OH956UntS8M\nhnPguDbw82P3Y8bGxuKOO+5Au3btmOgRiWVXH2uismFz+mL9ufV4csiT7Nw2a9cCs2eLxWsZsLKs\nDPMjI51u380UR3+j06oNvTtG9HqrHgWaAgwKddEHJSQA48czl03UatHfz++m8M+LbpsZTGPyY2Ji\n8PDDDzPTS0kR/3/rrWz0bLwNmy5swpND2LhZIAiioX/mGSZyRp7HpqoqPNOCY+dvJlq1oXfHiD6t\nPA2DwwajjefVBpOZ/9FF/vmDKhXuZJzXngSq19BL7QuFYidTt41Op8PRo0eZFgGvrcbX1L3I0b7Y\nl7sPfTv2Rc/gntIbBwDHj4sLsLewSbq3taoKYwIDEVVnkbK53Kw++huRVmvo1WpAowFcvS6SXHr1\nRimmqNVAfj4wfDhz6YMqFe5gbOj15/Tw9PeEby92Od0tlnKYTNlo355N4WoA2LVrF2677bbrdjQ6\nC5EYzMIwrTvWp6/HU0PYxbpjzRpxNM9oVrSyvBx/YxQ7b9famei0Zlqtoa+NuHF1xtQzZWfqNfRM\n/I8nTwIjRzLzqdaisdtx3mBgXh9WtU9Vr9tGSl8oFH8gOHg6PDzYha5u3boVsxjmKEhMBNq1cyyt\nuyN9obPosDtnNx4b1ITD31EMBjHLGqNNUnkmE9INBsyQWC6wti8qN1YyaFXrplUbendE3NS3EMsM\nF4ZVjgsMhA/j2GdX+OdZu21MJhP27dvHNKzSUbeNo2y7uA23d7sdIX6M6q5u2yamPIiMZCK3urwc\nT4WFoS2DURQRofTnUgatat20WkOfkeF6Q6+z6FCoKaw3NTET/6ML/fOs3TZ2jR260zoETb0+1M7Z\nvuB5E9TqI0zTHuzbtw8jRoxAaGgoEz2eBzZvdjw/mCN94TK3DQN4IqwuL8d8BjeNuLg46E7rYFfJ\nrhuptFpDf+GC6xdiU8tTMSRsyHULsUwwm4HTp28aQ6/cr0T729rDy5/dZmy1+jD8/W9FmzbsZglb\nt25lGm0THy/mBmOVT6lcX45Txacwsx+jGUdRkRgSxCifzz6lEuHe3hjGaGdtyfISdPq7i7eutwJa\nraF3x4i+sYVYyT76pCTxTsXYj15sNqPKamVeIEKxS4GO99efnsDZvqiu3s7UbWO1WvHHH3/gIYbp\nJH7/vXnZfpvqi+jz0Xig/wPwa+MnrWG1rF0LPPooICE6pi7/Ky3F3xnlFBk/eDyqt1cj4lnXJkRj\nybWlBAsKCuDh4YHAwEAEBAQgMDAQixcvdnu7WmWuG41GDFhxdcRNYnEiZvRlZ4iuIi4OmDSJuewh\ntRpTg4PhwTAmnXiCcrcS3Rd1Z6dJPKqrd+DWW99jpnno0CH0798fnTt3ZqJnswFbtohr5qxYc3YN\nltyxhI2YIAArVoghQQwoMpuRoNFgI6OpcvmqcoQ8EALvEBfliHIB15YSBACO46DRaJju82gurXJE\n766Im5PFJzGuS/1VnyT76I8eBVywc9AVbhvtaS28w73h273+sEpn+kKtPo62bbvA17eHxNb9SXR0\nNJ5gGAO5bx/Qu3fzqvE11hdp5WlQGBWY2mOq9MYBwOHDYmrrEWzCf38tK8PssDC0Y7CITwLhj6//\nQOd/sLnpuoP6SgkC4oKyIAh/UatEWq2hd7V/vlhbDJPNhF7BbGpuXoXFIrpubmNUCLoGInKJoVf8\n0bDbxlmqq2MQGsouBNJkMiE2NhaPPvooM81164Cnn2Ymh1WpqzBv+Dx4ejCKhvrlF2DBAibhQHZB\nwIqyMrzAyG2j3KeEZ4AnAka7t0CHszRUShAQR/Tdu3dH165d8eyzz0KhULi9fa3S0LvDP3+q+BTG\nRo1tcLomyUeflAT068e80MgFoxG+Hh7o6ctuQxPQtKFvbl8QCaiq2oqQEHaLprt378att96KSEYh\nhhoNsHdv07ltrqWhvrDyVmw4vwHzhs+T3DYAQFWVOOVgFDv/h0KBbj4+GMJobad0eSlm/otdZS9X\n01ApwZCQEJw+fRoFBQU4c+YMdDodnmJY1MVRWqWP/sIFYCqj2W9DnCw+iXFR7It1AxD98y5w2xxQ\nKjGN8WjeXGSGpdiCwLHsFo212iR4eQWhXTt2pcFYu222bAGmTQM6MprI7MzaiUGhg9ilPFi3Dpg5\nE2CUWfKnsjJmi7CmyyZoTmow8HfHp93cR2xuCPRh88tM1ZYSrFtJqpZ27drh1poER6Ghofjhhx8Q\nGRkJg8HALGGeQ0ipWtLcB26QKjFRUUR5ea79jAkrJtDB3IMNvn7kyBHnxadNI9q50/n3N8CdaWkU\nU1nJVLP4x2K6MOdCo+c0ty8uXXqT8vI+kNCqq9FoNBQYGEgKhYKZ5qRJRFu3Nv99DfXFvevvpTVp\nayS16QqCQNS/P9GxY0zk8oxG6nj8OBkZVeLKeS2HLr116aq+uFFsR318++235O/vT5GRkRQREUH+\n/v7k6+tLI0aMuO7c8vJy8vDwIK1WW69WQ98TEitMtTrXjVYLKJUuSd9+BStvRVp5GkZ3Ht30yc3F\nYgFOnWLun9fb7Tip1d7w/nkiQlUVW/98bGwsJk6ceKUikFQKCoDz58VCTSwo1ZXiRNEJzBrA6Duf\nOCEm4GF0Df2vtBRzIyLgy2AR1q61o3xNOTq/cvMswr7wwgvIzc1FWloazp49i7///e+4//77sW/f\nPiQlJSE7O7umZoICCxcuxJQpU+Ti4K7GHRE3Z8vPomdwTwS0bfjHdNpHf/o00Lcvsyl3LYfUaowO\nCECgFztvHm/koTmmQfD0xm8ezekLvT4NAId27YZKa1wdNm7ciNmzZzPT++03MTS9bdvmv7e+vlh3\ndh1mDZiFdt6MpvoMF2ENPI+V5eV4mVFIatmKMgTfFQyfLj43TT56Hx8fhIWFXXnUlhLs2LEj8vLy\ncPfddyMwMBBDhw6Fj48PNmzY4PY2tjofvTsibhKLEzE2aqxrxF3kn9+lUOA+Vg7lGpT7lQgYGYA2\nQex2BteO5lkt0ikUCsTHxzOr40kkur9XrWIiJ5bjS1uFlQ+sZCNYXQ3s2AF8+SUTud8qKjA+MJDJ\nAr5gF1D8XTEGbXJDEioX8uGHH1759xNPPMF07cdZHBrXchx3N8dxFzmOy+Y47l/1vP4kx3Fnax7x\nHMcNYd9UNpw751gWQSkkliQ2uRDrdBy9C+LniQi7XWDoq2OqEfpI0zljHO0LImIeVhkTE4Pp06cz\nm0onJ4v5bcY6eZ+/ti+OFRwDx3HsFvZ//RV48EFAYmZJQPw9visuxmtRUQwaBlRvr0bbzm0R5CZ+\nfQAAIABJREFUOFpcuJfz0bOjSUPPcZwHgB8ATAcwCMBsjuOuDXfIA3A7EQ0D8AmAX1g3lBWpqeyq\n/DTEyaKTrhnRW61iztuJ7HKvA8A5gwFtPTzQl2FYpWARoPhDgZCHGGVYhOi2EQQzAgLYrX389ttv\neJJRMWxAtKPz5rHLVPlj8o/4x8h/sJnB2O3A8uXAK69I1wJwQKVCG45jVle4+JtiRL3O5qYhczWO\nuG5GA8ghogIA4DguGsADAC7WnkBEiXXOTwRwQ66kCAKQlsasiE69VOgroDKr0C+kX6PnOeV/TEx0\niX++djTPMmZZdUiFdoPboW1k045qR/uiomI9wsLY1Ui9dOkSsrKycN999zHR0+vFbAIZGc5r1O2L\ncn059uXuw0/3/yS9cQCwcyfQpQuzkc53xcV4NSqKye+hTdLCUmJByIN/DgxuFh/9zYAjrpvOAIrq\nPC9G44Z8AYA9UhrlKi5fFvcYMfZQXMWpklMY03kMPDgXrPbu3QvcfTdz2V0KBe5lFHFSS9WWKofc\nNo5CxKOyciPCw9ltNlm9ejWefPJJtGFUuCU6Wkw/xCicHCtSVuDRgY+ivQ+jjXFLlwIvv8xEKtto\nxGmdDk+GhTHRK/qiCFELo+Dh1eriQ9wC017lOG4KgPkArvPj3wikprp2NA847rZxyv/oAkOvsNlw\nzmBgNv0GAMEmoDq2GiEPO+a2caQv1Opj8PYOQ7t2bFbSeZ7HmjVrMH/+fCZ6APDzz8Dzz0vTqO0L\nXuDx05mf8OLIF6U3DBDjPS9eBBhVzlpaUoLnIyOZhFQaMg1QH1Wj0/NX3yFlHz07HHHdlACom+cx\nqubYVXAcNxTAzwDuJiJVQ2Lz5s1D9+7dAQBBQUEYPnz4lSla7Q/rqufbt8dBDBN33eftObgHny34\njL1+eTnisrMBi6Wm9Wz0DymVmNSnD3w8PZm1d6htKHx7+yIxNxHIbfr8WhrTr6hYj9zcsdDr45j0\n56FDh+Dr6wulUunQ5zf1PC0NyM+PqwmpdL59aWlpmDx5Mnbl7IJ/qT80WRogUnr7sGwZ4u66Czhx\nQnL/DR4/HusrKvCz0Yi4wkLJeuGrw9H51c44fvr4Va/X7jRtjS6cuLg4rF69GgCu2EtJNLWjCoAn\ngEsAugHwBpAGYMA153QFkANgbBNaTe4ycyV33020fbvr9E02E/l/6k8as4a9+Jo1RA8/zFz2qYwM\n+rG4mKnmxecvUsGSAmZ6druJjh/vQGYzu3bOnj2bli5dykzvH/8gWrSImRxNXzed3U5YpZIoKIio\ntJSJ3H/y8mjBxYtMtIyXjXS8w3GyKq2NnvdX2w530dD3hKt3xhIRD+BlAPsBZACIJqJMjuNe4Diu\ndqL6AYAOAJZzHJfKcVyS9FsQe1ztukksTsTA0IEIbMu2GAgAMQEVY7eNTRCwV6lkGlZJPKF6WzVC\nZ7HzzyuVu+HvPwxt27JZ41er1di9ezezTVIGA7BxI/Dss0zkkKvMxZmyM+yKf//4o5jXhkHCNp3d\njuWlpXi7SxcGDQOKlhSh0/Od0CbYBVXYZK7gkI+eiPYSUT8i6kNEn9Uc+4mIfq7593NE1JGIbiWi\nW4jIBXv/pVFWJhaCYHR91suRy0cwpfsUh85tlv+R54H9+4Hp051rWAMcVqvR29cXXRhVFwIA9XE1\n2nZpC9+ejodqNtUXFRXrmS7CRkdH484770RHRje4TZuA8ePZXFtxcXFYmrQUzw5/Fj5eDH4Xkwn4\n/nvg7belawH4uawMU4OC0MdPeoUrS5kFldGViHqt/pBK2UfPjlazxF07mndl1tMj+UfYFYWoS0qK\nWHiUcUmsLVVVeJRR1EQtVb9XMR3N22xqqFQHERLCbpPUqlWrmC3CEokDZqmLsLXoLDqsPbsWr455\nlY3gmjXAqFFM8nJbBAFfFxXhHUbXYfHXxQifEw7v8JunglRDTJ48Gb6+vldKBg4YMODKa4cOHcKA\nAQPg7++PadOmobCw0O3ta3WG3lUYbUaklKVgQhfHinU3a4HJBdE2NkHA9upqzGKwQ7IWwSKgcnMl\nwp8Kb9b7GuuLyspoBAffiTZt2EQFpaamorS0FHfddRcTvRMnxCR5jELxcd7vPGb0m4HOgQzcVDwv\npjpgNJpfW16Oof7+uIXBLmJrpRVlK8vQ5c2Gp0E30yIsx3FYvnw5tFotdDodMjMzAYgpNmbNmoXF\nixdDqVRixIgReLw5RYQZIRt6RiQUJuCWyFvYJZ6qiwsMfZxajR4+PujOcDds9c5q+A/1h083Nq4g\nIkJZ2U/o1InRcBnADz/8gBdffBFejJK3ffMNsHAhwCDKEBa7BUuTluKf4/4pXQwAYmLEmSCDLJU8\nEZYUFeFdRqP5gk8LEP5UOHy6snMb/tUQXZ/LfuvWrRg8eDAefvhheHt7Y9GiRTh79iyys7Pd2rZW\nZehdmfrg8OXDDvvngWb4H1UqID2dedqDzVVVeDSUnYsFACrWViB8bvNG80DDfaHTJcNu1yA4+A6J\nLRNRKBTYunUrnnvuOSZ6ly+LOeZYheJvSN+AKGUUhoYzyMxJBHz+OfCvfzHxV26sqEB4mzaYyKCq\nmbnAjIp1Fej278Zzhd9sPvp3330XYWFhmDhxIo4ePQoAyMjIwLBhw66c4+fnh969eyNDyvZpJ2gV\nhl6jASorgT59XPcZR/IdX4htFgcPiiMyhgum9hq3zSMMDb210gr1MTVT/3xZ2c+IjHwOHKNdxitW\nrMDMmTMRyuh7L10qGnkW1fMEEvDlyS/x+CBG0/pDh8SF2BkzJEtZBQEf5udjcc+eTNId5H+cj05/\n79QifPO1LFmyBHl5eSgpKcFzzz2HmTNn4vLly9Dr9Wh/zc0xMDAQOp3Ore1rFWmK09KAIUPYTK/r\nQ2vRIqMqA+O6OJ5h0GH/486d7CpY1HBUo0E3Hx/0YOi2qdhQgZCZIfAKaP4lVV9f2O1aVFVtwahR\nmQxaJ+6EXb58ObZs2cJET6sV1zlTU5nIYe+lvfD29MYbs9+QLkYELFoEvPMOk8ILK8vK0NvXF5MY\n7J42XDRAEavA6JymA/Oa5aNnFWVRj/vFEUaNGnXl33PnzkV0dDR27doFf39/aLXaq87VaDRy4RFX\n4Gr//PGC4xjVaRSbcLi6WCyioX+YXRFsANhcWcl0NA8AFWsqEPFMBDu9ivUIDr4Dbduy0fzjjz8Q\nERGBkSNHMtFbuRK44w52gVBLEpbgzXFvsknYtnevuELMoAi1iefx34ICLO7RQ3q7AOT/Jx9d3uzC\ntEYBANFAs3gwZtCgQVfVkjUYDMjNzcUgBlFQzUE29AxwJqzSIf/jgQNi8nxWWbIgum22VVcz9c/r\nz+lhU9gQNMW5Ed+1fUFEKC39CZGRbBdhX2GUnpfnxdD0119nIofDlw+jVFeKxwc/Lt0vTQS8/z7w\n8cdMprDLSkowNjAQIwOlbwLUJmuhidc4XCbwZvHRazQa7N+/HxaLBTzPY/369Th+/DjuuecePPTQ\nQ8jIyMC2bdtgsVjw0UcfYfjw4ejbt69b29gqDH1KimsNfXMXYh1m82axJh1DjqjViGrblklFoFrK\n15YjfE44OA8202ed7jR4Xofg4GlM9DIzM5Geno5HHnmEid7GjUDnzs4XF6kLEeH9w+/jw0kfwsuD\ngSd12zbR2DOYBWrtdiwpKsJ/GYzmiQiXXr2E7h93h6efi3yofxE2mw3vv/8+wsLCEBoaimXLlmHH\njh3o1asXQkJCEBMTg/feew8dOnRAcnIys2pmzUJK/oTmPvAX5KtQKokCAoisjafScF7fqKSATwPI\nYrewFTabiYKDiRjnoXkiI4OWFhUx0+MtPCVEJJA+U89MMzNzHuXn/x8zvfnz59MiRolobDaiPn2I\nDh5kIke7s3fTwGUDyc7bpYvZ7UQDBxLt2iVdi4g+yMujuRcuMNEqW1dGySOTSeAFp97/V9iOv4KG\nvick5rpp8Yuxx4+LIy9GKcevY++lvbi92+3w9mQcQXDwoFjcllHRZUBMSbxHocByhuFHVZur4DfQ\nD+36s9k/YLGUoLp6B8aMyWGil5+fjx07duDSpUtM9DZuBCIigKkMNkATEf4T9x8smrQInh4MRrkb\nN4pFae65R7JUvsmE5SUlSGGwpmHX2ZH3rzwMihnEbNYn0zxavOvm6FGxGISr2HZxGx7q/1Cz39ek\n/9EFbpvfKipwf8eOCGZ01yMiFH9XjKiF0sq/1e2L4uJvER4+F23asMlD89lnn+GFF15AsJifWhJ2\nu+j6XrSITZBHbFYsrLwVswb+md7Bab+02Qx8+CGweDGTxr2Rm4vXoqLQlUFYb8EnBQi+MxjtxzYv\nBv9m8dHfDMiGXgJmuxn7c/djZr+ZbIWtViA2llmRCEA0yr+WlWEBgwyGtWgTtbApbOh4HxujbLOp\nUVa2El26MAgxBFBcXIxNmzbhdUarpuvXi+viUxgsxwgk4IMjH+DjyR+zqUb21VfA0KFMCscfUCpx\nVq/HmwyytBmzjShbUYae/9dTspaMBKT4fZr7gJv9bGo1kb+/6O52BbEXY2nSqknshXftIho/nqlk\nokZDvRMTSRCc85HWx/nHz1PhN4XM9PLzP6ULF55mpvfqq6/SG2+8wUTLZiPq1Yvo8GEmcrQ6dTWN\n/mU0m9+jsJCoY0eivDzJUhaep/6nTtGOqirJWoIgUNqdaVTwhfTaBO62HX8VDX1PyD76hklIEBP3\ntW26PrVTOOu2aRIXuG1+LSvD3yIimBXWNhebodqvQr+fGi+C7ig8b0ZJyfcYOvQAE73y8nKsW7eO\n2VbztWuBqCg2o3mNWYN3D72LbY9vY/N7vPmmWAuWQXTM0pISdPfxwQwGKZzLVpTBprQ1mIZYxn20\naNeNK902dsGOndk78WD/B516f4P+R6NRdNswCgUExGIRW6qq8EwEuw1NpctLET4nHF7tpY8V4uLi\nUFGxBv7+I+DvP5hB64CvvvoKTz31FCIZuKq0WjE0/fPPGTQMwMdHP8bdve/GmKgx173WbL/04cPA\nqVNiThuJFJjN+L+CAnzbu7fkG5C50IzL715G/1X9nS74XdsXJtNlSW2RaeEpEOLi2P1xXkt8YTy6\ntu+KbkGNJ2ZqNhs3AhMmiMNHRmyqqsKk9u0RyWhqw5t4lP1ahlsS2GxOEAQeRUVfol+/lUz08vPz\nsXLlSpw9e5aJ3iefiDVfxlxvl5tNZlUm1p5bi/MvnpcuZrMBr74KfP01IHFfhECEv128iH926YJ+\nEouKEBGyns9C1GtR8B8iLREQEY+sLEalu1oxLdbQ63RARgabP8762JYpzW1Tbx4PImDZMuDTT51v\n2DUIRPi2uBhf9GS3GFa+shyBYwPh10d6lSEA6N//Eioro9C+vfR0ugDw9ttvY+HChYhicLPMyhLT\nHZxnYJeJCK/ufRX/nvhvhPvXn+WzWfldvvhCDL99SLr78H+lpdDxPN5isABbtqIMtmobuvxLmtbk\nyZORn/+R5PbIoOUuxu7dSzRxomu0BUGgLl93ofMV59kKJyYS9exJxPPMJLdXVdGtp08zW4S16+2U\nEJlA2hQtEz2bTUMJCRGk1Z5honf06FHq2rUrGQwGyVqCIBaU//JLBg0jopgLMTRw2UCy2hns3ktN\nJQoNFRdiJXLJaKSOx49Tpl76pjfjJSPFh8ST7pxOspZSeYgSEiLIbC6VF2NdXRz8ZsWV/vkzZWfg\n28YXA0MHOq1Rry92+XLgxReZZBwExJv4JwUF+He3bswWYYu/L0b7ie0RcAub7HuFhUuQkzMEAQHS\niwXwPI/XXnsNn3/+OfwY1DTdtUvMOc8iRU6VoQov734ZP973I9p4NryPwSEfvcUCPP20GFIpcQQu\nEGH+xYt4r1s39G8nbdMbb+KR8UgGun3YTbLLxmIpx/r1j6F//7Vo25ZdSLCraKiUYEFBATw8PK4c\nDwwMxOLFi93evhbrujl6FPjIRbO+mAsxeKj/Q8yMJwCgulpchP36a2aSB1QqGHkeDzIqF2hT2VD8\ndTEz37zZXITS0h8RGfkjE73Vq1fDz8+PSak2vV6sHLVsGeAtcdMzEeH5P57HnKFzcHu32yW3Df/5\nj1hcYc4cyVKf1dQvXcjAzZXzSg78+vuh80vSdnMT8cjMnIOOHe9Bhw53Sm6XO6gtJVhfLWKO46DR\naNjai+YiZTrQ3AfcNP0yGIj8/IgYzESvw2wzU/gX4XShkk0OkCssWUL0zDNMJSempNBv5eXM9HLf\nyaWLCy4y07twYS7l5r7HREulUlFkZCQlJycz0VuwgGj+fCZStCp1FQ1ZPoTMNgYbOuLjiSIiiCor\nJUvtVSioU0ICFTPYaFK6spRO9T9FNp1NslZ29kJKTZ1KPP+nlrtsh7NMnjyZVqxYcd3x/Px84jiO\n7HbHchk19D0h0XXTIg39nj3M9xtdYd3ZdXTH2jvYivI8UY8eRKdOMZM8qlJRz5MnycbI328uNdPx\nDsfJVGRioqfVJlNCQgTZbGx8/U8++SS99NJLTLR27BB/Do1GutZl1WUKWRJCZ8vPShcrKyPq0oUo\nNlayVJ7RSGHx8XRMpZKspU3RUnxIPOkzpI+sioq+pVOnBpDVqrzq+M1g6MPCwig0NJRuu+02iouL\nIyLR0Ht4eFBUVBR16dKF5s+fT9XV1Q3qyIa+GcyZQ/Ttt67RHv3LaNpxcYdknSNHjvz5ZNs2ohEj\nxNU/RtyVlkY/l5Qw08v6exblvJHDRIvnzZSUNJRKS1cR0TV94QQbN26kfv36MVmALS8XB8zHj0uW\nIqvdShNXTqQl8Uscfk+DfWE2i6OXDz+U3C6D3U7DkpLoOwZZTI2XjZTQOYEqNlVI1qqq2k4JCZFk\nNF4moqv74kY39ElJSaTX68lqtdKaNWsoICCA8vLySK/X05kzZ4jneaqsrKRHHnmEpk+f3qCOqwx9\ni/PR6/ViUaavvmKvnVSShEpDJe7rcx87UbsdePddscGMfHh7FApcMpkwl9EGKfUxNapjqzEqfVTT\nJztAfv4i+Ph0R0TEM5K1iouLsXDhQuzatUvyAiwRsGCBWAf2NomRnkSEV/e8Cn9vf7wxTmLuHiJx\n52tYmOifl4BAhL9lZWGIvz9ekZgZ1Vptxbnp59D1X10R9miYJC2tNglZWQswZMhu+Pp2b/b7OUYJ\n0MjJXEHXlhLcuHEjdu/ejZdeegm33ioGGoSGhuKHH35AZGQkDAYD2klc/G4WUu4SzX3ADXfldeuI\n7r3XNdpzts6hLxK+YCv6009EU6YwG80b7HbqcfIk7VUomOjZ9XY62eskVe2QnvuEiEitjqeEhAiy\nWKSPAHmep2nTptF///tfBi0j+uwzcWJlYVBaYOmppTRw2UDSmBn4f5YtIxo0iEgrzc0lCAK9mJVF\nt6ekkNFBn3FD2PV2Sh6dTLnv5krSISLSaBIpPj6Uqqp2NniOO2wHS+655x5aunTpdcfLy8vJw8OD\ntA38lg19T8ium6u56y6iDRvY65bryinosyBSGNkYUCIi0umIIiOJTp9mJvlObi49kZHBTC/71Wy6\nMIfNwrPNpqOTJ3tRZeVWJnqfffYZjRs3jmw26QuAmzcTRUURsajJsu/SPor4MoJyldKNIMXEEIWH\nE126JFnqndxcGpmcTBqJ/WU32Cltehplzs+UvD9DrY6n+PhQqq7+o9HzbmRDr1arad++fWQ2m8lu\nt9Nvv/1G/v7+lJOTQ6dOnaKsrCwSBIGqq6vp8ccfp2nTpjWoJRt6BygtJQoKEqNuWPNx3Mf0XOxz\nzPSOHDlC9NFHRLNnM9NM1+koJD6eyhil61QdVVFCpwSyKqRv8BEEgS5eXEAXLjxz3WvO+Og3b95M\nnTt3poIC6ZkRExOJQkKIzjDYs5VekU6hS0LpWP4xp95/VV9s304UFkaUkiK5XZ/m59PAU6eoWmKp\nNavKSim3pdCFOReIt0pb6Fep4ig+PpQUin31vn6z+Oirqqpo1KhRFBgYSMHBwTRu3Dg6dOgQEYnr\nRz169CB/f3/q1KkTPfPMM1RR0fBsVjb0DvDVV0Tz5rHX1Vv01OmrTmwiJ2o4EhPDLLUsEREvCDT+\nzBn6kVHpQZvaxtRlU1j4JSUlDSab7XpXRnMNfUJCAoWGhlIKAwN4+bI4qWIQyEIppSkU/kU4bUzf\n6LTGlb7YuVM08hLDRQVBoPfz8qh3YiKVSBwAmMvMlDQsibJfzXa6JGAt5eW/UXx8CCmVhxo852Yx\n9CyRDb0DDB9OdKjh68Zp3tr/Fs3ZOoet6Pz5RIxypRMRfVlYSGPPnCGega+ft/CUOi2Vsl7KYtAy\nooqKaDpxIopMJunb9bOzsyk8PJz27NkjWevSJTHjxPffS5aipOIkCvsijLZkbJEutnWrmN5AYrit\nmefpqYwMGnvmDFVKXHjQZ+opsXciXf74siR3jSDY6dKlN+nkyR6k051z+H2yoZcNPRERpacTde4s\n1kdmSWpZKoV9EUYVeumLh1dYs4aob182gdpEFFtVRZEJCZRvkh7jLggCXZh7gc7NPEeCXfpNQ6U6\nSvHxoaTTSZ8N5eTkUI8ePejnn3+WrHX2LFGnTkT/+59kKYoviKfQJaEUe1HitEAQiD7+WFwskDiS\nV1qtNDk1lR5OT5e88Fq+vpziQ+KpdEWpJB2rtZrS0qZTaupUslobjiWvD9nQy4aeiMTB8dtvs9W0\n83Ya/cto+vXMr+xEz54lCgmhIytXMpFL1WopND6eEhndNPLez6Pk0clkN0i/Y2o0pyg+PoyUyoON\nnueI6yYxMZEiIiLop59+ktyu+HjRKxIdLU1HEARalrSMQpeE0t6cvdLE9HqiRx6hIwMGiItNEohT\nqajbiRP0z5wcsksYfdtNdrr4/EVK7JNIujRpScoqK7dRQkIk5eT886odr40hu26uO966DX1GhriY\nxiJioi4/nPqBJq6cSLzAKJukSiXWo/vtN8mbhIiISsxm6nLiBG1qZHHHUQRBoMIvC+lkr5NkqZAe\nX1hRsYni40MaDZmrpam+iI2NpZCQENq5s2mtxhAEouXLxWtlr0S7rLPoaPaW2TT0x6GUXZ0tTSwp\nSQyfnDePjuyrf2HSEcw8T29fukSRCQm0q5Hdl46gPKykU/1P0fnHzpNN43yUjsVSRRkZsykxsTep\nVM3bhSYb+uuOt15Dz/PihsFly9jqFmmKKGRJCLucNjYb0cyZRIy26eebTDQ4KYk+yc+XrMVbeLq4\n4CIlDUki42WjJC1BECg//1M6caILabWpkrQsFgstWrSIIiIi6JREf3VlJdGMGWKcfJbEpYdj+ceo\n/w/96W87/kZGq4T+0uuJXn9dDJ9cv17SXooDCgUNSUqiB9PTJfnjzaVmyngyg050PUGV2yqd9sfb\n7QYqKPiM4uNDKCfndbLbnQuFEwSBduzYIRv61m7oly8XDT3DFO5UrCmmvkv70lcnvmIjqNEQ3XMP\n0Z13MqlUflSlooiEBPqmsFByHLOl0kIpE1Po3APnyKaVFl9tNpdQevosOn16BJnN0tIvpKSk0LBh\nw+jee++lIglTNZ4n2rhR9Me/8460zVDFmmKavWU2RX0dRb+f/915IZuN6LffxIQ6c+YQSSjEfVan\no+lpadTr5EnaXFHh9PVgLjPTpbcv0fEOxyn3nVyy651z3dntJiop+R8lJHSm9PRZpNdnOqVjs9lo\nw4YNNHToUBo6dKhs6N1h6AHcDeAigGwA/2rgnO8B5ABIAzC8gXNY9gkVF4vTcIb7gyhflU+9vuvF\nbgdsfj7R4MFEL7xAVCeG2RnXjSAI9GNxMYXFx9N+iTtfBV6g8t/K6UTUCcr9d66kcDmet1Fh4Td0\n/HhHys39N9ntzRvl1u2L8vJyevPNNyksLIzWrFnjtOESBDG53fDhRKNGER1zLqydiEQD/+7Bd6nj\n5x3p/UPvk97iZPIum41o9WqiPn2Ibrut3hAxR64LQRDooFJJD6anU1h8PH1fVEQWJ0c6+gt6yvpH\nFh0PPk7Zr2STqcC5BX2jMY8uXXqb4uND6ezZe0ijSXJKp7CwkD799FPq2bMnDRkyhHbv3k2CIMiG\n3tWGHmIB8UsAugFoU2PI+19zzj0AdtX8ewyAxAa0mHWI3S56QhjkeLpCrjKXun3Tjb5L/E66mCAQ\n7dolBml/88110/JvvvmmWXIJajXdnpJCQ5KSKEfijjD1CTUlj0mm0yNOk+qY89kL7XYTlZWtoaSk\nIZSaOpUMBudSGH/zzTeUn59PL730EgUHB9NLL71EpU4uSOp0RKtWEY0bR9Svn7ix1Jl7hSAIdKLw\nBD0Z8yQFfxZMr+x+hfKUTu55SE8neust8VqYNIno8OEGG9XYdXHZaKSvCgup/6lTNDgpif5XUkI6\nJ3a5WiosVPR9ESWPTKaEyATKfTeXLOXNn+qYzSVUXPwDpaZOpePHO1JOzutkMDR/vaKoqIh+/vln\nmjZtGnXo0IFeeOEFSkxMvKovZEMvzdA7ktRsNIAcIioAAI7jogE8UDPCr+UBAGtrLPkpjuPacxwX\nTkQVDug3m8uXgblzxYIQ774rXU9r0eLbxG/x/anv8cnUT/D3kX93XowIOHxYTD6lVAKrVomVpa9B\nrVY3KSUQIVGrxeeFhUjV6/FR9+54OjwcXk5UoLJr7KjcXImKNRUwXTah56c9ET4nHJxH8xKpEREM\nhnRUVm5EWdkKBASMQM+e/4cOHe5tdmGF8vJybNu2Dd999x0+/vhjvPDCC8jMzER4eP31VBtCrRYL\nwe/cCWzdCkycCLz1FjBjBuDVjLR9ZrsZCYUJ2HZxG7Zf3A5/b388P+J5LLt3GYJ8gpohZAYSEoCD\nB4E9e8SiMk8/LT4f2HhVsrrXhYHncVqrxVGNBturq1FisWBmx474X9++uL19e4f72663Q3daB9UB\nFZT7lTDlmNBxRkf0WNwDwdOCwXk6pmOzqaDRJECjOQ61Og4mUw46drwPnTu/jA4dpsOEvKDZAAAO\n2ElEQVTT07GkcsXFxTh16hROnjyJffv2oaysDHfddReef/55zJw5Ez4+PgCAPXv2OKQn0zSO/Bl0\nBlBU53kxROPf2DklNceYGnoiYPVq4O23gXfeAV5/3fmqe0SEy+rLiLkQgy9Pfom7et2FUwtOoVeH\nXs0XEwQgJQXYt0+0NkolsGgR8PjjgKdns6QUNhtSdDrsqK7GtupqBHl5YUFkJH4fOBA+zdASLAJ0\nZ3TQnNBAE6+B+ogawVODEfXPKHS8tyM8vB3rOJ43w2jMhMFwDirVEahU++Hh4YuQkAdwyy3x8PPr\n65CO2WxGXl4ekpOTkZiYiMTERFy+fBn33XcfxowZg1WrVsHX19cBHSAzEzh3TnwkJIhF4MeNE++n\nn3wCRDpQeU5lUiFbkY0sRRZSylJwsvgkzleex+CwwXig3wM48PQBDAgd0LgIEVBeLo48zp8H0tLE\nR3o6MGQIcMcdwNKlwPjxjV4HdkFAidWKLKMRiVotXsjKQopejwsGA4b6++O29u3xXe/emNC+PTwb\nMe68iYcp1wTTJRNMOSYY0g3QndHBnG9GuyHt0OHODuj9TW8Ejg2ER5v6f38igs1WCZPpMszmyzAa\nL8JgOAe9Ph02WwUCAsYgKOh29Oz5Gdq3nwAPj/rLb1ksFhQXF6OoqAjZ2dnIzMzExYsXce7cOVit\nVowZMwZjxozBL7/8glGjRsGzmX8nNyLR0dH4+OOPUVhYiMjISKxevRoTJkzAoUOH8PLLL6OoqOjK\ntd61a1e3to0TZwWNnMBxswBMJ6Lna57PATCaiF6tc85OAP9HRCdqnh8E8DYRpVyjRR++8tZ1n3Ft\nCwii7RR4wMYDRgOg1QIajThCGzAAaHfd4IFAIIj/iYo8AMAOngg2mxVWwQIzb4XeaoDKpAKI0NGv\nA7oHd4e/d0CdxtS0iBfEhgCAYAfsPGCzAzabaHGMBsBgBAw6oI03EBIKhIcBoWGgmr9HAQBfIycQ\nwS4Q7CRg+46tuGvmAzDzBIsgwMDzUFltEAC0b+OFCG9vdPLyQfs2nPideNGugAjgAbKJ7bNbBZCF\nB1kIglmAXceD19vAGwleQR5oG94GniFt0DbCGx4+HIgEEHgAPEiwg4gHkQ0Cb4EAGwTBCLvNIP6f\n14G36eHZpj28PIPh7R0BL89O8PQMhCAIsPF2EC/AZrfDZrPBarPCYrHCYjbDaDTBaDRCrzdCrVHB\naDDB3z8QERERCA/vhPCISHTsEArAA6s2/IKnHlkAqw2wWgCrFbCYAaMZMBoBg4Gg1QJqLcFkFtCx\no4DQMEJYGI+ISAERnXh4ePCwCXbYeBtsvBUW3gqLzQSz3QyTzQSDTQ+dRQutWQOdSQuetyPCLwQh\nviHo7BuBbgFR6NwuHF7kITbAagXZbCCLGTCZAJMZpNeD9HpAq4Wg1YFXK0E+PhBCQsF3joI9qjP4\nqCjwEZ1g8faGDQSrXYBFIBh5O8yCAJ1dgN5ug94uQGOzQ2mxQ2ezIcjTE509vXFu7a947rnX0cWz\nDbp6eMPTDpBNgGDiwZvsEMwCeL0dvM4Om84Gu9YGm9ICu9IK3mZH20gvtI3yQpsoL7TpwsGnpwfa\ndAIE2MDzppqHEXa7Hna7BlarFna7GhZLNaxWBSwWBYh84OXVCZ6ekfD07AxPz67w8IgCUUdYLFYY\njUYYjUbodDpotVpotVqoVCpUVVWhqqoKlZWVUKvV6NSpE7p06YLevXtjwIAB6N+/PwYPHozu3bs7\nNBuZN28eVq9eXWs70JSt+is5cOAAnn/+eWzatAmjRo1CWVkZAMDb2xu9evXCypUrcf/99+P999/H\n8ePHcfLkyXp1GvqeNcedz2PelG8HwFgAe+s8fwfXLMgC+B+Ax+s8vwggvB4tkh/yQ37ID2ceNzLj\nx4+nlfVsgvz5559pwoQJV54bDAby9fWlrAZifJv4/i710Z8G0JvjuG4AygA8AWD2NefEAngJwO8c\nx40FoKZ6/PMk5Y4kIyPTauE4jv7qNjSEIAhITk7GzJkz0adPH1gsFjz44INYsmQJMjIyMGzYsCvn\n+vn5oXfv3sjIyEDfvvW7PV1hJ5s09ETEcxz3MoD9ECNwVhBRJsdxL4gv089EtJvjuHs5jrsEwADg\n+lLoMjIyMi2QiooK2Gw2xMTEICEhAV5eXpg5cyY++eQT6PV6hIVdXX0rMDAQOp3OrW10KCaBiPYC\n6HfNsZ+uef4yw3bJyMjIOEwcF8dEZzJNbvZ7aoMIXn311StG/Y033sAnn3yCSZMmQavVXnW+RqNB\nQECA5LY2B7fVjOU47m4A3+LPWcHn7vrsvxqO46Ighp+GQ1yf/YWIvuc4LhjA7xD3KOQDeIyINH9Z\nQ90Ex3EeAJIBFBPRzFbcD+0B/ApgMMTr4lmImxJbY1+8C2AOxBiKdIhegXb4sy8axRkDzYqgoCBE\nRUVddYzjOHAch0GDBl1ZUAYAg8GA3NxcDBo0qEndxv4uavrrWQB2AAuJaH9jWk4GJzaPmj/sHwBM\nBzAIwGyO4/q747NvEOwA3iCiQQDGAXip5vu/A+AgEfUDcBgAg10BNwULAVyo87y19sN3AHYT0QAA\nwyAGMbS6vqhZ/3sOwC1ENBTiAHQ2ru6LG5r58+dj6dKlqKqqgkqlwjfffIMZM2bgwQcfREZGBrZt\n2waLxYKPPvoIw4cPb9A/fw31Xgscxw0E8BiAARA3qy7nmgpjkrKS6+gDYuTOnjrPr4vcaU0PANsB\n3IE60UkAIgBc/Kvb5obvHgXgAIDJAGJrjrXGfggEkFvP8dbYF8E13zsYopGPrefvg25kbDYb/eMf\n/6CgoCCKjIyk1157jSw1iZUOHTpE/fv3Jz8/P5oyZUqj5S9RJ7qmoWvhWvsJYA+AMdRIH7vLdePI\npqtWAcdx3QEMB5AI8UesAAAiKuc4LqyRt7YUvgHwFoD2dY61xn7oAaCa47hVEEfzyQBeQyvsCyJS\ncRz3FYBCAEYA+4nooCt317PGy8sLy5Ytw7Jly657berUqcjMzHRGNqyBa6EzgLqB+LUbVBvELa4b\nGRGO4/wBbIHoU9NDjI+tyw0bQsYCjuPuA1BBRGkAGptqtuh+qMELwK0AlhHRrRCj1d5BK7smAIDj\nuJ4AXofoi+4EoB3HcU+hFXz3ZuJ0f7jL0JcAqLvnN6rmWKuB4zgviEZ+HRHtqDlcwXFceM3rEQAq\n/6r2uYkJAGZyHJcHYCOAqRzHrQNQ3sr6ARBntUVElFzzPAai4W9t1wQAjASQQERKIuIBbAMwHnX6\nopXS0LVQAqBLnfOatKfuMvRXNl1xHOcNcdNVrJs++0ZhJYALRPRdnWOxAObV/PsZADuufVNLgoje\nI6KuRNQT4jVwmIieBrATragfAKBmSl7EcVztqtw0ABloZddEDVkAxnIc51OzqDgN4mJ93b5ojTR0\nLcQCeILjOG+O43oA6A0gqTGhJnPdsKImvPI7/Ble+ZlbPvgGgOO4CQCOQQwbq93S/B7EH2cTxLtz\nAcTwqabTWrYAOI6bBOCfJIZXdkAr7AeO44ZBDK9sAyAPYkihJ1pnX7wF0ajxAFIBLAAQgD/7oq+7\nbNVfSU3wTBGADyEGbWxGPddCTXjl3wDY4EB4pdsMvYyMjIyzcBxHrcFWSU5e1gDyYqyMjIxMC0c2\n9DIyMjItHNnQy8jIyLRwZEMvIyMj08KRDb2MjIwMA6KjozFw4ED4+/ujT58+SEhIQEFBATw8PBAY\nGIiAgAAEBgZi8eLFbm+b27JXysjIyLRUDhw4gHffffe6UoJWqxUcx0Gj0ThczN0VyOGVMjIyNzw3\nenjlhAkTsGDBAsyff3XNpYKCAvTo0QM2m82hAuhyeKWMjIzMDUhtKcHKykr06dMHXbt2xSuvvAKL\nxQJANN7du3dH165d8eyzz0KhULi9jbKhl5GRkZHAtaUE09LSkJqaik8++QShoaE4ffo0CgoKcObM\nGeh0Ojz11FNub6PsupGRkbnhacp1ExfHxtsxeXLz7aFarUaHDh2wdu1azJkzBwCwdetWLF68GGfO\nnLnq3IqKCkRGRkKn06Fdu3bXabnKdSMvxsrIyNz0OGOgWdFQKcGG4DgOgiC4ullXIbtuZGRkZCRS\nXynB+++/H0lJScjOzgYRQaFQYOHChZgyZYrbi4PLhl5GRkZGIh988AFGjhyJvn37YtCgQRgxYgT+\n/e9/Iy8vD3fffTcCAwMxdOhQ+Pj4YMOGDW5vn+yjl5GRueG50cMrWSGHV8rIyMjIOIVs6GVkZGRa\nOLKhl5GRkWnhyIZeRkZGpoUjG3oZGRmZFo5s6GVkZGRaOLKhl5GRkWnhyIZeRkZGpoUjG3oZGRmZ\nFo5s6GVkZGQkUFsisLZcoJeXFxYuXHjl9UOHDmHAgAHw9/fHtGnTUFhY6PY2yoZeRkZGRgI6nQ5a\nrRZarRbl5eXw8/PDY489BgBQKBSYNWsWFi9eDKVSiREjRuDxxx93extlQy8jIyPDiC1btiAsLAwT\nJkwAIOalHzx4MB5++GF4e3tj0aJFOHv2LLKzs93aLtnQy8jIyDBi7dq1mDt37pXnGRkZGDZs2JXn\nfn5+6N27NzIyMtzaLtnQy8jIyDCgoKAAx44dwzPPPHPlmF6vR/v27a86LzAwEDqdzq1tkytMycjI\n3PQ0VtGpOUhJhbxu3Trcdttt6Nat25Vj/v7+0Gq1V52n0WjkwiMyMjIyzYWImDyksG7dOsybN++q\nY4MGDUJaWtqV5waDAbm5uRg0aJCkz2ousqGXkZGRkciJEydQWlqKRx555KrjDz30EDIyMrBt2zZY\nLBZ89NFHGD58OPr27evW9smGXkZGRkYia9euxaxZs9CuXburjoeEhCAmJgbvvfceOnTogOTkZERH\nR7u9fXIpQRkZmRseuZSgNOQRvYyMjEwLRzb0MjIyMi0c2dDLyMjItHBkQy8jIyPTwpENvYyMjEwL\nRzb0MjIyMi0c2dDLyMjItHDkXDcyMjI3PD4+PhUcx4X/1e1wNT4+PhWu0JU3TMnIyMi0cGTXjYyM\njEwLRzb0MjIyMi0c2dDLyMjItHBkQy8jIyPTwpENvYyMjEwL5/8BHEaM9UYxcWUAAAAASUVORK5C\nYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x110682710>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def org_f(nmin, nmax, x):\n",
    "    A = np.zeros((nmin, nmin))\n",
    "    alp = (nmax - nmin -1)/2.0\n",
    "    pvec = np.zeros(nmin)\n",
    "    qvec = np.zeros(2*nmin-1)\n",
    "    gamvec = np.zeros(nmin)\n",
    "    for l in range(1,nmin+1): #python vs matlab off by 1\n",
    "        pvec[l-1] = gammainc(alp + l, x/2.0) #gammainc argument switch python/matlab\n",
    "        gamvec[l-1] = gamma(alp + l)\n",
    "    \n",
    "    for l in range(2, 2*nmin):\n",
    "        qvec[l-1] = 2**(-(2*alp+l)) * gamma(2*alp+l) * gammainc(2*alp+l, x)\n",
    "    \n",
    "    for i in range(1, nmin+1):\n",
    "        b = pvec[i-1]**2 / 2.0\n",
    "        for j in range(i, nmin):\n",
    "            b = b - qvec[i+j-1] / (gamvec[i-1]*gamvec[j])\n",
    "            #print qvec[i+j-1] / (gamvec[i-1]*gamvec[j-1] + 1)\n",
    "            A[i-1, j] = pvec[i-1]* pvec[j] - 2*b\n",
    "    #print A\n",
    "    if nmin%2 != 0:\n",
    "        oddvec = np.zeros((nmin,1))\n",
    "        for i in range(1, nmin+1):\n",
    "            oddvec[i-1] = 2**(-(alp + nmin + 1)) / gamma(alp+nmin+1) * gammainc(alp+i, x/2.0)\n",
    "        A = np.concatenate((A, oddvec), axis=1 )\n",
    "        A = np.concatenate((A, np.zeros((1,nmin+1))), axis=0)\n",
    "    #print A\n",
    "    \n",
    "    A = A-A.T\n",
    "    return sqrt(np.linalg.det(A))*org_k(nmin, nmax)\n",
    "\n",
    "def f_fixed(nmin, nmax, x):\n",
    "    A = np.zeros((nmin, nmin), dtype=np.double)\n",
    "    alp = (nmax - nmin -1)/2.0\n",
    "    pvec = np.zeros(nmin)\n",
    "    qvec = np.zeros(2*nmin-1)\n",
    "    logqvec = np.zeros(2*nmin-1) #qvec[0] never gets used?\n",
    "    gamvec = np.zeros(nmin)\n",
    "    loggamvec = np.zeros(nmin)\n",
    "    \n",
    "    for l in range(1,nmin+1): #python vs matlab off by 1\n",
    "        pvec[l-1] = gammainc(alp + l, x/2.0)\n",
    "        loggamvec[l-1] = gammaln(alp + l)\n",
    "        gamvec[l-1] = gamma(alp + l)\n",
    "    \n",
    "    for l in range(2, 2*nmin):\n",
    "        logqvec[l-1] = (-(2*alp+l))*log(2) + gammaln(2*alp+l) + log(gammainc(2*alp+l, x))\n",
    "        qvec[l-1] = 2**(-(2*alp+l)) * gamma(2*alp+l) * gammainc(2*alp+l, x)\n",
    "    \n",
    "    for i in range(1, nmin+1):\n",
    "        b = pvec[i-1]**2 / 2.0\n",
    "        for j in range(i, nmin):\n",
    "            #b = b - qvec[i+j-1] / (gamvec[i-1]*gamvec[j-1] + 1)\n",
    "            if gamvec[i-1]*gamvec[j-1] > 100000: #ignore +1 in denom\n",
    "                #print('there')\n",
    "                b = b - exp( logqvec[i+j-1] - loggamvec[i-1] - loggamvec[j] )\n",
    "                #print exp( logqvec[i+j-1] - loggamvec[i-1] - loggamvec[j-1] ), qvec[i+j-1] / (gamvec[i-1]*gamvec[j-1] + 1.0)\n",
    "            else:\n",
    "                b = b - qvec[i+j-1] / (gamvec[i-1]*gamvec[j])\n",
    "            A[i-1, j] = pvec[i-1]* pvec[j] - 2*b\n",
    "\n",
    "    if nmin%2 != 0:\n",
    "        oddvec = np.zeros((nmin,1))\n",
    "        for i in range(1, nmin+1):\n",
    "            oddvec[i-1] = 2**(-(alp + nmin + 1)) / gamma(alp+nmin+1) * gammainc(alp+i,x/2.0)\n",
    "        A = np.concatenate((A, oddvec), axis=1 )\n",
    "        A = np.concatenate((A, np.zeros((1,nmin+1))), axis=0)\n",
    "    A = A-A.T\n",
    "    #print sqrt(np.linalg.det(A)), exp(log_k(nmin, nmax))\n",
    "    return sqrt(np.linalg.det(A))*exp(log_k(nmin, nmax))\n",
    "\n",
    "xs = np.linspace(0,100,100)\n",
    "nmaxes = [5, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70]\n",
    "for nmax in nmaxes:\n",
    "    ys = [org_f(5,nmax,x) for x in xs]\n",
    "    #print nmax, ys\n",
    "    plt.plot(xs, ys, label=nmax)\n",
    "plt.legend()\n",
    "plt.grid()      "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 310,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# function f = FfunW(nmin, nmax, x)\n",
    "# % This function calculates A defined \n",
    "# % by equation 5,6,7,8,9 (page 3)\n",
    "# % and then returns CDF of the largest eigenvalue of \n",
    "# % real Wishart matrix (equation 4 page 3) using KfunW \n",
    "# % to calculate K'.\n",
    "\n",
    "# % Algorithm : Using Algorithm 1 (page 5)\n",
    "\n",
    "# A = zeros(nmin, nmin);\n",
    "# alp = (nmax - nmin - 1)/2;\n",
    "# pvec = zeros(1,nmin);\n",
    "# qvec = zeros(1, 2*nmin-1);\n",
    "# gamvec = zeros(1, nmin);\n",
    "\n",
    "# for l = 1:nmin\n",
    "#    pvec(l) =  gammainc(x/2, alp+l);\n",
    "#    gamvec(l) = gamma(alp + l);\n",
    "# end\n",
    "# for l = 2:(2*nmin - 1)\n",
    "#    qvec(l) = (2^(-(2*alp+l)))*gamma(2*alp+l)*gammainc(x, 2*alp+l); \n",
    "# end\n",
    "\n",
    "# for i = 1:nmin\n",
    "#     b = (pvec(i)^2)/2;\n",
    "#     for j = i:(nmin-1)\n",
    "#         b = b - qvec(i+j)/(gamvec(i)*gamvec(j) + 1);\n",
    "#         A(i,j+1) = pvec(i)*pvec(j+1) - 2*b;\n",
    "#     end\n",
    "# end\n",
    "\n",
    "# if mod(nmin,2) ~= 0 \n",
    "#     oddvec = zeros(nmin,1);\n",
    "#     for i = 1:nmin\n",
    "#         oddvec(i) = (2^(-(alp+nmin+1)))/gamma(alp+nmin+1) ...\n",
    "#             * gammainc(x/2, alp + i);\n",
    "#     end\n",
    "#     A = [A,oddvec];\n",
    "#     A = [A; zeros(1,nmin+1)];\n",
    "# end\n",
    "# A = A - transpose(A);\n",
    "# f = sqrt(det(A))*KfunW(nmin,nmax);\n",
    "# end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 311,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])"
      ]
     },
     "execution_count": 311,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.zeros((1,10))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 312,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/ipykernel/__main__.py:47: RuntimeWarning: divide by zero encountered in log\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEACAYAAABI5zaHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXd4VGX2xz83kJBOCElICIRA6F2aWAGRpqgILpbVlbKK\niq6u/nYtK7a1sc0CCCgqSEepgjSR0KQGCBACgQQSUiad9DYz7++PS6gJmUlm7nsD9/M88+hMZu75\n5jI5c+f7nvccRQiBgYGBgcGNhYtsAQYGBgYGjsdI7gYGBgY3IEZyNzAwMLgBMZK7gYGBwQ2IkdwN\nDAwMbkCM5G5gYGBwA1JjclcU5VtFUdIVRTlyned8qSjKKUVRDiuK0tOxEg0MDAwM7MWWK/fvgWHV\n/VBRlBFAhBCiHTAJmOUgbQYGBgYGtaTG5C6E2AnkXucpDwE/XHjuXqCxoijNHCPPwMDAwKA2OMJz\nDwXOXXY/5cJjBgYGBgaSMBZUDQwMDG5AGjrgGClAy8vut7jw2DUoimI0sjEwMDCoBUIIxZ7n25rc\nlQu3qlgDTAaWKorSHzgvhEi/jkB79OmK9957j/fee0+2jFqjtf6sLDh6FI4cufTf48fB2xtatoQW\nLa68hYZCUBA0aaLe3N3l6gcoKi8iJjOGYxnHOJp+lGOZx0jITSAlPwWfRj609G1Jy8Yt1f/6tiTI\nKwh/D3+aejalqUdTmno2pYl7E1wbuDpff3ExxMSoJ/r0aUhMVG9nz6r/GM2bQ6tWEBwMgYHqyQ4M\nvHQLCAA/P/D1BU9PcLn0xd4R2kstFuJKSogtLuZEcTFnSkpILivjXFkZyWVluCoKLRo1okWjRjRz\nc6OpqysBF26V/+/fsCGNGzbEt0EDvBs0oKGLbebD9fRbiiyUnCmh9EwppQml6v8nlFJ6ppSylDIs\nBRZcg1xxC3a7dGvmRkP/hrg2caVhk4bqze/Sfxt4N8ClYe2NEbM5n5KSeEpKTlNScprw8LfsPkaN\nyV1RlEXAQKCpoihJwLuAGyCEEF8LIX5RFOU+RVFOA0XAeLtV1BPOnj0rW0KdcLZ+kwk2bID162HH\nDigqgu7doVs36NMHxo+Hrl2hcePaHd/Z+s1WM1GpUWw5s4V9Kfs4mnGUtII0OgZ0pGtQV7oFdWNo\nxFDa+relhW8LPFw97Dq+Q/Wnp8O+fRAdrSbzI0fURN6hg3rCO3SAESMgPFxN6M2bQ8Paf1G3R7sQ\ngrOlpfyen8/hwkJii4qILS4mpayMNh4edPL0pKOnJ3f5+dHyQjJv2agRPnXQZ4t+YREUnyqmKLqI\nwuhCCo8UUhRdREVWBe7h7ri3dse9jTserT3wG+CHe2t3GrVohKu/K4qLXRfNNlNRkU1BwSEKCw9R\nVHTsYjK3WArx8IjAw6MtHh5ta3XsGs+mEOIJG57zYq2iG9RrzGbYs0dN5uvXw5kzcO+9ak6ZOlXN\nKYpz/iYcghCC2KxYtiRsYcuZLWxL3EZL35YMbj2YJ7s/SbegbkT4R9DQxXlJx2by8mDbNtiyRb2l\npMCtt0KPHvDQQ/DOO2pCd3XVXFqZ1crBggJ+z8/n97w8fs/PB+AOX196+fgwMSSETp6eRHh44Grj\nlbYjKM8qJ297Hucjz5P1cxY7fHfgFuyGdw9vvHt4EzIxBO8e3riHuzsteV9ORUUOeXm/U1h4kMLC\nQxQUHMRszsXbuyfe3rfg53c3ISET8fBoi5tbCMoVfzz/tjueoqVNoiiKqM+2TGRkJAMHDpQto9Y4\nSn9iInz5JcydC2FhajIfMQL693dubnGU/mhTNHMOzmF57HLcGrhxb5t7Gdx6MPe0vodm3s6r4rVL\nvxAQFQUrVqjJ/Phx9QQPHqzeevWCBg2cpvVqrtaeWlbGmqwsVmZlsSsvjw6entzm68vtjRtzu68v\nrdzdr0pOzqciu4Lz289zfut5zkeepzSxlMZ3NMZvoB9HPI4w/OnhNPTV7oPaYikiL28nublbyM39\njZKSOHx9b8Xbuzc+Pr3w9r4FD48IFKXmDzxFUez23I3kbmAz+/bB//4HmzfDhAnw4ovq1Xl9IK80\nj8XHFvPtoW8xFZqY0HMCT3Z/krb+bTVPQtclLQ0WLIB586CkBB59FIYMgdtuu3YRQmPiiotZmZXF\nqqwsThYXc5+/P6MCAhjm7+9US+V6lKeXk7k8k4ylGRQeLqTx7Woy9xvoh3cvb1xctS0ILCo6Tmbm\nT+Tm/kpBwUF8fHrRpMlg/Pzuwdf3Vlxc3Gp13Nokd4QQmt3UcPWXrVu3ypZQJ2qj32wWYuVKIe68\nU4hWrYT43/+EyMtzuDSbsFe/1WoVOxJ3iKdXPi0af9JYjFk6Rqw/tV6YLWbnCKyBavWXlAixZIkQ\nI0YI4ecnxIQJQmzfLoTVqqm+qjhfUSH+l5QkWn31lQjZtUs8f/Kk2JSdLcosFk11tGrVSgA3/K1V\nq1ZV/v4Xcqdd+VYHZqKBXtm+HZ55Ri2geO01GD26TmtymnIk/QivbnyVpLwkJvWexL+G/IsgryDZ\nsq4kPx8++0z1uHr1gqefhh9/BC8v2cqIKy5mWkoKC9PTGe7vz2stWzL5tttwkfQtJzExsV5X2tmK\nI79FGraMwTWUlsKUKbBwIcyaBQ88oO+F0csxFZqY8tsU1sSt4d0B7/Js72f1sSB6OcXFMGMG/Pvf\nMHw4vPsuRETIVoVVCDbn5vJFcjJRBQU8ExLC86GhhDZqJE9TuZWMpRmE/CnkpknuVf2etbFldPau\nN5DN4cPw1FPQvr1aZRcYKFuRbZSaS/ls92f8d/d/GddzHCdfPImfu59sWVdSXg7ffgsffqgujm7d\nCl26yFaFEIJ12dn8PSEBN0Xh5RYtWNGlC+4aLthejTnPTOrXqSR/kYxnR09pOuo19vo4dblheO5S\nuZ7+igohPvpIiMBAIebP14Xdew3V6V92bJlo9Vkr8fCSh8Wp7FPairIFq1WI+fPF1uBgIYYNE2L/\nftmKLnKiqEiMiI4WHfbsEeuysoS1mn94rd77FXkVIv6NeLHDf4eIeTxG5EflCyEues43PFX9nlar\n1fDcDWrHqVOq3evpqVbftWxZ82v0QLmlnJd+eYntSduZO2ouA8MHypZ0LZmZamlRaiq88Qa8/LJs\nRQDkm818cPYsc00m3mrVihe7dsVNwxr0qxFCkLEog/i/x+M/1J8+B/vg3kpudZAeOHToEK+99lrt\nXmzvp0Fdbtwkn771iagoIYKChPjiCyE0LoCoE6n5qeL2b28Xo5aMEvml+bLlVM3GjUI0by7E668L\nUVYmW40QQgiL1Sq+S00VIbt2iQmxscKkA10FRwvEwbsPiv0994vzv5+v8jl6zx0DBgwQ7u7uwsfH\nR3h7e4uOHTvW6jiVv2dycrJ4+umnRXBwsJg5c2atrtyN5H4Tc+SIEM2aCbFihWwl9rE3ea9o8b8W\n4v3I94XFqsNPpNJSIV59VYgWLYTYskW2mosklpSI26OiRP+oKLFPVj3rZVTkVYhTfz0ldgbsFMkz\nkoXVXL0XqPfcMXDgQPHdd9/V+TiAmDJlivD39xdvvvmmyLvw71Sb5G60/LWDyMhI2RLqxOX6T5yA\nYcPgiy/g4YflabKHyMhIvj/0PSMXjWTGfTN4Z8A7uNiwu09TTpxQF0sTEtTV6Xvuufgjme+fTTk5\n9IuKYlRAALtuuYW+vr52vd7R2jN+ymBfx32Y88z0jelL6AuhKA3qSUlWNQgHVfMkJCRw6NAhPv74\nY3zt/He6HJ39ZRhowenT6qbHTz9VN0DWByosFXy590s+3fUp28Zt48EOD8qWdC3ffQd33QXPP6+2\nDWjaVLYirELwz7NnGX/iBEu7dOFvYWHSatVBLW089ZdTJLyeQJefutDx2464BdVu16beePPNNwkK\nCuKuu+5i27ZttT7OggULCAsLq7Meo879JuPsWRgwAN5+W92gVB8oNZfywOIHcGvgxsLRC/VX4gjq\nJ+WcOfDzz9Cpk2w1AORUVPBkbCwFFgvLOncmRGK9OkBZShkxY2Nw9Xel4w8dcW1ieyOi6uq/L/3c\nEQrVlj61Yf/+/XTu3Bk3NzcWL17Miy++SHR0NK1bt7brOI6sczeS+01EcrKa2P/6V7UvTH3AYrXw\n6E+P0sClAYtGL6KBi7za6yoRAt57D5YtUxt8NW8uWxEAUQUF/CEmhtEBAXzSpo2m3RirIndrLrF/\njCX0xVDC3gizuwtjTcldb4wYMYKRI0cyefJku17nyORu2DJ2UJ8997Q06N8/khdeqD+JXQjByxte\nJrskmx9G/cCO7TtkS7oSIeD112HVKrUdbw2JXav3zw8mE8OPHOFfbdrwn7ZtHZLYa6tdCEHSv5I4\n/vhxOv7QkVZvtdKkva5s9PBhZNS53wSYzeqi6ZAhao+Y+sKnOz9lR9IOto/bTqOGci2Fa7Ba1Zr1\n3bvVnab+/rIVAfBtWhrvnz3L9p496SS5R405z8yJcScoSy2j9/7euLe8MevW8/Ly2Lt3LwMGDKBh\nw4YsWbKEHTt28OWXX0rVZdgyNwFTp8KmTfDrr/WnR8y8w/N4N/Jdfp/4O8199GF1XMRigeeeU3us\n//JL7UdLOZjv0tJ49+xZfuvRg3aecrfsV2RXEH1vND79fGj3ZTtcGtXt24MeroSrIysri/vuu4+T\nJ0/SoEEDOnbsyIcffsg9l1VK2YrhuRvYzPHjqs++f786ca0+sOH0Bp5e9TSRT0fSKVAfi5MXMZth\n3Dh1EtLPP6sDYXXA92lpvKOTxF6eVU70vdH4D/WnzdQ2Dul0qOfk7kgMz10S9c1zN5vVuaX//Kea\n2OuD/v0p+3lq5VOsfHTlNYldun4h1FYCWVnqFbudid1Z+uempTHlzBm2ODGx26q9PLOc6HuiaXpf\nU4cldoPaYXjuNzD//S/4+MCkSbKV2MbpnNM8tOQh5jwwh9tb3i5bzrV8/jnExMDOneBh33BsZzHP\nZOLtM2fY0rMn7WVfsaeXc3jwYQIfDiT8g3AjsUvGsGVuUCrtmAMH6scovHJLOX2/6cuk3pN4oe8L\nsuVcS2QkPPYY7N2rmxM632TijYQEtvToQUfJi6dlpjKi74km6NEgwt8Nd/jxDVvGsGUMuGQLf/ih\nbvJQjXy0/SNaNW7F832ely3lWpKT4fHHYf583ZzQZRkZvJGQwK96SOypZRweeJigx52T2A1qh5Hc\n7UC652sj//mPWsDx7LNXPq5X/YdNh5l5YCazRs667ld5KfrLymDMGHjlFbWWtA44Sv+RwkImnzrF\nL927a1buWJ32sjQ1sQc/HUz4lHBNtBjYhuG532DExKhe+4ED9aPsscJSwbhV4/jP0P/or+QR4C9/\nURvc//3vspUAcL6igjExMXwWEUEPyZU6llILx0Ydo9kfm9HqTX18ozG4hOG530CYzXD77TBxYv1Z\nRP1g2wfsTdnL2sfX6m8Bbs4c9ZNy3z51ZVoyViEYdewYrdzdmdaunVQtQghOjD+BtcRK5yWdnf5v\nZ3juxgzVm5rp08HX91o7Rq9Em6KZvm86hyYd0l9i37cP3noLduzQRWIH+CQpiayKCn7SwdzV5M+T\nKYou4padt+jv384AMDx3u9CrZw1QWgr/+pd6oVnd35qe9FdYKhi/ejxT751KqG+oTa/RTH9GBjzy\nCHz9NXTo4LDD1kX/ppwcZqSk8GOXLlLG4V2uPWdzDuf+dY6uq7rSwEtnjdwMLmIk9xuEuXOhd2/o\n0UO2EtuYumsqzbybMa7nONlSruW55+CJJ2DUKNlKADhbUsJTsbEs6tSJUMlte4tPFxP7ZCydl3Q2\nZpxeYMaMGfTt2xd3d3cmTJhwxc+2bNlCp06d8Pb2ZvDgwSQlJWmmy/DcbwDMZmjXDhYtgttuk62m\nZo6mH+WeH+7h4LMHadlYZ9O4f/5Z7a525Ai4y09epRYLdx46xOPNmvGa5Mnl5gIzB/sfJPTFUEKf\nt+3blqPQs+e+atUqXFxc2LhxIyUlJXz33XcAZGdnExERwXfffcfIkSN5++232bFjB7t37672WIbn\nbnAFS5ao7QXqQ2I3W82MXz2eTwd/qr/EXlQEL72kLqTqILEDvHT6NG08PHi1RQupOoRVEPtULI3v\nbEzz53RY1SSRURe+4e3fv5+UlJSLj69YsYKuXbsyevRoAN577z0CAgKIi4ujffv2Ttdl2DJ2oCfP\nuhKrFT75RF37qwk96P9izxc09WzKhFsm1Pzkq3C6/g8+gDvugHvvdcrh7dW/KD2dnXl5fNuhg/RF\nyyXjlmDONtNuWjvpWuoLMTEx9LjMJ/X09KRt27bExMRoEt+4cq/nrFkDnp5Oy0cOJb8sn6m7phI5\nLlJ/CeLoUXUG6tGjspUAkF1RwaunT7OmWzd8Gsr9M83ZmEPOhhy6HO2Ci5s+rweV9x3zfhLvOs76\nKSwsJCgo6IrHfH19KSgocFiM62EkdzsYOHCgbAlXIAR8/LF61W5LrpStf9reaQyJGELnwM61er3T\n9Fut6lDrf/4TgoOdEwP79P8tPp5Hg4Lo5+vrND22YM4zc/KZkzy+6HHcmul3kLUjk7Kj8Pb2Jj8/\n/4rH8vLy8NGotFafH8MGNrFlCxQWwkMPyVZSM3mleXy+93Peufsd2VKu5fvv1VVpnWwQ2Jqby6+5\nuXxo53BlZ3D6tdP4j/DH/155k6bKy7Okxa4LXbp04fDhwxfvFxUVER8fTxeN9ikYyd0O9OBZX87H\nH8Obb4KtZc8y9X++53Pub3c/HQJqXzfuFP2ZmepXn1mzbD+RtcQW/aUWC5Pi4pjerp10OyZ7Qza5\nv+YS8e8Iae8dIQQnT06UEttWLBYLpaWlWCwWzGYzZWVlWCwWHn74YWJiYli5ciVlZWW8//779OzZ\nU5PFVDCSe71l9244c0btQqt3ckpymLZvGlPuniJbyrX8/e/wxz9Cz56ylQDwcVIS3by8eDAgQKoO\nc56ZuGfj6DCnAw195X3IpKbOpKwspeYnSuTDDz/E09OTqVOnsnDhQjw9Pfnoo48ICAhg+fLlvPXW\nW/j7+3PgwAGWLFmimS6b6twVRRkOfI76YfCtEGLqVT9vCiwAQoAGwH+FEHOrOI5R5+4gHnwQRoxQ\nrWK98/Zvb2MqNDHnwTmypVzJtm3w5JNq83sdtBg4XlTEgMOHOdynj/TNSicmnkBxVegwy3E7dO2l\nsPAY0dGDuOWWXXh5ddBtnbsj0bTOXVEUF2A6MBhIBfYrirJaCHHisqe9CBwWQoxQFCUAOKkoygIh\nhNkeMQa2ceSI2vVx2TLZSmomqziLmQdmEvVslGwpV1Jeru5E/eILXSR2qxA8e/Ik74eHS0/s2Ruy\nyd2SS98jfaVpsFhKiI19nDZt/oWnpzY2xo2GLbZMP+CUECJRCFEBLAGuXsIzAZV/IT5A9o2Y2PXi\nuX/6Kfz1r/bvs5Gh/z+//4c/dP4D4X7hdT6WQ/XPnAmtW8PDDzvumDVwPf1z0tKwAs81l7tBqDo7\nRuv3TkLC63h6diY4eJymcW8kbDHTQoFzl91PRk34l/MNsEVRlFTAG3jUMfIMriY+HjZvhtmzZSup\nmYyiDL6O+prDzx2u+claUlICU6fCunW6aHqfVlbG22fO8FuPHrhI1nP61dP43ye3OiY7ex1ZWWvo\n00eH3ULrEY5aKXkTiBZCDFIUJQLYrChKdyFE4dVPHDduHOHh4QD4+fnRs2fPi/W/lVcHer1f+ZhM\nPXPmwLhxA/Hx0b/+l756iQFiAGGNwxxyPIfpj46Gfv2IzMsDDc9HdfpnBgby55AQsg4cINKJ8Wu6\nv3rqapLXJTPp1KRrfj5w4EBN9JjNhXh6PkvnzkvZtSv6ip/fTERGRjJ37lyAi/nSXmpcUFUUpT/w\nnhBi+IX7bwDi8kVVRVF+AT4SQuy6cH8L8LoQ4sBVxzIWVOuA1QoREbBqlf67P6YVpNHlqy4cff6o\nzS19NaGkBNq2VRuE9eolWw1bc3P588mTHOvbF48G8trnmgvN7O+8nw7fdZB61R4X9yJCmOnQYdYV\nj+u5cZgj0XpA9n6graIorRRFcQMeA9Zc9ZxY4N4LIpoB7YEEe4TUB7T2Ha/m99/Bywu6d6/d67XU\n/+nOT3m6x9MOTewO0f/NN9Cnj5TEfrV+IQRvnTnDh61bS03sAOf+fY7GdzauNrFr8d4pKIgiM/Mn\n2rT52OmxbgZqtGWEEBZFUV4ENnGpFDJWUZRJ6o/F18AnwPeKokQDCvB3IUSOM4XfjCxcqJZk692G\nTMlPYf6R+RyffFy2lCspLVW99jVXX5vIYW12NkUWC49e1X9Ea0qTS0mZnkKfg32kaRDCQlzc87Rp\n8ymurvK+OdxIGP3c6wnl5RAaqpZAttL5LOJ/bPkHBeUFfDniS9lSrmTaNHU1WgfJ3SoEtxw4wD9b\nt5a+YSl2XCyNmjeizcdtpGlISZlFRsZCevbchlp9fSWGLWP0c79h2bgROnXSf2KvsFTw3eHv+O1P\nv8mWciWlpWoNqQ4SO8CyjAw8XFx4oGlTqToKDhaQuzGXfievLoDTjvLydM6efYcePbZUmdgNaodx\nJu1Apue+YIFqydQFLfSvPrma9k3b0ymwk8OPXSf9c+aoPnvv3g7TYy+V+iusVqacPctHbdpILfUT\nQhD/Wjzh74XX2GLAme+d+Pi/Exz8NN7e3ZwWw5lUN2YvMTERFxcXfH198fHxwdfXl48++kgzXcaV\nez0gP1+9cp85U7aSmpl1YBaTek+SLeNKKq/aV62SrQSAeSYTYY0aMbhJE6k6stdkU55ZTvBE57U5\nronz57dx/vxv9O0bK01DXQkNDWXKlCkXx+xdjqIo5OXlSfkQN5K7Hciqt125EgYMAP86rjM5W/+p\n7FMcST/CmE5jnHL8Wuv/9lu1MVgfeQuGoOovtVj4IDGRZZ1r19PeUVjLrcT/LZ5209rh0rDmL/DO\neO9YrRXExb1A27af07Cht8OPrxXVjdkD9duR1WqlgYRqKMOWqQdUVsnona+jvmZcz3E0aii3N8oV\nlJWpcwjffVe2EgBmpaZyi7c3/Rs3lqojdVYq7m3c8R8mrzIlOflz3N3DCAgYLU2Ds1EUhfDwcMLC\nwpgwYQLZ2dmaxTaSux3I8NzT0mD/fnjggbofy5n6y8xlzIuex7O9nTfwolb6K6/a+8prglXJ+i1b\n+DQpiX9KHsJRkVtB4oeJRPwnwubXOPq9U1p6jqSkqbRtO80xloWiOObmQAICAti/fz+JiYlERUVR\nUFDAHzW8SjNsGZ2zZAmMGgUeHrKVXJ/lscvpEdyDtv5tZUu5ROVV+4oVspUA8FNmJvd07053b7kW\nROKHiQSMDsC7qzwd8fH/R2joi3h6Ouj9osMySS8vL3pd2CwXGBjI9OnTCQkJoaioCC8vL6fHN5K7\nHcjw3BcuVNcCHYEz9c+Oms1L/V5y2vGhFvqXLYOOHXVx1Z5TUcHqFi3YXcs+IY6i+HQxpnkm+sXY\nV/royPdOfv4B8vJ20rHj9w47Zn1BURSsVqsmsQxbRsecOAGpqTBokGwl1+d45nHisuN4qIPOhrl+\n9RW8+KJsFQD8+9w5Hg4MpJ2np1QdCW8k0PK1llKHXZ858xatWr1NgwZyz4WjqG7M3r59+4iLi0MI\nQXZ2Ni+//DKDBg0yBmTrEa0994UL1TF6jlpod5b+r6O+ZuItE3Ft4OqU41dil/6DB9VPxpEjnabH\nVrIrKpidmsqQxESpOgqPFJK/K58Wr7Sw+7WOeu/k5v5GSUkCISF/dsjx9EB1Y/YSEhIYPnw4vr6+\ndO/eHXd3dxYtWqSZLsOW0SlCwKJF+p+2VFxRzPwj8zn47EHZUq7kq69g0iTHfTLWgdmpqYwKCCDI\nZJKqI/HjRFq81oIGHnLOiRCChIQ3ad36A1xcnHshoCXvvvsu71ZTjfWYxCHHRm8ZnbJ7N0yYoI73\n1HOjsLmH5/Lj8R9Z98Q62VIukZsLbdrAyZMguSlXmdVK6z172Ni9O90kLqQWnyzm0F2HuDXhVhp6\ny7mmy8xcxdmz714YwmGfaWD0lnFOy18DCdSXDpCzo2brb0fq3Llw333SEzvAkowMunp5SU3sAImf\nJBL6l1BpiV0IC2fO/IM2bT42+sdohHGW7UArz72iQrVjnnjCscd1tP5oUzTJ+cnc1+4+hx63OmzS\nb7WqlszkyU7XUxNCCP537hyvtlA9blm9iUrOlJC9NpvQF2vfW7+u2tPTF+Dq6o+/vzbvFQPDc9cl\nO3eq3R/byOvAahOzo2bzTK9naOiio7fRr7+qE01uu022En47fx6zEAyra9+IOpI0NYnmzzXH1U+O\nz221lnHmzLt06rTAmImqIYbnrkNefx3c3eH992UrqZ7C8kLCPgvT3xi9hx5SK2SeeUa2Eu47coQx\ngYFMDAmRpqE0uZQD3Q/QL64fbgFyyh+Tk78kJ2cT3buvrfUxDM/d8NxvCNavh+HDZau4Pj+f/Jn+\nLfrrK7EnJqpfexztZ9WC40VFHCwo4I+Sff9z/zlH8IRgaYndbC4gMfFj2rTRrtWtgYqR3O1AC880\nJUW99XPC7ARH6l8as5RHuzzqsOPZQo36v/4annpKtWUk83lyMs+HhuJ+WSmm1p57eXo56T+k0/K1\nlnU+Vm21Jyd/TpMmg/H21vlE9xsQHZmlBqD2bR8yRBfl2dWSX5bP1rNbmTtqrmwplygrU5uEbdsm\nWwmZ5eX8mJnJSWd8QtvBuc/OEfREEI1C5HTpLC/PIjn5C3r33isl/s2O4bnrjLFj1Sq+ceNkK6me\nBUcWsCxmGWse18fIOkDd8fX99+qMVMl8cPYs58rK+KZDB2kaKnIq2NtuL30O9cE9zF2Khvj4N7BY\n8mjfvu5TZgzP3fDc6zVms1rsMWyYbCXXZ1nMMsZ2GStbxpXMmAEvvCBbBaUWC1+lpPDXFvZv8Xck\nyV8kE/BwgLTEXlFxnrS0bwgLe0NKfC0pLy/nz3/+M+Hh4TRu3JhevXqxYcOGiz/fsmULnTp1wtvb\nm8GDB5OUlKSJLiO524GzPdO9eyEsDJxVXOEI/edLz7MtcRsPdniw7oLspFr9hw9DUpJjmt7XkYUZ\nGfTy8aEAMJvXAAAgAElEQVRzFb6/Vp67Oc9MyowUwt4Ic9gx7dWemvoVTZuOxN1d5xPdHYDZbCYs\nLIwdO3aQl5fHP//5T8aOHUtSUhLZ2dmMGTOGjz76iJycHHr37s2jj2qzVmV47jpiwwYYMUK2iuuz\n+sRq7ml9D76NfGVLuURlH5mGct/OlZuWvmgrt6d9ylcp+A/3x7OtnK6LFksxyclf0rPnFinxtcbT\n05N33nnn4v3777+f1q1bExUVRVZWFl27dmX0aHXa1HvvvUdAQABxcXG0b9/eqbqMK3c7cHY/9w0b\nnFsC6Qj9y44vY2xnOZZMlfoLCtTtvH+W32VwU24uDRSl2sHXWswDsJZbSZmWQtjfHXfVDvZpT0v7\nDl/f/nh5dXGohvpCeno6p06dokuXLsTExNCjx6VKIU9PT9q2bUtMTIzTdRhX7johIwNOndLFxspq\nyS3JZWfSTpaMWSJbyiVWrIC774bgYNlKLrYakLkLM3N5Jp4dPfHuLqeXjdVawblz/6FLl6WaxlUc\nZHmJOn4Am81mnnzyScaNG0f79u0pLCwk6Kq9Dr6+vhQUFNQpji0Yyd0OIiMjnXb1tXmzOpTDzYl7\nTeqqf+WJldzb5l58GmkzbOBqqtS/YIEudqMmlJQQVVDA6q5dq32OM98/lSR/kexQr70SW7VnZCzG\nw6MNvr63OlzD9ahrUnaIBiF48sknadSoEdOmTQPA29ub/Pz8K56Xl5enycAOw5bRCc62ZBzBsphl\nmm9cui4pKRAVpYuF1G/T0niyWbMrNi1pTf7efCrSKwh4IEBKfCGsJCVNvSkqZKpi4sSJZGVlsWLF\nChpceB906dKFw4cPX3xOUVER8fHxdOmigWUlhNDspoYzuBqLRYjAQCHOnpWtpHqyirKE7ye+orCs\nULaUS/z730JMmCBbhaiwWETIrl3iWKHccxPzRIxI+k+StPiZmavE/v29hNVqdfix9Z47Jk2aJG67\n7TZRVFR0xeOZmZnCz89PrFixQpSWloq//e1v4rbbbqv2ONX9nhcetyvfGlfuOuDgQQgIUDtB6pWV\nJ1YyLGIYXm7yt/ZfZMECtd2AZNbl5BDu7k4XiW0PylLLyPklh+CJctYehBAkJn5CWNibN13nx6Sk\nJL7++msOHz5Ms2bN8PHxwdfXl8WLFxMQEMDy5ct566238Pf358CBAyxZos2alZHc7cBZdcpaWTJ1\n0b80Zqn0jUtX6D96FHJy1MVUyXyTmsozNmxOcGade+rsVIIeD3JaW9+atJ8/H4nZnEtg4MNOia9n\nwsLCsFqtFBcXU1BQQEFBAfn5+Tz++OMA3HPPPcTGxlJUVMRvv/1GWJjj10SqwkjuOkDvfntmUSb7\nU/ZrNpTDJhYsUEdVuch9CyeXlvJ7fj5jJXZ/tJZZSZ2dSuhL8jp0JiV9SljY31EUHTdFuskwestI\nJjdXtWMyMtQe7npk9oHZRCZGsnjMYtlSVCwW9aRt3AhaLExdhw/OniWtvJyZTt6Qcj1MP5hIX5hO\nj41yOi8WFERx9OhD9O8fj4uLc5qUGb1ljN4y9Y5ff4U779RvYge5G5eqZNs2CAyUntgtQvBtWppN\nloyzEEKQ/GUyLf4ir5dNUtKntGz5qtMSu0HtMJK7HTjDM9Wy5UBt9KcXpnMw7SDD28r3jS7qnz9f\nFwupv+bmEuDqSi8ba5ad8f7J/z0fS54F/xHOHeVXnfbi4jjOn48kJORZp8Y3sB8juUtECP377ctj\nl3N/u/vxcPWQLUWluBhWrYILi1UysXUh1Zkkf5lM6EuhKC5yKlSSk78gJGQSDRvK2RFrUD02JXdF\nUYYrinJCUZQ4RVFer+Y5AxVFOaQoyjFFUbY6VqY+cPTuwmPHVDtGqz5TtdGvp/a+AwcOhDVr1DFV\nkpNqenk5W86f54lmzWx+jaPfP6XJpeRuziV4nPPLH6vSXlFxnoyMRYSGym+1bHAtNbYfUBTFBZgO\nDAZSgf2KoqwWQpy47DmNgRnAUCFEiqIocrbI1TMqr9r1WhacWZTJYdNhhkYMlS3lEgsWwJNPylbB\nPJOJhwMC8JXYiTL1q1SaPdWMhr5yNJhM3+HvP4JGjZpLiW9wfWy5cu8HnBJCJAohKoAlwENXPecJ\nYLkQIgVACJHlWJn6wNGeqdaWjL36159ez+A2g3FvqI/V3shVq9QB2A/LraUWQjCnFgupjnz/WEos\npM1JI/RFbcofr9YuhIWUlGm0aPGyJvEN7MeW5B4KnLvsfvKFxy6nPeCvKMpWRVH2K4oif7VL5xQX\nq8M5Bg2SraR61satZWS7kbJlXOK332DkSPCW6+9uO38eN0Whv6+8nvYZizPw6eeDZzs5Pduzs9fi\n6hqkeYMwA9tx1IJqQ6AXMAIYDkxRFEXuxAIn4EjPdO9e6NZN2zxlj/4KSwWbEzYzop1+pocM3LtX\nF1Uy36Sl8Uzz5nZvs3fk+yd1diqhL2i3aelq7cnJXxhX7Re43pi9xMREXFxc8PX1vdiW4KOPPtJE\nly1mXQpw+X7ZFhceu5xkIEsIUQqUKoqyHegBnL76YOPGjSM8PBwAPz8/evbsefGNU/nV72a4v2MH\ntGoVSWSkPvRcfX9n0k6aZTbjxIETBA8Mlq6HkyeJPHUKGjZE/akcPflmM+saNWJau3bSzkefJn0o\nTyvnSKMjKJGK9vH7+FNcfIKYmEBiYyM1i69XLh+z17JlS9atW8fYsWM5duwYoG5AysvLs+tiIDIy\nkrlz5wJczJd2U1NnMaABapJuBbgBh4FOVz2nI7D5wnM9gaNA5yqOVWN3NT2zdetWhx1r8GAhfv7Z\nYYezCXv0v7rhVfF+5PvOE2MvU6aIrWPGyFYhvjh3TjwRE1Or1zrq/XNy8klx5r0zDjmWrVyu/cSJ\nP4szZz7QNH59yx3du3cXK1asEGfPnhWKogiz2WzT66r7PXFGV0ghhAV4EdgExABLhBCxiqJMUhTl\n2QvPOQFsBI4Ae4CvhRDHa/dxc+NTUaHaMnfcIVtJ9aw9tZb7290vW4aKELBoEQwZIlsJ35tMTJRY\nhmkptpCxOIPgCXK6P5aXZ5GZ+RPNm0+SEr8+kJ6eTlxcHF0vDG5RFIXw8HDCwsKYMGEC2dnZmugw\nestIYN8+deTnkSOylVTNqexTDJg7gORXk3FRdLDP7dAheOQROH1aat3oscJCRhw9SmL//rhI0mGa\nZyLjxwy6r+0uJX5i4qeUlJykY8fvNY1bU2+ZSCXSIXEGioF1er3ZbGbEiBG0a9eOr776iqKiIk6e\nPEnPnj3Jzs7mhRdeoKCg4KInfzWO7C1jjNmTwPbtuuhUWy3rTq3j/nb36yOxA/z0k5rcJW8IWJCe\nzh+DgqQldoDUb1IJ+5s2LWOvxmqtIDV1Bl27rpES/3rUNSk7AlHFmD0vLy969eoFQGBgINOnTyck\nJISioiK8nNz/Xyd/vfWDq2t9a8uOHXDXXQ45lF3Yqn9t3FpGttdJCaQQ8OOP8MgjDjv/tcEqBAsz\nMnjSjh2pV1NX/UUxRZSeKcX/fuf2kamKyMhIsrJW4e7eGh+fWzSPXx+oasxeVSiKgtVqdboeI7lr\njNWq7sORkdxtIb8sn30p+xjcZrBsKSrHjkFZGfTpI1VG5PnzBLi60lVijX3anDSCxwfj0lDOn21y\n8heEhv5FSmy989xzz3HixAnWrFmDm9ulKff79u0jLi4OIQTZ2dm8/PLLDBo0yBiQrTccUZJ1/Dg0\naQLNJezYtkX/5vjN3N7ydrzddNII6jJLRmZJ3IL0dJ6qw1U71O39Yym1kL4gnZCJchZze/f2oazs\nHAEBozSPnZGRoXlMe7jemL2EhASGDx+Or68v3bt3x93dnUWLFmmiy/DcNUaWJWMra0/pyJIBNbnP\nmSNVQrHFwsqsLD5q3VqahqwVWXj38sajtZzunMnJXxIaOhkXF+1TxqxZszSPaQ+VY/aq47HHHtNQ\nzSWMK3c7cITnK3MxtSb9VmHll1O/6KcE8vhxyMuDW9Ut7rI89zVZWfTz8SGkUd2GUdRFf+rXqYQ8\nI+eqvbw8g19/XU5IyJ81j11WVsbMmTM1j3sjYCR3DRFC31fuB1IPEOAZQOsm8q5Qr2D5chgzRvqc\nVEdYMnWhOK6Y4hPFBDwop9lqWtq3NG48AFdX7RdylyxZQvfucso+6ztGnbuGnDmjblxKSZFe1Vcl\n72x9hzJzGVOHTJUtRaVHD5g+XeqnYUZ5Oe337iX5ttvwltTeN/7v8aBAxNQIzWMLYWHPngi6dl2O\nj09vjWMLbrnlFj799FNGjBhhzFA1Zqjql+3b1Tylx8QOOiuBjItTp4bffrtUGUsyMnggIEBaYreW\nWzHNMxHyZzmWTE7ORtzcAjVP7KDaWGVlZQwdqqN5AvUII7nbQV09X9mWzPX0p+SnkJiXyG0tb9NO\n0PWotGQuqxeW4bkvSE+vU2375dRGf9bqLLy6eElr7ZuaOpPmzZ+Xcu4/++wzXnnlFVwk23L1FeOs\naciOHfrdmfrLqV8YFjGMhhKqIaqksgRSIieLizlXVsZgPz9pGtK+TpO2kFpamkhe3u8EBWlf7XHq\n1Cl2797NUzpo8VxfMTx3jTCZoFMnyM6Wvj5YJQ8teYhHuzzKE92ekC0FEhLgttsgNfWKK3eteTsh\ngRKrlf9qNeT2KkoSSjh460H6n+tPA3ftz0NCwttYLAW0a/eF5rFffPFFGjdufLH3eU29ZW4UjN4y\n9ZCdO+HOO/WZ2EvNpUSejeT7h7RtBlUty5ero/QkJvbKdgMrunSRpsH0vYlmTzaTktit1nJMpm/p\n0eM3zWPn5uaycOFCYmJiNI99I6HDVKNf6uI7Vi6myqQ6/ZFnI+nerDv+HtqXulVJNZaMlr7vrrw8\nvFxc6OnAdgP26BcWgWmeieDxclr7ZmWtxsOjA15enQBtz/0333zDyJEjaS5jG/cNhJHcNUL2Yur1\n0NWs1MRE1ZYZMECqjAXp6TwVHGz3KD1Hkbs1F9cAV7y7y2kDkZo6k9DQ5zWPazabmT59Oq+88orm\nsevCU089RUhICI0bNyYiIuKKUXpbtmyhU6dOeHt7M3jwYJKSkrQRZe90j7rcqGfTVBzF+fNCeHsL\nUVYmW8m1WK1WEf55uDhiOiJbisr//ifExIlSJZSYzcJ/xw6RVFIiTUPMEzHi3JfnpMQuLIwVO3c2\nExaL9m/Yn376Sdx+++3XPK733BETEyNKLrxfTp48KZo1ayY2bNggsrKyROPGjcXy5ctFWVmZ+Nvf\n/ib69+9f7XGq+z1xxiQmg7qzaxf07QuXNYvTDadzTlNhqaBrUFfZUlR0UCWzLieHHt7etHR3lxK/\n4nwF2euyafaEnF2xaWmzCQmZgIuL9m/YadOm8dJLL2ket6507twZ9wvvFyEErq6uBAYGsmLFCrp2\n7cro0aNxc3PjvffeIzo6mri4OKdrMpK7HdTWd9RLCWRV+jfGb2RoxFBp9sMVpKTAiRNwzz1V/lgr\n33dRejp/dEK7AVv1Zy7NpMm9TXBt6upwDTVhsRRjMs0nJOTZKx7X4twfOXKEU6dOMWbMGKfHcgaT\nJ0/Gy8uLrl278o9//INevXoRExNDjx49Lj7H09OTtm3barJYbFTLaMCOHfD++7JVVM3G+I082e1J\n2TJUVqyABx6Q+hUnz2zm19xc5nToIE2Daa6JVm+3khI7I2MZvr634uERrnns6dOnM2nSJFxd7f9Q\ni4x0zMXJwIG1L7ecMWMG06dPZ9u2bTzyyCP06tWLwsJCgoKCrnier68vBQUFdZVaI0adu5MpKYGA\nAHUnvZOnatlNuaWcgH8FcOblMzT1bCpbDgwaBK++qiZ4ScwzmViRmcnqbt2kxC+KLSL6nmj6n+sv\nZShHVNSttGo1hYAAbRfYc3JyiIiIIDY2luDgayuE6lud+/PPP4+7uztCiIuLxJV069aNDz74gIcf\nfvia1xm9ZeoR+/ZBt276S+wAu5J20TGgoz4Se1YWHDwI994rVcbi9HQev+pKS0tMc000e6qZlMRe\nUHCQ8nITTZuO0Dz2d999x/33319lYq+PmM1mvLy86NKlC4cPH774eFFREfHx8XTRYP+EkdztoDa+\no55KIK/WvzF+I8MihskRczVr18KQIeBR/TAKZ/u+meXl7MnP54EA57TWrbGfvtlK+vx0abXtqamz\naN78WRTl2k1Tzjz3FouFGTNm1MuFVIDMzEyWLl1KUVERVquVjRs38uOPPzJq1CgefvhhYmJiWLly\nJWVlZbz//vv07NmT9u3bO12XkdydjMzhHDWxMX4jw9rqJLmvXKnuSpXIT5mZ3Ne0KV6Sdsbmbsql\nUVgjvDpp/zXPbM4nM/NHgoMnah573bp1BAYGcuuFoSz1DUVRmDlzJi1btqRp06ZMmTKF+fPn06dP\nHwICAli+fDlvvfUW/v7+HDhwgCVLlmijy/DcnYfFAv7+6p6cpjpwPi4nvTCdjjM6kvF/Gbg20L4q\n4wqKitShsmfPqgNmJXH3oUP8X8uWPOikK/eaiPlDDH6D/Qh9LlTz2CkpMzl//je6dPlR89hDhgzh\nT3/603WbhNU3z722GJ57PSEmBoKD9ZfYATYnbGZQ+CD5iR1g40Z1lJ7ExH6utJSYoiKG+ctpwVCR\nXUHO5hyCHtPe7xdCkJo6+5ryRy2IjY3l6NGjjB07VvPYNzpGcrcDe33H3bvV5oZ64XL9uvLbV62C\nUaNqfJozfd+lGRk8HBBAIyd2drue/vTF6TS9rymuftp/2BYUHMBiKaBJk8HVPsdZ53769Ok888wz\nNKrjfFqDazGSuxPZvVv6IKEqsQorm+I3MTRCBxNuKipg3Tp46CGpMhZnZPC4xDmppu9NBI+Ts5Ca\nlvY1ISHPoCjapoO8vDwWLVrEc889p2ncmwUjudvBwIED7Xq+3q7cK/VHm6Jp3KixPgZhb98ObdtC\naM0+s73n31biiotJLS9noJOHclSnv/BIIRUZFTQZrL0tpS6k/kRw8LjrPs8Z537u3LkMHTqUUBv+\n7Q3sx9ih6iSystQBHZ07y1ZyLfXRknEmizMyGBsYSANJLRhMc000e7oZSgPt46enL8LPbzCNGmn7\nrcFqtTJ9+nS+/14nMwRuQIwrdzuwx3fcswf69ZM6b+IaKvXrpgRSCDW521gC6QzfVwjBkowMTTYu\nVaXfWmElfWE6wU9rb8kIIUhLm03z5pNqfK6jz/3GjRvx9vbmjjvucOhxDS5hJHcnoVe/vbC8kAOp\nBxgYPlC2FIiKUrfuduwoTUJ0YSGlViu3+vpKiZ/zSw4e7TykDMAuKDiA2Zx33YVUZzF9+nReeukl\nfTSsu0Exkrsd2OM76s1vB1V/5NlI+jbvi7ebnCEQV7BypV2WjDN838UZGTwWFKRJkqlKv2muiZDx\ncgZg27OQ6shzf/r0afbt28fjjz/usGMaXIuR3J2A2Qz796ul23pj4+mN+qiSAbssGWdg1dCSqYry\nzHJyt+YS+IdAzWNfWkgdr3nsr776igkTJuBxnVYTBnXHSO52YKvvePQotGwpdU9OlURGRupnMTUu\nDnJz1SkmNuJo33d3fj7eDRrQTaOublfrz1iUQcADATT01b6uISNjsV0LqY4694WFhcybN4/nn9d+\nhJ8zqW7MXmJiIi4uLvj6+uLj44Ovr+8VI/iciVEt4wT06renFaSRV5ZHj+AeNT/Z2VRWyThx01BN\nVHaAlOX7muaaiPhPhOZxK3ektmnzqeaxFy5cyF133UV4eLjmsZ3Jm2++yTfffIO7uztxcXHcfffd\n9OnTh44dO6IoCnl5eZq/z4wrdzuw1XfUo98OkBeSx9CIobhovFmlSuz028Gxvq/ZauXHzEwe09CS\nuVx/weECKrIr8Bvk3Nr6qigoiMJsPk+TJra3V3bEuRdC1NsxejVR3Zi9yvtWq1VzTTb9lSuKMlxR\nlBOKosQpivL6dZ7XV1GUCkVRRjtOYv1Dr8l9U/wmfVgyaWlw8iQ4aVOSLWw9f54wd3faempfpQKQ\nPk8tf1RctP/WIGtH6rZt27BardxTzRjF+k5VY/ZAbfoVHh5OWFgYEyZMIDs7WxtBNU3QRv0AOA20\nAlyBw0DHap63BVgLjK7mWDWPEdcxW7durfE56elC+PkJYbE4X489lJvLhdczXsJUYJItRYhZs4R4\n4gm7X2bL+beVcbGx4n9JSQ47ni1U6reUWcTOwJ2i6FSRpvGFEKKiIk/s2OEnSkvT7HqdI8796NGj\nxYwZM2r12ppyB+CQW12xWq0iMjJSNG3aVOzbt08UFhaKqKgoYbFYREZGhnjkkUfEsGHD7P49Lzxe\nY76+/GbLR3c/4JQQIlEIUQEsAapqBPIS8BOQYesHy43I7t1qlYxEK7lK9qbsJdg7mGbe8vqnXKQW\nlowjKbVYWJ2VxaOSqmRy1ufg2cETz7baf2uwdyHVUSQlJbF161b+9Kc/OeX49ia+6m51RVEUBgwY\nwB/+8AcWL16Ml5cXvXr1wsXFhcDAQKZPn86mTZsoKipywG99fWxJQaHAucvuJ1947CKKojQHRgkh\nZgI37K4EW3xHvS6mbjy9kUfue0S2DMjLg99/h+HD7X6pozz39Tk59PT2prnGnQgr9ZvmymkSJi4s\npDZv/ozdr63ruZ81axZPPfUU3t462F+hAWazGc9qLD9FUTTx4B11ffk5cLkXf8Mm+JrQq9+umxLI\nX35RR1P5+EiTsEgPte1jta9tv7QjdYimcUtLS5kzZw6TJ0/WNK5WVDdm76GHHmLfvn3ExcUhhCA7\nO5uXX36ZQYMG4aPB+9+WUsgUIOyy+y0uPHY5fYAlilrrEwCMUBSlQgix5uqDjRs37mIZlJ+fHz17\n9rx4VVBZS6vX+59//vl19f76ayT79sGtt+pDb+X9rv26cjL7JFErolB6K3L1zJ7NwAsTdxx9/m25\nX2SxsMnNjdnt20t5/4SeDqXrg11p6NNQ8/irVr1Lo0aD6d/fxe7XX17nbm/8xMREevXqRWpqKqmp\nqbXWr1cqx+w9//zzCCFo164d8+fPp2/fvixZsoS33nqLzMxMfH19GTJkCIsWLarxmJGRkcydOxeg\n9mWjNnhQDbi0oOqGuqDa6TrP/56bdEF1/34hunbVRos9LDqySDyw6AGHLkjWiuJiIXx9hcjIqNXL\nHaF/XlqaeODIkTofpzZs3bpV7O+5X+RsydE8dnl5rtixw0+UlaXX6vW1PfdWq1X07t1brF27tlav\nr6S+5w5bqe73xBkLqkIIC/AisAmIAZYIIWIVRZmkKEpVc7lu2EGHNV1B6NZvj9/I8LbD5V8B/for\n9OoFgbWzJByhf7FES6a3X28qcirwG6h9bXt6+nz8/Yfj5la73722537v3r3k5uYyvBZrLAZ1w6Yd\nqkKIDUCHqx6bXc1zJzhAV71k924YqpO2LZUIIdgYv5Epd0+RLQVWrJDaSyazvJzdeXn81KWLlPiy\natuFEKSmzqJ9+680jQtq98fJkyfTQE+9r28SdFawp29q6q/x++/6W0w9kn4EbzdvIvwjHNYfpFZU\nVMDPP9cpuddV/4+ZmdzXtCleEhKNtdzK+rnrpfRtz8vbCVhp3PjuWh+jNufeZDKxbt06xo/XvjmZ\ngZHcHUZaGhQUQPv2spVcyYbTG/RRJbN9O7Rpo3ZUk8TijAyekFjb7t7SHY8I7TshpqbOonnz5zTv\nbfL1118zduxYmuitg95NgpHc7eB6vmNlCaTeZg9U+u0guepg5UoYXbeuFHXRn1RaSmxREUP9/euk\nobaY5poY+cpIzeOWl2eSk/MLzZrVbfOQvee+rKyMmTNn8pe//KVOcQ1qj5HcHYQe69sLywvZn7pf\n/tQlq1VN7hL99iUZGYwJDMTNRfu3fHnGhdr2R7SvbTeZ5hIQMApXV22vnpcuXUq3bt3oIml9w8BI\n7nZxPd9Rj3771jNb6Rfa7+LUJWme+7594OcHHTrU/NzrUBf9iy6095WB6QcTgQ8HsvPgTk3jCmG9\nsCP1uTofy55zL4Tg888/569//Wud41bSqlUrFEW54W+tWrVy2DkzkrsDKC+Hw4fVgdh6Qjd+uwMs\nmbpwvKiIzIoK7vLTvgRRCIHpWxPBE7VfSM3N3ULDhj74+Gj7xtyxYwfFxcUMG+a4997Zs2fr1DNm\n69atDus/48zb2bNnHXbOFOGAZjk2B1MUoWU8rdi7F557Dg4dkq3kStp+2ZYVj66ge7Pu8kQIoa4y\nL1sGt9wiRcKUM2cotlj4b9u2msfO25XHiYkn6BfbT/MFzWPHxuDvP5TmzSdpGnf06NEMGTLkhpu2\nJBNFURBC2PUGMq7cHYAe/fbTOacpriimW1A3uUKOHVPLIHv2lBJeCHFx4pIM0uakEfLnEM0Te1lZ\nKufPbyUo6AlN4yYkJLB9+3andX80sB0judtBdb6jHv32ykHYlycVKZ57pSXjgORWG/0HCgpwURR6\nS2hUZs43k7kyk+A/qZaMluc/Le1bgoIepWFDx/zetmqfNm0aEydOxEujubS2InWPhySM5F5HhICd\nO+HOO2UruZIN8RsulkBKZcUKqX57ZQdIGXNSM5Zk0GRwE9yC3DSNa7WaSUv7xiELqfaQn5/PvHnz\nbtjuj/UNw3OvI/HxMGAAnDunnxr3MnMZgf8O5MzLZ2jq2VSekPh4uOMOSEkBCbtCLULQcvdufuvR\ng44SriSj+kUR/n44TUdo+2+QlfUzSUkf06vXbk3jfvHFF/z+++8sXbpU07g3A7Xx3G3qLWNQPZVX\n7XpJ7AC7zu2iU2AnuYkdVEvmoYekJHaArbm5BLu5SUnshUcKKU8rx3+o9pumUlO/IiRE20VUi8XC\nl19+yYIFCzSNa1A9hi1jB1X5dnq0ZDaernowh+a+o4NLIO3VPy89naeDtS9BBEj7No3g8cEoDbRd\n8ygujqOg4CBBQY859Lg1aV+7di0BAQH079/foXEdheG5G9jNjh1w112yVVzJ5S0HpJGWBrGxMGiQ\nlPKanw0AACAASURBVPAFZjM/Z2VJ6SVjKbWQvjCd4Anaf7CkpEwnJOQZGjRw1zTuZ599xiuvvCJl\nbcOgagzPvQ5kZkK7dpCdLc15uIa0gjQ6f9WZzL9l0tBFous2cybs2gWSvqZ/n5bGqqwsVnfTvhQ0\nfUk6pm9N9NjcQ9O4ZnM+e/aE07fvURo1Cq35BQ7i0KFDPPDAA5w5cwZXV1fN4t5MGHXuGrNrl1oC\nqZfEDrApfhP3trlXbmIH6btS55lM8iyZC7XtWmMyzaVJk6GaJnZQF1InT55sJHadYSR3O7jat9Oj\nJbMhvvqWA5r5jjk56rZdB24/B9v1nykpIaa4mJFNtV9QLkkooSi6iIBRAdf8zJnnXwgrKSnTaNHC\nOV0Yq9NuMplYvXo1zz5b1VA2/WB47gZ2obfFVIvVwub4zfL7yaxZA4MHg6SNLD+kp/NYUJCUDpCm\n700E/TEIl0baxs7J2UCDBo3x9dV2N9306dN57LHHaCrhg9Tg+hieey0pKoJmzVTf3UP7+QtVsi9l\nH+NXjyfmhRi5QoYNg4kTYexYzUMLIWi7dy9LO3emj6+vprGtZit7wvfQfX13vLt5axo7Ono4zZo9\nQXCwdtv+8/PzadOmDXv37iUiIkKzuDcjhueuIXv3Qo8e+knsAOtPrZd/1Z6Robb4Han9YAqAnXl5\nuLu4SGk3kLsxl0ahjTRP7EVFJygsPExQ0KOaxp09ezZDhgwxErtOMZK7HVzu2+nNkgFYE7eGBzs8\nWO3PNfEdf/wR7r8fPD0dfmhb9FcupMooyUv79voLqc46/ykp02je/FlcXBo55fhwrfbS0lI+++wz\nXn/9dafFdCSG525gM3pL7ufyzpF4PpE7wySLWrwYHn9cSuhii4UVWVk82ayZ5rFLz5VyPvI8QY9p\nW1dfUXGejIxFmveRmT9/Pj169KCnpG6fBjVjeO61wGwGf384e1b9rx6YsW8Ge1P28sPDP8gTkZgI\nvXtDaiq4adssC9RpSz+YTGzooW19OUDCWwlYiiy0+6KdpnHPnfuMgoL9dO68SLOYFouFjh078u23\n33L33XdrFvdmxvDcNSI6GsLC9JPYAVafXM1DHR6SK2LJEhgzRkpiB3m17ZZSC2lz0gidrG19uRAW\nUlKmExqq7RDqFStWEBgYyF16qwM2uAIjudtBpW+nN0smrzSPPcl7GNb2+oupTvcdnWzJXE9/SlkZ\n+wsKGBVwbX25s8lcmol3L288219/ncHR5z87+xdcXQNo3Nj5/VwqtQsh+OSTT3jjjTfqVasBw3M3\nsAm9bV5af3o9d7W66+IgbCnExqp1oZJOzIL0dMYEBuKh8XZhIQTJ05Jp8VILTeMCpKR86bRNS9Wx\nefNmysrKGCmpGsrAdgzP3U6EgJAQtRTSgYPK68Tjyx9nUPggnu0tcZfgO+9AYSH873+ahxZC0GX/\nfr5u3547NR6Cnbc7j9inYrk17lYUF+2uZIuKjhMdPZj+/RNxcdHOBrvnnnsYP348Tz31lGYxDQzP\nXRPi48HVVfXc9UC5pZwNpzfwQPsH5IkQQmqVzP6CAsqtVu5o3Fjz2CnTUgidHKppYgdITv6M5s2f\n0zSx7927l/j4eB57zLHthA2cg5Hc7SAyMvKiJaMXu3Hb2W10DOhIiE/Njaqc5jtGRakJvk8f5xz/\nAtXpn2cy8ScJte1laWXkrM8heLxti7iOOv+lpclkZi4nNPRFhxzPFiIjI5k6dSr/93//Vy8bhN2M\nnrsxiclO9LaYqosqmcqrdgmfeKUWC0szMjjQu7fmsVNnpxL4aCCuftomu+Tk/xIcPAFXV+36uSQm\nJrJr1y5j0lI9wvDc7aRDB3UTZvfuspWoXnPY52FsenITnQI7yRFhtaoe1aZN0Lmz5uHnmUwsychg\nvcb/INZyK3ta7aH75u54d9VuIbu8PJN9+zrQt+8xGjVqrlnc8ePHExERwdtvv61ZTINLGDNUnUx6\nuto6pWtX2UpUDqYdxKOhBx0DOsoTsWMHNG0qJbELIZiWnMz7rVtrHjtzeSaenTw1TewAyclfEBg4\nVtPEnpSUxOrVq4mPj9cspkHdMTx3O5g9O5LbbwcJnWSrpNKSsdVrdorvqOFC6tX69xUUkGM2M1zC\nbrKUaSmEvmTfpqW6nn+zOY/U1FmEhf29Tsexl/fff5/777+fJk2aaBrXkdyMnrtO0lT94OhRHfrt\nHSX67eXl8NNPIKl6YnpKCi80b04Djb3+gqgCylLKaPqAtj3MU1Jm4u8/HA+PNprFPHnyJGvWrOHR\nR7XtOGlQdwzP3Q769oXPPtNHgj+Te4Zb59xK2mtpNHCRNOdv3Tr46CP4/XfNQ2eUl9Nh3z7ib70V\nf42rN2LHxeLVyYuw17Wrh7VYitmzpw09evyKt7d2vuDYsWPp1asXb7zxhmYxDa7FqHN3IoWFcPy4\n06v9bGbNyTWMbD9SXmIHqbXtc9LSGBMQoHliL88sJ2tVFsETte1hk5b2Hb6+/TVN7AcPHmTnzp28\n9NJLmsU0cBw2JXdFUYYrinJCUZQ4RVGuaeCsKMoTiqJEX7jtVBRF+5HzTmbvXmjTJhJ3d9lKVGpT\nAulQ37G4GNau1XTaUqV+s9XKzNRUJodq26gLIO2bNAJHB+IWYP/modqef6u1nHPn/k2rVm/W6vW1\n5R//+Af/+Mc/8PLyqveedX3XXxtqrJZRFMUFmA4MBlKB/YqirBZCnLjsaQnA3UKIPEVRhgPfAM7v\nZqQhkZHQTScfWTklORxIPcCQiCHyRCxbBrffrs4a1Jg12dmENWrELRpPW7IUW0iZlkL3jdqWXaan\nL8LDox2+vrdqFnP79u2cOHGC1atXaxbTwMEIIa57Q03S6y+7/wbw+nWe7wecq+Znor7Sq5cQ27bJ\nVqEyP3q+eHDxg3JF9OsnxJo1UkLfc+iQWGQyaR733OfnxNFRRzWNabWaxZ49HUROzhYNY1rFHXfc\nIebNm6dZTIPrcyF31pivL7/ZYsuEAucuu5984bHq+DOw3t4PGT2TlgZnzqgXqnpA+q7UqCgwmeC+\n+zQPfbyoiOPFxYwJDNQ0rqXUQtK/kmg1RdtucZmZK2jY0A8/v0GaxVy/fj25ubn88Y9/1CymgeNx\n6CYmRVEGAeOBautJxo0bR3h4OAB+fn707NmTgQMHApd8Mb3dj48fyNChMH3659L1llvK2Ry/mRn3\nzfj/9s48Pooq2+PfGwIEEgQFIWgSAsgia9h5z2UYl8FxfYw74ojzdEQHHd8bR9HRJwiDCuKORAQB\nUQQRFQV0Ipo4iUFkSUIWQiBISEJCNrLQ0El313l/VAcjsiSd7ttpqO/nU590Varr/Oqmcvr2ueee\n2+T3v/qql/S//z488AAJiYla7//VV18loXNn7r/0UtoEBWlt/+LFxWRHZVNXXcc4PNfflPaPj48n\nJ+dJbr75ZZRSWu7XMAyeeuopZs2aRWKDv2/DmLW//x892Q80/QkJCSxduhTgmL9sMqfr2mOGZb5q\nsH/CsAwwBNgN9D7FtXz+9cUXTJgg8t57IvHx8f6WIut2rZNLFl/i0Xu9ov/QIZGOHUX8EBZZt3Gj\ndEpMlAK7Xatdl90lyZHJUrW5qlnXaWr7l5Z+IT/+OEgMw9Usu01h5cqVMmrUKDEM4xfHW8Kz3xwC\nXT8ehGUa49xbAXuAHkAbIBW4+LhzotyOfexprqWhGbyL3W76spISfysxufWjW2XBlgX+E/DqqyJ3\n3OEX02/k58utGRna7Ra+XSip41O12nS5HLJ58wApLV2rzWZdXZ306dNHvv76a202LRqHJ879tGEZ\nEXEppaYCcZipk4tFZKdS6gG3wYXAM8B5wFvKnAvvEJHRnn2XaFkkJsLFF4PmEO8JKTtSRlxuHAtv\nWOgfASIQGwtvv+0H08KbhYUs7NdPq13DYbD/+f1c/IHewmzFxe/SuvX5dO6sr07/smXLiIiI4Mor\nr9Rm08J3NCrPXUS+EpF+ItJHRF5wH3vb7dgRkftFpLOIDBeRYWeKYwdzEuZ115mvG8bt/MGK9BVc\n3/d6OoV4ttpQs/UnJECrVn5ZSu+bQ4eoTUnhMs0Lchx8/yAhvUPo+J/Nt9vY9nc6a9i371kuumie\nthr1NpuNGTNmMHv27BPa9Pez31wCXb8nWDNUT0ND5+5v3k15l3tj7vWfgLfeggcf9Evd9jcKC5nQ\npYvWBTkMp0HeP/OIfiZam02A/Pw5nHvuVXTooK9G/cyZM7n00ksZO/aMmp5yVmPVljkFu3fDuHFQ\nUOD/lZdSilKYsGoCe/+6lyDlh8/koiKzrG9eHpxzjlbTqTU1/D49nT1jxhCqcQHs4uXFFC0qYth3\nw7TZtNsL2Lp1KCNHphASoqd2TWZmJuPGjSM9PZ3wcL1lFSwah1XP3cusX2+mcvvbscPPvXa/OHaA\nRYvg9tu1O3aAZ/ft44nISK2OXVxC3qw8+r7VV5tNgJ9+epoLLnhAm2MXER588EGmT59uOfYzDCss\ncwqOD8n4K25nd9r5MONDJsdMbtZ1PNbvdMLChWZIRjNbqqvZVlPDlAsu0Nr+JatLaN2lNZ2u8Gx8\n40ScTn9NTQoVFV8RFaWvAuOyZcs4evQoU6ZMOeV5gR6zDnT9nmD13E9CTY1ZLOzTT/2tBNZmryUm\nPIYenfTOjjzGunUQGQlDh2o3/cxPP/GPHj0I0dlrN4S8mXn0ntdbW4xfRMjNfYzo6GcJDtbz7ai8\nvJxp06axfv16WmlsXws9WDH3k/Dpp7Bggbk0qL8Z//547hl6DxMHT/STgPFw990waZJWs99XVXFX\nVhY5Y8bQJkjfl8yS1SXkz81n+Obh2px7efl6cnMfY+TIdIKC9PS57r//ftq1a8frr7+uxZ6F51gx\ndy/SUrJk9lftZ+uBrXx2+2f+EbBnD6SkgB+qAz7z0088Ex2t1bG7jrjY+/he+r7TV5tjNwwnubl/\np1evudoce3JyMhs2bCArK0uLPQv9WDH3EyACGzb8ui6WP+J276W9x+0Db6dd63bNvpZH+mNj4d57\n0V3IPv7QIfbb7fyxQUlhHe2fNzuPDqM7cN5V3l+X9WT6i4oW0aZNOJ076+lNOBwOpkyZwssvv0zH\nRs4bCPSYdaDr9wSr534CUlIgLAz69PGvDkMMlqQuYeXNK/0joKYGli0zBx80IiI889NPPBsdTWuN\nvfYjOUc4EHuAUWmjtNl0OqvIy5vB4MHrtX1TeP311wkPD+c2jQutWOjHirmfgJkzoaLCXC/VnyTs\nS+DhLx9mx5QdWifvHGP6dMjNheXLtZr9V0UFj+7ZQ8aoUdoWvxYRdlyzg/N+dx6Rf4vUYhMgO/te\nlGpLv36xWuzl5+czbNgwNm3aRB9/914sGo0Vc/cS69fDrFn+VmHmtv8p5k/+cewlJfDGG7B1q1az\n9b32GdHR2hw7QOmaUmoLa7nwEX1L95WVraWy8t+MHJmmxZ6IMHXqVB5++GHLsZ8FWDH34ygthZ07\n4fLLf/07nXG7KnsVn+/6nElDvJeh0iT9s2aZ2TE9e3rNfmNYV15OrWFwywkqtfmq/Z2HneT+Ty59\n5/clqLXv/iUa6q+rKyEnZwr9+y8jODjMZzYbMn/+fAoKCpg2rel59IEesw50/Z5g9dyP48sv4cor\noU3T1z/2KqsyV3FFzys4P9QP5Sj37oUPPjA/5TRiiPB/+/YxIzqaII299ryZeXQa14lOv/HehKVT\nISLk5EyhW7c/0qnTSde18SpbtmzhueeeY9OmTbRt21aLTQs/09Qawc3ZCIB67rfdJrJokb9ViIxd\nNFa+2PWFf4zfdZfI9OnazX5cUiIjtmz51UIRvuRw5mFJ6pIk9iJ9C4AUFb0nP/44SFwuPTYrKiok\nOjpa1qxZo8WehffBg3ru1oBqAxwO6NoVsrKge3f/6dh2YBs3rryRvEfzCNaU93yM1FS45hqzalqH\nDtrM1jidDN26ldi+ffnded5PQzwRIkLaFWl0mdCFiEcitNi02/PZtm0EQ4bE0aFDjM/tiQg33XQT\nvXv35hV/ZwhYeIwnA6pWzL0BiYnQq9fJHbuuuN20b6bx9GVPe92xN0r/k0/CP/6h1bED/D03l3Gd\nOp3SsXu7/UtWluCsdHLBQxd49bonIz7+W7Kz7yUi4lEtjh1g3rx5lJSU8OKLLzbrOoEesw50/Z5g\nxdwb8MYb5nwdfxKXG8f+qv3cN/w+/cYTEmDXLu2zUeMqKthQUUH6KI355dVOch/LZeDHAwkK1tPH\nKStbS8eONiIjH9diLykpiblz57Jlyxba+HsQyUI/TY3jNGejBcfcd+0SOf98kcOH/afBZbgkJjZG\nPs78WL9xwxAZPVrkgw+0mj1UVyeRycnydXm5NpuGYUjGbRmS/UC2Nps2W7YkJnYWmy1Hi72SkhKJ\niIiQdevWabFn4VvwxRqqZwvz5pkVbUND/adhRfoKQoJD+MPFf9Bv/NNPoa4O7rhDq9n/yc3l+s6d\nuUpTnB0gf24+9r12YhL1hEYMo5adO/9Iz54zaN/e9/nlLpeLSZMmMWnSJK5rCQWSLPyCFXMHDh6E\n1ath6tRTn+fLuJ3daefpb59mzlVzfDZp6aT6nU546il4/nnQON1/XVkZ/66sZE6vXo063xvtXxFX\nQcErBQz8ZCCtQnxf5lbEYOfOuwkJiSInR88i2zNmzODo0aPMnDnTa9cM9Jh1oOv3BMu5A2++aS4y\ndIJ5M9p4a8tbDOk2hMt66F98mqVLzVHk8eO1mSx3OHggJ4cl/fsTFqznC+TRvUfZefdOBqwcQEik\n7wuhiQh79jxKXV0J/fsvR2lYRWvOnDmsXLmSjz76iGBN7WrRMjnrUyFtNoiOhuRk/xUKq7RX0veN\nviRMTmDA+QP0Gi8shJEjzUHU0aO1mZ2YlUW3Nm145aKLtNhz2Vxs/8/tdP/v7trSHvPyXqCkZAUx\nMf+mdWvfT5B6+eWXWbBgAQkJCVx4ob4yCha+x6ot4wHvvmuWGvBnqY0Xkl7gxn436nfsdXVw663w\n8MNaHfua0lK21dSQOnKkFnsiwq77dhEWE8aFD+txekVFSykqepthw77X4thfe+015s+fbzl2i59p\n6ghsczZaWLaMwyESHS2yaVPjzo+Pj/e6hvyqfDnvxfOkoKrA69c+nl/pnzpV5IYbRFwun9uu52Bt\nrYR//71sqqxs8ns9bf/9L+2XLcO3iPOI06P3N5WysvWSlNRNbLZfZuP44vkREXnzzTelZ8+ekpeX\n55Pri/hOuy4CXT9WtkzTWLMGIiJg7Fj/aXg2/ln+PPzPXHiO5t7WBx+YhXS2btU2iFrjdDIhI4M/\nhYcztpGLRDSXio0V5L9kLpnXqp3vB1CrqzeTnT2ZQYM+p337fj63Fxsby9y5c0lISCAqKsrn9iwC\niKZ+GjRnowX13A1DZMQIkbVr/ach42CGdJ3bVQ4dPaTXcFqaSJcu5k9NHHY65bLt2+W+7Gxxaaod\nU7OjRpK6JUnFtxVa7Nls2ZKU1E3KyvTkli9cuFAiIyMlNzdXiz0L/4HVc288CQnmYOr11/vHvojw\nxMYnmHbJNDqF6KlGCEBlJdx8s7kSyZAhWkwecbm4IT2d3u3a8XbfvloqPh765hBZd2Zx0esXce5v\nz/W5verqzWRm3kKvXs/7fLk8EeG1115j3rx5xMfH06uRqaQWZxdnbSrk3Lnwt781LSLhzVzZWf+e\nRUF1AQ+Neshr1zwdCd9+C5MnmymPk7xXJ/5U2F0u/isjgwvbtmVRv37NcuyNbf/i94rJmpjFwNUD\n6XZHt9O/oZkcOPAO6ek30KfPm3TvfvL6Fd54fsrLy5kwYQLLly/nu+++4yJN2UaBnice6Po94azs\nuWdkmOukfvKJf+zHbo1lWdoykv6URNtgjbW1P/zQnLH10UdazNUaBn/IzKRz69Ys6dfP5ysriQh5\n/8yjeHExMQkxhF7s2+nGLpedPXsepqoqmWHDEn0eY09MTOSuu+7illtuYdWqVVZddotTclbmuU+e\nDH37mpMydbM6czWP/utREu9NpNe5Gr9Ox8WZN/7jj+Yoso+pMwxuycykjVKsHDCAYB8P2hoOg90P\n7aZmew2D1w2mbXffOj67PZ/MzJsJCelBv37vEhzsuyqaLpeL2bNnM3/+fBYvXmyVFDgLsfLcG0Fq\nKnz+ubnus2427t3I1C+nEjcpTp9jF4EFC+DZZ39OD/IxDsPgjqwsgoAPNTh2Z42TzFszUa0UMd/F\nEBzm28f60KF4du6cSETE/xIZ+ZhP17g9cOAAd911FwDbtm2zctgtGs1ZFXNPSzPXoXjnHTjXgzG2\n5sTtthRuYeKaiXx868cMDR/q8XWahM0Gd98NsbGQnEyCYfjc5NbqakZv344AqwYOpLUXHfuJ2r96\nczUpl6YQEhXCoLWDfOrYXS4b+/bNICvrTvr3X05U1N+b5Nib8vw4HA4WL17M8OHDueKKK9i4caNf\nHXugx6wDXb8nnDU99x07zHHEN94wk0V0kl2WzY0rb2TRjYv01Y7Ztcu80REj4IcfoH17s9SAjzjs\ndPLMvn18ePAgc3v3ZlK3bj7t0drz7Ox9ci+V31XS8589Cb8n3Gf2XC47RUVvs3//C3TseDkjRmwm\nJKSHT2zV1dWxbNkyZs+eTe/evVm7di1jxozxiS2LM5ym5k42Z8NPee7p6SLh4SIrV+q3nV+VLz1e\n6SFLUpboM7p6tZnH/vbbZkK/j1lXViZRyclyT1aWlNbW+tSWo8ohudNyJfG8RNn77F5x1Dh8Zsvl\nqpPCwnckOTlSduy4XqqrU3xmq7a2VmJjY6VHjx5y9dVXS1JSks9sWQQeeJDnfsY794wMke7dRVas\n0GvX6XLK+2nvS49XesicpDl6jNbViTz6qFlTYetWn5srstvltowM6b1pk2ys8O1EIZfDJQULCuT7\n8O9l5+SdYi/w3eLShuGS4uIP5IcfLpKUlCuksjLZZ7YqKytlwYIFEhUVJePHj5fkZN/ZsghcfObc\ngWuAbCAHeOIk57wO7AZSgZiTnOP7VmhAZqbp2L21uFBj6lMYhiGfZH0iA+cPlLGLxso3e7/xjvFT\nkZ8vMmOGSFSUWSvmJKsaeaO+hsswJKmyUqbm5EiXpCR5MjdXjjh9U7PFMAyp3l4tuf/IlU29N8k7\nMe9I9fZqn9hyueqkvPxr2bXrIfn++wtl27axUlGx0as26tu/sLBQFixYIOPHj5cOHTrITTfdJJsa\nW+DITwR6bZZA1++Jcz9tzF2ZRajfBK4EDgBblFJrRSS7wTm/B3qLSB+l1BggFvBjxRbIzoarr4Y5\nc2DiRO9cMzU1lXHjxp3wdyJCXG4cT8c/jcPl4IWrXuC6Ptf5Lu7sdMKGDbBwoVmv+M47zbK9MSdf\nXehU+k+FiJB6+DAflpSwqqSEsFatuLNrVzYPH06vdu2acRMnsGUI1ZurKV1TStknZRAE5998PgNW\nDGDT95voMMx7KYcul42Kiq8oK/uM8vL1tGvXhy5dJjB06EZCQ/t7zY5hGOzcuZOXXnqJadOmkZOT\nw7XXXst9993H6tWr6aB5MXJP8PTZaSkEun5PaMyA6mhgt4jkASilVgI3Yfbk67kJeA9ARDYrpToq\npbqJyEFvCz4VlZXwzTfw1Vfw2Wfm0nnenIhZWVn5i30RobCmkK0HtvLKD69w8PBBnvvtc9wy4BaC\nfLEwQ1kZZGbC11/DkiVmIfr774dVqxq1PuDx+k+G0zDItdvJstnYVlPDx6Wl1IpwZ9eufDF4MIND\nQ73yoSUi1B2ow5ZpM7cMGxX/qiD4nGDOv/l8Bn06iNAhP9uq2lDlsa26ujJstgxstvRf/DznnDF0\n6TKBXr1eoG3b5mejVFVVkZ6eTlpaGjt27CAtLY2MjAy6du1K586def7557n88ssDbsHqxj47LZVA\n1+8JjXHuFwL5DfYLMB3+qc4pdB/zqXM3DDNv/auvzAKHqalw6aVmuuPjj3uvRrshBrY6GyW2Et7f\n8T6pxanHtlZBrYgJj2Hy0MncPfRugoM8TEASMVMXq6uhqgpKSyEry3Tm9VttLQwcaJaxjIszXzcR\np2FQ6XRyqMFW4XCw++hRsmw2so4cYffRo1zQpg0DQkMZFBrK0v79GXPOOU1y6EadgaPcYW5l5uYs\nd+Ioc2DfZzedeZaNoLZBhA4MJXRgKB1GdiDyschGzSwVEQyjFsM4gst1GIejjLq6EhyOkl/8rK0t\nwGbLwDDshIYOIixsMGFhQ+jWbSJhYTEEB5+8OqWIYLfbOXz48C+20tJSioqKKC4upqio6NjrAwcO\ncOjQIQYNGsSQIUMYOnQokyZNYsiQIXTs2JHp06dz1VVXNboNLSyaw2lnqCqlbgbGi8if3fuTgNEi\n8kiDc74AnheRZPf+RuBxEdl+3LVk1EA9pV59wd6CI/SKaG/ueNpxbewEXXWcEQ97yg3N7d1/mF6R\nYceuX39F1XBT6te3drxm+fmF1O9Lg+P1z1SQgiC39CBQQaZR1UqhgoFWCvMLjrjfJu63Gscu6A43\nIiLk5BylT59WiLgQMdzL1gUBrVCqDUFBbVCqDWD+/HlrD7TGMAz3h4KBYRg4nU4cDsexra6u7tjr\neqfepk0bQkNDCQsLIywsjNDQULp06UL37t2PbeHh4cdeR0ZG0qrViUsLT548maVLl57ir9VyCWTt\nEPj6PZmh2hjnPhaYLiLXuPenYQb3X2xwTiwQLyKr3PvZwG+OD8sopfxfe8DCwsIiAGmqc29MDGEL\ncJFSqgdQBNwB3HncOZ8DfwFWuT8MKk8Ub2+qOAsLCwsLzzitcxcRl1JqKhCH+R14sYjsVEo9YP5a\nForIBqXUtUqpPYANOHndUwsLCwsLn6O1KqSFhYWFhR60FQ5TSl2jlMpWSuUopZ7QZddTlFKLlVIH\nlVI7Ghw7VykVp5TapZT6l1KqRY4OK6UilFLfKqUylVLpSqlH3McDRX9bpdRmpVSK+x5mu48H6krC\nrwAAA2tJREFUhP56lFJBSqntSqnP3fsBo18ptU8pleb+G/zoPhZI+jsqpVYrpXa6n6ExgaBfKdXX\n3ebb3T+rlFKPeKJdi3NvMBFqPDAQuFMp5b1ZIr5hCabehkwDNopIP+Bb4EntqhqHE/hfERkI/Afw\nF3d7B4R+EakFfisiw4AhwBVKqUsIEP0N+CuQ1WA/kPQbwDgRGSYi9anPgaT/NWCDiFwMDMWcl9Pi\n9YtIjrvNhwMjMMPcn+KJ9qZOafVkw5yt+mWD/WmcpIxBS9qAHsCOBvvZQDf363Ag298aG3kfnwFX\nBaJ+oD3wIzAgkPQDEcDXwDjg80B7foCfgM7HHQsI/cA5QO4JjgeE/gZ6fwckeqpdV1jmRBOhAnHV\nga7izgISkWKgq5/1nBalVDQQA/yA+XAEhH53SCMFKAYSRCSLANIPvAL8nV/OEggk/QJ8rZTaopS6\nz30sUPT3BMqUUkvc4Y2FypzsECj667kdWOF+3WTtZ9ViHT6gRY9GK6XCgI+Bv4rIYU4xHamlISKG\nmGGZCOAypdQ4AkS/Uuo64KCIpHLq6W4tUr+bS8QMDVyLGda7jABpf8wswOHAfPc92DCjBYGiH6VU\na+BGYLX7UJO163LuhUBUg/0I97FA46BSqhuAUiocKPGznpOilArGdOzLRWSt+3DA6K9HRKqBDcBI\nAkf/JcCNSqm9wIeYYwbLgeIA0Y+IFLl/lmKG9UYTOO1fAOSLyFb3/hpMZx8o+gF+D2wTkTL3fpO1\n63LuxyZCKXN++B2YE59aOvWz8uv5HJjsfn0PsPb4N7Qg3gWyROS1BscCQr9Sqkt9NoBSqh1wNZBC\ngOgXkadEJEpEemE+69+KyN3AFwSAfqVUe/e3PpRSoZix33QCp/0PAvlKqb7uQ1cCmQSIfjd3YnYM\n6mm6do2DA9cAuzBrvk/z92BFI/SuwCxxXAvsx5yYdS6w0X0fcUAnf+s8ifZLABdmbf0UYLu7/c8L\nEP2D3ZpTgDTgMffxgNB/3L38hp8HVANCP2bMuv7ZSa//fw0U/W6tQzE7lanAJ0DHQNGPmURQCnRo\ncKzJ2q1JTBYWFhZnINaAqoWFhcUZiOXcLSwsLM5ALOduYWFhcQZiOXcLCwuLMxDLuVtYWFicgVjO\n3cLCwuIMxHLuFhYWFmcglnO3sLCwOAP5fyJbRiQtlAE0AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1106ba650>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#every thing with mpmath\n",
    "#mp.det seems to have a bug :(\n",
    "# import mpmath as mp\n",
    "# mp.dps = 1500\n",
    "\n",
    "# def mp_t(m, a):\n",
    "#     t = mp.pi**(m*(m-1)/4.0)\n",
    "#     for i in range(1, m+1): #matlab is inclusive\n",
    "#         t = t*mp.gamma(a - (i-1)/2.0)\n",
    "#     return t\n",
    "\n",
    "# def mp_k(nmin, nmax):\n",
    "#     nmat = nmin if nmin%2==0 else nmin + 1\n",
    "#     alp = (nmax - nmin - 1)/2.0\n",
    "#     K = (mp.pi**(nmin**2/2.0)/2**(nmin*nmax/2.0)) #that's gonna overflow real quick\n",
    "#     K = K/mp_t(nmin, nmax/2.0) #I can see a lot of gamma cancellation here too\n",
    "#     K = K/mp_t(nmin, nmin/2.0)\n",
    "#     Kd = K*2**(alp*nmat+nmat*(nmat+1)/2.0)\n",
    "#     for k in range(1, nmat+1):\n",
    "#         Kd = Kd*mp.gamma(alp+k)\n",
    "#     return Kd\n",
    "\n",
    "# def mp_log_t(m, a):\n",
    "#     lt = mp.log(mp.pi)*(m*(m-1)/4.0)\n",
    "#     #need to take care of the sign?\n",
    "#     return lt + sum(mp.loggamma(a - (i-1)/2.0) for i in range(1, m+1))\n",
    "\n",
    "# def mp_log_k(nmin, nmax):\n",
    "#     nmat = nmin if nmin%2==0 else nmin + 1\n",
    "#     alp = (nmax - nmin - 1)/2.0\n",
    "#     logK = nmin**2/2.0 * mp.log(mp.pi) - nmin*nmax/2.0*mp.log(2.0)\n",
    "#     logK = logK - mp_log_t(nmin, nmax/2.0)\n",
    "#     logK = logK - mp_log_t(nmin, nmin/2.0)\n",
    "#     logKd = logK + (alp*nmat+nmat*(nmat+1)/2.0)*mp.log(2.0)\n",
    "#     for k in range(1, nmat+1):\n",
    "#         logKd = logKd + mp.loggamma(alp+k)\n",
    "#     return logKd\n",
    "\n",
    "# def mp_f(nmin, nmax, x):\n",
    "    \n",
    "#     A = mp.matrix(nmin, nmin) if nmin%2==0 else mp.matrix(nmin+1, nmin+1)\n",
    "\n",
    "#     alp = (nmax - nmin -1)/2.0\n",
    "#     pvec = mp.matrix(1, nmin)\n",
    "#     qvec = mp.matrix(1, 2*nmin-1)\n",
    "#     logqvec = mp.matrix(1, 2*nmin-1)\n",
    "#     gamvec = mp.matrix(1, nmin)\n",
    "#     for l in range(1,nmin+1): #python vs matlab off by 1\n",
    "#         pvec[l-1] = mp.gammainc(alp + l, 0, x/2.0, regularized=True) #gammainc argument switch python/matlab\n",
    "#         gamvec[l-1] = mp.loggamma(alp + l)\n",
    "    \n",
    "#     for l in range(2, 2*nmin):\n",
    "#         #qvec[l-1] = 2**(-(2*alp+l)) * mp.gamma(2*alp+l) * mp.gammainc(2*alp+l, 0, x, regularized=True)\n",
    "#         logqvec[l-1] = -(2*alp+l)*mp.log(2) + mp.loggamma(2*alp +l) + mp.log(mp.gammainc(2*alp+l, 0, x, regularized=True))\n",
    "    \n",
    "#     #print pvec\n",
    "#     for i in range(1, nmin+1):\n",
    "#         b = pvec[i-1]**2 / 2.0\n",
    "#         for j in range(i, nmin):\n",
    "#             b = b - mp.exp(logqvec[i+j-1] - gamvec[i-1] - gamvec[j])\n",
    "#             #print b\n",
    "#             A[i-1, j] = pvec[i-1]* pvec[j] - 2*b\n",
    "#     #print A\n",
    "#     if nmin%2 != 0:\n",
    "#         for i in range(1, nmin+1):\n",
    "#             A[i-1, nmin] = 2**(-(alp + nmin + 1)) / mp.gamma(alp+nmin+1) * mp.gammainc(alp+i, 0, x/2.0, regularized=True)\n",
    "#     #print A\n",
    "#     A = A-A.T\n",
    "#     #print mp.sqrt(mp.det(A)), mp.exp(mp_log_k(nmin, nmax))\n",
    "    \n",
    "#     return mp.sqrt(mp.det(A))*mp.exp(mp_log_k(nmin, nmax))\n",
    "\n",
    "# org_v = org_f(5, 5, 10)\n",
    "\n",
    "# mp_v = mp_f(5, 5, 10)\n",
    "# xs = np.linspace(0,70,50)\n",
    "# nmaxes = [5, 10, 15, 20, 25, 30, 35]\n",
    "# for nmax in nmaxes:\n",
    "#     ys = [f_fixed(5,nmax,x) for x in xs]\n",
    "#     #print nmax, ys\n",
    "#     plt.plot(xs, ys, label=nmax)\n",
    "# plt.legend()\n",
    "# plt.grid()\n",
    "#print mp_v, 'mpv'\n",
    "#print 'compare', org_v, mp_v\n",
    "# a = 10\n",
    "# ms = np.arange(10)\n",
    "# org_ts = [org_t(m,a) for m in ms]\n",
    "# log_ts = [log_t(m,a) for m in ms]\n",
    "# mp_ts = [mp_t(m,a) for m in ms]\n",
    "\n",
    "# print org_k(10,50), mp_k(10,50)\n",
    "\n",
    "# print org_ts, mp_ts"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 252,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.8092513943330656e+75"
      ]
     },
     "execution_count": 252,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "2**250.0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}