{
"cells": [
{
"cell_type": "markdown",
"source": [
"# Tutorial"
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"DFTK is a Julia package for playing with plane-wave\n",
"density-functional theory algorithms. In its basic formulation it\n",
"solves periodic Kohn-Sham equations.\n",
"\n",
"This document provides an overview of the structure of the code\n",
"and how to access basic information about calculations.\n",
"Basic familiarity with the concepts of plane-wave density functional theory\n",
"is assumed throughout. Feel free to take a look at the\n",
"[density-functional theory chapter](https://juliamolsim.github.io/DFTK.jl/dev/#density-functional-theory)\n",
"for some introductory material on the topic.\n",
"\n",
"!!! note \"Convergence parameters in the documentation\"\n",
" We use rough parameters in order to be able\n",
" to automatically generate this documentation very quickly.\n",
" Therefore results are far from converged.\n",
" Tighter thresholds and larger grids should be used for more realistic results.\n",
"\n",
"For our discussion we will use the classic example of\n",
"computing the LDA ground state of the\n",
"[silicon crystal](https://www.materialsproject.org/materials/mp-149).\n",
"Performing such a calculation roughly proceeds in three steps."
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "3Γ3 Matrix{Unitful.Quantity{Float64, π, Unitful.FreeUnits{(β«,), π, nothing}}}:\n 0.0 β« 2.7155 β« 2.7155 β«\n 2.7155 β« 0.0 β« 2.7155 β«\n 2.7155 β« 2.7155 β« 0.0 β«"
},
"metadata": {},
"execution_count": 1
}
],
"cell_type": "code",
"source": [
"using DFTK\n",
"using Plots\n",
"using Unitful\n",
"using UnitfulAtomic\n",
"\n",
"# 1. Define lattice and atomic positions\n",
"a = 5.431u\"angstrom\" # Silicon lattice constant\n",
"lattice = a / 2 * [[0 1 1.]; # Silicon lattice vectors\n",
" [1 0 1.]; # specified column by column\n",
" [1 1 0.]]"
],
"metadata": {},
"execution_count": 1
},
{
"cell_type": "markdown",
"source": [
"By default, all numbers passed as arguments are assumed to be in atomic\n",
"units. Quantities such as temperature, energy cutoffs, lattice vectors, and\n",
"the k-point grid spacing can optionally be annotated with Unitful units,\n",
"which are automatically converted to the atomic units used internally. For\n",
"more details, see the [Unitful package\n",
"documentation](https://painterqubits.github.io/Unitful.jl/stable/) and the\n",
"[UnitfulAtomic.jl package](https://github.com/sostock/UnitfulAtomic.jl)."
],
"metadata": {}
},
{
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"n Energy Eβ-Eβββ Οout-Οin Diag\n",
"--- --------------- --------- -------- ----\n",
" 1 -7.905203413380 NaN 1.95e-01 4.0 \n",
" 2 -7.909822292766 -4.62e-03 2.98e-02 1.0 \n",
" 3 -7.910051958765 -2.30e-04 3.05e-03 3.1 \n",
" 4 -7.910052832343 -8.74e-07 4.73e-04 2.4 \n",
" 5 -7.910052852374 -2.00e-08 3.89e-05 2.0 \n",
" 6 -7.910052854035 -1.66e-09 4.83e-06 3.1 \n"
]
}
],
"cell_type": "code",
"source": [
"# Load HGH pseudopotential for Silicon\n",
"Si = ElementPsp(:Si, psp=load_psp(\"hgh/lda/Si-q4\"))\n",
"\n",
"# Specify type and positions of atoms\n",
"atoms = [Si => [ones(3)/8, -ones(3)/8]]\n",
"\n",
"# 2. Select model and basis\n",
"model = model_LDA(lattice, atoms)\n",
"kgrid = [4, 4, 4] # k-point grid (Regular Monkhorst-Pack grid)\n",
"Ecut = 7 # kinetic energy cutoff\n",
"# Ecut = 190.5u\"eV\" # Could also use eV or other energy-compatible units\n",
"basis = PlaneWaveBasis(model, Ecut; kgrid=kgrid)\n",
"\n",
"# 3. Run the SCF procedure to obtain the ground state\n",
"scfres = self_consistent_field(basis, tol=1e-8);"
],
"metadata": {},
"execution_count": 2
},
{
"cell_type": "markdown",
"source": [
"That's it! Now you can get various quantities from the result of the SCF.\n",
"For instance, the different components of the energy:"
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "Energy breakdown:\n Kinetic 3.0795154 \n AtomicLocal -2.1806377\n AtomicNonlocal 1.7339395 \n Ewald -8.3979253\n PspCorrection -0.2946254\n Hartree 0.5417570 \n Xc -2.3920764\n\n total -7.910052854035\n"
},
"metadata": {},
"execution_count": 3
}
],
"cell_type": "code",
"source": [
"scfres.energies"
],
"metadata": {},
"execution_count": 3
},
{
"cell_type": "markdown",
"source": [
"Eigenvalues:"
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "7Γ10 Matrix{Float64}:\n -0.170182 -0.131801 -0.0883282 β¦ -0.0562607 -0.11494 -0.0700168\n 0.201343 0.0909036 0.0122924 0.0111144 0.0420613 0.0176495\n 0.249295 0.174773 0.176137 0.132964 0.220114 0.11233\n 0.249295 0.231428 0.20237 0.161043 0.220114 0.190456\n 0.350986 0.360027 0.340134 0.291807 0.320727 0.327376\n 0.36997 0.395897 0.389483 β¦ 0.331814 0.38819 0.460233\n 0.36997 0.401676 0.412474 0.565537 0.38819 0.462717"
},
"metadata": {},
"execution_count": 4
}
],
"cell_type": "code",
"source": [
"hcat(scfres.eigenvalues...)"
],
"metadata": {},
"execution_count": 4
},
{
"cell_type": "markdown",
"source": [
"`eigenvalues` is an array (indexed by kpoints) of arrays (indexed by\n",
"eigenvalue number). The \"splatting\" operation `...` calls `hcat`\n",
"with all the inner arrays as arguments, which collects them into a\n",
"matrix.\n",
"\n",
"The resulting matrix is 7 (number of computed eigenvalues) by 8\n",
"(number of kpoints). There are 7 eigenvalues per kpoint because\n",
"there are 4 occupied states in the system (4 valence electrons per\n",
"silicon atom, two atoms per unit cell, and paired spins), and the\n",
"eigensolver gives itself some breathing room by computing some extra\n",
"states (see `n_ep_extra` argument to `self_consistent_field`).\n",
"\n",
"We can check the occupations ..."
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "7Γ10 Matrix{Float64}:\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0"
},
"metadata": {},
"execution_count": 5
}
],
"cell_type": "code",
"source": [
"hcat(scfres.occupation...)"
],
"metadata": {},
"execution_count": 5
},
{
"cell_type": "markdown",
"source": [
"... and density, where we use that the density objects in DFTK are\n",
"indexed as Ο[iΟ, ix, iy, iz], i.e. first in the spin component and then\n",
"in the 3-dimensional real-space grid."
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "Plot{Plots.GRBackend() n=1}",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdeUBU1f4A8O+5M8OO7DsIioCCgoiIJG64W5Zm+jR9VuaaqWmm+Vps8ZWmPTPrZ2ZqWanly9Q0FwQXENwX3AA3iG1A9h1m5p7fH9MjYlHAuXO37+cv53YYDqc58733nPM9h1BKASGEEJIrhu8KIIQQQnzCQIgQQkjWMBAihBCSNQyECCGEZA0DIUIIIVnDQIgQQkjWMBAihBCSNQyECCGEZA0DIUIIIVnDQIgQQkjWpBYIMzMz8/PzW1lYp9NxWhlJwkZrB2y0dsBGawdstPaRWiBcs2bNTz/91MrCVVVVnFZGkrDR2gEbrR2w0doBG619pBYIEUIIoTbBQIgQQkjWMBAihBCSNQyECCGEZA0DIUIIIVnDQIgQQkjWMBAihBCSNQyECMlUdnZ2v1Hjnpw4DZPPkMxhIERIjtLL6azP9yQ6Dj1Sav/Gj/H51XxXCCH+YCBESEYKa+H/brJRv2n77NPa9R7hk/LfLtV3Ct3DA3ZrRh/R/nCHrdDwXUWEjE7JdwUQQpyr0UFMNvv9bXokix3gRhZ2Z57xZkwYf5hwvFGBeac1A9zIND99AX5rjZCRYCBESLJYCol59Ps77M/32CA7Ms2P2TJAZa1qpqSZAsZ0ZMZ0hOJaxW9/sF+nsHMSdE96MdP8mCEehBi95ggZE4eBMD09PS4uztbW9qmnnjIxMWlaoLKy8sCBAzU1NSNGjHB1da2/npCQkJqa6ujoOHLkSFNTUwAoKCi4cuVKfYGePXs6OjpyV3OExO5GMf3+DvtdGmtvSqb5MbcmqFzNW/WDdqYwzY+Z5sdkVtI99+nSc7qCGhjnQ170Z0IdMCAiaeIqEMbHxz/99NMTJ05MTU1du3btyZMnVaq/3YiWlZVFRER07tzZxcXl9ddfP3nyZFBQEADMmDFD/7PJyclvvvnmmTNnbGxskpKSpk2b1rt3b/3P/vvf/8ZAiFBTGRV01126NY2t1cEkX3LyKaW/TTujl5clWdidLOzO3Cimu++z44/pzBUwoTN5wY/pZI0REUkKV4Hw/ffff+eddxYvXqzVanv27Ll3794JEyY0LLBt2zZXV9cDBw4QQpYuXbpq1arvv/++vLx869att27dCggIoJR27dr18OHD//jHPwAgKCgoJiaGo9oiJGqFtfDLfXb7bTalhI7vxGzpr+jnarDxzCA7EmSneDcUEvPo7vts3/1avw5kQidmShfG0cxAvwMhXnEyG15dXR0XFzdu3DgAUCqVY8aM+f333xuV+f3338eOHUsIAYBnn3324MGDAGBqampra1taWgoAdXV1VVVV9UOmFRUV+/bti4+Pr6mp4aLOCIlOtRZ232fHHNV23qU5lk2XhTC5U1SbohRRhouC9RgCUa5kfaTij0mqZSHMxQLqv1sz5qh2+222SmvoX4aQcXHyRJiTk0Mp9fDw0L/08PC4ePFiozLZ2dnu7u71BYqLi6uqqiwsLPbs2TNt2rSuXbumpKTMnz9/4MCB+jKEkB9++CElJaWqqurQoUP+/v7N/mqtVnv06NHy8vKGF2fMmGFnZ9e0sEaj0WhwtXjbYKO1gwEbbc3n/3flRuozc948VOF2IAt6O8KUzrA9CqxULACATqfh+IhyBmCkG4x0g9I6+C2T/nxPt/iMbqQHhGtuHPtuw/gnh02ZON4gvwg/ae2AjdaUUqkkj7oz5CQQsiwLAPW/m2EYna5x76SUNiyg/ymWZVevXt2rV6/JkyffuHFj48aNEyZM6NSp05NPPjlmzBj9T82YMWPJkiX79+9v6bdXVVUVFRU1vKLRaJpWAAB0Ol2z19FDYKO1g6Ea7fbt26t3HCoLfi5+1YZ/vbfyk95gb0L/9yse/+3bxkoBk31gsg9kV8HuDOat2e9VDpiX8O7rzz3zlFJpgC8W/KS1AzZaUwqFgp9AqH/Uy8/P1z8UqtXq+oe/em5ubvn5+fp/q9Vqa2trKyurU6dOnT17Nj8/Xz+gmpycvGnTplWrVukjJQAQQsaPHz937twW/x6lcuzYsfPnz29NPTUajZkZznK0DTZaOxiq0Xw6+9ZWVdkmbvp63Qdjgk0f/w0NwtcM3rQH5T9Hvr/hX0pXPysrK4O8LX7S2gEbrX04mSO0tLSMiIg4dOgQAFBKjxw5Eh0dDQBarVY/agoA0dHRhw8f1pc/fPjwkCFDAMDExESr1dbV1emvV1ZWNs27OH/+vI+PDxfVRkjgvs8w6/dpXP71xDGjR/Fdl8aWzJ+jvpZo99pPMdmU77og1DZcrRp96623XnzxRbVaff369dLS0kmTJgHA7du3AwMDi4uLbW1tZ86c+cUXX8yaNcvNzW39+vX6oNi7d++goKDRo0ePHz/++vXrJ0+e/PTTTwFg3rx5Wq3W29v75s2b+/fvf8i4KEJSVVIHH17WHRmlVKkEmr1gaapaHc4uOqO7Mk6pxF1pkHhw9WkdM2bMkSNHdDpd//79z5w5Y25uDgDu7u7ffvutpaUlADg5OV2+fLlr164qlSoxMbFv374AoFQqT548OWvWrMLCwuDg4JSUlC5dugDAvHnzunfvXl1d3a9fv5s3bw4aNIijaiMkWB9e1o31YYLtBRoF9cb5MG4WsDWN5bsiCLUB0Q9USsaCBQv8/PxaOUdYXl5ubW3NdZUkBhutHR6/0e6W0b77tdfGt3aDGB5dLaIjD2lTJqhsmtlOqg3wk9YO2Gjtg+MXCInA62fZpcEK4UdBAAixJ091ZD66gmsXkWhgIERI6I7n0mtFdEGQaHrrv3srtqWxt0slNdqEJEw0XQsheWIpLDmrWxvBmCr4rkqrOZvDoh6K5RdwphCJAwZChATtm1S2gwrG+Yisq77eg7laSI9hKgUSA5H1LoRkpVwD719i10aI52Hwf0wY+CiceeOcToehEAmefANhbW3tyZMn9Rt8IyRMKy/rRnuRMEdBp0y0ZEInpoMKvruNA6RI6OQbCMdOffn5jTHBA0bwXRGEmnevnG5NYz/sLb7HwXr/6at45wJbjrtAI2GTbyCsqqqusXDJLtME79F+dIW9V44jOEhYlp5jF/cQR8pES8IcyXBPsuoqplIgQZNvINy/Y+um0W4ZCfu3DVDkVdN++7VB/9WuvsrmVGFERPw7nUcvFdBF3UXfQ1eFK75JZdPxRhMJmOi7WbvZ2Ng899xzHh4eYY5kfaQi+3nVpihFThUN/VUb9Zt2/XU2v5rvKiK5Yim8lqRb3YcxE/Gw6J9czGF+oGLZeZwpRMIl30DYSLMHcEf9pv06hS3DGQ5kXNvSWBMFPNdJIt1zSTBz7gGNV+NDIRIoifQ0AzJVwJiOzPZBitznVctCmGPZ1GuHZsxR7fbbbKWW78ohGajQwIpL7Gd9FaJcKtocMwV81JtZmKRjMRQiQcJA2CJzJYzpyPw8RJExWTWhE7P7Puv+o2ZirO63P9g6HOZBnPnoim64Bwl3kkwcBACY5MtYquD7O9hzkBBhIHw0WxOY5sf8Nlx5b5LqqY7k8xus+4+aaSd0v/3BarFfI4PKrKRb0tiVvaXWMQnA+r6Kty/gsAoSIqn1N045mMI0PyZmlPLys8owR7L6Kuvzk3Zhki5BTe/cvbtz566Kigq+64jEbclZdn6gwt1CUo+Der0cySA38gmmUqDH8+DBgx9++FGtVhvwPTEQtoeXJVnYnUkYozw2SmFnCtNP1HYdPO6Fn26+OG8x31VDInY6j57Np6/3kGyv/Dic+b9bbGYlThWi9hs6btLLh9SDn55owPeUbJczjq625L1eipSJJo7mjK5EXcWIOfkZ8YqlsOiM7uNwxlzJd1U442lJXg1UvHkOZxRQ+1USc1qUbW5mZsD3lG6fMyKGYW6dOf7Zkav7VZE6CtJZ7YeM6Ps7rILAJF+J35suDWa6/VeboKZRrthPUJuVa6B65q7Prc+9MPLfBnxbifc6o7Gzs3t/0iBbM8WWVLzbRW1WpYUVF9n1kdK/iTJXwsrezJJzeCgFao/3Lume6mw+Z1y0ubkhh98wEBrSF08o3rmoK6zlux5IbFZd1Q1wJX2klTLRkildGCWBHZhKgdroVgn94Q67koNt6DEQGlKQHXmuE/PBJVwXh9ogq5JuvCXBlImWEIC1EYrl5zGVArXN4jO6t3oqnAw5OfgnufQ9o1kZpvjpHptchAM/qLWWnmNfDVR0tJLF46BeX2cS5UrWJuNDIWqtvRnsHxUwtxsnMQsDoYHZmcK7oYpFZ/ChELXKmXwar6ZLgmXXEz/pw3xxU/dHBd4yoker1cHSc+y6SIWKm44iu+5nBLO7MYU18N/7eLeLHoECLDmr+zicsZTf8m1PSzKnG/PORewm6NHWJLMh9mS4B1ejJhgIDU9B4Mt+itfP4hQIeoQf77BaClO6yLQbLg9RnMil5x7gQyF6mKxKuv6Gbk0fDruJTHsg1/q5kCdcyJpkHCBFLarWwtsX2LV9pJ8y0RILJbzfi3ktCVMp0MO8cY6dF8j4WHPYUTAQcmVNH+bLm3gwN2rRJ8nsEy5E5nnl0/wYLYWf7+EAKWre6TyamEffCOb2iGoMhFzxtCQLgxRv4G5SqDnZlfSLm7qPw+XeARkC6yMVb5xlq3AeATWhozDvtO7TCM4n0eXeDzn1RjBztYgezcaHQtTY8vPs3G6Mt5xSJloS6Uz6OpN11/GWETW26RbrYAbPdeI8TmEg5JCpAj7pwyxK0mmwj6MGLhXQYzks16M9IrK6D/PZdV1OFd4yor8U18IHl3Xr+hqjm2Ag5NZYb6ajFWy8hZEQ/WXhGd3H4QprFd/1EIxO1mRmAPP2Bewm6C9vX9T9ozMTbG+MURMMhJz7T1/Fv6/oHtTwXQ8kDLvuspUa+KdcUyZa8q+eiqPZ9EIBPhQiAIAbxfS/99l3exlp1AR7I+e62ZKpXZi3L2AqBYIaHSy/wH4WqWBwcvDvrFTwfi9mIaZSIAAAeDVR92GYwsHUSL8OA6ExvNdLcTCTnsfEYdlbm8z2cSID5J0y0ZKX/JlaHfyCWzLJ3q67bGkdvBxgvPCEgdAYrFXwQRje7cpdXjWsv6H7SDanTLQVQ2BthGLpObYGR09krFoLb55nP4tUGHOnCeyTRvKiH8MC/IhnsMnY8vO6WV0Z3w74ONiiQW4k1IGsx1QKGfvoqq6/q7FHTTAQGglD4MsnFMvOsWUavquC+HC5kB7OYpeFYMrEI6yJYNZe0+VW8V0PxId75fSrW+wqo280gYHQeMIcyXBP8vEVHPeRo9eSdB+GKTpgysSjdLYmL/kzK/B0a1l6/Qz7eg+Fh6WxR00wEBrVqnDFllQ2rRTnCmXk7t27sz/eXKjOetEfu1urvBOq2Jd066113+Tl5fFdF2Q8sTn0WjFd1J2HboI906hczGFpiOL1s3i3KyNDxz2/U21T/fV0+R4z0UbWKqjdNOWLu2bjX5jNd12QkWhYWJCoW9eXMeVj9gADobEtCGJul8LvmfhQKBfmNnbKOwldvFz5roiYuNvbwO3Tzm7ufFcEGckXN1lvaxjTkZ+QJL+DsflmwsAXTyjmntYNcVfycu+DjGzyZ/tvXr3844vhfFdETK4kHBv09ZXnR4fxXRFkDPnV8PEVXfwY3uIRV+E3LS2tf//+dnZ2vXv3vnjxYtMCtbW1M2fOdHR09Pb2/uqrr+qvf/PNN/7+/o6OjuHh4bGxsfqLOp1u8eLFzs7OHh4eq1ev5qjORjPUgwTakfU3cI24LJzIY8ZH9WAYHH1pAxMTk6f69jiu5rseyCiWn9e95M8E2PA2ecBV53z++eejo6PVavX06dPHjx+v0zWeFfvkk09u3bp1+/bt/fv3v/XWWxcuXACAq1evLliwYOfOnQUFBfPmzXv22Wdra2sBYPPmzceOHbt27dqpU6c2bNhw+PBhjqptNOv6MmuScbt96avRwfkH9AknvOlps4HOulg8wkwGLhbQQ1nsv3ryOT7GSSC8cuVKamrq8uXLTU1N586dy7Ls0aNHG5XZsmXL8uXL7ezsQkJCJk+evHXrVgC4f/9+x44dw8LCAOC5554rKyt78OABAHzzzTeLFy92cXHx9fWdOXPmli1buKi2MXW2JjMDmOXn8ftR4k7n0R72xFqJX+ht1sOOFtXSrEpsOimjAK8m6laFK2xM+KwGJ4EwLS3N39/fzMwMAAghPXr0SE1NbVigpqYmIyMjODhY/zI4OFhfYOjQoUqlcs2aNfHx8YsWLRo/frynpycApKamNi0sdm+HKk7m0ng19nMpi8tho91xtWh7EICBbszxXOwgUvZdGqtlYSrfh7FwMjlZVFRkZWVV/9LGxqawsLBRAQCwtrauL1BQUAAAlpaWU6dO3bBhw8GDBzMyMtauXQsAdXV1FRUV9YU7dOigL9ysmpqaBQsWLFiwoOHFK1eu+Pr6Ni1cWVlJCJ9fUu90ZxYkKk8OqxPRWQS8N5q4xGSq3g/RYaO1Q2VlZZSD8kgGGeeq5bsuoiGuT1qFlrx1wWRHlKaqksNj6iwsLB45Q89JILS3ty8vL69/WVpa6ujo2LCAg4MDIaSsrMzW1hYASkpKnJycAGDXrl3btm1LSUmxsLBIT0/v0aNH9+7dAwICrKysysrK6t9NX7hZZmZmn3/++fz581tTT0ppw4BtfNO7w/Z07c85FjOMuM/6Y+K90USkXAMpZZrB3qaaKh02WltRSkd1Mvv0ls7KyozvuoiGuLrn+2d1o7xgYEcLvivCzdCov7//7du3a2pqAIBSev36dX9//4YFTE1Nvby8rl+/rn95/fp1Pz8/ALhx40avXr0sLCwAwMfHx93d/ebNm/o3bFpYAgjAZ5GKdy7oSur4rgriwIlcNsKZmGGSTHv52xCGwG3ciUmK7pTRb2+zH4YJontwEgh79uzp7++/evVqjUazefNmSumwYcMA4NixY8uXL9eXefnll1etWlVaWnr9+vWdO3dOnz4dAPr06XPs2LEbN24AwG+//ZaZmRkaGqovvG7duoKCgvv372/evFlfWBpCHcgz3sz7uLOiFMXl0Gh30TzrC9MgNxKH04RStCBJ96+eCjf+nwYBuEuo37Fjx8svv/yf//zH19d3z549SqUSAIqKiu7evasvsHTp0oyMjE6dOllZWX344Yfh4eEA8PTTTy9duvSZZ54pLi7u2LHjjz/+6OPjAwCzZ89OS0vr1q2bUqlcsGDB6NGjOao2Lz4KVwT+VzPdn+lhL5rBfdQacTl0c38MhI8l2p0czKSzu/JdD2RQ+zPY9HJ4NVAovYNQKqm7rQULFvj5+bVyjrC8vLx+DQ6/Ntxg92awsaNFsNGPcBpN4PKrIWC3puCfKgXBRmsPfaPlVNGQPdq8KSoRLSjjkSg+aXUs9PhFuz5SMdJTKP9ThRKQZe6VQKagBvakY1qhdBzPZQe4MbjR9mNytyD2puRasaTu12VubTIbZEeEEwUBA6FAKAh8+YRi0Rm2CheKS0VcDo12E1BXF68h7iQuBwOhRGRX0k+v6VYb/ejdhxNWbeQsypVEOJFPr+FDoUTE5tAhHhgIDSDancRmY7+QiKXn2HmBjB9/24o2CwOhgHzal/n8hi6jAm9+Re+PClquoUF2wurtIjXYnUnIoxoMheKXmEfj1XRZiCBSJhrCQCggXpbk1UDFsnPY40UvNodGu+PyDsNwMIVO1uRCAd4gihtL4bUzujURjKXwFgViIBSWpcHMuQf0BCZOiVxcDh2CW4wazhB3gidRiN3mVNZUARM7CzHoCLFOcmauhDV9mGnrfvl623atFlfOiNXxXIp7bRtQtDsTl4MjJWJVUVGxav1X//rh+BdPCHQZNQZCwXHKTMiN3T7/u1M/7tzFd11Qe6SUUBUDna2F2eVFaYAruVBAq/HOUJzeWrn6reM5VT//y6kml++6NA8DoeA4OztblmXR3Ftenh581wW1RywmThialQqC7cnpPBwdFaVOHb3g3vkOUCPYDcGFN2spe127dk2OPxq0s6JPlA/fdUHtEZdDx3fCQGhg0e7keC471ENwCw7RIw2dPNNNMzDlZQ/BBkJ8IhSijq6OYf6eeGavGLEUTqnZQfhEaGjR7kwsptWLU0w2fSqsi2CjIGAgFKxhHkwMJhGL0OVC6mxG3C0wEBrYEy7kVjHFA8vEKCabHSbszSUwEArUMA8Sg+vFRSgON5ThhgkDEc4kXo13hyJTx8JpNR0s7PPIBF05OQtzJLlVNLeK73qgNorNYTFxgiM4OipGiXm0qy2xN+W7Hg+FgVCgFAQGuTPHMHdKVOpYSMqjA12xW3EiGnffFiHhj4sCBkIhG+ZBYrKw24vJ2Xzqb0PshH3zK15hjiS7kqqr+a4HaouYbDrMQ+iBRuj1k7MRHiQmm8VIKCKxOSxOEHJHQaC/K3MCh0nEo7gWUktoX2ehdwoMhMLlY00sVeR6EYZC0YjLodHCXhQgdtHuJA534hWPYzlsf1diKvjkT+y0goZrR0WkSguXC2k/F6Hf/IoaHtIrLqIYFwUMhAI3zINgNqFYxKtpmCMR4BEzUhJoRyo1NL0cY6E4HMumwzxFcGuIgVDQot2ZxDxaq+O7HqgV4nJYHBflGgEY7M7g6Kgo3CmjtTroZouBED0eWxPoZkcS87HbiwDutW0cmEQhFjHZdJgHEUWXwEAodMM8SEwWjo4KXUkd3C6lfQS/Ok4ChnmQY7iaWgz0gZDvWrQKBkKhG+bB4HoZ4Tuew/ZzJSbYn7jnbUUslCSlBDuFoOkonMwVzWSBOGopZ5HO5E4ZLajhux7ooWIxccKIot0J7rUmcOcfUC9L4mbBdz1aB7uu0KkYiHIlx3NxdFTQ4nCC0IhwmlD4RDQuChgIRQFHRwUutwryq2lPB9F0e7Eb4s6cyGV12CcELCabHeYpmvgimorK2TAPchQ3HRWw2Bx2kDvDYBw0FhdzcLcglwuxUwhUuQauFNIo8WwugYFQBLrZEpbC7VLs9gKF46LGh1vMCNmJXLaPE7EQz+YSGAjFYQjutSZgx3PxMF5ji3Yncbj7tlCJZWe1emKqq5zhpqOCdbeM1ukgwAYDoVENcsNNl4RLXCtlAAOhWAz3YI7nshq8Axae2Bw6BI+kNzobE+hqS84+wLtDwcmupAU1Ils7hoFQHBzNoJM1uVCA3V5w4nJoNAZCPgzB0VFBOppNh4ht7RgGQtHA0VEBogAnctnBGAj5EO3O4HoZARLduChgIBSRYR4MHskkNNeLaAcT4m0lsm4vDf1cyJVCWqHhux6oAQoQm8OKbrIAA6Fo9HclVwtpaR3f9UAN4IkTPLJQQi9HkpCHD4UCcqWQ2pkQH2uRdQoMhKJhpoAIZ3IS91oTEpwg5Fe0O4PThIIixnFRwEAoLrjXmqDoKCTksYNxr23+DMHdtwUmJpvFQIi4hetlBOX8A9rRkjiZ8V0PGYtwInfxbBbBqNHB2Xw60E18YUV8NZazng6kuI6ml2MsFIS4HNxQhmdKBvq5kFNqHB0VhHg1DbYnNiZ816PtMBCKCdFPiuRiIBSEuBzRnDsqYZhEIRwx2ay4dlarJ8pKyxmOjgpEjQ7OPaD9XfGJkGd4SK9wiHSlDABwuD34xYsXDx486Ojo+Pzzz9va2jYtkJeXt2PHjqqqqmeffbZbt24AUFpa+tNPPzUsExUVFRgYmJGRceTIkfqLo0aN8vLy4q7mQjbCgyw7p2OpQlwbN0hPYh7tbkc6qPiuh+yF2JPCGppdST0ssUvwqaAG0stpuJMo/y9w9UR48ODBYcOGAcDJkycjIyOrq6sbFSgoKOjVq9eNGzdqamoiIyPPnTsHADU1NRf/JzExcfbs2bm5uQCQnJy8YsWK+v9UWlrKUbWFz8OSOJqRK3gSG9/iclicIBQChsBAN+Y4zhfw7Wg2O8iNUYlzkJGrJ8KVK1euWrVq1qxZlNKIiIiffvrpxRdfbFhg8+bNoaGh33zzDQCYmJisWrVqz549Li4umzZt0hf46aefTp48OXjwYP1LX1/f+v8kc/rR0V6O+C3Mp7gc+u9wBd+1QAB/HslEp3bhux7yJt5xUeDoibCysvLMmTOjR48GAELIyJEjY2NjG5WJjY0dNWqU/t+jR49uWmDr1q0vvfQSw/xZw6Kioq+//nr37t2FhYVc1FlEhnkQ3GuNX+UauF5M+4pzFEh6cJpQCGLFHAg5eSLUj2c6OzvrX7q6uiYmJjYt4+Liov+3i4tLWVlZRUWFlZWV/kpWVtbx48e//vpr/UtTU9POnTtfv349JSVl7ty5R44cCQsLa/ZXa7XavXv3ZmRkNLy4aNEiBweHpoVrampUKvFN8vS1gykPFEUVNbwcAC3SRjOsmCzS24EQbU2NtlXlsdHaofWN5m0KOlZx40GNrzXXlRI6vj5pKaVAQOFlWlsjvJxOExOT+geqlnDyVar/rZT+eY+m0+kUisaDSAzDsOyfjzX6fzQss2XLlujoaG9vb/3L4cOHDx8+XP/vJUuWvPnmmzExMS39dgsLC3t7+4ZXVCpV0wrof2Oz1wXOVgHB9nCmUDHMnYffLtJGM6xT+TDYDVrfDtho7dCmRhvsBifzFf7NrMmTF74+acfzYKh7G3qEMRHy6OdUTgKhm5sbISQ3N9fHxwcAcnNz3dzcGpVxd3fXPzgCQE5Ojp2dnbm5uf4lpXT79u0ff/xxs28+ZMiQ3bt3t/SrlUrl8OHD58+f35p6qlQqkd6nD/dkj6vpaG8ePnbibTQDOq7Wfh2lUKlaOxCEjdYObWq0oZ7s75n0lSAhfhEbE1+ftDi19p9dGJVIl8pwNEdobm4+ePDgvXv3AoBWqz1w4IB+vrCmpubChQs6nQ4ARo0atW/fPv1T4759++rnCwEgNja2pKTkmWeeqb+i1f41AnXs2LGAgAAuqi0imE3Io8JayKigYbhYSUiGupO4HJbFPsEHLQsJairqTXe5mmVasWLF2LFjb58uIA4AACAASURBVN++nZKSYmlpOXbsWAC4f/9+eHh4cXGxra3tSy+99NVXXz3zzDMuLi6//vrriRMn6n92y5YtU6dONTU1rb8yffr04uJib2/vW7duXbt2rWFOoTz1cSKZlVRdDa7mfFdFfmKz2YFujFLEvV6CPCyJvSm5XkyD7fEGxdiS8mmXDuLedJerQDhgwICLFy/GxsZGR0c/9dRTJiYmANCxY8ejR4/qV8TY2NhcuHDhwIEDVVVVH3zwQcOx01deeaXRM98nn3xy+vTpgoKCoUOHRkdHd+jQgaNqi4WCwEA3JjabndIFv4+NLQ7PIBQk/dpRDITGJ9ITJxricN1hp06dZsyY0fCKpaWlPstez8rKatKkSU1/sH///o2uuLq6jh8/notKipd+dHQK5k4ZXWwOnReI9x+CE+1Ott9mF3XH/zXGdjSbfizynFr80IiVPhDinIiRZVbSMg3tjo8dwhPtzsSrqQYzbI2rpA5uFdMnnMXdIzAQilWXDsRUAbdKMBQaVWw2jXbHfV6FyMEUOlmTiwXYI4wqLod9woWYivuBEAOhmA31IDFZ2O2NCicIhQy3mDG+mGwq0qOXGhL9HyBnuNea8R3PpdHuGAgFKtqdicvBHmFUot5itB4GQhEb6s7Eq2mtju96yEZKCWUI+HYQfbeXqgGu5PwDWt26fe/Q40svpxWSmDLHQChidqYQYEvOPsCxICOJzaFD8XFQwKxV0MOeJOZjjzCSo9l0uIcUpswxEIobjo4aU1wOjosKXbQ7wdFRo4nJpsM8pdAjMBCK2zAP5iiulzEKlsIpNTsIV8oI2xB3Jg7XyxiFjsLxHHaomHdWqyeFv0HOnnAhKSW0qJbvesjAlULqbEY8LDEQCtoTLuRmMS2t47seMnCxgLpZEDcLvuthCBgIxc2EgX6u5DiOBXEvFsdFxcCEgT7O5JQaewTnpLFeVA8DoegN82DwJAojiMthMRCKQjSOjhpFTDYrgQxCPYn8GXI2zIMcwUDIMQ0LiXl0oBv2FxEY4k4wEHKtUguXCmh/V4ncGmLHFr3udkTDwr1y7PkcOptP/W2IvemjSyLehTmSzEqaV813PSTtZC4NcyRWUjltGgOhFAxxx3N6uYUThCKiINDflTmRi9OEHJLSuChgIJQGPLCea3E5bLQklonLRLQbjo5yKyabDpdEBqEe9m0pGObBxOWwOuz43KjSwuVCGiWV6RA5GOKBu29zKKeK5lTRUAfp9AgMhFLgYg6eluQC7rXGjXg1DXUglhweYo0MLMiOVGpoOk6ccyMmmw51ZxTSiYMYCKUCR0e5g+OiokMABrkzx3OxR3BCShmEeti9JWKYB4ObjnIkLocOwZUyYoNJFByhAHE57FAMhEiABriSS4W0TMN3PSSnpA7SSmkfZ0l1ezmIdiexOSxGQoO7VkQtlaSTtaR6BAZCibBQQrgjiVdjxzew4znsEy7EBDuK2HS2JmYKklKCPcLApDcuChgIpWSYJ46OGl5cDsUJQpGKxtFRDsRksxgIkXAN9yAxeCSToeEEoXhhIDS4OhaS8uggye01KLW/R85CHciDGppZiT3fYPKqIa+a9pRQvpSsDHFnjudifq0hxatpkB2xk9xegxgIpYMhMNidicUkCsM5ls0OdGMYjIPi5GIObubkSiH2CIOR5LgoYCCUGMwmNKw43GJU5IZ44OioIcVkUyltMVpPgn+SnA3zIMdwybjhxOXiBKG46ZMo+K6FRBTWwt0yGiHFVCIMhJLibUVsTEhyEYZCA7hbRmt1tKutBLu9fAxyYxLzaK2O73pIwrFsdoAro5Ji0JDi3yRvODpqKHE5dAgmToicrQkE2JBzuA2vIUgyg1AP+7nUDPMgMVk4FmQAOC4qDUM8cHTUMGJzMBAikRjiziTl02ot3/UQOQpwIocdjIFQ/KLdGVwv8/jSSqmWBanOFGAglBprFfSwJ6fzsOc/lutF1EpFvK2k2e1lJcqFXCmklXhr+HiOZtPhEn0cBAyEkjTMg+Bea4+jpKTk3bVfdi8+x3dFkAFYKMEn/fC/1m2pra3luy4iJuEJQsBAKEnDPBhcL/M4Xl36zr7bFbGfvFJWVsZ3XdDjunr16p1f1v9fTPIXX23muy5ipWXhVK6UT+WU7B8mZxFO5H45za/mux6i1cnbi9xNslaBmZkZ33VBj8ve3t6s6gFkX/P28uC7LmJ19gHtZE2czfmuB2cwEEqQkoEBbkwcrpRrr3+++obb1H/fuZxkYmLCd13Q4/Ly8rqWdMLkxU1jnhnHd13ESqo7q9XDQChNmE34OJLy6YDQbubm0r0BlhkPJ7uALp1w09F2k+rOavWk/LfJ2TAPchQDYXsl5dFIKe4jJWeRLiQpH3tEe5Rr4FoR7eci5R7xiEBYWVmZlpZWWlpqnNogQwmwIUoGUkux57dHUj6NlHS3l6FIZwyE7RSXw0Y6E3Ml3/XgUouB8Ntvv/X397eysgoICLC1te3YseNnn31GKX6SRGOIOzmK5/S2XbkG7pbRYHsMhJLS15kkYXJtu8Rk02GeEh87bP7PW7NmzUsvvWRvb79mzZrt27evW7cuICBg0aJFr7/+upHrh9oNpwnb59wDGupITCTe8WWnSwdSx9IsPLa67Y5KOoNQr5nH3bq6upUrV86ePfurr76qv/jaa6+tXLnyvffeW758uZOTUyvfvbS01MbG5iEFqqqqVCqVSqVqzbvV1NQwDIML+VppqAczO0GjYRWS3C2eOzhBKFV9nJgz+fS5Tvg/tw0yKmhpnfQHSJr5jnzw4EFZWdm8efMaXX/llVd0Ol1GRkZr3vf8+fN+fn5+fn6enp7Hjh1rWqCsrOzJJ5/08PBwdHR877339Be3b99u/3eZmZkAUFtbO2nSJFdXVycnp4ULF+IIbWs4mIJvB3IW50XaKCmfxUAoSbheph1isulQd0by/aGZQGhlZUUIKSwsbHS9oKBAoVB4e3s/8k0ppdOmTXvttdfy8/M3bNgwZcqUppsb/fvf/waABw8e3Lx5c9OmTSdPngSAf/zjH3f/Z9myZZ06dfLy8gKADRs2ZGZmqtXq+/fv//7777/88kv7/lq5wb3W2ooCnHsgzaNHUaQzOYOBsI2kvbNavWYCoY2NzZAhQ+bPn3/37t36izk5ObNmzZo3b15rxkXPnz+vVqtnzZoFAOPGjevQocPhw4cbldm+ffuiRYuUSqWHh8eUKVO2b98OAKampnb/8/PPP7/88sv6wt99992rr75qZmZmb28/Y8YMfWH0SLjXWlvdLqVWKuJuIf2eL0PhTuRqIR7S2wYsheM57BAZBMLml8Ru27ZtyJAh/v7+ISEhrq6uBQUFV69eNTExcXJymjhxor7MnDlzoqOjm/3xu3fvdunSpX7mr2vXrg1jKgBUVVWp1epu3brVF/jhhx8aFkhOTr558+bkyZP1L+/du9ew8EMCIaW0sLDw3r17DS96eXm1chpSYqJcyI1iWlwLdqZ8V0UkkvJxglCyLJXgb0OuFOITf2tdLqSOZsTLUvrN1Xwg9PT0vHjx4vfff3/ixImsrCyWZXv06AEA9+/fry/zkP2IS0tLLSws6l9aW1sXFxc3LFBSUgIA9WWsrKz0V+p98803zz77rJ2dHQDU1dVVVVU9pHBDdXV1X3zxxXfffdfw4oEDBzp27Ni0cEVFRUvvIxnhDqrD9+ue8jDYAKm0G+1UlqqXDS0vN/CZPdJuNI5w0Wi97VUnMjWB5pI9k8mwjXbgnnKQM5SXi3vbYgsLC4VC8fAyLSZJWllZzZ07d+7cue34xU5OTg1z8EtKSpydnRsWcHR0JISUlpbqQ11JSUnDEde6urqdO3fu2rVL/9LExMTGxqY+7jYq3IipqemKFSvmz5/fyqpaW1u3sqRIjfJmE4ro5K6P+By0iYQb7UKxdm4PhbW14W+BJdxo3DF4ow3wYH/7g1pbS3nzPAM22qkC7WJuuoPQcLKyPjAwMDU1VX9vwrLs5cuXu3fv3rCAiYmJn5/fpUuX9C8vXrzYsMCePXusrKwGDx5cfyUoKOjixYvNFkYPN8yDHLqejSextQam0kseptW3Xl5R2dn7Bf1dZdEdOAmE3bp169Onz5tvvpmVlfXBBx/Y2toOGjQIAPbv3z9jxgx9mdmzZ3/wwQdpaWlHjx7dvXt3/XUA2LJly/Tp0xnmr7rNmTNnzZo1165dS0xM3Lx5s34ZDmqNpL3bMz57sVNI3/Lycr7rInSYSi95mFbfSmlpaV0jBtateyb5bALfdTEGrvaP27lz56JFiwYNGhQYGLh//35CCAAoFApT0z+XbSxcuLC4uHjcuHFWVlZbt24NCgrSXy8rK2MY5oUXXmj4blOnTs3NzX3++edNTEzWrFkzYMAAjqotPddT0mjg0PIbB0tKSnB07uHO5NO+TrK4/5UzTKtvjaysrGqHLtTK+fbde1FRUXxXh3NEYsnpCxYs8PPza+UcYXl5ueRjQ2lp6dQPvq51DTz6xpMGeUMJN9pTR7QvBzDjfAz/SCjhRuMOR4328VW2sIaujTDkrLlwGKrRKKV+i78d4VT92Rsz5bDkHoeBJM7Gxuajt5ak+4zguyJCRwHOYiq9DOAxFK1RoyPqkKlrlr4ihygIGAjloLs9KaqhuVV810PYMJVeJjCtvjXOPaDd7YmFpI9eaggDofQRgEgXkpiHe609DKbSy0R9Wj3fFRG0hDzaX05HcmIglIUoFyYBV40/FB46IR+4+/YjJajZKHkkTuhhIJSFKFcSr8ae/zBJ+bQvBkJ56IvThA+lo3Amnz7hIqPoIKM/Vc56O5LUUlqu4bseQqVPpQ9xwEAoC5GYVv9QyUXUzYI4mfFdDyPCQCgLpgro5YBn0LQIU+llBdPqHy5BTWWyoUw97PpyEeVK4tW4XqZ5OEEoN/q0er5rIVDxaiqrCULAQCgfUS5MAk4TtuAMnkovM5EuOEDSotN5NEpOS0YBA6F89HMl5wtoHT4TNoGp9DKEafUtuVtGCYFOMjhxoiEMhHLRQQVdOpBLBdj5G8NUehnCtPqWxMtvghAwEMpKlAvBbMKmMJVehjCtviUJ8hsXBQyEshLlSnCasClcKSNPmFbfrAT5rZQBDISyMsCViVezLPb9v8NUennCtPqmHtRAXjXtbie77oCBUEbcLMDWhKSUYuf/C6bSyxam1TcVr2afcCEK+fUGDITy0h9HR/8OU+llC9Pqm0pQ0yhXOXYGOf7NcobThI3gBKGcYVp9I3I7dKIeBkJ5iXIh8Tgc1ACm0ssZptU3VKmFWyW0t5McuwMGQnkJsCVVWpqJw0EAgKn0sodp9Q0l5dGeDsRMwXc9+ICBUF4IQD8X5jSOjgIAptLLHqbVN5SQx8pzXBQwEMoQptXXw1R6mcO0+oZku1IGMBDKEB7SWw9XyiBMq9fTsnD+gXwTajEQyk4vB5JeTotq+a6HAGAqPcK0er1LhdTHmtib8l0PnmAglB0lA+FOuFgOU+kRAKbV/48MD+NtCAOhHOEhvYCp9AgAMK3+f+LVctxrux5+DchRlAuD62VwghDpYVo9BUjMZ/vhEyGSlUgXcqWQ1sh71Tim0iM9TKtPLaEWSuJlKd/ugIFQjiyV0M2WnH8g386PqfSoHqbVx6tlurNaPQyEMiXzbEJMpUf1MK0+IU+OZxA2hIFQpqJcSYKM18tgKj2qh2n18jyMtyEMhDLV35VJzKM6ufb9pDzMIER/kXNafW4VlGloN1tZdwcMhDLlZAYu5uR6sUw7Pz4RoobknFZ/Ss32c2Fk3hkwEMqXbA/pxVR61Iic0+pxXBQwEMpZlKtM18tgKj1qRM5p9bI9jLch/DKQrygXme6+jan0qCl5ptWXaeBOGQ11lHt3wEAoX74dCKVwv1x2nR9T6VFT8kyrP62m4Tg6goFQ5mQ4Ooqp9KhZ8kyrT8hj5bzXdj0MhLIW5SK79TKYSo+aJc+0ejkfxtsQNoGsyfCQXkycQM2SYVp9rQ4uFWJCLQAGQpkLtie5VfRBDd/1MCJMpUctkVta/YUCGmBDrFV810MAMBDKmoJAX2eSmCejvdbwiRC1RG5p9TI/jLchDIRyF+XKyGd0FFPp0UPILa0+Xs3K+TDehpTcvfXBgwcPHjxoa2v7yiuveHp6Ni2Qmpr6zTff1NTUTJgwYcCAAfXX09LSvv3226KiIn9//7lz55qbm9+5c+eXX36pLzBhwoTOnTtzV3NZiXIhy87L5YkQU+nRQ9Sn1XvK4GQ+CpCUTzf3x84AwN0T4Y4dO2bMmBEaGlpZWdm3b9+ysrJGBbKysiIjI83Nzbt16/bMM8/Exsbqr8fHx0dERGi12t69e6enp5eXlwPArVu3NmzYUPw/Go2Go2rLUIQzuVFMK7V818MoMJUePZx80uqvF1EHM+JmwXc9hIGrJ8LVq1d/+umnzz//PABcvXr1xx9/nDt3bsMCmzZtGj58+AcffAAAtbW1a9asGTJkCKV09uzZn3zyycyZMxu9oY+Pz6pVqziqrZyZKSDYnpzNp9Hu0o8QZ/LZlwPwFhi1SJ9W/1wnvuvBPTyMtyFOvhTKy8uTk5Ojo6P1L6OjoxMSEhqVSUhIGDx4cKMCmZmZqampUVFRn3766fr163Nzc+vLq9Xqd999d/369ffv3+eiznImk0N6MZUePZJ80urxMN6GOHki1AcwR0dH/UsnJ6fjx483KqNWq52cnOoLVFZWlpWVpaenm5qaTp06dcqUKampqcHBwZcuXfLy8urQocOgQYMsLCwuXbr09ttv79+/vz6INlJXV7dly5YTJ040vLhmzRpXV9emhaurqxUKxeP9rVIQbkc2pjJVAa3KJRZvo90uA0ulyhaqq6qM/avF22g84qXRgizhaqGquLzKVJz/u1rfaPG5qje7aYzfF4zPzMyMYR7xyMdJIDQ1NQUAjUajVCoBoK6uztzcvGmZuro6/b/r6uoIISYmJkqlsrq6evXq1UOHDgWA3NzcjRs3fvTRRwMHDhw4cKC+cOfOnd999934+Phmf7VCoejVq9fIkSMbXnR0dDQzM2taWKPRNHtdbgZ7wsuJOqWJStmKAQLxNtrlLBrpTHmpvHgbjUe8NJoZgL8Nm1JpKtKRg1Y2WkYF1VC2u7MsPpOPjILAUSB0dXVVKBRZWVl+fn4AkJWV5e7u3qiMh4dHVlaW/t+ZmZn6WKVfXKr/KQDw9/fPzs5u9IPh4eHbtm1r6VcrFIrQ0NCJEye2pp4Mw7SmjSTPwRw6WrHJJaR3KzahF2+jncnXRbrwU3nxNhqP+Gq0SBd6toBEinPjsVY22ul8tr8rfib/wklDmJqajho1aseOHQBQXV29d+/esWPHAkBFRcVvv/2m1WoBYOzYsbt379bpdACwc+dOfYGOHTuGh4frn/ZYlk1ISOjRowcA1C86pZT+97//DQkJ4aLaciaHQ3oxlR61hhzS6hPUFDMIG+Jq1egHH3wwcuTIK1eu3L17t2vXrqNGjQKAzMzMp59+uri42NbWdurUqdu2bYuMjHR2dr569eqpU6f0P7h27doJEyYcPHjw7t27JiYmr7zyCgDMmjXrzp07HTt2TElJ0el0Bw4c4KjashXlSvak09e6810PzmAqPWqlSGfyltQzaxPy6Kyu+Dj4F64CYWhoaEpKyunTpx0cHCIiIvTP4L6+vtevX+/QoQMAmJubnzp16vTp01VVVf3797eystL/4IABA27evJmYmOjq6hoWFqb/wW+//fby5ct5eXnu7u6hoaEqFe6OZ2ADXMlrSToKCqkGCkylR60k+bT6olrIrKDB9tL869qHw51l7OzsnnrqqYZXTExMgoKC/vrdSmX9EpiGHBwcxowZ0/CKmZlZZGQkR/VEAOBpScyV5HYp9beRZvfAVHrUevq0+uc6SfMDk6Bm+zqT1qyMkw9sDPSn/pI+pBdPpUetJ+3T6hPy8AzCxrA50J8kfEgvptKjNpF2Wj0eOtEUBkL0Jwkf0oun0qM2kfBp9dVaSC6ifZywL/wNBkL0pyA7UlxLc6W40wQmTqA2kfBp9Wcf0B72xILDxSGihIEQ/YkAPOHCnJbiIb14Kj1qK6lmE+K4aLMwEKK/SHV0FJ8IUVtFukgzEOJhvM3CQIj+IsljKDCVHrWDJE+r11E4+4BGuuDXfmPYIugvvZ3I7VJaWsd3PQwKU+lRO+jT6rMrJRULrxZSDwviJIutttsGvx7QX0wY6OUotQyqMzguitpFeqfVx+MEYQswEKK/iXIhCdJaL5OUh6n0qD2kN02Ih/G2BAMh+psoV0ZKafWYSo/aTXpp9afzcKVM8zAQor95woVcKKB1UnkmxFR61G4SS6u/U0YZQnyssS80AwMh+psOKvDrQC4WSORGGBMnULtJLK0+Xk0H4LhoCzAQosakdEgvptKjxyGltHo8jPchMBCixqIkdAwFPhGixyGl9TIJebhktEUYCFFj/V2ZBDXLir/7Yyo9ekySSat/UAMPamiQHfaF5mEgRI25moOdKblVIvr+j6n06DFJJq3+VC77hDNhMA62AL8kUDOkcUgvptKjxyeNtHo8jPfhsGlQM6RxSC+m0qPHJ41pQjx04uEwEKJmSOAYCkylRwYhgbT6Cg2klNIwR+wLLcJAiJrhb0OqdTRTzFMjmEqPDEICafVJ+TTUgZgp+K6HgGEgRM0gAFEu4t5rDRMnkEFIIK0+Qc3iuOjDYSBEzRP76Cim0iNDEXtafbyaRuEZhA+FrYOaJ/ZDevGJEBmKqNfLaFi4UEAjcU+Zh8JAiJoX6kAyymlRLd/1aBdMpUcGJOq0+ksFtLM1sTXhux7ChoEQNU/JQB9nkijO/o+p9MiARJ1WH487q7UCflWgFkW5MCI9pBdT6ZFhiTetPkGNh/E+GgZC1KIo0R5Dgan0yLBEOk1IARLz2H44QfgoGAhRiyKdydUiWq3lux5thKn0yOBEmlafUkKtVMTTEvvCI2AgRC2yUEKgLTkvtkN6MZUeGZxI0+rjcWe11sFAiB5GjIf0YuIEMjiRptXjYbythIEQPUyUKxHdehlMpUdcEGNaPR7G20oYCNHDRLkwiXlUJ6ruj0+EiAuiWy+TXUnLNTTAFvvCo2EgRA/jaAZuFuRakWj6P6bSI45EOhNxZVDod1bDntAaGAjRI4jrkN7zmEqPuNGlA6nViSmtPiEPMwhbC78w0COI65BeHBdF3BFXWj0extt6GAjRI0S5klNq0ayXwVR6xB0RTROW1sHdctoT5whaBwMheoTO1kRByL1yEfR/TKVHnBJRWv3pPNrHCecIWgvbCT1aP5GMjmIqPeKUiNLq8TDeNsFAiB5NLIf0JuVjBiHikKUS/ESSVh+fh4fxtgG2FHo0sRzSm5SHK2UQt0QxOlqrgyuFOEfQBhwGwurq6rS0tIqKipYKUErv3LmTn5/f9D9lZGTcu3ePZf+2RuP+/fs5OTmGryh6lGB7oq6i+dV81+NRcMko4poo1sucf0C72hBrFd/1EA+uAmFMTIy3t/fEiRO9vb1/+umnpgXy8vLCwsKefPLJHj16zJo1i9I/P1t37tzp1atXnz59Ro4c2bt3b/3FsrKyAQMGDBkypHfv3hMmTNBoNBxVGzWLIRDpQk4Le681TKVHRiCKtHo8jLetOAmEOp1u5syZn3/++ZUrV3799dc5c+ZUVlY2KrNy5cquXbumpqampqbGxMQcOnQIALRa7dNPPz1u3Li8vLy0tLSYmBh94XXr1llYWNy5c+fu3btpaWk7duzgotroIaJcGCGPjlZXV0+et9T6wNtsXQ3fdUFS1qHmQf7Wha+9t7r+3l2AEtQsptK3CSeBMCkpqbKycsKECQAwYMAAd3f333//vVGZHTt2zJ07FwBsbW0nT568c+dOAIiNja2oqPjXv/5VVVUFAA4ODvWFZ82axTCMubn5iy++qC+MjKm/sNfL7N2790g2FFRqfvvtN77rgqRs3f99XevcbcuR85cvX+a7Ls1jKSTl0364UqYtOGmsjIyMzp07KxQK/csuXbqkp6c3LFBeXl5UVNSlS5f6AhkZGQCQkpLi6ek5YsSIwMBAZ2fnb7/9Vl8gMzOzaeFmsSybmZl58e9qa2sN/BfKTx8ncquEVgh1TDqkVzi9ddwpMyEsLIzvuiApGxk90Pb8t9oSta+vL991ad71YupkRlzM+a6HqCi5eNOKigozM7P6lxYWFuXl5Q0L6F/Wl7GwsCgrKwOAwsLCM2fOHDp0aMSIEefOnRs4cGBUVJS3t3dNTU3Tws3SarU7d+48duxYw4s7duzw9PRsWriyspIQHEBore42qhN/VPW2FGKjxVR7DP/i1K4oDcMwD1mfxRf8pLWDMBstrFfo3bPHnjhqkVio7a8Q4iftWK6yrwOpqMA5gj9ZWFgwzCMe+TgJhM7OzsXFxfUvi4qKXF1dGxZwcnJiGKa4uNjOzq5hAVdXVy8vrxEjRgBAnz59goKCEhMTu3TpYmdnV/+GTd+tIRMTk6VLl86fP7819aSUWllZtf3vk6lBHrqLpaqBzjqhNZqOwpe3tdsGKDp0EOhtMH7S2kHIjbYkhP08jRnVmZPvz8dBKT1fYjLCk1hZmT26NPofToZGQ0JCUlNTS0pKAECr1V64cCE0NLRhAZVK1b179zNnzuhfnjlzRl+gV69eFRUVOt2fOzeUlJRYW1sDQGhoaH3hpKSkRu+GjCPKlRHmIb2/3GftTaEfnsSNjGWaH3OjGC4LMrP+NO613Xac3NF07tx5xIgRs2bNWrJkybZt23x9fSMjIwFg586du3bt2rdvHwDMnz//nXfe8fDwSE9P37dvn37muW/fvgEBAYsXL54+ffqePXtqa2uHDh2qLzxnzpzAwMDS0tKtW7fGxcVxUW30cP1cyOQ4qhFeKPwkmV3RC5cGIONRMTA/iFmbzP44WMF3Xf7mj0qiYalvBwyE9W1FsAAAF7hJREFUbcPVo/327dtXrFjxxhtvBAQE7N+/X3/Rzc2t/mFuxowZGo3m/ffft7W1PXDgQOfOnfXX9+3b99577y1cuDAgICAhIUH/RPjMM8+Ul5evWbPGxMRk586dvXr14qja6CFsTcDHmlwrYQba8F2VBo5l00otPOmFgRAZ1ZxuTKddmrtljKCiTmIB098V+0KbESFnw7TDggUL/Pz8WjlHWF5erg+0qJXmJeq8TOveDBPQVNywQ9qpXZgX/ATd+fGT1g7Cb7Tl53UVGtjwhIAeCqcfr+npbLIgSNDdQYCwvVAb9Lau3HfoWGlpKd8V+dPVIppSApN98WOMePBad8WOu+wDwSzPVKvVx2KP97UXapKTgOE3CGqDzxdMOnsmKWrkM3xX5E+rrrKLujN46BrihYs5PNeJ+fKmII5l0mq1oQNGZJ05uOH9pXzXRXzwKwS1gYKyxNSiqEYQC2buldPYbHZmV/wMI94s6cFsvMVWavmuBwAAVGhYxtRCpxVGbUQFv0RQGxzbu+udoZ3gld1VAuhra5PZ2d0Y3GIf8cjPhgxwZbak8n9rmFerUC3ct+65nlu/+JTvuogPBkLUBra2tq9PfTrK1/E/13ju+fnVsOseOy9QQOsUkDwtC2HWJrO8pxW9eZ6d39/7xefGNNzVC7USBkLUZqvCmc+u63Kr+KzD5zd0k30ZVwEtX0Uy1duR+NnAT/f4jISXC2lsDrskGO8L2wkDIWqzTtbkJX9mxSXe1ghUamFzKruoO356kSAsC1F8kszymIi25KzuwzAFThO0G36VoPZ4K1SxP4O9VsRP3990i412Z7oIKZEZydlwD6Ji4FAmP93h13S2oAZe9Mcv8/bDtkPtYWsCb4cqXjvDw0OhhoX1N9g3euBHFwnIGz2Y1Vf56Q5vnmfXRCgUeFv4GPDbBLXTnK6MugqOZBn7LvjHO2xXG+jliP0eCciEzkxOFSTmGbs7/N9NtksHGO6B3eGxYCBE7aRk4ONwZtEZndaIqwQowNpr7LIQXBSAhEVBYFF3Zk2yUZfMlNTBx1d1q/tgd3hcGAhR+z3tzXhYwrY043X+3zJYcwVEu+P9LxKc6QHM2QfszRLjPRR+eFk3zofpbofd4XFhIESPZU0fxbsXdWXG2t3wk2R2WQh+aJEQmSnglW6KT431UHi/nG6/za7ohY+DBoDfKeix9HQgwz2ZNcnGWCaQoKa5VTDOBz+0SKDmBTJ7M9jMSmM8FC49xy7uocBUWoPA7xT0uD7qzWy8yf5RwXnnX52sWxbC4Oo4JFh2pvCCH/P5dc4fCs/k07P5dCEet2Qg2I7ocXlYkrmBzDsXue38t0rohQf0n13wE4sEbXEPZlsaW1LH4a+gAEvO6j4OZyy4OlhddvBrBRnAsmDFsWx6oYDDh8JVV9mF3RXm2PORsHlakqc6MhtvcXhfuOsuW6XFYzgNCZsSGYCVCt7rxSw5y9VMYVYlPfgHO6cbflyRCLwZwnx+XVfNzQktdSy8c5H9LFLB4ByB4eA3CzKM6QFMUS3sy+DkRvjTa+z0AMbWhIv3RsjAutqScCdm+x1O+sK6a2ywPRngimHQkDAQIsNQEPisr+L1s2ydobt/US1sv83Ox3UBSDze6sl8cpXVGXqu4EENfHpNtyoc+4KBYYMig4l2J34dYJOhZ0e+vMmO82G8LPEWGIlGhDNxs4A96QbuC+9d0k3pwvjbYF8wMAyEyJD+01ex8oquqNZgb1ijg423dItxi20kNstCmI+uGPJsptRSuvse+3YoZtAbHn6/IEPqZkvGejMfXzHYqpktqWxfZybQFm+Bkcg81ZHRshCXY7BQuOSsbnlPhYOpod4P/QUDITKwD8MU395m75QZoP/rKKy7zr4RjJ9SJD4EYEmwwc5mOpFLbxbDK7hwmhvYrMjAnM3hte6Kty4YYHbk53uspyVEOuPjIBKl532ZtFK4+Nj5tSyFJWd1n/RhTHFYlBsYCJHhLe7OnM2nCerH7f9rr7HLgrHrI7FSMbDQEGczbb/NmjDwbCf8uuYKtiwyPHMlfNibWXLusVaPH8miGhZGeuHjIBKxWV2Z47mPNVNQrYUVl9i1EbjJLocwECJOTO3C6Fj4+V7774VXX9UtC8bdM5C4WSphVlfmP9fa3xE+SWb7uZAnXLArcAgDIeIEAVgboXjzPFvTrrUC5x/QO2UwsTN+PpHozQ9U7LrHqqvb87N51fDlTd1HvbEjcAvbF3FloBsJsScbbrTnXnh1MrskmFHhxxOJn7M5TOrMfHGjPbeE/zqvmx7A+Fjj4yC38JsGcejTCGbtNV1BTdt+6m4ZPZXLTvfHDyeSiDeCma9T2HJN234quYgezGSXh+B6Mc7hdw3ikG8HMqkz8/6ltt0Lr05mXwlkrFQcVQohY+tkTaLdmc0pbRsdeeOs7t1QhQ3uNc89DISIW+/1Uvx8n71Z0tpVc3nV8Mt9dl4g3gUjSVnek1l3vQ1b0v+eSf+ohJld8SvaGLCVEbfsTGFpsGL5+dZ+Aay7rpvShXEy47RSCBlbiD3pZgs7Wnc2k47CsnO6tREKnCY3DmxmxLn5QczNYnos+9EPhWUa2JLKLu6OH0skQctCFKuT2dbsw705hXUygycxidZY8BsHcc6EgY/CmTfO6R75FfDVLXaEJ66RQ9I0xJ3YmMCBPx7xUFiugQ8vs2sicHbAeDAQImOY0ImxUsEPDx0XqtXB5zfY1/HEJSRdr/dgPrr6iED48RXdCE8S5oi3g8aDXzrISNb2Ubx1ga3UtljghztssD2EOmD/R5L1rA9TVAun81ocG8mqpF+nsO+H4TezUWFzIyOJcCb9XMi6FvaaogD/wS22kdQpCLzeg1nd8kPhm+fZV4MYL0u8HTQqDITIeFaFM+tv6HKrmvlPe9NZSxUMdMP+jyTuRT/mYgG9XtzMQ+HlQhqXw77eA28HjQ0DITIeH2vykj/z7sVm8uvXJLPLQ/DTiKTPVAGvBjFrmzubaWGSbmVvhTVuJWF0+NWDjOqdUMXvmfRa0d9uh0/m0qJaeMYbP41IFuYFMgcz2T8q/tYL9qSzpXXwgh/2Ah4oOXpfSumWLVv279/v6Oi4ePHi7t27Ny2TlJT0xRdfVFdXT5o0aeLEifqL7733Xm5urv7f/v7+r7/+OgAkJyd/+eWX9T84f/78Zt8QCZ+1Cpb3ZF47o4sd/ddnb3Wy7o1gBo9cQjLRQQUv+jGfXWf/0/fPUVANC2+eZ794QoGnDvKCq7uPjRs3rlq1avbs2f7+/gMHDiwoKGhU4M6dOyNGjOjXr9+0adMWLFjw66+/6q///PPPtra2YWFhYWFh/v7++osZGRmxsbFh/2NjY8NRtZERzOnK5FXD4aw/b4evFdErhXRKF7wRRjLyWnfm29tsYe2fL7+8yfp3gOEeGAb5wdUT4bp169asWfPkk08++eSTp06d+u677/TPdvU2btw4fvz4V155BQDUavW6devGjRun/09PP/10v379Gr2hq6vrrFmzOKotMiYlAx+HM4vP6IY+q1QysDqZXdRdYYbrA5CceFiScd7Mxpvs26FMcS2suvq3MRJkZJzchpeUlNy5cycqKkr/sl+/fufPn29U5sKFC/UFoqKiLly4UP+fNmzYMGfOnI0bN9bV1dVfzMzMnDVr1rJly5q+FRKdMR0ZD0vYmsZmVtLDmews3FkYyc/SEObzG7pKLXx4WfesDxNkh4+DvOHkHkStVhNC7O3t9S8dHBzUanWjMnl5eQ0LVFdXl5SU2Nrajh071sfHh1L6zTff7Nq1Ky4uTqFQODs7v/DCC76+vikpKYMHD/7+++/rHx8bqaur27BhQ/1Aq97GjRvd3d2bFq6oqHjcP1V+DNVoK4LI+IOlezS3/xkaztTWlNc++kfECz9p7SD5RnNnINwG/rnx+EmzsPNjTcvLW3tCy0NIvtHawcLCQqF4xIgTJ4HQ0tKSUlpbW2thYQEA1dXVVlZWTStXU/Pnga3V1dWEEH3hjz76SH9x8uTJ3t7ep0+fHjBgQEREREREhP66o6PjqlWrWgqESqVy+PDhY8eObXjRy8tL/+ZNWVtbt/OPlDGDNFooU1n+yagj7iG+uaes+7/7+G8ocPhJawfJN1rR5ilnqhwcc1d0nplkqPeUfKNxgZNA6OrqamJicv/+/aCgIABIT0/38vJqVMbLyys9PV3/7/T0dP2PNCzQoUMHNze3po+SXbt2zc/Pb+lXMwzj5+c3dOjQx/8rEKe0Wq2litHauetqKvmuC0L8UGqqFA4h5rmJfFdE7jiZm1GpVGPHjt2yZQsAlJSU7NmzR58dUVJSsmnTJv3M38SJE3fs2FFdXQ0AW7du1RcoKysrLi7Wv0lsbOy9e/d69+4NAH/88QelFABqa2u3bNnSt29fLqqNjMnGxubk3h3fTQld99F7fNcFIX7s37F165Nupw//+uiiiEtcrVNauXLl8OHDT58+nZWVNXr06EGDBgFAbm7unDlz/vGPf5iYmEycOHH37t1BQUG2trYajSYmJgYA0tPT+/Xr5+vrSyn9448/vvzyy86dOwPAihUrjh496uXlde/ePX9//59//pmjaiNjCgwMDAwM5LsWCPHGzs5u2rR/8l0LxFkg9PPzS0tLu3Hjhq2trY+Pj/5iQEBAUVGRPgtQpVLt3bv39u3bNTU1QUFBDMMAQHBwcG5u7t27dxUKha+vr7m5uf4Ht23blp6eXlhY6OLi4unpyVGdEUIIyRCHmSsqlapnz54NrzAMY2dn1/CKn59fo5+ysrIKCQlp+m4+Pj71ARUhhBAyFFnnbx06dKh+5SpqpX379rHsI04WRQ1ptdrffvuN71qITHV19eHDh/muhciUlpbGxcXxXQtRknUgXLBgQUlJCd+1EJkZM2ZoNBq+ayEmtbW1M2fO5LsWIlNYWLho0SK+ayEyGRkZb7/9Nt+1ECVZB0KEEEIIAyFCCCFZk9o2r05OTps3b96/f39rCms0mkmTJqlUeA5mG5iZmY0ePVq/yhe1BsuyKpVq2LBhfFdETOrq6mpqarDR2qSqqio/Px8brZGvvvrK19f34WWIPlFdMjQazcmTJ/muBUIIIUGIiIh45LZzUguECCGEUJvgABdC/9/evYY09cZxAJ8ptpzsag2zTU27Z6hsJ3TmrYaUZmFYgkRqLyzEpoYWmphpoGgXS9R8EwqVealXYfiqvOUFY84L5jS11Gpmc9PNy6bP/8UBiU37v9s57vw+r55n/V58eZj7xtlhBwBAaVCEAAAAKA2KEAAAAKVBEQIAAKA0qhRheXm5RCIJCgqqr6/fcODTp0+nT5/GMOzOnTv4g6JAfX19UFCQv79/WVmZ+b8ODQ3JZLKAgIDAwMD8/Hz4sToajba8vJyZmYlh2JkzZ7q6uv4xmZ2dffXqVYsFI7OlpaWMjAyxWBwREdHT07PhzNjYWEJCgkgkkkql+JNqKE6n06WlpYnF4sjISIVCYT6wtLSUk5MTEBAQEBBw7949+Ez7H4gCamtr9+zZ09bW1tjYyOFwOjo6TAZmZmbYbHZlZWVvb6+/v39WVhYhOUmls7OTzWa/e/eura1NIBDU1NSYDDx//vzu3bstLS0fPnzw9va+fv06ITlJ5datWxKJpLe3t6KigsPhzM7Objj29u1bT09PHo9n4XjklJKSEhwcrFAonj59yuPxNBqNycDPnz93796dmZnZ3d3d1NTU3NxMSE5SSUxMDAsLUygUxcXFfD5fr9ebDNy+fVssFsvl8s+fP/v4+OTk5BARc8ugRBEGBgaWlpbi6/T09CtXrpgMPHjwQCqV4uvW1lYnJyeDwWDJhCQUFxeXlpaGr8vKyk6cOPGP4Tdv3ri7u1skF3ktLy9zudz29nZ8GxoaWlJSYj6mVqu9vLzq6uqgCBFCer2eyWT29PTgW4lEUlFRYTKTlpZ28eJFi0cjL41G4+DgMDg4iG99fX2rqqpMZqRS6frbr6ioKCIiwqIRtxpKXBrt6+sTi8X4GsMw8ysJCoUCwzB8LRaLZ2dnp6enLRqRfBQKxfqhHT9+fMPLL+vkcvn+/fstkou8pqam1Gq1SCTCtxiG9fb2mo/JZLLU1FQ+n2/ZdCQ1MTGh1+t9fHzw7YaH1tnZiWGYTCY7d+7cw4cPjUajxWOSC/7E1kOHDuHbDT/ToqKiXr9+PTw8PDQ0VFtbe+HCBYvH3EqsvwhXVlbUajWbzca3HA7n169fJjMqlWp9wN7ensFgmM9QDX65GF9zOByNRrPZt4Dd3d2PHz8uLCy0YDoyUqlUDAZj/Rf7OByOSqUymWlsbPz+/XtcXJylw5GVSqVisVg2Njb4dsND+/btW1FREYZhqampr169unnzpsVjksvfn1e0TQ4tNjaWwWD4+fn5+/s7OTlFR0dbNuMWY/1FaG9v7+DgoNfr8e3CwsLf7yEci8XS6XT4em1tTa/Xm89QDZPJ/PvQ6HT69u3bzcf6+/sjIyOrqqo2fJwypbBYrMXFxfWHNZq/03Q6XXJycmFh4dzc3Pz8PEJIrVavrq4SEZYs/v7To23y58lkMqOjo2NjY4ODg4uLi6urqy2bkXRYLNb63yZtk0O7fPnyvn37VCrVzMyMi4tLQkKCZTNuMdZfhDQazdXVValU4mulUmn+pHs3N7eRkRF8PTo6amdn5+LiYsmEJOTm5mZyaOv/bV83NDQUFhb26NGj8+fPWzwg6bi4uNjY2IyNjeHbkZERk3eaSqX68+dPWFiYh4dHTEyMWq328PAYHR0lICtpCAQCg8EwOTmJb80PjUajubu7c7lcfM3j8XQ6HcWvjgqFwrm5ud+/f+PbDQ/t48ePly5dsrW1tbW1jYmJaW5utnTKrYXoLyktIT8/PyQkxGAwzM/PHz58+MWLFwih1dXV7OzsqakphFBfXx+Lxfr69StCKDk5Gb6ZRwi9fPny4MGDWq3WaDSePHkyLy8Pf724uLi/vx8hNDw8LBAIqqurCY1JLlFRUTKZDCE0OjrKZDLx2xkmJyezs7NNJpubm+FmGVx4eHhGRgZC6MuXL46OjiMjIwih8fHx3NxcfKCurs7Ly0ur1SKE0tPTQ0JCCExLEqGhofiNoH19fQwGY3JyEiGkVCrv37+PD/j5+aWnp6+tra2trclkMji0f6NEES4sLISHh+/atYvL5cbHxxuNRoSQwWCg0+nrt6sVFBSw2WyBQODj4zMxMUFoXlIwGo0JCQlcLpfP54eHhy8sLOCve3p6NjQ0IITy8vI4f3F2diY0LymMjY0dO3ZMKBSy2eyioiL8xa6uLjqdbjLZ3t7u4eFh8YBkpFQqjxw54urqymaz1290bG1tZTKZ+Hp1dTUpKWnnzp179+4ViURKpZK4sGQxODh44MABNzc3Dofz7Nkz/MWmpiY+n4+v5XL50aNHXV1dBQKBt7f3wMAAcWG3AAo9fWJ2dtbOzo7FYm02sLi4ODc35+zsbMlUJKfVag0GA4/HIzrIVjI9Pc3lcul0OtFBtpL/PbSFhYWVlZX1a6QAIfTjxw8ej7fhl/e4ubk5Gxubf3zoARyFihAAAAAwR4mbZQAAAIDNQBECAACgNChCAAAAlAZFCAAAgNKgCAEAAFAaFCEAAABKgyIEAABAaVCEAAAAKA2KEAAAAKVBEQIAAKA0KEIArEdiYqJQKBwfH8e3KysrgYGBGIYtLi4SmgsAUoPfGgXAeszPz4tEIjab3dLSYm9vn5qaWl5e3t7e7uvrS3Q0AMgLihAAq9LT0yORSG7cuBEUFHT27NmSkpLk5GSiQwFAalCEAFibJ0+epKSkODo6njp1qqGhwcbGhuhEAJAaFCEA1kaj0QiFQq1Wq1AovLy8iI4DANnBzTIAWJtr165t27ZNIBAkJSUZjUai4wBAdlCEAFiVysrKmpqasrKy+vr6jo6O3NxcohMBQHZwaRQA6zEwMIBhWHx8fGlpKY1GKygoyMrKev/+vVQqJToaAOQFRQiAldDpdBiG2dradnZ27tixg0ajIYQiIyO7urrkcrmzszPRAQEgKShCAAAAlAbfEQIAAKA0KEIAAACUBkUIAACA0qAIAQAAUBoUIQAAAEqDIgQAAEBpUIQAAAAo7T810QU8ucl7ZwAAAABJRU5ErkJggg==",
"text/html": [
"\n",
"\n"
],
"image/svg+xml": [
"\n",
"\n"
]
},
"metadata": {},
"execution_count": 6
}
],
"cell_type": "code",
"source": [
"rvecs = collect(r_vectors(basis))[:, 1, 1] # slice along the x axis\n",
"x = [r[1] for r in rvecs] # only keep the x coordinate\n",
"plot(x, scfres.Ο[1, :, 1, 1], label=\"\", xlabel=\"x\", ylabel=\"Ο\", marker=2)"
],
"metadata": {},
"execution_count": 6
},
{
"cell_type": "markdown",
"source": [
"We can also perform various postprocessing steps:\n",
"for instance compute a band structure"
],
"metadata": {}
},
{
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Computing bands along kpath:\n",
" Ξ -> X -> W -> K -> Ξ -> L -> U -> W -> L -> K and U -> X\n",
"\rDiagonalising Hamiltonian kblocks: 6%|β | ETA: 0:00:02\u001b[K\rDiagonalising Hamiltonian kblocks: 13%|βββ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 24%|ββββ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 35%|ββββββ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 44%|ββββββββ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 54%|βββββββββ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 63%|βββββββββββ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 74%|ββββββββββββ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 85%|ββββββββββββββ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 94%|ββββββββββββββββ| ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 100%|ββββββββββββββββ| Time: 0:00:01\u001b[K\n"
]
},
{
"output_type": "execute_result",
"data": {
"text/plain": "Plot{Plots.GRBackend() n=210}",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOydZ0ATWReG34QmIBZQAUHBBhYEu2JBcUUQRUXFLhaUoutiXeywlrX3gii2FRsWsKEIlhXBXgDFVURFUEFBkN6S+/2Aj2YgbWaSQJ5fzOTOOSdhMif33lNYhBDIkSNHjhw5tRW2pA2QI0eOHDlyJIncEcqRI0eOnFqN3BHKkSNHjpxajdwRypEjR46cWo3cEcqRI0eOnFqN3BHKkSNHjpxajSQd4cGDBzds2CBBA+TIkSNHjhyWpPII79y5M3bsWFVV1fj4eIkYIEeOHDly5EBSM8KcnBx3d3cvLy+JaJcjR44cOXJKkYwjXLZs2bRp01q2bCkR7XLkyJEjR04pisyrfPDgwcOHD7du3RoSElL9yJSUFHd39/z8/PInhw8fPmHCBAF1/fyJnByWri4B8OEDKyMDZmaEw8GrVyxTU6muLcfhcBQUFJjXm5KCFy/YRkYcNTVWo0bCXRsTkzhqVKMuXb5t29ZER6fOrwPi4rhaWuwGDagxVRAKCvD2LcvEpOR/nZeH9+9Z7dtT/6/PycHnz6w2bQSSXLwfwWKx+I7MykJSEqt1a6m+V8vD5XLZbOF+XqenIy2N1aJFyXv8+hUAdHWptCoxkaWszGnShIdh795luLvnfP2qdvVqrq6utlBiv35lARxdXfaTJ6xu3ci9e6xOnUjduhQZLRhpacjMZDVvLugdItqD5dUrVps2RFm55DAqimViQoT8P1dHfj5iY7kmJrRMzBQUFPjek0w7wry8PGdnZz8/P0H+GRcvXvz27dusWbPKnzQ2NhbwH0kIJk1i9e6NlSsJgOHD2YWFiI3l5uXB0pKdmMhl+JYVipycHA0NDeb17t/P2ruXNWIEt3NnxdmzhXv+Hjt2ND6ek5AQGhx8QUtLw8yMdO6Mzp3RqRPR10dKSsrw4Vdmzeq7cGFrmoz/lS9fYGPD/vCBW6cOALx/Dzs7dlwcl/LfGFFRWLiQ/eABV5DBhYWFXC5XRUWF78i7d1ne3ggKkg1HSAjJzc2tK+T36sIF1t27OH685D36+LBYLKxeTeVbXrDgflbW9YAADzU1tffv8eIF6/lzPH+OFy9Y2dlxubmeLFbzoKAezs7ThRK7ezdLXb1o5UoFJyf2vn3cc+dYf/3FunKFq6ZGoe18uHqVFRICPz9BPy4RHiwcDqys2M+ecfX0ACAtDb/9xv76laukJKyxVfLhQ5aFxdbgYOtevXpRJhQAwOVyORwO/x9nhFkiIyNVVFRatmzZsmVLHR0dRUXFli1bfvnyhedgX1/f6dOni6xryxbSowcpKCCEkJQUUqcOWbCg5CVLS3LlisiCmSAjI0Miei0sSN++3G3b8l1chL42Ozvb29s3PDycyyWxscTfnyxZQqytSZMmpHFjoqzcEVimoGBAvdHVMngwOXGi7LBLF3LjBvVa8vJI3bokO1ugwQUFBXl5eYKM9PIiy5aJZRiTcLnczMxMYa8aPpycPFl22LUrCQuj0qqkpCQWSw9wU1f/vUED0qwZsbMjq1aRgADy8SMhhFy9GnTs2InCwkJhJVtbk7NncwghBw6QYcMIl0ucncmgQSQnh0r7qycmhrRuLcR4ER4sd++SLl3KDk+fJsOHCyuDD+3aWQCubHbT1NRUaiVzOJyCYh9QLUzPCDt06PC1eO0DuHnz5vz58588eVK/fn3KFUVHY/Nm3L+P4p8toaFQU8OQISWvWlsjOBhDh1KuVuaJiYGLC0xMuGfPCn2tmpqaq6tT8d+tW6N1azg4lLz05Qvat2cVFOSzWEXUGSsQbm7YuhUTJ5YcTp+OI0dgZUWxFhUVmJjg8WP070+l2OfPMXkylQKljfx8/PsvDh8uOUxJQVwcevakUgWLxWKxQEiRtnbyw4f4dcHf1nYIr+v4ExWFDh04ABwd4eWFmBh4e2PqVIwciYsXUYfH5gD1tG2L1FSkpPB4X1Rx9WqFR2VwMKytKVbBYrEBBUCZ/1B6YDpYRkFBoeH/qVu3LpvNbtiwobCbCnzJycHYsdi2DS1alJy5ehU5OejXr+TQxgbXr1OrsybA5SI1FWPGcDt2JC9fgsLMmqZN8f79bQ8PxejoUMqECoadHRIS8Px5yeHEibh2DWlp1Cvq3RsRERTLfPYMXbpQLFOquHULZmbQ0io5DA7GwIGgcM0NQJMmTcLDz61cqRMdfZxCb5Gaitxc6OkRACoqmD0b27aBzcbRo9DSgr09KsY20AWLhS5d8PgxjSquXClzhITQ4ggfP77avXsXd/ermpqaFIsWDEkm1FtYWISFhdEhed489OxZNgkAcP06evZE6aaMqSlychAXR4dyGebmTSgooGNH1KtHGjbEhw9UCtfU1Fy+fHnbtm2pFCoACgpwcsLBg6VmwNoaJ09Sr8jcHPfvUykwJQUZGWU/5mokDMw2APTq1Wvx4sVqlO7dRUbC1BSlAU+zZ+PiRXz5AgUFHD+O+vUxYQIKCylUWCU9e+LRI7qEf/qE79/RvXvJYXQ01NXRqhXFWtTU1BwcJhUUtKdYrsBI0hGqqak1b96ccrEBAbh1C7t2lZ2JjUV2NkaPLjvDYsHaWj4prExAQFm0nqkpoqMlag11zJqFM2eQkVFyWLw6SjnFM0IKp9HF00EBYktlmKCgCrON0FBaHCEdREejY8eyw4YNMXlyyWOn2BdyOJg4EUX0bwX06EGjI7x8Gba2KF2zu34dNja0KGrVihsbS4tkQahptUY/f4abG/z8UK9e2cmQELBYlf9/xduEcsoTEQEzs5K/O3ZEVJREraEOHR389hv8/EoOrazw/TsiIynW0rQpNDRA4Ze5xq+LvnoFDgcdOpQcPn+OBg1gYCBRmwSmkiMEsHAhfH3x8ycAKCnh7Fnk5mLiRM7jx08L6Zwb9uyJhw+p/AVWHmam7JA7QgrhcuHoCHd3VArBPX8eampo06bCSSsr3L3L0Dq+rBAXV3aXd+xYc2aEANzcsHdvycOCzcbUqTh6lHot1G4TPn+Ozp0pkyaFXLmC4cPLDul7yNJBVBRMTSucadYM1tY4dKjkUFkZZ88iJGSshcWBfv1G0meJjg7U1PD+PfWSc3MRHl4WWZadjSdPKA4HK6V5c+7XrxJ7INcoR7h+PYqK8OefFU5yOHjwAMOGVR7csCE6dMC9e4xZJ+1kZyM7uyzO09S05swIAVhagsUq+3dPnw4/P+q/dYaGL729d5TGRYvJs2c13BEyNtugHC4Xr1+XzWVL8fDA9u0oKCg5VFWFru6P/Pzunz+n02oPTduEISHo3h2lQf23b6N7d9CUfq2oiGbNKI5LEJya4wifPMGuXTh+HJVypR89goIC7O15XCLfJizP+fOoUwfa/6+tYWyMhATk5EjUJkpxdoa3d8nfLVqgY0dcvkyxigMHJj96VN/e3ll8URkZ+PoVxsbiS5JS0tIQGVk2vcjMxPPnsLCQqE0C8+4dmjSpsP9SjKkp2rbFmTNlZ0JC/EaNKurUiYborHL06IGHD6kXy/AvldatqdxZEIqa4AiTk5M7dbLp23fw6tVffw2+CQpCfj4GDOBxoTyJojxXr1bYnlFUhJERYmIkZxDVTJ2K4GAkJ5cc0hEyo63dgM0ObdCgmZhysrOzra2nKiiM/v49iRLDpJDr1zFgAFRVSw5v3oS5OZisySIO0dGV10VLWbwYmzeX7djp6ekdOeJ6/75BYiKN9tAxIySkQigT6IyUKaZNG7x7R6P8aqgJjvDOnTuvXvUpLLSsW/fWr68GBKB9e/AsKtStG759w6dPtFsoEzx9ih49KpypYduE9etj9Oiy3O0xY/DgARISqFTx/HnInDnLDA33iSnn7t27z583zsy0PHcugBLDpBDZXRcFr0iZUgYPhqIibtwoO6OhgQkTcOAAjfZ064bo6LIlWUp48QIqKjAyKjn88AHZ2TAxoVJFJdq0kc8IxaB798EsVnTPns9sfvm5kp2N2NgKiRPlYbMxaFCFW7Y2k5BQIXIBNc4RApgzB/v3g8MBAFVVjBlTFkpKCUpKSsuWdfD3L0vVEI2ePXuxWDEGBhdsbWXHOQgDh4PgYNjalp25cUOWHGFUVJWOEMDChdi8ucKZP/7AgQPIy6PLHjU1tGxJ8ab+1asVHgjXrsHGht5kHvnSqFj880/DadP8IyLOapUWqPg/t29DQaHy87088iSKYhITUVhY4cGEGhcvA8DMDE2bIiio5LB4dZTauHMdHVha4tQpsYQ8fdrQyCjow4dbNbVV2f37aNYM+volh2/foqAA7SWWTi00v4aMlmfcOLx7h6dPy860aYNOnSBC2ULBoXx1lPkpu3xpVHRycrB/PxYt4v3qhQtQVKzut9uQIQgNZagAhDRz8iTq1SvbsCnG1JT6ZDuJ4+ZWFjLTqxeUlKiPHHZxKVMhGj4+cHOjyBqppNJD9vr1sjrA0k92NpKS0LrqHiqKinB3x5YtFU7OnYsdO2i0itq0+u/f8d9/ZTUpCwrw778YNIgy+TwxNERSEo3z5mqQeUd48CAsLMoWsitx7RoGDKhuOt+4MVq2pCXgSrYIDeURoKirCxYLSTUrXGPcODx7VrYCQ0fIjJUVcnJEL/+YlIRbtyoUCKx5yPoGYbt2laPTK+HsjFu3KhRxtLVFVhYePKDLKmoDR4OCMGgQShsQhoejXTv8suJGMQoKaN5cMhkUsu0ICwuxfTsWLuT96pcv+PGD/wPFxka+OoroaPTty+N8TaovU4yKChwd4etbcujoiIAATnJyNoUqWCw4OcHHR8TLDx7EuHE8QvNrDMeOXf/0KbQ0Mis/H+HhsLSUqE3CUE2kTCnq6nByqlDokcWCqyt276bLqg4dkJj4Njb2m/iiEhMT9+79p3//1NIzjP1SkVS8jGw7wpMn0aZNlU1brl0DIfwb7sizCQnBt2+8Q4pqXrwMAFdXHD1asgKjoJBaUNDLyGjg9eshFKpwcsL586L0uOBycegQKvairlFcvnzZ1dU3N3fn7ds3i8+EhcHEBBLqOiAKgjhCAO7u8PNDSkrZmZkzcf06Pn+mxaoLFy7k58/p0cPm48ePYorq12/0kydphw87lp6RO0LphRBs2QIPjyoHnDqFZs34T+d790ZcHL5R8ENKVineWjA35/FSjXSELVuiSxecOwcAnz59UlBompk5Jjz8Kb/rhKBRIwweLEqPi6tXoatbk0uM/vxZr6AgVU0trbSRvWyti4JfpEwp2tqwt8f+/WVnNDQwblxZIxRq+fLlG9AmK0vzx4+fYopKS6ujrPxBU7PkH5SUhE+fKudW0UTr1pKJl5FhR3j5MhQV8dtvvF8lBA8ewM6OvxxFRQwcWKuTKM6eRZMmvHdSa17gaDGlITOdO3des2ZokyYZxsau1KpwccE+4fMJfXzg4kKtIVJEdjY2buzv6bnv0aMjPf+/kkN3mjblvHwp0IwQwOLF2LsXubllZ+bOxf79tFTUnD3b6eDBvt27Lw0IMOM/umr27EGzZtcuXRp17VpJalFwMAYN4rMnShUSSyUUsvE9o/j6+k6fPr2qV/v0If7+VV4bFUWUlcmDBwIpOnCATJ4svH10kpGRwZiuHj2IlVXZYVFRUXZ2dvHfOTlETY0UFlKjiMk3VT1FRcTQkDx9WnJ4+TLp1IlwuVSq4HKJkREJD+fxUkFBQV5e3q/n4+OJlhb5/2cve3C53MzMzGoGODpW/qIlJpJGjQiHQ69hhLp7LyGB6OgIIdbG5tuMGSe+fv1aembQIOLnR4ktPPj2jejrk9DQKgdUb/CDB0Rbm7x7V+HkhAnk0CGK7KuaYsPevyfNm1MplsPhFBQU8B0mqzPCsDAkJWHUqCoHnD8PFgvdugkkbcgQBAeDy6XKOhnjzRsMHMj7JVVV6Ovj7VtmDaIfBQXMnFlW7GPYMCgrIzCQShUsFpydhQuZOXAAU6bITJkxYdm5E5GRlT+Q4GBYWZW1u5N+qimuxpPY2AlHjiQPGDC29MzcuTSGzDRujOPHMWNGhb1JAfn2DQ4O8Pau0HeXy0VoKP9IC6po3hzfvlWYQzOD7NyAFdm4ER4e1c3Wz59H9+6CTuf19dGkCZ49o8o6WaKgABkZGDu2ygE1dXXU2RnnzpV0jwOwciVWr6Y4uX7aNFy+jB8/BBpcVIQjRzBzJpUGSA/372P9ely4UNnNy+IGoYDrosU0bqyqqPglM1Ol9MywYfj+ncZWugMGYOJEODoKdzNzOJg8GdOmVe5P8OQJdHTQTNzquYKioABDQ1paSlWPTDrC6Gg8fYopU6ocUFCAN28wfrwQMm1scO2a+KbJHleuQEkJ1RQwqZHxMgAaN4aVFY4fLzmkY1KopQVbW/zzj0CDAwJgZMSjs08NIDkZ48bhyJHKtxmHg5s3MXiwhMwSCQFDRku5dcv/2LEhBQUXS2NA2GzMno09e+iwroQ1a5CeLpyKpUvB5cLTs/J55n+pSCReRiYd4caNmD8fdepUOSAsDITwbr1UFbW21tqlS2WVrnhSUx0hADc37NtX9sOZjkmhiwt8fASSWVPDZIqKMG4cnJx41I559AjNmkFXVxJmiYqAIaOlqKqqTpgwcMkStVmzym4DJydcuUJjqQpFRfj5Ye1aPH8u0PiLF+Hvj9OneSyhMe8IJRIvI3uOMCEBwcF8HhnHj0NLC02bCiHWwgIvXwq6ilWTePQIXbtWN6CmLo0CsLCAoiLu3i05pGNS2K8fFBQQFsZnWFwcoqOF++kmK3h4QFkZK1bweEnm1kULCxEXh3bthL5w3jxkZeHYsZLDBg0wdixdeRTFtGyJPXswbhwyM/mMjI0t2SZo1KjySxkZePmyrNAaM8gdoUBs2oRZs8qaJvMkJKTKtIqqUFFB3764eVMc02SSjx8rFLv6lZYtkZaGdHo7bEuMSnVB6ZgUzprFP2TG2xszZkBFhc8w0Siotj1PTk6Ovn43NTWj6zTUlQgMxIULOHWK9269zDnC169haFjdWlRVKCjg0CEsWVLWDvOPP7B/P8WNkyrh4IC+fTF/fnVjsrNhb49163jHFd64gb596botq0IiPShkzBGmpuLUKcydW92YtDQkJ2PaNKGF18LV0e/fkZeHkSOrG8NioX17vHzJlE3M4uiI69cj/fyucTgc0DMpnDoV165VV7EhPx/Hj8PJiUqlxXz48EFFxVBVtcOYMScXL8bUqRg2DL17o0MHGBigUSPUrYt69aI/f26Xm7vewcGf2vrjb97AxQXnzvEuapGWhpgY9O5NpUa6EXaDsDymppg6tcwttW8PY2NcuECVabzZvRvh4dUVdnBzQ5cuVYZoSeSXikR6UMiYI9yxA2PH8tlUuHgRLBYsLIQWXtywntrZgJTj7w91dT7TawCmpjV2mzA5+V1+vquTk/+ePSULVZRPChs0wPDh1YXM+PujS5fquhmIzI0bdwoKenK57hcvvjx5Eg8eIDkZ6uro2BHDh2PePOzahXPn2tevH6WouEJLa1r//tDUhIsLBYWWsrIwahTWr69y4f3GDQwYwPRsQ0zEcYQAvLzw+DEuXSo5pDWPohh1dfj7Y+FC3pWsd+xAVFR1yxUSaRLZvDlSUhjPoKAyd5FqKiXUZ2WRJk1IbCyfq377jZiYiKixVSsSFSXitdTCTO65nR0xM6t8snxCfTG7dxM3NwrUSU9CfSkfPnzQ1OzGYjls2lSWM9yjB7lwgUotERGkVauyhP1KCfV9+pCAACrVFcPhkDZtihQVF+nrm79+/VqQSzIzyYoVRF+fsFhEX59s3ChcqntpQj2XSxwciItLdYOnTyd79wohXEwoufeGDCGXLokl9vZtYmBAii8qKiItWpBHj8S3iw/bt5Pu3UlxWnmpwRERRFubxMVVedWrV8TQkHbbSin/SbZrR6KjqRFbAxPq9+/HwIH8fzg/eiR60EHxpLD28OKFQGtTNa8HRSmGhoaPHp0eP37e27fTS09SPik0N0e9erh9m8dLMTH48IHPNq1odO+OpCSFz583JyREtG3bVpBL6tbFmjVISEBsLAYMgJcXlJVhbo7bt/H48WMfHx+uYFUntm5FfDx27qxuTGiojG0QQviQ0V8ZMAADB2LlSgBQUICrK/bupcS06nB3h64uvLzKzlSV0FIeCTaJlEC8DDVulx7KzwgLCkjz5mU1sari/XvCZlcuESQ4ly+T334T8VpqYWbypKhIgoMrn/x1RvjjB6lXj4IKZFI4IywmI4Po6ZGHD8vOUD4p3LuXjB1b8nf5GeGcOcTTk0pFxfTtS9TUSGKiuHKOHydmZoTN/gp0AEZ37Djz+/cqBxfPCMPDibY2ef++OrFRUaR1a3FtEwrx773UVB7fAhHEpqYSXd2S2ns/fpCGDUlSkpim8ef7d6KvT0JCSEZGRmEhsbAgf/3F55LBg0lgIO2GlVL+k1y4kGzcSI3YmjYjPH4c7dvzr8p/6BDq1q1QIkgoLC3x+DGyskS8XLZ49QocTpXF1crTsCHq1UN8PP02SQgNDaxejXnzaMwpnDIFISH4+rXCyZwcnD6NGTMo01LMiBF49AhPn0JPT1xRkyfjxQv8++8XFosLNI6J0WvcGCoqMDCArS3+/rvyUsGrV6njxuHYMbRoUZ1YmSu0jf8XV6umy7eAaGpi2za4uqKwEA0bYvTosu6Y9NGoEfz8MHUqYmJS/vwTGhq8E1pKyc3FgwcYMIB2w3jCfE69bDhCQrB1a3Udl0q5eBG9eomuSF0dPXrwXsKqeZw5A01NKCoKNLgGZxMWM20aOJyy+DrKw0c1NDB6dFkmWTEnT6JvXzRvTpkWADNmICgI9+5BsNVQgejbt4u//+olSxpmZPyZnY1TpzBkCL59w9at6NQJCgpo3BidO39ms/V69RrJ4Szmu+Ypc4kTEL64WjWMH48WLbBpEwD8/jv270dRETWSq6F/f/z82atXr6Hbt1sfO8anvuu//6JzZ/5hdDTB/NKobDjCgACoq/P/ecLh4M2b6kqvCULtSaK4cwft2ws6uAYHjhbDZmPHDixejIyMkjOUTwpdXHDgQIXa7pRXk1m8GP/8g6AgdO9OpVgAY8aMWb/+bzU1NTU1jBqF/fvx5AlSU8HlIiICrq7Iz/8XaArM//o1okMHODhg1SqcPo0XL0p6IJeSnY3HjyU22xAZMUNGK7F7N3bswOvXMDODhsbZ/v1/f0fzJOjFC2RnfwG2Au/atv05cCBmz8bevbh1q/JCBST9S0XuCHmzZQuWLuU/LCwMHE51LSkEwcYGQUFiSZAVYmKEeBjV4EJrpZibY9AgbNxYckj5pLBbN2hqIiSk5PDFC3z7RmWlzS1bsG0b/vmHuV4BxfTsiTVrEBMzUU+vSFHR88iR2adPw8EBbDYCAuDoCE1NtG6N4cPh4YFNm+IMDMwLC/ukpSUwaqXYiB8pU57mzbFyJVxdkZ2dk5i4PiLC2tn5l0KfFPH0KUaOxNChGDrUjc2eY2LSNDKy/rJlaNcOL19izRqYmUFTE336YNYsbN0KK6sFe/Z0+/hxO0328EVfH2lpyM5mUCU1O5L0UBwsc+sWMTYWKIx70iRiYECB3ubNyZs3FMgRB7rjSoqKCIvFO0b512AZQkhUFGnXTlylUhssU0pSEmncmLx9W3JIeZ9CHx9ib18SLDNzJvn7b8okHzpE2GyyZw9lAkWgqn6EhYXk7VsSEEDWryd9+vzDYnkqKy8+d+4ck7aJee9xOERDg6SlUSmWwyG9exMfH263bjZqar2aNj3w5Ys4NvLg4UMybBjR1ye7dpHcXEKqNjglhdy9S3x8yLx5RFnZDIht3bo3xdZUSyXD2rcnkZEUiBUwWEYGHKG1taBtIZs2JTNnUqB35kyycycFcsSBbp9x/TpRVOT9Ek9HWFBAVFVJTo5YSqXfERJC1q8nI0aUHVIbPpqZSTQ1ycePhcnJeRSGC54/T9hssno1NdJEhm9jXkJIRkbGuHGzJ0/+I0fMm0lIxLz3YmN5J9WJKTY6mjRqRBITSU5Ozo4dpEULyn6CR0SQYcNI8+Zkx44KX1tBDPbzO9ev35jQ0DvUmCIYlQwbMYKcP0+B2BriCEeMWKWvT/Lz+Q/OziZsNnn8mAK9Z88SW1sK5IgD3T7D1ZU0a8b7JZ6OkBBiaso/faV6ZMIR5ucTIyNy7VrJIeWTQldX4ulZtH174bhx1Ai8dYsoKBB3d2qkiYMgjlBSiHnvXbhA7OyoF0sIWbGCjBlT8vexY0Rbm9y/L6Kod+/epaSkhIeTYcOIgQHZsaNkFlgeqf0OVjJs0SKyYQMFYmtI+sSLF/3mz4eyMv+RJ05ASUnQlvTVM3gw7t2TQJdkJomIQKdOwl1S4wNHi1FWxpYtWLAAhYUADTuFbm7w8fm8b98P8cNkMjIygoK+WFtj4kTs2EGFcXKqgNoNwvKsWIHg4FmNGvX097/o6IgDBzBiBEJDhZZz/nxgly7OTZsOmjAhfuRIxMbC3V2U+uBSAsPxMtLuCBMStjZrJlDD3JMnKYvpqlcPpqb8W+fINO/fC53IVRviZYqxs4OhYVnJD2rDR+PjLycl/RYba/v+/WFx5Dx9+rRhQ5OhQ4e3aXNSwN6/ckSG2pDR8nA4OWx2ZGqq9759pwEMH47z5zF5Ms6cEc48T8+vmZkdlJQaXbiQ5uQEJSVarGUMhntQSLsjVFRsHx8vUFTx06cYPZoyvTU7iSIrC9nZQn9ctccRAti2DWvXlrROpXZSGBkZyWKZAmYvXrwQR87jx4+53K7A4DZtLlNjmZyqKc6mpwM1NbU//hirqbmiRYuSzhTF/eD+/BNbt/K/PCEBLi4YNAgTJzrt3dv13LlFXbsKudQjlTDdg4KCVVja8PX17dNnYO6v69y/cPv2fyzW65QUylQ/fkyMjL6Ur4zMMLQu5R85QlRVq3y1qj3CxETSuLFYeqV2f4In8+eTWbNK/t627XGdOgMnTvydI1QVal5wOGz7evYAACAASURBVJzJk2fa208W5MauBi8vwmbvMjcf+b2aomfMUlP3CHNyiJoaKSykWGx5isOVy0fKfPxIjI2Jh0eV+9MpKcTDg2hpEQ8Pkp4uqCKp/Q5WMozLJWpqRHxja8geoZGRQR1+69wXLlwYOHASITNCQk5TpffJk2Pv3k1q0cI8k2+DZxnk6lUYGgp9lZ4euFwKGvTICl5euHoVjx8DwMOHR/PyVl269D4xMVFMsWw2+/DhfadO+fK9savh2zesXYtNm+ZGRAQ0+rWzuBxKefkSbdsKWoNJNLS1sWgRFi0qO2NggIgI3L2LGTMq153JzsbGjWjbFmlpePkSGzZIrAQMfbBYaNkS798zpE7aHaEgPHv2jBBdQOv58+dUyXz16i0hNhkZauk1sTX7s2fo0UOUC2twG4pfqVevrADpH39MNjRcWVTUNC5OX9J2AcDgwWjVCgsXStqO2gGFxdWqYd48vH1boZpHcfmFpCSMHl0Su1dYiAMHYGSEp09x/z58fKCjQ7thkoLJeJma4AgVFFaxWOMtLZVXrVpFlcy1a/+0t2eZmS1r1qwZVTKlh8REDB8uyoW1JHC0lOnTUVSEU6fQu3evDx/uXrp0cPJkdoKkK6IcO4aXL3FNoBgyORRAX6RMeZSVsXs33N2Rn192Ul0dFy9CWTlDW9tWR8eiVauYwEAEBcHfn5ZOzlIFk/EyMu8ICwqwaZPy4sWTbt06r66uTpXY+vXr79u3+PVrWwaK4TJMfDwKC2FrK8q1NTJe5suXL5GRkTxfKi5A6uFR0pDEygrz5mHECOTkMGphefLy4OaG2bP5tHeQQyH05U5UwsoK7dtXzoRRVoab27OCAsPkZEcbm+tBQTAzY8IYicNkvIzMO0IXFygr4++/qZesrQ1Dw5ItoprEqVOoX1/EBKOa5wgTEhJMTYf17bty587jPAeYm2PAAGzYUHK4eDE6daK4UrZQjBwJdXXs2iUxA2ohL18yMSMsZvt2bNmCz58rnOzTx9zOjmNhcd3TcxxDdkgB8qVRQfn2DcePY+tWKCjQIn/w4BqYRHHzJoyMRLzWxAT//cdEyxjGiIrKTUtTz8nRX7LkZ5cu8PBASEjlWgqbN8PHp+w7uW8f3rzBtm3MG4t79xASgnPnJKC61vLlCwDmtuJatoSra+WWcyoqKmfP+vz77zk98TtMyg5yRygo48dDTw8zZ9Ilv0ZmE758iX79RLxWXR26upy3b7n8h8oCjx/D2dnIy2vTiROW6emuBw9CUxObNkFXF1ZW2LgRT5+Cy4WODqZOfde37xRPzy0A6tRBYCC2b5fALt3IkbCzQ//+TOutzbx4UWRqSl0vLgFYtgzh4fj3XyZ1SiNNmyIzEwyF7YubpkEnxUW3q3o1MpKw2eTWLRoNyM8n9euT1FQaVVQFTek+XC5hs0lERHVjqsojJIS8evVKVbWblla3Dx8+iKBdqnKYQkOJjk5ZTdHyZGWRkBDi4UG6diWNGxMHB9K27XggUFGxX+kbDw8nOjokNlYU1cXdJ4S9auZMoqrKo3qk9FDz8ggfPHhYr163+vV7JlVRH52mW9rfn5iY8M5cFBOp+g6Wh6dhHTuSZ8/EEltD8girYcIEdO4MS0saVSgro18/3LxJowqGefAAAHr1EvHyJ0+eFhRYZWb2ErMqisTx88OUKbh4kXedOXV1DBqEDRvw5AmePsWQIUhNzQVWFRWpbdigUVxrrXdveHrCzq6skS+txMbi8GH4+spw9UhZ5N69h5mZowoLjV+/fs2kXgcHNG2KgweZ1CmNMBYvI6uOMCgI//2HU6doV2RtjRs3aNfCGOfOoUkTsFgiXj5mzGgrK+jra9mKFnUqHezZg6VLERwsUDJls2aYPh3Xrq1q27ZxvXrdgoO1pk8vqcft6or+/eHoWKHpPE3Y2KBbN0ycSLsiOeWZOdOxYcMf48a1tbCwYFj19u3w8kJKCsNqpQvGtgll1RE6OWHYMLRpQ7uiwYNx/TrtWhgjLEysQHA1NbW9e/8uKvJSFqQhiPRBCLy8sG8fIiKEiwPs2rXL69ehL1+uJQQvXsDWtmQiuHs30tOxZg1N9pawfTs+fcKlS/RqkfMrqqr18/I279u3lM1m+lHZvj3Gj8fKlQyrlS4YSyWUSUe4axdSUsBMxX0jIygrg9l1ERp5+xa//SaWhBYtkJ6OtDSKDGIQDgeurggKwt27EK1MQrNmCAlBSgrYbPTti8REKCnhzBkcPoyzZ6k29//8+AEPDyxbBm1tulTIqYrXr9GihcSWo9eswaVLNTCDS3DkS6NVwuFg2TLMns1ceT0rqxoSO5qXh4wMjB0rlhAWCyYmsldfJj8fEyYgLg43b0Kc2pxt2iAoCJGR6N0bFhZ4/Rra2ggMxO+/4+VL6swtx9Ch0NPDX3/RIlxO9URGMpRKz5PyRf5qJ/Kl0SqZOxcAo1lcNSaJ4soVKCmJUm67EmZmMuYIs7JgZwclJVy7Bg0NcaWZmiIwEAEBmDgRlpa4exedO2PbNgwfTv2OzrlzePQIV69SLFaOgERFSbiMy/Tp4HBw8qQkbZAgurrIzsbPn7QrkjFHmJkJX19s2EBXBj1PfvsNERE1oWH95csiLglWwtQUVZQkk0aSk2FhAWNjHD9OWbfSXr1w8iQOHoSXFxwccPo0Jk3CqFEYNy7/8eNnXIqCZ4qKMGMGpk9H+/aUyJMjNJGREnaEbDb27sWffzIUnCxtsFho1QpxcbQrkjFHOH48GjXC778zqrS4Yf29e4wqpYNHj9C1KwVyzMxkwxEGBgZ17z66c+crgwdj925QG+7w2284cACrV+PgQXh4wMsLGzfi4cNBFhZ7HBwoqMCWnp5uZ5eioIADB8QXJkdEGKsyWg1du8LaGuvWSdgMScHM6qgsOcL373H9Onx9JaB68OCakETx8SOGDqVAjqkpXr8Gh0OBKFr5/fdVT55sy8vzLK0USi0jRmDDBsydi7NnERiIBQugoZGVn989LOxbXp5Yku/evaul1en6dZtp044xHq4op4SvX8HhoGlTSdsBLFv2c9s2W0PDfjExMZK2hWmYCRyVpS/ZmDEwMRGxbYKYyPo2IZfL7dDBKi+vZ14eBb2L1dWhq8tcGUCR0dOzV1Ud5egoUscpwXB0hIcHJk3CmTOIioK+voe5+fUePXaYmyM+XnSxO3e+4nJ7A/1TUmp9oS3JERWFTp0kbQQAIDHxuZKSQXy848WLNSiXSzCYCRyVGUf477+IjGQig54n3bohKQliNyeXGK9fv379OgPYuGbNdkoESv/q6OfPiItb/vbt0x07PGlVNHs2pk3DmDE4eRIvX/4dETHz/Xs3NzeYm+POHaGlcbkYNQqBga5t22r17/9+925J1PaWA0AKNghLMTc3t7GBouKNoUPHS9oWppEvjVZg8mQMGiSxqAE2GwMHIjRUMtrFp02bNoqKLVisZY6OwygRKP0deletgosL9BnpJ798OYYOxciR6NixuYbG3o8fu8TFYd8+jB8v3Er+588wMMCNG7h9m/X69e47dwIaNGhAm9Vy+CDZ3InyqKioXLjgPWXK2eBgKVioZRa5Iyzj2DEkJUk4hlimV0eVlZU1NU/Pnh2ybh01lSqkfEb433+4fBmLFjGncf16dOkCNnvLggWWz57Nz8vD3LlYuBDbt8PFpaQkW/WcPIkWLaChgaQkMF7PSw4PJJ47UYlZs+DjU+tyCnV0kJdHewaFDDhCQuDujunToaUlSTMGD0ZoKBNVJemAy8X373B0VKdKoJQ7Qg8PLFmChg2Z08hiYe9eREWNX726wdixTjt34sIFnD0LLS3ExWHgQCQnV3e5iwumTMGcOYiJQd26TBktp2ry8xEXh3btJG1HOczNoaFRG9sztWpF+zahDDjCP/9Efj727pWwGXp60NHBs2cSNkM0bt0Ci4Xu3SkTaGCArCwprQh87x6iojBnDtN62WwYGtZTU7v/4YNup054+BDXr2PyZERGghB064anT3lclZwMQ0McP44bN7Cdmg1cORQQE4PWraGiImk7KjJjRm1sScHA6qi0O8KiIqVdu+DpSVkqtDjI7uqovz/09ERvOvErxYXWoqMpE0ghS5ZgzRrJPMKePw+5d29eRsb+Q4cQHY3WrREair170bIlcnIwYAD8/CqMDwqCoSEUFPDpk7g1YOVQi/REypTH0RFBQfj+XdJ2MIvcESIszLl+fSxZImk7AACDB8uqIwwPR7duFMuUzkJrFy4gM1Ni7YpUVFQ6derEYrG6doWPD96/x6BB+PtvPHwIBwc0aQJnZ8ycWbLA7uoKOzuMG4e4uMrlT7OysqysJvTqNfzz588SeSNypCGV/lfq18fw4ZV/TtV4GEgllHZH+PHjcgcHClLfKMHCApGRTBS+o5z37zGc6mw6KSy0xuFg5Ups3EhxERmRadAAzs548QInT4IQpKejVSscPXqpTp2Wdep08/V9dfYsDhzAx4+4dw8nT2LzZri7Y+RIdOkSdvNms0ePrE+elPdekgzSOSMEMGsWDhyoXSEzDKQSSscDozqGhoRskbQNJdSpA3Nz3L4taTuEJDER+fkYPZpisVIYL3P4MLS1eTedlyylE0Q3N7DZ7oA94MBizfn9d2hooH9/LF2KK1fw7RtatcLUqfDx6d2583+amoGHDg1JSJC09bUS6ZwRAujbFwoKNaHio+CcP7/54cNurq5L6VOhSJ9oitirqTmnoABS0gi2eJtw5EhJ2yEMfn6oV4/6WEQTE/z3H4qKoCgdN1FODv76CwEBkrajaurXx+zZYLM93NzWAewJE2Zu2AAdHZ7z1/pPn14CsH07evfGpUvo3Jlxc2sxnz+DxYKOjqTtqIKZM3HwIPr1k7QdTHHx4nku91xgoP3+/etpUiHtM8IxY3o2bTp34EB8+yZpUwAA1tay17A+OJiWKHB1dejr4+1b6iWLxvbt6NuXyshYmnB1dX3+/HJExOl//lnZtCmfVdz587F7N6ytceUKU/bJASIjpaW4Gk+mTMHly/jxQ9J2MMXevWvr1Vvo7LyKPhXS7gg1NBTOn8egQejdG9JQb7Z9e3A4DDVNporoaAwcSItk6VkdTUnBjh1Ys0bSdghGhw4dunTpIuDgkSNx5QpcXCSfQVR7kLZU+kpoaWHo0FoUMjNkyKDJk883amRPnwppd4QAWCx4eWHtWlhaSkWHUtlqWF9YiB8/MGkSLcKlp9DaunUYPx5t2kjaDnro0QP37mHPHri7y2pJB9lCeoqrVcWsWdi/X9JGMEjHjvQma8mAIyxm/HgEBsLZGbt2SdgS2comvHIFiop01WiVkhnhx484fhzLl0vaDjpp0QIREYiMxNixNaFHtJQjtSGjpfTvDwAREZK2gylMTeWO8P+Ym+PePRw4IGjxRpqwssK//yI/X2IGCMWFC2jenC7hUpJBsWIF5s6V3tAGqmjYEMHBqFMHAwfWupRqJsnPR3w82raVtB38qFVVZjp2xKtXNC6HyJIjBNCiBR48wNevGDoU6emSsaFhQ7Rrh/v3JaNdWB48QM+edAlv3hy5uRKOY4qKQmgoFiyQpA2MoaKC48dhbQ1zc7x5I2lraigvX6JNG2kJU6+G6dMRGIi0NEnbwQgaGmjcGHFxdMmXMUcIoG5dBASgSxf06CGxZ4G1tcw0rI+Px6hRdAlnsWhfu+fLn39i5UpoaEjSBiYp3jJftgz9++Pw4Re3bt0itSq5mn6kf120GC0t2NhIuCcPk9D6qJE9RwhAQQEbNuDPPzFggGTS22Wl1lpsLIqKMIyaFoS8kWyhtTt3EBuLWbMkZoCkmDEDXl7Rs2bNtbM78M8/teZZyAjSHylTSnGVmVoCraF5MukIi5k5E/7+mDgR3t5Mq+7VC/Hx0pLaWA3Hj0NTk97y0xLcJiQEixZh/XoZWMWiAwsLhXr1OLm5hR8+KEjalhqFlOdOlMfSErm5ePhQ0nYwgnxGWCX9+iEsDNu3/9DXd+jXb3QKU22BFBTQvz9CQpjRJjqhoTAxoVeFBANHz5wBiwUHB8lolzjt27cPD/fdsmWBj8842UpslXKio2VmRshiwcmptoTM1LQZYVZWlp+f3+zZs52cnHx9fQvFCwBt3RoeHiFfv5rev98tmMH1SplIooiJweDB9KowMcHbt0wH8RJCQkJue3i82rCByt5SMkf79u0XLOizZg3Lzq62BE3QTUIClJXRpImk7RAYJycEBCAjQ9J20E+bNvj6FVlZtAjn7QgfP368YcOGiRMnDh482NbWdtq0ad7e3u8o+tl58+bN48ePd+jQoV+/fjt37pwxY4aYAu3sLE1MIlisf1u3Zq6lW3G8jDSHKeTkICOD9oZEqqpo3pzpqKUDB44OH34wKcmpWTOpqfAmOWbOhK0tRo5EQYGkTZF9ZGiDsJhGjTBwYK0ImVFQQNu2ePWKFuEVHGFhYaG3t3eHDh169OixYsWK8PDw79+/JyQk3LhxY/bs2UZGRpaWlpcuidsXxs7OLjg4eM6cOdOmTTt69OiZM2cKxPsGN2nSJDLy2rJl1729mUslMzBAgwZSkUVXFefPQ1kZhoa0K2J+dbSoiJWfT9TUGFUqzWzeDC0tuLpK2g7ZR1ZCRstTe6rM0JdWX+YIo6OjTUxMVq5caWlpefPmzfT09Pj4+OfPn0dHR3/58iU5Ofn8+fPa2trjxo2ztLRMFyOJj12uzPC3b98aNGigREX7+cWLceMGo09kKU+iuHQJLVsyoYj5Qmv5+VP79nUODz9sZGTEqGJphc3GiROIicF6uqrz1xZkKFKmFCsrZGfj6VNJ20E/9MXLlHXQSUhIcHFxcXNzU1VV/XVckyZN7O3t7e3tk5OTN27cmJ6e3qBBAzF1Z2Zmzp8/39PTk1X1Pk94ePioinlwQ4cOHTduHM/B8+YpLV2q4O+fJ6ZhAmJhobh7t9Ls2bQUvMrOzq7mYxGEhw/VBg7kZGUJXQKHw+Hk5+dzBa7iYGysuH+/UlYW/89B/DcFIDcXW7eqnzvXs3lzbhZNOwY0U1hYyOVyxdwd/5UTJ1gDB6ppa+ePHVtErWQBIYTk5ORIRDVfBLz3nj9XW7gwLytL0JufkltafCZNUt63j7VzJ/8vu5QY/CuCGNa6tUJAgLIgj5pSuFyuigBx8yxJZePm5OTY2toaGxvv37+/qvd/6NCh06dPz6qYJmZsbNyxY0ee4wsK0KED++BB7oABlNvLg+xs6OmxExO5lLf6A5CZmakhXpa4igr70iWutbXQFxY7QjWBVx4TEtCrF/vzZ/7PDvHfFIAtW1hPnuD0aSnenuVHsSMU5PspLDExsLJinz3L7d2bctn8IYRkZ2fXpeP7IDaC3Hu5uWjShP3jB1fwJSpKbmnxSUqCiQn7wwcuX1ukxOBfEcSw79/Rvj37+3chKq1xuVxCCN9Fxwo9Vbdt2zZ06FBjY2PB1YhGXl6evb29gYGBt7d39b8CmjVrNnbsWAHF1qlTXHSDff8+E8GEGhro3h1hYeyhQ6kXzmaz2dW3qquWp0/B5cLKShQZhBChtBsYoLAQ37+ztbX5jBTzTQHIzsb27bhxA2y2NP6qFZDiD0HMj4InJiY4ehQODuywMLRuTbl4Pgh45xBCmJ+UCGLYq1cwNoaKihD/F/FvaUpo2hQDBsDfn823uISABjP/PxLEMG1tKCnh61e2np4QkjkcDn/t5Q927drVtm3bbt26HThwIDMzUwhVwlBQUODg4FCvXr1Dhw5Rfg9NmoT8fFy+TK3UKpHaJIpTp9C4MXO94zt2ZGibcPduWFqiihUBOQBgbY01ayC1CRXv3r1r0cK8RQvzT58+SdqWysjiBmEps2Zh167XX79+FV+Uo+M8HZ0eq1ZtEV8U5dAUkVDBD4WFhW3YsOHHjx8uLi5NmjQZO3ZsaGgo5Wun/v7+V65cefTokbGxcatWrVq1akXJP68YNhtr18LDA0WM7JJIbbzMnTuMtthmptBaVhZ27sTKlbQrknWkOaEiIuJ+QsLgjx/7Dh/+eNs2vH4taYPKIXO5E+VJSzsbEzO/Q4ehHz58EFNUUNCdb9/8NmwImDsXQUGQqm1fJhxhs2bNPDw83r17FxYW5ujoGBQUZGVlZWBgsGTJko8fP1KlcsSIEXFxcbdv3w75P40bN6ZKOIChQ6Gry1D7ZlNTZGRA7BuPet68gY0Nc+qYKbS2cycGDUK7drQrqgFIbUJF27b2ysrp48Zx1q2zjY+HrS10deHoiLNn8fOnhG2TxdyJUlJTfygoGOTkNBB/MU9Pz9PIyPPo0Q2tW2P7djRtCisrbNwoFYGpdAWOkqpJT0/fv3+/ubk5AAUFBRsbG39//2rGU46vr+/06dNFuPDhQ6KvT3JyKLeIB1OmEB8f6sVmZGSIfG1aGmGxyNevIl5eVFSUnZ0t1CVPnhAzM/7DxHlTP3+Sxo3Jf/+JLECKKCgoyMvLo1tLTg7p2ZMsWPAmLCysOGSAbrhcbmZmZvVjRowg27dXOBMXR3x8yLBhRF2ddO1KPD3JkyckNzcvMDAwNjaWKtv43ntcLmnYkHz/TrFYxigsLDx06HTDhneq/47wNfjePWJgQPLzy85kZ5OQEPLHH8TAgBgaEmdn4u9P0tPJq1evLl++XFhYSIX5gn6Sz54RU1MhxHI4nIKCAr7DqnOEpbx+/fr3338v3s8TwgSxEdkREkJGjiRbt1JrDm/8/MioUdSLFecL5u1N1NREVy2CI8zLI2pqFb48PBHnTXl5kWnTRL5aumDGERJCIiJiFRR6qKk5HDhwlAF1fB3ho0dEX5/k5vJ+NTubXLlC5swhrVoRdfUliorzGzUyo8rT8JXz8SPR06NeLMOsXs3na8LX4P79yZEjVb766hXZsoUMGkTq1v2ipNRVRcV11aqNohgqvGHFCPioKUVAR8g/vuvu3bubNm06duwYl8tt1aoVDZNSWtiwARs2MBEvkJZ2JjCwm43NFNo1CcyVK2jThlGNKiowNMR//9El/+dP7NmDpUvpkl9T0dLi1qunkJurcvu2ZDILK7FyJZYvR506vF9VU8PQodizB+/ewcFBmc3OzsoijMVkyvS6aCnu7rhyBbGxIl4eGoqkJEyp+mHWvj0WLkRICKKjFdTUOAUF2VlZFJRDERyaHjVV3mSJiYkbN240MjLq37//qVOnbGxsQkJCYkX+gBnH2Bh2dti6lXZF589f4HL33L//OjeXlsx6EXj+HP36Ma2U1kJr27ZhxAjIy8gIi5GR0Z07+48enfLkyfQlSyRszL17ePMGApYW9vVdFRAwVlv72oMH6jTbVYJMR8qUUq8e5szBpk0iXr5yJby8oCBAXy9DwyaRkYEuLk7R0e4iKhMVWrYJK80Qc3Nz/f39hw0bpqCgAKBr1647duxITU0VdCJKKeIsjRJCEhOJlhZJTKTQIh48fvy0VavRnTrto1asyEsuXC5hs8m9e6KrFmFplBCyfj1ZtIjPGNHeVFoaadyYxMWJcKmUwtjSaClJSaRLFzJnDuFwaNRS/dLogAHkqJALtNeukTZtqlxKFQq+996YMeTUKerFMk96OmncmLx/z/vVagy+dImYmAh3hxQWEjMzQknoiOCf5Nq1xMNDULGiLI0uWLBAW1t77NixT548WbBgwevXr588eeLu7q6pqUm1/2UCPT1Mn4516+jV0q1bl6dPz3386Pb9O72KBCQsDADMzZnWS9+McPNm2NszVDe1pqKtjX//RWwsxoxBHkMlCCtw4wa+fMGkScJdZWMDExNs3kyPTRWpGUujAOrXx6xZQk8KCYGnJ9auhVBL0YqK8PHB/PmMRvzSkUFR4U1fvHixR48e/v7+nz592rRpU9u2bSnWxjjLluHcORr3roqpXx92djhxgl4tAnL6NJo2Fe5upgRTU7x4Qb3Y1FT4+EDiy3o1gLp1cfkylJRgayuB9nWenli9WpQKDzt3YudOvKW53VZ2Nj5/ZnpnnT4WLIC/P+Ljhbjk7Fmw2Rg+XGhdPXvCxgZeXkJfKDK0O8Jnz56FhIQ4ODhQ0g5CGmjYEPPnM/FPmjkTBw7QrkUQwsLQpYsE9BYXPUpKoljs5s0YNw4tWlAstnairIyTJ2FsjIED8e0bc3ovXUJWFhwcRLm2WTMsWYK5c6m2qSLR0WjXjrlKTHSjpYWZM4WIkOBw4OWFdetErEy5aRNOn8bz56JcKwLNmyM7G6mpVMqs4Ajr169f+vfHjx+XLVs2YsQIW1vb4jOnTp26du0alcoZYd48hIfj4UN6tVhYAMCDB/RqEYS4ONjZSUZ1x44Ur46mpODQIXmwKJUoKMDbGw4OMDcHRZ22+VC85rZuneirFPPm4ds3nD1LqVkVkeniajxZvBgnT+LzZ4EGnzgBTU2IUKC/GE1NrFsHFxcI3LFGLFgsmJhQHC/D+9588OCBmZmZr6/v58+fo/4/C01ISJg/fz6VyhlBVRUrVjCxtjZtGnx9addSPUlJyMvD6NGS0U55obWNGzFhAvT1qZQpB4CHBxYuRP/+tKxmV8LfH4qKYv04Y2AjqsZsEJbSqBGmThVoUlhYiNWrsXatWOqmT4e6OnOrYpR36OXtCN3c3Hr27BkXF7dlS1nd1aFDh7558yY5OZlK/YwwcyaSkxEaSq+WqVNx4YIEdl/Kc+IE6tZFw4aS0U5tobWkJBw+DA8PygTKKc/s2dizB0OGlERX0QSHg7/+En3NrZQePTBkCDw9KTLrF2pG7kQlFi/GsWPgW8j5yBG0bAkxW9exWNizBytX8ldHCZRnUPBwhD9+/Hjx4sVff/2loaFRvhOHgYEBgC9fvlCpnxEUFLB6NRYtonfmrq0NS0ucOUOjCr5cuwYJRjhROyPcuBGOjhCq34ocobC3x8mTGD0a27a9vHbtmiDdaoTFzw+amhg8mAJRGzfC3x/PnlEgqhKE4OXLGtjSREcHkyZhxxLtiQAAIABJREFU+/bqxhQUYMMGrF5NgboOHTBjBkO/XClvd8PDERYUFABQV6+cx5qamgpAUTY3lEePhqoqzp2jV8vMmTh0iF4V1RMVBUtLiWlv3x6xsTHJyenii/r6Ff/8g8WLxZckpzosLXHo0MfFi53s7c9u2+ZNrfDCQqxZI+6aWymlG1GU++snT5LV1N5raVEsVhpYuhSHD6OazC5vb3TsiF69qFHn6Yl793DrFjXSqsHUFK9eUTmx4eEImzRpoq2tfenSJQDlZ4QnTpyoW7cuA2176YDFwoYNWLqU3sY01tZISmKiDwNPioqQmoqJEyWjHcDevd6FhUtMTPqnp4vrC9evx4wZaNqUErvkVEenTgoNGnCKinJ9fJTevKFS8uHDaNVK3DW38kybBg0NHDxImUAAsbGxVla2qalOly9fpVKudKCri7FjsWMH71dzc7FlC/76izJ1amrYuxdubsjPp0wmTzQ00KgR3r+nTCAPR8hmsxcuXLhmzRpPT8+3b99yudwXL14sXbrUy8vL3d1dWVmZMuXM0r8/WrfG4cM0qmCzMW0avSqq4fp1sNmS3PP//DkZ6JiTUydHvA5mX7/i1Cn5dJAhmjVr9ujRmdu35yxe7NynD5YsoeYplpeHdeuoWXMrhcXC7t0Ub0Slp6fn5TUGjJKSGEwoYZClS7F/P1JSeLy0axd696a4cemQIWjbFlvo7+lLcTYhz3ozXC53yZIl5bMJ2Wy2s7OzILVqKETMEmu/cvHiG0XFgebmI/k2ixGZT5+IlhYFHaBEKN00bRpp0UJcvUTUEmuEkOzsbE/PQ40aRVR1mwjypqKiogYN8lu4kJEeWhKC+RJrAvLlCxkzhrRpQ27dEuXy8iXWtm8nI0ZQaVspHh5k8mShr6rq3ktLI3XrXtuz52RRUZEIxkhhibVfmTWLrFpV8nepwZmZREeHxMRQry4+XpSaiMJ+kitWEE9P/sPE6j7BYrHWr18fHx9/4sSJLVu2HDly5O3btz4+PrKeaP/y5eWiosnPnzd9SluLyWbN0K0bAgJoEl8d9++jRw8J6C1FTU3Ny2tGu3bmV66IKCE9Pd3Scmpo6KvUVOrWa+QIjK4uzp7F1q2YNg2OjvjxQ0Q52dnYtImuQharViE8HDdvUiPtyBHY29vMmTNBQZBS07LJ8uXYt69yK56tWzF4MC1trps3x6JFmD2besnloTZwtLocV11d3YkTJy5cuHDatGky1ICpGiZNGmNmdoHLzWrenEaP4eQkmZCZDx8wYoQE9FbCzQ3eokZdKCkpZWdzVFU/NW1an/9oOfRgZ4eoKDRsCFNT/POPKBJ274aFBcVrbqVQuBFFCHx84OZGhVlSjIEBhg3D7t1lZ9LTsW8fVq2iS+OCBUhKwvnzdMkH5amEpXPD5ORkAVc+f/z4kZWVJchIMaF8abSYRYuIszPlUsvIzyfa2kTM3trCLhTExREWi4i0olkZkZdGi8nPJzo6vFvJ831TN24QQ8OUsLCHzHRUlxRSuzRaifBw0qEDGTaMxMcLNL54afTnT9K4MS1rbuUZMYKsXSvEeJ733o0bxMxMLDNkYmmUEPLuHdHSIunpJQYvWULvM5AQEhZGmjYl6emCjhf2kywqIurqhO8el9BLozdu3DA2Nj5w4MDPqus3JCQkeHp6tmzZUhazCUtZvhwXL+LlS7rkKytjyhSmQ2b8/NCgAdTUGFXKE2VlzJgBHx+hLywowNy52LNHq2/fHiwxE7DlUEHv3nj+HH37onNnrrn58n79xrwToCzbtm2wtaVlza08u3dj507ExYklxNub9hU8KaFVKwwZgn37AOD7d/j6YsUKejX27Qtra4qjpcqjoIC2bRETQ5G4UpdYUFCwa9cuLS0tVVVVOzu7devWnT17NjQ0NDg4+OTJk0uXLu3Xrx+bzTY2Nr548aJQrltkaJoREkJ27iRDhtAhuIT//iM6OkSc0CJhfx/17Uv69BFdXXnEnBESQuLjiZYWj+lp9W9qzRpiby+OWplBVmaEpVy9GqOkZA+cb9PGY+dOEh7OOxyMy+UmJGQx1jxy40ZibS3o4F/vvc+fiaYmEXNGJyszQkLI27ekSRPy5UvmvHnE3Z0JjSkpREeHPHsm0GDRwgMPHuQzRsAZYeWo0czMTB8fnz59+lRKnK9bt+6wYcMCAwNFi60SDfocYUEBMTIiwcF0yC6hb18SGCj65cLeFg0bEi8v0dWVR3xHSAgZPpwcOlT5ZDVvKj6eNGpUZTfRGobMOcK8vLzevUfo6PT08AibPZt0707U1IiZGXFyIvv3kydPSn7zxcTEWFtfdnVlyKrCQtK69f0+feaEh9/nO/jXe2/VKvLHH+LaIEOOkBAydmyhtfXBBg3SPn9mSOO+fTna2p7Ll28oLCysfqQIn+TWrfz/gyI6wlKysrIeP3589erVGzduREVF8X0bdECfIySEXLhATE0JfW79yBFiZyf65ULdFnl5hMUib9+Krq48lDjCa9dI586VT1bzpoTd8pFpZM4R/kp+Pnn8mHh7kxkziKkpUVcnnTp9ZLH0gbF9+kxkzAwdne7AfX397nxHVrr3CguJnh55+VJcA2TLERoY9AWmsNltGHueHzp0hMVaoqIy6+rVq9WPFOGTDAkhAwbwGSNW+gQAdXX1bt262draWllZdezYUUYrq1WDvT0aNcKxY3TJHzcOERFISKBLfnkCA6GkJF1tRa2tkZWFx48FGhwcjFevsHAhzTbJoQ5lZXTrBldXHDqEyEh8+4ZFi7IBFaBpXh6v5G16GDiwp7q6R15eT2EjSAMCYGSEDh3oMUtaKSzMBVoDOQW0VtgqR9eunRs3Di0qevbjR3vKhVNY3JjxRubSxJYt8PREdjYtwlVV4eAgYvS5sAQEwNCQCUWCw2Jh1iyB8ijy8/HHH9i1C3Xq0G+WHHpQU8OkSe1Pnlzr7Jx/6xadnQMrcuLE7qSkoMGDdzs6ghAhLvT2rvlZE7/y5MnlqVPjr107rMZUWJ2ZmdmnT/du3XqwYIEh5fGJjRtDSUnQnovVU6sdYefOsLCgsRpQcQ1uBppVPnoEc3PatQjLjBkIDOTfSHrDBpiaYsgQRmySQyfjxo3bunVTvXr1mFRat676oUP49Al//y3oJa9f47//MHIknWZJJbq6urt37xhMSTcQgVFRUbGwUNyyBcOHV1f+WzSoKrRWqx0hgI0bsXs3EhNpEd61Kxo0YKIWe2KixJrxVoOWFuzscPRodWPi47FnDzZvZsgkOTWSOnUQEICDB3H6tEDjvb0xaxZkvEyWjOHoiPHjMWoUxfW4qUqrr+2OUF8fM2fSVQsKjFSZiY5GURGsrenVIhrFVWaqmRP//jsWLZK6dV05MoeODgID8ccfePSIz8icHJw+jVmzGDFLTjnWrYOeHlxcqJRJVaG12u4IASxfjqAgWhp+Apg8GdevU78gUJ4TJ9CoEaSzKUivXmjQAKGhvF+9dAlv3mDePGZtklND6dQJR49i9Gg+EWp+fujXD/r6TJkl5/+wWDhyBK9fY9MmymRS1aFX7gihoYHly+lq+lO/PuzscOIELcKLuXULpqY0yhcTV1feITO5uZg/H7t2QUWFcZvk1FBsbfHHHxgxoroIuNpQXFRqUVVFYCD27kVgIDUCO3TAu3cUdJmVO0IAcHFBcjKuX6dF+MyZOHCAFskAUlJSXr26MWgQzX0wxWDiRNy7h/j4yufXr0f37rCxkYRNcmouixeje3dMmcJ7QT4iAhkZGDiQcbPk/B9dXQQGYtYsQXOrqkdFBQYGEL+htNwRAoCiIv7+GwsWoKiIeuEWFgDw4AH1kgkhBgY9cnKCT5yQgq4TVaCmhkmT4Otb4WRcHLy9sXWrhGySU6PZswc/f/KupVlcXJQtf+xJlM6d4eODUaOoiVKkJHBUfkeUMHw49PRw5AgtwqdNq+wJKIHL5ebmKgDqCQmfqJdOHcVp1+WXL9zd4eEBPT3J2SSn5qKkhLNnce4cDh6scD4lBVevYupUCZklpxyjRsHNDSNGICdHXFGUxMuUOcLU1NRXr16JK0+W2bwZq1YhI4N6yVOn4sIF6iXn5iqw2dc0NR+cOrWDYtGU0rYt2rUra1YcEIAPH+DuLlGb5NRoNDVx+TJWrMCdO2UnDx2CvT00NSVmlZzyLFuGDh0gbCWEX7lxY+6OHd19fI6LI6TMEV66dGns2LHFf/fs2fMZTWGUUkynTrCyoiW/XlsblpY4c4Zisa6u0NBo/f37jSFDGM2QFYHSbr25uVi4ELt3y7O45NCLsTH8/DB+PIo7RxECX195mIx04euL1FSxstcIIdHR/+bn/+Pjc0ocS8ocYd26dbP/H2sVHx+fl5cnjlwZZf167NuHTzQsNBZXmaGQlBScOYMNG2Rjw2PkSMTFISaGvXYtzM3l0QpymMDKCuvWwc4O6ekICVHU1ES3bpK2SU45lJVx7hxOnoSfn4gSWCzWypW/16v359Chy8SxpKyUtpmZWUJCgouLi5mZWW5u7sWLF1/yqg3n7Owsjj4ppzjf09OT+s1Ca2u4uSEyEmZm1AicMQOampCV/4aiIkaOfDt27OafP+e/ekV9+V05cnji5ITISFhY+KWkfF2yZAagJWmL5FRASwuXLqF//69Hjy5xcrKZMGGCsBLmz3du3tx51y7x7CjfimL37t3/Y+++42re/ziAv85p74koUklmSyShEJfIyAgJZXftfe11jYt70bVFZGYLJbOIiqzscZOZFWlonfP9/XH8Ujn7fM/3nFOf5x/3ke/5fD+fd91T7/P9TAsLC/HLy5tcj2ESJCeHqlWLunGD/pqnTHnn5bX+8ePH4hQWfijJq1eUmhq1bx9NkZVHyzFMvzIwaACEsVh2tNescirBMUx8cbncnJwcRUdRUXr6SxbLBlhQr14b2itXrWOYKGUNuEYNd2ARi1U3MzNTitt5h2qlpfF5SZpjmMaOHfvu3bvc3Nxq1arFxsZm8SNb2lUB+vqYOxdTp9Jfc1xccHy8Zvv2En/k+VVQEGrVQv/+stfEHBMTQzY7XkNDeZc8EpWSsbEBm13AYiXVqGGs6FgI/qytzVisGxSlV1KiJ8Xt6uoICZFpuTaf8SU9Pb158+Y1adLEhB9emXfv3sXExEjfrHIbMeLHTGt62dhYaGtfzc83krGep09x+XLFqeHK7/nzxK1bu7x/T8fOgAQhNmNj44yM6zt3Drxy5YSiYyH4S04+FRHRe/jw+MBAfel25R45Env3Sn+mHv+JFmPHjrUUusjr5s2bkyZNkrJNpaemhsDAhF69PLy9+xYXF9NV7dGj4RcujK9WLVbG6aMDB8LeXkl32RZCXV29b9++xsbkUznBNEtLy55V8NQlldKrV6/Nm81MTTF6tDS3W1nB01Pc40d+pQozDhXhyZPjxcVzrl8vyPh1czBpsdlsDw/nvXs1x4+XfmLq9eu4eVMuy/MJgiAUiM3Gnj149AhLlkhz+6hR2LxZ2qalvK+ymz59RLNmW9lsx1u37OituVkzjBsncC9EkYYMgYsLPD3pDYogCELxdHRw9CjCwxEp+fr4zp3x6RNSU6VplyRC/ho0aHDjxrGEhD/HjmXRvt/OrFlgsfDPPxLfePEinjyR1z5wBEEQCmdhgVOnMHWqxEeas9kYMQKbNknTKEmEwri4YPVq+Pvj61c6q2WzERmJlStx+7ZkN4aEwNMTTZvSGQxBEIRSadQIUVEIDJT4WInhw3H4ML58kbhFkghFGDQIXbogIAAcDp3V1q6NlSsRGIjv38W95cQJvHqFiAg6wyAIglBCXl5YuhS+vvjwQYK7qlVDp07SnP9KEqFoq1eDxcK8eTRXGxQEJyfMmCFu+dGj0aULbGxoDoMgCEIJBQdjwAB06ybZCRWjRmHDBok38iaJUDQ1NezdiwMH6N81e8MGnDgh1oLFrVvx4YPqrR0kCIKQ2uLFaNAAQ4YgJ0fcFYLt2oHFQmKiZA1JmQjbt29//vx56e5VRaamOHwYY8eC3jM5jI0RGYkRI/D+vbBiXC5mzkS/fhC1/x1BEETlwWJh61ZcuTLdwsK3Z89hYt41fLjEU2b4J8LFixdHR0dXuLhjx47GjRvzvtbR0RG+4r7ycXLC5s3o3RufPtFZbZs2CApCcLCwZ/m//0ZODtavp7NdgiAI5aelBV3dxPz8ZXFxt3JyxLolOBinT0s2uMgnEXI4nOXLl3M4HABZWVmlZ1C0b9/+wYMHr1+/lqD6ysXfH/36oX9/lJTQWe2SJfj0SWC3Z3ExFi5ESAj+v70dQRBEFXLgQFi/fvt8fdc5OiIpSXR5Y2P06IGdOyVogk8i/PTpU35+vq2tLYCYmBh/f3/e9Vq1arFYrMzMTAmqr3SWLYOWFmbOpLNODQ3s3InZs/nPFV60CMXF+PtvOlskCIJQFW5urgcOhB061Prvv9GrF9auFT0XhrfLjPiblggcIywpKQGQk5OT8//HUd4BK4aGhuLWXRmx2di9G8eP07yqvWFDLFyIgQNRVFTuekEBVq/G+PHQ1aWzOYIgCJXTqxeSkhAVhZ498fmzsJItW8LICOLPY+GTCKtVq2ZiYrJ9+/b8/Pz9+/fn5uYmJiYCiIyM1NPTq127tsThVy4mJoiOxsyZuH6dzmrHjIGlJRYuLHdxyhSwWFi8mM6GCIIgVJS1NeLj4eKCZs1w+bKwkhLtMsMnEbLZ7ClTpqxfv97AwODLly9//PFHp06d3N3dp0yZMnLkSB0dHcmDr2waNMDmzejbV7LxWOFYLISHIyICly79uJKbi23bMGcOtLRoa4UgCEKlqatjwQJs24YBA7BggcCtTgYNwqVLePNGvDr5Xp09e7aLi0t6enr//v2NjIxYLNatW7cGDBgwbtw4aYOvbHr2xM2b8PfHhQvQ1KSnzmrVEB6OwYNx5w7U1TFmDHR0MG0aPZUTBEFUGj4+SE3F4MGwtV2mpnb677//6NnTt2wBfX0EBGD7dtasWaJrY1GSLsFnUHh4eGJi4vbt2xUdCH9cLnr2hLU1wsLorDY0FN++YfHi3Pr19deuRWgonZWLg8PhFBYW6sphWDInJ8fAwID2alVOcXExl8vVqnRP+hRF5eXl6evrKzoQPuT03lO5t7TSBixdYCUllJGRc37+SUPDkbNnx3h7w9UV6v9/vrt7F9264cmTYm1tDeH1kJ1lpMc7PSsm5r6jY/Dmzbvoqvbvv3H8+CBb27ZqauFjxtBVK0EQRGWjrs4aMaJX3boB48ePzspCaChMTdGxIxYswLlzsLfn5OSMXrr0juh6GIi1EjMwgLn54uTkkZMmTe3Wra+lJQ0DqNrayM29CDwqLGzMYom7mQJBEEQVtGbNgjVrFpT+8/NnXL6MS5cwdSqePbuSl/cqI8NBZCXkiVBWQ4Z0MTefYWxs6+ioM2GC9EfP81y8iO7dAbCBZurqLHpCJAiCqBrMzNCzJ9aswe3bSEqyU1d/2KhRmsi7SCKU1ZgxQzIzk96+PfTwIUxM4OoKPz+JT0kuKsLBg2jZEqNGoUMHPH+eMm/egHfvpDprmSAIggCaNLF68uR8SIi9yJIkEdJATU0NQPXqWLAA6enw8UGPHmjdGr9s18pHdjbWrkW9eli79sfmMhMmwNa25tSpU83NzeUeOkEQROVlbW1tbGwsspiIRFhYWBgfH5+bm0tTVJWfgQEmTMB//2HkSEyfDjc37NrFf6XL8+eYMAE2NkhNxalTuHIFfn5gkd5QgiAIZolIhO/fv/f29n769Ckz0VQampoYPBj372P+fGzYgAYNsHYtCgrw8OFDDodz5Qr69UPr1tDRQVoadu1C06aKjpggCKKqIrNG5YjNhp8f/Pxw7hz++guzZo0E8ijqg43N2YkTsXMnyC49BEEQCkcSIRN8fODjg4YNMx49GmxmtvrePdIFShAEoSzIZBnmxMVtW7IkMz4+kmRBgiAI5SHiibB27dq5ublko21a1K5de/bsKYqOgiAIgihHRCJksVh6enrMhEIQBEEQzCNdowRBEESVRhIhQRAEUaWRREgQBEFUaSQREgRBEFWaWImwqKjo2bNnJSUl8o6GIAiCIBjGPxFOmzbtn3/+4X197969unXr2tvb16xZ89q1awzGRhAEQRByxycRcrncTZs22dra8v45ZcoUDQ2NHTt2uLq6jhw5kqIoZiMkCIIgCDnis47wy5cvubm5DRo0AJCdnX3x4sV///136NChXl5etra2b968sbKykrHVe/fubdy4MS8vr2/fvl27dpWxNoIgCIKQGp8nQt4zH5vNBnD+/Pni4uJOnToBsLS0BJCZmSljk69fv27Tpk2tWrU6dOgQHBx84sQJGSskCIIgCKnxeSI0MzMzMTGJjo6ePHlyRESEg4ND3bp1Abx58waAiYmJjE1u2bKlY8eOs2fPBpCXl7dq1aru3bvLWCdBEARBSIfPEyGLxZoyZcrUqVNr1KgRHR09btw43vULFy4YGhpaW1vL2GRSUpK3tzfvay8vr+TkZDLuSBAEQSgK/71GZ8+e3bhx45s3b7q6uvbs2ZN38fv373PnzlVXl/XkpszMTDMzM97X5ubmRUVFnz9/Njc351s4MTHR39+/7JWuXbsGBATIGIPyy8vLYynolAoOh1NYWMjlcmmvWYHflFIpLi7mcrnFxcWKDoRmFEXl5+crOgr+5PTeU7m3tNIGLKfAuFyulpaWyGICs1rPnj1LUyDP2LFjaYgL0NbWLioq4n1dWFgIQMjpFnXq1Onfv3/ZKw4ODrq6urREosw4HI6ivk0Oh6OmpiaP1hX4TSkVXiIU5/dTtVAURVGUcv4vltN7T+Xe0kobsJwC43K54vQ4CkyE+fn5R48evX//PofDWbFiBYDr168bGBjwZpPKwsrK6uXLl7yvX758aWJiIuSAi9q1a/fr10/GFlURm83mzVdiHkVRcmpdgd+UUuH9ECrfj0J+7xzZkbc0j9IGLL/AOByO6Nb5Xk1PT3d0dAwKCtqyZcuePXt4F48ePRoUFCR7WL179z5w4ADvWTAyMrJPnz6y10kQBEEQ0uGfCIcPH66jo/PkyZNDhw6VXuzdu3dqauqXL19kbLJv3741a9Z0dXX19vaOjY3lTR8lCIIgCIXg0zWak5Nz6dKl2NjYevXq8ZZM8NjZ2VEU9fr1axlXUGhqasbGxt66dSs/P7958+aVb6SEIAiCUCH8EyGXy+Utny/r+/fvAGjZepvFYrm6uspeD0EQBEHIiE/XaPXq1Q0NDS9fvgyg7HzW06dPa2ho2NvbMxcdsGfPlStXrjDZIkHI1YMHD7S0bLW168XFxSk6FoIgAL5PhOrq6sHBwX/88YeJiYmxsTGA/Pz8I0eOTJs2bdCgQfr6+kzGV1Q0eciQ8c+f32SyUYKgV0kJnj/HvXt4+BBLlw6gqH6Aub//vNzcTooOjSAIAcsnli9fnp6eHhAQwJvPamRkVFJS0qZNm9KzmRh0QkdnJuONEoQ0jh49tXv3iZkzxxgaOt+/jwcPfiS/p09haYkmTdCwITw8gi9cWAUY5+cftbVFeDjatVN03ARRtbGELDaMj4+Pi4v79OmTvr6+t7d3165dGV6AEh4evmBBszdvnA8dQvntZSq/nJwcAwMDhTTN21lGHotbFfhNMYDLhbGxY07OZjZ7rq3tucaN0agReP9t2BDa2gCQmwtHR3z5UqitjYMHtVaswOnTqFsXu3fDw0PR34DMKIrKy8tjuNNITHJ676ncW1ppA5ZTYFwul8PhaGhoCC8mbL80Ly8vLy8vWqOSWMuWe0+dch4wAC9fokYNxcZCEALFxWH6dKirOxsbTxo50n/FCv7FQkNhaoru3TWKiqiEBERHIyMDQ4bA0xONG2PvXjRtymzcBEEIWkeoPAwMPvXsCSMjtGql6FAIgp8HD+Dnh3HjMHs2Pn/e9f59wooV0/mW3LULqanIzsaAAdxevbiHDwOAtTUuXcKdO9DSgpMTPDzw/22XCIJgCP9EaGFhwRKA4fgAzJwJDQ28fYvgYOYbJwiB3rzBqFFo1w6tWyMtDX37gsWCpqYm38LPn2PaNMycCXV1uLpSrVtz373D8+c/Xm3aFDdu4OpVZGfDxgYdO2LOnFVTp86gZbUSQRDC8e8anTNnTl5eXuk/v379mpCQcP/+/UmTJjEV2E+OjnB0hIMDwsLQrRt692Y+BIIoJzcXq1bh338xfDiePIGRkYjyxcUIDMTChUhORmAgALDZ6N4dR45g2rSfxVq2xIMHOH0agwbFnDt3DjDJzZ24adO/cvxOCIIQlAj5HjQxbty4Bw8eyDke/mbMwOjRGDYMAweSwUJCkYqLsWMHFixA69ZITYWYp3POnAlLSwwfDisrXL3642Lv3pg7t1wi5PH1xY4dRb16ZVBU0cOHnnRGTxAEPxKMEU6YMOHgwYNlN11jjLc3zMzw22+wtiaDhYQCPHv2rHPnwQMGrG7cGAcPIiYGUVHiZsHYWBw6hC1bEBuLevVga/vjert2+O8//iOCPXr0uHYtYvbseVeujF22jLbvgiAIviRIhGpqahRFffr0SX7RCDFtGlaswLVrePsWISEKCYGouubMWXvmjP+hQ7GLFr06exZOTuLe+P49hg3Dzp0wM8OePT/6RXnU1eHnhyNH+N/o7u6+ZIn3hg2YMweRkbLGTxCEEGIlQg6H8/Tp0/Hjx2tpadWvX1/eMfHVowfy83H7Ng4dws6dAv98EATtiorw7Jmvjs6fzZpp9+5tIf6NXC6CgjB6NLy9kZeHM2fQt2+5Ar17gzd3VJBRozBrFoYORXS0VKETBCEG/mOEFhYW79+/r3BRW1v7n3/+EXKavFyxWJg6FStWIC4Ow4eDrCwkmMHhYMgQWFl1uXats4aGZLOm//oLhYWYNQsADh9G69YwNy9XoGNHDB6Mt29Rq5bAShYvxufP6NULFy63s2jnAAAgAElEQVSgbVspvgOCIEQQa9aompparVq1vL29awn5fZW/wEAsWIDUVGzejAsX0KrVz9nnBCEPFIXRo/HpE6KjIWkWvH4da9YgJQVqagCwZw+fLn0NDXTpgmPHEBoqrKoNG5CVBR8f3LgBR0eJoiAIQjQJZo0qnIYGJk7EX3/hwAEkJcHKCsOGITxc0WERldfUqbh/H3FxPzZIE192Nvr3x+bNqFMHAD58wPXrOHqUT8nevREWJiIRAti/Hx06wN0d9+7Bzk6yYAiCEE7Zd5apYORIxMfj6VOYmeHwYUREkMFCQl5mzcKFCzh1ClLsnRkaiq5d0aPHj3/u2YMePcB399bOnXHrFj58EF3n+fNo2hTOzvhl1IIgCJn8fCJMTEzcsWOHyBu2bdsmz3hE0NXF6NFYtQqbN8PX98dg4atXqF5dgUERldDSpYiOxqVLMDGR+N7wcNy9i5SUn1f27MHy5fwLa2ujUyccP44RI0TXnJSEhg3RqBHS02FoKHFgBEHw9TMRvn//PiEhQYGhiGn8eNSvj/nzUavWj8FCDw8yWEjQaf167NiB+HiYmUl246dPn44cSZo7t8OlSz9nlT19infvhJ211Ls3wsPFSoRsNtLSYGeHhg2Rng4Bu7kRBCGZn4nQ39/fXxXOOjI1xaBBWLsWvA3+eYOFAQFvli4tsCODJ4TMdu3CX3/h0iVhMzkFad7c9+XLjg0b7m3YcG/ZCvv3/zFlhq+uXTFiBLKyYGoquglNTTx8CBsbNGqEJ0/A7MFoBFE5qeSv0dSpCA/H168AYGaG8eOPRUV1t7f33bNnj6JDI1Tb0aOYOROxsbCxkeb2798BqNesyS178cCBcuvof6Wriw4dcOKEuK3o6+PhQ3z8iBo11vXqNejbt2/SxEoQxP+JSISfP3/+rzxmwhLOygrdumHjxh//1NRMBRwoquWdO3cUGheh2uLiMGYMTp1Cw4ZS1tC+/cngYNdjx35OZb56FWpqcHUVcaPIlfUVmJtjxowdnz7FHTtmFhLyu1TBEgTxA//lE/n5+VOmTNm3b192dnaFl4ScaM+kmTPRrh0mToSODubPn//wYeiRI26dO5ODmggpJSYiKAhHjsDFRcoaiosRF1f95s0eeno/L+7Zg0GDRN/r54fQUGRniz7IopS3twObfY/L/fzq1VDJgyUI4if+T4Tjx4/ftWvXpEmTvLy8/P39w8LCOnfurK+vHxYWxnB8gjRoAHd3REQAgLq6+qFDW1xdR86eraHgsAgVlJqaunLlCX//kgMH4CnDYQ8XLqB+/R8LB3mKi3HwIAYOFH2vgQHatsXJkxI016pVq5cvE7dt256aOmrpUomjJQiiFP9EePDgwZUrV86fP9/GxqZBgwZjx46NiYkJDQ09ePAgw/EJMWsWVq5E6cGlq1cjJQWfPys0JkLVPH36tGPH8TNmRPfoscXbW6aqDh+ueFhmbCwcHMQdbpS0dxSApaXlsGENt27F3LnYu1d0eYIg+OKTCD9+/Pjt2zdvb28AGhoapXutTZw4MSEhISMjg8n4hGjRArVrozQ1e3nB3JzP6W4EIYS6ulZuboGOTpabG7/l7mLjcHD8OHr1KnexwnETwvXogfPnkZMjcdPBwZg1C0FBOHtW4nsJggDfRKinpwegsLAQQM2aNdPT03nXtbS0AHxWpmeuGTOwfDlKRy3HjcP+/QoNiFA1J07UcXU9cuHCjJEjh8pSz8WLsLb+edYggJwcxMaiTx9xazA2RqtWiImRpvXFizFwILp2hYJOziYI1cYnEerq6tra2t69exdAmzZt4uLijh49mpGRMWPGDC0trXr16jEepEC+vlBTQ2zsj3/OmgUOB5s3KzQmQnW8fInFi7Fjh7W7ewsZq/q1X/TIEXh7VzxuQjgpekdLRUbC3R0tWoi1WxtBEGXxHyMMDQ19/vw5AB8fnw4dOvj7+9etW3f79u1LliwxVLKdnXhnM/Gw2ejaFeREb0JM48Zh8mTpF0uU4nJx/HjFRChRvyhPz544cwb5+VKGER8PS0s0bYqCAilrIIiqif/yiSlTppR+feLEieTk5PT0dGdn50aNGjEVmLgCAjB/Pq5eRatWABAWhtq1cfOm6JVbRBW3fz/S03HoEA1VJSTAwgJl+0revcONGzh+XLJ6zM3RrBnOnKk41igm3gZs1tZwdsajR9LUQBBVE/8nQt4A4Y8SbLaHh8fAgQOVMAsCUFPDpElYufLHPy0t0agRyuRxguAjOxvTpmHTJmjQseLm8OGKY4H79qFnT0hxiLUsvaMANDWRloa3b4VtbUoQRAX8E2HTpk0HDx587tw5JVk+L5yPT0Z0tIeVlQdv45vly3H5MsqcK0wQFU2eDH//H70IMuJyceQIKmzTK0W/KE/v3jh1Sqa+TXNz3LiBq1cRECB9JQRRpfBPhL17946Nje3YsWP9+vWXLFny8uVLhsOSyO3bySyWV2Zmh6SkJADdusHQEH/8oeiwCGUVH49z57BkCT21Xb0Kc3M0aPDzyqNHyMyEdKsSa9RA06Y4f16mkOrXR1wcDh/GzJky1UMQVQT/RLhs2bJ3796dPXvWxcVl8eLFNjY2rVu33rJlS25uLsPxiaNbt24+PgXGxjk9/n8Q6vDh2LlTsUERSqqwEKNHY906GBjQU+Gv80V378bAgcKOmxBOxt5RHi8vRERg5Ups3MjNl3r6DUFUDQI33VZTU/Px8YmKinrz5s2aNWsKCgpGjRplaWnJZHBi0tXVPXFiTVHR2pKSH5s8LlmC/HxERio2LkIZLVwIR8efZ8fLiKJw5Ei5AUKKwr59UvaL8vTpg+PHUVQka2yDBmHChJehoS0NDBzDwjbIWh1BVF6ij2EyNzcPCgoaMWKEtbW10h74oqGB5s1x9eqPf2pqon172vq+iEojLQ1bt+Kff2irMDkZOjooO43s6lVoacHZWfo6LS1Rvz4uXZI5OKBp0/NALS73940bE2mojiAqKWGJkMPhnD59un///jVr1hwzZoy1tfWOHTsYi0xSbdsiIeHnP8PC8PQpHj9WXECEkuFyMWoUli+X5sRdQQ4fLjcnJSsra/nyMwMGyLqOj5beUQBDhgzx8mKbmZ1//Hjt+PE0VEgQlRL/RPjw4cMFCxbY2dl17do1MTFxwoQJT548iY+PHzp0KLPhScDLq1wirF8ftraYPFlxARFKJiwMmpoICaGzziNHyg0Qurl1OXXqwrlzw2Sstk8fHD36c0N5qbHZ7EuXjnz6dPLgQfONG6Wcv0MQlR7/BfXt2rXLzs728/MLCwvz9fVVk3rcn0EtWyItDXl5KD0NbvFiDB6MoiJoaio0MkIJvHyJJUuQkAAWi7Y6b9yAmhocHX9eKSjgslh66uqyZjAzs5yCghGOjvlnz26kZWDe3x8pKWjdGg4OuHULujJtME4QlQ3/J8KwsLAPHz5ERUX5+fmpRBYEoK0NJyckJf28MmAAtLWxaJHiYiKUxrhxmDiRht3Uyjp8GH37lrsSGnrCx6fx8ePhAu4Q15UrV4qK6jx61Pno0RMyVlXKxQXp6fj2DbVr49UrumoliMqAfyLs27evAV2zyxlUoXcUQFAQNm1SUDSE0jhwAP/9R/8RXYcOVVw48fhxzf79e+vr68tYc6tWrWxtH+nqHu3WrYuMVZVVvToyMlC7Nuztf84sIwiCf9cogIcPH0ZHR7969aqo/DzuzUp8uEPbthV33F65Elu2IDoafn4KiolQtOxsTJ2KAwdo7iG/fRscDlxcyl1MTqZnDbuRkdHlyyfq1YO1NQ21laWpidu3ERSEtm0RHo4hQ2iunyBUEf9EuHbt2ilTprBYLAsLC03VGWFr1Qo3bqCgANraP67o6aFVK8yZQxJh1TV5Mnr1omc3tbIOHULfvuVGHL9+xfv35baYkYWZGXR18fIl/bkQQGQkGjZESAhSU7FuHf31E4Rq4ZMIKYqaP39+t27dtm/fbmpqynxMUtPXR6NGuH4dbdr8vLhuHVxd8fo1rKwUFxmhCA8fPly8+ND58/5PnjSmvfLDhyvuXpScjGbNpN9Q5lfOzrhzRy6JEMCsWahfHwMG4O5detYsEoTq4jNG+OHDh+zs7FmzZqlWFuT5dZjQ2RlWVmQdRVXUrduwffucWKxhRkY018ybn9y8ebmLKSlwd6ezFScn3LlDZ4UV9OmDlBRcv87V1h5lbu6SnJwsx8YIQonxSYSmpqaGhoZKu4mMcBWW1fNMn47jx8HhKCIgQnE0Ne20tHY6OdnSXjNvvmiFlRgpKRVTo4zknQgBuLjgn38OFxZ++vx5ckDAStE3EERlxCcRamhozJkzZ8mSJcq5xbZwbdogKQnFxeUujh0LNTWsWqWgmAhF+PIFWVmR0dF/nzmzl/bKf50vCuD6dfqfCG/fprNCvnr29NLRucNi/fnq1fhatXDqlNxbJAhlw3+yzIcPHx4/fmxvb+/p6WlmZlb2JWWeNQrAyAi2trh5s+KfpL59sWYNZsxQUFgE4+bPR79+6NiR/hG2x4+RnY2WLctdfPECbDbo3ZTe3h6Zmfj2DYaGdFZbQfXq1fPznwHIzcWYMejeHbVqYc8etG0rx0YJQqnwT4QJCQm6uroAbt26xWw8NOANE1ZIhH//jerVER8PLy8FhUUw6PFj7NmDBw/kUnlUFHr3Brt8Z0pKSsXUKDs1NTRujLQ0eHrSXDNf+vqIjMSqVQgOhrc3GjTAkSO0TYIlCGXGf0F9cnLycwEYjk8KfIcJzczg6orp0xUREMG46dMxcyZq1JBL5b8eQAggJQUtWtDfFm/iKJNq1MDp03jwAJqaaNwYvr7IymI0AIJgnuhjmFRO27a4coXP1JjVq3H9Oj5+VERMBIPi45GWBjkdtvD0Kd6/57MqUU6JkIH5Mnw1aIDbt3HqFNLSUL06goJw+/ajAwcOKCAUgpA/gYnw5cuXc+bM8ff39/X15V3Zu3dvTEwMU4FJz9wctWrh7t2K19u2RfXq5KGwkuNyMXUqli+HlpZc6udNk6mwWJDDwe3bcHOjvzlm5ssI0rkzXr3Cli04ceKdi0ufAQN2tWgxRgWn0BGECAK7Rh0dHTdv3pyRkXH3/ynl9evXkyZNYjA26f26mpBn7Fju7t3RBw8eZDwigiGRkVBXr7gXNo349oumpaF2bblMaXF0xP37Cl75ExKCuLiXLBaLohxSU6sZGEBHB7a26NYNy5aV+8RZUlLy7NkzxUVKiPaYnNHKD/9EOGbMGDc3t+fPn69evbr0YteuXR8/fvz+/XumYpMe32FCAOnpY0tK9gYELD1xgrZN/Qnl8f075s3DqlV0nrVU1osXePmy3L5FPLQvpS9lYAALCzx9KpfKxefu7r5t26RRo4q/fJmanY3wcLRrh7dvsXIlnJ2hpgZzczg6vtTUrOPi0s/RsbWCwyUEqFatafPmQSwWfSdTVxZ8EmFWVtatW7cWLVpkaGjIKvMXxdraGsDbt2+Zi05aXl64fBkUVfG6uroakENRbDa7Eg6OEitXwsNDjnMsDx6Evz/Uf5lqTftS+rKYny/DV0hIyKZNYYaGhoaGGDgQ4eG4eRNZWeBykZKCCRMAXAMsgBkPHpiq5m4cldynT/j8OQdYAZgdO6YCzzNM4pMPeMdN/HqUzOfPnwGo//pnQPnUrAkjIz6z5zdtWrdiRRc2exfQTRFxEXL0/j3CwvDnn3Jsgm+/KOT5RAjFzZcRX7NmmDsXd+8G2NqqsVi7WrRY0qgRtmwBl6voyAgAAJeLXbvg6Ah7+9ks1mw9vd6TJtXo148cS/kTn0RYvXr16tWrR0dHAyj7RLh37159fX0HBwfmopMB32FCFos1ffrvXl6NaT+ajlC42bMREgI7O3nV//o1nj9Hu3YVr+flIT0dTZvKq13lT4Slnj27vmRJXIMGjidOICICLVsiJUXRMVV5qanw9MS2bYiNhZHRiD17UrS0Fpw/j0aN4OKCBQtQ/py9KopPImSz2VOmTFm0aNHChQufPn3K5XLv3r07e/bs+fPnjx8/XlVOZRI0TAhg2zY8fozUVGYDIuTp7l2cPIlZs+RV/+XLiS4ubXR1Q1isihNXbtyAoyM0NOTVtGInjkpq0KCS48dRrx4SEzF2LHr0wODBZM2SYmRlYcIE9OiBMWMQH4+SEnz+DF/fkt69cfAgFixAUhJu3EDTpoiNVXSsCkfxw+Vyp0+frlHml5vFYg0bNqyoqIhveTnZtm1bcHCwdPemp1M1awp81cmJatlSyqiY8e3bN0U1XVJSkpeXJ4+a5fdNdepE/fuvnOqmKIrq3/934KKRkW9GRkaFl1asoCZOlKy2oqKigoIC8cubmFDv30vWhEJwudycnJw+fajNm39c+fKFGj+eqlaNWrOGKilRZGxyeu8p8PdUCA6H2rmTsrCgxo+nvn79cXHoUOqvv6hv374lJVF2dhSX++P6iROUrS3VrRuVnq6gcCmKkttPksPhiJO2+M8ZYbFYK1asSE9P3717919//bVt27bHjx9v27ZNQ36fe+lWty40NQVOt9u4EcnJePmS2ZgI+Th9GhkZGDlSjk2MHh3CZv/p5+dQu3btCi9dvy6XpfRlOTryWRertIYNQ3j4j6+NjbF2LRIScOoUmjdHYqJCI6vsuFwugNRUtGqF8HDExWHtWvDOIPvyBcePY+hQAHB3h77+zw4zPz88fAgfH7i7Y8ECFBQoKHqFEjbzxdLSMjAwkLFQaMfrHbW35/OShwesrREaipMnGQ+LoBWHgxkzsGqVHDsnAbx/7+rrezYyks9LyclYvlyOTeP/E0d9fOTbCl06dcLo0bh9G87OP640aIC4OERHIzAQbdpwGjXaYWamPnz4YDJ5m0ZJScn+/uPy87V0dI6vWmU6cGC5RUTh4ejRA9WqIScHAIKDER7+c+NlTU1MmIA+ffDHH2jSBKNH3/369dSwYf1tbGwU8J0oAv9EmJGRweG3iNfExMTExETOIdGGlwiHDeP/6sqV6N9f7lv7E/K2ZQvMzdFNzrOAY2LQpQuf65mZyM+HLf0nHpbj5ISLF+XbBI3YbAwdih07sHZtuet+fmjXDoMGHZk9O1lbu8DSslrXrl0VFGMldPbs5czMQZqa16OiHnbqVG4JEUVhyxbs3v3zSlAQFizAly8o++fc0hK7duHcOXTpMrSkZNbx48PT0s4zFb6C8f9E5u7ubsePqampsbFxcHAwbymFkmvbVtifjz59YGICFdkqh+AvJweLF8v9pEmKwpkz6NyZz0vJyXB3l9f6/VKqNV8GwPDh2LePTyebvj7mzLExNr5dWPiwuLiOIkKrtMzNhxobP/r9d3sfH48KL8XGQl+/XAe+qSk6d8a+fXzq8fGBq2stLa39BQX8OtMqKf6JcNmyZSYmJh07dly3bt3+/fv//vvvFi1aWFlZbd68OSQkJCoqyt/fn+FApVC/PjgcZGQILDBzJnbvRkkJgzERtFq6FJ07o1kz+baSmgpjY/6PfXJdSl+qcWM8e4bCQrk3RBcrK7i64uhRPi+5ubk9e3Zmy5ZzM2Y0Jevu6fL9O5YvNz9xYsPq1fN+7XDesAHjxlW8ZdgwCDpb9tq1EzdurFRX38R3LKBy4juFxs/Pb/LkyWWvcLlcf3//cePGURTF23o7NTWVllk9Qsgya5QnIIDatUtYAT09avZsWVqQFzJrVKRXrygzM+rlSxqr5G/RImrKFP4v+fhQp05JXKGks0YpimrShLp5U+KGGMabNcr7OiqK6tBBWOFRo6h+/ZiIiqdyzxqdO5cKDOT/UkYGZW5Olf5ClwbM5VJ2dpSQv+J371LVq1P379MbqUBKN2s0Jyfn5MmTISEhZS+yWKzg4OD9+/cD6Ny5s5mZ2cOHD5lJ1bIQspqQZ+RIrFvHVDQEraZPx9ix+GUWJ/0EDRBSFFJT5T5llEdJNloTX48euHcPQrbgXrsW//2HsDAGY6qkXr/Gxo1YupT/qxs3YsgQ6OpWvM5iYejQn/N7f9W0Kf7+G/7+qAoP7nwSYX5+PkVRv44Cfvr0KS8vj/e1kZGRWoWjaJRS27aIjxdWYMUKFBZi40amAiLocP36jRo1Whw61GHUqGx5t5WVhQcP0JrfPtKPH8PUFObm8g4BUKn9ZXg0NTFoECIiBBbQ0sLhw1i6lKypkNXUqRg3DnX4DbkWFiIiQuDKouBg7N+P/HyBNQcGok0b+S5MUhJ8EmGNGjUcHBwmTZr0ssw6u/v378+fP79t27YAcnJyXr9+/euCKiXUuDGysyFkn3ANDfTujUWLGIyJkNnJk+c/fBimqWn14sUv+8nSLSYG7dvzP90wOZmhx0GoYCIEMGIEtm8XNgZfpw62bcPAgWTrGeldvYqrVzFlCv9Xo6Lg7Iz69fm/amkJDw8cPiys/rAwPHuG9etljVPJ8Z8ss3379mfPntWrV8/Nza1z585OTk5NmzYtKSlZt24dgBs3brRr185NHueQ0o3FgqcnLl8WVmb9enz4AHIukwqpUWOIkVFSUJBlC/knIkH9omBkKX0pZ2fcusXnQBVl5uAAGxsR23d17YrAQAwYoOAzF1UUl4uJE/HXX9DT419g40aEhgqroezuB3xpa+PwYSxZUskf3PknwlatWt2/f3/GjBlWVlZfv361t7f/888/7927Z29vD6Bdu3axsbFacjoCnG4ihwlNTNC+Pcg23Kri+3f89ZdFdPSOjRuXyrt/nsvF2bP47Tf+rzL5RFitGnR08Po1Q83RReTfWQBLloDNxuLFjARUuUREQF0dAQH8X719G2/ewNdXWA1+fnjyBI8eCStjbY2tWzFwID59kj5UJSdwZxkrK6vFleK92bYttm0TUWbbNtjY4Pp1JqbCEzJavhyennxOx5WHlBRYWAgcfXn4EC4uTITBw1tNqAojEj8FBGDaNLx9i1qCz4Jls7F7N9zc4O4u8OGb+FVODubNw5EjApexrl+P0aMh/LOiujqCgrBzJ5YtE1asWzdcuYL+/XHmjIgKVVTl3+LIyQlv3+LDB2FlrK3h4oKxY5mKiZDW69fYsEHgBDnaxcQI/EB96xYcHKCjw1AkUMGJowD09NCnD0QuR6teHVFRCA7GixdMRFU5/PknOnUS2Cfx9SsOH0b5uf/8DRuGiAgUF4tujs2W73mfCvTziTAmJmb58uUjR44MDAzs1atXVlYW3xvihc/CVD5qamjVComJ6NVLWLFNm+DujvR0VJnd9VTS1KkYOxbW1gw1d/o0Vq/m/5JcD+Ply8lJxLwG5TRsGAYOxPTpIvbfadkS06YhIAAJCfynJhFl/fcfwsOFbcW+Ywe6dkWNGqKrql8f9evj1Cn07CmsmJoaIiPRvDlatOC/y5JK+/lEqKGhoaenxztfQldXV08AxYUqPZHDhACaN4eNjYiBZUKxeBPkpk5lqLkPH/DsGTwq7lf1AzN7ypSlchut8bRoUe6sAyEmT0adOsz9/1VpU6diyhTUrMn/VYrC5s0YM0bc2sQZygVQowYOHMDQoZXxwV0ei/npIvvOMjzXrlEuLqKLHTtGqan9PL5LscjOMhVwOJSbG7V/P43hiBARQfXtK/BVe3vpN92QYmcZiqJKSig9PSo7W8pGGVB2Z5my1qyhgoLEquHbN6pBAxG7QUmnMu0sc+ECZWNDff8usEBcHOXkxP8lvgHn50uwSdNff1Hu7lRhoViFxad0O8tUPm5u+O8/fP0qoliPHjA1xfjxjMRESGj7dmhqol8/5loUsnAiKwvv36NBA+aCAaCmhkaNcO8eo43SIigI0dH48kV0SQMDREVh6lQ8kPsCUVXF4WDiRKxaBW1tgWU2bMDvv0tQp44O+vXDrl1iFZ46FVZWlW2avcBEePv27ZCQEHd391atWvGuhIeHR0VFyd5kTk7Orl27Ro8eHRwcvGXLlmKRo7QyU1dHixZirYOZPRv79okeNyYYlpOD+fOxZo3cz3koxeHg3DmBCyeuX4ebG5g/TU8V58tA6FkHv6pSO3tJYetWmJpCyKkHr17h8mUMHChZtbzeUS5XdEkWCzt2IC6u3LlOqo7/r3JcXJy7u/vly5eNjY1L95f5/v37jBkzKJnX9F64cGH//v1OTk7e3t7//vtvcHCwjBWKQ5xhQgATJkBbG3Pnyj8gQhKLFsHXl9ExuatXUbeuwEn/TK4gLEsV95fhEXLWwa8CA1G37tmaNVt27TqEK87f5irj2zcsWiTi3LHNmxEUJHCJvSDNmsHYGJcuiVXYwAB79hQNG9bXysrzypWrkrWklPgnwnHjxvXo0ePBgwd//PFH6UUfH58XL168e/dOxib9/PxOnz49ZsyYIUOGREREREVFFRUVyVinSCI3HS01alTl309ItTx/jogIprfBE9IvCmb3lClLdRNhhw7Iy8PNm+KW19A4kJ+/KiHh1adKvIpbcgsWoHt3YeeOFRcjIgKjR0tTuZhTZnj09NI1NUvevJkVHs7vtC1VwycRfvz48cmTJ9OnT9fQ0GCV6YribS4qeyIse1xWZmamsbGxpqamjHWK5O6O+/eRkyO65LJlKCoiuVCJTJqEGTMETpCTk9OnlTEROjri3j2x+q+UjcizDipYuPD3Jk2Wcrltnz2rLs+4VMmzZ4iMxIIFwsocOoSGDeHgIE39gwbh9Glxt4+xt7cPCKhXrdqmDx+GSdOYkuGzswyv8/PX0x0zMzMB6Ii3hPjjx4+/9mno6+uXXYDx7du3SZMmLVy4UEg9iYmJFQ4B7tq1a4CgPYWEcnbWuXixqH170Xsa9uqlvXCh2pAheVK0Qpe8vDwWYwNi5XE4nMLCQnl0SUnxTV26pHbvntaOHfm5ubSHI1BmJuvVK93GjfP4NpqRwWazdYyM+L8qjuLiYi6XK8XoOJsNc3Pdu3cL6tVTxmRIUVS+4LMM+vVjtWypO39+3q9HAv2qfn37a9eiYmPV+/enLtRycUUAACAASURBVF/ONzOTdURGTr9QjP2eHjsWvXRp9oQJgfr6LCFvvLAwnXHjinNzBe50LiRgNTV07qy9YwdnzBix3plr1sxfvhwdO+quWVM4fLisEyvk9JPkcrni7AbKJxFWq1atTp06Bw4ccHV1LRvZli1bTE1N6wvaybw8Hx+fX/s0Jk+ePOX/26Tn5+f7+fm1b99+tNDH+Dp16vTv37/sFQcHB11xfpN+4e3NSknR7tZN9G/U5s0wM2PFxuoKGZGWNw6HI923SUvTampq8mhd0m+qpAR//MFes4YyMWH0R3HxIuu332BgwL/RtDRWy5aQ5efDS4TS7dbr7Mx6/FjH0VEZt9/mzUQX9JOpVw8eHoiN1Rs0SNzg/f2RnIyQEN3YWErGnb3k9AvFzO9pSkrKmDG78/LsdXV36+qO4lvmw4cPgwbNvXu3Zo8ec7W1BfaxCQ94xAhMnKg+ZYqGmIHp6uLAAXh6ajVvriFoxa2Y5PST5HK54sxr4ZMIWSzWnDlzRo0a9fXrV1tb2+Li4tjY2KioqIiIiBUrVqirC9yetKw7Qocyvn//3qNHD1tb2/Xr1wv/FFC7du1+NE2Z9/LCokVgs0V/6DAyQvv23NDQaFvb2q6urrS0Lik2m/3rQzkzKIqSU+uSVrtpE2rUQPfuLIDRh+PYWPTqJfCtwusXFeeNJAjvhyDdT9jZGXfvIiBAMb0Fwol85wwfjrVrMXiwBMEvW4YOHbB8OUvGKWxK8paWTs2aNYuL3+rp5TRo0FlQc1u3Rl644KSllXLrVoqnp6egqoQH3K4dCgqQmsoWf2Ja/foID0dgIPvGDZnO5pTfT5IjzskmghYYrl271sjIqLSYtrb2vHnzOByO7CscCwsLu3Xr1rdv3+LiYuEl6VpQz5OXRxkYUPn5YhVu3ToAGMJm133y5AldAUiELKj//JmqUYNKS5NHIMIUFVEmJtT79wILtG5NXbggYxPSLKjnOXaM6tpVptblR9CC+lLFxVStWtTDh5JVm5lJWVlRZ87IFJtKL6h/9YoyMfn67NkbIWWSk1PU1Zvb2LTNysoSUkxkwEuXUqNGSRzh5MmUry8lS35Q0gX148ePf/Pmzblz5/bv33/q1Kk3b94sXLiQlowdFRV18uTJ5ORkBwcHOzs7Ozs72SfgiENXF02aICVFrMIcTi6gyeXqFxQUyDkugr+5cxEQgCZNmG73yhXUr4/qAqZocDi4c0fYtD15U9GN1nh4Zx0IObaerxo1EBmJoUNV7xQqumzbhsBAIzs7wUd4AHXrNtfXT3r6NN7ExESWtoKDERUFSce/V6xATg6WL5elZUUS1s+pp6fXoUMH2pvs0aPH8+fPy16pVq0a7a3wxVtE4eUluuTZswdmzly8cePe06ebNm0q/8iI8h48wKFDitleRPjCibQ01KkDQ0MGAyrP2hr5+fj0SaZuKAUKCUGrVin+/uotWkgw6ODtjbFj0bcv4uMh/znmyoXLRUQEjh0TUSw+Hm3bsmU/I8nCAi1bZs2de23JEm/xN5dWV0dUFNzc0Lw5OnaUNQbmKWAUysDAwLY8MccdZSfmsnoAenp6YWHLZ89uOm+eWIsuCHpNmoT582FmpoCmT58WdpapopbSl2Kx0LSpqq4mBPDiRVx29vwOHaZduXJFohv/+AMWFiizsLmqOH0aNWvC2VlEMTE/4ovjzh2/sLDkLl2CJLrLwgKRkRgyBG/e0BMGk6rEXqOlWrdGcnJ+fr7AucUVLFwIU1NG97ckKIrauzfn9WuMHKmA1l+9wsePwno+FbWCsCzVXVYPgKIoDQ214mI1SsI9qlgsbN+OY8dw6JCcQlNSW7eK9btAYyLU1aW4XI2iIomX6LRrh9BQ9O0L+W+RQrOqlQhv3kwoLPS2sXHnrYkUx8GDiIvD5ctyjYv4oaioqEGDNoMHd+7YcRdT3QTlnDyJzp2FbSKq8CdCqHgi/O233zZu/MPCYnGbNm0kvdfEBEeOIDQUjx7JIzRl9O4dEhNFfxb//BmvX4t+ahTTtWsnHB2bTpki6jxlfmbPRvXqmD2bnkgYU7US4Y0btzmcbjk59SoMUgrRujU6dCAPhQzJzs5+9YrL4Yz48EG8SU10Ez5AmJuL9HQofMxYRbfeLhUU5Pntm/uHD9Lc6+SEhQvRrx8EL9yvVLZuRf/+ojcOjY+HpydkHyDkMTc39/fvmZpqIMW9vC25Dx9WsUOkq1YiHDUquEsXrp1d29IjNcRx9Ci+fMG8efKLi/jh8eNqWlqTBg16uGrVLOZbLyxEQoKwof4bN+DsDA1xVxvLS5MmePpU9XqfSrHZ8PAQ6zQYvsaMgaurYnrOGcblYvt2DBNjCzMa+0V5WreWvhvMxAQHDmD0aDx+TGdIclW1EqGBgcGGDQs+fPhdor189PSwdCmWLcPnz/ILjUB2NoKCEBnZNzJyRS1B5z7IU3w8mjYVNkMnJUXx/aIAtLRgY4OHDxUdhwxat4aEc2XKWb8et29LsHOpioqNhYUFXFxEl7x0ieZE2LIl7tyR/rG7eXMsWKBKD+5VKxECqFMHGhpIT5fsrsmTUauWsDPACNmNGwdfX3TrprAAhPeLAkhJYfQoKCFUejUhZHvgAKCnhyNHMGuWBMdZqKItWzBihOhiWVl48QL0boHFW3V9/br0Nfz+O5ycVOac8yqXCAG0bIlr1yS+6/hxXLmCuDg5BEQAR47g2jWsWKHIGIQvnACQkgJ3d6aiEUql58sAaNECDx/KtDCpfn2EhcHfvzAmJrGwsJC+0JTFu3dISBBrdkJCAlq1Au2Ty2R8agewcSOuXcPSpY/SJX3yYFxVTIQeHkhKkvguZ2d064YBAyDzycRERW/eIDQUe/ZAX19hMaSnIycHTk4CC7x7h/x82NgwGJNgqj5fRksLzs5ITpapkn79UFDg37PnwbZtK2FfzbZt6N8fBmJMWKF9gJBH9kSop4fJk8/NnRvq6jrgtnL3YFTRRCjFEyGAgwdRUIBp0+gOqGqjKIwYgXHjFDz8dvIkfH0hZOw4ORnu7sIKMMnZWbW7RkHH31kApqZ5xcX2jx/nloi7Nlg18KbJDB8uVuH4eHh70x9D69a4dg3ibFgthJFRtpaWRXa22fXr32iKSy6qYiJs1gyPHkkziqupiTVrsGYNGNkbtapYswbZ2Zg5U8FhCB8gvHnz1qhRnd69myLpMnA5qV4dmpqqvfcmLYnw/Pl9Gzboenjs79ULlWlX4DNnUL26WMN+2dl49ozmAUIec3PUqoW7d2WqpHdv/337BixaNGnu3LZSzxNmQFVMhFpaaNQIqanS3DtiBGxt0aMH3TFVVQ8eYOlSRETQtgRKOt+/IzERPj4CC6xcGf7hw6z//nv46tUrBuMSphLMl7l+HZIfTlxOzZo1R48OPnGipo4OfH0l3ipaaW3dKtY0GQAJCfDwkNf+q23ayPphhcVi9ejhN2eOz/796N0bZ8/SFBndqmIihAy9owCio5GaiqNHaQ2oSiosxMCBWLkS9vYKjuTCBTRrhjLHjlU0evQANntOy5Y1LC0tGYxLGFWfL2NkhLp1cesWDVVpaGDfPtStC19ffFPqHjixZGbi0iUEBIhVWE4DhDy0PLXzeHvj0CEEBorePVwhqm4ilGK+DI+DAwICMGSIrL3nxJw5sLHB0KGKjkOMhROmpp729ldiY3eoKfbRtQxVT4Sg44GjlJoawsPh4oL27VV+vW94OAICxJomA/knQjFPKRCztpgYhIYq426xVTQRtmyJq1elvz0yEhSF0FD6Aqp6Ll/Gnj3YvFnRcQBcLjcm5rvwRHjtGjw8mApIPKo+cRS0PnAAYLGwZg3atYOXlwqP4lMUduwQt1/02zc8fgw3N3kFY2MDDQ2IvR+laM2a4dw5TJqEHTtoq5MWVTQRWltDXR0vXkh5u5oatmzBtm347z86o6o6eJvIbNsm8PxbxmRlZVlbe2RkeL95c05IsaQkZVlBWMrBAW/eqPaoWJs2uHyZzvVILBZWrkRQENq3V9WZRGfOwMRE3MkvV66gRQtoackxHk9Pmo8caNQI585h/nysW0dntTKqookQgLu79MOEAAYMQKNG6NmTvoCqktBQdO8uYvU6MzIyMrKyanK5AVeuCNtFIylJ6Z4I1dTQsCHu3VN0HDKwtIS+Pv07Us6YgZAQtGlD56MMY8SfJgM594vy0PvUzuPggMuXERaG1atprllqVTcRyjJMyHPqFB48wK5dNAVUZRw+jFu3FLyJTCkXF5caNTp16fJh8uTRgsp8/YrXr9G4MZNxiVZQUPDx47BBg3q9VtFnHwC0DhOWNW0aZsyAlxfu36e/cvnJzMTFi+jfX9zytG8x+is5/Q+ytkZCAnbsUPy6KZ4qnQhleSIEUKcOhg3DmDEqfA4A83ibyOzcCR0dRYcCAPjwAVlZoUeOLDcxMRFUJikJbm70b2Elo2vXrn34oPP8edc9ew4qOhbpybjpqBCjR2PFCnTqJOtKOCZt346+fWFoKFbh3Fw8eCD3bSiaNMHHj5DuzCzhatbEhQuIjcX06fRXLqmqmwibNcODB7Jujr5pE9TUPjZsOOv48eM0xVWZffqUNXgwJk1Slq2rAZw4gS5dRIyyKGG/KIBmzZrVrZuho7O7Z08l6GKWljx63koFBuKff9CxI5KTkZWVJa9maMKbJiP+8VKJiXBzg7a2PGMC2Gy0bCmv/0fVq+PSJVy+jDFj8OmTIv8HVd1EqK2NRo1k3b2exYKGRs///qvWq9eEXJWetyB/rVr1sLPzT0joqlR71B0/Lnp7hKQktGzJSDSSMDQ0vHo1Wl39kr29g6JjkV7DhsjJkePEln79sH59gYeHtY1N56Cg3+XVDB3OnoWeHpo1E7c8AwOEPHL9sGJsjDNncOjQCDu7zjo6dYsU1L1WdRMh6OgdBWBhoQkkUJTx3bvynLyl+m7dekBRKzmcNKVZiYe8PCQkoHNnYWUoSlmOIfyVsTFMTCQ+U0ypsFho1Ur6Q3rF4eT0CtCjqLnnz9Oxel9utm7FmDESlGdggJCHN7lXfgwNoaZ2k6IWFhYavXv3Xo4tCVbVE6GM82UA3L17LiKiZ9eu5729NWTcTb8S43JhZ/e3icmsefMk+V2Xs9hYeHjA2FhYmUePYGqKGjWYiklCTk6qNAbGl1wfOADY29sPG/absXF006YRcmxGNu/f4/x5DBggbvn8fNy7x1BHRfPmePRIpjOzRDpyJMzael2tWuMXLKjN5cqxIUGqdCKU7mDCCtTU1IYMGXLypJmfH9q0kfVkmcpq/XoYGvqlp59dsOAPRcfykzj9okq4lL4sR0eVX1Yvp3mJZW3d+s/Tp1tevaqvnPt7Adi+HX36iDtNBkBiIpydGZpxRsuZWcK1atUqLS3myZNhGRkYMQLM58IqnQjr1gWLhYwMemo7fBgdOqBNG5mOda6UXr7EokUIDwdbmd5uHA5iYuDnJ6KYEi6lL8vRUeWfCF1d8fw5vn6VbytaWti4EWPHyr0hKfj5Bc+b16JOnSPi3yKno5cEYeDDCgBdXZw8iefPMWoU08e+KtNfJkWQcVl9BTExaN8enp4kF5YzahSmTEHDhoqOo7z4eNjawspKRDHlnDJaqhJ0jWpowM2Nzl9DQby84OenLAvXSuXn5yck3CspWX/xogTLYBibKcMjv1UuFejq4tQpPH6MiROZaK5UVU+EtAwTlhUbi3bt0KYN0tLorFZ17dyJd+8wZYqi4/iFOP2iOTlIT4ejIyMBSaVePWRmqvyRC8w8cABYsQIxMTgnbDc9punq6trYBNaps3jZMnF/Sb5/x+3bjM5k9vRESgpDC6b19BAdjaQkTJrERHM8JBHS/1H0zBm0bYvmzUkuxKdP+OMPhIdDQ0PRofwiOlr0DnnJyXBxUcbgS6mpoVEj1d5oDQw+cBgaYuNGjByJvDwmmhNHbi4yMiampJxo0ULczbOvXYOTE/T05BpXOUZGsLNj7vxLIyPExSExkbkP0FU9ETZrhvv38f07zdXGxaFNGzRvrmI7PNEuNBTBwRIsjWLMrVs/UohwSt4vylMJekc9PHDzJkNHzPv6wsMDCxcy0ZY4Dh6El5dk05IZ7hflYezDCg8vF8bHg5llx1U9EerooGFDWZfV83X2LFq3/pFoq6aTJ3HnDubOVXQc/Bw/jl69RBdTzqX0FVSC+TL6+mjQAKmpDDW3bh1271aWgfydOzFkiGS3KCoRMtN9XcrYGLGxOHMGCxciX8Y9wESp6okQ8ukd5Tl3Dp6eaNYMDx/KpX5llp2N0FBs2yb3LaCkc+yY6AFCikJyslJPGeWpBIkQ8l+1XZaZGVavxrBhit8l+MUL3L8v2TEshYVITVVARwVvHJfhyZzm5jh3DmvXTq9evXOvXsPl1xBJhPTPlynr/Hm0agUXlyqXCydPRs+eaNNG0XHw8+IFMjNFP+o9ewY9PdSqxUhMMnByQloa03+haMfwA8eAAbCzU/wRKDt3YsAAyQ4UTEpCkybinl9PI0tLGBjg0SOm261eHcbGiXl5K2Jjb758Ka9WSCKkZ1m9EBcuwMMDrq4KeA8pysWLOH8ef/6p6DgEOHYM3btD5E5vSr6UvpSxMYyNVXujNQBt2+LqVUZXUq9fj/XrFTlyQVGIjJS4X5SxndV+xdjk3goOHFjXr9+ewYPXurhgxQpwOPQ3QRIhbGzA5UJ+nzUAXLyIZs3g5PS6e/eRl5kccVaE/HyMGIFNmxTwoVVM4iycgNIvpS+rEvSOVquGatUYTUu1amHhQgwbJpc/rOK4fBna2hJPJVPIACEP88OEPM2bNztw4N/Nm9ukpODsWTRvTv+sDpIIAfk/FAK4fBkcTu/o6Jbe3sHKM3VbHmbPRuvWInayVqDPn3H7Njp0EF1SJaaM8jg5qfxGa2B8XiKAkSNhYICwMEYbLbVzJ4YOleyWoiKkpsLTUy7xiMT8/6AK7Oxw9iwmToSvL2bOpHOaMUmEgJyHCXlYLJiafgfOcbnVDA3RsCH+/lvxY/X0unQpoWZNjw0b+ixdWqzoWAQ6eRI+PqKn8OTl4ckTODszEpPMKsETIRTxwMFiYeXKL9Om/Va3bpsHDx4w2fT37zh+HIGBkt2VnAwHBwm2JKVXgwbIy8OrV4ppnYfFwuDBuH0bb9+iSROcP09PtSQRAvKcOFpWRkby1q3tP38+lZqKZs2wYMGPxRsrVqCkRO6tM2DXruOZmbO1tArz82nav1UOxOwXvX4djo6SzWJQoMqRCNu0QUIC041mZd1SV6/38uXA06fPMtnu4cNwd0fNmpLdxfAWoxUwcGaWmCwssGsX1q7FsGFwclrl4NDu1Kk4WSokiRAA3Nxw/77c1/Pq6OgMHz7c1NTU2Rm7d+PbN5w9i7p1sXAhtLXh4oI9e1R77l9Q0EgNjfBBg5zs7OwUHQt/37/jwgV07Sq6pAr1iwKwt8e7dyq/0Vq9euBy8eIFo416enp26QI1tYRevQKYbFeK5YNQ6AAhj6KGCfnq2hV37lBPnkQ+eRIxcOA/8+fj9Gl8/ixNVSQRAoCODhwc5LKsXrj27RETg/x8xMTAwgIhIdDSQsuWcHUNNDFxDAwMYTog2SQmOowceXTDhiUsFkvRsfAXFwc3N5iYiC6pEkvpS/F2yakEWzd4ejL9d1ZLS+vIkfW//bbv2jULxhp98wY3b4o++aSC4mKkpChsgJBHURNHBTEyYg0Z4lu7dkBo6HAWCxs3on591KqFfv2wdi2uXIGf39DoaNHjXiQR/sBM76ggHTsiJgbfv2PTJhQV4datSxzOqb17lWlvYDHs3y/xmAfDxOwXBZCcrEqJEJWld1RRDxyBgdizh7nmdu5Ev34SnyZ4/Trs7UWcIy1vrq5IT1euo6w2bVr28mXSsmW9FyxAdDQ+fsTZs+jcGffvY/DgtNOnP9+4IXqonyTCHxiYLyMSm42QENy8CXV1DaAfMDA+XsEhiS81FQUFSp08OBycOiXWx/D0dLDZqF1b/jHRp3IkQib3lymrZ0+kpCAzk6HmpFg+CCXoFwWgrg43N1y9quAwhGCz0bgxQkKwZQsuX7ZgsR5wuaL7+kgi/KFlSyX6v1tc/CItbX9w8PIOHXD6tKKjEc+ePQgMhLL2iQJAYiKsrFC3ruiS166hVSu5x0OvSnBUPQAnJ7x9i48fmW5XRwfduuHAASbaSkoClyvNElVlSIRQvt5RITZvrtar17M5c0Qv1SSJ8AdbW3C5Cp4ZXJa1tfX27QgJQffuiI5WdDSicLmIisLAgYqOQyjx+0VVaCl9KWfnyrDRmpoazWdli4+x3lHe8kFJPzIWF1NJSSWKHSDkUfhqQjE9fozNm7FmDUtLjMnfJBH+pKjfQCG2bMGECejZE3v3KjoUoc6fR61acHBQdBxCiXMAIY9qTRnlMTaGkRHTUy7lQVHDhD4+ePsWjx/Lt5WCAhw8iEGDJLvr48ePNjat8vM9nj1LkU9cEvDwwK1bDJ2ZJTWKwpgxWLgQlpZilSeJ8CdlGCb81erVmDYNQUGIjFR0KILx+kWVWVoaiovFOmv++3c8eABXV/nHRLfK0TuqqAcONht9+8r9E+eJE3B1lXj4+dGjR58/25eU9IuPV/z4jZ4eGjbEjRuKjkOorVuRn4+RI8UtTxLhT4qdOCrE8uWYNw9Dh2L9ekWHwk9BAaKjEcDoKiyJHTsm1gGEAFJT0aiRxDP6lEElOKEXgLs70tIg5+Pn+OP1jsq1e1m65YOenp4mJg06dXo3cuRQ+mOSnKLmNInp/XvMmYPNm8EWO7+RRPhT8+a4d09JH/nnz8eyZRg/HmvXKjqUXxw7hubNYcHcKixpSDRAqHL9ojyVY+Kojg4cHZGiiC5ANzdoaSE5WV71v3+PpCRxP5CV9e4du7h41smTfxsrdvHE/ynVsvpfTZiAkSPh5CTBLSQR/qSjg/r1ceuWouMQYPp0rF2LyZOxZImiQ/lfe2ca0NS1LeCVgUHKrIRZVAYFwYlB8TKKUC3OCBSr4oxefW3pq7bWa7XW2We9KsUBQbhQpFIRUKwCEYIIWAUpIlZrlXmQRJE5QHLej1NzEREznOQkcX+/knP2XnuRgZW99hpeR/79onV1UFUlbCayYqXS90c5XKNAajjGxx9LMWQmPh4WLgQNDZEnpqWBvz/Q6VLQSSzc3aGggLSuHUNz9SoUF8O2baLNQobwNeTWO4qzaROcPAk7dsCOHWSr8ornzyE/X9ggFLK4eBHmzBH2/4jiGkIbG2hogPZ2svWQGBI3HMuXw/nz0CudovHipQ+CKP4M2WBgAIaG8ljJqLMTNm6EEydEPtpAhvA15DNepj9r18KZM7B7N4SHk60KAAAkJcHs2fLbehBH+P8jNTXQ2wujR0tZIelAo4GtLZSXk62HxLi5QVEROZXoLSzAygoyJSrgPDjFxdDaCm5uIk98+RKKimDmTOJVkgT5TKL49ltwdxfntUKG8DWmTZOL2upDs3Il/PwzHD8O//M/ZKuiCH7Rly/h9m3w9RVqsCKm0vdHObyj+vpgbk7aHyKlhELx0gcB4NdfwdNT7n5ryuExYVkZJCTAoUPizEWG8DUsLYHHg9pasvV4F4sXwy+/wIkTsHx5L4kt76uq4PFj8PMja32huHwZvL3hgw+EGqyIqfT9UY54GSC1fMnHH8Ovv0JbG5Eye3rg559h2TJx5sqbXxRH3gJH+XwIC4ODB8HAQJzpyBAORA7T6gdlwQK4dAni4909PQ/Y288gRYf4eAgKAhUVUhYXFpH+jyhuyCiOcmRQAKkbDn19cHOD1FQiZWZkgK0tjBkj8sTeXsjMFKpxmIyxtAQMg6dPydbjFceOgZqamD81ABnCN5H/Y0IBvr59FEo9hi14+pScavDy326Cy4XsbGH/j3C5UFYGju8uTCi/4DtCRS+0BgBOTl1ZWecfPHhAyuqEe0fFSx8EgJwcsLWV09wke/vSI0cu9PT0kK0IVFby9uyBkyfFr3WMDOFA5DxwtD90Oj06epeRUSGPl97RIevV8XYTcu5IZDJhwgRhvSV374KNDWhqSlknaaKvD9raylBo7YcfdrS2lri5hbSS0W54wQK4fZuwZhQcDuTlweLF4syVT78oANTX1xcUrDtxgvXdd4fJ1WT27GXjxk1zdIwcN058IcgQDsTZGcrKgMslWw/hWLlyRUNDtI6OmRhZuhIi/+0mHj16tH37flfXP4Qcr+h+URzl8I5qaX2gotLU14fRycieU1eHuXOJaUbR3d29cydr9uwuMaJdMAzS02HePALUIBwVFRU1tR4+v1lLS7jjd+mAYVhh4e9c7jE2O0MSOcgQDkRDQ67T6gfl0iVgMiEtTXYrKkS7CX//lSUlVomJK4Ucr7gZhP1RjsDRgwe379q1esKEbA0x8s+JgCjvqLd34IkT6b/9Js62rrgYtLTktJa9gYFBWdkVc/PPfHzIDF6nUCjOztttbE5FRe2RRA4yhIOgQN5RHBcXCA6GpUullQj8JtnZCtBuQlt7lIrKeWvrkUKOVxpDqAQ7QiqV+vnnbuXlBmw2OQr4+BDTjKK5uZPPt6DRusSYK7d+URwTE5OQkGnp6WQ6hfr64N69wMuXYydPfncb+iFAhnAQFCheRkB8PKioyK7ytfynDwKAmdlPu3bty8pKEmYwXpPFykraSkkd5XCNAoCaGvj4QIZEHi/xoVIhKIiAZhR2dueWLNHNzU0WY25qqlwbQgCYP5/g8FpRycqCMWPA2lpSOcgQDsK0aQq2IwQAGg0uXIC0NLh+XeprdXXB5cvy3m6ioQFu3IBNmywpwh1jFhaCq6tcH3kKiY0N1NcrQ6E1AJg/X6YO/wFI3oyiqgoKCxknTy43Ej3us7ISmpvBxUX81WXA1KnQ0gKPHpGmQGIiMQc0yBAOgpUV9PQoQFr9ALy9wd8fAgKAz5fuQmlpCtBu4swZWi74AwAAGblJREFUCAkRIQRU0VPpBdBoMG6cPNaBFIO5c+H6dXJaMgGAoyOoq0vkHDp6FFatEjMOOSUF5s8XoZEQKVAoMGcOpKeTs3pnJ1y+DIGBBIiS75eZPFxcFM87CgAXLkBfn5gZS8Ij/35RPh9iYmDNGhGmKEfIKM7EicoQLwMAurrg5ARZWaQpIEkzirY2iI+HjRvFnC7nB4QCSNy1p6XBtGlgaEiAKGQIB0cRjwkBQEUFEhMhMVFarl0Mw+LiLubkXJZ9toZIZGaCgQFMnizs+N5euHsXnJ2lqZMMUY54GRxyvaNLlvDi43/+9VdxTHFUFPj6wkhhQ7Veg8OBsjKYQU7BKNHw8YGKCsJyLkXi3DnCAteRIRwchQscFTB3Lnh5wbx5UikvkpJyMSwsuacnuqgom3jpxBEVBWvXijD+999h9GjQ1paaQrJFmQzhokVw6RI5nSgAIDs7rr09Oyjo0J07d0SayOPBjz/CZ5+Jue6lSzBzJqirizldlqiogJ8fXL4s63WfP4e8PMI2zcgQDo69fVdx8X9u3FDAXSFARsbffbkIR09vRF9fjYZGg76+PvHSCaKxEXJy4OOPRZiiTH5ReOUaVYJCawBgagoWFqT1hDE0NNDQeMrlcnR0dESamJICpqbinzoril8Uh5Rd+/nzMHs2YT9ekSEcnH37DvT0VMyb91l9fT3ZuoiMujqcPQunThF/UFRZ6TF1aty9e8lTpkwhWDRxxMRAYKBobWuUI4NQAF5oraqKbD0IgkTv6Pz5c3/7LdLC4nJ1tWgR+keOiN8xtKsLcnLksdD22/D3hxs3CO7X8U6IihfFQYZwcExMGGpqD/n87mGitjqWD4KCwMWF4O8Sjwf798O+fWPMzc2JlEsoGAZnz4rmFwWlM4SgXN7RBQvg4kXSVre1tdmxw/jbb0WYgtcpFbs0WlYWODqCnp6Y02WPlha4usK1a7JbsaYGKirgww8JE4gM4eD87//+81//2uXnl6unQJ/H17l2Ddhs2LKFMIH/+Q+Ym4OHB2ECpUF2NnzwATg5iTDl9Olf6uv/bWYm87Ll0kRpAkcBwMEBVFTI/HNCQoDDESFD9//+Dz7/HGg0MZdTLL8ojox37YmJsHgxqKoSJhAZwrcSEOBQUqKoVhAAtLUhIgIOH4Y//yRAGr4dFOl3MSlERcH69SKMv3fv3pdfnunufn7kSKTUlCIBBwfl2RECwLx5ZMaO0miwbZuwH/6qKsjOhhUrxFyLz4crV2DuXDGnk8XChXDliuxKPBLrFwVkCIdg7Fh48QKePSNbDwlYswYcHOCjjwgQFR8Ppqbg6UmAKOnBZgOTKdo3hMFgYFiThkauvb18F04VEaUptIZDbhIFACxZAmw25OS8e+Tx47B6tfhBHAUFYGQEo0eLOZ0sDA3Bxgby8mSxVkUFvHgBbm5EykSG8K1QKODsDLdvk62HZGRlQVUVfPedREJ4PNi3D3bsIEgnqRETAwsXivY/yNDQcNy4gp9+Sl60SC673YjL2LFQV6ckhdYA4B//gPp6Mvuh02jwzTewffs7hrW1QVwcbNok/kKK6BfFkdmPlYQECAkhuOYOMoRDMXUq3LpFthKSYWAABw7A7t0SVYxLSAATE3nfDmIYREeLHCbT3g5//DHsww+JqE4hT9BoMHaskhRaAwAqFfz94dIlMnX45BNgsyE3d6gx0dEwc6aYSfQ46emwYIH400kkIABSUqSetINhkJREfAM4ZAiHYupUhawvM4DwcLC0hNmzxZyuKNvBnBxQUxM5c+vmTXByUozMZVFB3lFiodFg61b417/eOoDHg4gI+Pxz8ZeoqIDubpg4UXwJJGJtDVpaUFIi3VVu3gR1deJfImQIh2LaNLh9W+o1rGXA1at95eXBNNrYvXsPizr3p59gxAjw8pKCWoQSFQXr1ok8i8WS952u2ChTBgUA+PlBSQmQ1Z4QZ+lSaG4GFmvwu6mpYGQkUen2tDRYsECBW6AsWCD1HyuJibB0KfFikSEciuHDYfhwMpuMEEVrawWFUsXnn963L1OkiTwe7N0L338vJb0Ig8OBq1fFKQWu3IZQaTIo4FV7witXyNQB3xS+7cRdkiR6HMU9IMSRdnvCvj5ISYGQEOIlI0P4DpTDO+rg4DBhgraq6iYud5eJiQgbhcREGDECvL2lqRwRxMbCggUi5yB3dkJZmZJ0X3oT3DWqHIXWcEj3jgLAsmVQVzfIpvDOHairk8iM1dfDo0fg7i6JdiQzdSq8eEFMvtagXLsGVlZSCakl0xCeOnVq//79Qwxgs9lscl0hJMXL8Hi8dEJ7fFEolNLSTC73Hocz1cYGJk+GVaveOri6uvr27dsAwOPBnj2waxdhaly7dq1TOs3lYmJEDpMBgMJCmDQJNDSkoNCQPHz4sLy8XNqrDB8OmppQXS3tdf5LY2NjQUGB9OTPnQtMpjjtCdlsdh5Bof00Gnz99X+/FHl5efj/qMOH4bPPgE4XX/KlS+DvDyoqRGj5dq5fv/7y5UspCadQwN9fzPaELBaLw+EMPUaM9MG6urpCYfonYCSRnZ3NYDBGjhw5xJjQ0FAbGxuZqTQoRUXY5MmyXrSrq0tdXV168uPjMXV1jMHA7twZ9G58UFAQPmz6dCLXNTMzq6qqIlIihmEYlpuL2dpifL7IE7dvx7ZtI1ydd7N79+7NmzfLYKHZs7G0NBms8zdXrlzx9fWV6hIzZojzF+Xl5U0n7qPc14dZW2MsFoZhmKur640bN2prseHDsZYWicTOno0lJxOi4FCMHz++rKxMevKvXMHc3cWZ6OzsXFhYOMSA9nZMVxdrahJNbEpKyrx58945jJwdYWdnZ3h4+M6dO0lZXSQmT4ZHj6BDqcpv/X3mP2ECuLjAsmWDe88I3w5Kj6goCAsTJ8RAiQ8IcZQscBTkwzs6YFMIAEePwooVIGKDitdob4eCAiKLZ5LFjBlw755U2hOmpsL06cBgEC8ZyHKNfv3112vWrBmtCOUTVFXB3l7qMcGyR1MTsrLg/HlISQFDQ/jtt4EDkpJAXx98fMhQThRaWuDXX2HZMpEncrlQUqJstbYHoGSBowCwaBGkp5PWnlBAaChUV/9dSKWzkxYTI2nXsytX4B//EK1linyipgYffggZGcRLJrysWn8omBQO0zs7O1tbWweuRKEYGhoCQGFh4RdffJGfn5+VlRUWFlb19m4x33zzzZEjRzReP8NhMBgmJiaE6zwEVVUhGho1Bgb5MluRz+ffvXvX0dFRBmvxeMN+/31PR8eo6dM/ptF6AIDD4bS0tPB4Rw0M8vT0Sglcq7S01M7OTpXAWrkAz59Pef7c2crqlKgTOzvNKiuX2dntI1AZIWloaOjr65NBE4/OTrPKyuV2dnulvRBOa2trQ0PD2LHSLVZXXr7d0vLMsGENwk9pa2urqamxs7MjUI2mJu+uLtPOzm90dee3t/vZ2PwoibTq6sVqas8NDYUu7C0u9+7ds7S01JDmwTib7drWZjV6dLxIs+7fv29hYaGpqTnoXQyjl5bunzBhO43WJZLYlpaWSZMmRUdHDz1MKoYwJibm2zcq1NLp9MrKyq6uLhcXl8TERAcHh6tXrw5tCAEgMjKy/fU6UUZGRjI2hAgEAoFQUKZOnar1rr22VAzhEJSVlbm4uJiamgJAZ2cnm80eOXJkfn6+sbGxLNVAIBAIBAJH1oaQx+MJvKZMJjM8PLysrExHR4dKbAlVBAKBQCCEQ4K0F7Gg0WiCVreamppUKlVxO98iEAgEQgmQ9Y6wPwLXKFkKIBAIBAJBpiFEIBAIBILL5XZ1denq6pKlgJyezEVHRwf0Y/PmzWRrJCNOnTp19OhRwdPz58/LoOzAH3/8sWLFCsHTqKioDRs2CJ4ePXo0MTFREvlfffXVjRs38Me1tbVLly7NHbqrm/Jy5MiR/h/s3bt3k60RMVy8eHHocomk0N3dHRQU9OzZM8GVU6dOnT17VhKZx44dE3wdent7Q0JCUl/VmW5vbw8KCiK9KuSg8Pn8gICAuro62S9dU1MTHBwseHrhwoU1a9a0tLT0H5OdnR0aGjpgIoZhQUFBNTU1giuxsbEnT54UZlE+n79x40ZWv5qw27Ztu3jx4tvGy6khLC0t7ejoWP6KOXPmkK2RjJg7d+7333+fn58PAHV1dRs3bvzoo4+kvaiZmVliYuLDhw/xpydPnjx//nzjq+IQx48fV5esX19WVhaeJPP06VMvLy8zMzMv+e/qJB1u3bpFp9MFH+yZM2eSrRExPHz48ObNm2RrMZC+vr7k5OSOfnWhiouLS0slSo3l8/kxMTH445KSkszMzPj4vxPmCgoKbty4MXz4cEnkSwk+n5+SktLW1ib7pV++fJmcnIw/jo+P37Bhw5o1a4Tc/CUnJ/dPSS8tLS0uLhZmIpVK9fPzCw0Nxf/kc+fOJSYmDvF1k3WwjPCMGjVqvkK3JBELExOTQ4cOrVmzpqSkZPXq1evWrXNxcZH2opqamk5OTiwWa+zYsS9fvnz27NnixYtZLFZwcHBtbe3Tp089PDwkX+XBgwd+fn7h4eFffPGF5NIUFxsbm/fwg600eHp6btu2jcvlqqmp5eTkbNiwISYmhs/nU6nU3NxcLy8viuK2E5QyP/744969e5lMpoODgwyWmz9/fmJi4ubNm7///vvw8PCkpKQhsgnldEf4PrNy5UpLS0s3N7e6uro36xJICS8vL9yNkJeX5+7u7uHhgT/Nzc11cHAYMWKEhPJv377t5eW1e/fu99wKIhSdiRMnDhs2DG/PwmKxZs2aNWbMmPv37wMAbghJ1k9e2blz5w8//JCfny8bK4gTERGRmprq5+cXGBg49FuDDKE8smrVqrt3727evFlNTU02K3p6eubk5AAAi8Xy9PT09PTEj/FYLBYh3+0TJ05MmjRp+fLlkotCIEiESqW6ubmxWKy+vr6SkhInJyf8V2NnZ2dxcTEyhIOCYdjevXt37twp4/rSBgYGc+fOLSsr27p169AjkSGUOzo6OrZs2RIYGLh3796uLtEK64mNu7s7h8P5888/c3NzPTw8zMzMenp6mpqaiPqRu3///tra2nXr1vH5fMmlIRDCgJfp6P+R4/F4krsuPT09WSxWSUnJ+PHj1dXVcUN48+bN4cOHS7vOqoJCoVCSkpL++c9/ZohYjZtCoUjy9t2+fTslJcXX1/edIYfIEModW7ZscXBwOH/+/JgxY3bs2CGbRTU0NBwdHdPT0+vq6saNGwcA7u7uSUlJT548cSeiZzaDwbh+/XpBQUFYWBiyhQjZoKGhoaur29Dw3/Lc9fX1eH1HSfD29r5582ZmZiZ+du7m5nbz5s2cnBxvb28JJSsxixYtiouLW7JkSXZ2tpBTKBSKsbGx2G8fl8tdtWrVgQMHEhIS0tPTMzMzhxiMDKF8kZOTk5ycfOrUKQA4c+ZMdHS0zILxvL29Dx8+7Obmhv/m8vDwOHDgwIQJE4iKgjM0NMRt4fr161H2qvKBYVjvK/pI75P0ilmzZh07dqy3txcA7t69y2Kx/Pz8JJQ5YcIEdXX1yMhI3BBqamqamZnFxsZ6yn1zy76+PsF7JPvv4KJFi2JiYgICAoS3hbNmzYqIiOByuQBQUVFx7dq1D4Xu2bhjxw5DQ8PVq1ePGDHi2LFj69atGyJoVk4NoaqqKrHNehSC7u7uL7/88vjx43i/KhMTk4MHD3711Vey+bcyY8aM7u5uQYSxp6dnd3e38B+7IdDW1sbfTUNDw6ysrFu3bsksCEjeUFVVVVFRIVsL4lFRUWEymdqvsLa2Jlujv/n3v//d1dVlYmJiaWn50UcfRURETJw4UUKZVCp11qxZPB7P1dUVv+Lr69vd3S3PO0IKhaKuru7s7Cx4j5KSkmSzdP+ymgEBAVFRUStXrvztzQ6og3Hw4EE6nW5qamplZeXt7X3o0KFpwnUQvXfv3rlz52JiYvCf9UFBQdOnTx8i1RVVlkEgEEpOT09PS0sLQ0rdzRESk5GRcfr06bS0tEHvyuDtk9MdIQKBQBCFqqoqsoKKiwzePmQIEQgEAkEmw4cPHz9+PIkKINcoAoFAIN5r0I4QgUAgEO81yBAiEAgE4r0GGUIEAoFAvNcgQ4hAIBCI9xpkCBEIBALxXoMMIQKBeI3KysrTp08P6CGOQCgxyBAiEIjXKC0tDQsL61/sGIFQbpAhRCCkBZfLbW5uFm9uT09PY2MjIX24urq6Ghsbe3p6Br3L5XKHuDvo4O7u7rcNwNWWn6LbCIQwIEOIQAiLu7v7unXr8McYhllZWRkaGnZ0dOBXtm/fPnbsWLxCRXR0NN6sjsFgaGtrL126VOBpXLt2rb29/YBeVE5OTp988gn++MWLFytWrNDV1TU2NtbR0QkKCnr+/Pmg+nh5ec2bN6//lY6ODlNTU0H3rsrKynnz5mlraxsbG+vp6W3atAkv5I9TVVW1cOFCHR0dY2NjdXV1V1fXhoaGhISEZcuWAcC0adP09fX19fUrKioAoL29ffXq1Xp6erhWwcHBbDZbIIrBYOzfvz88PBxX+9atW+K9wggEOWAIBEI4Pv30UxMTE/xxaWkpAKioqFy9ehW/MmXKlODgYPzxnj17Tpw4UVRUdP/+/aioKD09vcDAQPzW5cuXAeD69esCsYWFhQCQnJyMYRiXy3VycjIzM0tISLh///4vv/xiZmbm7u7O5/Pf1Ofw4cNUKrWqqkpwJTY2FgBKS0sxDGtubjYzM3NwcEhPT79///6ZM2d0dHRWrVqFj2xubjY3NzcyMoqNjS0vL2exWFu2bHny5EldXd13330HANHR0VlZWVlZWW1tbRiG+fv7q6mpHT16tKysLDo6WktLy9nZua+vD5dGp9MZDIavr29GRkZOTk5dXR1hLzoCIX2QIUQghAWvjl9RUYFh2OHDh21tbb29vTdv3oxhGJvNplKpp0+fHnRiREQEjUbr6urCMKy3t9fY2Dg0NFRwd/369fr6+t3d3RiGnT17FgCKiooEd/GGogUFBW+KbWpqUlFR2bNnj+CKt7f3lClT8Mdbt2794IMP+tskXI2mpiYMw7Zt20alUouLi98Ue/HiRcGfiXPnzh0A6L9QTEwMAKSmpuJP6XS6hYUFl8t9yyuHQMg1dHL2oQiEAuLt7U2n05lMpq2tLZPJ9PHxMTIySklJAYCcnBw+n+/j4yMYXFRUlJWVhR+Y1dfX83i8p0+f2tra0un0pUuXRkZGRkREaGpq9vT0JCcnL1myRE1NDQAyMzP19PTa2toEzUt7e3spFEp5ebmg+50ABoMxe/bs2NjYrVu3UiiUqqoqFot15MgR/O61a9fGjBlTUVGB+zYBQFVVlcfjPXjwgMFgZGVlTZo0acqUKcL84fj2Nzg4WHAlODh49erVeXl58+fPx6/4+/u/hz1EEcoBMoQIhLBoaWk5Ojoymcz169fn5+evXbvW2Nj422+/ZbPZTCbTwsJizJgx+MiVK1cmJCTMmDHDysoKN2wA8PLlS/zuihUrDh06dOHChdDQ0NTUVA6HExoait9qampqa2sLCgrqv66uru6LFy8GVSk0NDQgIKCoqMjV1TU+Pp5Op4eEhOC3nj179uzZswGi9PT08LM9Npttb28v5B9eXV0NAMbGxoIrGhoaOjo6HA5HcAXvJo1AKCLIECIQIuDj4xMZGZmfn9/R0eHp6Yk3+87NzWUymTNnzsTH/PXXX7GxsSdPngwLC8Ov/Pzzz+fOnRMIsbOzc3Z2jouLCw0NjYuLGz9+vJOTE35LR0fHyMiopqZGSH3mzJljYGAQFxfn6uqakJCAP8VvaWtr29jYMJnMQSfq6uo2NTUJuYqGhgYAsNnskSNH4le4XG5ra6uOjo5gDN4KHIFQRFDUKAIhAj4+Pi0tLYcOHXJyctLT06PRaB4eHnFxcX/++afAL1pZWQkAjo6OgllXrlwZICc0NDQ3N7eoqCgzM3PFihWC6x4eHrW1tXj4jDCoqqqGhIQkJSUxmcyHDx8KdpYA4OnpWVRU9Dab6unpeffu3cePH795S1NTEwD6Z27gXlk8zAcnIyODz+e/6a1FIBQSsg8pEQhFoqura9iwYQCwdetW/MqxY8cAgEKhNDQ04Ffq6upUVVUDAwObm5s5HM6+ffu0tbUBoLCwUCCHw+GoqamNGjWKTqfX19cLrre2tlpbW5ubm6ekpDx//pzD4RQWFn766adPnz59m0olJSUAMGrUKAaD0dPTI7heWVmpp6c3ceLE69evt7W1NTY2MpnM0NBQHo+HYVh1dbWenp6dnV1OTk57e3t1dfWJEydqamowDKupqaHT6atWrcrPz79z505nZyefz8ezKdLS0lpbW7Ozs01NTa2trfEAHwzD6HT6rl27iHqREQgZgwwhAiEa+M6PyWTiT8vLywHA3t6+/5jTp0/j7kQAmDRpEh5j2d8QYhi2ePFiAPD39x8gv6amxt/fn0r921tDpVLd3NwaGxuHUGnixIkA8MUXXwy4/vvvv/fftKmqqs6aNUuQiVFSUjJp0iTBXXNzc0EmRmRk5OjRo1VUVACgrKwMw7CGhgaB7xcApk6d+vjxY8FCyBAiFBrUoR6BkAqtra2PHz/W0tKytrYWYzqHw/nrr780NDTMzc37H8WJQX19fW1trZaWloWFhcA8C3jy5AmbzdbX17e0tHznOV9tbW19fb2BgcHo0aMlUQmBkCuQIUQgEAjEew0KlkEgEAjEew0yhAgEAoF4r0GGEIFAIBDvNcgQIhAIBOK9BhlCBAKBQLzXIEOIQCAQiPea/wfpGKZEL34G4gAAAABJRU5ErkJggg==",
"text/html": [
"\n",
"\n"
],
"image/svg+xml": [
"\n",
"\n"
]
},
"metadata": {},
"execution_count": 7
}
],
"cell_type": "code",
"source": [
"plot_bandstructure(scfres, kline_density=5)"
],
"metadata": {},
"execution_count": 7
},
{
"cell_type": "markdown",
"source": [
"or get the cartesian forces (in Hartree / Bohr)"
],
"metadata": {}
},
{
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"β Warning: Forces for shifted k-Grids not tested\n",
"β @ DFTK /home/runner/work/DFTK.jl/DFTK.jl/src/terms/terms.jl:72\n"
]
},
{
"output_type": "execute_result",
"data": {
"text/plain": "2-element Vector{StaticArrays.SVector{3, Float64}}:\n [0.00044880407371949554, 0.00044880481864267277, 0.0004488052606748624]\n [-0.00044880481864318646, -0.00044880407371887104, -0.00044880526067524775]"
},
"metadata": {},
"execution_count": 8
}
],
"cell_type": "code",
"source": [
"compute_forces_cart(scfres)[1] # Select silicon forces"
],
"metadata": {},
"execution_count": 8
},
{
"cell_type": "markdown",
"source": [
"The `[1]` extracts the forces for the first kind of atoms,\n",
"i.e. `Si` (silicon) in the setup of the `atoms` list of step 1 above.\n",
"As expected, they are almost zero in this highly symmetric configuration."
],
"metadata": {}
},
{
"cell_type": "markdown",
"source": [
"## Where to go from here\n",
"Take a look at the\n",
"[example index](https://juliamolsim.github.io/DFTK.jl/dev/#example-index-1)\n",
"to continue exploring DFTK."
],
"metadata": {}
}
],
"nbformat_minor": 3,
"metadata": {
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.6.1"
},
"kernelspec": {
"name": "julia-1.6",
"display_name": "Julia 1.6.1",
"language": "julia"
}
},
"nbformat": 4
}