{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Elastic network analysis\n", "\n", "Here we use a Gaussian network model to characterise conformational states of a trajectory.\n", "\n", "**Last executed:** Feb 06, 2020 with MDAnalysis 0.20.2-dev0\n", "\n", "**Last updated:** January 2020\n", "\n", "**Minimum version of MDAnalysis:** 0.17.0\n", "\n", "**Packages required:**\n", " \n", "* MDAnalysis (Michaud-Agrawal *et al.*, 2011, Gowers *et al.*, 2016)\n", "* MDAnalysisTests\n", "\n", "**Optional packages for visualisation:**\n", "\n", "* [matplotlib](https://matplotlib.org)\n", "\n", "
\n", " \n", "**Note**\n", "\n", "The elastic network analysis follows the approach of (Hall *et al.*, 2007). Please cite them when using the ``MDAnalysis.analysis.gnm`` module in published work.\n", "\n", "
" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import MDAnalysis as mda\n", "from MDAnalysis.tests.datafiles import PSF, DCD, DCD2\n", "from MDAnalysis.analysis import gnm\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Loading files" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The test files we will be working with here feature adenylate kinase (AdK), a phosophotransferase enzyme. (Beckstein *et al.*, 2009)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "u1 = mda.Universe(PSF, DCD)\n", "u2 = mda.Universe(PSF, DCD2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Using a Gaussian network model" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Using a Gaussian network model to represent a molecule as an elastic network, we can characterise the concerted motions of a protein, and the dominance of these motions, over a trajectory. The analysis is applied to the atoms in the `selection`. If two atoms are within the `cutoff` distance (default: 7 ångström), they are considered to be bound by a spring. This analysis is reasonably robust to the choice of cutoff (between 5-9 Å), but the singular value decomposition may not converge with a lower cutoff." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "nma1 = gnm.GNMAnalysis(u1,\n", " selection='name CA',\n", " cutoff=7.0)\n", "nma1.run()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The output is saved in `nma1.results`: the time in picoseconds, the first eigenvalue, and the first eigenvector, associated with each frame." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "98" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "len(nma1.results)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "nma2 = gnm.GNMAnalysis(u2,\n", " selection='name CA',\n", " cutoff=7.0)\n", "nma2.run()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Unlike normal mode analysis, Gaussian network model analysis uses only a single eigenvalue to represent the rotation and translation of each frame. The motion with the lowest positive eigenvalue represents the dominant motion of a structure. The frequency of this motion is the square root of the eigenvalue.\n", "\n", "Plotting the probability distribution of the frequency for the first eigenvector can highlight variation in the probability distribution, which can indicate trajectories in different states.\n", "\n", "Below, we plot the distribution of eigenvalues. The dominant conformation state is represented by the peak at 0.06." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEHCAYAAACp9y31AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAXo0lEQVR4nO3de7BlZXnn8e+Pbgg0yEVBJFxsmBBMj9fmQLSIV8aMIxEwQ9AEk45Sdqa8V5KRljiRmoqpdkqDeEnGtjHVGC25eAHFhAAqMVMzQLcQERAh2GpzkVYxLUhoG575Y69jTvpyzjqXtXefs76fql1nve+67OflNM9597ve/a5UFZKk/thj1AFIkobLxC9JPWPil6SeMfFLUs+Y+CWpZ0z8ktQzi7u8eJIDgbXA04ECXgfcAVwMLAU2AmdW1YOTXefggw+upUuXdhmqJC04GzZs+EFVHbJ9fbqcx59kHfDVqlqbZC9gCXAu8KOqWp1kFXBQVZ0z2XXGxsZq/fr1ncUpSQtRkg1VNbZ9fWdDPUkOAF4AXAhQVVur6sfAacC65rB1wOldxSBJ2lGXY/xHA5uBv05yU5K1SfYFDq2q+5pj7gcO3dnJSVYmWZ9k/ebNmzsMU5L6pcvEvxhYDvxVVT0HeBhYNfGAGowz7XSsqarWVNVYVY0dcsgOQ1SSpBnq8ubuJmBTVV3flC9jkPi/n+SwqrovyWHAAx3G0MrSVVcO9f02rj5lqO8nSRN11uOvqvuB7yU5rqk6GbgNuAJY0dStAC7vKgZJ0o46nc4JvBn4RDOj527gtQz+2FyS5GzgO8CZHccgSZqg08RfVTcDO0wlYtD7lySNgN/claSeMfFLUs+Y+CWpZ0z8ktQzJn5J6hkTvyT1jIlfknrGxC9JPWPil6SeMfFLUs+Y+CWpZ0z8ktQzJn5J6hkTvyT1jIlfknrGxC9JPWPil6SeaZX4kzyj60AkScPRtsf/l0luSPKGJAd0GpEkqVOtEn9VPR84CzgS2JDkk0le2mlkkqROtB7jr6o7gXcC5wAvBD6Q5JtJfrOr4CRJc6/tGP8zk5wP3A68BHhFVf1Ks31+h/FJkubY4pbHfRBYC5xbVY+MV1bVvUne2UlkkqROtB3qOQX45HjST7JHkiUAVfXxyU5MsijJTUm+0JSPTnJ9kruSXJxkr9k0QJI0PW0T/zXAPhPKS5q6Nt7KYIho3HuA86vql4AHgbNbXkeSNAfaJv69q+qh8UKzvWSqk5IcweDTwtqmHAb3BS5rDlkHnD6dgCVJs9M28T+cZPl4IcnxwCOTHD/u/cDbgceb8pOAH1fVtqa8CTh8ZycmWZlkfZL1mzdvbhmmJGkqbW/uvg24NMm9QICnAK+a7IQkvwE8UFUbkrxouoFV1RpgDcDY2FhN93xJ0s61SvxVdWOSpwHHNVV3VNXPpjjtJODUJC8H9gb2By4ADkyyuOn1HwHcM7PQJUkzMZ1F2k4AngksB347ye9NdnBVvaOqjqiqpcCrgS9V1VnAl4EzmsNWAJdPO2pJ0oy16vEn+TjwH4Cbgcea6gIumsF7ngN8KsmfATcBF87gGpKkGWo7xj8GLKuqGY21V9VXgK8023cDJ87kOpKk2Ws71PMNBjd0JUnzXNse/8HAbUluAB4dr6yqUzuJSpLUmbaJ/7wug5AkDU/b6ZzXJXkqcGxVXdOs07Oo29AkSV1ouyzz6xkss/CRpupw4HNdBSVJ6k7bm7tvZPCFrC3w84eyPLmroCRJ3Wmb+B+tqq3jhSSLGczjlyTNM20T/3VJzgX2aZ61eynw+e7CkiR1pW3iXwVsBm4B/gD4IoPn70qS5pm2s3oeBz7avCRJ81jbtXq+zU7G9KvqmDmPSJLUqems1TNub+C3gCfOfTiSpK61GuOvqh9OeN1TVe9n8EhFSdI803aoZ/mE4h4MPgG0/bQgSdqNtE3e75uwvQ3YCJw559FIkjrXdlbPi7sORJI0HG2Hev5wsv1V9RdzE44kqWvTmdVzAnBFU34FcANwZxdBSZK60zbxHwEsr6qfACQ5D7iyql7TVWCSpG60XbLhUGDrhPLWpk6SNM+07fFfBNyQ5LNN+XRgXTchSZK61HZWz7uT/C3w/KbqtVV1U3dhSZK60naoB2AJsKWqLgA2JTm6o5gkSR1q++jFdwHnAO9oqvYE/maKc45M8uUktyW5Nclbm/onJrk6yZ3Nz4Nm0wBJ0vS07fG/EjgVeBigqu4FnjDFOduAP6qqZcBzgTcmWcZgbf9rq+pY4NqmLEkakraJf2tVFc3SzEn2neqEqrqvqr7WbP8EuJ3BQ9pP499uDK9jcKNYkjQkbRP/JUk+AhyY5PXANUzjoSxJlgLPAa4HDq2q+5pd97OLaaFJViZZn2T95s2b276VJGkKbWf1vLd51u4W4DjgT6vq6jbnJtkP+DTwtqrakmTidSvJTh/aXlVrgDUAY2NjPthdkubIlIk/ySLgmmahtlbJfsK5ezJI+p+oqs801d9PclhV3ZfkMOCB6QYtSZq5KYd6quox4PEkB0znwhl07S8Ebt9uEbcrgBXN9grg8ulcV5I0O22/ufsQcEuSq2lm9gBU1VsmOeck4Heb825u6s4FVjO4Z3A28B1c11+Shqpt4v9M82qtqv4RyC52nzyda0mS5s6kiT/JUVX13apyXR5JWiCmGuP/3PhGkk93HIskaQimSvwTh2qO6TIQSdJwTJX4axfbkqR5aqqbu89KsoVBz3+fZpumXFW1f6fRLVBLV105tPfauPqUob2XpPlh0sRfVYuGFYgkaTimsx6/JGkBMPFLUs+Y+CWpZ0z8ktQzJn5J6hkTvyT1jIlfknrGxC9JPWPil6SeMfFLUs+Y+CWpZ0z8ktQzbR+9KLUyzJVHYfirjw67fcPmaq79YI9fknrGHr/mtYXeA5e6YOKXNDILfWhwd+VQjyT1zEgSf5KXJbkjyV1JVo0iBknqq6En/iSLgA8D/wVYBvx2kmXDjkOS+moUY/wnAndV1d0AST4FnAbcNoJYFjxvfmo6Fvq/F+8pDIwi8R8OfG9CeRPwq9sflGQlsLIpPpTkjiHENpcOBn4w6iBGwHb3Rx/bDNNod97TcSRTe+rOKnfbWT1VtQZYM+o4ZirJ+qoaG3Ucw2a7+6OPbYaF0e5R3Ny9BzhyQvmIpk6SNASjSPw3AscmOTrJXsCrgStGEIck9dLQh3qqaluSNwFXAYuAj1XVrcOOYwjm7TDVLNnu/uhjm2EBtDtVNeoYJElD5Dd3JalnTPyS1DMmfknqGRO/JPWMiV+SesbEL0k9Y+KXpJ4x8UtSz5j4JalnTPyS1DMmfknqGRO/JPWMiV+SesbEL0k90+l6/EkOBNYCTwcKeB1wB3AxsBTYCJxZVQ9Odp2DDz64li5d2mWokrTgbNiw4QdVdcj29Z2ux59kHfDVqlrbPG1rCXAu8KOqWp1kFXBQVZ0z2XXGxsZq/fr1ncUpSQtRkg07ez5wZ0M9SQ4AXgBcCFBVW6vqx8BpwLrmsHXA6V3FIEnaUZdj/EcDm4G/TnJTkrVJ9gUOrar7mmPuBw7d2clJViZZn2T95s2bOwxTkvqly8S/GFgO/FVVPQd4GFg18YAajDPtdKypqtZU1VhVjR1yyA5DVJKkGeoy8W8CNlXV9U35MgZ/CL6f5DCA5ucDHcYgSdpOZ4m/qu4HvpfkuKbqZOA24ApgRVO3Ari8qxgkSTvqdDon8GbgE82MnruB1zL4Y3NJkrOB7wBndhyDJGmCThN/Vd0M7DCViEHvX5I0Al33+Edu6aorR/K+G1efMpL3laSpuGSDJPWMiV+SesbEL0k9Y+KXpJ4x8UtSz5j4JalnTPyS1DMmfknqGRO/JPWMiV+SesbEL0k9Y+KXpJ4x8UtSz5j4JalnTPyS1DOtEn+SZ3QdiCRpONr2+P8yyQ1J3pDkgE4jkiR1qlXir6rnA2cBRwIbknwyyUs7jUyS1InWY/xVdSfwTuAc4IXAB5J8M8lvdhWcJGnutR3jf2aS84HbgZcAr6iqX2m2z+8wPknSHGv7sPUPAmuBc6vqkfHKqro3yTs7iUyS1Im2if8U4JGqegwgyR7A3lX106r6eGfRSZLmXNsx/muAfSaUlzR1U0qyKMlNSb7QlI9Ocn2Su5JcnGSv6YUsSZqNtol/76p6aLzQbC9pee5bGdwbGPce4Pyq+iXgQeDslteRJM2Bton/4STLxwtJjgcemeT48eOOYDBMtLYph8EN4cuaQ9YBp08nYEnS7LQd438bcGmSe4EATwFe1eK89wNvB57QlJ8E/LiqtjXlTcDhOzsxyUpgJcBRRx3VMkxJ0lRaJf6qujHJ04Djmqo7qupnk52T5DeAB6pqQ5IXTTewqloDrAEYGxur6Z4vSdq5tj1+gBOApc05y5NQVRdNcvxJwKlJXg7sDewPXAAcmGRx0+s/ArhnRpFLkmak7Re4Pg68F/g1Bn8ATgDGJjunqt5RVUdU1VLg1cCXquos4MvAGc1hK4DLZxa6JGkm2vb4x4BlVTUXQy7nAJ9K8mfATcCFc3BNSVJLbRP/Nxjc0L1vJm9SVV8BvtJs3w2cOJPrSJJmr23iPxi4LckNwKPjlVV1aidRSZI60zbxn9dlEJKk4Wk7nfO6JE8Fjq2qa5IsARZ1G5okqQttZ/W8nsG3bT/SVB0OfK6roCRJ3Wm7ZMMbGczL3wI/fyjLk7sKSpLUnbaJ/9Gq2jpeSLIY8Nu0kjQPtU381yU5F9inedbupcDnuwtLktSVtol/FbAZuAX4A+CLDJ6/K0maZ9rO6nkc+GjzkiTNY60Sf5Jvs5Mx/ao6Zs4jkiR1ajpr9YzbG/gt4IlzH44kqWutxvir6ocTXvdU1fsZPFlLkjTPtB3qWT6huAeDTwDTWctfkrSbaJu83zdhexuwEThzzqORJHWu7ayeF3cdiCRpONoO9fzhZPur6i/mJhxJUtemM6vnBOCKpvwK4Abgzi6CkiR1p23iPwJYXlU/AUhyHnBlVb2mq8AkSd1ou2TDocDWCeWtTZ0kaZ5p2+O/CLghyWeb8unAum5CkiR1qe2snncn+Vvg+U3Va6vqpu7Cmv+WrrpyZO+9cbXfrZO0a22HegCWAFuq6gJgU5KjO4pJktShto9efBdwDvCOpmpP4G+6CkqS1J22Pf5XAqcCDwNU1b3AEyY7IcmRSb6c5LYktyZ5a1P/xCRXJ7mz+XnQbBogSZqetol/a1UVzdLMSfZtcc424I+qahnwXOCNSZYxeKjLtVV1LHBtU5YkDUnbxH9Jko8AByZ5PXANUzyUparuq6qvNds/AW4HDgdO499mBK1jMENIkjQkbWf1vLd51u4W4DjgT6vq6rZvkmQp8BzgeuDQqrqv2XU/u/g+QJKVwEqAo446qu1bSZKmMGXiT7IIuKZZqK11sp9w/n7Ap4G3VdWWJD/fV1WVZIcnezX71gBrAMbGxnZ6jCRp+qYc6qmqx4DHkxww3Ysn2ZNB0v9EVX2mqf5+ksOa/YcBD0z3upKkmWv7zd2HgFuSXE0zswegqt6yqxMy6NpfCNy+3eqdVwArgNXNz8unG7QkaebaJv7PNK/pOAn4XQZ/MG5u6s5lkPAvSXI28B18oIskDdWkiT/JUVX13aqa9ro8VfWPQHax++TpXk+SNDemGuP/3PhGkk93HIskaQimSvwTe+zHdBmIJGk4pkr8tYttSdI8NdXN3Wcl2cKg579Ps01Trqrav9PoJElzbtLEX1WLhhWIJGk4prMevyRpATDxS1LPmPglqWdM/JLUMyZ+SeoZE78k9YyJX5J6xsQvST3Tdllmabe2dNWVI3nfjatPGcn7SrNhj1+SesbEL0k9Y+KXpJ4x8UtSz5j4JalnnNWzADnDRdJk7PFLUs/Y49ecGdUnjVEaZZv9hKWZsscvST0zkh5/kpcBFwCLgLVVtXoUcUiaPj/lzH9DT/xJFgEfBl4KbAJuTHJFVd027FgkzS99G07s6g/dKIZ6TgTuqqq7q2or8CngtBHEIUm9NIqhnsOB700obwJ+dfuDkqwEVjbFh5LcMYTY5tLBwA9GHcQI2O4hyXuG+W475e+6Y3PwO37qzip321k9VbUGWDPqOGYqyfqqGht1HMNmu/ujj22GhdHuUQz13AMcOaF8RFMnSRqCUST+G4FjkxydZC/g1cAVI4hDknpp6EM9VbUtyZuAqxhM5/xYVd067DiGYN4OU82S7e6PPrYZFkC7U1WjjkGSNER+c1eSesbEL0k9Y+KfgSQvS3JHkruSrNrJ/l9IcnGz//okSyfse2aS/5vk1iS3JNl7mLHPxkzbnWTPJOua9t6e5B3Djn2mWrT5BUm+lmRbkjO227ciyZ3Na8Xwop69mbY7ybMn/Pv+epJXDTfy2ZnN77vZv3+STUk+NJyIZ6iqfE3jxeCG9D8DxwB7Af8ELNvumDcA/7vZfjVwcbO9GPg68Kym/CRg0ajbNIR2/w7wqWZ7CbARWDrqNs1Rm5cCzwQuAs6YUP9E4O7m50HN9kGjbtMQ2v3LwLHN9i8C9wEHjrpNXbd7wv4LgE8CHxp1eyZ72eOfvjZLTpwGrGu2LwNOThLg14GvV9U/AVTVD6vqsSHFPVuzaXcB+yZZDOwDbAW2DCfsWZmyzVW1saq+Djy+3bn/Gbi6qn5UVQ8CVwMvG0bQc2DG7a6qb1XVnc32vcADwCHDCXvWZvP7JsnxwKHA3w8j2Nkw8U/fzpacOHxXx1TVNuBfGPTufxmoJFc1HxffPoR458ps2n0Z8DCD3t93gfdW1Y+6DngOtGlzF+eO2pzEnuREBj3nf56juLo243Yn2QN4H/DHHcQ153bbJRsWqMXArwEnAD8Frk2yoaquHW1YnTsReIzBR/+DgK8muaaq7h5tWOpKksOAjwMrqmqH3vEC9Abgi1W1afAhd/dmj3/62iw58fNjmuGNA4AfMuhB/ENV/aCqfgp8EVjeecRzYzbt/h3g76rqZ1X1APB/gPmw1slslheZz0uTzCr2JPsDVwJ/UlX/b45j69Js2v084E1JNgLvBX4vyW77nBET//S1WXLiCmB8FscZwJdqcOfnKuAZSZY0ifGFwHx5DsFs2v1d4CUASfYFngt8cyhRz85slhe5Cvj1JAclOYjB/Z2rOopzrs243c3xnwUuqqrLOoyxCzNud1WdVVVHVdVSBsM9F1XVDrOCdhujvrs8H1/Ay4FvMRi7/JOm7n8CpzbbewOXAncBNwDHTDj3NcCtwDeA/zXqtgyj3cB+Tf2tDP7Q/fdRt2UO23wCg09yDzP4dHPrhHNf1/y3uAt47ajbMox2N/++fwbcPOH17FG3Zxi/7wnX+H1281k9LtkgST3jUI8k9YyJX5J6xsQvST1j4peknjHxS1LPmPi1ICV5LMnNE16rmvq1SZbtBvE9NOoY1F9O59SClOShqtpv1HHsyu4enxY2e/zqlSRfSTLWbJ+d5FtJbkjy0fE11JMckuTTSW5sXic19ecl+VhzjbuTvKWpX53kjRPe47wkf5xkvyTXNgvy3ZJk+9VMSfKiJF+YUP5Qkt9vto9Pcl2SDc3Cfod1+h9HvWHi10K1z3ZDPf/ugSBJfhH4HwyWjzgJeNqE3RcA51fVCcB/BdZO2Pc0Bksunwi8K8mewMXAmROOObOp+1fglVW1HHgx8L60XMGrue4HGaz5fjzwMeDd7ZouTc7VObVQPVJVz55k/4nAddUsD53kUgbLZgP8J2DZhBy9f5LxYZkrq+pR4NEkDwCHVtVNSZ7c/DE5BHiwqr7XJO8/T/ICBuu3H85gvfb7W8R/HPB04OomjkUMlrWWZs3EL+1oD+C5VfWvEyubBPzohKrH+Lf/hy5lsDDdUxj09gHOYvCH4Piq+lmzcuP2j9rcxr//5D2+PwzWgXnerFoi7YRDPeqrG4EXNqtnLmYwpDPu74E3jxeSTPbJYdzFDFZzPIPBHwEYLEv9QJP0Xww8dSfnfYfBp4tfSHIgcHJTfwdwSJLnNTHsmeQ/tm+etGv2+LVQ7ZPk5gnlv6sJy+RW1T1J/pzBKqI/YrBM9L80u98CfDjJ1xn8P/IPwH+b7M2q6tYkTwDuqarxIZlPAJ9Pcguwnp0sRd0MCV3CYLXWbwM3NfVbm4d5fyDJAU0c72ewwqk0K07nVG8l2a+qHmp6/J8FPlZVnx11XFLXHOpRn53XfCoY721/bsTxSENhj1+SesYevyT1jIlfknrGxC9JPWPil6SeMfFLUs/8f7+lIrPqnChnAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "eigenvalues1 = [res[1] for res in nma1.results]\n", "eigenvalues2 = [res[1] for res in nma2.results]\n", "\n", "histfig, histax = plt.subplots(nrows=2, sharex=True, sharey=True)\n", "histax[0].hist(eigenvalues1)\n", "histax[1].hist(eigenvalues2)\n", "\n", "histax[1].set_xlabel('Eigenvalue')\n", "histax[0].set_ylabel('Frequency')\n", "histax[1].set_ylabel('Frequency');" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "When we plot how the eigenvalue varies with time, we can see that the simulation transitions into the dominant conformation and stays there in both trajectories." ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "scrolled": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOzdd3jb1fX48feVZMnb8U48Qpy9FyETwgwkjAYoIxAg0DJL2lJGf3RACf22FMpehbB3gEAhQCBNGFmQvUP2tJ3hvYds6/7+uJIt27ItJ1Hi2Of1PH4kfYZ0xfDxueNcpbVGCCGEaMhyohsghBCibZIAIYQQwicJEEIIIXySACGEEMInCRBCCCF8sp3oBhwrcXFxulu3bie6GUIIcVJZvXp1jtY63te5dhMgunXrxqpVq050M4QQ4qSilNrX1DnpYhJCCOGTBAghhBA+SYAQQgjhU7sZgxBCiJZUVVWRkZFBRUXFiW7KcRccHExKSgpBQUF+3yMBQgjRYWRkZBAREUG3bt1QSp3o5hw3Wmtyc3PJyMggLS3N7/uki0kI0WFUVFQQGxvboYIDgFKK2NjYVmdOEiCEEB1KRwsOHkfyvSVAlBfAD49C5uoT3RIhhGhTJEAoBT/8E/YuPdEtEUK0c1arlaFDhzJgwACGDBnCE088gcvlqj2/YsUKxo8fT58+fRg2bBg333wzZWVlvPnmm8THxzNs2DB69erFBRdcwI8//hjw9sogdXAUOKKgMP1Et0QI0c6FhISwbt06ALKysrj22mspKipixowZHD58mCuvvJJZs2YxZswYAGbPnk1xcTEAV199Nc8//zwA33//PZdffjnff/89/fr1C1h7JYMAiEqBAgkQQojjJyEhgZkzZ/L888+jteaFF15g2rRptcEB4IorriAxMbHRvWeffTa33norM2fODGgbJYMA6JQKhRknuhVCiONoxheb+flA0TF9z/5JkfztkgF+X9+9e3dqamrIyspi06ZNTJs2ze97hw8fzssvv3wkzfSbZBAAUalQuP9Et0IIIfymtQ74Z0gGAaaLqaIQKoogOPJEt0YIcRy05i/9QNm9ezdWq5WEhAQGDBjA6tWrmTx5sl/3rl27NqDjDyAZhNEp1TxKN5MQ4jjJzs7m9ttvZ/r06SilmD59Om+99RbLly+vvebTTz/l8OHDje5duHAhM2fO5JZbbgloGyWDAIjqah4L0yGx/4ltixCi3SovL2fo0KFUVVVhs9m4/vrrufvuuwFITExk1qxZ3HvvvWRlZWGxWBg/fjwTJ04E4MMPP2TJkiWUlZWRlpbGJ598EvAMQgIEmC4mgAIZhxBCBE5NTU2z58eMGcPixYsbHb/xxhu58cYbA9SqpkkXE0B4Iljt0sUkhBBeJEAAWCwQmdzkYrkal+Z3H6xl9b7849wwIYQ4cSRAeDSzWC6npJI56w+wdGfOcW6UEEKcOBIgPDp1bbKLKbu4EoDSyurj2SIhhDihJEB4RKVC8UGodjY65QkQJRIghBAdiAQIj6gUQENRZqNTkkEIITqigAYIpdREpdQ2pdROpdT9Ps6PV0qtUUpVK6Wu8HE+UimVoZR6PpDtBJpdLJdd4skgmp+iJoQQzQlUue/77ruPvn37MnjwYC677DIKCgqOSXsDFiCUUlbgBWAS0B+4RinVcBXafuBG4P0m3ubvwKJAtbGeKE+AaDxQXdfFVHVcmiKEaJ885b43b97M/Pnz+frrr5kxYwZAbbnvRx99lG3btrF27VomTpxYr9z32rVr2bFjB/fffz+XX345W7ZsAWDChAls2rSJDRs20Lt3bx555JFj0t5AZhAjgZ1a691aaycwC6hXZERrvVdrvQFwNbxZKXUqkAj8L4BtrBOZbB59zGSq62KSDEIIcWwcy3Lf559/PjabWfc8evRoMjKOzZquQK6kTga8f9tmAKP8uVEpZQGeAK4DzmvmuluBWwG6du16xA0FICjYLJhrJoOQMQgh2pGv74dDG4/te3YeBJP+5fflgSj3/frrr3P11Vf7/T7NaauD1L8B5mqtmw2DWuuZWusRWusR8fHxR/+pUSm+A0SJzGISQrQtvsp9/+Mf/8BmszF16tRj8hmBzCAygVSv1ynuY/4YA5yhlPoNEA7YlVIlWutGA93HVFSqz78oJIMQoh1qxV/6gXIsy32/+eabfPnll3z77bcopY5J+wKZQawEeiml0pRSdmAKMMefG7XWU7XWXbXW3YB7gbcDHhygbmc5r8hc7qyhpLIau81CqbMGlyvwm3QIIdq/Y1nu+5tvvuGxxx5jzpw5hIaGHrM2BiyD0FpXK6WmA/MAK/C61nqzUuphYJXWeo5S6jTgv0A0cIlSaobW+sTt4hGVCjWVUJoN4QmAKbMB0C02lO2HSyirqiHcIUVwhRCtF6hy39OnT6eyspIJEyYAZqD6pZdeOur2BvQ3ndZ6LjC3wbEHvZ6vxHQ9NfcebwJvBqB5jXmmuhak1waIrGJPgAhj++ESSiurJUAIIY5IoMp979y582ib5lNbHaQ+MTo1XgvhGX9IiwsDZKBaCNFxSIDwVptB1G0c5JnB1M0dIGSgWgjRUUiA8BbSCUKiIX9P7aHs4kosCrrGmIEfySCEOLn5mh7aERzJ95YA0VB0GuTtrn2ZXVxJTJiDyOAgAEoqJEAIcbIKDg4mNze3wwUJrTW5ubkEBwe36j4ZbW0opjtkrKx9mV1cSVy4nTCHFYBSpwQIIU5WKSkpZGRkkJ2dfaKbctwFBweTktLsnKBGJEB4qXFp1hRFMaIwHVXtBJud7JJK4iMctTOXpKKrECevoKAg0tLSTnQzThrSxeRl9up0Zu0MQmlX7UymnGITIMLcAUIGqYUQHYUECLeSymr+PW87+7RZ/0DebrTWZLsDRKjdilISIIQQHYcECLf//LCTnJJK4rqalYlV2bsoKq/GWeMiPtyBUoowu01mMQkhOgwJEEBGfhmvLN7D5KFJXHbGMEq1g7z0rWSXVAAQH+EAIMxhlQxCCNFhyCA18Ng321DAHyf2xWGzsF8nEpK1s7bMhidAhDtssmmQEKLD6PAZxO7sEuasP8AtZ3QnuVMIceEOcoKSsBftqy2zkeAVIKSLSQjRUXT4DKJ7fDgf3TaGAUmRtcdqOqURm7ua7KJyAOLDzeKSMIdNupiEEB1Gh88gAEamxdROYwUI79ITB1Vs27Edu9VCZIg5FyYZhBCiA5EA4UNS9/4AHNr7M/ERjtrdmcIdNllJLYToMDp8F5MvXdJMgEhyHaLIPf4AZhaT1GISQnQUEiB8UFGpVGOjmzpEbrh3gJBZTEKIjkO6mHyxWCkJSaarOlw7xRUg3G7DWePCWe06gY0TQojjQwJEE1RMGt0aBAipxySE6EgkQDQhPKk33a3ZjOseA/n74KUz6Fn4EyCbBgkhOgYJEE2wxnYnRJcxKjIX3v0lHNpAYskWQPaEEEJ0DDJI3ZSY7ubx7UuhLBeUlVDMwjnpYhJCdASSQTQl2r2pSMkhuOJ1CInGoU2AkE2DhBAdgQSIpkR3g5SRcMkz0O9isIfhcPmXQWit+X+zN7B0Z85xaKgQQgSGdDE1xWaHm+fXvbaHE1RTBrQ8SF1UXs2Hq9IJdVgZ1zMukK0UQoiAkQzCX/YwgqpNgGgpg0jPN9flljgD3iwhhAgUCRD+sodhdQeIlsptZOSbrqi8UgkQQoiTlwQIfznCsVSVYrdaKGlhmmuGJ4OQACGEOIlJgPCXPRycJX5tO1qXQVQej5YJIURASIDwlz0MnKV+Fezz7mLSWh+P1gkhxDEnAcJf9jCoLPFr21FPF1NVjaZIyoMLIU5SAQ0QSqmJSqltSqmdSqn7fZwfr5Rao5SqVkpd4XV8qFLqJ6XUZqXUBqXU1YFsp1/s4VBTSaS9+VlMWmsy88uJCgkCZKBaCHHyCliAUEpZgReASUB/4BqlVP8Gl+0HbgTeb3C8DLhBaz0AmAg8rZTqFKi2+sUeDkCMvbrZAFFUXk1xZTWDU6IAyC2RcQghxMkpkBnESGCn1nq31toJzAIme1+gtd6rtd4AuBoc36613uF+fgDIAuID2NaW2cMAiLVVNdvF5FkDMSTFxDOZySSEOFkFMkAkA+lerzPcx1pFKTUSsAO7fJy7VSm1Sim1Kjs7+4gb6hd3gIi2VTQ7SO0ZoPZkENLFJIQ4WbXpQWqlVBfgHeAmrXWjbdy01jO11iO01iPi4wOcYLi7mKJszma7mDwD1ENSTQYhAUIIcbIKZIDIBFK9Xqe4j/lFKRUJfAX8RWu97Bi3rfUc7gBhdVLirG5y+mpGfjnhDhsJEQ7C7FZyZAxCCHGSCmSAWAn0UkqlKaXswBRgjj83uq//L/C21np2ANvoP3cXU4SlEq2hzOm7mykjv5yU6BCUUsSGOySDEEKctAIWILTW1cB0YB6wBfhIa71ZKfWwUuoXAEqp05RSGcCVwMtKqc3u268CxgM3KqXWuX+GBqqtfnF3MYUrkxE01c2UkV9GSnQIADFhdgkQQoiTVkDLfWut5wJzGxx70Ov5SkzXU8P73gXeDWTbWs2dQYQrz6ZB1SQ0uMSzBmJ091gAYsPsHCysOJ6tFEKIY6ZND1K3Ke4AEYong2jcxeRZA+GdQeRKPSYhxElKAoS/3F1MIdpkBL7WQnjWQHgChGcMQuoxCSFORhIg/GWxgi2EYN30tqOeNRAp0aGA6WKqqtEUt1C7SQgh2iIJEK1hD8Phcu8q52NPiIwGGURMmB2APNlZTghxEpIA0Rr2MOw1dYPUDXnWQHgK9cWEmwAh4xBCiJORBIjWcEQQ5M4gvtuSxcHC8nqnPVNclVIAxIU5ANmbWghxcpIA0Rr2MGzVZdx+Zg8Wbs/mzMd+4IHPNrFqbx4VVTXuRXKhtZd7MghZCyGEOBkFdB1Eu2MPg8pi7p/Ul6mjuvLiD7v4YMV+3lm2D5tF4dK6dg0EmEFqkIquQoiTkwSI1rCHQdFBAFJjQnnk8kHce35v1uwvYO3+fLYcLGLSwM61lwcHWQm1W6WLSQhxUpIA0Rr2cHCW1jsUG+5gQv9EJvRP9HlLbLidPBmkFkKchGQMojXs4eAsadUtMWGOel1Mq/bmkV0sAUMI0fb5FSCUUolKqdeUUl+7X/dXSv06sE1rg+xhjTKIlsSG2Wu7mHJLKpkycxmvLt4diNYJIcQx5W8G8SamKmuS+/V24K5ANKhNs4dDTSXUVPl9i3dF1683HaLapcmSDEIIcRLwN0DEaa0/wr13tLuUd9P7brZX7oJ9relmMmMQph7TlxsOADLtVQhxcvA3QJQqpWIBDaCUGg0UBqxVbZV7V7nWdDPFhtlx1rjYnVPK8j15AOSXSYAQQrR9/s5iuhuzG1wPpdRSIB64ImCtaqtqMwj/A0SMezX1u8v2oTUMSo6qm/ZadACCO4E9tJl3EEKIE8OvDEJrvQY4ExgL3AYM0FpvCGTD2iS7J4NoRReTe7Hc7FUZ9O0cwci0GJNBaA0vnQFLngxES4UQ4qj5lUEopW5ocGi4Ugqt9dsBaFPb5ckgKls3BgFQXFnNxYO7oJSizFlDRf4BgstyIHtbIFoqhBBHzd8uptO8ngcD5wJrgI4ZIFrVxWSvfX7x4CR+3JULQMnBHQQDFGY0umd9egFdooJJiAw+mtYKIcRR8StAaK1/6/1aKdUJmBWQFrVl9gjz2KpBajMGMTA5km5xYWw9VAxARdZOc0FRZr3rndUurnllGZOHJvPI5YOOvs1CCHGEjrTURimQdiwbclI4gmmuIXYrZ/SK45fDU4C6jKImx71YruQwVFeCzQSSDRkFlDlr2JXVuhXbQghxrPk7BvEF7imumIHt/sBHgWpUm3UEAQLgnV+Pqn0eE2Y2E7IW7q27oCgTYroDsGy36YLandO6FdtCCHGs+ZtBPO71vBrYp7Vu3Hne3h3BGERDnmmvjqJ9YAuG6goo9A4QZq1ETkklxRVVRAQHHV2bhRDiCPk7zXWh18/SDhkcACxWCAptXQZRlgc/vQA5ZswhKiQIpSC8LB1S3ZmFe6DaWe1i9b58kqLM4PS+3LJj2nwhhGiNZgOEUqpYKVXk46dYKVV0vBrZpvhbsK80FxbMgKcHwbw/w9KnALBaFKkhVYRUF0K308217gCxMbOA8qoarjotFZBuJiHEidVsF5PWOuJ4NeSk4U+AcLng1XMgfx8MvBzy9kDmmtrT/YNzoAxI6Adh8VBkAoSne+mqEak8vWAHeyVACCFOoFbtB6GUSlBKdfX8BKpRbZo9vOWFcmW5kL8XJsyAK16HXhMgawtUmimuvYNyzHXRaRCZXJtBLNudS5/ECJI6hZAUFcweCRBCiBPI3/0gfqGU2gHsARYCe4GvA9iutsse1vIYRLGp2kq0eyZw8ghAw4F1AKRZs9znu0FUChRmUFVjxh9Gd48x18SHSYAQQpxQ/mYQfwdGA9u11mmYldTLAtaqtszHtqONFLkDRGSyeUwebh4zVwOQog+RQydTHTYqFQoz2Ohe/zCqeywA3WIlQAghTix/A0SV1joXsCilLFrr74ERAWxX2+XPGIRndXSke3+lsDiTLbgDRGL1Qfa6EnG5NEQlg7OENdv3ATAyzZ1BxIVRWF5FvtfeEQVlTiqqGm/D8adPN/D+8v1H972EEKIBfwNEgVIqHFgEvKeUegazmrrj8Wdf6qIDoKwQnlB3LPnU2gAR48xkn06guKLadDEBe3Zto3diOHHhZp1EWpxZc+GZyVTj0lz07BJmfPFzvY/KLq7kgxXpfLf18LH4dkIIUcvfADEZM+/mD8A3wC7gkpZuUkpNVEptU0rtVErd7+P8eKXUGqVUtVLqigbnpimldrh/pvnZzsDzZwyi6ABEdDHrJjySR5jMIm8P4ZWH2edKJLe00nQxAQUHd9dmD1AXIDwzmVbvyyezoJxvNh2kusZVe92i7dkAso2pEOKY8zdA3AZ00VpXa63f0lo/6+5yapJSygq8AEzClOa4RinVv8Fl+4Ebgfcb3BsD/A0YBYwE/qaUivazrYHl8GcMIrOue8kj+VTzuPlTAPbqRLMvhHucIro6i35dImsvT40JxWpRteMQX286CEB+WRVr9hfUXveDO0BkS4AQQhxj/gaICOB/SqnFSqnpSqlEP+4ZCezUWu/WWjsx1V8ne1+gtd7r3njI1eDeC4D5Wus8rXU+MB+Y6GdbA8seBjVOqG5m29CiA40DRJfBYLHBxtkA7NeJ5JVWQXgiLksQScpMcfUIslpIjQ5hT24pWmvmbTrEqLQYgqyKBVtMd1KNS9dmENnFlWZMQwghjhF/S23M0FoPAO4EugALlVILWrgtGUj3ep3hPuYPv+5VSt2qlFqllFqVnZ3t51sfJc+uclVNZBFaQ9HBuhlMHkEhkDgAsswYwj6dQF5pJVgslNgTSFI59Eqsvy6xW1wYe7JL2ZhZyIHCCq44NYUxPeKY//NhtNasS8+nsLyKEadEU+3SFJRXHetvK4TowFq1UA7IAg4BuUBCC9cGnNZ6ptZ6hNZ6RHx8/PH50JZ2lasoNMEjskvjc+5uJu2IJJ8Ik0EA2ZY4TrEVEBVSvzBfWlwYe3NLmbvxEFaLYkL/RCb0S2BPTim7skv5YVs2FgWXu0uJZxVXHJvvKIQQ+L9Q7jdKqR+Ab4FY4Bat9eAWbssEUr1ep7iP+eNo7g2sliq61q6BSGp8Ltk9MzgmDYfNasYggH3VMaRYGg/ppMWFUeas4eNV6YzpHkunUDvn9jO9e/N/PswP27I59ZRoeiaYrKbhOMSi7dlsyCho9L5CCOEPfzOIVOAurfUArfVDWuufW7wDVgK9lFJpSik7MAWY4+fnzQPOV0pFuwenz3cfO/Fa2lWu4SI5b+4MQkWnERNmJ7fESY1Ls70iitiaHHDVX+PgmcmUW+pk4sDOACR1CmFgciQfr05nY2YhZ/VJICHCTI3NKqofIP706UbueHcNldWN104IIURL/B2D+BOwUSmV5G8tJq11NTAd84t9C/CR1nqzUuphpdQvAJRSpymlMoArgZeVUpvd9+ZhVm+vdP887D524rW0aVDDRXLe4npDVFdIGUFMmJ38MifpeWWk18RgocbsLuelW6z5LKXg/P518wLO65fI7mwToM7sHU+8O0Bkl9QFiKoaFwcLy8ksKOe9ZbKITgjRev7uKDcdeAg4TN2MIw00282ktZ4LzG1w7EGv5ysx3Ue+7n0deN2f9h1XLQaIA4CC8M6Nz1ks8Ls1YLERs2UFeaVOth8uJlOb8hoUZtQLLEmdQrDbLAxOjiIhMrj2+Hn9Enl6wQ4SIhwMSIpEKUWY3VovgzhUWIFLg91m4YXvd3LVaamEO450h1khREfkbxfTXUAfdxfTIPdPS2MQ7ZNnFlNTXUzFB8wKapvd93lrEChFTJi9NkAc0HHmXGF6/Ustir9e1I/7LuhT7/iApEh6JoQzaWBnlFIAJEQG18sg0vPNZkN3ndeL3FInry/Z08ovKoTo6Pz9kzIdKAxkQ04aDk+AaCaD8NW91EB0qJ38UifbD5eY8YpKzNajDdwwplujY0opvvzt6dgsqvZYfLiDrKK6WUwZ+eUAXDI4iXX7C3hl0W6uG30KMWFNBC4hhGjA3wxiN/CDUupPSqm7PT+BbFib5eliOrjBrHloqOiA7wHqBmLC7BRXVrPpQCEpXTpDSDRkb/W7GcFBVmzWun998ZGOehlERn45FgWdo4K594I+lDireWXxbr/fXwgh/A0Q+zGrme2YVdWen47HHg4Dfwmr34APppitRb0VZZo6TC3w/CW/O7uUXonh0Ot82PoV1BzZYrf4cAfZRd4BoowuUSEEWS30ToxgWGon1u2XKa9CCP/51cWktZ4BoJQK1VqXBbZJbZxS8MvXIGUkzH8AXhoHU2dD54Fm8VxFoV9dTN5dPX0SIyDsMtjwIez+wexA15LyArAFQ5AZvE6IdFBcWU25s4YQu5WM/HKSo0NqL0+ODpU1EUKIVvF3odwYpdTPwFb36yFKqRcD2rK2TCkYfTvcvACqK2DRY+Z4sSmo508XU3RoXYDonRgBPc4BRxRs/q9/bXjjQpj3p9qX8e4y4Z7Fcpn55aR4BYikTsEcLKyQek1CCL/528X0NKaAXi6A1no9MD5QjTppdBkCg6+GbV9DeX7zayAaiA03AUIpzEpomwP6XgRbvoTqFiqzulyQsw12zK8dB/FMg80uqahdA5ESHVp7S3KnEJzVLnJLnfDz57Dr+yP4wkKIjsTvWkxa6/QGh2R5LsCQKaa66+b/miJ94PcsJoBTYkIJDnLvGzHwcqgsbPmXd1kuuKrNtNgCsxOdJ4PIKqqsXQPhnUF0iTLPDxSUw3f/gDm/M4FGCCGa4G+ASFdKjQW0UipIKXUvZnW06DIU4vvC+lmtyiCiQ01hvt7eFVy7nwXBnWr3jGhS8YG653sWA2YMAsxqas8aiIZdTAAHC8uhshgK98PeRS22UwjRcfkbIG7HlPpOxhTNG+p+LZQyWUT6cti3FEJiTGnvFtisFsb3juc8rxIaWIOg3yWwdS5UNVOZtfhQ3fO9SwCICbVjtSiyiipr10CkNuhiAsgsqDABAmDNO35+SSFER+RvLaYcrfVUrXWi1jpBa31dSzvKdSiDrgIU7PrOrwFqj7d/NZKrRqTWPzjwcnAWw85mttvwFARMOc0ECK2xWBRx4XayiivqrYHwiAoJIiTIyoH8UvP+Fhts+cKMnQghhA/+zmJ61sfP35VSk1u+uwOISobuZ5rnvvaBaI1u4003045mitcWHwIUDLwCijIgfy8A8REOsosrycgvo3NkMEFeC+mUUiR1CiYv313zcMDlUFNZu8OdEEI05G8XUzCmW2mH+2cwpsjer5VSTweobSeXIdeYRz/GH5pltZntSQ9vbvqa4gMQFg89zjav3d1MCRHBZBWbLibvGUweSZ1CKCxwZwzdxkHnQbBWupmEEL75GyAGA2drrZ/TWj8HnAf0BS7D7NUg+l4MYQnml+7RShgAWVuanmVUfMhkKnG9TaDYawaq48NNBtFwDYRHUlQIJUXuDMIRAcOuh4PrTdkQIYRowN8AEQ2Ee70OA2K01jWYMnPCEQ5/2Awjfn3075U4AKrKIL+JCqxFB005D6Wg2+m14xAJkQ5ySirdayB8BIhOIVSVumsuOiJh0JVgtcO694++zUKIdsffAPEYsE4p9YZS6k1gLfBvpVQY0Mxoagdjs5tf2kcrsb95bKqbqfhgXb2nbqeb6bX5e4iPcODSuNdA+OpiCiZMuWdH2cMhNAaShsPhTUffZiFEu+PvLKbXgLHAZ8B/gdO11q9qrUu11vcFsoEdUnw/QEGWj51dqyuhLMcrQJxhHvcuqd16FGgygwjHTIHF4V5/ERYHZW1jsz4hRNvSbIBQSvV1Pw4HumD2hUgHOruPiUCwh0JsD99/2Xu2JfXMlorrDaGxkL68dutRaCqDCCFCuWstegJEaIxZmS2EEA20VM31HuAW4Akf5zRwzjFvkTAS+vsOEJ5yHp4MQimz9qIkm4QIs+6h4RoIjy5RwY0ziNBYEyC0PjbdY0KIdqPZAKG1vsX9ePbxaY6olTjQLGRzltZtUgR1FWO995wIjYHy/NoMonNkMHZb4+QwOMhKgt1pQrt3gHBVQWURBEcF6MsIIU5GLXUx/dHr+ZUNzv0zUI0SuAeqNWQ12GXOV4AIiYbyPIKDrEQE23x2L9W+raOKSuUwZT3ABAiQbiYhRCMtDVJP8Xr+pwbnJh7jtghviQPMY1aDmUzFB83U1NCYumMhMbUlM/p3iWRo105Nvm2cvZJSvAawQ+PMowxUCyEaaGkMQjXx3NdrcSx16gZBYY2nuhYdhIjO9ccLQqJNgHC5+OCW0c0OJURbKylyhRCtNUopySCEEE1qKYPQTTz39VocSxYLJPRrHCCKD0JEg3IeIdGgXVBZhMWizC/+JkRZKijSwRRVVJsDnkxEAoQQooGWAsQQpVSRUqoYGOx+7nl9DGpKiGYlDjABQnvF4mJ3BuHN80u+vOVuonDKKdEhZl8IqMsgSnP8atKDn29ixhfN1IkSQrQbzQYIrbVVax2ptY7QWtvczz2vg45XIzusxBbA/gAAACAASURBVAHml75n7QO46zD5yCDAr9LdwbqMEkLMznJgZjNZgvzKIMqdNXy4Mp2fdkm2IURH4PeWo+IE8AxUe9ZDVBSBs6RxBhHiySBaDhD2mlJKCDEbB4EZy/CshWjBkp05VHr2tRZCtHsSINqyBE9NJnfJDc9Ocr7GIADKWg4QVmcJZd4ZBPhdbuPbLSaTyS914nLJEJQQ7Z0EiLYsNAaiu8H2b8xrz17UTY5BtBAgtEZVFoMjon6A8KPchsulWbAlC4uCapemqKLK/+8hhDgpSYBo60bdbva63r+8LoNoOAYR7F730NIgdXWlWTUdHElOiVeVdncX0+p9+by+ZA+fr8vkx505lFZW116yPqOAnJJKxveOB5BuJiE6gJbWQYgTbfgNsOjfsORJSB1ljjXMIKw2cES1nEFUFgOgHBHklXplAKGxUJbDQ3M2szGzsPZwr4Rw5kw/nRC7lQVbDmO1KH45PIUftmWTV+qkR/yx+IJCiLYqoBmEUmqiUmqbUmqnUup+H+cdSqkP3eeXK6W6uY8HKaXeUkptVEptUUo1XMXdcdjDYNQdpptp57cmEHjXZvII6dTyOEJlEQDW4EjyvTOA0FgoLyCvqIxLhiSx4O7xPH7lEHZml/Dwl2b849stWZzWLZq0OPPZuSWSQQjR3gUsQCilrMALwCSgP3CNUqp/g8t+DeRrrXsCTwGPuo9fCTi01oOAU4HbPMGjQxp5M9gjYN+SxtmDh2c1dXPcGYQtNIq8Mifas74iNBbQVJflk9QpmJ4JEVxxagq3je/BByv2M3PRLrYeKua8fonEhtsByC2VjQSFaO8CmUGMBHZqrXdrrZ3ALGByg2smA2+5n88GzlVmGbAGwpRSNiAEcAJFAWxr2xYSDaf9yjyP7OL7mtCYlscg3AHCERaFs9pFmbPGfa9ZLBfuKiQm1F57+T3n92ZIaif+OdcUDJzQP5GYMHM+TzIIIdq9QAaIZMzmQh4Z7mM+r9FaVwOFQCwmWJQCB4H9wONa60a//ZRStyqlVimlVmVnZx/7b9CWjL4TrA6ISvF9vhUZRHC4GdTO83QzuQNEDMW1AQAgyGrh2SlDCXfY6JkQzimxYThsViIcNhmkFqIDaKuD1COBGiAJiAYWK6UWaK13e1+ktZ4JzAQYMWJE+56YH5EIN35Zv8y3t5AYP8YgTIAIjYgGcsgvc5IaE1oXIFRxbReSxymxYbx38yhs1rr6TjHhdgkQQnQAgQwQmUCq1+sU9zFf12S4u5OigFzgWuAbrXUVkKWUWgqMAHbTkaWObPpcSDRUFIKrBixW39c4TYAIjzIBomEGEa2KiQlzNLptSGr98uGxYXbyZAxCiHYvkF1MK4FeSqk0pZQds7fEnAbXzAGmuZ9fAXynzcjpftzbmSqlwoDRQIOdc0Q9oTGANkGiKe4MIqqTWViXX9awi6mI2DC7z1u9xYQ5ZBaTEB1AwAKEe0xhOjAP2AJ8pLXerJR6WCn1C/dlrwGxSqmdwN2AZyrsC0C4UmozJtC8obXeEKi2tgv+FOyrLAZlJToiEqBuLURQME5rKDGq/hhEU2LDpItJiI4goGMQWuu5wNwGxx70el6BmdLa8L4SX8dFM/wp2OcusxEREoTVoup1E5VZo4izlBBqb6J7yktsuJ38UjNNtrm9J4QQJzcptdFe1Bbsa2agurIYHJFYLIroUHu91dTFligSbKV+/cKPCbObekzl1S1eK4Q4eUmAaC/8KdjnziAAYsKC6q2mLlARxKpivz7KM9MpRwaqhWjXJEC0F7VjEM1lEEW1ASI61E5eWV2AyHGFE42fAcI90ylPxiGEaNckQLQXwVGAakUGYa+XQWTXhBOpm5kB5cUzkC0zmYRo3yRAtBcWqwkS3mMQX/wefny+7rVXgIgOs9dNcwUOOMMIdpVDVUWLHyX1mIToGCRAtCfe5TaqnbDufdj8ad157wwi1E5+WRUul6aiqobDNe4KsS3VcwKpxyREByEBoj0JjakLEIc3QY0TsraCy2WOVZbUyyBqXJriimryy5zkaXO8dme5T2+DeX/x+TFSj0mIjkECRHsSEl2XARxYYx6rSqEw3ZTgqCoFh1kkFxMWBEBemZPcEif53gHiwDrYMAu2zW34CbWkHpMQ7Z8EiPYkxCuDyFxbdzx7a22ZDe9ZTGBmIuWVOsnDK0AsfcY8z98LVV57V3uRekxCtH8SINqTkGgo8wSI1XVblGZt8QoQ4YDXOIInQGiTWZCxGn7+DGJ7gnZB7k6fHyX1mIRo/yRAtCehMVBZaAr25WyD7meb8uD1AkT9DCK/1EluqZNCwtAoWPkqWGww6TFzffY2nx8l9ZiEaP8kQLQnnsVyuxeav/6Th0NCP8huHCBqM4gyp1kPYbGZfa1rKmHIFOh2Oiir6Z7ywbsekxCifZIA0Z54AsTOBeYxaTjE94Ps7VBRYI65B6lD7VYcNkttBhEdGoQKjQUUjP092BwQ073JACH1mIRo/9rqjnLiSHgquu76DqJSITweEvpCdTkc2mjOuTMIpRQxYXbySp0UVVSZjKLrGDhlHMT1NNfG92m6i8lrsVxUaFBAv5YQ4sSQANGeeDKIwnTo595yI6G/ecxYaR7dAQLMOER+mZPCcneAmOy16hogvi9s+9osurPV3yfCU48pt9RJ9/hj/k2EEG2AdDG1J6HRdc+Th5vH+D7mMX2FefQKEJ4MIrfUWfsLv574vqBrIG9Xo1NSj0mI9k8CRHsS4hUgktwBwhEBUV3rFtDZw2svMfWYqsgrdRId5qObyBNcfIxDeLqYWqzomrurrntLCHFSkQDRnjiiQLn/lSYNrTue0Nc82sNNUT+3mNAgsosr3V1MPjKIuF6A8jkOUZdBNLNYrvgwvHY+vH0p1FQ1fZ0Qok2SANGeWCwQ3Alie7nLf7vFuwOEV/cSmAyipLIarc26hkaCQiC6m88MosV6TC4XfHY7lOWYH8/MKiHESUMCRHsT3xd6nFP/mGegukGAiPEKCjG+AoTn/ZqYyRQTbm+6i2nZi7DrO9YM+DOltmhq1rzrV/Ob9L+/mmykidIfQohjTwJEezPtC5j4SP1jCU1kEKF1QcFnBgFmHCJ3J9Q0Xu8QG2bnQIGPX9gH1sGChyhNm8iUdYP4oGI0rm3fsGXX3tZ8kzo/fw4/Pgfpy+G7/zuy9xBCtJoEiPbGaqs3zgBAXB9A1RughgYZRHgzGUSN0xTua+DM3gms2pfPx6vS65/47u/o0BjuqbgZh81Kz/NvJYhqPnrzKWZ8sZmvNx4kPa/Mv1XYhRkw53eQNAyG3wA/vQD7fmz5PiHEUZMA0RHYQyGuN4TVX7DgnUHEhDaTQYDPcYg7z+7BmO6x/PWzTfx8oMgcrKmG/cvYm3Au3+xxcu8FfThr/DlUJwzkxrCfeOenfdzx3hrOeOx7rnzpJyqqapput6sGPr0VXNXwy9fggkcg+hT47A6zt4UQIqAkQHQU13wA59fvnvHOIKKb6mKK620efQQIm9XCc9cOo1NoEHe8t5rC8iqzUZGzhFf2JTIoOYrrRp9irh02lVMqtrF5elc+v3Mcf5zYh1X78nnw801Nt3np07BvKVz4OMT2MJVoL/0P5O+DBX9r3ff30UXmj1V787j/kw24XFJzSnQ8EiA6itgeENml3iHP2ofIYBtB1ib+U3CEm3UUB9aCjy6huHAHL04dTmZ+OWc//gPPvfkOAN+Xdef/Lh2I1aLMhYOuBIsNx9o3GZIcyW/O6sn0s3vy0aoMZq3Y3/hzc3bAD/+C/pNN8UCPU8bCaTfDqjeg6KB/333BDHhmCJRk+3e9l6cX7GDWynS2HS5u9b1CnOwkQHRgDpuVcIeN2HAfayC89b4Atn4Jn95SVxXWy6mnxPDSdadydp8Exjl2kmtL5PeXncWQ1E51F4XHm1/2K1+BpwfBgof4w+AqzugZy4NzNvP9tixW7c3jf5sPsXZfLnzxezPNdtK/Qan6Hzj6DrPCe+07LX/J7G1mA6SiDJh7rx//VOqk55WxdFcOAMt257bqXiHaA6nF1MFFhwU1PcXVY9JjEJEI3//TZBJXvQOJ/etdcl7/RM7rlwBPbIV+ZzJlZNfG7zP5RehzIWz4CJY+i3XJU7wZ04s3HcP4vzfT2aWTAbja+j3Dgpayd9wjdItIbPw+sT2g+1mw+i04457Gg/Le5v3ZDM4Puw6WvWBmRPWfDIDWmvKqGkLtvv83mL06AzBdcT/tyuWmcWnN/3MSop2RDKKDG9AligFJkc1fZLHA+Pvghjkmg3j/aqiqaHxd/l4oOQRdR/t+n6BgGHQFTP0I7tkGFz2JNbIzv6r5mG8d97E+9UmWnLePh0NmsZr+nP1tKte/tpzlu3Mbz3ga8SuTFXgvwNvyJXxxFxQfMq+3/8+cP/OPMGEGdB4MX90DpSYbmPHFz4z913ccLGw8Vdfl0sxencHpPeM4t28Cy/fktX4cotoJmWtad48QbYhqLxu+jBgxQq9atepEN6P92/U9vHMpTPg7jPtd/XPrPjCrp+/4qVGG0aziw7DhQ1j9pikMaHVQfvNC3t5u55XFu8kpcTLilGjuntCbsT3jzD01VfDUAFNz6tpZkL4S3rzQTMl1RMK5D8LylwFt2mOzm5pQM8+C1NEcihjA5+syqMBOp+Q+TLv4XDOlN9gEyyU7crjuteU8d80wqmpc3P3Rer763ekMSHKvUK+qMLv2dRnS9Pf65s9mweDv1pi9NYRog5RSq7XWI3ydkwxCtE6Ps6HnebD4cSjLq39u/0+mxIentIe/IhJNsPntarjxK7jhM0K69OO2M3uw+I/nMOMXAzhYWMF1ry2v7fbBGmS6jXbMM5VqZ11LWXAib/V/1ayZmHsv5O6AC/5JUbViQ0YBdB4EE/6OzlxNzKbXmWabz29tnzPt8L/gtQnw7FDzVz/w4ap0OoUGcf6AREZ3jwVg2W6v7ztnOrw83gQmX/L3woqZgDb7cwhxEpIAIVpvwt9NV9Oix+sf378MUkebLqkjoZTZ6vSUsbWHQuxWpo3txvy7xzO2Rxz3fryed5ftMyeHTzMzq968iBpnGVcU/p6/rQllwzlvweWvwrl/g17n89f/buIXzy/lrllryR98My+OW0LvirdYctUmqv90kGlhL/Jq0DVQlguHNlJQ5mTe5kNcOjQZh81KUqcQTokN5add7oHq3T/Axo/N8/kP+JzdxXf/Z7ZxDU+EXd/zyeoMzn9qITnNFTcULcoqrqC6xnWim9FhBDRAKKUmKqW2KaV2KqXu93HeoZT60H1+uVKqm9e5wUqpn5RSm5VSG5VSwYFsq2iFxP4wdKr5CzlvjzlWmmu6XJoafzhKoXYbr04bwbl9E/jrZ5t44fuduKK6Qs/z0DVVPBh0N1nB3Qi1W3lveToMvhLOuJtDRZXM3XiQAUmRfLnhIBOeWsiz3+7gwkGdOa9/InZHMLddfgEzi08HYN68z5kycxnOahdXjUit/fwx3WNZsSeXGmeFGceITqPm/EdM1rRtbv3GHlhrAsiYO6HPJGp2LeSBT9ey/XAJHzVcdd4G5Jc6ue/j9Wxv41N59+WWcsaj33Ptq8spKJN9SI6HgAUIpZQVeAGYBPQHrlFKNeyY/jWQr7XuCTwFPOq+1wa8C9yutR4AnAVIvei25Oy/mG6eT242XTzpy83xrmMC9pHBQVb+c92pXDS4C/+et43L/vMjW0Y/ysxeL/JeXl+evGook4cm8fn6TLNoD3hn2V5qtOY/U09lzvTTSYwMJsRu5W+XDKh937E94jh92EDSXfHofcuIDA7ib5f0p7/X4P3o7rEUVVSTM+/fkLuT92J/y5gF3amO6QXz/1a3EE9r+N8DEBoL435PUfIZWKuKOTNsP0NSO/H+8v1tbtHdrJXpfLw6gxteW0FGftmJbk6Tnl6wAw2s21/AFS/91Kbb2loVVTXUtLH/LiCwGcRIYKfWerfW2gnMAiY3uGYy8Jb7+WzgXKWUAs4HNmit1wNorXO11s3UZBDHXWQXuPhpyNtt+u/nTAer3fT/B5DdZuH5a4bx1NVDyMwv46LXtvLIxihuG9+d8b3jmTrqFCqqXPx3TQYVVTW8v3w/5/VLpGtsKP2TIvli+uks/X/nkBhZPyH95+WDiOw9jgui9vPRbaMbTWkd3T2W7uoAsWufY3vsufxlU2eyymr4LO5WM9ax6nUz1jDnt7B3MZz5/6gOCufuFVHUaMWMgYe5+fQ0MvLLWbij9Qv2AkVrzcer0umVEE6ps5obXl/R8iZQJ8D2w8V8ti6Tm8Z14+1fjySrqILLXvyRrYeKTlibtNYs351Let7RBSqtNde8sowLnl7ku/jlCRTIAJEMeOfTGe5jPq/RWlcDhUAs0BvQSql5Sqk1Sqk/+voApdStSqlVSqlV2dlt53+6DmPI1XDXRlPCw2IzaxOCAt8TqJTismEpfHv3WVw/+hQmDezMPeebmlEDk6MYktqJd5fv5/N1meSXVXHTuG6191osijBH43UPwUFWonqfjio+aPb09qY1nXd9zBeOByhz2bg+8zJ+OTyFiwZ14aFtXalOGQ1f3wfvXAabPsE5cAqvV5zFOU8sZMFeJwUxg0nI+pELBnQmLtzBe54xlDZg1b58dueUcuv47rw27TQy88u56Y0VlDmPrDTJ0XK5NCv25PHYN1tZva9uUsAT/9tGuN3G7eN7MLp7LLPvGItVKW56YyWHi+qmXOeUVHLbO6uYu9HPVfZHQGvNou3Z/PI/P3L1zGVMfmEpmzILj/j9Fu/IYe3+AnZnl3DlSz+xO7vt1Blrq4PUNuB0YKr78TKl1LkNL9Jaz9Raj9Baj4iPj294WhwPjnAY+1u4ewtM+eC4fnRUaBAzJg/kP9edit1W95/y1FFd2ZlVwiNfb6Vv5wjGuGchtSh1lHncv7zuWGmuWfcxZzqHw/pwYeU/GNSvH4/+chC/ObsHJZU1fNj5Hhh5K0z5gDkXLGHI+st4+OudJEY6ePn6U4kdPBEyV2OvKuTq01L4bmsWmW3kL8WPVqYTZrdy0eAujEyL4blrhrEhs5C/f/mzX/dnFpQz8elFR/0LWWvNSwt3Mf7f33PVyz/x4g+7uPKln3h6wXbW7M9n3ubD3HxG99qaYb0TI3j9xtMoKq/i12+tpMxZTWZBOVe99BPzNh/mno/WszPr2P+i3X64mKtnLuOG11dwqLCCv1zYj5AgK9e8soxVe/NafgMfXl60i8RIB7PvGEtFVQ1XvvTTUQWcYymQASITSPV6neI+5vMa97hDFJCLyTYWaa1ztNZlwFxgeADbKo6WxWpKjbcBlwxOIjLYRoE7e1ANS3U0JaG/WXWd7hUg5j8Au7+Hif9CT/uCi8eP4flrh2GzWhiQFMVZfeJ5Yq2i/Lx/MbdqGHd9spVBKVF89bvT+fj2sVwwoDN0Pxu0C/YsYsppXdHgu/5UC8qdNbyyaDdTX112TPrfSyqr+WrjQS4ZklS7mvz8AZ25dXx3PliRzvyfD9deu/lAIde/trzRL65nF+xg66Fi/vDhOtanF7T4mS6X5qsNB2vHiDy+3ZLFv77eSnKnEJ68aggr/nwulw5N5ukFO5gycxkxYXZ+fUb9br/+SZE8d+0wfj5QxG3vrOaK//xIdkklL1w7nOAgC7/9YC1VS541kwaOUkVVDf+et5ULn1nMjsPF/P3SgXx/31ncMr47H90+hrhwB9e/toKvNhz0r4y928aMQpbuzOVX49IY3jWaj28fg8Nm4YHPN7XqfQIlkAFiJdBLKZWmlLIDU4A5Da6ZA0xzP78C+E6bfyrzgEFKqVB34DgT8O9PGtHhhditTB19Cl2igpk8tGGvZjOsNkg+tS5AFGaasiCn3gij76BHQiT3T+pLcFBdaY/fnNWTvFIn985ez+8+WMvwrtG8edNpdQvqAFJGgD0Cdn1HakwoZ/dJ4IPl+zm09mv44Bp4vI8pTuhFa82Wg0Us353Lwu3ZvLJoN2c89j3/mLuFZbvzuGvWuqOe7jl3w0HKnDVcOSK13vG7J/SmX5dI7v9kA9nFlSzdmcPVLy9j8Y4c7pu9gSr35+7NKWX2mgwuG5ZMXLiDW95exaFCHyvs3Vwuzf2fbuDO99fwp0831Dv+xPztdIsN5b2bR3H58BQSIoN58uqhPH31UELtVu45vzfhProGz+mbyAMX92fxjhyqajQf3jqGiwZ34fErh1B4cDdBCx6Aj244qvLwWmuufWUZL3y/i8lDk/n2HtO16bCZ/w6SO4Xw4W2jSYsL48731zD5haX8sC3Lr1/wLy/aRYTDxjWjTGma7vHhzL5jLC9OHe7/HzYBFLA/+bTW1Uqp6Zhf9lbgda31ZqXUw8AqrfUc4DXgHaXUTiAPE0TQWucrpZ7EBBkNzNVafxWotor2577z+/D7c3vV+2Xul9RRZhFgZYlZBa1dZrpqE0amxTDilGi+2nCQQclRvH7TaY1rO1mDIG28Kf2x9BketW+hsHoRnT/PpMTaiVCbC8vnd8JNX4PFitaah7/8mTeW7q33NmN7xPKf64ZzoKCc389ax3Pf7eQPE3q37vt5VFfy849f0COuH8O7dqp3ymGz8syUoVz83BJueH0FO7OKGRireDrkdabtn8RrS5K4/cwePPvtDoKsij9d2Je8Uie/fPFHbnl7FR/dNoYQe/1/7jUuzf/7ZAOzV2cwJCWKuRsP8e2Ww5zbL5GvNx1iy8Einr56KLYGVYUvHZbM5KFJzf6yvHFsN7pEBTMwOYqU6FAAzu2XSHnvDNgPFOyHb2fAhf+uvefzdZkMSIqkZ0JEE+9aZ+uhYtbsL+DPF/bl1vE9fF6TEBHMnOnj+HRtJs8s2MGNb6xkTPdYHp48gF6Jvj9jf24Zczce5Jbx3YkMDqo9ntQpxDzJ3wt7Fps9WVJHttjOQAhon4DWei6me8j72INezyuAK5u4913MVFchWs1iUQQ3V8SvKV1HmaCw61tT+mPAZRDdrdlbHrykP28s3cuDF/ev9z96Pb0vgG1fwfwHiQ+NIyq1N59abuKBXb24yLqcx9JfRC9/CTXmTmYu2s0bS/dy7aiuXDyoC44gC1Ehdnom1O0IuHB7Ns99t4NxPeMYmRbTqq+ot8+j5LN7eahsP/MH/Aulzq476SyDJU/Se8yd3D+xLw9/+TOju8fw+qn7Cf1iAf+MtfCrBcn0SYzgs3WZ3HxGdxIigkmICOaZKcO45Z1V/OrNlbw6bUTtZIDK6hr+9OlGPl2TyV3n9eI3Z/Xk4ucW8+DnmzktLYanFmynV0I4lwxJ8tnelv6SVkoxcWCXRscnBW/iEPFsjBjHhBUzzb/LU8aSVVTBXR+uY8Qp0Xx8+1gf71jf3I0HsSi4fHhKs9fZrBauGpHKpUOT+WDFfp6cv51Jzyzm16enccmQJMqraiiprCavxElWcSU/bMvCalH8qmERyNVvwZIn63ZxjEyG3284IV24UotJCG/lBfBoN4joDMUH4bZFzddb8perxnQjRXSGkLq/2NPzyvjzpxu4cd/9nGH7ma/GzeYPC4q5aHAXnpsyDIvF9y/HkspqLnp2MZaqcj6/e0LTganed8un9MNbCNs7n12uLsRbSwnufTb2a96uu2bNO2bK8sR/oUfdzvI9eQzr2gnHV7+Hte+grXbOqn6ejKoIHDYLi/94dr1y8Z+tzeSej9czNLUTb9x0Gvtyyrj34/VsO1zM3RN687tzewFmI6YrXvqJfl0i2XKwiJlX9eL8sN2QdgbYw47sn7G3aic8lsbqqAncdOAXrIt7CIvVCrcv5e3VWTz4+WYAZt8+hhHdmg6wWmvOfXIhiRHBfHCrexHojvkQlQIJ/ZptQm5JJY9+s5WPVmX4PB/usHHb+O781v3PBDDla57sb7KGoVNNWftv7odrZkGfSb4/aOe3YAuGbuOabU9TmqvF1DZGFYVoK0I6mf/xs36GHuccm+AAZhA/oXGNqtSYUN761Sje//YRKpZcTtKi+xib9hRPXjWkyeAA5pfLK+Mr6D53Kgs/vJtzpz3g87pth4p5+MvNZBWU8I+SvzFUb+EJdR1dL7qbK7KeQ234CKrKzd4bAJtmm8ft81Cj76itQ8WehZA4CHV4I8/0XMulm8/gpnHdGu0lcumwZBw2M0B88bNLyCwoJy7czhs3nsbZfRNqrxvRLYapI1MoWTWLhyJWMvKr9abQ4ul3w3mt3C3Ql/0/gbOEqMGTKEp38OOABzn9x1/Du79kc/lU0uJSKChz8tGCHxnRawv0PMeMPzX853e4mN3ZpXXrYvL3mXGj2J5wx9JmS83Hhjt47IohTBvbjfS8MsIcNsIcNqJD7SREOHxOt2bde1BdDpNfgM4DTVHKJU+ZDbJ8BYiKQvj8TgiLg1sXHXmZmyZIgBCiodRRJkCMu+u4fJzForhuwhgy1AOMWvInXjmzrHYAtDm9C5eCcnHunscp+SGK8LPqV9fdcrCIqa8ux4LmmbA3GMkm5vV+iOt/cQcJEcGw6xew+g2zwK/vRaaq7p5FZibXvqVmHMYRbsqpFOw3mzft+B9DDn7CmzfcyZjevruEJg3qwswgC9PfX8svhiTx0CUDiAptkOEUpDOj4M/Y7EuotCehBt8KGatg/Sw456/N7/Hhj53zwWqn+2kXkrh4Oe8cTuT0yS/gmvcAj5RP5+fEiwi2l5CWvhAyNKx9G36z3Ozf7mXuBtO9NHFAZ3Ng8ePgqoLsLbD5v6Z8fQsGJEXVn7TQFFcNrHgFThlnggO4i1Jeb7qcCtKhU/0JBXz7d1Pe/ur3jnlwgLa7DkKIE2fU7XDeDDOwfBylnHkTOCIJ2zLbvxv2LqUicThfu0YS/sMDsPTZ2sKBmw8Ucu0ry7BbLfxv1HrGFc2FM+7lgql/MMEBoNsZENwJtnxhXm/+rxl/Oe8h89f8noXmuOex+5kw+g5UaRZnOZc0G8TO6ZvIHniI8AAADg5JREFUhr+dz1NXDyWqJtf8Fbz4SfOX8I/Pw3/GYju0Dn7xHI57f4YL/gGjb4fiA3WfdzR2LICuY7AER3DBgM4s3J5N2YApzB43h9dqJjEgZx7dyzfyqp7Mu8l/hYL97J/zf0x8ehEPf/EzWmu01ny18SCj0mKJj3CYQLn2PbPlbXw/syXuEe517rvN86FgH4y8pf7xU91FKde8Xf94+gpY+SqMug1SGmc/x4JkEEI0lNDXZ3dQwAWFwIBLYeMncNETzffFVxbDwfUEn/4HliRchGv9n7ho/gPote+yrctk/rIpiSnWLfwuYTMhPy4xA7Rn/6X+e1iDzA5/274yffYbP4bEQaZK7oIZsON/JrPYvRDCO0Ncb/dPHzPDa8iUxtvBerGVHDDbva5+C2oaVLHtOhYu+0/9CQC9J5ly8es+MN17TXG5wFkCzlJTLbfhX84F6eYv/GFTAZg0sAtv/7SPH7Zl89nWUg5F38bNv3kRiy2Ew9/s5u2f9tIr9jyGbnyZYEsvXl9aTLXLxTUju7Iru5QbPd1Lix43FQPOuNf88fDRDeaf2dBrmm4rmJLwhzeZf59BIabbspOPHRdXzISILtD34vrHO5milKx9B878f2awuub/t3fvUVbW9R7H35+5wHBHBqVBHEFmSpDDIFAOgomKeEnEky0BtRTtqonSbWlrtU5ZrWMezzHLolQoPRVaSEp0jkbokZKGSyBCQIUIAqEICUpiXvieP77PMHuGZ+S2Z/aw9/e11l57nmc/+9m/Hz/W831+l+f3e8uX5e3a22tcLSQCRAhtyeCJfqe49tcw+LLmj9u0yDsw+47ksz0HcM7yG9h23EhGvfa/nLz9dh4BeBvYU+UXlVFT05sgBoyDFT/zi8+WpV5zKmnn6378dZ5fjJ9fAFXnNASD2k/D3Kmw8Lsw4rONz7v3He80XXY//OUx31czyX+/a2/Y84pf2Hv03z89pWUw6FIPEG+8um/xJsA7b5fO8GCzaxM++p19a3zQP2Mk1rp5/l51LuBDkcs7teOnizZSt34H142uQh2OAeDjZ/TjgT9sYOrfP8ITZXU8XPlLbiv/d+79/QaeWLutoXlpx3OwYqbfrXetgJPH+W8/9S1vZipuZpDAsz+HX37Ka2b1VAyDPgwjb/RzAGxf5yPn6ifBbGr4ZHjwcm8K69DD14jfthom/gzaH3io7uGKABFCW1I5ArpVelv8uwWIDU/73ewJp1HRrgNXnt6fry0Qndufxq0jihnf/TmK+54OvQa9610+/c+G0k4wLxl9PuhSf68eC2vm+B3y69uh35kN36mZ5AFs3le8eWrct/0CvuZX/nrtb9CxJ9R+xqcgybxbru8Mb07N5R4IVj8KQz/q/SDzb/Wg+fYefyq9ZoLXNAAW3eMrHFaNgYGXwDEn+tKz3U7wkUBAcZEYe0ovZi72ObYu/JeGIbG9u3fgwU/WUt65PWXrd8P/fIEvD13LOzaAGU8/T+1JPbx5afbtPhllfb9UUZFfzGdO9H+7YZOhZ3Xjf+tlD8CcKb7Gyfi7PXj+81VY9bA3ta38BZRXe+B8YycUlXrtLU31edCltwdm8LIffq3X8FpQBIgQ2pKiIr8A/u4/vfOxy3vSj9v4NFQM2dcMNeWcao7r0p6La3pzXNcyYMzB/V5pGbx3rPc/VI5o6AStHuvv82/195MyAkRpB7hilgexx2+BacmzBCVlHnBqbvPmopJ2h5Z38KfOy6v8br2yFh66Erb/BYZcDrXX77+U7fs/AUvu9eafzPXJh1/T6GJ9/qAKZi7eRN/yjgyoaHzHvW+Ya49rYPlP0GM385XrF1PdqzOD+3SDjQvh2Qf9jr9Lr4Yvvvd8bw6q+76/uld6E12HYwDzEUlVY2DCTxoHxt6nwhmf9+dstizz4dT/2OG1k8zzZyougcvu96aqiho47pTWmRgznoMIoY3Zvg7uHuaz5J5+w/6fv/k63FYJI66Dc2898t9bNRtmTYYL72jcQfrDM2HrM94cNGVZ+nd3v+wXup7VfjFs3zn9uEOx4A544us+mqqkPVw6vXETUpp33oJdm72Td9cWqD4XOjcMq33rnb2M/o//4/LTKrn+rKrmz7N1Bdxzltdext3la4//YKR32l9Xl94vtPMFb45b91sfBvvGTn+e5n0XwCXf9zy0YfEcRAhHk55VPiZ/xUPpAWLzEh9qeeLhPRi1nwEXw7jvwOAJjfdXj/UAkVl7aKrzsXDmF7OTjno1E+Gp272J6LIH/KG0AykuhR79/JWitLiIp744muJ3ebYE8LvzEdd5/8rgCX7h37EOPvpI84MGulfC+6/1V56JYa4htEU1k+CllTD3c96pm2nj06Ci7C3vWlziQymbNlkMuAiQN6W0pm594MZn4JrHDy44HKSS4qKDmwBv9C1+0Z/9KR+FNeTKA9dg8lQEiBDaomFX+wihpTNg2kh4oa7hs40LffRL2UE8fHUkKmpg6iqfR6q1de3d/MigltauE3zoTtj1gi8de943cpOONiACRAhtUXGpPzx29VwfzjrjPF+xbs1cb2LKVvPSgWTxDv6oUj3G+yAmzUw6nQtT9EGE0Jb1HQWfWQh1P4Cl0+Ehf/ir1QJEIRt2da5TkHMRIEJo69p38Y7gUTf5A1KbFvuDayG0sAgQIRwtikt9yoxT/jXXKQkFIvogQgghpIoAEUIIIVUEiBBCCKkiQIQQQkgVASKEEEKqCBAhhBBSRYAIIYSQKgJECCGEVHmzHoSkl4GNR3CKnsD2LCXnaBD5zW+Fll8ovDxnK78nmtmxaR/kTYA4UpKWNrdoRj6K/Oa3QssvFF6eWyO/0cQUQgghVQSIEEIIqSJANLgn1wloZZHf/FZo+YXCy3OL5zf6IEIIIaSKGkQIIYRUESBCCCGkKvgAIel8SX+WtE7SzblOT7ZJOkHSk5JWS/qTpBuT/T0kzZP01+Q9rxbelVQsabmkucl2P0mLknJ+SFK7XKcxmyR1lzRL0lpJaySNyOcyljQ1+f+8StJMSWX5VsaSZkjaJmlVxr7UMpX7TpL3ZyUNzUYaCjpASCoGvgdcAAwEJkkamNtUZd3bwOfNbCBQC1yf5PFmYL6ZVQPzk+18ciOwJmP7W8CdZlYFvAJcm5NUtZy7gMfM7GSgBs97XpaxpOOBKcBwMxsEFAMTyb8y/jFwfpN9zZXpBUB18vokMC0bCSjoAAF8AFhnZuvN7E3gQWB8jtOUVWa21cyWJX+/hl84jsfzeX9y2P3AJblJYfZJ6gN8CLgv2RZwNjArOSTf8tsN+CAwHcDM3jSzneRxGePLJXeQVAJ0BLaSZ2VsZguAvzfZ3VyZjgceMFcHdJdUcaRpKPQAcTywKWN7c7IvL0nqC5wKLAJ6mdnW5KMXgV45SlZL+DbwJWBvsl0O7DSzt5PtfCvnfsDLwI+SZrX7JHUiT8vYzLYAdwAv4IFhF/BH8ruM6zVXpi1yLSv0AFEwJHUGHgZuMrNXMz8zH+ucF+OdJV0EbDOzP+Y6La2oBBgKTDOzU4F/0KQ5Kc/K+Bj8jrkf0BvoxP5NMXmvNcq00APEFuCEjO0+yb68IqkUDw4/NbPZye6X6qugyfu2XKUvy0YCF0vagDcZno23z3dPmiMg/8p5M7DZzBYl27PwgJGvZTwGeN7MXjazt4DZeLnncxnXa65MW+RaVugBYglQnYx+aId3dM3JcZqyKml/nw6sMbP/yvhoDnBV8vdVwKOtnbaWYGa3mFkfM+uLl+cTZnYF8CTwkeSwvMkvgJm9CGyS9L5k1znAavK0jPGmpVpJHZP/3/X5zdsyztBcmc4BPpaMZqoFdmU0RR22gn+SWtKFeJt1MTDDzL6Z4yRllaRRwO+AlTS0yX8Z74f4OVCJT5N+mZk17RA7qkkaDXzBzC6SdBJeo+gBLAeuNLN/5jJ92SRpCN4p3w5YD0zGbwDzsowlfQ2YgI/SWw58HG9zz5syljQTGI1P6/0S8G/AI6SUaRIo78ab2l4HJpvZ0iNOQ6EHiBBCCOkKvYkphBBCMyJAhBBCSBUBIoQQQqoIECGEEFJFgAghhJAqAkQIgKRySc8krxclbcnYXthCv3mqpOmH+d3f5tPsrKFtimGuITQh6avAbjO7o4V/5xfAN8xsxWF89yqgT749txPalqhBhHAAknYn76MlPSXpUUnrJd0m6QpJiyWtlNQ/Oe5YSQ9LWpK8RqacswswuD44SPqqpP+W9Idkrv9PJPsrJC1IajKrJJ2RnGIOMKlV/gFCwSo58CEhhAw1wAB8Gub1wH1m9gH5Qkw3ADfhcz/daWa/l1QJPJ58J9NwYFWTfYPxNTs6Acsl/RoPAo+b2TeT9Us6ApjZK5LaSyo3sx0tktNQ8CJAhHBoltTPcSPpOeA3yf6VwFnJ32OAgT77AQBdJXU2s90Z56nAp+jO9KiZ7QH2SHoSX69kCTAjmXDxETN7JuP4bfhsphEgQouIJqYQDk3m3D57M7b30nDDVQTUmtmQ5HV8k+AAsAcoa7KvaYegJYvGfBCfmfPHkj6W8XlZcp4QWkQEiBCy7zd4cxOwbyK9ptYAVU32jU/WVi7HJ2lbIulE4CUzuxefjG9ock4B7wE2ZD31ISQiQISQfVOA4cni8auBTzc9wMzWAt2Szup6z+JTVtcBXzezv+GBYoWk5fjspXclxw4D6jJWUAsh62KYawg5Imkq8JqZ3XeoQ2sl3QXMMbP5LZnGUNiiBhFC7kyjcZ/GoVgVwSG0tKhBhBBCSBU1iBBCCKkiQIQQQkgVASKEEEKqCBAhhBBSRYAIIYSQ6v8Bu/tKHe27KsUAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "time1 = [res[0] for res in nma1.results]\n", "time2 = [res[0] for res in nma2.results]\n", "linefig, lineax = plt.subplots()\n", "plt.plot(time1, eigenvalues1, label='DCD')\n", "plt.plot(time2, eigenvalues2, label='DCD2')\n", "lineax.set_xlabel('Time (ps)')\n", "lineax.set_ylabel('Eigenvalue')\n", "plt.legend();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "DCD and DCD2 appear to be in similar conformation states." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Using a Gaussian network model with only close contacts" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `MDAnalysis.analysis.gnm.closeContactGNMAnalysis` class provides a version of the analysis where the Kirchhoff contact matrix is generated from close contacts between individual atoms in different residues, whereas the `GNMAnalysis` class generates it directly from all the atoms. In this close contacts class, you can weight the contact matrix by the number of atoms in the residues." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "nma_close = gnm.closeContactGNMAnalysis(u1,\n", " selection='name CA',\n", " cutoff=7.0,\n", " weights='size')\n", "nma_close.run()" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAATpUlEQVR4nO3dfbRldX3f8fcHxgoIomZGQtFxEhc+UKMGRmOCRpA2QaeKVqMhT2Ko07RETWtaJ640spJlMlmpDzEiFQ0BjRqhgsGgKNIASZUgRJQnUZeMCqIToxYxCgLf/rH3bS7jfdh35u5z7pnf+7XWXfecvfc5+zN3nfOZ39nnnN9OVSFJasc+0w4gSZosi1+SGmPxS1JjLH5JaozFL0mNWTftAEOsX7++Nm3aNO0YkjRTrr766q9X1YZdl89E8W/atImrrrpq2jEkaaYk+eJCyz3UI0mNsfglqTEWvyQ1xuKXpMZY/JLUGItfkhpj8UtSYyx+SWqMxS9JjZmJb+5qNmzaduHU9r1j+5ap7VuaNY74JakxFr8kNcbil6TGWPyS1BiLX5IaY/FLUmMsfklqjMUvSY2x+CWpMRa/JDXG4pekxlj8ktQYi1+SGmPxS1JjLH5JaozFL0mNsfglqTEWvyQ1xuKXpMZY/JLUGItfkhpj8UtSYyx+SWrMaMWf5OFJ/jrJDUmuT/KKfvlDklyc5HP97wePlUGS9IPGHPHfDbyyqo4AngKckuQIYBtwSVUdDlzSX5ckTchoxV9Vt1XV3/eXvw3cCBwGnACc3W92NvDcsTJIkn7QuknsJMkm4MeBvwMOqarb+lVfBQ5Z5DZbga0AGzduHD+kZtqmbRdOZb87tm+Zyn6lPTH6m7tJDgTeB/xGVd0+f11VFVAL3a6qzqiqzVW1ecOGDWPHlKRmjFr8Se5HV/rvqqrz+sVfS3Jov/5QYOeYGSRJ9zXmp3oC/ClwY1W9ft6qC4AX95dfDPzlWBkkST9ozGP8RwO/DFyb5Jp+2auB7cA5SU4Gvgi8cMQMkqRdjFb8VfW3QBZZfdxY+5UkLc1v7kpSYyx+SWqMxS9JjbH4JakxFr8kNcbil6TGWPyS1BiLX5IaY/FLUmMsfklqjMUvSY2x+CWpMRa/JDXG4pekxlj8ktQYi1+SGmPxS1JjLH5JaozFL0mNsfglqTEWvyQ1xuKXpMZY/JLUGItfkhpj8UtSYyx+SWqMxS9JjbH4JakxFr8kNcbil6TGWPyS1BiLX5IaY/FLUmMsfklqjMUvSY2x+CWpMaMVf5Izk+xMct28ZacmuTXJNf3Ps8bavyRpYWOO+M8Cjl9g+Ruq6on9zwdH3L8kaQGjFX9VXQ58Y6z7lyTtnnVT2OevJ/kV4CrglVX1zYU2SrIV2AqwcePGCcaThtu07cKp7XvH9i1T27dm26Tf3D0deCTwROA24HWLbVhVZ1TV5qravGHDhknlk6S93qDiT/Jjq7GzqvpaVd1TVfcCbwOevBr3K0kabuiI/y1Jrkzyn5IcvLs7S3LovKvPA65bbFtJ0jgGHeOvqqclORz4VeDqJFcCf1ZVFy92myTvAY4B1ie5BXgNcEySJwIF7AD+w57FlySt1OA3d6vqc0l+m+5N2TcBP54kwKur6rwFtj9xgbv5091OKklaFUOP8T8+yRuAG4FnAM+uqsf2l98wYj5J0iobOuL/E+DtdKP7784trKqv9K8CJEkzYmjxbwG+W1X3ACTZB9ivqv6pqt45WjpJ0qob+qmejwL7z7t+QL9MkjRjhhb/flV1x9yV/vIB40SSJI1paPF/J8mRc1eSHAV8d4ntJUlr1NBj/L8BnJvkK0CAHwZeNFoqSdJohn6B6xNJHgM8ul90U1V9f7xYkqSxrGR2zicBm/rbHJmEqnrHKKkkSaMZVPxJ3kk3q+Y1wD394gIsfkmaMUNH/JuBI6qqxgwjSRrf0E/1XEf3hq4kacYNHfGvB27oZ+W8c25hVT1nlFSSpNEMLf5TxwwhSZqcoR/nvCzJI4DDq+qjSQ4A9h03miRpDEOnZX4p8L+At/aLDgPeP1YoSdJ4hr65ewpwNHA7dCdlAR46VihJ0niGFv+dVXXX3JUk6+g+xy9JmjFDi/+yJK8G9k/yb4BzgQ+MF0uSNJahxb8N+AfgWroTpH8Q8MxbkjSDhn6q517gbf2PJGmGDZ2r52YWOKZfVT+66okkSaNayVw9c/YDfg54yOrHkSSNbdAx/qr6x3k/t1bVG+lOwC5JmjFDD/UcOe/qPnSvAFYyl78kaY0YWt6vm3f5bmAH8MJVTyNJGt3QT/UcO3YQSdJkDD3U81+WWl9Vr1+dOJKksa3kUz1PAi7orz8buBL43BihJEnjGVr8DwOOrKpvAyQ5Fbiwqn5prGCSpHEMLf5DgLvmXb+rX6Y1aNO2C6cdQdIaNrT43wFcmeT8/vpzgbPHiSRJGtPQT/W8NsmHgKf1i15SVZ8cL5YkaSxDZ+cEOAC4var+GLglyY+MlEmSNKKhp158DfAq4Lf6RfcD/nysUJKk8Qwd8T8PeA7wHYCq+gpw0FihJEnjGVr8d1VV0U/NnOQBy90gyZlJdia5bt6yhyS5OMnn+t8P3r3YkqTdNbT4z0nyVuBBSV4KfJTlT8pyFnD8Lsu2AZdU1eHAJf11SdIEDf1Uz//oz7V7O/Bo4Heq6uJlbnN5kk27LD4BOKa/fDZwKd17B5KkCVm2+JPsC3y0n6htybIf4JCquq2//FWW+BJYkq3AVoCNGzfu4W4lSXOWPdRTVfcA9yY5eDV3PP89g0XWn1FVm6tq84YNG1Zz15LUtKHf3L0DuDbJxfSf7AGoqpevcH9fS3JoVd2W5FBg5wpvL0naQ0OL/7z+Z09dALwY2N7//stVuE9J0gosWfxJNlbVl6pqxfPyJHkP3Ru565PcAryGrvDPSXIy8EU8i5ckTdxyI/73A0cCJHlfVT1/6B1X1YmLrDpu6H1Iklbfcm/uZt7lHx0ziCRpMpYr/lrksiRpRi13qOcJSW6nG/nv31+mv15V9cBR00mSVt2SxV9V+04qiKSVmdaZ1nZs3zKV/Wr1rGQ+fknSXsDil6TGWPyS1BiLX5IaY/FLUmMsfklqjMUvSY2x+CWpMRa/JDXG4pekxlj8ktQYi1+SGmPxS1JjLH5JaozFL0mNsfglqTEWvyQ1xuKXpMZY/JLUGItfkhpj8UtSYyx+SWqMxS9JjbH4JakxFr8kNcbil6TGWPyS1BiLX5IaY/FLUmMsfklqjMUvSY2x+CWpMeumsdMkO4BvA/cAd1fV5mnkkKQWTaX4e8dW1denuH9JapKHeiSpMdMa8RfwkSQFvLWqzth1gyRbga0AGzdu3O0dbdp24W7fdk/t2L5lavuWxuJzavZNa8T/1Ko6EngmcEqSn951g6o6o6o2V9XmDRs2TD6hJO2lplL8VXVr/3sncD7w5GnkkKQWTbz4kzwgyUFzl4GfAa6bdA5JatU0jvEfApyfZG7/766qi6aQQ5KaNPHir6ovAE+Y9H4lSR0/zilJjbH4JakxFr8kNWaaUzZI0kzY27605ohfkhpj8UtSYyx+SWqMxS9JjbH4JakxFr8kNcbil6TGWPyS1Bi/wDWiaX7pQ9ob+ZxaHY74JakxFr8kNcbil6TGWPyS1BiLX5IaY/FLUmMsfklqjMUvSY2x+CWpMRa/JDXG4pekxlj8ktQYi1+SGmPxS1JjLH5JaozFL0mNsfglqTEWvyQ1xuKXpMZY/JLUGItfkhpj8UtSYyx+SWrMVIo/yfFJbkry+STbppFBklo18eJPsi9wGvBM4AjgxCRHTDqHJLVqGiP+JwOfr6ovVNVdwF8AJ0whhyQ1ad0U9nkY8OV5128BfmLXjZJsBbb2V+9IctPA+18PfH2PEk7erGU27/hmLfOs5YUZyZw//P8XdyfvIxZaOI3iH6SqzgDOWOntklxVVZtHiDSaWcts3vHNWuZZywuzl3k1807jUM+twMPnXX9Yv0ySNAHTKP5PAIcn+ZEk/wL4eeCCKeSQpCZN/FBPVd2d5NeBDwP7AmdW1fWruIsVHx5aA2Yts3nHN2uZZy0vzF7mVcubqlqt+5IkzQC/uStJjbH4JakxM138SfZN8skkf7XENs9PUkmm/rGt5fImeWGSG5Jcn+Tdk863kKUyJ9mY5K/79Z9O8qxpZJyXZ0eSa5Nck+SqBdYnyZv6qUI+neTIaeTcJdNymX+xz3ptko8lecI0cs7Ls2Teeds9KcndSV4wyXwL5Fg2b5Jj+vXXJ7ls0hkXyLPcY+LgJB9I8qk+80tWuo81+zn+gV4B3Ag8cKGVSQ7qt/m7SYZawqJ5kxwO/BZwdFV9M8lDJx1uEUv9jX8bOKeqTu+n3fggsGmC2RZybFUt9iWXZwKH9z8/AZzOAl8enIKlMt8MPL1/TDyT7g2+aWdeKu/ctCx/CHxkcpGWtGjeJA8C3gIcX1VfWkPPu6X+xqcAN1TVs5NsAG5K8q5+JoRBZnbEn+RhwBbg7Uts9nt0D8DvTSTUEgbkfSlwWlV9E6Cqdk4q22IGZC7++T+Eg4GvTCLXHjgBeEd1rgAelOTQaYdaSlV9bO4xAVxB972Xte5lwPuAqT+GB/gF4Lyq+hKsjefdAAUclCTAgcA3gLtXcgczW/zAG4H/Bty70Mr+ZfzDq+rCiaZa3JJ5gUcBj0ryf5JckeT4yUVb1HKZTwV+KcktdKP9l00o12IK+EiSq/spP3a10HQhh00k2eKWyzzfycCHJpBpKUvmTXIY8Dy6V1NrwXJ/30cBD05yab/Nr0w430KWy/xm4LF0A61rgVdU1WLP0QXN5KGeJP8W2FlVVyc5ZoH1+wCvB06acLQFLZe3t47uEMQxdKO6y5P8WFV9azIp72tg5hOBs6rqdUl+Enhnkset9EG4ip5aVbf2L9cvTvKZqrp8SlmGGpQ5ybF0xf/UiSe8r+XyvhF4VVXd2w1Ip265vOuAo4DjgP2Bjye5oqo+O42wveUy/yxwDfAM4JH9Nn9TVbcP3cGsjviPBp6TZAfd7J7PSPLn89YfBDwOuLTf5inABVN8g3e5vNCNPi+oqu9X1c3AZ+n+I5iWIZlPBs4BqKqPA/vRTSQ1FVV1a/97J3A+3Uyw86256UIGZCbJ4+kOt51QVf842YT3NSDvZuAv+sfNC4C3JHnuREPOMyDvLcCHq+o7/TH1y4GpvoE+IPNL6A5PVVV9nu59oMesdCcz/UM3Qv6rZba5FNg87axL5QWOB87uL6+nOyTxQ9POu0zmDwEn9ZfnXnpmShkfABw07/LH6N6wm7/Nlj5z6AYDV0757zok80bg88BPrYHHwbJ5d9n+LOAFazlv/7i9hG7kfwBwHfC4NZ75dODU/vIhdIOX9SvZz0we6llMkt8FrqqqmZj7Z5e8HwZ+JskNwD3Af60pj+4WskvmVwJvS/Kf6Y5LnlT9o3EKDgHO7w8vrAPeXVUXJfk1gKr6n3TvQzyLrkj/iW7kNE1DMv8O8EN0I2eAu2t6M0oOybuWLJu3qm5MchHwabr3st5eVddNLfGwv/HvAWcluZZuEPOqWuJTVgtxygZJasysHuOXJO0mi1+SGmPxS1JjLH5JaozFL0kjSHJmkp1JVuVTQkkuSvKt7DJhYjqvTfLZJDcmefly92Xxa6+U5J5+dsO5n2398rf3E8pNO98d086g0Z1F9/2c1fJHwC8vsPwkui8mPqaqHkv3hcsl7VWf45fm+W5VPXHXhVX176cRRu2pqsuTbJq/LMkjgdOADXTfJXlpVX1m4P1dssj0Kf8R+IXqp0qpARPNOeJXU/rJuDb3l0/uXx5fmeRtSd7cL9+Q5H1JPtH/HN0vP7V/+X5pki/MvaROsj3JKfP2cWqS30xyYJJLkvx9uvnVT1ggzzHzX7oneXOSk/rLRyW5rJ+s68NZ4zOJapAzgJdV1VHAb9JNCb2nHgm8KMlVST6Ubor3JTni195q/yTXzLv+B1X13rkrSf4l8N+BI4FvA/8b+FS/+o+BN1TV3ybZSPet6sf26x4DHEs3H9RNSU4H3ks3Odlp/TYvpJtI63vA86rq9iTrgSuSXDDk281J7gf8Cd38PP+Q5EXAa4FfXekfQmtDkgOBnwLOnTeB3f37df8O+N0FbnZrVf3sMnd9f+B7VbW5v58zgactdQOLX3urBQ/1zPNk4LKq+gZAknPppugF+NfAEfOenA/sn7QAF1bVncCdSXYCh1TVJ5M8tP/PZAPwzar6cl/ev5/kp+mmAziM7iv5Xx2Q/9F0Ew1e3OfYF7ht0L9ca9U+wLcWOQR5HnDebt7vLfNuez7wZ8vdwOKXftA+wFOq6j4n8OkL+M55i+7hn59D59LNRvnDdK8AAH6R7j+Co6rq+/2Mlfvtsq+7ue8h17n1Aa6vqp/co3+J1oz+ld/NSX6uqs5N94B6fFV9atkbL+39dK9CbwaeTjez75I8xq9WfQJ4epIHJ1kHPH/euo8w76QySZZ65TDnvcDP05X/uf2yg+nOafD9dPPpP2KB232R7tXF/dOdBvC4fvlNwIZ05zkgyf2S/Kvh/zxNW5L3AB8HHp3kliQn0w0GTk7yKeB6urPCDb2/v6F7bB3X39/cIaDtwPP7Sdv+AFj2AwyO+LW32vUY/0VVtW3uSnUnuvh94Eq6U9d9Bvi//eqXA6cl+TTdc+Ry4NeW2llVXZ/uHM+3VtXcIZl3AR/on5BX9fvY9XZfTnIO3XTANwOf7Jffle5E5W9KcnCf4410ZaEZUFUnLrJqtz7iWVULHrev7mRNW1ZyX87OqWYlObCq7uhH/OcDZ1bV+dPOJY3NQz1q2an9q4K50fb7p5xHmghH/JLUGEf8ktQYi1+SGmPxS1JjLH5JaozFL0mN+X89LmisZvRN6wAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "eigenvalues_close = [res[1] for res in nma_close.results]\n", "\n", "plt.hist(eigenvalues_close)\n", "plt.xlabel('Eigenvalue')\n", "plt.ylabel('Frequency');" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAERCAYAAACU1LsdAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOy9eZgk11nm+56MyIhca19671Yvau1rWxKyjS3skfeFGWAwGAZf20KD8YULA1wzwIAZZpgBZsxybeNr4w2wMcY2WMbyhhctlqUW1uaW1N3qvau71q6q3CMy4swfJ07kicgTmZFVmVlLn9/z9NO1REZGRmWe73zf+y2EUgqFQqFQKAAgsdYXoFAoFIr1gzIKCoVCofBRRkGhUCgUPsooKBQKhcJHGQWFQqFQ+CijoFAoFAqfDWkUCCF/RQiZIYQ806Xz3U8IWSSE3Bf6OSGE/AEh5Cgh5FlCyP/djedTKBSK9cqGNAoAPgbg1V083x8B+BnJz38OwE4AV1FKrwbw6S4+p0KhUKw7NqRRoJR+B8CC+DNCyD5vx/84IeQBQshVHZzvGwAKkl/9RwDvpZS63nEzq7luhUKhWO9sSKMQwYcAvJtSeiuA/wTg/V045z4A/54QcpgQ8mVCyIEunFOhUCjWLfpaX0A3IITkANwJ4O8JIfzHpve7fwvgvZKHnaeUvqrNqU0AVUrpIe88fwXgpd25aoVCoVh/bAqjAObxLFJKbwr/glL6OQCfW+F5zwmP/TyAj67wPAqFQrEh2BThI0rpMoCThJAfB/ysoRu7cOovALjL+/plAI524ZwKhUKxbiEbsUsqIeRTAF4OYAzANID/AuBfAHwAwFYASQCfppTKwkay8z0A4CoAOQDzAN5OKf0KIWQIwN8A2AWgCOBeSumT3X01CoVCsX7YkEZBoVAoFL1hU4SPFAqFQtEdNpzQPDY2Rvfs2bPWl6FQKBQbiscff3yOUjre7rgNZxT27NmDw4cPr/VlKBQKxYaCEHI6znEqfKRQKBQKH2UUFAqFQuGjjIJCoVAofJRRUCgUCoVPT4VmQsgpsO6jDoA6pfRQ6PeDAP4arDhMB/DHlFLVSkKhUCjWiH5kH91FKZ2L+N27AByhlL6BEDIO4HlCyN9QSq0+XJdCoVAoQqx1+IgCyBPW2jQHNiOhvraXpFAoFJcvvTYKFMBXvcE390h+/xcArgYwBeBpAL/EB9qIEELu8WYaHJ6dne3tFSsUinXFwy/M4di0bAaWohf02ii8hFJ6C4DXAHgXIeSHQ79/FYAnAGwDcBOAvyCEDIRPQin9EKX0EKX00Ph424I8hUKxifjNzz2NP/zyc2t9GZcNPTUKlNLz3v8zYPMIbgsd8jYAn6OM4wBOgnUrVSgUCgBA1XZx5MLyWl/GZUPPjAIhJEsIyfOvAdwN4JnQYWcAvMI7ZhLAQQAnenVNCoVi42E5Li4sVbFQUvkn/aCXnsIkgAcJIU8CeBTAlyil9xNC7iWE3Osd8/sA7iSEPA3gGwB+o0WmkkKhuAyx60xmfFZ5C32hZymplNITAJqmn1FKPyh8PQXmQSgUissMq+7iM4fP4i237YKWINHHOcwo/GBqCS/eP9avy7tsWeuUVIVCcZnytSPT+K0vPIPHT19qeZztGYUjU8pT6AfKKCgUijXhuYtskV8o1SKPcVwK1xsOqcTm/qCMgkKhWBOeu8hqDxbLduQxlqcnpJIJvDBbQtV2+nJtlzPKKCgUijXhec8oXGplFLzQ0fXbB+G4FEdVEVvPUUZBoVD0nVKtjjMLZQDAYjk61ZTrCTftHAIA/EDpCj1HGQWFQtF3xB3/pRhGYe94DnlTV2JzH1BGQaFQ9B0eOsqZektNwa4zldnQErh624ASm/uAMgoKhaLvPHexgIyh4ZqtA62FZs9TSOoJXLN1AM9eWIbL05EiWChZ+K/3HUGpphourwRlFBQKRd95/mIBBybzGMkaLcNHPPvI0Aiu2TaAsuXg1Hyp5bkffmEOH37wJP7usbNdvebLBWUUFApFX6GU4vnpAq6azGM4m2yZfcQ1BcPzFID29QrFKvMQPvHdU229CkUzyigoFIq+MlusYaFk4eCWPAbTBpYqFiiVL97cKCS1BK6czCOpkbZic9ELG52aL+Pbx9T8lU5RRkGhUPQVLjJftSWP4UwStkNRsuRFaZZgFAw9gf0T+baeQsHzFMbzJj7+8KnuXfhlgjIKCoWir3CjcHBLHsMZAwBwKaIttu0wDyKpsaVqx3AaF5eqLc9frNWRM3X89O278K3nZ3FyrrUGoQiijIJCoegrz18sYCxnYjRnYiiTBAAsVeS6gu0LzWypGkonI4/lFKo2cqaOn7p9F5IawSe+e6pr1345oIyCQqHoK89PF3DVljwAYIh7ChEZSI2UVNZaezCGUSjW6sildEzkU3jt9Vvx2cPnVHpqByijoFAo+gbvX3TQMwrDnqcQlYHkZx9xTyGTRNly/FRVGYUqCx8BwFtu24VCrY4HjqnZXXFRRkGhUPSNMwtlVG3XNwrcU4jqf8QXf64pDKZbh5sA5inkU8wobB9KAwCWq629C0UDZRQUCkXfeN6boXBwkhkFvshHVTVzodnQPaPgGZGlSnTBW1HwFPj/KnwUH2UUFApF3+CdUa8YzwJgi33O1CM1BbFOAWhvRIBG9hEAZEwNAFCOSHlVNNOzGc0AQAg5BaAAwAFQp5QekhzzcgDvA5AEMEcpfVkvr0mhUKwdc0ULpp5A3mwsPUOZZOQi3wgfMaF5KE74qFpHPsWOM3UNSY34BW2K9vTUKHjcRSmVqjyEkCEA7wfwakrpGULIRB+uR6FQrBFzhRrGciYIIf7PhjPR/Y+skKfAU1ijjIjrUhQtln3EyRg6ysooxGatw0c/BeBzlNIzAEApnVnj61EoFD1ktljDWM4I/KyVpxDOPmonNJdtB5Qi4IlkDS2yYlrRTK+NAgXwVULI44SQeyS/vxLAMCHkW94xP9vj61EoFGvIXNHCWM4M/GwoY0RmH9mOCz1BkEgwzyKfSoIQYDHCKBS8LCPRU8iauhKaO6DX4aOXUErPe2GhrxFCnqOUfif0/LcCeAWANIDvEkIeoZQeFU/iGZR7AGDXrl09vmSFQtEr5os13LB9MPCz4Ux0p1TboX7oCAC0BEHe1LEcYRR4h9Sc4ClkTF15Ch3QU0+BUnre+38GwOcB3BY65ByAr1BKS57u8B0AN0rO8yFK6SFK6aHx8fFeXrJCoegRrksxX7Iwlg+Fj9JJLFdtOJI211bd9UVm//gWnkXB8whETyFnaspT6ICeGQVCSJYQkudfA7gbwDOhw/4RwEsIITohJAPgdgDP9uqaFArF2rFYYQu/LHxEKaS7f8tx/RoFTqtWF9xTEDWFjKHCR53Qy/DRJIDPe1kGOoC/pZTeTwi5FwAopR+klD5LCLkfwFMAXAAfppSGDYdCodgEzBVrANBkFIazvNWFheFs0Iuw664vMnOGMslITaEo8RSY0KyMQlx6ZhQopScgDwV9MPT9HwH4o15dh0KhWB/MFZhRGG3KPuJN8ZoXettxkQx5CgPpJM5fqkifQ6YpZE0d5ZrSFOKy1impCoXiMmHW8xTGw+EjP820WScIC838+KjwEdcU8mbS/1nW1JWn0AHKKCgUir4wX2SLflP4yB+0I9cUwkZhMM3CR7IRnr6nEAgf6ajaLupOdGdVRQNlFBQKRV+YK9agJ4hfgMYZbjFTwaq7MJqyj5JwXPkIz2LNRsbQoCUaj8ny/ke2CiHFQRkFhULRF+aKNYzmDL8QjZNP6UgQeZWyLfEUhtLR7bbFZnicrOqU2hHKKCgUir4gq2YGgITnPcg8BVuSkjrQotXFcjXY9wgAMgbzFEpKbI6FMgqrxHUpTswW1/oyFIp1D/MUmo0CwJviyTQFidDM5zpLji9W64EaBYBpCoDyFOKijMIq+coPLuKV/+vbmF6urvWlKBTrGtYh1ZD+jjXFk3gKdbnQDMg9BT6fWcQPH6kMpFgoo7BKzi9W4FJgZrm21peiUKxbKKWYK1lN6agc1rpCrikYerPQDMib4olT1zi+0KzCR7FQRmGV8NJ8NQNWoYimUKvDqrtSTQGIbp8dlZIKyGcqMKE5mN2kPIXOUEZhlXAXNqpro0KhaFQzh5vhcaIG7cjaXKSTGgwtIQ0fFao28uHwka8pKE8hDsoorJIl5SkoNiG24+LJs4tdO99cROEaZyidRNlyUKsHF27LoU1tLgghGEgnmyqgKaXSlFQ+p1kJzfFQRmGVNDwF9YZTbB6+/MxFvPn9D3UtgYI3wxvNRhgFrxFeOKPIdpo9BYCFm8KeQsV24FJEewoqfBQLZRRWifIUFJuR+WINlMqrjFeC3yE1MnzEO6U2G4XwPAXAa3UROlbW4gJgg3lSyYTyFGKijMIqUZqCYiPw7z7wMD760MnYx5e9FhLdWkjnihYIAUYyESmpaXmrC1lFMzu+2VPwB+yYzc2fc2r6WmyUUVglS17YaLmqdiGK9QmlFE+eXcTR6ULsx1S8BbTYJXF2rljDSMaALlngASHNVNj9uy6VdkkFgEFJtlKBD9hJNRuFjKGjrDyFWCijsAoopY2UVOUpKNYpVdtF3aX+7j8OPP7eNU+hUIsUmYFGmqkYhrVd1tU03OaCHx/+zDVmKSSbjs+aetcM3GZHGYVVULVdWF47XqUpKNYrBe+9WenAKDQ8hW6Fj2qRegIApJIsQ6gmdDK1HdYaWyo0pw0UanXYQjvsYo29Tln4KGtoKCuhORbKKKwC0RCo7CPFeoWHNisdtI7uhaYQlXkEAGmvaZ14jXadLfhyoZkt/KK30DJ8ZKo5zXFRRmEVcKErndSUp6BYt/DdfieeAt9VdxJyasVcsXX4KOWFiCpWY+fPvYBwnQLQGOEpis3FlkKzFltofmG2iIePz8U6djOijMIq4G/InSNppSko1i1++GgFnkI3wkdlq46y5bQMH+laAkmNBK6x5nsKck0BCPY/ikpJBZjQHNdT+OOvPI9f/NT3Yx27GempUSCEnCKEPE0IeYIQcrjFcS8ihNQJIT/Wy+vpNrzQZudwBiXLUeP+FOsSHlbpzFPoXvgoagxnmFRSQzWgKbDPkykTmjPNnVKLtTpSyYTUiOQ6CB89P13AQsmStua+HOiHp3AXpfQmSukh2S8JIRqA/wHgq324lq7S8BQyABofPoViPbEyT4G9l7vhKcx6hWtRHVI56SajwITmVp6CuHAXJM3wOBlDQ9lypHOdRWp1B6fnywCAU/OllsduVtZD+OjdAP4BwMxaX0incKOwYzgNQGUgKdYnfLPSiT7QTU/Bb4bXzigYWlBodqLDR0OSmQrFal0qMgMsJbXuUj8kFcWJ2RIclxkOZRR6AwXwVULI44SQe8K/JIRsB/CjAD7Q6iSEkHsIIYcJIYdnZ2d7dKmd0zAKzFNQGUiK9chKso8qvlFYvdDMm+GNRgzY4aSTWiDEZTnR2UcDkvbZhaotFZkBlpIKtDeMYoHfqblyy2M3K702Ci+hlN4C4DUA3kUI+eHQ798H4DcopS3NN6X0Q5TSQ5TSQ+Pj47261o5ZqtjIm7rft0V5Cor1CA8fWXXX3wW3wy9e60Juv98Mr41RMJMaqsJOnqekyuoUkloCOVNv0hSijELGjDeS8+h0AXqCYCxn4vRl6inI72CXoJSe9/6fIYR8HsBtAL4jHHIIwKcJIQAwBuC1hJA6pfQLvbyubrFcsTGQTvq7FpWBpFiPiFpXxXYiF06O61JUbbYgdyN8NFuoYTCdhKlrLY9LJxOoyjwFidAMeE3xhPbZhWrd1/fC5GIO2jk6XcSesSzGcyZOXqZGoWeeAiEkSwjJ868B3A3gGfEYSukVlNI9lNI9AD4L4Bc2ikEAmKcwKBoF5Sko1iEF4X0ZJwNJDDO1Cx9dWKoExGEZU4sVbBtKt33edFKuKcg8BYAZBVFoLtbqyEd5Cka8mQrHpgu4cjKHPWNZX3C+3Ohl+GgSwIOEkCcBPArgS5TS+wkh9xJC7u3h8/YN3yikeHWl0hQU64+ApxDDKPDdtJ4gbbOP3vDnD+L933qh5THnFyvYPpRq+7xhodmqR2cfAc0zFYq1urRGARA8hRZGrmo7OL1QxoGJPPaMZlha6mXo/fcsfEQpPQHgRsnPPxhx/M/16lp6xVLFxr7xHHKmjgRRnoJifRIOH7WDG46xnIm5Yg2UUngh3gCOSzFXtPCD80stzze1WMFtV4y0fd6ULq9TMPTm5waYp3BspgjAm7pWbaEpeIN2WvU/Oj5TBKXAlZN56J64fXq+hBt2DLW99s3EekhJ3bAsVWwMpHV/PKDSFBTrkULV9r3ZOE3heIbOxIDZMo2Tj848MRcdey/W6liu1mOFj1KG3ChEeQrDWQPTy1XYjotanXWCjfIUst5IzladUnnm0cEtOewZzQIATrZ4bZsVZRRWAQ8fAcBAKqlmKijWJYVqHRMDLHwTx1PghoMXm0WlcXIx+sxCGVaE4biwWAEAbB2MET4Kp6S2aHMBAK+8egKFah1f+cFFoRmevHgta7Y3ikeni0hqBLtHs9g9ygTry1FXUEZhhVRtB7W62zAKaV15CpuA84sVvPZPH8BMoTuzidcDhWodE3m2wLcThYGGERj3HhMlznID47gUZxbkO+rznlHY3oHQzKuO23kKL7tyArtGMvj4w6d87SNKaOZzmltpJMemC9g7lkNSSyCV1LB1MIVTylNQxIXrB0FPQRmFjc7RiwUcubCMY9PFtb6UrlC1HViOi0nPU4hT1eyHjzyjELWQigbmhVn54jm1yIxrrOwjQ4NLG6moFp+nEJGSqiUIfvaHduOxU5fw6Ml5APIOqQCQSiaQIEC5VfhopoADkzn/+z2j2cuyqlkZhRXCvYIB0Sio7KMND1+QutUyeq3hYRW+wMfJPuIhlrE2nkLQKMiN6NRiBVqC+M/fCt74joel2qWkAsCP37oT6aSGD3gZUFGaAiEEWUOPNHClWh1nFyq4cjLv/2zPWAanVPhIEReeqhYIHylPYcPD49jdntJVqNr4X187Ghl77xW8RqEzTcELH+XiewonojyFpQom82bkbGYRPmiHn7fVkB3OYCaJN9+83V+8WxXmZU098u963MtiujLkKVyOaanKKKyQJqOQUtlHmwG+O+2kzXQcHjo+jz/7xjE8dmqhq+dtx0o8hYqQfQRE5/bzHb2hJ3CihacQJ3QEME1BfH7bcUEICxO14j/cudv/OqohHgBkWgza4ZlHBwRPYbeXgXS5tbtQRmGFNHsKSTVTYRPQ8BS6axT89M2IxbNXcKPAReM4r4sbAd7VNKo1BN/RH5zM44XZkrQt9dRitXOj4J235rhIaglpjYTIVVsGcMdeVgfR0lNoMWjn2EwRhpbAbqFNxhVjzChcbiEkZRRWCC+vb3gK7M2oZipsbHxPoYOOonHgxiZKkO0VBSEhwtQT8bKP7DpMPYEBL72zXfbRtdsGsFSxsVCyAr93XYoLS/E9hVQyHD6iMGOEnQDg1151EG++aRuGM9FN97KmFik0H5suYO94NhDm2uUZiMstA6mnDfE2M0ueqDwgeAoAy0oazrbuBqlYv/CMl25rClzAjhJke4U4zD4TaiMRRcVykDE0P7c/Wmhmr+nabQMAWBHbqDAzYa5Yg+3QWC0ugIZR4NdoO25kM7wwt+4ewa27W1dNZw0dF5bkqcZzRaupliJteGmpKnykiMNSxUbW0Pwcar6rUhlIG5tehY98T2Gmv0aBJz/kU0mkk1rs8FHG0GHoCRhaIrIKmC/e13hGIfzazvuFazHDR2Gh2XFbisyd0kpoZn2Tmgvfdo9mLjtPQRmFFSJWMwNQnVI3Cb0SmrlRmFqqdqUddVy4p5Az9aaGc1FU7LrfVTRjapHXW/POtW88x8Tm0OLJd+WdC828TsGNLFxbCVlTizRwhYi+SVeMZZWmoIjHkjdLgTOQ5p1SlVHYyNg9qlMQU1H72U+HL3ZagjCjENtTYAt01tDbCs0ZQ8cVo9kmEX2qg2pmgBWYiee1HdqyRqFTskYrT8GWZi7tG89hoWRhennzVLi3I9YdJ4RMEkI+Qgj5svf9NYSQt/f20tY3y2FPIaU8hc1Az8JHQlZaP3WFQrWx2GWSeuyUVB7KyZnRGTsV24GWIEhqBPsmsk0i+vnFCrKG5m+Y2hHOPrLqTlc9hYypo2w5cEPT52zHRdV2pS0yeHfXR07Mr/r5//grz+MdH39s1efpNXHv+McAfAXANu/7owB+uRcXtFGIDB8pTWFDY/nZR10Wmusu9ARBgvRXVygIw+xTMcNHZbvu9wrKmlrLOoWUzlJG947lmhrj8RqFdimlnFSTpkAjW1ysBH9Oc+geFHmITeIpXLN1ADlTx/dOrr6+5Imzi/j+mcVVn6fXxL3jY5TSzwBwAYBSWgewOfoArJCwUcgampqpsAnoVfioVneRNjTsGsn0NS21ULP9zqHpZCJem4taw1PImtGtIaq242cM7R3Peo3xGvH3qcUqtsYMHbHray5e67bQDADl0Ovhr0+mKehaAof2DON7HXgKH33oJJ4427z4zxVrWChbsedkrxVxjUKJEDIKgAIAIeQOAK0na2xywkZBzVTYHPCdbteFZseFqSewbzzX5/BRw1PIGHrsNheZmOEjbhT2jbP2EOJru7AUb+IaJ6kloCeIED7qvtAMNLftENN2Zdx+xShemC1htlBr+xwLJQvvve8IPvW9M02/mytaoBS4VLYkj1w/xL3jvwLgnwDsI4Q8BOATAN7ds6ta51h1FxXbCRgFQM1U2AzYfp1C94VmQ0tg30QOJ+ZKq9otnl+s4O8Pn411LDMK7H2aipmSWrbq/qSyTIsq4Jrt+uLw3nFW/ct7IFVtB3NFC9tipqNyUkkt0BCvu+EjPlMhFD7yPQX5LIbbvWrpOC1KvvvCPChlXoGI41IslNjP5oubwChQSv8VwMsA3Ang5wFcSyl9qpcXtp7x22ZnQkZBzVTY8PSqS6pVZwvcvvEsrLqL85cqKz7XZw+fw6999qlY3kxAaA5NNosi6ClE9wsSw0f5VBITeRPPTLEAQqfpqJxUsqF7dD8lVT5ToVhjn9moDqvXbx9ExtBihZAePD4HAJgNGYXFsgW+D5gvtvc41pJYaQGEkJ8N/egWQggopZ/owTWte8J9jzhqpsLGpxE+6r7QbHjhI4CFWXaNZto8Sg5vXVGo2X7sP4plIXzEitfqkTOX+XXWXdpISfXCR7LHVGzH1wEA4PU3bMNfPXQSP3brjJ9K2qlRSBuJQJuLbqak8tcUTksVazlkJLUEbt09HEtsfogbhVCoaV5oATJX2gSeAoAXCf9eCuB3Abyx3YMIIacIIU8TQp4ghByW/P6nCSFPecc8TAi5sYNrXzO4URhISYyCyj7a0PhCszABrBtYTrNRWCm8bqDYJlRZqzuw6q7/Pg0PsZHBvY+Mn32kR85pFj0FAPj1Vx/EVVvy+E+fedIXWrd1oCkAwZGcnbS5iEPOb9shDx+16rB6+xUjeO5iAZdaLOhn5ss4s1DGUCaJuWItkPo6JxiJuRjaxFoSN3z0buHfOwHcAiDX7nEed1FKb6KUHpL87iSAl1FKrwfw+wA+FPOca8pSaMAOR81U2Phwo0ApIgfWrwSuKQxnDYxkjVUZBV6V22q0JNAsoPJdfdWKfl3c4IhCMyDvf1Sx3YBRSCU1/PlbbkbJquN/f+0oAGBLjNnMIqmkhmpdDB91L/soE/Faim08BQC4fe8oAODRFrrCQy8wL+F112+F7dDAHAbRO5gvbQKjIKEE4IrVPjml9GFK6SXv20cA7FjtOfvBcqvwkdIUNjRinn03dQUePgKAfeNZvDCz8rRUvqi18xSajIKfpx/9OP6axZRU9pzN96JmO77QzDkwmcfvvP5a1F2K8bwJU28d3gqTCnkK3Qwf5TzvJ6yRFGt1JEjDEMq4YccgTD2B752INgoPHp/D5ICJOzwDIuoK3Dsw9cS6F5rjagpfhJeOCmZIrgHwmRgPpQC+SgihAP6SUtrKE3g7gC9HPP89AO4BgF27dsW55J4SqSkIMxXiTJpSrD94l1SAxZ5HutTxtua4GDTY+2XfeA5fOzK94nNxD6FdphvXHnhWDV/0WgnU/Hd+8ZohT+MEmsNHnLfcthP/eubSijKs0kkNi17Kpu3QnqSkhjduvBVIqyI7U9dwy65hfO+kXGx2XYqHj8/hrqsm/IFGM8s1f7znfKkGPUGwZzSLuc1gFAD8sfB1HcBpSum5GI97CaX0PCFkAsDXCCHPUUq/Ez6IEHIXmFF4iewknjH5EAAcOnRozSs/wrMUOOJMBdU+m+G4FB9/+BTectuutqLoekD0FLpZq8DDRwAzCp8uncWlkrWi94nvKXQYPuILeCsPKBw+8gu+JMJ7WGjmEELwxz++MnkwndRwoUd1CrqWQD6l+0aHI6bttuL2vSP4028ca6pRAoAjF5ZxqWzjJfvH/IFGs8VGv6S5goXRnIHxvLk5wkeU0m8L/x6KaRBAKT3v/T8D4PMAbgsfQwi5AcCHAbyJUrr6BiOr5Le/8Ay+/PSFlscslC3kTL0ph1p1Sm3meyfm8d77juCBY7NrfSmxsL0iM6Db4SPHP+++CS+nf25lukIjfNT6fVbw22Y3UlIBtExLrUSEj+SegtsUPlotaaNRp2B1uU4BAIYzBhZDnkKxZrfUEzgHJ/OgtNHoT4RnHb1YNAqCoDxfqmE0a2I0Z6z78FHLO04IKRBCliX/CoSQ5TaPzRJC8vxrAHcDeCZ0zC4AnwPwM5TSo6t7KauHUoq/e+wsvvn8TMvjFkqWNKxwuc1U+PADJ/D2j7Vu8HXM6/PT7UlmvcJ2XH8X2FWjICxwfKLXuRXWKsQVmnl4yc8+CjWck8FfMzcGURk7lFJU6/Lw0WpIJROoeJlfTFPontAMAMOZJC6Vw0ahHlmjINLKa3rw+BwOTOQwOZBCztSRSiYCRmG2yDyF0ay57usUWhoFSmmeUjog+ZenlA60OfckgAcJIU8CeBTAlyil9xNC7iWE3Osd8zsARgG8PypttZ8Ua3VYjhvZAIwTZRRGcuxnv/bZJ8k0ye4AACAASURBVPGZx87GKhTayPzrmUttqzyPe0Zho9wLq+5iyCtK7GZTPDF8xGP87d5nUXBPodDD8BE3IDwOH87YqdVdUIoeGAUNVcuB41JQiq6GjwBgKGM0hY+KEbMUwjTaZATvn1V38dipBbx4/xgAFj4bz5tBT6FYw3jOxFjeQMlyut5GpZt0dMcJIROEkF38X6tjKaUnKKU3ev+upZT+gffzD1JKP+h9/Q5K6bCXshqVtto3+IzZdjuw+aKFUYlRuGnHEP7nj90AAPj1f3gKL/7Df8FzF1s6VBuapYqNolfYFAU3Cuv5QyDSK09B7PgZtdDGwXGpv9Nvn33EhebOw0fhlNTwZ6LmhXi6bRTSXkUzbzfSzToFABjKJJt6DxU69BTCf7eFkoWq7eLAZCNLfzxn+tlHlFLMFWsYzRkYy7LQUrgNxnoi7jyFNxJCjoHVFXwbwClEZAptZHhWQLv5vFGeQiJB8BOHduLLv/RS/O07bsdC2cJXnll5lsl6Z7Fsw6XNKX4ix2d5+Kh7Of+9xKq7GEyzv22vUlJ5Zk+7nb4McXEuxEhJzRianwnnp6S2eF3h8BEvYgsvhNwwyYTm1ZBOaqi71P8MdttTGM4YWAyHj6p16SyFMNmIeyHrsip6CmXLQdV2MZZjmgIQrHBeb8S9478P4A4ARymlVwB4BVhdwaai4SlEf2gopcwo5KKzRgghuHP/GLYOpHB6Ew/95qm5hQjBc6ls+x+MjaIpWIKn0PXsI88oJBIEWSN6zGUrxMe0zz4KThPLJNnXrV5X2aqDEPiiOJ/THDb83NvottDMPQ+uh3RbaB7KJFGo1lEXqrqLtbjho5UZBe4VjOZMjOaYp7CedYW4d9z2MoMShJAEpfSbANY01NMLeBfDVh/Wgqc7yMJHYXaPZnHqsjAK8vt1fLbgf13bIEbBdmjXw0eUUiY0C7vebIuW1K0IGIUYnoKYapky2PO3E5ozSS2Qs5+VzGnmVcdd1xSMYC1B94Vm9rnlGUiOS1G2nJjhIy/sF3pf8HuTFYzCRD6FS2UbVt31IxBjOcNfN9ZzBlJco7BICMkB+A6AvyGE/ClYVfOmgrt0rcJHC94fc9SLDbZiz1gGpzfp0G/Hpb4xiPIUjgsTxjaCp+C4FI5L/fGR3WqKx3sNibveXMTwmv/+5Wfxzy1SovljDC0RS2gWPQVDSyBB2nkKjt8OgiMzYH7qag/CR0Ajrbv7QjMzklxsbjVgJ4ypa0hqJLanALBUVO4VjOVMjHmewtw6rlWIe8ffBKAM4P8BcD+AFwC8oVcXtVZw693KLeeGo1X4iLN7NIv5krUp6xbEqtCoytrjM0UYegKTA+aGEJp53yNT17yOot25Zl4QZ+pBT0H2Pvv0o2fxj0+cjzwXz1iaGDD9ls9RsPBRcBBUu0E7bJZCcKGXGTBeS2B2u06BG4VKbzSFIc9T4Gmp/HWFm1tGIZsvUZIZhVyjVoF7CqM5A2lDQ9bQNoWn8PMAtlJK65TSj1NK/2w9FJp1G64pVG03EHOUHRMnfLTHa418ZhN6C2Kzr8jw0UwRe8eyyBo6ql1sLtcr+I4+qRFkDK1plu+Kz1uXewrhxYVSimKtjqnFKqLgi9iWgVTM8FFwB9xu0E7Zaq5Szpq6n6rK6Vn4yDMyhR55CsO+pxDUw+KEjwDv7xYxpCcr8RRmCw1PgUcXRnPru1Yh7h3Pg/UweoAQ8ouEkMleXtRaIWYERGXUcN0hTk+c3aOscnUz6gqLAaMQET6aLeLAZD7Q5Gw9YwuLd9ro3jX74SMt7CkEz1/28vNlFbMcbkgmB1MoVFunAy9X637rFU67QTsVywksbvwx4Wut9jh8xDcdht4bTYGnpcbpkCoi01dahY9mCjXMFWsYSDU6IIzmjHXd/yhum4vfo5ReC+BdALYC+DYh5Os9vbI1YEGI80XpCvOl+JoCr1zdjLpCO0+hYjk4d6mC/eM5pJKJDVG8Ji7eGUNrm5os40tPXcBvfv7p4HmlnkLz4sLv43zJijRIfMe+ZSAVOeegcT67qacPH7QTRSkifBQedt9zobnKhebunj+sKXBdJq6nkDGaQ2klr8uqmInFU09nCzXMlSyM5RvrxWjW3Ph1CgIzAC4CmAcw0f3LWVvmi5b/h43KDFkoWkgntVjN3bKmjvG82ZSW2s3hLWvFUhtP4YXZIigF9k/k2K57nRmFzz5+Dj/zke8FfmbXvYIpLYG0oa9IU/jO0Vl88YmpwM9kRkGmKYj3cWpJ7i2I4SP2GPn71Kq7qNXdpvz7lKG1rBmpRIWPIjSFngnNvqbQXU8hZ+rQE6ShKfCq75ieQs5sfl+Uak5Tl1VT1zCUSTKjUKj5RWsAy0La8HUKhJBfIIR8C8A3wNpSvJNSekMvL6zfUEoxX7Kwc5jt7qNqFeYjCtei2DOawamQp/BzH30M7/3ikZVf7DpgydtpESJfmPgQmf0TOaST8WYD95PvvjCP774QlMV8TUFPILPCkFe17vi7aA7fzYvho1xKYhSE76NmOJdqdWgJ4u9Eo5IioqaJZbw2ElGUJeEjmdDM703XG+KFs4+6XKdACMFQJulrCsUOPQVZ+KgQ0SZjPMdqFeZLFsbyjTVjNGdgoWQFJrOtJ+Le8Z0AftlrV/G7lNKNvaJJKFlsdCEP+UR5CvMlC2MxMo84u0ezAU9huWrjgWOzePr84uoueI3hnsJE3pQaheMzRSQIS8sVh7F3ysxytW2DwhWdt1BF3aWBhAKrHg4fdX7NFYu1aLDF88pSUg0dVt0NHCfexyhdoVRzkDU0PywUJTb7O+Bw+MjQ2gzZqTd5wVlTQ8kKjiftndDMzsfvRTeH7HDE/kcdawoR4aOwIQW8ArYi0xRGA56CCccNTmZbT8TVFN4D4GlCyLa4vY82Grz+YGcbo7BQqnXsKUwv1/w47uFTC3Ap1rXQFIelio2MoWEka0rDR8dnitg9moWpa36Ts5XwkQdP4m0ffayl+LoSppdZho8Yk7f9xZusOOTFHyN6RlHhIyD4PguEjyJeL6++5R5AISItdTkiq6adgM6L10Sypg4npF/wv6fZ5Z18I3zUm+wjgHdKDWoKvIVFO7Ky8JEl7500njdxYbGCxbLte3YAGlXN67RWIW746BcBTAP4GoAvef/u6+F19R1eTOJ7ChFi3ELRwkgMkZnDM5DOLLAQEh/nN7vOh3e3Y7HMBo3kU7q0TuHYTNEfUp+O8BR+8W//FX9/+GzL5+GZW60KulbC9DK7/4HF209JXbnQzBvFVe1mD8QMpaQCQe+Af53UCM5FGYUq25XKHh84riaPladbhMVcr9leU/Ga0dwUr1pnsxRaTStbCbzqupF91CtPoaEp5EwdiUS815ExtSZPIapNxnjOxNQS23zwojUAGPM2lXxj6LoUx2cKTY9fK+Le8V8GcNALH13v/dtUmgL3FLhRkGkKXHcY7SB8tIenpc4xo/DIiXnv/PUNkaYZBZ8+NZDSmxYm23Fxaq6E/ROeURAGp4h849kZfPmZiy2f58wCWxy/+FT3jELVdvxFR6yfsAPho5UJzS09BSGThu8sxc0H9xT2T+Sjw0ferpR7Cp2GjzItPKBq3QGlzbOK/TkCwmciahTnauFV142K5u4aHQAYSouagt2ku7RCFvYrVutST2M8HxSXOdxT4BlIH3nwJO7+399ZNxvFuEbhLIClXl7IWsOL0naNRoePypaDWt3tKHzEz3d6voRC1cbT55ewdZBljqzntLR2cKOQTyWbwken58uouxQHPKPANQUxJs1bQD9/MXqHRCnF2YUyMoaGJ88udq0IcGa5cd9rMk9hFXUKVZlRkGgK8vARa0Z3YCIXWcDGd6VRLa39c9Uiwkctitf4z5tTUpvnNMuylLoBIQSppNZTTWE4azTqFGI2w+PIDGQrTYEjegp+p9SixcbVfvcUXNpYg9aauHf8BIBvEULeQwj5Ff6vlxfWb3iK2PahNAhBU1420PijdWIUBtNJjGQNnF4o4/CpS3Ap8PobtgKA3299I9IwCs2eAs882ucbBfY2E2PSPDRzfrESWfx2qczmNfzki5h8dd/TU9LjOmW60FhwZWEeQ2PZR3WXBmY2x4HvwittNIWcZGALz2LZPpzGhaWKNDulVGO7Ur7YR2YfRQioqaSGWt2VnrsxS6G59xEQ9GpY+Kg3M7dFw9ULTWEok0St7qJiOeyed+Ap+IN2hHtRrDVXjgNBozAqGIXhjAFCWKfUbx+d8Sfwtet62y/i3vEzYHqCAVbdzP9tGuaLNaSTGrKmjkyyuYITEFrgdjhsffdoBqfnS3jk5DySGsHd124BsHF0hb968KQ/g5YjGoXwoJ0ZT8Td5nlE/hhIYYcq7laPTstnFXMd5s59o7h51xDue7I7ISQuMgMIpI/ywS68ojl8zXGottAUZJ6CGP5ZrtoYSCWxbSgN26HSTUOpxlJGTV1jTfEiwkeFqJRUPmin3vy6+KIf9hR4X6AlYQ5Br8JHQDCjqdspqUCwqjkqnTSKsIfHW5NwYyESFT7SEgQjGQNzJQuf+O5pcFlmJV1ze0Gsu0Ep/T0AIIRkKKWbrzwXwcE5UW2N/b5HufhCM8B0hUdPLqBYc3DjjiFft9gI4SNKKf7kq8/jZQfH/XGDABOahzIsfMTbD/MPDBfQ+P3kRkFciMT7+/zFAm7dPdz03DyVd9doBq+/YRt+/74jeGG2IWCvlGkhfCSGeeyA0OyFCew6BhGvWZp4PtFTqEnaXOQiwkf5lI4dQ2kAbIbzpFekxmHhDnY/mUGWe1mFah1JjTRlB4mDdsIeATfU4ZRU/ndcECaWMaPQ/QU7/Py90hQA9h4u1urYNpRq84gGYaNQtV24FPLwkbdOGHqiyfCM5Ux8/8winru4jFddswX3/+DiujEKcbOPfogQcgTAc973NxJC3t/TK+szooCckzQA48cAK/MUppYqeOb8Eu7YO+p/yFbiKZRqdbzp/3sIT5/rj8SzXKmjZDk4L8S4rbqLiu34ngIQzIKZK9YwnEk2TfyK9hTkusJZz1PYOZzB667fCkLQFW9hRvAUapIdPW+IF77OOEg1hYiGeEAwZFD0jMI2zyiExWZKaSB+nZOE7sRzhatsAbnXxuE/C4umvlEohYyC3rvwEadXdQoAa3URdz4zpzF9zRuJGpHlBTCPREsQjGWNpr/DaM7AsxeWkSAE73jpFYFzrTVx7/j7ALwKrL0FKKVPAvjhXl3UWtCJp9CJpgAwT4FSJq7esXcUSS2BkayxIk/h2EwRT55dxPdOdrdJba3u4Le/8EwgtAI02i1cEBYonrnDhWYgmGM/V6wFhDXTWzzE3TO/v4QgUmw+s1DGeN5E2tCwZTCFF+0ZwX1PrV5XEF9jrS4XhFcSPrIdF3UvVi8zCuHW2UDIU6ixXkV85xo2CrU6O79vFEw9MvtI1vcIEAy0JAOJX0s4fJQxNJh6ImQU3FitXlYC90CSGul6yisADGfZfeGaVc6M7wlmQ6K7bMAOJ5EgGMsZgb5HHB5teNW1k9jreb4bylMAAEppOKG87aeFEHKKEPI0IeQJQshhye8JIeTPCCHHCSFPEUJuiXs93WZeqDrMGJq/ExBZKFkw9UTTh6Ydu70MJD1BcMvuIQAsxrgST4EvFDNtHvufP/803v2p78c+79PnlvDJR07j688GZ0qLz8cX0KUKWxwGBE9BrFWYLwbTdvniIcbZuSd2cDKP56cL0n5QZxbKfqgNAF573RYcmyn6HsRKmV6u+d6eeE12qCEe0JmnIBoCeUpq4+OW1BIw9ETTzOV8Skc+xVJ9z4eMQrhvf87UIwftRGXVtPQUbHn4iBCCkawRMAqVHoaPuKbQC5EZaGgKC6Uau0+dpKTy7CPv/Strmy2ydyznp6WL8Pffz9yxJ3Ki21oR926cJYTcCYASQpIAfgnAszEfexeldC7id68BcMD7dzuAD3j/95Vw/UHO1HFhqTklcL5oYVTiCraDvylu3Dnkx3HH8+aKqpr5Ih3e0Yd56PhcRwsaF3XPhBZccbc6vVTDrtGM7ykMZQw/vh32FK7bPuh/72sKAU+BfX3zrmF86tEzmCtaAWEOAM4uVHDbFSP+91dvHQAAnJgr+ZXnK2G6UMWu0QzmS5Z08U4GjEL83Zu4+xYXXctxoCdIU4FUPtRTSJx/sG0o3eQp8HvGF6B8So9MXY3KqmnlKZQjwkcAmoxCP8JHvTIKfNwqD4nGbYYHNDKzuIFuN7ntg2+9FZpEF3n9DVuhJQju2DsCQgiSGtlw4aN7wdpmbwdwHsBN3ver5U0APkEZjwAYIoRs7cJ5OyJcf8BK2WXho1qsiWthhjJJHJjI4TXXbfF/NpYzV+gpsDdyK6NQtR2cWShjplCLHf446xWJhXfhopbAd67y8FHQUxDDR/LsI3b8LbuY5xTWFay6i6mlSmDx3zPmVYevcj7F7HLN90BkQrOhJ5COMeQ+TNVqeB2VkFYhq8wVw5SU0kDIZ8dw2k9V5DQWIM37Xz69DWB/j/AsBaC1p8CvRRYWkhqFHoWP+PP3yiikkmyy3tlL7L2+Ek+BZyfKpq6JDGaS0t8d2jOC3379Nf4GUxayppTiYw+d7PtAnri9j+YopT9NKZ2klE5QSt8ac/IaBRvO8zgh5B7J77eDFcZxznk/C0AIuYcQcpgQcnh2djbOJXdEWCvImvKU1PmSFWuOQhhCCL72Ky/DO1661//ZeG5lPdX9cM5y9GNPzJbA09DDO/8o+AeEGwfx+XgGCH9uXg0qE5qrtoNCrR5IweNhhkqEpwAAz4V0hfOLFVCKQPhoIm8ilUw0dZ3thFKtjkKtjt3cKAh1CJaXkqonViY0i9lVYQ8kyijw91mt7sJ2aGtPwQqGKvKpZMsuqbLFiO90ZZ5CJaJ4DZAZBbdnngI/r9GDzCPOcCaJc95no5OK5lSSVVzHDR/FRdZo79ylCn73i0faVv13m1ivhBDyZ5IfLwE4TCn9xxYPfQml9DwhZALA1wghz1FKv9PpRVJKPwTgQwBw6NChrvebDdcfZCVzWAG2A96/ynRIzljeRNlyIqsho+DCbytP4ZjQR+X0fAkHt7QvKWkVPrpm2yCePLuIC0tBT2EonfTzyHn4iGdoiZ5CShI+4h+qXSMZjGYNHA0ZBX4dolEghGBPqOtsp3AtZpcX0gtnHxlawptl7BmFDpriibvvcEWzLIsmZ2p+Silv68BDGduH0liu1gPeAxeVg9lHNiilTSHNqFg59xRkxm65asPUE9Id+kjWwKWQp5A2epuS2ou+R5yhjIGznifWSfYRISSwgLcLH8VFNp6Vf85W0oNrNcS96ymwkNEx798NAHYAeDsh5H1RD6KUnvf+nwHweQC3hQ45D9aWm7PD+1lfCdcfZE023NwJVX0udDhLoRXjof4nceHho5LlRO4Sj88U/YKYuJ4C3zUtVexAS9+pxQr2jWUxljP8UBL//UA6iayhISHMVJjzFl2xlqMhNAueguXA8MTWKz2xWURmFPj37TyFrx+ZjhSjuTHdOpiCoSVCxWuNHX0j+6gR3nnHxw/jX56bRhRRQnMtwlNgCwE7rhDqVdRIS20Y/3D6Y87UYTvN09fCoSgR3nBO5inIdB3OSMZAoVZHre74WVY98xR6rCkALAOJf+478RSAYKjHDx91eI7mczYnt/D3xEp6cK2GuHf9BjDB+M8ppX8O4JUArgLwowDulj2AEJIlhOT5195xz4QO+ycAP+tlId0BYIlS2t12mDEI1x+EMwwAtgus2M6KNAUZY8Jg77hUbQdzxRr2jrNdbpS3cGy6iCvGssin9FhGoVZ3cGG5ioOTzKPgC2rdcTFdqGHbUDoQzlgs28ibOrQESxnMmbrvKXAjFwwfyVNSedbFwS15HJ0uBFovnF0ow9ATmAgtUnvGsjizUI4cUOK4FP/xbx7HXz9yWvp7fs8mB0yYeqJJU+ChMr94zftAzhRq+Pqz03jE63IrIyA0xwwf+YKlbxQa4SMgKPSH0x/zEa0ueCiqZfhIsvsMpxKL8Pf9Ytn271kv21wAvTUKvFYBQEcpqUBjvgTQ0BbC7cY7RTaJj3uP/Z5aGPeuDwMQ4yZZACOUUgdA1Ko2CeBBQsiTAB4F8CVK6f2EkHsJIfd6x/wzWF+l4wD+fwC/0OkL6AZhTSHDU8QEy817n3dauBbFSjyFi15G1M07WRw+0ijMFHBgIue112hvFKYWq6AUuHP/KADgnKcvzBRqcFyKbUNpbB1M+QvUcsXGYKbxQWJN8dgber4oCR95C2JFEGJLtUZF7cEteZQtJ5CCeWa+jJ3D6aaMnd2jGVh1FxcjXvtCyYLt0MjW51yLmRhIwUwGu7eKi7eWIDD0hB8SOjFb8l97FOK5ms4rDR81FoKwp7BjmBkF8Z6E49e+6BmqVYiaugaIQnNzT6fZQq2lpwCwvy9/bb0Smv06hV6Gj9KN92+nu/ywMc8aWuzW21HIwkf8PdHvbspx7/r/BPAEIeSjhJCPAfg+gD/yPICvyx5AKT1BKb3R+3ctpfQPvJ9/kFL6Qe9rSil9F6V0n9eOu6mWoR/MF2uB+gNZtWnDcHQuNMvg4/k68RS4nnCTl7EjE5utuotT82UcmMhj10gmlqfAPYMX72NtLPhjuBHYNpTyPQVKqd/3iCPOVJj1PYXGfdK1RFOopmw1PIUrPQ9FLGIL1yhw/FbkEboCv5+yRQ9ghjSd1JA3daSSiabiNXF3Kk5fOznHni+qghhohIxMPRH0FJxWQjM3Cp6m4C1Q4zkTSY0EjIKfktrifSpeo8wo+MZOGj5q4SkIVc2+p9CjRdvXFHoqNIueQodGwQiGj1YrMgPy7KNlX1NYh0aBUvoRAHcC+AKYNvASSumHKaUlSumv9fIC+wEbsWk20sNCucj8GKDzauYoRrMmEgSY7aBWgceXb97JjILMUzg1X4LjUhyYzGHXSBbnLpWbtJEw3Ahcu30Ag+mkn4HEF6TtQ2lsH0qjZDlYrtaxGDIKA0L77PmihayhNaU1ppKJwI6nJPTeuXKSOaFcV+Ats2VGYbffilxu7LhRipoJPV2oYXLA9Fs0y4RmTkbo1nnC6/y6HNHRVXzO4YzRnH0U4SmULQeuS5sW8kSCYOtgMAOpZNWRSib89iGydGBA7JAqD4uwQTvBxzguxULJwnhEeFTsf8RfW+8qmvshNAueQqdGQdCCihFT1zpFll68Lj0FQshV3v+3ANgKlj56FsCWtaw+7jZhAVnWKrgRFumOUdASrEq0E0+Bt5rYP5FD1tACjd04x7yOo/u98JHtUD9rKIqzl8owtAQm8ynsHEkLnoInynqaAvtZBUsVO/ChEttnzxVr0oaBqaQWzD4SNIV8KontQ2l8/8wlACxuXajVpQVqWwfTSGok0ijwvkZRcdjp5SomvCZzqWSzpiAuRGwkJ3td3FOQTZnj8OccyiRjpaTmhPeZn30kiMPbh9IBwTycZtpIBw4aKn+WQsRiJxu0M1+qwaWIDh9xo1Cs+Y/dqMVrQMNTyBgatA5DP0xTaHgKq808apwzOHOEvyf6nX3U7tX8KoB3AvgTye8ogB/p+hWtAc1GoVlTWPA0hW55CgALsXSiKUwtVTCWM5BKapgcSAXmAnCOzRRACLBvPOfXE5yZL2PHcHQF8NmFMnZ48ftdIxm/ZmBqsYKBFBvowgcDcaMQDh8dneEpqTWp4QzPPC7W6hjONq7pTTdtw/u/9QL+4fFz/sQ2maegJQh2jmQi01K5pxC1u5pZruL6HczTMnWtqXV2MHykS8JH7TWFkawR2L1bjivdTWaF8E/B3903jts7nsV9T13wU07DoYqVhI8A+aCduUKzFiQy5M0AWCjbDU1hQwvN7P27kgW9WVPoTviIz8Hm97VQXYfhI0rpO73/75L82xQGAWi0r+DImpXNlywYWnML3NUwnu+sqvn8YtXfsU8MmIFun5xjM0XsGskgldT8RfV0G13h7EIFO7xjdw5ncG6BDXi5sFTxn2/7UEP4XCrbGEjLhea5giX1FNJhT8FyAvfyV/7NlfihvaP4zc8/jX9+hiWg8al1YfaMZiPTUrnOIvMUKKWYXq5h0tsNM08hGD4SWzWnPU3Bdlzfe1qutPcUBtPJ5uwjyQLX2HzU/b7+4q5133gOSxXb17P4gB1O1KCdcCZTmJRkTjM3plGegpYgGEonsVCqCeGj3vY+6kWHVA7PPlpJ6EdMJS52SVOQGXg/fLSeso8IIb8ufP3jod/9t15dVL+ZL9UCHoBUaC4yb6KbXRvDVc2UUmljOM6FxQq2DbLFeXIgJQ0fHZ8u+mMwtw2xUEs7sfnspTJ2jbDz7hzJwHJczBRqOL9Y9Y3BmCd8npgtwXJcDKUb94uHjyilkWIlG8kZnLwmVs7qWgJ/8VM3YzRr4C+/fYJdS4R3s3s0gzPzJem9aqUpFGp1VGzHn1GQ0rUmoVkM82S8kZxnF9h40bGc2cZTYE3i2Ezq9uGjRkqp49UVBBcXPrnuBS/zKRw+4l83aQptCqrG82ZTQ0VeXxLlKQDeGMuS7S9SZq/CR0ajS2qvGPY8BVktRzt4+M1xWZZbp3UOMmQ6pp+Sup48BQA/KXz9ntDvXt3la1kTzi6UUbXdwO5W1gzt4nIVkwPdyTzijHmeAl/cPvjtE3jFn3xbeiylFFOLFWz12iozo1ANLIx1x8WJuSL2T7BsHi1BsGM403K28XLVxmLZ9hdgHsc/s1DG1GLDU0gkCLYMpnDkwjIAhMJHbNBOsVbHQlkuVqaSCVRFobnmNO2wRnMm/vJnDsHQExjNGpE7sD2jWZQsR9pQcLaFp8A9q4kB7ikEU1JtafZR3Q8d3bhjELW6GzAkInwaWVg/iSpeExcC2QSwvV6vJz7elN2zxkLMqo+bG6lxwxW1C94+nG7qwOrXl0R4CgBLx54XPIVeT17rh6bQSTM8jqgFku2OKQAAIABJREFUhf8mKyXbwlNYV+EjACTia9n3Gw7Xpfh/P/cUMobmz00GGh9Wsf/R+cUKtnu5491iPGeiVnf9cZZ/++hpnJgrSXejfNgN37lP5NljxXDG6YUybIf6ngLAFvnTC9FtIc6GKof5/89dXMZSxfaNAgBsG0zjWalRYPfrzEIZlMon06WTjfi941JUbEfaY+f6HYP4y7feit949VWR17zLz0Bqfl2tNAXuWXFPIVy8Fg7zpJM6Kpbj1yjc6GV9RaWl8mH26VB4xnLcpgloQGMhKFTr0jm/24fSMPUEXpjhRiEYqiCEsP5Hoesp1Oow9ETkTn77UBoLJSuw6ZkteONoW2QUjXieAs/Y6lX2UboP2UcD6SQIWZmm4Bc21hymKXQxfCTqmDwldV2Fj8DEZNnXsu83HJ985DQeOj6P33rdNYFMl4TXEE3sYDm1WPEX5G4h1io8eW7JTwW9KGnbzWsU+CLNF7YZQWzmmUcHJhtGYfcIK2DjHgWlNLAY+NPNvNe/bSgFQoDveZW74qjC7UNpf0EMZx8BDTFWFoJIG42Fkj9/lEB311UT+IkX7ZT+DhBrFZo9oFbZR41qZs8oSDyFcPiobDs4MVfCcCaJnV6ILcoo8GH2qWQC1brr3/NWKakA9xSa21IkEgR7x3O+p1CQZLqI1eScYrXecgfsF8YJXVjnijWM5VuHR0eyBuZLlpB9tDHnKQDMix5MJ1ekKXDP4FLZguW4K/I2os4ZHs8KrL/eRzcSQpYJIQUAN3hf8++v78P1dY0jU8v4sQ88jH984jysuosTs0X89y8/i5cfHMdbbmtegMT22azvvtt1ozCeY4vTXNHCF59sTBSTzXLg+eo8C4gvbKKucNxrhCfOMN49mkGhWvczkT7y4Em86L9+3U9T5YaIGwVT17B1IOVPdhNfs+g1hOsUAOCUZxRGpeGjRvaR37d/hR+m7UNpaAnS1EK7VGPelJ4gUk2B36sJQWiuBVJSw9lHTGg+OVfE3vGc/zqjqporluO3ZXZcCtsRjIIsJTXVCEOIsxRE9o5nccK7r7L0x6j89lZxbv43PSeEkGaLNb/KPoqRrIFL5YZR2MieAgD83huvxc/duafjx/G/AddleiE0U0p9TaFqu5FtXXpBu+wjjVI6QCnNU0p172v+fecKzRoyX6phtljDL336Cdz5h9/A//Wxx2DqGv7Hv7tBujvKGo322XxHta1HnsL0chVfeuoCrvK6mUo9BaGQDICvb4gFbMdmitg+lA68SXcJGkGhauMvvnkcJcvxxdyzl8oYSOmBRX7HSMaP128VXvNWwWuQhY9OtPAUxDh7o4fPyhYVQ09g+1C6yVPgmVw7htOwHerPR+DMFKrIm7p/f1JJLdBMjmUfBesUrLqL4zO8l5S8WIxTqzOh2e8K64XLoiqaxYVguVqXip77xnM4u1BGxXJQtpp1GNmc5nbTxHh6csBTKFgtRWaAxeEdl/oZXj2rU/DnKfQ2Qv2mm7YHhkHFhYeP+GevWxXNQOOzwftXcY+8nyGk3pridcRLD4zjm7/6cnzsbS/CTTuHcO5SBf/tR6/3d9xhxFxkv7K3B5oCANz/g4u4uFz15y1MSYrNppaqSGrE/+BO5D1PQQgfPX+x4Of4c/z4+0IZH3/4FBbLNm7ZNYRPPXoGM4Uqaycx2tyJFAASBH76JhDyFEK9j4CGpyDbcaaFUA33FDKryO9mfZ2CngLXE3ZKBugALF11fEDsyaTBcly/4luWfQQwT443GASiaxW4puAbBYtlqDguhaE1L6CmnoCWICh6LbJlQ3H2jWfhUuDIhSUAzTHw8PQ2AG2H0U/kWSbZuabwUWujwD3AqcUKDC2x6n4/UZh6AjtH0tIxlusB31PwjEJ3iteCngL3Eia9z3k/xebLxigALEb78oMT+PB/eBGe/f1X43WCuBxG7EvDd+k7hlY+AlLGcMaAliC4/5mLSCUTeM11WzCWMyM9hS2DKf+DmDY0DKR0f9d2camK5y4WAuMrgcYC/4PzS/jQd07glVdP4E9+4ibYjouPPHASZxfKTamf/DFbBlJ+SwWg4aUkCJATFnRRU0hqBANpeSO2is0qNn1PYRXhh92jzS20+b3grTDCu6vlqu2HgADA9Bqv8Wwi23ED/XbSwmvcO5b1azOiWl1U643sI4C5/f58ZomnwHrza7hUtlGru9KQDw8FPnmWGQWZpyDrrtkq1ZK30OCbnbrjellj7T0FgG1aejWfGWD35YFf/xH8+KFoXWkt4R4uD0d2xSgYwYJZnkDCM+Wi2rb0gsvKKIi0E7FygqZw7lIFOVOXLnarIZEgGM0yl/wVV00ia+rYNpSK1BR4jQKHp6UCwNeOsOlMr7p2MnBMxtAxnjfxsYdPYblaxy+/8kpcMZbFG27chk8+chpnFypN7SS4oBoOl3E9YyCdDOwS+WJ2qWxjNGtKw3GpZMKPs/uewio+THtGs1iq2FgsN9JSZz2vyR+1GWqKx8IvDUPEhVLuwYTDR2I75L3juaYpc2HE7COAGaVWRgFgXhb/G8oWct4m/enz3CgEDWk+pUu7pLYTP7cPpXHe64a7ULJAaet0VAD+1MGpxUrP0lE3AtwwdzN8pGsJpJIJv30G90Z5JEN5CusAln3kaQpe5lE3C9c4PBzEU2K3DKSkvYqmhGpmjmgUvvKDaewdz/o1CiK7RzKo1V3cfc2kH0N91137UbYcWI7bZBR2+ZlIwefLp9j4TbHtMMCyiPitkYnMQHCmQrErnkJzBtJMoQYtQfzrDnsKpVo9ELIKT4STZR8BACHM+8h5rzOq/1HVZtlHvPiqajuoOezcUUYha2q+JyrzFDKGjm2DKTx5bhFA8640ZyY71hSAYK0CF0zbegpZ9nefK1rKKIA1VwS64ynw8zTCR+x/rh32MwNJGYUIxD/Q+UuVQGpmNxnPm8gaGu66agIA242HPQXHpbi4XG26hom8ienlGpbKNh45MY9XXbtF+hxcM/jlV17p/+zKyTxecx07fmdIK+HhJJmwvn0oHRCZAebx8A9GlFjJxcOa7TRSUlfxYeK1GOIYz9kC67uUjZhDXLGdgCHiC1vNSx8NZx/xa942mEYqqfmvMyr7iFc0cwFW9BTMCM80a+r+bIioxWXfRM6vlQjfs3xKh+U0CurY1LX2VbY7htOYXq6hVnf8wrXxfOu+XuJ88vRlbBS4B+lrCl2oaAaCOib3FLZ4nkI/q5q7Gw/ZRIh/oKmlCm7ZPdST5/nFH9mP+WLNX6C2erUAYkuDmULVH3YjMjGQwkyhim88N426S3H3NZNN5weAd750L+7cN4Zrtg0Efv6rd1+JYq2Om3YGX9t43sTPv2wv3nBjs+by1jt2Q+YwDXj9j6I8BTGk0pgLsPK3366RDPKmjmemlvAT3kTXmUINE/lUwyuxwp6CEwhZmXpjR285zWEe7lXwEI74OmVUbE9oFsaPtgsf5UzdTxeO0gH2jefwwLE5//jw4wEmLps5JuY7Lm07TYzrQxcWq36mWbvso7Sh+f2ieqkprHd4HRPPdst1oSEeEJzT0NAU+h8+UkYhgqyXo17w2kBs77LIzHnRnqAwzOP2F5eqfiaRP+ymSVMwYTsUn37sLCYHTNy4Q264rt46gKu3DjT9fP9EHp98++1NPyeE4D2vuVp6rrfesVv6c3E4jAwxfMQ9hdXkuScSBNdsG8AzXqwdYJ7ClsGUdCY04A32kXgKVdvxawoMrTl8xNtNAHygULOnQClttLnQBaMgMTYislbYYUSj1CQ0C1XRo0JvpjjhI4CFRmdj9D3ijGZNnL/MNQWA1zHxOoXu3IvgJL6gp1BWQvPawz98x70WA91OR42CvwlEXYFfwxVjwRQ9LkI9enIBd1+zpWcpgnHgWT2R4SNh916yHBhaYtXFSddtH8SRC8t+SinzFMyAV8JxXSZwi5oCzz4Ss4TE3HhuFMT7Lg4UErEcFy5lho4bpYDQ3CJ8JJ5bhliMGF6A+CaCp5cW+CjONqE5nkl3/lIFc8UaMoYWK5zHdYXL3ShwY2zqiUCG3mrImg0dc7lqQ0sQ3/OWzdTuFcooRMA/IEe9aWDbe6QphNnqeQOirvDshQIyhtY0X0Bs0Hf3tfLQUb/gu9x2QnPVdpng24Xd1bXbBlC1WXU6mxwWMgqCy80NRCD7yNcUHL/QTZwLfMVYFv/5tVfjzTdvD7xOWfiIZzCZesIPrbRLSQXieQqiUQiHj7g3yavZ27XN5mwZZO1MznmeQlTL7DB8HO3lHD4CGhuGbnRI5QQ1BaYL8U1MPzWFy/sv2wL+4Xv+oucp9Ch8FGZykH3oxFqFIxeWcXBLvskT4AVsAykdd+wd7cv1RcE/HO2E5qqnKXRjMAnPpHpmagnzxcbksJSX/SN6CjzVL5B9pDcMlWxHTwjBO394r997H2DpuLLwkTiiUjRK7YyCaKSiQj6TAywZIUGaBd7xvIl8Ssdx3h9JMqxHhqEnsGUghXOXyi1nM4fhc0cuZ6EZaGwau5GOygmGj+oYSCX9+7ypwkeEEI0Q8n1CyH2S3+0ihHzT+/1ThJDX9vp64sJ3AsdmCkhqxO+X02tMXcNYzvDDR5RSPHdhWaoJTAyYIAR4xdWTPW0eFgcukkZ7Co2FumzVuxKH3TuWRSqZwDPnlxtplfmU/0EKjv9s9hTE4rV2sX9OtKfQGFEptrmotdUU2H1LJ7XIvyEhrDFe1tSb0qIJITgwkfObIRZr8TQFgNcqsPBRu3RUDi9gU+Ejzyh0SWQGgp7CcoXN10glEyBk83kKvwTg2Yjf/RaAz1BKbwab3fD+PlxPLHJC+GjrYLqv8fotQlrq1FIVy9W61CiYuoY//cmb8Sv/5sqm3/WbdkKzuFCXQrH9laJrCVy9lYnN4uSwcP0BEOEpBIRmrinEMwrhAT9ikzg/q8lqrynkzHhhiKu35iPv7f4JoZOqZ7Ci9AkRXqswW6j5fbjawY3+5W4U+KaxW+mogGcULAeu20grJoRIx6f2kp4aBULIDgCvA/DhiEMoAL7aDQKYijiu7zSqFmtd747ajq2DaT989OwUm19wzdbmojQAeOON26QD7vvNtdsGsXc8GznDWhR/y7XueAoAcN22QRyZWm4M0MmbSGps+IwYPvI7swbCR43Yv11vzj6SMeANFAp/SBtzixP+B7lab4SlZPMUgMb7rJ1R+M3XXo2Pve026e/2T+QwV7SwWLZih48A5ilcWKriUtmOHT5SngKD399ujuflG4Sy7QRasvBuvf2i157C+wD8OgA34ve/C+CthJBzAP4ZwLtlBxFC7iGEHCaEHJ6dne3JhYYRY4Xd7o7ajq2DKT8NlQ+1Obil2VNYT7zuhq34l199eWQmBs/d59lH3fAUAOC67QMo1Oo4fOoSgMaMYTaHuPG24255JkJotrzK42Tb8JG8Uyp37/k5U8lETE3BW1za7OyHMkbkzOqG2FxsjOKMEz4aTvuZW/GFZm4ULm85sheagt8Ur8pbqXuhRUPbHNlHhJDXA5ihlD7e4rC3APgYpXQHgNcC+CQhpOmaKKUfopQeopQeGh8f79EVBxF3sv1KR+VsGUxhuVpHqVbHsxeXWYuFLr751gIxd58NoO/OTvPabUxs/tbRWQykdH9R5g34ODJPwdTFLKF4nkJUp1TeJjstPH+cOgWeOirrkBqXA15rk2OeUUglE7E0JtEDjuspjCihGUCjRUuuSx4vO5fQSr1i+73WMkl907TOfjGANxJCTgH4NIAfIYT8deiYtwP4DABQSr8LIAVgrIfXFBtxB7Cjz54CL1K7uFzFsxcK/pyFjUxSI9ASBFXb9YTm7hi5KyfzSGqkKa0ybQTnJPuegmCMdC3hD+Sx/cW7tXYU1Sm12uQpaB3VKawmtXH7UBqpZALHZ4ooVO221cycHUJ33M49hcvcKPQgfMQ3LIWqjaIV9BQ2RfiIUvoeSukOSukeMBH5Xyilbw0ddgbAKwCAEHI1mFHoT3yoDeKOci08BQA4MVvCqfmSVGTeaPA4O29z0S2jYOgJXDnJjCZP0QWYZyJmbDRmOAQXs5Q356FRvBbPUwg3xfOFZsEoxKlT8I1CzIVcRiJBsHcs5xmF9n2POKKnEDf7aMdwGi/eP4pbdg2v6Fo3C70MH00v10Bpw3sMz/zuNX0PDBJC3ksIeaP37a8CeCch5EkAnwLwczSc1rFGaAnix037LzSzxe1bz8+AUmwKowCwhbJUq6NiO02L82q4zgshibvdlKFJ6xTCH2I2T1n0FNoLzUDzSM6G0NzQFDppc7HaIqj9EzlfU4h7rrSh+XUHccNHqaSGv3nHHbh+R+cTyzYTPLzcXaGZneuil46+VkJzXwLVlNJvAfiW9/XvCD8/AhZmWpfkTB1V2/J37v2Ct6/45nMzAIBrNo1RSGChxJqvdTO/+7rtA/i7wwjUkqSTiaCmUHOQIM1ZQKauoWa7/uLdzlMY8DWF1p5C2jNKtXYpqSluFFY33fbARA7/9OQUBtJJDGfin2v7cBq1utuzecubFf7+7Wr4yDM0PB2dG/d0aIPTay7vFII2ZE09kPfeL1JJtoObWmIzhXf0OXzVK9JJDfOeUehGmwvOtdubPQUu9HJKVt2b+xDUDEzPU2gX++dEZR/x5+IFcSld87ukGloichZHztTxX95wDX5UaKWxEngG0tHpQkcL1f7x3KZ5f/WTXA/CR/yc3Chw/Yp5Cv3LPtrYKS09Jmvofl52v9kymMJ8ycJVW/M9Ge6zFqQNDfPF7g4mAVj46M03bcPLD04EnqtyKegpyAxRStdQE7uktgkfscwe0iw02w6I4Inw8JVVd9ue820vvqL1C4wBNwqOSzvyOn779df0tYXCZoH3juqm3pj1jUJw6FLG0Ddf+Gij8rYX71kzt3rrYAo/mFrGVeu8PqETUrqGk0U2LKZbdQoAW8jf95M3B58r2awpyEJWfD5A3IpmQgjykk6pVdtBStd8A57yw1LOqrvBxmH3aBZ6gqDu0o70ieGsgctbMl4Ze8dzeOQ9r/BDvd0gY2ggRPAUhOyjfs5oVkahBWs5OJx3S90sIjPAds+FLozijEM4fFS2IjyFpMaK1ySts6OQ9T+q2E5gA5E2Eg1PoQ99qQw9gd2jGbwwW9rwNS0bhW4aBIBtOLKGLszs5nUKGmyHwnbcvvQ4U5rCOoWL21dHtLfYiKSFKthMjxeucBpfeD4zx9SZpxC3IR7AdnCy7KOU8Fi/eC1G+Khb8BBSN9s5K/pL1tT8UKZYpwD0b/qaMgrrlLsOTuA1121pGqG5kRGrYHvuKXgxfZ7hXLYc6XOmhMUbAJKJ9h+JKE8hJZzfL15z+m8UutmkTdFfuK6QSjaGUKWFFjH9QBmFdco12wbwgbfeClPfPKmCYhZXrz2FVFKDS+F7ACWrLn3OVFLz6xT0BInVDXcg1TxToebNZxbPSynLUupH+AhotLtQ4aONS6NupZEskPE9hf5kICmjoOgbolHo1rDzKPxW3V5TvHItylNoCM1xd/SRnkLIKACs8rlfnsJNO4dgaAnsHcu1P1ixLuHJEGIvrHSSfd2v8JHaUij6RlCI7X34CGCL9SCSzFOQagosJdWqxxfxWPZRc5dU8fzcKC1X7Nh9hVbLnv/T3r3H2HGWdxz//vbsrnfX8S22mzh2nITEKVg0JI6FTNOmTkBqUiK7UhElTUSEWqJKqLmUNiL9g9JQpCIhIBUoUnBSaIHQNqDEAomkDXHT0jrYxkAcm7bgXEhwYpMaOxfb6/U+/WNmzs6ePed4z+7OrD3n95FWPpfZc2Y8s/PM+z7vvM+SuTx912/PesElm7q5bVoKZd3A5qPHSpOdKPtrPYVfPefrN0Qk9Q+a1XBIbl4bZfhETPpkOn8wKZuYTTsNaaI5l0jPHh8+crxlLYUiOCCc3rJZV7Mb1yAXFJxTsKrJTpQzeTdz6+8a+0M6NjLKidFo2lIY6K0xPDLKsZETkz55Z1dxr+VaC0cbuo+yoHToyPHScgp2+ms2a65HH1llZSfKmZz3qOV35ZrcY7UUmo8+guQEP5l7FCA/U+pYsrkxKGSPR0ajtJyCnf6yRHO+nOpYq9eJZquY7EQ5U6U42xlXE7peda35Hc2QnOAn3X00MLGmwpEmo48yDgo2WXPnTEw0Zy1ctxSscrKr95mc4qLld+W6j5pVXcvUWwrHJj9KqNlMqa1yCnDySfbMMvWgkMsp+D4Fq6ysJGcpLYX+tNTmyIl6LYVmuYwsj/Dq0ZEOEs3jaypExISWQn50lVsKNllZonneuJaCcwpWUWW2FPKJ5jeOnbylcLiDhPC8hpZCVjNhTpNEMzgo2OQ1SzT31ZKZeR0UrHLqOYUSZp7N5xSyO0GbVXvLunle7eAms7GaCsfr35H/zuRzHRSsc3ObJJph4gSPRfLRaqUZG5I6S6OPmiWae8dGCXU++igJNlkpzsH+5kFhjnMKNklvOXs+bz573oTZkZOaCuWMPvIdzVaa7Eq6jLl5spP9keHRsZxCk5bCnHxCeJJX9H21Hgb7avWWQnanactEs1sKNklnLxjg27ddOeH1wRLrNPtotdKM5RSK7z7q6RH9vUlNgyyn0DQo5CYc7ORu4HkDvRw+kgSbbFTIuDxCrYdsbj0HBZuuxqngi1T40SqpJmmnpG+2eP+9knZLelrSV4teH5s9WYJ5ukXqJyvrhx1rKbRONENnQ0fPWTjIM68kVeSOjmT1mcc+S1L9sz0k1aZrqMSWQhndR7cCe4AJhQEkrQLuBK6IiIOSfqVxGauOBYN93HPDGta9aXEp35ddXfX39jDQ10OtybTYU+3muWzlQh743vMcPzHK0SYthez5G8Mn6K/Q9Oc2Owb7a/UcVtEKvYSRtAJ4N7CpxSIfBD4fEQcBImJ/ketjs+/aX1vGorn9pXxXVmjn9WPN6zPD+JZCJ91Hl5+3iKPHR9mz73C9pTDQEBTqLQV3H9k0DfXXOFKRegqfBe4ARlu8fzFwsaTvStoq6ZpmC0m6WdJ2SdsPHDhQ1LpaxWTVz1rVZwbGTYLXaVAA2PHcQY6kNRsaWwpZK8RBwaZrqL/39J86W9J1wP6I2NFmsV5gFbAeuB74gqSFjQtFxL0RsTYi1i5durSQ9bXqGezrqc99NJmWQicn72ULBjlnwQA7njtYHz+e74rKf7ZzCjZdg/3VSDRfAWyQ9CzwNeBqSV9uWOYFYHNEHI+IZ4D/IQkSZtOW/SG9MXyi5YinvtpYrqF/kvcpZNact4jvP3ewfgXXLKcAlFpPwaopy0+VobCjNSLujIgVEXE+8D7gOxFxY8NiD5G0EpC0hKQ7aW9R62TdZTDtPnp9eKTpjWuZgfSk3WmBmjUrF/HzQ0d59hfJKKQ5jUGh3zkFmxlDaX4sIk6+8DSVfrRKukvShvTpI8ArknYDjwN/HhGvlL1OVk31nMKx1i0FGDuZd3ryzvIK3/1pcsg2thSyeyAcFGy6BvtrRIzdPV+kUu5ojogtwJb08Udzrwfwp+mP2Ywa7KtxdPgEwz2jbQv7TLWlsPqc+Qz09bBn32FqPZowTUa9peCcgk3TUF82U+pI4fXNfbRaZWVDUtuNPoKxhHBfh1f0fbUeLlmRjIsY6O1BGh8UsmDjloJNV5mFdny0WmXVcwptRh/BWPfRVCauy7qQml29OadgMyU7lsqYKdUT4lllDfTV6n2w7Wo4ZENJ+3o7G30ESbIZxs+hlP9+cPeRTV+ZhXZ8tFpl5a/e21V7mzPFnALAmpULJ3xXZsBDUm2GDDoomE1ffjRQ+5bC1K/oF58xhwuWzJ1w41ryuc4p2Myo1xw/XvxUF+4+ssrKB4V2LYWs9kKniebMh666iGMjE6/gfv3CJVx3yeHSZoW16ioz0eygYJU10D/ZlkJ6RT/Fvv/3XL6i6euXnruQz/3Bmil9plmecwpmM2BcS6HdzWu+ycxOcfXysg4KZlM3LqfQbpqLvqknms3KMJSrOV40/xVYZQ32jx3e7VoK9ZvXOpwQz6wsWd7L3Udm05C/d6Dd1ABzPHTUTnE9PWLlmUOlHKNONFtljbtPYTI3r7n7yE5hT9xxVSnf478Cq6zxOYWTJ5odFMwcFKzCsqDQ26O2w01XL5vPxWedwZkl1Y42O5W5+8gqK+s+GuqvTZjBNO8dFy7m0dt/q6zVMjuluaVglZUl5dpVXTOz8RwUrLIkMdhXa1t1zczGc1CwShvsr7mlYNYBBwWrNLcUzDpTeFCQVJO0U9I32yzze5JC0tqi18e6y0BfT9vJ8MxsvDL+Wm4F9gDzm70paV66zJMlrIt1mVveuYpFQx5qajZZhbYUJK0A3g1sarPYx4FPAkeLXBfrThsvXc6VFy+d7dUwO20U3X30WeAOYLTZm5LWAOdGxLcKXg8zM5uEwoKCpOuA/RGxo8X7PcCngQ9P4rNulrRd0vYDBw7M8JqamVmmyJbCFcAGSc8CXwOulvTl3PvzgLcCW9Jl1gGbmyWbI+LeiFgbEWuXLnVXgJlZUQoLChFxZ0SsiIjzgfcB34mIG3PvH4qIJRFxfrrMVmBDRGwvap3MzKy90u9TkHSXpA1lf6+ZmZ1cKQO4I2ILsCV9/NEWy6wvY13MzKw139FsZmZ1DgpmZlaniJjtdeiIpAPAc1P89SXAL2ZwdU4n3vbu1K3b3q3bDa23/byIOOnwzdMuKEyHpO0R0ZXzK3nbve3dpFu3G6a/7e4+MjOzOgcFMzOr67agcO9sr8As8rZ3p27d9m7dbpjmtndVTsHMzNrrtpaCmZm14aBgZmZ1XRMUJF0j6b8l/UTSR2Z7fYok6VxJj0vaLelpSbemr58p6V8k/W/676LZXtciNJaAlXSBpCfTff+PkipZik3SQkkPSvqxpD2S3tFF+/z29FjfJekBSQNV3e+S7pe0X9Ku3GtN97MSf5v+H/worWHTVlcEBUk14PO3OUA1AAAFAElEQVTAtcBq4HpJq2d3rQo1Anw4IlaTTEn+oXR7PwI8FhGrgMfS51WUlYDNfBL4TERcBBwE/nBW1qp4dwPfjog3A28j+T+o/D6XtBy4BVgbEW8FaiQzM1d1v38RuKbhtVb7+VpgVfpzM3DPyT68K4IC8HbgJxGxNyKGSeo7bJzldSpMROyLiO+nj18lOTksJ9nmL6WLfQn43dlZw+I0loCVJOBq4MF0kapu9wLgSuA+gIgYjohf0gX7PNULDErqBYaAfVR0v0fEE8D/Nbzcaj9vBP4+EluBhZKWtfv8bgkKy4Gf5Z6/kL5WeZLOBy4DngTOioh96VsvAWfN0moVqbEE7GLglxExkj6v6r6/ADgA/F3adbZJ0ly6YJ9HxIvAp4DnSYLBIWAH3bHfM632c8fnvm4JCl1J0hnA14HbIuJw/r1IxiJXajzyyUrAVlwvsAa4JyIuA16noauoivscIO0/30gSGM8B5jKxe6VrTHc/d0tQeBE4N/d8RfpaZUnqIwkIX4mIb6Qvv5w1HdN/98/W+hVkQglYkn72hWm3AlR3378AvBART6bPHyQJElXf5wDvAp6JiAMRcRz4Bsmx0A37PdNqP3d87uuWoLANWJWORugnSUJtnuV1Kkzaj34fsCciPp17azNwU/r4JuDhstetSC1KwN4APA68J12sctsNEBEvAT+T9KvpS+8EdlPxfZ56HlgnaSg99rNtr/x+z2m1nzcD709HIa0DDuW6mZrqmjuaJf0OSX9zDbg/Ij4xy6tUGEm/Afw78BRjfet/QZJX+CdgJcn04++NiMaEVSVIWg/8WURcJ+lNJC2HM4GdwI0RcWw2168Iki4lSbD3A3uBD5Bc+FV+n0v6K+D3SUbe7QT+iKTvvHL7XdIDwHqSKbJfBv4SeIgm+zkNkp8j6U57A/hARGxv+/ndEhTMzOzkuqX7yMzMJsFBwczM6hwUzMyszkHBzMzqHBTMzKzOQcG6iqTFkn6Q/rwk6cXc8/8s6Dsvk3TfFH/3X6s6s6mdmjwk1bqWpI8Br0XEpwr+nn8G/joifjiF370JWFHl+2rs1OKWgllK0mvpv+sl/ZukhyXtlfQ3km6Q9D1JT0m6MF1uqaSvS9qW/lzR5DPnAZdkAUHSxyT9g6T/Sue+/2D6+jJJT6Qtll2SfjP9iM3A9aX8B5iRTKJlZhO9DXgLyRTFe4FNEfF2JQWL/gS4jWRepc9ExH9IWgk8kv5O3lpgV8Nrl5DUuZgL7JT0LZIT/yMR8Ym0/scQQEQclDRH0uKIeKWQLTXLcVAwa25bNkeMpJ8Cj6avPwVclT5+F7A6mUkAgPmSzoiI13Kfs4xkSuu8hyPiCHBE0uMk9T62AfenExk+FBE/yC2/n2T2TwcFK5y7j8yay8+RM5p7PsrYxVQPsC4iLk1/ljcEBIAjwEDDa42JvEgLp1xJMoPlFyW9P/f+QPo5ZoVzUDCbukdJupKA+oR0jfYAFzW8tjGtIbyYZGKzbZLOA16OiC+QTGq3Jv1MAWcDz8742ps14aBgNnW3AGvTgui7gT9uXCAifgwsSBPOmR+RTOu8Ffh4RPycJDj8UNJOktk+706XvRzYmqsgZlYoD0k1K5ik24FXI2JTp8NgJd0NbI6Ix4pcR7OMWwpmxbuH8TmKTuxyQLAyuaVgZmZ1bimYmVmdg4KZmdU5KJiZWZ2DgpmZ1TkomJlZ3f8DuZutgW984gYAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "time_close = [res[0] for res in nma_close.results]\n", "ax = plt.plot(time_close, eigenvalues_close)\n", "plt.xlabel('Time (ps)')\n", "plt.ylabel('Eigenvalue');" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## References\n", "\n", "[1] Oliver Beckstein, Elizabeth J. Denning, Juan R. Perilla, and Thomas B. Woolf.\n", "Zipping and Unzipping of Adenylate Kinase: Atomistic Insights into the Ensemble of OpenClosed Transitions.\n", "Journal of Molecular Biology, 394(1):160–176, November 2009.\n", "00107.\n", "URL: https://linkinghub.elsevier.com/retrieve/pii/S0022283609011164, doi:10.1016/j.jmb.2009.09.009.\n", "\n", "[2] Richard J. Gowers, Max Linke, Jonathan Barnoud, Tyler J. E. Reddy, Manuel N. Melo, Sean L. Seyler, Jan Domański, David L. Dotson, Sébastien Buchoux, Ian M. Kenney, and Oliver Beckstein.\n", "MDAnalysis: A Python Package for the Rapid Analysis of Molecular Dynamics Simulations.\n", "Proceedings of the 15th Python in Science Conference, pages 98–105, 2016.\n", "00152.\n", "URL: https://conference.scipy.org/proceedings/scipy2016/oliver_beckstein.html, doi:10.25080/Majora-629e541a-00e.\n", "\n", "[3] Benjamin A. Hall, Samantha L. Kaye, Andy Pang, Rafael Perera, and Philip C. Biggin.\n", "Characterization of Protein Conformational States by Normal-Mode Frequencies.\n", "Journal of the American Chemical Society, 129(37):11394–11401, September 2007.\n", "00020.\n", "URL: https://doi.org/10.1021/ja071797y, doi:10.1021/ja071797y.\n", "\n", "[4] Naveen Michaud-Agrawal, Elizabeth J. Denning, Thomas B. Woolf, and Oliver Beckstein.\n", "MDAnalysis: A toolkit for the analysis of molecular dynamics simulations.\n", "Journal of Computational Chemistry, 32(10):2319–2327, July 2011.\n", "00778.\n", "URL: http://doi.wiley.com/10.1002/jcc.21787, doi:10.1002/jcc.21787." ] } ], "metadata": { "kernelspec": { "display_name": "Python (mda0170)", "language": "python", "name": "mda0170" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.7" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": false, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false } }, "nbformat": 4, "nbformat_minor": 2 }