{
"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](https://juliamolsim.github.io/DFTK.jl/dev/#density-functional-theory)\n",
"chapter for some introductory lectures and 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": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"n Energy Eₙ-Eₙ₋₁ ρout-ρin Diag\n",
"--- --------------- --------- -------- ----\n",
" 1 -7.905260123614 NaN 1.95e-01 4.2 \n",
" 2 -7.909881622689 -4.62e-03 2.99e-02 1.0 \n",
" 3 -7.910078962844 -1.97e-04 2.98e-03 1.3 \n",
" 4 -7.910108570943 -2.96e-05 1.24e-03 2.4 \n",
" 5 -7.910109230188 -6.59e-07 5.88e-04 1.0 \n",
" 6 -7.910109405857 -1.76e-07 1.97e-05 1.0 \n",
" 7 -7.910109414090 -8.23e-09 6.80e-06 2.9 \n"
]
}
],
"cell_type": "code",
"source": [
"using DFTK\n",
"using Plots\n",
"\n",
"# 1. Define lattice and atomic positions\n",
"a = 10.26 # Silicon lattice constant in Bohr\n",
"lattice = a / 2 * [[0 1 1.]; # Silicon lattice vectors\n",
" [1 0 1.]; # specified column by column\n",
" [1 1 0.]]\n",
"\n",
"# 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 in Hartree\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": 1
},
{
"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.0807686 \n AtomicLocal -2.1784368\n AtomicNonlocal 1.7340995 \n Ewald -8.4004648\n PspCorrection -0.2948928\n Hartree 0.5413872 \n Xc -2.3925705\n\n total -7.910109414090\n"
},
"metadata": {},
"execution_count": 2
}
],
"cell_type": "code",
"source": [
"scfres.energies"
],
"metadata": {},
"execution_count": 2
},
{
"cell_type": "markdown",
"source": [
"Eigenvalues:"
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "7×10 Array{Float64,2}:\n -0.170032 -0.131619 -0.0881013 … -0.0560199 -0.114743 -0.0697767\n 0.20171 0.0912103 0.0125268 0.011408 0.0423192 0.017919\n 0.249665 0.175087 0.176467 0.133227 0.220467 0.112584\n 0.249665 0.231789 0.202704 0.16134 0.220467 0.190779\n 0.35146 0.360507 0.340562 0.292133 0.321217 0.32774\n 0.37038 0.396309 0.389909 … 0.33217 0.388572 0.460941\n 0.37038 0.402093 0.412871 0.56613 0.388572 0.463185"
},
"metadata": {},
"execution_count": 3
}
],
"cell_type": "code",
"source": [
"hcat(scfres.eigenvalues...)"
],
"metadata": {},
"execution_count": 3
},
{
"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 Array{Float64,2}:\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": 4
}
],
"cell_type": "code",
"source": [
"hcat(scfres.occupation...)"
],
"metadata": {},
"execution_count": 4
},
{
"cell_type": "markdown",
"source": [
"And density:"
],
"metadata": {}
},
{
"outputs": [
{
"output_type": "execute_result",
"data": {
"text/plain": "Plot{Plots.GRBackend() n=1}",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdZ2BUVdoH8HPunfTeJxUIvYQBkkASkqAQEvXFgoJlVVhEURcVBBFQVFAUCyqirK5rW1BBEV3FRocEEkoooYcWSsqk9zqZc94PYWNEShJm5rb/75NMxskzN/e5z72nUs45AQAA0CpB6gAAAACkhEIIAACahkIIAACahkIIAACahkIIAACahkIIAACahkIIAACahkIIAACahkIIAACahkIIAACaZslCmJ2d3dDQ0P73m81mC/52wPG0IM45Vh+0IMaY1CGoCpLdsixZCO++++7s7Oz2v7+urs6Cvx1wPC2oubm5qalJ6ijUAyenZeF4WhaaRgEAQNNQCAEAQNNQCAEAQNNQCAEAQNNQCAEAQNNQCAEAQNNQCAEAQNNQCAEU7OzZs6WlpVJHAaBsKIQASvXftT8bbv17r5hRZ86ckToWAAVDIQRQnl1FfPpO88Qfc2v8+lU4BoxaXbxgHztRiTXhADpDJ3UAANBeRyv4t2fYytO8iZF7wmn66w9lrF2p9x/TLW7oilPshl/MXvZ0fDid2FPo5kalDhZAMVAIAeTuQi3/PoevzmHnasidXemnCWK8vqXOif0nTWx5z+vR4mtRYnohX53DYn9q7uFOx3cT7usu+DtJGDiAMqAQAshUWSP5+TxbcYrtL+G3hAqzDcItoYJ45Sc9gZJ4PY3Xi28OFdfnsdVn+Ev7TNF+9MEewtiugpudDUMHUBQUQgB5qWwiP55jq3NYagFPCRGe6i+khAj2HenNdxDJrWHCrWGkvlncmM9WnORPppsSA+mEnsLtXTr2UQBaoFu8ePHgwYNHjRrV8u/i4uK1a9dWV1dHRUUNHz6cEHLw4MG2myuNHTtWp0P5BLCwBjPZkMdWn+G/XGBxAXR8N2HljYLr9T3GOekuVsTyRnHtebb8JHt8u/mWUGF8OL05RNChIgIQQgjRlZeXP/bYY/fff//8+fNPnjwZExMzduzY0NDQ+++//8EHH3zllVe++uqrn376KSIiouV/uPXWW1EIASzFzElGIV9xiq3OYf086fhuwjsxdr6OFv4tXg5kQk9hQk8hr5Z/l8PfyGIPpZrv7Co82EMYrqcYVwMaRznnx44di46OzsvL++CDDzIyMn7++WdCyIYNGx544IHCwsLZs2cLgrBo0aJrfpbBYFi+fLnBYGjn766urnZzc7uu8KENHE8LMplMjDEHBwfr/Yq9JXz5SfbtGebtQFuqVKCz9X7bpc7V8FWn+RcnWV0zuaMLndhTGOJrxYJYU1Pj6upqvc/XGiS7ZekIIX379g0KCtq+fXtISEhRURFjTBCE/Pz80NDQljfl5OR8+eWXYWFh8fHxgoD2FIAOY4yt+HqVo4P9gKS7VuewL09xB4GMD6epY3Q9PSR4JOviSmcb6GyDcKScr85h4zeZHUUyPpyOC25M/fHrXuHdkpJG2T4qAElcbOQMDAzMy8t7+OGHs7Ky+vTpExgYWFJS8t///pcQ4uTkVF1dvWnTprS0ND8/v40bN7q4uFz2sxobG5ctWxYQEND2xZSUlOjo6Cu9397e3qJfR9NwPC2o5YnQgh/4/fffP/rx76b6Wr/9zn+/PenbEby/58UfNTZa8Pd0WA9nMrc/md2PpBfTb3PIsCffra+pdD//n/RvA7t3726p39LY2Ghnh3GrFoNkbz97e3t6reb/i4WQc04p3bFjx3fffffKK68EBwe/9957r7zyyvLly+fPn9/ynqamptjY2KVLl86dO9eqcQOoj7tfkCn/pKtg+mmsn2Gg7JaAESiJ9+fx/sRwInjW+1vrqqvcPTykDgrARijnnBDSq1evd95557PPPouIiFiwYAEhpKyszMfH5/z5860NpISQuXPn5ubmrlix4rKfhT5CaeF4WpDF+whn7DQXFuS9E+dwSZOJDJ08dfq+PZ5PR3vd38Ni/SDoI7QsJLtlCYSQEydO5Ofnx8fHOzs7l5WVtfyg5T9cXFxaKiUhhHO+Y8cOC7aWAGjE8Qq+/CR7+6Yw+VdBQkjPHt3fH+U9Zw+rMUkdCoBN6F588cVVq1bNmDHD09PziSeeSElJaW5uDgwM/OqrryZPnuzt7T106NAhQ4a4u7unpqaWl5c/9dRTUscMoDAzdpnnDRb1ylntLNaf3hBIF2WZX40SpY4FwOp0rq6u//znP5OSkgghMTExR44c+f3336urqz/88MORI0cSQpYtW5aZmVlXV/fss8+OGTMGPbQAHbL2PDtbTab2U9hw67eGiRFrTH/vKUgyqBXAlmhry+f1Qx+htHA8LchSfYRNjESsaX4vVrwpRHnlZFEW213EfxhtgYdC9BFaFpLdshR2lwqgLO8eYn09qRKrICFkZoRwtIKvy5XdGFcAy0IhBLCWwnry9iHzW0OVmmX2AnlrqPD0TrPJkjMqAWRHqSkKIH/P7jY/3FvZfWy3dRG6uJIPj6ESgpqhEAJYxd4SvjGPzxmk+FGX78SIrx4wFzdIHQeA1aAQAlgeJ2RquvmNoYK78pcV6+tJ/9ZdeCHTLHUgANYiWSEsLi7+1yeftt3pEEA1/nOCCYRYcGUWaS2IFNee55klGDUD6iRZoo598JHZ+3Sx/3d3fbNUIQBYRbWJPJ/J3osVFdw3+GfudmT+EGHmTjMqIaiSZIUwUB/gemFPk51r0Nem+7eYfzrHGtH0AqqwcL85JYRG+6mmDhJCyOTeQrWJfHsGo2ZAhSTba/6bzz5MT0+Pjo6uJXY/n2f/zmYPpZpvCRXGh9OUEMFeJU1KoDmnq/hnJ9jBO5XfN/hnAiUfxIn3bjaPCRNcJLtsAFiFZAVHEASDweDg4ODtQCb0FNYm67LvtksKpkuPsMCvTBO2mteeZ5i9BIozLcM8xyDacq95m4kLoPF6+mYWmm5AbWT05OXjQCb0FDbcrMu6UxfpS9/I+qMiNqMighJsyOPZleQJpS0r2n5vDRX+eYydrUZfIaiKHDM2xIVOGyBsv1W3b+z/KuLXpke3m7cbLbcuKoClmRh5Mt38XqzooPipg1cU7EKn9Rdn7cadKaiKHAthqzDXixVxz+26fp50zh5z2MrmaRmoiCBHS4+wHu7kllBVjZH5q2cGCntL+IY8pCCoh6wLYauubhcr4rqbRS8HMmW7uesqVESQkaJ68kaW+e0Y9T4M/o+jSBYPE57eaUaHBaiGMgphq36edP4Q8eg43a8popcDmZxmDv/mYkWUOjTQtLl7zH/vJfRW8rKi7XdnVyHYmXx8HJUQVEJhhbBVfy86f4iYPV73c7Lo5UAmpZr7f9ec+OQbfn0ip86aJ3V0oC37Svhvuex55S8r2n7vxooL9ptLsAAp2BBjbORtdwf3j/ruh58s+8lKLYStWiriibt1nySKWRvWlDz03bc//ix1UKAhnJAnM8yvRoke9lKHYkP9POnd3YT5+zCVAmynvLw881R+/s2vrfjuR8t+suILYQtKSKw//fLdV7r/8rTHnfMYGkrBVr46xRrNZGJPlaRS+70cKX6Xww6WIdnARgRXHz749sGH//3aczMs/MmW/Thp3Xpzyslt//WPu+OLk+i9AFuoMZE5e9iSGFHQROfgn3g5kJeGiFPTsQAp2MgLe80Tps7Y9/vq/v37W/aTVVUICSGUkGXDxbl7zKWNUocCGvDaAfOoIBqv114ZJIQQMqWPUN1Evs/BfSdY3aEy/l0OeznSKj3xaiuEhJDBPvSursJLe9F7AdZ1ppr/O5stilZhErWTSMkHceKMXawOe8iANXFCnkg3vxol+jhY5fPVmcOvRonf5bA9xWizASuasZM9EyEGOWv0cbBFvJ4O86OLD+GhEKzoixOsrplM6mWtgqXOQujlQF6NEqemmzFqBqxkUz4/XM6nD1BnBnXI2zHC+0fM52qQbGAVVSYyL5MtG27FnnjVpvGkXoJIyXKMmgEraGZkeob53RhBxcuKtl+oC53aT5iNBUjBOp7fY76tCx1qzQ0+VVsIBUqWxYnPZZrLMWoGLO2Do0zvTG4NU236dNTsgeKuYr6tAA+FYGGHyvi3OewV64yRaaXmTB7iS2/vIryIUTNgUaWNZFGW+V0NLCvafk468ka0MH2nGXMpwIJax8j4Olr3F6m5EBJCXosSv8thB0qRnWAxc/eY7+8uDPDS9BiZv7o7XPC0J59ko4EULOY/J1hdM3nIamNkWqm8EHo5kJcjMecXLGZ/KV97jr04BI+Dl/FBnPjiXkzhBcuoMpF5e607RqaVygshIWRyb4FxjJoBy5ieYX4lSvTU0rKi7dffi97ZVXgZC5CCJczLNI8Jte4YmVbqL4QCJR/EiXP3mCuapA4FFG7laVZtsuJkJhV4NUr85gw7hAVI4focLuffnmELo2zU9KKJlI70pbeFYdQMXJf6ZjJ3D1sSK4roHLwybwfy/CBx+k7kGnReyxiZVyKtPkamlSYKISHktWhx9RmMmoHOey3LPDyAJmp1WdH2+0c/obiB/PccOiOgk1acZDUmMrm37cqTVgqhtwNZECk+gVEz0Cnna/hHxzS9rGj7iZQsiRFn7mQNeCyEjqsykecy2bI4m+7ooqHEfri30MzJCoyagY57eiebPkAMc8XjYLuMDKIGH/oOFiCFjnsh0/x/oXSYv01zTUOFsGWtmTkYNQMdtDmf7y/lMyM0lCzX7+1hwtuHzOexACl0xJFyvuoMe9VWY2RaaSu3I33pmDBhPoZ3Q7uZOXl6p/ntYYIjpg52RDc3+nhf4blMPBRCe7WMkXnZhmNkWmmrEBJCXo8WvznNsjC8G9rnw2PM15GM7aq5TLl+zw0Stxt5mhG5Bu3y5SlWbSIP23CMTCvNpbe3A5kfKU7dgVEzcG1ljWThfvMSLCvaKc468lq0MC0Du6HBtVWZyNw97IM4aaYnaa4QEkIe6S2YGPnyFBpt4BrmZZrvDhcivDFGppPu6y642ZHPTyDX4Bpe3Gu+JZTG2HaMTCstFkKBkmXDxdm7MWoGruZoBVlzlr2EZUWvAyVkSaw4L9NcacLNBFzRkXK+8rQEY2RaabEQEkKifOn/hQoLMGoGrmDlN6vveW3FSwbi4yB1KAo32Icm+1TfM/9fW7dtkzoWkKOWMTILhoh+Nh8j00qjhZAQsiha/Po0FkWEy9i4cePjH6w5vT/DLWu11LGoAVu7aEdu/Z2PPlNcXCx1LCA7X59i1SbySB8pi5F2C6GvI1kwBGvNwGUEBQWZi3Kci4526xImdSxqYOjdzfF0Gms2ubi4SB0LyEuViTy7my2TaIxMK52Uv1xqU/oIn59gX59i9/fQ7g0B/FXffv08Zv/6Q2JTXF8UQgt45ql/6PqO+LU+zNnZWepYQF7m7zXfbPN1ZP5K0wWgZdTMs7tZJUbNQBsHy7i7l09keIDUgajH+KiuO8sdG9EpD20cKedfnWaLoqUfjya4ubmNHTu2oqKCEGI0Gm+77TZvb29fX9/x48eXl5cTQvLz85OSklxcXMLCwtasWSN1wBYW5UtvDqUYNQNtbczjowKlDkJdPOx4X0+aUYSOCPjDk+nm+ZKOkWkllJSUUErnzZtHCJk9e7a9vb3RaLxw4UJlZeXLL79MCJk2bVq3bt0qKytXrFjx97//vaioSOqYLeyNoeLXp9lBjJqB/9mUz0aiEFpaUjDdlI8JhXDRV6dYpYlMkXSMTCvBwcFhzpw5K1asaG5uzs3NvfHGG+3t7Z2cnEaMGHHhwoXy8vIff/zx+eef1+l0I0aMiIuL+/rrr6WO2cJ8HMiLg8WpGDUDhBBCmhjZYeQj9FLHoTpJwcLGPCQZEEJItYnMkW4dmb/SEUIGDBhQVVVVWFg4ffr05557LiAgwGQyffnllx999FFOTo6jo2PXrl1b3t2/f/+TJ09e6bM45/X19bW1tW1fdHBw0OnkPiTnsb7Cf06ylafZ37rL4vYEJJRRyPt4Um8HzvD0YlHDA+jRcl7eSLwwNVPz5u8zJwfTWKnHyLTSEUKcnZ1FUSwrK+vTp4+Hh8eyZctMJlNQUFB4ePiJEyfajnh2d3c/fvz4lT6rqqrqxhtvpPRP323hwoWPPPLIZd9fU1NjoW9hAW8OEu7bbj/Cq87dTqk3rbI6nsr161ldoh+prq5njDU1YRiVZdTW1rpwHu1j99uZpltDcItxvRSd7Mer6IqT9jtTmqqrbXGxbSlwV3+PjhBSXV1tNpv9/PxuvfXWBx54YNq0aYSQBQsWTJ48+a233qqqqmp9d0VFhb+//5U+y8PD48cffzQYDO0P0c3Nrf1vtqoRbuTmc+Z3TzotHib9EKZOk8/xVK7U4uZXo0Q3Nx1jzMEBDy+WQSl1dXVNCWPp5fxvfRWcYvKh3GSfndo8f4gQ7mcvdSB/EAghBw8e9Pb29vPzy87Ojo6ObvnB0KFDs7Ozu3XrZjabW5tDs7KyevfuLVmwVvbmMPHLUxg1o2mVTeRwOZdPi43KJAVTdBNq3NenWUUTebSvvDqhhNLS0gULFjz00EOiKA4fPnzp0qUVFRUlJSXLli2Lj493d3e/++67X3jhhYqKiu+//37//v333Xef1DFbi48DeWEw1prRtC0FLM6fOsm9U1upBvnQskZ+DtvWa1W1iTy7my2JlcsYmVbCgAEDunTpsmDBAkLIp59+KoqiwWCIjo729/dfunQpIWTJkiWEkN69e7/66qtr1qzx9vaWOGRreryv0GAmq06jD0OjNuXxUcHyuldVE0rIyCBhcz4KoUYt2GdODqaJepmVQUJ0BQUFrf8ICgr66quvLnmHt7f3qlWrbBuVZARK3osV79rYfEuo4CGjFmywkY35/OsbUQitKCmYbsrnk3pJHQfY3NEK/sUJduguO6kDuQzk/KVi/WlKiLBwP9aa0Zy8Wl7SwA3YhteaRgfTjXnYsl6Lnkw3z48UA2W53CwK4WW8ES1+cRI7NGnO+jw+KkgQUAetqYsrdbOjSC6tWXWaFTeQx+SxjsxfyTQsafk7kZeGYK0ZzdmUz0cFoQxaHcaOak21iczazT6IE3VyLThyjUtqj/cVqk3k2zMYNaMVnJDN+SwpGIXQ6kYFYdFRbXl5nzlJlmNkWqEQXp5IybI4ccZOVmWSOhSwicNl3FlHu7nJN1dVY1SQsN3IsSWTRhyr4J+fYK9FyXoVBRTCK4oLoKODKUbNaMTGfJ6EdlGb8HIgvTzormK0jmrCE+nml4bIdIxMKxTCq3lzqPjZgbJlq35W9Mp+0B4b89AuajstY0eljgKsy2g0zv3818Kqxsdlto7MX8k9Pmn5OxG7D8dNW7kr6Y57pY4FrMjEyHYjvzEI6WAj2JJJ9RhjkTfc9MaqjQE/zZTtGJlWsg9Qal72hAm6hmY0kKrZziLe04P6YIVtW4kLoIfLeSX29lC1+mYm6Oy87RVwx4NCeA3bf/shYciAyW+vlDoQsKKNeWw02kVtyFEkMf50awFaR1VLEIQxi3+eeFPsio/ekzqWa0MhvAZvb+/J945Nq1DqjifQHhvz+Si0i9rWqCBhExYdVbUdjUEzH7zd0dFR6kCuDcl/bckhwqZ8htn1alVtIofK+PAAPBHa1GhMq1e1k5W8kZF+nspIKxTCa9M7kVAXugejvVVqSz6LwdZLNjfIh5Y08Au1SCt1Wp/HU5TT3YBC2C4pIXRdLjJWndAuKgmBkhuDhE14KFSpdbk8JQSFUF1SQoR1uejYV6eNeRwzCCUxKoiim1CVmhhJNTIFbe2pmEClFa+nRyt4WaPUcYCl5dfxono+2AeFUALJwXQDtmRSo/RC3sdTSfORUAjbxV4gwwPoFqwUrDob8vhIbL0kka5u1EVHj5SjFKrN+lymoA5CgkLYfikhwjr0Z6jOJrSLSgpbMqnSujyeHKKk4qKkWKWVEkJ/u4CMVZvNBSiEUkoKxpZMalPSQM5U8aF+SkorFML26u1B7QRyvAK1UD2OlHM7gYRj6yXpjAwS0ozchFKoIuty2cggwU5RtUVRwUotKZiuRzOOimzI48l4HJSUjwPp4U53FiGt1ENZEydaoBB2QEowxSQKNdmUz0ZhD0KpoXVUTTghG/OVt3IvCmEHjArGztrq0Yytl+RhVBC2ZFKPg2XczY52U1p3A64CHeBpTwZ40+2FSFo12FnEw92pnwIWBFa5BD09WMarTFLHAZagxHZRgkLYUcnBwnq0jqqCEhtwVMlRJMP86TZsyaQK63JZiqImTrRQXsTSwqKjqrExD0uMygVaR9WhrpnsKeYj9Mq7v8SFoGOi/WhuLS+okzoOuD7VJnKwjMdj6yV5wLR6ddhawKN8qaud1HF0HAphx4iUjAwSNuShGUfZthawoX7YekkuhvjQogaeiy2ZFE6h7aIEhbAT0DqqApvyuYKWxlc9gZIbAoXN2IlC4RQ6UoagEHZCSghdn8ewZr6ibczjSZhBKCfYkknpztXwiiZuUOZGLiiEHRbiQv0c6YFSJK1SGetJfh0f4qvIjFWr0cF0Yx5HUinXulyeHKzUfVxQCDsjJYRiJwrl2pDHRgUJokJTVqW6u1MHkRzDWr6KtT6PJyuzXZSgEHYONqxXtI15HCuryVBSMN2A3ndlMnOyJZ8pdz6SUuOW1gg93V+CtTCUanM+tl6SI3QTKteuIt7FlQY6Sx1HZ6EQdoaTjgzzp1uxUrACHa3gIiU93FEIZScpWEg1MmzJpETr85hCx4u2QCHspOQQAVsyKdHGPI6V1eTJx4GEu9HdxUgr5VmXq7At6S+h4NCllRKM2YSKtAntojKGJWaUqLyRHC3ncf4KTisUwk4a4E3rmsnpKiStkjQzklrARiq2S1/1koIF7E2oOBvzWWIgdRCljuM64IrQSZSQ5BBsWK8wu4t5NzdsvSRfiXp6oBTD0BRmXS5X6MpqrZQdvbSSg+l6tI4qyoY8tIvKmqNIov1oagHSSkk25PEUhacVCmHnjQ4WthZgkJuSbFLyVCeNGBWE1lElOVbBBUp6eqAQapWvI+nhTncW4e5VGWqbyf5SHq/AzdI0BeNllGVdruIfBwkK4XVKCaFYYkYpthbwaF/qgq2X5C3SlxbUYctPxViXq+wZhC1QCK9LCmYTKsemPJaErZdkT6TkBrSOKkSDmaQXchUMw1b8F5BWrD89WcmLG6SOA9oBI2WUYlQQWkeVIc3IB/pQD3up47huujlz5gQHB0+ePNnZ2fnEiRNbt25t++Nx48adPHkyKyur9ZVJkybZ2dnZOky50glkRKCwKY/d2x23FLJWWE/y6ngktl5SgqQgunA/nggVYH0uS1ZFK4tOr9evX7/+m2++SUtLKy0t3bt3b8sPzp49u3379nvuuef7779fv3790KFDW16fOHGidNHKUcuWTPd2lzoOuKoNeezGQGy9pAw9Pai9SI5V8L6e+IPJ2ro8/kmCKgrh9OnTp06d2q1bt61bt954442xsbEtP5g6dWpAQICHhwch5Kabblq0aJGkccrXTSH05X1mTnCNlbVN+dh6SUlaWkdRCOXMWE/yalXSyiIQQuzs7EaOHLlly5bWVxsaGlauXDl58uSWfx4+fHjx4sWrV69uaEBv2KW6uVFXO3q4DF0asrYJa20rShK2ZJK93y+w0cEqaWW5OJbc39+/oKCg9dU1a9Z4enomJia2/KikpMRoNK5evXrevHkZGRne3t6X/az6+vqXX37Zx8en7Yt33XXXiBEjLvv+hoYGdXQ3jgoUfj3HejpL3KuhmuNpcdmVhBAx1KGx/TdyJpOJMcY5rsWW0dDQoNN1YOZKgi99bIdQXWeyU0PDm+XJIdl/Oy+MDOQNDXJfEM/e3l4QrnEaXTw1m5ub2x7Wzz777OGHH6aUEkJmzpzZ8iJjbMSIEUuXLp0/f/5lP0sUxd69ewcHB7d9Ua/XX+kPZmdnJ/nf0iJuCuHLjvFnBkr8XVRzPC1uWzEfHUw6enAYYzieltLRk1NvR7q6soNVdjF+1gtKwSRPdsbJ1kL2xlDBzk7uj4QthezqLhbC/Pz8gQMHtvx3Tk5OamrqF198cclbBUGIj48/ffr0lT7L3t7+nnvuMRgM7YxPFEVRVPKK5f8zOpRMTDU1ctFZ0snaqjmeFre5wHxPOBXFDjxcMMYopTieltKJk3N0MNlUQIfr8Uh4GZIn+/4S7ufIu7qrZH0KgRBSUVGxadOmMWPGtLz0+eefJycnh4aGtvyzsbGx5T+ampo2bNjQr18/SQKVMxcdGexLU41oRpMjMyepRnaj8uf8as2oYGFTHiZRyNS6XK6CBWVaCVOmTElISLjlllsGDRpECGGM/ec//2kdJkMIGThw4B133DFp0qR+/frpdLonn3xSumjlKzlYWI+11mRpdzEPc6F6J6njgA5K1NP9pbxG7j1QGrU+lyl966W2dPHx8X/7299ax7M0Nzd///33rc2khJD169fv27evqqrqkUceiY2NbU97qwalhNAJW1EI5WgjFpRRJmcdifSlqUZ+Syj+fPJSbSIHSnlCgHr+LroJEya0/be9vX1kZGTbV7p06dKlSxfbRqU8g31oSSM/X8PDXNVzcqjDpnw2x4CuPkVq2bD+llD8+eRlcz6L8adOKukfJARrjVqKQElSkLABCyTKTG0z2VfCE7D1kjIlBVHklAypYEv6S6jqy0grJYRiJwq52VbAI7H1kmJF+dG8WmzJJDvr81Q1UoagEFpQcoiwMY+ZUQrlBFvSK5pIyYhAYTO2ZJKTU1W8wUz6eaEQwuXonUiIC80sRiWUEYyUUbpRWGtNZtbn8pRgtY2ZRCG0pJadKKSOAi4qqifna3iUKhYF1qykYPQ4yMu6PJ6srnZRgkJoWSkhwjrMJpSNjfnshkBBh3NcyXp7UDuBZFeiFsqCiZHUAqaCLekvobbvI60EPT1SzssbpY4DCCGEbEK7qCqMDMSG9XKRXsh7eVA/R6njsDQUQkuyF8jwALqlAA+FsrC5AIVQDUYFoxGj0UwAACAASURBVBDKxfo8prLxoi1QCC0sOVhYl4uklV52JTcx0ttDhUmrNUlBwtYC1ozbSxlYl8uTg1VYNVT4laSVEkJRCOVgYx5PxuOgKvg7kTBXmlmCtJJYSQM5XcWH+aswrVAILayPJxUo+valtymfjwpSYcZqUxJaR2VgfR67MUhQ5VbJavxOUksKpuvxUCgpMyfb1Di2TbOSgoRNmFYvtfXq2nqpLVwpLC85mGIShbQyi3mICw10ljoOsJARgXRvCa9tljoODeOEbMjjo1Xa3YBCaHlJwUKakTeapY5DwzbmY7yoqjjryBBfmoa9r6VzqIy72JFwN3WmFQqh5Xnak/5edEchklYyG/OwxKjajAoSNmLDeumsy+Up6r25xMXCKlJChPVIWonUNZO9JTwRWy+pC8bLSGudurakv4Rqv5i0kjGJQjqpRj7Yh7raSR0HWNRQP3q+lhvrpY5Dk+qaye5ifkOgam8uUQitYqgfvVCDfdSksSmPJalxzq/GiZQk6oUtGDsqhW0FPMpXzTeXuF5YhUjJjUHCRiStFDBSRq2SgrElkzTW5bFk9baLEhRC68ESM5IoaSA51dh6SZ2SgugGdBNKYZ16ZxC2QCG0lptC6PpcxpC2trUhj90QqM7FL6CPJ+WcnMCyTbaVW8tLG7jBG4UQOi7Ehfo60qwyJK1NYWU1dcNOFLb3ey5PCREEVWcVCqEVYeyo7W1CB6GqjQpCN6GtrctV4Zb0l0AhtKKUYGxYb1MnK7mJkb6eKk9aLRsdLGwpYGaUQlsxc7IlnyWpfXkKlX89ad0QSPeX8hqT1HFoxsZ81a6FCC0CnEiIC92LLZlsZXcxD3VV/7K9KIRW5KQj0X50awGS1kbQQagFSUHoJrSddbnq3JL+EiiE1pWM1lFbMXOyNR9LjKrfqGAsOmo761W6Jf0l1P8NpZUSQtfh7tUm9pbwQGf1t+HADYF0D7ZksomKJnKknA8PwBMhXJ8Ib1prImeqUQutbmMexotqgouODPah27Elk/VtzGMJeuogSh2H9aEQWhclZDQ2rLeJTflYYlQrkoKxYb0trM/jKt5xoi1NfElppYTQ9WgdtbIGM9lTzBOw9ZI2YLyMbaxX+8pqrVAIrW50sLAln5lw/2pNqUY+yIe6q3d1fGhrqB/NqebFDVLHoWrHKzghpJcHCiFYgq8j6e5OdxXhBtaKsCW9pugEkhgobEbrqDWpfqHttnDtsIWUELoOA76tCSNltGYUWketbF2eJmYQtkAhtIWUEAHjZayEc77w3Q+PLH+5j0O11LGA7fRj51e9Ofv7H9dKHYg6NTGSXshvDNRKgdDK95RWrD89UclL0KVhBbt37379m40mIn78yadSxwK2s/jlF2qCoyY/80J9fb3UsahQmpEP8KJeDlLHYSsohLZgJ5ARgRjwbRXdunWjpWddjv4yLHKw1LGA7YyMj3XY+oGbh5ejo6PUsajQulyVb0l/CZ3UAWhFcjBdl8vvCZc6DtXx9/fv9fqOtyKbR3Z1lToWsJ1np/2Dxtx3yuROqVb6sWxpXS7/OF5DhVBDX1VaN4fS9Xkc/YQWV9tMsqvF4aGogpozqrvnjkKklOUZ60luLY/y09AdBgqhjXRzo44iOVKOvLWwXUV8kI8mVoGCSxh86IUaXtYodRyqsy6XJQULoobqIAqhDaWEYK01y0sv5HEaWBQY/kqkJMoPM3QtTzsLyrRCIbSd5GCKLZksLr2IoRBq1vAAml6EnLIkxsnGfKa1Da5RCG1nZJCws4jXYfsYy+GE7CriMf44jTUq1l9IRzehRe0v5b4ONNQFhRCsw82ODPKhadg+xnKOlnNvB6p3kjoOkEhsAN1TzJvxTGg5mlpZrZVQUlLC2LXPo5KSksZG9Epfr5QQYT3WWrMcdBBqnKc9CXOlB8twc2kx6/OYRrZeaksYPHhwjx499u7dSwiZOnWqdxsRERGEkFOnTg0cONBgMOj1+o8//ljqgJUtOYSuw3gZy0kvQiHUuuEBNB3jZSyk2kT2lWhxOzPhwoULM2bMeOihhwghy5YtK/uf5OTk0aNHE0KmTZt288035+XlpaWlzZw588KFC1LHrGBDfGhxA79Qi7y1jPRCHuevuaSFtmL9KboJLWVLPovxp87aW2dFIIQ88sgjp0+fzsrKan21rKzsxx9/nDhxYnFx8bp166ZNm0YIGTBgwA033LBy5UrJglU+gZJRQcIGrJpvCSUNpLCe9/NCIdS04QEohBazTjNb0l9CIIQ4ODh07do1Jyen9dUvv/wyIiLCYDCcO3fOxcUlKCio5fXevXu3fdslGGO5ubln/qyqqsra30FZMJvQUtILWaw/1dS0X/irHh603sxz0cpiCRqcQdji4jOwq6trRUVF66uff/75o48+SgiprKx0dnZufd3FxSUvL+9Kn1VbWztlyhQ7uz9tEz537tz777//Su/X4DqB8V706Qy7yup6i1/BtXY8t+aKkV60psYqm3qYTCbGmMlkssaHa1Btba31Pjza227zufo7wzQ0DM0ayX6ultaY7LrYNdTUWPaDJebs7CwI13jMvVgIy8vLAwICWv5779692dnZ9957LyHE39+/srKy9d0VFRX+/v5X+iw3N7cffvjBYDC0Mz7Ouaur5paIdHUloa7N2Q0uQy29lJ/WjmdmefNLQ0RXV6tsPtBSCB0cNLMPjfVZ7+RMDGb7q3QTXDW0zp41kj31Ars5lLtp6RrSSiCElJWVnT17tn///i0vffbZZ+PGjfP09CSEdOvWTRCEw4cPt/xoz549LUNJ4XokuJS+v+pX7KN2PUyM7C/l0VpaFxiuJM6fYvXt61RaWrr8+19HeGv0oiScOnXq6aefTklJCQsLI4Q0NDSsXLmyZRApIcTV1XXChAmzZs06ceLERx99dOLEiXvuuUfSgNVgzYwxK39PHfvgw1IHomD7SngPd+pud+13gupF+dGj5Viz6boMHfV/uzJ2fPrcZKkDkYYwduxYURQ///zzln8fO3Zs7NixI0aMaH3HW2+91bt377vuuuvHH39ct26dm5ubRKGqhz3lZtGh0WSWOhAFwwxCaOUokghvuqcYD4Wd12AyU529yDXUz9qW7tChQ23/PXjw4E8//bTtKy4uLkuWLLFtVCq3a+PPo97fPvPvz0odiILtKOR3dkUhhItaptWPCMQp0UkPvf/TgT3pK6fOkjoQaWhxyojk9Hr92LFj91RhiczO21mEqfTwh1h/ml6o0acZi9hv8pt0z1hNjbZrC4VQGgl6IRWrb3dWTjXnnHR1QyGEi+L1QkYhR0Z1DuMkvYjHB2i3HGj3m0treADNLOEN6CXsFKy1DZcIcCLu9jS7AqWwM7LKeKAT9ddwExUKoTRc7UhfT5qJ7v1OwUgZ+Cusvt1pqQU8UdvdqyiEkknUU7SOdg7W2oa/wurbnZZWyBO1t+NEWyiEkknU0zQjuvc7rNpETlXxwb6azlv4q+F6FMLO4ISkGVk8CiFIIkEvpBdic+0O21XEh/hQe5y58GcRXjS/jpdYZelZNTtewd3saKgLCiFIwcuBdHWjB7C5dgftwEgZuByBkqF+dCe6CTso1aj1dlGCQiitRD1NLUDedkxGEUMhhMuKC6AZRWhj6Zg0oxa3pL8ECqGUEvQ0DeNlOoJxsquIx/jjvIXLiAsQsPp2R6XhiRCFUFqJgUKqkTFkbrsdLucBTtTPKjsvgeLF+NO9JbwJz4TtllPNTYx3d0chBOnonYivIz2KWcDthqn0cBXudiTcjWaVIqHaK9XIbwhEFUAhlBq6CTsEhRCubngAJlF0ADoIW6AQSixBT9OQt+2GNWXg6mKxvkxHYMhoCxRCiSXq6bYC9Gm0S1E9KWvkfTyQt3BFcf50OwagtY+xnpQ18L6eSCgUQql1daN2Aj1VhdS9th2FLNafCkhbuLLu7pQTfr4GCXVt2wpYYqCAhCIohHKARUfbKb2Qx2l4pxhop1h/Ad2E7YEOwla4rEgPswnbCR2E0B7oJmwndBC2QiGUHgaOtkejmRws49FYaxuuBQNH26OskVyo4YN8kFCEoBDKQW9PWtvML9Qida9mbwnv40Fd7aSOA2RviA/NruQ1JqnjkLc0I4sNoCLqICEEhVAOKCEJegGto1eHtbahnRxEYvCmu7Hr9VWlGXmCHtf/i3AgZAHdhNeUgQ5CaLc4tI5eCzoI20IhlIVEPd2GbsKryijEphPQXnEBNB3bUFxZjYkcr+BR6HH/HxRCWRjoTYvqeVG91HHI1akqLgpa3zsU2i/OX9hZxLGc/ZVsL+RRvtRBlDoO2UAhlAWBktgAmmbEPezlpRfyeDwOQrv5OxEfB3oMy9lfQZqRJQYiof6AQigXCXoBi45eCToIoaMwieIqUjFS5s9wLOQCswmvAkNGoaMwrf5KGswkq5TH+COh/oBCKBdRvvRMNa9okjoO+akykbPV3OCNvIUOiPPHE+Hl7SziA7yoi07qOOQEhVAudAKJ9qM7kLp/kVHII32pHU5V6Ij+XrSkgRc3SB2H/KQaOToIL4Gri4wk6AWMl/mr9EI2HO2i0EECJUP9aEYhEupSaUaGDsJL4HDIyAh0E15OehGPxaYT0HGxAQK6CS/RzMieYo47y0vg+iIjMf70cDmvbZY6Djkxc5JZzGPRsQ8dh4Gjf5VZwru7U097qeOQGRRCGXEQySAfuhP3sG0cKuNBztTbQeo4QIFi/On+Ut6ExtE2sLLaZaEQykuiHtPq/wQTJ6DTXHSkpzvdV4I7yz+kGRk24/0rFEJ5SdAL6CZsC1Pp4Xpg9e22GCcZhTwePe5/gSMiL8MDaGYJbzRLHYds4IkQrkccptW3cbCMBzhRfyep45AfFEJ5cbUjfTxpJhpzCCGEGOtJdRPv5YFCCJ0U50+3o6/hfzCD8EpQCGUnUU9TsTchIYSQ7UYWFyAgcaHTurpRO4HmVCOhCLm4GS/y6TJQCGUnAeNl/gcdhHD9YrHWGiGEEI6RMleGQig7CXoho5CbkbnoIARLwOrbLbIruIsOm3peHgqh7Hg7kFBXeqBU66lb30wOl2ETbbhemFbfAh2EV4FCKEfoJiSEZJbwAd7UGWvkw/UZ7ENPVfEqk9RxSA0dhFeBQihHCXqapvlCmF7I47CyGlw3O4EM8aG7Nd86ijVlrgKFUI5GBAqpBYxpO3PTi3gsOgjBEuICtL7B2dlqbmK8hzsS6vJQCOVI70S8HemxCu2mLicko5BhrW2wiLgAmlGk6ZHYqUY+IhBX+yvSTZw4cejQoY8++qhOpyOEmEymL774Yvv27e7u7uPHj09MTNy6deuuXbta/4enn37a3h5Ll1tdSzdhfy+NVoKTldzFjoZghBtYQlyA8OBWs5kTUasnVJqRJ6B95cqEW265Zfny5c888wwhhDF2xx13LF++PD4+vm/fvsePHyeE/Pbbb7///rvUcWqOxrsJd6CDECzHx4HoneiRcu0mFIaMXp3unnvuiYyMNBgM8+fP/+2337Kzs48dO2ZnZ9f2TTExMbNnz5YqRG1K1NO5e8yEiFIHIg1MpQfLall9e6C3Fk8qYz0paeD9PLX43dtJIIT06NHD29v7wIEDmzZtGjdu3Jo1a+bNm7dmzRrOL95A7dq1a/bs2cuWLSsvL5c0Wg3p5kbtBHq6SqP3sDuMKIRgSVrehiK1gCXqBaxVeBUXZ2n5+fkVFBTk5OT8/vvv5eXlgwYNeumll1JTU997773w8HBRFH18fH799deFCxfu3bs3KCjosp9VW1s7ffp0d3f3ti9OnDjxpptuuuz76+vrRVGjTzztEesrbjjfHBje3k5+1RzPShM9XyN2d2iqq5MsBpPJxBgzm7EPiGXU1dUJgpSDNYZ40EUHdHV1jRLGYEEdSvbNueIwb66a795Rjo6O1zz3LhbCxsZGR0dHe3v78PDwf/3rX4SQhIQEg8Hw6quvPvrooy3vmTlzZkpKyvvvv79o0aLLfpa9vf2YMWO6du3a9sWIiAhHR8fLvt9kMl3pR0AIuSGY7yziU/q19/KhmuO5pYQM9TO7Okv5XURRZIw5ODhIGIOaNDc3S3tyGhxJZZO5gjvqVbEJUYeSPb2EPdSHOjpq9JGwPXdgOkKI2WzOz88PCQkJDQ01mS4uwNC9e3fGWGFhoaura+u7Bw8enJeXd6XPsrOzS0pKMhgM7Y9P2ptEmRsRyN86aG7/IVLN8dxZbI7XS/xdWn67Oo6nHMjh5IzxZ7uKydiuavibtv94ljWSc9XmIb46qQ+/rAmEkJ9//tnDw2PIkCHjxo3bsWNHU1MTIWTLli1eXl5dunQpKytreWtVVdXatWsHDx4sZbxa0seT1jTzC7Wa69jYYeSx2EQbLC02QNBgN2GakcUFUB3y6ap0d9xxx44dOz7++GNRFEePHh0VFTVo0KB+/fqlpaV9/PHHOp1uyJAhgYGBXl5ee/bsiY2N/cc//iF1zFpBCUnQC2lG/rfuGmrTaGZkbwkf5qehrwy2MTyAPp+puWn1aUaeoEcZvAbd448//tFHH+n1ekIIpfTrr7/OysoqKyv76KOPfH19CSGHDx/Oysqqq6sLDw/v3r271AFrS0IATTPyv2npqGeV8S5u1At9c2Bpw/zowTLeYCaOahhS1l6pRv7OMBTCa9ClpKRc8tIlnXyurq7Dhw+3YUjwh8RA+vFxbd3DYq1tsBInHenjQfeVaGhmTo2JHKvgUWhfuRbcKciawZsW1POieqnjsCGstQ3Wo7XVt3cU8ihfqqkn4M5BIZQ1gZLhAXR7oYYeCtML+XAUQrCOuACaoaX9mNKMbARWVmsHFEK5axkvI3UUNpJXyxvM2CwGrCU+gO7Q0m1lKkbKtA+Okdxparf67YU8zh/nJFhLsAt1EOgpbaxc2GgmB0p5DHrc2wEXHbmL9KWnKnlFk9Rx2ATW2gZr086iozuL+AAv6qKTOg4lQCGUOzuBDPXXSupirW2wNu0UQmy91H4ohAqQoBfSjOrv2KhvJscreaQvUhesKC6ApmtjvEyakaGDsJ1wmBRAI92Eu4p5hBeGeoN1DfKm56p5pdr7GpoZ2V2MAdjthUKoADH+9FAZr2uWOg4rw8QJsAGdQCJ96U61PxTuLeHhbtTTXuo4FAKFUAEcRWLwUX/qphcyTKUHG4gLoOlqn0SBDsIOQSFUhkQ9TVV1NyEnZFcxj8XcCbC+uABB9d2EaUaegNvKdsN1RxlUP63+eAX3sKeBzlLHARoQF0B3FfFm9d5YMk7SCzFSpgNwpJQhPoDuKeZN6k1drLUNNuNpT0Jd6OFy1d5ZHirn/k7U30nqOJQDhVAZXO1Ibw+aWaza1E0vxFrbYDvqXn07tYAn6pFNHYBCqBiJgWqeRJFehCGjYDvqXn07zcgTUAg7AoVQMRL0VK3T6ksbSUEd7++F1AUbUfcT4fZChkLYISiEipGoF9ILuVmNyZteyIb5URGZC7bSy4PWmnherQrTKbuSO4o0zBXp1AEohIrh7UBCXGhWqQpTN6OQxwXgVATboYTE+AuqbB1FB2En4OqjJImBdJsauwnTsekE2JxauwlT0UHYcSiESpIQQNU3m9DEyL4SPtQPqQs2pdZuwlQjngg7DIVQSUYECmlGprLc3V/Kw92oBxZFBNuK9qWHVbeE79lq3mTmPT1QCDsGhVBJAp2Jpz09qq6JwOmFfDhuYMHmnHRkgDfdW6KqbEo18hGBuKp3GA6ZwqhvS6b0Qh6LNWVACnGq2/IaMwg7B4VQYRL0ausmzMBUepCI+jbpRQdh56AQKkxiIN1aoJ5p9edquInxbm5IXZBAy35MqqmExnpS0oCFKToDhVBhwt2ojtLTVSpJ3vRCHo818kEiQc7U1Y6erFRJNqUWsAS9IKAOdhyuQcqToKJuQnQQgrTi/NUziSIN7aKdhUKoPImB6ukmxFrbIC01TavHrvSdhkKoPKoZOFrbTLIr+GAfpC5IJi6A7lBFNpU3krPVfJA3sqkzUAiVp48nrTbxC8pfL3hnER/sSx1EqeMADRvoTfPqeFmj1HFctzQjiwugOlzROwWHTXkoIQl6Ybvyb2OxKz1ITqQkypfuVH7raJqRJ2DcWWfhwCmSOhYdzShiWGsbJDc8gGYUKX5KEmYQXg8UQkVSwW71nJBdRTzWH2cgSCw2QFB6N2GNiRyt4FFYub6zcBlSJIM3za/jRfVSx3EdjpRzH0fq7yR1HKB5cQE0s4Q3K/mZML2IR/lSR3S3dxYKoSIJtGX+k4JzN70QEydAFtztSFc3mlWm4IfCNCNDu+j1QCFUqgS9oOhuQkylB/lQ+urbqQUYKXNdcOyUSundhJhKD/Kh6NW3G81kXymPwW3ldUAhVKooX3qiklc2SR1Hp5Q0kJIG3tcTqQuyEBeg4CfCXcV8gBd1tZM6DiVDIVQqO4EM9VNq9u4oZDH+FKsDg0z0cKdNZqUuUpFagIkT1wuFUMES9EKaUZHjZdILMXEC5CXGX1DobWWakaGD8Drh8CmYcrsJ0UEIcqPQ1bebGdlVzOPxRHh9UAgVLMaPZpXxumap4+igJkaySvlQTP4FOVHo6tv7Snk3N+ppL3UcCodCqGBOOmLwpruKFZa9e0t4Tw/07YO8RPnSYxW8xiR1HB2EldUsAoVQ2RL1NLVAYYUQa22DDDmIZKA3zSxRWDalFvAEFMLrJhw7doyxPw24qKqqOnjwYElJSesrTU1Nhw8fLi0ttXl4cA0JeiFVaeNlMoo41toGGRqutEkUjJP0QoyUsQBh/PjxBoOhoKCg5d8LFiwICQl54IEH+vXrt3r1akLIgQMHevToMWHChF69er366quSRguXitfTPcW8SVGlML0Qm06AHMUG0HRFLVt4qJz7OdEALNh73YTDhw9HRka+8sorhJAVK1YsX778yJEjBw8eNBqNKSkphJAZM2Y89thj+/bt27t37xtvvHHq1CmpY4Y/uNmRXh40UyHdhA0NDYYbbimaF5ezN03qWAAudXTtZ789Nuzuh/4hdSDthRmEliIQQqZMmbJy5UpCyIcffjh79uzAwMDq6mpBENzd3Y1G47Zt26ZMmUII6dq16+jRo7/55huJQ4Y/S9QrZhJFTk7O6Tp78+inv/3pN6ljAbjUmtXfskmfrd+SKnUg7bW9EB2ElqEjhHTv3r2ioqKysvL48eOZmZmvv/56Y2Njz549V65cef78eXd3d19f35Z3d+/e/dy5c1f6LLPZfPDgwbq6urYvduvWTa/XW/U7aFyCnn6SzeYYFNBP0L1XHzFs4ND8n599/zWpYwG41NLX5k96YYHTXTOlDqS90ozsraE6qaNQAx0hxMnJiRBSWVlZWVl55syZY8eO6XS6+++//9lnn500aZKjo2Pru52cnC5cuHClz6qvr3/jjTecnZ3bvjh16tQ777zzsu+vqamxzJfQtkg3OsloX1FVXV8r9+O5IkeMmvjcf28wEUKqq6ulDudqTCYTY6ypSZlrucpPbW0t53JvtxhkiNj5w38MvzpuP19r8JJ1Z2FNTc3JampH7b14jbwzSXrOzs6ieI2tGnWEkLKyMkEQAgMD/fz87rvvPgcHB0LIhAkTpk6dOmfOnPLycs45pZQQUlpaepXHO1dX1+XLlxsMhvaH6Obm1v43w2W5ERLi2nymybWnq6yPJ+NkaXbzR/Gim5vjtd8ttZZC2JILcP0opa6urlJHcW1uhMyIYB+cclw1Uu673O7Nc74hiMs55RVEIIRkZGRERETY2dlFRkaWl5e3/KCsrMzDwyM8PNzJyWnfvn0tL2ZkZAwZMkSyYOEKFNFN+F0O83QgNwSiSwNk7bG+wtYCdqJS7gmVZkQHocUI69evf/7556dNm0YImTFjxrvvvrt27doNGza8/PLLDz30kJOT05QpU5544olt27bNnz+/sLBw3LhxUscMl0rUU/lv0vvmQTZvkNzvsgFcdOTxvuKbB2XdNEoISTXyESiEFiK89957zz///KRJkwgho0aN+uijjz777LOPPvpo3rx5Tz75JCFk4cKFt91226JFi86dO7dly5aWDkWQlRGBQpqRybkS/nqBmxj5vzDkLSjAk/2FH86yczXyTancOtJo5j09kFCWQS3Yg20wGDrUR1hdXY0Gbkvp+W3z18MbooNl2g2T8HPz1L7Cvd0VMLS1BfoILaumpkYRfYStnt1tbmJkSYxM2zD+fbhuU5GD/DsylUIxFya4ukQ93Xy27trvk0KqkRfUkXHdcLKBYsyMEL88yYrqpY7jCradr0MHoQXh2qQSm5+9+aWHxr7+7gdSB3IZiw6Y5xgEHc41UI4AJ3J3uLD0iFnqQC5j3N8fXTXzrt8WT5c6EPXAxUkNGhsbK0sKG0fP/GXzdqljudS+En64nDzYA2caKMysgcK/jrNK+U0l3b5rr/mOV47u3S11IOqBy5MaODg4fPL2wj7GLcEPvCJ1LJdalMVmRggO6MsApenmRm8KEf55TF7DRxknXpOWxJxZteoTOTb/KBSW51GJO2+/NTbxxqjfHPeW8EhfuXQeZFfybQXs80RswguK9NwgYeQvzdP6C86yuVL+5yTz6hPz++MGd4w0tBw8EaqHq46/NESYlmGWz6Dv17PYk/1FbEYPCtXXk8b4C59my+WhsMZEXtjLlsSKcrnVVQsUQlV5uLfQYCbf58giby/U8rXn2BP9cI6Bgj03SHjrIJPJlp+vHTAnBdGhfqiDFoaLlKoIlLw1TJy1mzXIYLDbm1ns4T6CF2bigZJF+9HenuSrU9JXwpxq/u9s9lo0LtqWh2OqNjcG0kE+dMlhifO2qJ58dZpN649BMqB4cw3iG1lM8qWbnt3Npg8Qg5zxOGh5KIQqtHiY8PYhc4Gk0+vfPWz+W3ch0Pna7wSQuZFB1MeRrDkr5c3ljkK+u5g/PQBXbKvAYVWhcDc6ubfwfKZkzaOVTeSTQkRaBQAAEylJREFUbPZMBM4uUInZA4VX90v2TMg4mZ5hfmuojAavqgwuVer0/CBxXS7fUyxN5r5/hP1fqNDVDW04oBK3dhEESn6/IE1CfXaC2YtkfDgu19aCI6tObnbk5Uhh+k4JplLUNZMPjppnDcSpBepBCXlmoPDKfglaWapN5KW9bEkMpkxYEa5WqjWpl9BkJqvP2Lpj4+PjLF4v9PdC2oKq3BMulDQS22/8uXC/+aYQGo0pE9aEQqhaAiXvxYrP7GJ1zbb7pSZGlhxmcww4r0BtREqeiRAWZdn0ofBMNf/sBFsYhdHX1oULlprFBdBh/vSdQ7Z7KFx+kvXxJFGyWeMNwIIm9RKOlJO9JbZ7KHxmF5sZIWL0tbWhEKrcOzHCe0fMF2ptkbpmTt48yOYacPcK6mQnkGn9hTcP2ujOcksBP1DKp2PKhPXhEKtcqAud0keYl2mL1P0uh3k7kBGBeBwE1Xq0r7C1gB2tsPqdJeNk1i7z4mGCI24srQ+FUP3mGsRN+Xy39adSvHmQzRuMrAU1c9GRqf3Et63/UPjvbOYokrFdcYm2BRxl9XO1IwsjhelW3pXilwvczMgtoXgcBJWb1l/46Tw7V2PFfKo2kZf3sfewy4StoBBqwoSegpmTVaeteBv7epZ57iABeQuq52FPHuolWHUM2oJ95v8Lo/LZWFT1UAg1QaBkSaw4ezertc5Uim0FvKCOjOuG0wk0YUaEuOIUM9Zb5cNPV/H/nGQvR6KXwXZw5dKKWH8aF0AXW6dvY1GWea5BQDsOaESAE7mvu/D+EavMKZyxiz07UNQ7WeOz4fJQCDVk8TDhg6Pm85bu2zhQyo+Wkwd74lwCDZk9UPj4OKtosvDHbs7nh8v4U/2RTTaFw60hIS708b7Cc5aeSrHwAJsZIdjjVAItCXOlt4QK/zxqyWwyczJ9p/mdGMEBzaK2hauXtsw1iGlGvt1y6yUer+CpBWxyb5xIoDlzDcKSw2YL9rv/6xjzdSC3d0E22RqOuLY46cjCKGH6TrOltlZ7PYtNGyC62lnm0wAUpI8nTdALn2Zb5qGwoom8vN+8JBYPgxJAIdScB3oIjiL52hJTKS7U8l8usKn9cBaBRs0bLCw+yJosUQoX7DOP7SoM9MaQMwngEqY5lJDFw8S5eywwleKNLPZIb8HT3hJhASjQYB/a15N8eep6K+GpKv7lKTZ/CB4HpYFCqEUx/nREIH3z+jaUKawnX59mT/ZH6oKmzR0kvnaAXee6TdMzzHMNYgCmTEgEhVCj3ogW/nnsupaJevew+YEeAjaIAY27IZDqncianM4/FG7K58crCboYJIRDr1HBLvSJfuKcPZ3M3som8mk2mxmB8weAzDGIrx3o5PizZkamZ5iXxIiYMiEhXMi069mBQkYhT+vUVIqlR9itYUIXV3TsA5AxYVQUyG8XOpNKHx5jAU5kTBhSSUoohNrlpCOvRQvTMjo8laKumSw7an5mIE4egItmRQiv7O9wp3t5I3n1gPndGDwMSgzXMk27r7vgYkdWdHDM27+Os0S90M8T97AAF40PF8oaSWoH21de2mce102IwJQJqaEQahol5L0Yce4ec5Wpvf+LiZElh9lsA84cgD+IlMwaKCw60IGHwmMV/JszmDIhC7icad0QXzo6WHij3VMpvjjB+nkS7JQGcImJPYWjFSSzpL0PhTN2mp8fJPo6WjUoaBcUQiBvDhU/yWZnq6+dwGZO3jrE5g7CPSzApewEMn2A8EZWuzoafrnAz9aQx/viCiwL+DMACXAiT/YTn9197QT+9gzzcSCJejwOAlzGo32EVCM7Un6Ne8pmRmbvNr8zTLTDBVge8HcAQgh5ZqCwp4RvK7haAnNC3jrIXhiMx0GAy3PWkSf6iYsPXeOe8oOjLNCZ3ByKG0q5QCEEQghxFMnr0dfYleLn84xxZC/A1TzVX/j5/NU6GsoayaIs8xJMmZATFEK46J5wwd2OfHHyijezb2Sx5wYJKIMAV+FhTyb3Ft6+8kPhC3vN94YL/b2QSTKCQgh/WBIrzss0VzZd5kdbCrixntzVDScMwDU8PUD86jQrqLvMj45W8O9y2IuYMiEzuK7BHwb70JtChEWXm0qx6ID5uUGCiLtYgGsJcCJ/6y4sPXKZPJqx0/ziYNHHwfZBwdWgEMKfLIoWP81mJyv/1MOxv5QfryAP9MDZAtAusw3Cv4+z8sY/vfjTOXa+hkzpgzySHd3NN98cHBz83HPPhYeHNzY2PvXUU60/S05Ovuuuu77//vt169a1vrh06VIHB9zPqFaAE5kRIc7Zw9Yk/dF688p+9sxAwR75C9A+oS50TJiw7CibN/hi2jQxMms3ey8WUybkSJg1a5a/v/+NN97Y2NhoMpk+/vhjg8EQGRkZGRkZEhJCCNm1a9f58+cj/0cQ8GdUuZkRwsEyviHv4kPhsQq+3cgm98bfHaAD5hiEpUfMNf9bvPD9I6yXB7kpBL0LcqQbOXLkyJEjf/nllx9++GHMmDGEkIcfftje3r7tmwYNGjRlyhSJIgRbsxfI69HCrF3mvWN1IiWLDrDpA0QXndRhAShKH0+aGCh8ks2mDxCKG8gbWea0W5FFMnXxNn/48OGZmZkt//3kk08+/vjj33zzDecXnwk2b948ceLEF1988cKFC9KECbZ1VzfB04F8ms3O1/BfL7B/YO9sgI6bN0h4+xBrYuSFTPODPYXeHngclKmLdyg+Pj45OTmiKE6bNi0yMrKsrGzOnDlpaWkffPDBwIEDg4KCfH19N2/ePGDAgH379nXv3v2yn1VTU/Pggw86Ozu3ffEf//jH2LFjr/R+y34ZjbPs8VwYQe/8vuDL2iMT4keKjQ3Vjdf+X9TEZDIxxpqaLjeVBDqutra29cZaO7rbkx66httf35rpHb//brfqdizn2064eLafs7OzKF5jvgptOTtnzpxZV1f34Ycftv7gwIEDQ4YMqaysdHNza33xzjvvDA8PX7x48WU/a8CAAbNmzerRo0fbF8PDwwMDAy/7/urq6rYfDtfJssezuro6IGJ4fXj8w0N8/r34FUt9rFK0FEKMC7OUmpoaV1dXqaOQQP8bxhz1iAzI2Wg8uMOCH4uLp2VdfCLMycmJiopq+4OePXtyzktLS9se7p49exqNxit9liiKgwYNMhgMVooVbMzVnpopc7FDew5AJ3k6iHaUezqgc0HWBELIuXPnNm7cOH78+Pz8/Pr6ekII53zp0qXBwcFhYWGnT59ueWtBQcG3334bGxsrZbxgK25ubpkb1/531h2LX3lR6lgAlGrdmq++mxydsf4nqQOBq9ENGzbszJkzs2bN6tmz5/Lly5966qmuXbtWVFSIorhq1SpBEG666SbOuZub26lTp+69995HHnlE6pjBRsLCwsLCwqSOAkDBXF1db7v1VqmjgGvQffzxx4GBgf7+/oSQCRMm3HrrrefPn3d3dw8LC2vpYMzOzj5z5kxtbW2XLl08PT2lDhgAAMCSdJd06Xl5eXl5ebV9RRCES8a/AAAAqIaUXbg1NTUaHFFtJYyx2tpaqaNQj6ampoaGBqmjUI+6ujqz+TKLUEPnVFdXSx2CqkhZCCMjI/Pz8yUMQE2OHDkyatQoqaNQj/fff3/BggVSR6Eed99999atW6WOQiVMJlOXLl2kjkJVMKgXAAA0DYUQAAA0zZKLwJrN5hMnTrT//SaT6ejRoyUlJRaMQbNOnjzZ0NCQlZUldSAqUVBQUFZWhuNpKTU1NadPn/b19ZU6EDVobm7mnOPkbKfevXs7Ojpe/T3UgsNV3n777c8///yaq7q1qqqqcnNzoxQLl1iA2Wyuq6vDqkuW0tjYyDm/Zv5AO9XW1jo4OOh02H7BMiorKz08PKSOQhm+/fbb3r17X/09liyEAAAAioM+QgAA0DQUQgAA0DQUQgAA0DQUQgAA0DRx/vz5tvlNjY2N27dvP3/+fFBQ0CUjSznn+/btO3TokJ+fH8bptVN9ff22bduMRmNQUJAg/OmG5vz587t37y4rK9Pr9Zf8CK7kyJEj+/btc3Nzu+z+sYyxnJwcSinOz/Ywm827du06fvx4YGCgnZ3dJT+tq6tLTU3Nycnx9PR0cnKSJEJlqa2t3bZtW3FxcVBQ0CXD7Ovr63ft2nXy5El3d3dnZ2epIlQ8bhPFxcV9+/aNj4+PjY01GAzl5eWtP2KMjRs3rlevXmPGjPH399+3b59tQlK0c+fOhYWFJSUlDR48OCEhob6+vvVHjz32mL+/f0pKSt++ffv3719QUCBhnEoxa9askJCQ22+/3dfX95dffvnrG95//31BEF599VXbx6Y4DQ0NI0aMMBgMycnJoaGhOTk5bX+6devWgICA2NjY5OTkuLg4iWJUklOnTgUHBycnJxsMhpEjR7ZM7Glx8uTJkJCQkSNH3nXXXd7e3mvXrpUwTkWzUSF88cUXb7vtNs45YywlJWXRokWtP9qwYUNoaGhVVRXnfOHChbfccottQlK0xx9/fPLkyZxzk8kUHR39ySeftP4oMzOzqamJc242m0ePHj1z5kzJolSIU6dOubi45Obmcs6//vrrPn36MMbavuHs2bMRERFJSUkohO3xxRdfDBo0qOUkfPzxxx955JHWH1VWVvr7+3/zzTct/2yZGA5XN2nSpKlTp3LOGxsbDQbDl19+2fqjaf/f3v2FNPWGcQB/d5ZnC+XgcSjnKEnlaksnnkzTSRRdDC8KzWluOBV2o3ihidCNIBgkqDUqFINEsS6i3SwIigZBQjSCQGOolbML3Ryo8z8Z29DzuzjwctCfFV6c02HP52rj2cWXl/flZe95z/veumWz2YTPjx49KikpkSei8km0bubxeOrr6xFCKpWqrq7O4/GIS5WVlcKb4A6H4+3btzs7O9KkUi7cnseOHbPZbOL2vHDhgrAYRRCE0Wjc3NyULaVCvHz58vLly1lZWQghq9U6Pz8/MzODqzzPNzc3u1yu/10yBQd5PB6bzSZ0QofDIe6cr1+/ZhjGarVOTU2trq7+/eEbiQwPdpIkb968KW5PrVar0WiEzyRJwrr9kUk0EQaDQXzXeXZ2digUwqVQKIRLJ06cQAgtLi5Kk0qhYrHY8vLyYe2JhcNht9ttt9ulTac8oVAIn+Wv0WgYhgkGg7g6PDzMsqzFYpEpnfLsG+yrq6v4grC5uTmNRsNxXFtbm8Fg6Ovrky+mMmxtbW1ubh422G/fvr22tlZbW9vU1PT06dOBgQGZYiqeRCceRaNR/MycJEnxTW/iEkEQarU6Go1Kk0qhhEUn3GgajebgzXnb29tWq7WxsRHuZvqjWCwmPvqLJEncA8Ph8L1793w+n0zRFCkWi4kHO0IoGo0mJycjhLa3t798+eL3+41G4/fv3zmOq6qqOnv2rJxx/21CVzxssM/Ozn779q2+vp6iKJ/P9+nTp/z8fHmCKpxEEyHLsvhw7UgkwrIsLjEMg0sbGxvxeFxcBQelpKSkpKREIhHhD/TKysq+Fvv58+f169cLCgr6+/tlyqgkDMN8/foVfxW3p8vlysjIcLlcCKHp6em1tTWWZZ1OpzxBFUI8oldWVrRaLU3TwleWZfV6vdFoRAgZDAa9Xj85OQkT4W/odLqkpKRIJJKeno4ODPbOzs6WlpaOjg6EUGlpqcVicTqdcKDrEUi0NGo2m/G1nO/fvy8rK8OlsrIyXBofHzcYDGlpadKkUi5xe46Pj4vb89evXxUVFadOnXr8+DEcaP43zGbzhw8fhPvT/X7/7u5uXl6eUCovL6+oqKBpmqbppKSk48ePw5PCP9o3os1mM+6Hly5dWl5ejsViCKF4PL60tJSRkSFXTkUgCKK0tPSwwR6NRvFzQa1WG4/HeTg7+mik2ZPz+fNniqIGBgYePHhAUZTf7+d5/sqVK/fv39/a2srKympvb3/27Fl2drZ4AyQ4jNfrpWl6eHj47t27qamp8/PzPM/n5uaOjY01NjbqdLqenp7e3t7e3t4XL17IHfZft7e3V1hYWFdX9/z5c47jOjs7eZ4fHR01mUzin924cQN2jf6NYDBI0/SdO3dGRkZ0Ot2bN294nnc6nS0tLTzPWywWYXuX3W4vLi6GjaN/9OrVK51ONzIy0t3dnZaWFg6HeZ4/ffq02+0eHBxkGObJkydut7uwsLChoUHusEol0Qv1mZmZV69efffu3fr6+sOHDzmOQwgRBHH+/PmcnJyamprJycmZmZnW1laHwyFBHqXLyckpLi72er3xeHxoaEiv1yOE1Gr1xYsXKYoymUz4PfrU1FSTySRr2H+dSqWqra0NBAITExM1NTUdHR0qlYogiMzMzKKiIvwzkiTz8/OF5WjwGxRFVVZWfvz48cePH93d3eXl5QghgiDOnTt35syZ6urqhYUFn8/Hcdzg4CDe9AgOYzAYOI7zer0IoaGhoZMnTyKE1Gp1SUnJtWvX8vLyfD5fIBCoqqrq6uqCjbhHA9cwAQAASGhw/hYAAICEBhMhAACAhAYTIQAAgIQGEyEAAICEBhMhAACAhAYTIQAAgIQGEyEAAICE9h+kQDZuOF1uIQAAAABJRU5ErkJggg==",
"text/html": [
"\n",
"\n"
],
"image/svg+xml": [
"\n",
"\n"
]
},
"metadata": {},
"execution_count": 5
}
],
"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.ρ.real[:, 1, 1], label=\"\", xlabel=\"x\", ylabel=\"ρ\", marker=2)"
],
"metadata": {},
"execution_count": 5
},
{
"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: 4%|▋ | ETA: 0:00:03\u001b[K\rDiagonalising Hamiltonian kblocks: 9%|█▌ | ETA: 0:00:02\u001b[K\rDiagonalising Hamiltonian kblocks: 19%|███ | ETA: 0:00:02\u001b[K\rDiagonalising Hamiltonian kblocks: 24%|███▉ | ETA: 0:00:02\u001b[K\rDiagonalising Hamiltonian kblocks: 33%|█████▍ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 41%|██████▌ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 50%|████████ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 59%|█████████▌ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 65%|██████████▍ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 74%|███████████▉ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 83%|█████████████▍ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 93%|██████████████▉ | 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+gvaeTAAAgAElEQVR4nOyddUBU2RfHvzNDtwICigUGBtjdAajYqGDhCgjmusbPREUwMNddA10DW7DFQjFWsBOwxSQEVAyaGWbu7w9YEASmXszAfP7ivbnvnAPMe+fdc889h0MIgQoVKlSoUFFZ4bJtgAoVKlSoUMEmKkeoQoUKFSoqNSpHqEKFChUqKjUqR6hChQoVKio1KkeoQoUKFSoqNSpHqEKFChUqKjWSOsLk5OR58+ZFRETQao0KFSpUqFDBMGoSjpsyZUpERISJiUnXrl1pNUiFChUqVKhgEolmhIcPH9bU1Gzbti3d1qhQoUKFChUMI94Rpqam+vr6/vnnnwxYo0KFChUqVDCM+NDotGnT5s6da2ZmJnZkcHDw/v37udxiztXd3d3JyUlyg1684NjYFFR9i4riGBmhTh2Sf75BA8JVhuQeoVDI4/HY0p6QgMREjqWlSEuLY2ws9eUvX8YlJJi+efPZzs7MyEirXj2iVvw7kpND3r+HjQ2HKoMlJykJOjowNCw4TE9HWhqnRg26agS+fMmpV49I+J8khBBCuJJ9QRMTOTo6pEoVucxjHpFIxOFwOByp//Xv3nHMzYm2dsFhYiJHX58YGFBsXnY2EhNF9eqV/BdkZODNG05KSkp8vIaVVU6LFlpGRkbSCn/3jlO1ap6hIS82lmNtTX78QFYWjd+9csjMxOfPnPynooTI9kTKy8O7d5z69YsUvXjBadiQSP//F0NiIkdHR1ilCl0Pdx6PJ/bGFOMIQ0NDv3z54ubmJok+Pz8/T0/PunXr/nyyadOmkv8Pdu/mbNjAuXdPpK6OnBx07codMoTs308AuLpyd+8WtWoloSQ2ycrK0tfXZ0v7vHmc6Gj06yeqWlVtwQKpb9Tr1yMXLtwjEqX36vXvkydqCQlo0gTNm5OWLdG8ObG1havrvMhIvdevJ5qYmNBhfzls2cIRCrFqVcEvFRHBWbeOc+2aiCZ1rq7cPXtELVpINFggEAiFQnV1dUkGz57NGToUrq5KVuZXIBDweDw1NUkTCwpxcuIeOSKysys49PLiTJ6MgQMp/vUPHrw2Y8bKjRtHDB7s/ugR59EjPHqEqChOXByaNIGZ2etr1+bq6vKvXz9iLP0b4vjx3EWLch0cNNzdufPmET09jBnDOXNG1LIltb+EeKKjMW8e98YNKb72sj2RIiM58+dzbt0qUJSQAHt7bmKiiHJHOHjwoh8/oi9d+quE76AEkUgkFArFv6GScvH09DQ1NbWysrKystLR0TE2Np44cWJZg21tbaOiosoXWA6xscTUlERHFxxeuECMjcnBgwWHM2aQ5ctlls0oaWlpLGqvV4+MGSPatSvXxUVGCUlJSXw+P//n9HRy/TrZuJGMH0+aNycaGrsAe2C+k9NwyiyWmPfvSbVqJDu74FAgIBYW5PlzutR5epJNmyQdzOfzswstE4e1NY1m00d2dnbhF0NyXr4k1asTkahQCDEwIN+/U2wbIURHpx5wHGhqYEC6diV//EH27CExMUQgKBiQmpqakZEhg2SRiBgYkA8f0gkhwcGkSxdCCDl3jpiZkQcPKLNfQjIziZ4eyc2V4hLZnkizZxNf36LDf/4hY8bIIEYMt27d4nDsgDVNm3anXjohQqFQki+tGEeYmZn59T/s7e39/PzS09PLGiyPI+TzSdu2ZOvWojOzZhFtbfLpU8FhWBjp2lU22UzDriPU1CTBwcI7d7IbN6Ze+J07D7ncGlxuY3//ldRLl4C+fcm+fUWHs2aR+fPp0rVzJxk9WtLBkjvCHz+Ivj7Jy5PdMLaQzRGuX0+8vIoO6buRGzbsxuV20dGxKXS6VPHuHbG0LLiv8/KItTW5eZMQQo4dIxYW5PFjitWJpVkzcu+eFONleyI1alRMy7BhxW49qkhOTubxagOdx42bQr10iR2hmAmjjo5Olf9QV1fX1tbW09OjbNb6E3PnokYNeHsXnQkNRZ06MDUtOOzWDVFR+PGDDuUVh8RE8PlwciINGojev0dODsXy27ZtkZT08NGjQz4+8ygWLRkTJ2Lr1qLD8eOxdy+EQlp0deyIW7eoF/vwIWxtwd4iMtOcO4d+/YoOL16EoyMtil68+Pfq1WXp6U8pj93FxKAwrsvj4fffsX49AAwdir//hqMjnj2jWGP5tG2Lu3fpVfHuHVJTURj4FQpx9Sp69aJekZmZWXz8U23tw0FBm6iXLjFSrE/u3LlzwoQJdBhx4QKOHsX27UVnUlMRF4chQ4rOaGmhQwdcvUqH/opDSAj09aGtDQ0NWFnhxQvqVVStWrV+/frUy5UMJyckJODRo4LDJk1QowYuXqRFV8OGSEtDcjLFYh89AvMLS2yRkYG7d9GzZ9GZCxfg4ECXuhYtWkiYryQVjx/D1rbo0NMTERF4/RoAhg3DihXo3ZuWe60s2rTBvXv0qsh/fSn8W969i1q1YGFBiy5zcx1dXfOUFFqES4gUXxpzc3PDwow96vj0Ce7u2L0bPy9gX7oEbW307VtspKMjLlygXH+FIjwchU7K1haPH7NqDQ3wePDwwI4dRWfGj0dQEC26OBy0bUv9pPDRI0iYgFMBCA9H+/YoTNRITERKivK9B5RwhDo68PTEhg0Fh+PGYflyODri3TuG7GFgRlhiHk/r6wsAa2sSG0ujfLGwvB2BELi7w9Oz2DsjgDNnkJuLdu2KnXR0RFgYk9YpH48fo1Ongp8rpCMEMGECgoORllZwOHIkwsPx5Qstujp0UDlCufj1eWpvD6XYBPUzJRwhgN9/x6FDRd+68eOxaBF69MDz51k5lC9I/ELTpoiPL7oFKCc7G9evw96+6MyFC3QFtPOxthblz7DZguWv5Lp1+PwZPj4lz1+4gE6dUCIXvXFjiERg98VBwUlOLoon29lVTEdobo6ePXHoUMGhoSH69Ss6pBbKlwmzs/H2LRo3plKmwkIIzp9naIGQPnJz8f49bGyKnTQzw5AhCAwsOuPpiWHDomxtu9Su3fHly5e0msTjoXlz3L9Pl/wrV9CyJQo3W37/jmfP0LEjXeoAWFlVYkf44AHWrEFISEmH9/o1MjPh7FzKJQ4OqklhmURHgxAU1oK1tUVMDKsG0cbEidiypehw/PhiwVIKsbH5fvfu2tOnKYvIx8SgYUNoalIlT6GJjoaWVlGsXiTClSvF5hlKwfPnsLKChkbJ8//7H7ZsQXZ20ZmWLZ8Dnb5+bUm3IwTN0dES8/hLl9C1K71fWmtrUkkdYWYmRo/Ghg2oU6fkR+HhIKT0kLRqmbAcQkJgbFwUd6pZE1lZSE1l1SZ66NkTubm4fbvgsFcvZGQgKop6RUuXrhAIMHasb1JSEiUCK1Vc9OxZDBhQdHj/PszNUb06ewbJxK9x0XwaNkSrVti/v+jM8OHDZs+upa5u27q1FLW0ZIPWfJmwMEYXCAFYW4sq4xrhyJFTzM071akTPnJkKZ+eOAFDQ1hZlfJR796IjERuLt0GKiUREWjSpOiQw0HTphUzOsrhwMuraB8Fh4MxY2hJmbG1baitHSYS8Q2oKAh2//5DX9+hKSkbxA+tEPy6QKh0cVGU7QgBzJ6NtWsh+q/Gi7q6ekDA7NGjpwcF0b45hr4Z4bNnyMsr9iRhIKBdr17lS5b58ePH+fPRGRmb8vL2/vqpSISbN0vmixZiZIQmTXDjBr0WKikvXqB792JnKmq+DAB3d5w5UzTfHT8eBw9S/4Y0darH9et71NVvArryS5s/f11Kyv+uXz+UkZEhvzQF5+tXPH2Kn5u2KeMCIYpvIixB9+6oUgVnzhQ7OX06AgMhENBrVd264PPx8SP1ks+eRf/+RYfPn4PLRYMG1Cv6GQMDoqUFFndQsOAIDQ0N9fS616gxa8mSib9++uABgNIXCPNRRUdLRSTC169wcSl2sgI7QiMj9O+PPXsKDuvUga0tQkOpV9SiRY3OnTVDQigQ1bv3MB7v944dG9JUlUKhCAtD9+5FC0tpaYiJKUppViLKmREC+OMPrF1b7EzTprCyouWrWII2bWiZFDJWAKEE9euzmQjJgiO8dw9c7rJ376506VLKbXHuHAQCdOtW5uV9+qjyZUohPBw8XsnctgrsCPFfygz5r3QzfRsKvbzwzz8UyHn9eoif371z50oJhFQ8SjxPL19Gp04obEChLHz9isxM1KxZ5oDhw/HxY9FydT5Tp2IT/WVS6FgmTEvDw4fo0aPoDAMLhPnUqwcW82VYcISrV+N//yuZKVrIiRNo1AjlvDG3bo3ERFpiAkrNiROlpCHY2uLJk6I1jApG+/YwNMSVKwWHw4bhzh3Ex1OvqG9fJCcXlbORjfR0HDuG336jxiQFRyjEhQvo06foDGPPU2qJiYGtLcqp2cbjYdq0goprhQwditevac/ZpmOZ8OJFdOoEHZ2Cw9xc3LxZcpM3TVQuR/j2LSIiMH586Z/m5OD58/LiogB4PPTqhfBwOqxTYm7fRrNmJU8aGqJqVbx/z4I9zODlVbSXS1sbzs7FsviogsuFh0exEoAysG8fevdWvpxJ2bhzBzVqoFatojPh4Uq5QFh+XDQfT09cu1bsIa6mBi+vYjt86KBdO9y/T/Frbol5fGQkmjaF9N0bZaFyhUZXrcKkSWVO+CIjweNh0CAxQlTLhL/y5k3pD5qKuq0+n9GjcfUqEhMLDvOjo4SGNn+enggORnq67BJ27ChWU75iU+J5GhsLPh+NGrFnkKxI4gh1deHhgb/+KnZy4kQcPoyvX+kzDVWrwtgYr15RJpAQXLhQLFGRyfymSjQjTEnBsWOYMqXMAceOgccT/81zdER4eIWN+MlARgYyMzFsWCkfVeBt9QD09ODqil27Cg47dACff+2PPwLT5XFZpWFhgW7dIHPKzK1bSEsrtvRSgTl58uyWLb83bVq0qbxEmFSJePy4zJTRn5k+HQcPFqvzZ2qKfv2wezdtlgEA6tSJWLdufy4V2dLZ2dljxvjy+QG1a+cVnmQyoF2JZoQbNmDUqKLOSr9y/jy6dCkvIp9PjRowNy/IL1UB4PhxaGrCzKyUjyp2vgyAKVOwfXtBJ6akpKQvX+Zu2vR5wYIAyhXJkzKzbRu8vZWvxqYMiESi8ePnffvmtHHj3MKTSrpASAiePSu2o64szMwweHCxBmEApk7Fli00vq+/evXq9m2foKD7a9dull9aSMjhkJCctLS3F/6LtiUnIyEBrVvLL1sijIygqcnaDgrmbs30dGzfjunTyxyQmoqkJIwaJZE0VXT0Z86dK6VATz4V3hE2boxatXD2LADo6ekZGOQQctvYuJa466TG0RGfP+PhQ6kv/P4doaEYN45yixSR7Gwun2+lr+/fp0/BFkI+H5GRtLSyo5t372BkJOkK2f/+h82bi1Vca98eVavi/HmarIO+vr62dhohr7Kzq8kvjc9vDlwyNn7Y+L9KuBcuoHdvRhtnshgdZc4Rbt2KPn1gbV3mgPyuchK+Oaoc4c88eIA2bUr/yMYGcXHF7s+Kx8SJBSkz+vr6sbE3PD3//vKF+uW4/JQZGSaFe/agb19Uo+BhpQRMnoyhQ08lJV1YsmRm/pkbN9CoEapWZdcuWZAwLpqPjQ1atiyZqzVlCjZTMFsrHQsLi+fPL23Zsj44eIycnSji4+Hr2+zChRtxcbfr1q2bf5L5Agj16rEWHWXIEQoE2LgRs2aVN+bQIZibS/q86NoVjx/j+3dKrFN64uOLFYP4GTU11KuH58+ZNYhZhg/Ho0cF75K6urorVtQLDi7KoKEQT08cPix1ysz27ZUlTWbzZkRFYds26OoWFeJR0oIy+G/vhOTMno1ly27cvl20p8HVFY8egb4S3CYmJt7eNvb2kKdjukCAkSMxaxZ69dJUU1PLPykS4dIlpiukV/wZ4d69aNpUTLnhyMgyK6v9iqYmOnYs2kNWmfn4EXx+serGJajw0VFNTbi5FTWgMDHBuHFYs4Z6Rebm6N4dwcFSXBIZCaEQXbpQb4yicfcu/P1x/HjRLrR8lHSBEJKljP5Mbm5YYuJKe/uF169fzz+jqQl395Jrh5Tz55948QL79sl4+R9/wMQEM2cWO/noEUxNy6skQAf161doR0gI1q/H3LnljXnzpqAfheSooqP5BAdDX7+8mh0V3hEC8PZGUFBRrdE5c7BvHy1VF7y9pdsflp8mIzb/S9n5+hWurggMLLn28fkzPnxA27YsmSUf0jpCLperrS3KzhYJBEX/78mTsXcvjU10AWhp4eBBzJoly1aKQ4cQHo49e0p+RVl5fangM8JTp6CrW17VNAAnT4LDQYcOUoh1dKRxIVqJCA8vavlWKhV7K2E+1tZo3hzHjxccmplh7FisW0e9IgcHpKdLmrGcmoqzZ+HmRr0ZCoVIhNGjMXJkUVPoQi5cQI8e+C/epkzk5iIuTrpi0w4ODpcuLenZc82VK0XFI2vUQI8eOHCAegt/pkkTLFqE0aPB50tx1ZMn+P13HD4MQ8OSH7ES0GZxBwUTjnDNGsybJ2bMkSOwtS2z7lqp2NhATQ0vXshjWkVAbC3jir2VsJDClJl85s7Fnj349IliLRwO3N0lTZkJCsLAgUqZJyIVvr7IzYWfXykfKe8C4dOnqFevlH685dOuXbu9e1tu24bo6KKTU6Zg0yZa6jz8zNSpqF4dS5ZIOj4jAyNGYO1aNG9e8qP0dDx6xEI8v0oVqKtTf89KAu2OMCICX75g8ODyxohEiIoSU1mtVOztVdFRpKRg6NDyBlhags/H589MGcQSAwbg3Ts8eVJwaGGBkSNpmRS6u+PIEfHBLkIqRZpMeDh27cKBA6Xk2ROCS5fQuzcbZsmNtHHRQszNsWwZvL2LdhDmz4mvXqXQulLgcLBjB/btw+XL4gcTgvHj0bNn6bt6rl5F+/Yl13qZga3oKO2OcNUq/O9/YrYSP3gAobBkCyFJUC0TRkVBJBL/7ta0aZGHqKioqcHDA9u2FZ2ZOxc7dlD/gmlujp49ceiQmGFXrxakdLHI6dNna9fuMHv2Iprkx8XBzQ3BwbCwKOXT6GgYGuK/bHwlQ2ZHCGDCBOjqFsuRmTSJiX4UpqbYvx/jxxercVMq69YhLq7M10QW85vYio7S6wgfP8bDhxgzRsyw/fuhp1d6S/ry6d0bN25U8E1y5XP4MExMxJcsqSTRUW9vHDxYtL3B0hIjRmADDQ3hfy72XRb5aTI0kZaW1qpVn9at+6SlpQH48AHh4di+HYsWwdMTAwagfXs0aIBBg/6Mi/Nbty743r08sTKlJTcXzs6YPx+dO5c+QHnzRSHlJsIScDgIDMTSpUhIKDjj5obISCbK33fvDldXjBtXXiT21i2sXYuQkKJukSVgMaBdMWeEq1djxgxoaYkZdv58mfdS+RgYwM4O/+UqV0YiIvBfIYjyqAyJowBMTATAoFq12pw/X9CdZMECbNsm/gVZWuztkZmJ+/fLHPD5My5elC4LWioWL1788GHDBw8aGhmd4XBQty6cnDBzJv75B//+i6QkGBmhY0fUrl2fw5nF41m2a6dmYYFFi6RLpiifadNgbY3ffy9zgPIuEEL6TYQlaNAAkyYVFdLS0cHYscXCFfSxfDlSU8t8UUtJwYgR2L27zFpU794hO1uiwnJ0UAEdYXw8zp+Hl5eYYTk5ePtW9sy6Sh4dffFColLOlcQRJiYmCoW5378vOnjwbP6ZmjXh7Iy//6ZYEYdTMgxbgp074exMY/+aqCgXIILLPRMYaJGeDpEIfD7S05GSgtevcf8+wsKwezfevQv88uXfvLyrnz7ByQl//w0dHXTujIgIeQ04cAAREeUlDeW/KHTtKq8iVvjyBTk5qFFDLiELFuDFC5w8WXA4ZQp27UJOjvzWiUFdHQcOwNe3lCBQXh5cXODlVV4N9LAwODiwtuGHtcRRQh22trZRUVGFh7//TubOFX9VaCjh8Uh6uoxK794lTZvKeC1NpKWlMaMoL49wOOTp01/P52VmZhY3iejpEaGQAqUCgSArK4sCQfQwb94KY2PX+fNjC8+8f09MTMjXrxQrSk4mVaqQHz8KDvl8fnZ2dv7PIhGpV4/cvUuxxkKmTiU8Hrl0KadQo+ScP0/atydcLjEwIGPGkNRUKa7Nzs7m8/mEkJgYUq1aKV+8nzlzhvTqJa118kLVrXf5MunalQLhERGkRg3y7VvBYb9+JChIXtskZNcu0qQJyb9TCy2fPZvY25O8vPIuHDSIHDpEv30/IRKJMjIy8n/++pUYGFApXCgU5n9py4cuR5iaSoyNSWKi+KuGDSN168quVCgk1aqRuDjZJVAOY44wLIyoqZVy/ldHSAipU4fExpYyWFoU3BESQmJjiYkJSUoqOjN+PFmyhHpFw4aRLVsKfv7ZEYaFkWbNqFeXT0AA4XLJ4cNyCUlPJzNmEDMzwuGQRo1Ily6zLCzaXLhwofyr8h1hWhqxsSH79olRMW0aWbVKLiNlgKpbb8MGMmUKNcK9vMjUqQU/nztHWrSQzzJpGDWqQHW+5adOkdq1yefP5V0iEBAjI/LpEyP2/cfPjpAQUrUqlQZI6AjpCo1u2oQhQyTqxx0RIddCApdbeRvWnzghRcfzShIdBVCvHsaNw6KfMiUXLsSmTdRXpvX2Lr161rZtmDyZYl35bN2KBQuwfTuGD5dLjp4e1q9HcjKuXoWeXmpk5L9JSUv69ft7zBjs2lXekioh+O039OkjPgNOqRcI5UkZLcHq1Th5EjdvAkCfPsjMxK1b1EgWy5YtOHsWp04BwOvX8PLC0aMwMSnvkps3Ub9+eZ3yGICV6CgtjjArC5s3lyxeVyqpqfjyBZ6ecqmrtMuEt2+jWTNJB9vZVYrE0XyWLMG5c7h3r+DQ2hpOTtTnr/fqhexs3L1b7GRSEv79FyNHUqwLwJEjmDIFK1bA3Z0ymd264fbtKnp6aRzOLEvL9vfuYdo0mJpCSwt166J/f6xeXVS4a8WKFS4u/378iFWrxIj98AE/fsiedck6FDpCQ0OsWwdvbwgE4HAY2kdRqDo4GJMm/XB1nTJ4cLq/v/jmgorw+sJOvgxlU9CfQqMbN5KhQyW6ZPNmoqlJRCK59CYlkapVxQS+mYSx0KieHtm0qZTzpYZGg4OJszMFShU/NJrPjh2kQ4eir1ZsLKlWjVD+nwkIIO7uhPwUGvXzIxMnUqyFEHLpEuHxyKxZ1EvO5+flxsxMcvgw8fQkLVoQIyPC4RAej2houALtgFYuLgG5uWKkbdtGxo6ly9RyoOTWEwqJnl7R6i8lwgcOJMuXE0LIt2+katVicXu64XAsAW+gVUKC+MFt2pCICPptKk6J0OiSJWTRIsqEsxYaFQqxYQNmz5ZocHAw7OzkzVAyN0fNmkWv/5WEjAxkZkpRjqfyhEbzGT8eQmHRtvd69WBvT/3LuIcHjh8vCrqKRNi5U66eOKXy6BH69sWoUVi7lmLJhWj9tMlJRwfDh2P7djx8iG/fIBLhyhUYGLwHagK6YWGaRkawsYGzMxYuxIEDePgQWVnFpJ07l876xEJm3r6FiQkMDKiUuXFjQY8IIyPUqLG4ceOOISEnxV8mH5mZWLMGhIiApkBK69YwMkL79vDwwJo1OHsWb98Wlb8BkJiYExsraN+ebrvEwEpXQuodYUgIataUtHz2w4cYNowCpZUwOnr8ODQ1YW4u6fgGDZCQUPKBVYHhcrFhA+bNQ2ZmwZmFC7Fhg9TdBMvHxAQODjh4sODw3DmYm6NlSypVxMejUyfY22PvXirFSkXXroiPv9q69UdHR93v3//48QPHj2PkSGhr48wZeHjAxARWVnBywpw5sLbuHxra7+FD2pw2zcizlb4satXCwoWYOBEiEYmPP/XtW8iSJTRuKszIwOrVsLbG/fvYuHFr7dobd+5cmpSEd++wfj06dEByMjZtQs+e0NdHq1YYPRrTpkU1bNg5J6d9fPxb+gyTBHaaMVE2Bf0vNNqiBTl/XqLxsbGEyyXJyRSovnyZdOhAgRxKYCY06uJCbGxK/6jU0CghpHlzcu+evHqVJTSaz6hRxMen6NDFhaxeTbGKS5eInV1BaLR/f7JrF5XCP38mBgakeXMqZcpM4faJX8nLI69fk9BQsnKlUFPTFjjTseNghs0jFN16vr5k4ULqhQuFpEMHEhRE5s5dUb16ZwODcxcvyiOvdDIyyIYNpHp10r8/efiw4GQ5lqenk3v3yN69xMlpD5f7Py0tjzNnzlBvVrmUCI1++UIMDSkTzs72ic2b39jZSbrmN3s2qVqVGtW5ucTQULpNUfTBjCOsV6/MZZiyHOHYsRQ8ppXLESYkEBMT8vZtweHTp6RaNdk3rZaKSEQaNCCRkYLY2BxjY1LaH15G0tOJiQmxtqZmA6j8lOMIf+b06fPjxs149eoVAyaVgJJbz9mZBAfTIjw6mpiZkZQUQgi5do2YmpIjR+QUWcTPLvDRo2IfSWJ5dnb2nDn+S5euzWM826KEIySEVK0qZpuH5LCzRhgUVG3OHEnX/M6eFdM/SHI0NNCli0Rl1ysM8fHo31+6SyrbMiGAGjUwfTr+97+Cw8aN0bWrpE2UJITDgacnduzg7trFGzWKspr9eXlo3Bjq6njyRHwtWYWif/8+u3evr19+k0wFhsKU0RLY2WHcOMyaBQBduyIsDL//jp07ZRcoEokyMzMzM/HXX6hfH5cu4exZnD5dSmclsWhpaa1a5bN48Szer21EGIeFZUJq3C4hhJDatZ1MTaMl8L6EECIUEjU1cvQoZdo3biQeHpRJkwcGZoSJiYTDKXPyUdaM8Px50ru3vKqVa0ZICMnJIfXqkcLN4k+eEAsLQu1v8Pkz0dYO0NMLfvRIQIlAoZA0alSsco0iIOGMkEXkv/WysoiODin1t6Tkvs7KIsbGm6ysHMPDrxJCXr4kdeqQdetkEfX169datVrr6rY1Mjo/ahR59qzMkYzlscvGrzPCUaPI3r3UCM4f6K4AACAASURBVGdhRpiR0S0zcxQhEpX1vXoVhGDgQMq0V6p8meBg6OtLPfmwsyvWL7SSoKmJVaswYwby8gCgSRO0a4ft26lUcfTo1uzssIyM0J07JcuWLhczsxY8XpPY2I3Pn1Ocu6hCLE+fokED6TqES4WWFhEItr19u2nOnNUAGjRAZCR27MDChdLJEQiwZk18QoJ5VtaIESMeHDiARo1oMZgVmM+XodIRamnd19eHhDPrnTtRuzaVX7j69aGpiWfPKBOoyISHQ4bIU/XqEInYaQDNLkOHwtKyqEb2kiVYtYrK7l3VqlXjchM5nJcWFmZyiuLz+Z8/ZwC7NTS2mskrTIXU0BcXzYfD4Qwe3Etbe1Tt2gW1eSwtceMGrl7FpEnFNjOUBSE4cgRNmuDePbtZswZOmfItIICeOkbsodyh0UaNGt24cUPCwRYWZNIkCpUTQsikSTIGGaiFgUBE9erk99/L/LSs0CghpHt3cumSXKqVLjSaT36aTOEKfLdu0V27LoiOjqZKfmRk5KlTp+SXk5BAOJxt2tptd+8WV82TcSpDaHTGjDLziim8r+PiiIkJef++6ExGBnFwIIMHk/LrqN+4QTp1Iq1bS3cXK11o9NYt0qYNNcJZCI2qqanp6upK4nobNOialNRKS2szhdoBcLkHFi/uvn37AWrFKiApKRgyRJYLK0mH3l9p3BjDh8PXt+AwNtYrIqLrwIHimoRJTLt27Ryo6ELbty8aNfLKyrozbpy4ap4qaICOTYS/UrMmpk7FnDlFZ3R1cfo01NTg5FT6VtdnzzBiBEaNgpsb7txBr160G8kizJcbZSEdLTEx8fXr78Dafft2UCv5xIn1mZmhS5eup1asohEVBUJk7PRWCRNHC/H3x9GjBe8BrVo10tH5OyensSTBKMbYswfPnuHMGbbtqMTI2Y9XcubMwb17uHq16IyGBoKDYW2NXr2K1T2Pj4e3N3r0QKtWePkSXl5KlkUsA8bG4HCo76ddDiz8RU1NawCjuNxZ06ZRXJl46tTxmpq97e3HUytW0Th8GMbGMt4MldkRVqmCRYswYwYAhIYGvXq1q3HjXYsXs23Wf+TkYNIkTJ2KunXZNqWy8vkz8vKk6OgiD9raWLsW06ZBICg6yeNh2zb06IGuXXO3bTt669azefPQogWqVMGrV5g7F5qaTNimCDBcepsFRzhtGkdXd15W1r3Fi+eIHy0N8+dP9fG5a2w8lVqxikZEBBo3lvFaW1s8fw6hkFKDlIeJE/HlC44fB4AaNcyOHEFwMA4fZtssAICzM3R1sWED23ZUYqKjGe2YMXQoatYs2cmLw8GqVahSZdmkSbc7dx79+fPXJ08QEABDQ+YMUwQYjo4y7Qi/f0dQEAICoKlJy7ZNR0dcvEiHYAXixQv06CHjtbq6MDPDmzeUGqQwiESi8eNnNmvmEBVV+jYRHg8bNmDGjIKaq8bGOHYM06axP0u+eRNhYTh6lGUzKjl0p4z+yvr1WLYMnz+XPN+9u46OzueqVUV//qkmeT3hikQFnxGOHg0TE0yZQpf8Vq2QnIyEBLrks45QiK9fMWKE7BIqcHT0w4cPR47ExsTMHTZs1+HD+Pq1lDE9eqBNG6z/bx25WTNs2YJBgxhdkPiVoUPRty+6dWPTBhXMO8JGjTB6NHx8Sp7395974YLX48cXDCrrTtKK7AjfvkVYGAIDaVTB5aJnT1y6RKMKdrl0CTyeXJtn7ewqpiMkBJs31+Jyq1pa+g4e7HrgAKys0L49lizBrVvFosF//om//kJMTIGfdHbGsGFwdS3Ycc8806bhxw9FidBWZph3hACWLsXZsyXbO3O53E6dOplXzskgAMa3EjLqCIcPR6NGGDyYXi329ggPp1cFi5w8Ke9ivo1NzrVrt/PYeurTg1AIb29ERPDevdsTHx+5dm2HU6eQmorNm8HhYNo0VKmCAQPwzz/48AE1a8LY2LNt2xHOzgXbJwICoKmJefNYsDw+HoGB+PtvyoqUqpCNL1++PX36oEkTpvXq68PfH3/8AUKYVq3IVNg1whs3EBWF/ftpV+ToiPBwiWo0KCO3bqFZM7kk+Pk5RUbud3CgOGWXRXJzMXw44uJw9SqMjYvO83ho1Qq+vrh/H7GxGDECERFo2xZNmuDNm39zcxecPn07/+nD5WL/fpw6haAgpo13cEDTptT38lUhFVlZWU2adM/NXbdixQrmtf/2G0QiJp6NSoSJCTgcpKYypI45Rzh6NLp1k6UsurRYWsLUFFFRtCtihTdvIOe+bT4/Uyhskpr6gyKLWCY9HU5O0NBAaCjKKedgZoaxY7F/P5KSsHcvdHS0gSCRaKGbW0H+epUqOH0a8+fj3j3GbMeuXYiNxalTzGlUUSoCgSAnhwvUS039xrx2DgebN2PuXPyoIDclNVhbM7dMyJAj3L8fCQk4wFTJFweHipk7mpGBzEw4O8sl5N9/D9eqpb5y5T6KjGKTb9/g6IiaNbF/PzQ0JLqEy0WrVrh6de+oUdWrVbN58QJ9+yItDQBsbLBtG4YPZ6gca0YGJk/G7NmoXZsJdSrKwdDQ0MUlqH//en/+6cuKAa1aoW9fLF/OinIFhcnoKBOOkBD8/jtGjoSFBQPaAMDevmI6wuPHoakp75+xVq1avXt7JiQofUXnpCR064auXbFrF9TUpLu2ZcsWBw6sunat2ceP4HLRq1eB8xs0CG5uGDoUfIl6qMjF4MEwMEBAAO2KVEjCp0/N3dzcJCkSSRMrV2L3brx8yZZ+hYPJxFEmHKGfHzIzKe6GWj7du+P+fWRkMKeRGc6dQ506FMipABVH375F164YOxYBAZI2gv6V+vVx7hxiYmBnhw4dCl4/ly6FuTlmzqTQ2FIID8eVKwVb+1UoAgzvpv+VatUwfz6mTWPTBoWiQjlCPh8BAZg9G9radKsqQkcHbdvi2jXmNDLDgwdo04YCOc2aKbcjfPIE3btj7tyi1vMy06wZTpzA6dNwdUWPHnjwABwOgoLw778U9ywsgasrBg9G5840qlAhOWlp+PIF1tYsmzFtGpKTERrKshkKQoUKjXp7Q0MD/v506ylJhYyOxsejf38K5OQ7QiVN146IQM+eWL8enp7UCOzQAQcOYMcO/O9/6NsXYWHQ18fx41i4EJGR1Kgogbc3cnMRHEyLcBUyEBODpk3ZL2atplZQ+Sgnh2VLFAEmtxLS+5//+hX79mHNGha+YRUvXyYxEXw+NY6walXo6SEujgJRDHP2LIYNw/79GDaMSrH29ggMxLp12LoV7u44eBANGmD/fowciZcv0wmlrwyxsdixA//8I2l2jwoGyA+PKwI9e6J5c6xbx7YdCoCpKYDS60NRjpQ5BlLi6orq1eFFWdM3KWjeHN+/Iy4OtWqxoJ0OgoOhr0/Ztms7O0RHK02+okAg+OefPU+e6J44MTI0FO3bU69i6FCkpGDOHBw6hHHj8OoVfH3Rvn2wre3ftWvj8eMrWlpalCjq2xctWmDUKEqEqaAGxXGEADZsQMuWGDNGaW5P+rC2Rmws2rWjXRGNM7WXL3H5Mnbtok9DeXA46NWrQk0KL11C/fqUSVOuZcL9+w/NmPHwn3+OrVr1Lx1eMJ9Jk+DujmnTcO4cTp7E9OnQ07srELh/+KD28iUFO3tXr/5LTc3m7duZZ8/KL0wFlSiUI6xZE6NGvenYcfyKFX+zbQvLMJYvQ6MjHDECdnbo3Zs+DWJwcKg4tdaeP38eEeHToEHpTRVkQLkSR83N6wiFj4yM3nboQG+zuAUL4OgIDw+cPYuoKHz5MmfEiEcDB3oNGFDjwQN5hfv5HRIKjwCXgRQqjFVBDYTg6VMWqoyWQ3Lynx8/Dlqz5mRiYiLbtrAJY/kydDnCK1fw5AkOHqRJvEQ4OODSpQrSe6916wFZWe2PHRtKlUDlmhHevNnVxeXUmzdXGjRoQLeu1athZwc3N5w8ibi40BMnYm/dWrdyZW6/fnLtdhg7FllZvlyui7W1gZmZ0m/irEi8eYMqVWBkxLYdPzFihKO+/nI1Ne1q1aqxbQubMFZchi5HOG4c7O3lapIgP+bmsLTE/fts2kAdusB2HnU9HBs2RHx8QVs+BefTJ2zbhoCAakaMPKs4HAQGwtQU48ejXbs3QmHP5GTN1NRvFy5g1izMmyd1GduUFNSujWPHcPlyH6HwWWwsPamoKmQlJkbe+r2U4+w84M2bO3l5ZzMy1Nm2hU3q11dmR7hvH5KSFKKGbIXJHXVxuaij88eDB2eoEqimhoYN8fQpVfJoZOlSuLkxmvTE5WLfPggEyMhYsGyZ9ubNy4ODzX//HXv34tYtuLpK8QJx7hxq14aaGuLiZG+nrIJWFGqBsBBTU66jI8tBNdapVw+vXjGhiBZHOG0afvsNJiZ0yJaOCrNMGB1t1qVLj0aUTrHzE0cVnHfvcOQI5s9nWq+6Oo4exbNnT1atOnL27L6ICKGnJ1xcYGsLDQ107IgPH8QLmTwZAwbAxQVv3ijE7aCiVBTTEQKYMIHRglwKSLVqEImY2EFBvSNcsAA5OdiyhXLBstClC6KjK0JN91ev0LcvxTKVokPvwoWYPr1YcyXG0NZGw4YhP374h4V9Xr8+wdkZT5+Cx8O//6JpU7Rrh4iIMq9NS0Pz5ti5EyEh2LOHQaNVSI/COsKePZGVxWg7FAWEmcRRih1hXh53/Xr4+CjKZmEtLXTogKtX2bZDPr5/R2YmXFwoFqv4M8LoaERE4I8/WDNg/nx3W9uVvXrZ3L9fq1YtzJuHceNw6hRiY2FmhqFDsXt3KVfdvo0aNZCSglevKN74r4Jy0tORnIx69di2ozQ4HLi701vqT/FRSke4eHFNHR34+FArVS4qQMP6w4ehrQ1zc4rFNmuG6GiFLrQ2Zw58fMrrMkg3LVo0j4m5GBa27vBhztOnsLLC8OHw9oa7Ozw8IBJh5kxMnlwsM9nXF506oUsXJCaWsiH6RwWITlQsHj9G48agMA2NWtzdcfQo0tPZtoM9mNlBQaUjzMkxDAvT/1vB9oA6Oip9vsyZM7SUAzY1hbY2EhKol0wJ167h9Wu4u7Ntx3+Ym2PuXMTGIiAAly9j0SJ07gw7O+zahRYtEB397saNe927Y9kybN6Mc+dKlhUUCuHjs6VOnf5Nm3YXSZt4qoI2FDBl9GfMzNC9O0JC2LaDPZiZEVJZYu3Dh6EiUetRox4w2fheLE2bIicHb9/CyoptU2Tl4UMMGkSLZDs7xMSgZk1ahMsDIZg3DytXKkqMvRAuF717o3dvJCdjzx5s3w5zc7x48U+bNssAYw2Nlffu9dHUxPnziItDfDzi4vD+PeLikJQEDudObu7UrKxVX75kV6vG3jxXxU/ExCjWVvpfmTABvr6UlZhXOurVQ2AgH6D3QUClx+LxbIBvn5hp7y0NvXsr8aSQECQnY8QIWoQr7DLh8ePg8zF8ONt2lE3+BPHVK+zYAU3N9UAroKtAsKJTJwwbhr/+woMH0NRE797w98fVq0hPx6tX/h4eD7p0WerkpJuczPYvoAKAws8IATg6IiUFUVFs28ESO3YsvHevm6vrFHrVEOowMHDicBru2iWiUCYlHDhAhgxhTl1aWhqF0m7eJFwuEQqluCQvLy8zM1OSkfv2EVdX6ewRCARZWVnSXSMlAgGxsSEXLtCqhEquXbvG4ZgC5suXb5ZkvL8/qVuXPHtGt10Uk52dzefz2baiPKS99UQiYmhIvnyhRTiF+PqSqVNlv5xFyyVBJBJlZGSU9amVVXvgobl5C9mEC4VCSb60VM4Ia9X6EBp6ctUqzvTpUlffoBUHB1y9CoGAbTtkIjgYZmZ09bFSzBlhUBDMzODgwLYdEtO1a9dv32I/fny4YMFkScb7+GDVKvTqhX//pdkyFeXy/j0MDNjZnCMVHh44dEg56kBRzt69642Nt8ycuZZWLVQ+XzkcTs2auTdvIiYGw4cr0L/NxATW1rh7l207ZOLaNbRoQZfwRo3w/j2ys+mSLwM5OfD3R0AA23ZIiY6OTpUqVSQfP3w4jhzB6NE4cIA+o1SIITpaQXcQlsDSEu3a4dgxtu1gg06dOri4bNfU7EmrFuonGlWr4sIF6OujUyfEx1MuXkaUt8TM69dwcqJLuLo6GjTAs2d0yZeBv/5C+/a0dBxUNDp1Qng4fHzg68u2KZUVxV8gLGTChMq7odDWlvbSH7RE3DQ0sHs3fvsNHTtC/uY1lGBvjwsX2DZCer59Q1YWvTkj+YmjCsL371i/HkuXsm0HUzRujFu3cOYMPDyUNXSv1Ch+ymgh/fvjzRvFemdlDGV1hPlMn44NG9C3L06dok+JpHTqhOfPmahZRy2HDkFHB6amNKpQKEe4ciWGDGG5aQnDmJvj2jV8+gQnJ6SlsW1NJUOJZoRqahg3jrU+5+xiZ4enT+nNO6F3w5+zMy5exLRp7Ad/NDTQpQuuXGHZDGk5e5b24k+K4wg/fsTOnYpVlogZdHVx8iTq10fnzopb36DikZmJxETUr8+2HRIzYQL27kVuLtt2MI6+PoyN8e4djSpo3/nevDlu3cLp0/D0ZDn4o4y11h49Qrdu9Kpo3lxREkeXLIGXFywt2baDDXg8bN4MDw906ICoKORWwqcd4zx5gkaNoEZlTRF6qVsXdnYKEWBjHro7BDBRAqZGDVy7hs+f0bcvvn9nQGHpODgo2TIhIUhJoWsrfSHVqkFNDYmJ9GoRy6tXOHUKs2ezbAa7TJ+O1avRtu0Qc/NuK1ZsZNucCk50tNLERQuptCkztrb0Bq4YqoWmp4fjx9G6Ndq2xcuXzOgsiY0NOBzWtMtAZCQAdOhAuyJFiI7On485c1C1KstmsM6IEUJt7Tffvy8LCVG2OL6y8fix0mTKFDJkCGJiGGrarlDQnS/DXFFQHg8BAZgzB926ISwsM4aNR69y1VoLCYGFBV1b6X+GdUd47x7u3sVkiTajV3B4PN727X6dOx9LSVkhSe9fFTKjjDNCDQ2MGYOgILbtYJyK4wjz8fREUJDAyalzp07LZs/2Y1i7cu0mjIxEy5ZMKGLdEc6bh6VLoaPDpg2Kw4gRgyMjA+fPb9SvX0XoKa2YEIInT5RvRgjAywtBQZVus03DhoiPp7FIi3hH+O3bt0OHDi1ZsmTjxo3JVJQK7tZNYGhIMjNbvn6dIr80qejdG9euKU3a1Zs36N+fCUUsOsLo6JjWrd2ePdvv5saOAQrL9Ono2RMuLsjLY9uUikhcHHR0YGLCth3S07AhrK1x7hzbdjCLmhrq18fz53TJF+8IXV1dQ0JCeDze3bt3GzZs+EzuLZ06OjqXLwf17m1sbLxSTlHSUqUKGjXC7dsMq5WFT5+Qnc1Qf/PGjfH2LXJymNBVgqlTlz54MCEzcy0hlewVVwI2bIC6uipiTAsxMcpRXK1UKmfKDK2Jo+Id4dGjR0+ePLl48eJ9+/bZ29vv3LlTfq0tWrQ4fHjC6dMGzKeuODgoxzLhoUPQ1WUoeURDA9bWNL5tlYONjZOa2sy2bW3U1dVZUK/Y8Hg4eBB37kDRml1XAJRxgbCQ4cNx5w7i4ti2g1loXSYU7wj19fULf+ZwOJqampQoNjLCH3+wUExLWRzh+fNo0IA5daxERwlBVJT70aO3L10KZlq3kqCvj9BQrFqF0FC2TalYKGPKaCHa2nBxwe7dbNvBLLTuoJBiN2lkZGR4ePjatWW2w8jOzl65cqVp8YJggwYN6ty5c6njJ02Cra3m7duC5s2Za9rUvDlev9ZMSOCbmBCaVOTk5Mg/v4mK0nR1FebkyLJAJBQKc3NzudLkmzZurPbwIVxcxKvLy8sTCAQcDkcGw0oQGsrLy+M5OAhycpQ+LioQCIRCIR2Szcxw9Ch38GB1c3O+nR1dX1pJyMnJ4fF4NP2alCD5rRcVpTl3Lj8nR4q/JyX3NVWMG8cZOlRj1qxcHk/8YIWy/FcIIflfrfKHNWjAiYnRyMmRLsVDJBKJlQzJHeGTJ09GjBixe/fu2rVrlzWGx+PVqFGjevXqP580MjIqyw49PcyZI1q6VD00lLl8AB4PXbuSa9fURoygy/vyeDxJ/vTlIBLh82eOiwtkliOtDXZ22LCBx+OJfy7kN9KU8xcEQAhWrFDz9RWqqckrShEQiUSA7P+v8mndGhs3ioYO1bh+XVijBmu+kPcfbBkgFgnNy85GQgLHxoYr1a+iUL978+YwN8fVq2qOjuK/Dwpl+a8QQiSxsGZNCIWc1FRetWpSCOdwOJK8tUvkCF++fOno6Lh+/frBgweXM0xDQ8PNza2ZNKH3SZPw99+4fl29Rw/JL5IXR0dcvozRo+n6Zqirq8v5/nXlCjgctG+vJtu8i8vlikQiqWxo2RIxMZDkkvxvlfwvmMePQ00NgwbJ+DsqIFwul7737hEj8Po1nJ3VIiKgq0uTEjEIhUIejyfJ75ibm0vVGopUSHjrRUWhYUPo6Ej3z5L/vqYWLy8EBalJklguoeV5eXkcDod5l0kIkdDCpk3x4oV6jRpSCBeJRJLEMMRHz16/ft27d29fX9+RI0dKoV8y1NWxdCnmzQNh8DW3b19cvMioRmkJCUH16mDSQ+Sro2J3jEQQAj8/LF3K6O+o7CxYgNat4eZGbxl++enWzbl69W7Llv3FtiFlotQpo4WMHIkrVxLOn79HibTHjx9bWravVavdBwWu40Bfvox4Rzh+/PisrKzDhw/b29vb29v7Ut1IYuRICASM5gLUrQstLYXu7HXjBlq1YlqprS1z1bfzp4P9+jGkrsKwaRO+f8eCBWzbUTYikejRo1dfv67w97/i7Y3QUGRmsm3TL1QMR/j9e0JOzkBn54DAQAqaM92//+DTJ8ePH9v06fN4+XI8eqSIUwX6HKH40Oi6devSfuqTVk2qAK0EcLnw98fs2XByYq4SfH7uaJMmDKmTlrdvMWsW00rzE0cdHWlXRAj8/bFihWo6KDXq6jhyBB06oE4dTJzItjWlweFwjY2X2dkdWbRo5YcP2L4dY8eibVv07o2BAxWl02R0NAYMYNsIueHz+VpaWj9+mH79miG/tPr1R+jpLffwqD5okGNYGMaMwZcvcHTEgAFwdISBgfwaKMDODjt20COaUIetrW1UVJRs13bvTnbvptAWMRw/Tvr0oUt4WlqaPJd//Eg4HPLjh+wS8vLyMjMzpb0qKIiMGSN+mEAgyMrKksWs/zhyhLRqRUQieWQoHHw+Pzs7mxldb96Q6tXJsWOfv3z5wozGfLKzs/l8fvljQkJI69bF/rmZmSQ0lHh5kRo1iJUV8fIioaEkN5cQQr58+ZKXl0ehhRLeesbGJCWFLuFMcu/evVatTmzeLCh/mCSW29uT7duLnXnzhmzbRvr3J7q6pFMnEhBAnj4t+ChFhj9f2YhEooyMDElGpqcTHR0i1VdGKBSK/dISQhTFEd6+TWrXJjk5FJpTHj9+EAMDIt/zvEzkvGHWriX6+nIZIJsjfPCA2NmJHyanIxSJSLNm5Nw5mQUoKEw6QkLIP//c53JbV63a+uHDh4wpFesI8/JI48bkwoXSPxWJyP37ZOlS0rYtMTIirVrt1dfvWqdOGznfq35GklsvLo5YWNAlnHnu3ye1ahW8WJSFWMuvXSPW1qSs/216Ojl+nHh6kurVSYMGxMysT5Uq9tOmLZLV5JJI7ggJIVZW5NUrKYRL6AiZLrpdFu3awc4OW7cypM7AANWq7fHyWp1FXxlXWQkLQ8OGLOht0gSvX4PPp1fLsWNQV0efPvRqqfCYmsarqzf+9s3mypV4tm0p4tAhGBvDwaH0TzkctGqFxYtx5w5evoSubnR6+vCUFPVv374xaWTFWCAspFUr2Nhg/365hCxejMWLUVbapp4ehgzB9u1ISMChQ8jISPj2bfK1a4/kUikrNG2rVxRHCGDFCqxciZ+WI2nkzp07CQknDx5M27KFgopx1BITg549WdCrqYk6degttJa/Oujvr1odlJeBAwdu3txrxgyHtWv7K0j9ZaEQy5bBT7KOMtWq4dixeRMnpqqr/+/Hj+riL6CO6OgK5QgB+PpixQrZi7NfuoSkJIwaJX4kh4OWLXHuXGCXLve1tf9iJZuGpnwZBXKETZvCwQEbNjChq2bNmgYGH4DLtraNmdAnMXl5+PwZLi7saG/WjN5Ca0ePQkdHNR2kAC6X6+Hhtm7d2NBQrqcnbRkE0rB7Nywt0b27pONNTEwCA5f4+w+eOJHRBMXHjyuaI+zQATVrIljWMoWLF8PPT4pExa5dO1+9uozPtwoJkVGjPFR8RwjAzw9//41Pn2hXVL169Q8fbrZpc14g6EW7Mmm4fBlcLlq0YEc7rdX8CMGyZaB6901lp00bREZi9WqW/7ACAZYvx5IlUl84ZQrS02V/iMtAxZsRAvDxwfLlsmwwPXsW6ekYPly6q3g8bNuGWbNY6JdJUw8KxXKEdepg5EisWsWELi0tLW9vIyp6aVDJkSOwtGQtckjrjPDIEejqMrE9o7JhbY2ICISGYto01vba79iBRo3QpYvUF+Y/UmfPxvfvNJj1C7m5+PABNjZM6GKSXr1gYoKjR6W7ihAsXQp/f0hTlriANm3g5MTC61f9+vj4kfrNqYrlCAEsXow9e/DuHRO6RoxARASSkpjQJSE3b6JNG9a009eDQiRSTQdpxNwcERGIjcWwYcjOZlp7Tg5WrpRlOpgPk4/UJ09Qvz40NJjQxTALFsDPT7o3oZMnIRRi0CAZNQYEIDgYUVEyXi4bPB4aNqS+HIrCOUJTU0yahGXLmNClqwtnZ+zbx4QuCXn7FgMHsqbd0hICAVJSqJd85Aj09MrMJ1QhP3p6CA2Flhb69WNodlXI1q1o1Qpt28ouISAAISF4RH8eYgVLGf2Zvn2hoyNFiS6RqGA6KHP8qWpV+PvD25vpOAQdy4QK5wgBzJmD8+cZKoHm4YGdOxWlmFBCAvh8DUrulQAAIABJREFU2V/QKIGOZUKRCMuXs9B7srKhoYEDB9C2LTp3RjxTuyqys7F2rbzzuapVsWwZE4/UCuwIASxcCD8/SZ9mhw9DQwN9+8ql0cMDamoICpJLiLRUFkeor4+ZM7FoERO62rWDtjauX2dCl1gOHICBAcvVjOhYJjx8GPr6sLenWKyKX+FwsGoVJkxAx440tvP+mY0b0bkzBd3e3d2hro5dFFTNLI+K7QgHDgSHg7Aw8SOFQvj5YeVKedMROBxs3YqFC/Hli1xypIKOl3VFdIQApk7F/fu4dYsJXb/9BgVJmblwgf1lfMq/ZKrpIPNMn46AANjb0/6Gl5GBP/+UfXXwZ/IfqT4+9D5SK7Yj5HAwb55EWzkPHICxMXpRkTJvawtXV0YLwVciR6ilhUWLMG8eE7rGjsWpU0wvq5TK48fUfDXlgfIZYUgIDAzQuzeVMlWIZfRoHDwIJ6dzFhadZ86UbJe79Pz5J+ztKSulbWuLkSNpfKQmJoLDgbk5XfIVAWdnZGTg8uXyxggE8POjMg/D3x/nzzM0bwFgYUF9zzgFdYQAxo/H58+4eJF2RcbGsLcHK5tDf0YgQGoqXF1ZNqNJE7x48frbN2rSk0UirFghabURFdTSsycsLbckJx/YuvV4Vhb1i28/fmDjRixcSKVMPz8aH6lXriQ0aZJKi2iFgcvF3Lnw9y9vzO7dsLJCt26UKdXXx+rV8PaWvbqNtDRtSnHkX3EdIY8HPz8sWMBEJkt+ygy7hIWBx4OtLctmbN26VSD43camU3p6uvzSgoNhYMD+NLfS4u8/sXbtUbVrD27WjBsRQbHwtWsxaBDFdXHpe6SGh1+aNGnk3bv2r169oli0gjFyJJKSUNa/m8/HypXUv5uOHIlq1bBlC8Viy4LybfWK6wgBODsXNGCjG3t7fPrEXFvaUjl+HDVrsmlAPm/fxhPSITNTKyND3iZnQmFBZVEVbDF0aP/37288f+67fj1Gj4a3N2U7kVNTERhISxiTpkdqfHxidnZTDsfyEwOVq1iFx8PcuWVGPrdvR5MmaN+eer2BgVi2DB8/Ui/5VyhPHFVoR8jhYM6cZHf3ce7uswUCAX2KuFz89hvTScAluHlTrp1YVLFy5fwZM4wNDVebmVnII4fP569eHWNsTFgpIK6iBAMGFDw47Oxw9SoFAlevhosL6talQNSv0PFIbdt2tIFBu1273Dt37kylXIVk7Fi8fo0bN0qez8lBQAAWL6ZFaf368PTEnDm0CC8B9fkyUnR2Eoc8/QjLYtmytRxOoI6OR0REBLWSSxAXR0xNCSUd5WTrW6apSQ4dokA7kbUf4c+0aUPOny/9Iwn7EdrYdObxxvbsOVkeM5QIhvsRyszZs6RWLeLlRdLTpb62sB/hp0/ExITEx1NvXiELFpBRo6S+qpxbb+pU4usrl0mK2Y+wLAIDSf/+BT8XWr5uHRk6lEalWVmkbl1y+bJ0V0nVjzCfjAyiqytRh14l60dYFv379zY3383nv69bl97Vs5o10aIFTp2iVUmZvH0LPh8DBrCj/VcmTkRgoFwSkpLSRaLOfD4jgRIVEtOvX8GrtJ0dLl2SUcjy5Rg7FpaWFNpVEh8f3L6NK1eokZaVhZAQeHpSI00pcHfH48e4f7/oTGYm1qyhazqYj7Y2/vwTkyYhN5dGLQB0dWFhgdhYygQquiNs1qzZx4+3XV0v7dhhRLcuFlNmDh2CoSF0ddnR/isjR+LmTXz4IOPlqalQUzs4axb/2LFtlNqlggIMDbFtGwID4ekJNzdI2xY3KQkHDmDuXHqM+w9tbWzaRNkj9cABdO6MGjUoEKUsaGhg5kysXFl0ZuNG9OhBQemD8snPn/rzT3q1gOplQkV3hPmsXIlNm5CQQK+WwYMRHY03b+jVUioXL6KxIjVG1NbG6NGyvxbMmYPffmu8Zs3UatWqUWqXCspwdERMDKpUga0tAgJuLFu27otkW9n9/eHhATMzug1E375o2BDr11Mgats2TJpEgRzlwssLd+4UeIuMDPz1F73TwUI2bcK6dbQ3TqA2cVQ5HKGlJTw8aK9OoqGBUaOwdy+9WkrlyROF23I+cSJ27IAMKUr37yMsjKEKeSrkwcAAf/2FbdsyfXymLF6sYW+/4N498PnlXRIXxwkJwaxZDFm4aRPWr5f3kXr3Lr5/r4x7eLS0MH06VqwAgHXr4OjIUOGqWrUwfTpmzqRXS2WcEQKYPx+nT+PJE3q1eHhg1y4IhfRqKUFuLr59Y38rfQlsbGBjI/WiqUiEKVOwejUMDekxSwXVODpqWlqq6epe5HCsvbxQpQratsWUKQgKwpMnJe+FZcu4U6fC1JQh22rVwh9/YMYMki1Hc6nAQEyaJEvLvQrA5Mn49188fMjbsoXRd9M5c/DiBY4dyyK0bQOvpI7QyAgLFtC+MtG0KapXR3g4vVpKcPYs1NQoq1NFITKkzPzzT8HEWoWyoKam9urVzejovx4+nPvoEb5/R1AQmjXD1atwcUGVKujcGdOnw8fnuoFB8z17nIYOfc+keZMmZZw/38nCovvx46dluPz7d5w6BTc3yu1SDnR1UbXqhB492nK5/7O2Zk6vhgaGDTvh4tKzYcPOWVlZdKiwtkZyMqgo+wEokSMEMGkSYmNp91LMp8ycOIFatRjVKCFDhuDFCyn6YX39Cl9fbNokb0l7FQyjoaFhZWWV/7O6Opo0gZcX9u7F06eIj4efH6pXx549WwWCKYTYnDy5h0nbUlOTtLT0f/z4LTT0tgyXBwWhf3/mprAKyKtXYYTc+PTpMMN6P326JRR6xMVp0VTBgMdDo0aUdetTJkeoro4VK/C//9HbtGzUKFy5gs+faVRRgtu30a4dc+okR10d48djxw5Jx8+fD1dX2tPSVDCJoSF69sTcuTh9erae3mZj40sTJ05k0oD69esvXz6iXbu3Dx78IUOxox07KmOazM94eIxUV28+cCDTOQh+frPc3F5aWIw9fLgOTSqo3FYv1TbG8qFjQ/2vdOpE9u6lV8Vvv5F162S/XNqNtxoa5OhR2dX9ivwb6guJiyMmJuRnYWVtqL9/n1hYkG/fKFGrfCjLhnp5KNxQzwre3mTwYCIUljemxK13+TKxs6PMAOXaUP8zLFr+8SOpVYucPFneGBk21Oezfj2ZNk3MmAqyof5X1q3DggWgJ+xcgIcHtm+nUf7PxMZCIICTE0PqpKVmTXTogMPiwioiEaZOxcqVMKJ9t6eKSsrGjfjxQ7qMj8BATJ5Mm0EqJMDCAidOwNub+g6CoHQHhfI5wnbt0K4d/v6bRhWdO4PDwW1ZliSk5uBBVKkCLS0mdMnGpEniU2Z27QKXW3lTElQwQH79/SNHsH+/ROOTk3HlCkaOpNksFeJo2RKbNmHQIKSkUCzZzo4y/6p8jhDA6tVYt476P+vPjB/PUMpMeDiaNGFCkcw4OuLzZzx4UOaAb9/g44MNG1Q5MiroxdgYoaGYPVuihoXbt8PVFQYG9JulQhzDhsHNDUOHUlx6zdQUamrUFGdXSkdoZYVRo7B8OY0q3Nxw7Bhlubnl8OyZwm2lLwGXCy8vbCu7VpqPD4YNQ5s2DNqkorJiY4M9ezB8OOLiyhsmFGLnzspVXFTB8fVF7drw8qJYLFXRUaV0hAB8fXH4MOhrsWlmhu7daW+FmJuL79+VIHrj6YmjR0svSvnoEY4epb3ojwoVhTg6YtYsDByIcpJIz56FpSVatGDQLBXlwuFg5068eIHVq6kUS1XiqLI6wipVMGMG5s+nUQUDGwp37YpXV/9avz69WuTHxAR9+uDAgZLnCcEff2DFChgbs2GWisrKjBlo3x5jx5a5kyq/mowKhUJbGydPYvNmKpv8UFVfRlkdIYDp0/HgAa5fp0t+nz6Ij6dsw+av7N9/cOrU4QJBlyd0F46jgkmTsGULStRL2rMHWVkYP54lm1RUYvKTSH18Svno7Vs8fAhnZ8ZtUiEOypNIVY4QWlpYtgyzZpV8OlMFjwc3N+zaRYtwACdPnheJWhFS/QpVXdfopEsX8HiIjCw6k5YGHx9s3lxJqziqYBd1dRw9iqNHS0ki3bYNv/2m0JnYlZmWLbF5M2VJpE2a4NUr5OXJK0e5n2GjR0MoxNGjdMn38MDevXQ1mbS2XsHh1LCySh83bhwtCqjG2xtbtxYd+vigf3+0bcueQSoqN1WrlpJEyudjzx5VmoxC4+yMceOoSSJNS0vhcr29vZeJ5Ks3ptyOkMPB2rWYO5cuX1W3Lpo2xZkz1EsWChEYWNPdfcGbN7cNlaRTw7hxuHix4D3u6VOEhGDZMrZtUlG5+TWJ9OhRNG8OxV93r+QsWYLatTFhgrxyNm7cmZnZMiTkyYNyNnhJgHI7QgDdu6NxY6mbJEgOTSkzy5YhN5eJPs4Uoq+PoUOxZw+XEEydCj8/mJiwbZOKSk9+Emm/fkhLA1RpMkpCfhLpy5dYtUouOX36dDcy2gkkNWjQQB45Su8IAaxdixUr8PUrLcKdnXHvnphNS9LC52P1akyeDH19KsUywJQp2LaNe/AgLy1NFX1S8X/27juu5v2PA/jrnE47qaxIURGiZbTL3kS2ECIRl2vP62ZcK1zr2rIpbkbIHmW0jJBkZoZco61xzvn9cfyM1Jnfcc7p83z8Ho/f7Xu+38/nTc55n89WFhMnwtMTAwciNZX7/LnyblhI/Eg0iXT9eoUmkXp6umdkxFeqFPPunUL9auqQCBs2hK8vliyhpXAdHfTvj52UnjwzaxYEgq8nR6sWQ8P0V6+sAwPdu3WL1NBgOxqC+L+1a5GQMNjV1SU7uyePx3Y0hHRq1sSePdm9etU3NGy8YIGcbUMdHe7QoYruDq0OiRDAggXYvh1PntBSeL160fPnu/r6UtMCysvDunWYMgW6upSUx6iTJ08Khc2BYdeu7WY7FoL4TlMTubnXgAM5ObfYjoWQQUnJDcBUIFi+Y8dpuQsZNQo7d+LLF/nDUJNEWL06xo3D3Lm0FH7y5P6SknXnz6dmi0YhFDN5Mng82TbRVx4BAQGNGr2vUSNs9WoVbM8Sam3FihmGhv2GDevFdiCEDFq3bu3mZqKnt4nHW5+XJ2chVlZwcsKhQ/KHoSaJEMC0aYiNLQoLu1ZI9RTSefN+a9hwDp/f7eNHRXfwzcpCWBjmzoWWFiWhMU1HR+f27TPp6VdtbW3ZjoUgfjJ27KhXr+K2bVOpGWgEcOXK0czMI61aNRwwAHy+nIUEBYnbD1ki9UmEurrQ1PQdPTrC25vir4Surs73758KCZk1eLD8vyeRMWOgp4fJkymKjCAIQi2sW4eCAsyYIefjPj548gT37sn5uPokQgA6OnklJTbPnpW/F68CJk+GtjZCQ+Uv4f17HDyIRYtAppkQBEH8SFMThw7h1CmsXy/P4zwehg3D1q1y1q5WifD8+f2hobo8Xnh0NPWFc7nYuRN//42kJDlLGDkSxsZkkRNBEEQZDA0RFYWFC+Xcw2TUKOzejfx8eZ5Vq0RYs2bNyZMD/v23ZkAAHj+mvvzatbFyJYYORUGBzM++fo0TJ7ByJTm9liAIomyWljh6FCNHyrMrt4UFXFzk3HFTrRKhiJsb5s6Fjw+omONZ2qBBcHTE9OkyPzh0KExNMXgw9SERBEGojRYtsHYtevbE27cyPyv3lBk1TIQAgoPh4YFhw2g5mGLDBkRFQabe16dPcekS1q6lPhiCIAg107cvhg9Ht26QdUFF16549QrJyTLXqJ6JEMC6dXj9muLTkEUqV8bu3Rg5EpmZ0j4yeDDq1IGvL/XBEARBqJ85c9CkCfz9yz17uUwaGhg+XJ4pM2qbCLW1ceQI1q3DyZPUF+7lBT8/BAVJdXNKChISFFrjQhAEUaFwONiyBVlZmDlTtgdHjsT+/TI3JdU2EQKoWRMRERg+nJat1/76C+np2L5d8p1DhqBhQ7RrR30MBEEQ6kpTEwcOICxsnbGx89Spf0n5VO3a8PREeLhsdalzIgTg7o45c9Crl8xfECTS1sa+fZg2DY8eibstMRF379JykBNBEIR6MzGBtvauz59Pbtp0UPo+UjmmzKh5IgQwbhxatIC/P/UTZ2xtMWsWhg0Tt93M0KFwcICrK8VVEwRBVAQrVsxu1GhA7dpT27fHmzdSPdKpE96/x82bMtSi/okQwLp1ePUKy5dTX/Lvv6NyZSxeXParFy7g0SOEhVFfL0EQREXQv3+P1NSzd+8O8vJCs2Y4c0byI1wuRo7E5s0y1FIhEqGODo4cwZo1OHWK4pJFI7rr1iEhoYxXR4yAmxscHCiulCAIokLR0EBICPbuRUAAZsyQvOfziBE4cECGpeQVIhECqFkT4eG0TJwxM8P69Rg8GLk/b3F6+DBeviTNQYIgCGq0bo2bN5GcjHbtkJEh7k5TU7Rujf37pS25oiRCAB4emDmTlokzvXrB1RXTpv10cdw4tGmD+vUprosgCKLCql4dJ0+iZ080by5haVxQEDZskLbYCpQIAYwfj+bNERhIfcmiBYvf9ordtQvv3sm/FTpBEARRJg4HEyZg/36MGoUJE3D//qMPHz78elv79sjPl/aMhIqVCAH88w8eP8bKlRQXK9puJjAQ795BKMTkyfDxgYUFxbUQBEEQAFq2xM2buHgxwt5+bP36bZ89e13qBg4HI0Zg82apTjmocIlQRweRkQgNxenTFJfs6YmhQxEUhOXLsz5/lqFVThAEQciqWjUMGvQOsM3JMW7SJKt9e/z1F65dQ3Hx1xuGD8fhw8L//isWWwxQARMhAHNzbNr0qVu3NmZmLsly7M9avnnzBFFR5gsW9DcymlSjBoUFEwRBEKVNnDh606amJ0/OefvWdvp05ORg2jRUq4b27RESgvj4nLy8Fps3p0gsh8dArErI2DhFS8syI6P5xInnIyIcq1enptj7958JhRrAxA8fpgBUd78SBEEQP9DS0goI8Bf9d7t2X3ey/PwZsbG4eBGTJsUXFWkLBDYSy6mILUIAbm5ufn5VWrdOrFNncKNG+O03PH+uUIFPn2LCBLRtW5fLLQamurs3pihSgiAIQgZGRvDxwd9/4/Zt92rVsurUkXzIbwVtEfJ4vC1bvh7RlJmJ9evRrBnatMGff6KxjCnsxg2sXo0TJzB4MJKTuebmr1+8eGFB5skQBEGwSl9f/82bO0VFRRLvrKAtwh9Vr46QEKSnw8MDHTuie3fExUl+SiDAsWPw9ES/fmjWDC9eYPVqmJsDgLGxMd0xEwRBEBJxOBweT3J7jyTCrypVwoQJePIEffti6FB4euLYsbL36c7NxebNsLXFvHkYNQoPHmDCBOjrMx4xQRAEQQWSCH+irQ1/f6SlYfp0zJsHR0fs2oWSEmRnZwN4+xYhIbC0xLFjCAvD9evw94cU3zYIgiAI5UU+xcvA5aJ7d3TrhuPHsWQJJk/e9OXLbh6Pp6Fxwd+fe/066tRhO0SCIAiCIqRFWC4OB9274+pV2NjE5+b+VlSUc/du/sqVJAsSBEGoFZIIJdu/f15AQNKOHX/UrGnAdiwEQRAExUjXqGQWFhbbttFwqi9BEAShBEiLkCAIgqjQSCIkCIIgKjSSCAmCIIgKjSRCgiAIokIjiZAgCIKo0EgiJAiCICo0kggJgiCICo0kQoIgCKJCI4mQIAiCqNCk2lnm4cOHV65cqVWrVocOHbhckjsJgiAI9SE5qx05csTNzS0pKWn27Nl9+vRhICaCIAiCYIzkFuHs2bPXrl3r5+eXm5tbr169a9euubu7MxAZQRAEQTBAQoswPT39wYMHPXv2BGBgYNCxY8cTJ04wEhhBEARBMEFCi/D169fGxsZ6enqiH83MzDIyMsq7uaioaNeuXbVq1frxYsuWLR0cHBQPVIUUFxcXFxezGACfz6cvhpKSkuLiYh6PnFvyk+LiYj6fr6GhwXYgNCouLhYIBGxHIQ6tbz3W39dyU/LIhUIhfREKBAIOhyPxNgkfZ3w+/8dSuFxuSUmJmJtfv35dVFT040UHBwc+ny8xDnXC5/PZ/SPz/0/lClddFeGvRfSnU+bpcrT+ClT396vkkQuFQvoiFAgE0nw9lZAIa9Wq9enTp+LiYk1NTQDv3r0r1eD7ka6u7syZMyta++9XxcXFOjo6LAYg+vpCUwwlJSUaGhrs/gGVkIaGBp/PV/u/Fg0NDdFHgXKi9a3H+vtabkoeuVAoFAgENEUoEAikSbESvtxZWVnVqlXr3LlzAEpKSs6fP9+6dWtqAiQIgiAIJSChRaihoTFr1qzAwMAJEyZcvny5evXqHTp0YCYygiAIgmCA5O7+oKCgnTt3ZmVlde7c+fz588o8QkAQBEEQspJq7l/btm3btm1LdygEQRAEwTwqm3f5+fkvXrygsECCUBXz588PDg5mOwqCIORB5WqwjIxBPj6jCgrSlXmGEkFQKzcXPj6DLl5MBqpevWr74EEq2xERBCEbKhNhUVFLICIvL48kQkJt5OXl6enpfVtNm5+P+/dx797X/6Wm4t07AF8AA6DKo0cD9+7FoEGsRkwQhIyo7Bo1MloBHJoxowqFZRIEi3btiqhTp52ZmdfUqV969EC9eqhaFYGBOH0aRkYIDMTZs8jJQWBgJIdTzOHw/fyGDx0KU1NERbEdOkEQUqOyRVir1stKlay2bYOvL7p0obBggmDB27dYtCjhw4dAbe2dXO6HoUPN7OxgZYVS+1ScOYO9e+HufsPQUNitG3fdOowYAV9f1K2L3btBNqgnCOVH8VqIiRPfGBmhTx98/kxtwQTBnPx8LF0KOzu0aDF10KD769ePWrrUrFcv1K9fOgu+f4+AANjaYsgQQc+e/MhIGBkhMhLPn8PcHF5esLPD3bss/TEIgpAOxYnQ2zvbzAyVKsHLi9qCCYIJAgEOHkSTJrhxA/Hx2L275p49oQEBZQ/6CYUICICfH27fRs+egm7dBGfOIC8PAGrXxqVLuH8fBgZwcICbG8h8aoJQWhQnQg4HU6agYUM8eICJE6ktmyDode4cmjXDihXYtQsHDsDaWsL9q1bh3TvUr482bWBiAhMTobMzzpz5foONDeLiEBeHrCxYWsLV9ba1datJk2bT+qcgCEJW1G8T4+eH588REoI1a3D+POXFEwT17t9Hv34IDsasWYiLg6en5Efu3sXSpYiIQETE92mivXsjMrL0nS4uSE1FVBSSkuY+fTpu1aqDmZmZFP8BCIJQAPWJUFMTEyaIOovg44PsbMprIAjKZGQgKAitWqFZM9y9i759IcXhZcjLQ79+WLkSOjq4efP71DBfX5w4gS9fynika1d062bM4SwUCmslJlal8s9AEIRiaNk4dNQoXLqExYthYoI2beiogSAUEhl5zM2t15AhJ+ztoaODBw8wfTq0taV9fMIEuLrCzw/798PXF7q6X6/XqAE7O1y4UPZTR4/ueP365Jgxl3r25F65QsGfgiAIStByzri+PkaNwt9/49o1WFtj1iwsWkRHPQQhp3HjQt6+PXznTq/797taWMj2bGQkLl3CrVsAsHcvli376dVevRAZWe7yoZo1a65fjw8f0KYNbt5EkyZyRU8QBKXoOkri999x4AB4PGzahKVLERdHUz0EIbNnz5CX16NyZd8xY3rKmgVfvcLYsdi7F5UqIS0Nb96gVaufbujTB1FRKCkRV0hEBFxd4eyM589ljZ0gCOrRlQirVMGAAVi3DsOHo1MndOyIggKaqiIIGWRmolMnzJs39/PnG8uXz5HpWYEAQ4ZgyhS4uADAvn0YOLD0ysLatWFpidhYCUXFxqJBAzg6khW3BME+Gg8XnDIFGzciKwtRUdDVBTnHiWDd58/o1An+/nKu7Zk/H5qamDTp64/h4WVvK1rm3NFfJSXB0BCNGpU9uYYgCMbQmAgtLdGhA7ZuhYYGLl9GUhLmz6evNoKQID8f3brB2xuzZsnz+JUr2LQJO3dCdDR1fDw0NNC0aRl39u2LQ4cgEEgokMfD/fvg89G4seSbCYKgD73Hzc+YgVWrUFQEGxusXYv5879OMSAIhhUVoVcv2Njg77/lefzzZwwZgs2bUbPm1yt798LPr+ybraxQvbpU4+J6ekhNxX//lZ1QCYJgBr2J0MEBtrbYuxcARo9Gy5Zo0wZFRbTWSRCl8fnw84OBAbZskWqZ4K/GjEGPHuje/euPJSU4eBADBpR7v5S9owCqVkVKCh4+JPvUEwRr6E2EAKZPx9KlX3t+Tp+GhgY6daK7ToL4TiCAvz+KirB/f+mJLVLauhUpKVi8+PuVs2dhaYn69ct9pHdv/PsvhEKpyjc3x9WrOHcO/v7yhEcQhIJoT4Rt2sDYGMeOAQCPh5gYxMZi+XK6qyUIABAKMXYs3rzBgQPQ1JT58cTEJGvr1r/9FhweLvi2ah7Avn0STt9t3Bj6+rhxQ9qKnJxw+jT27fs+E4cgCMbQnggBTJnyfUF948YIDcW0aRF167aLI6sLCZrNnImbN3H0KHR05Hl81apdT5/+qan50sDg5beL+fk4fhx9+0p4VrSyXnqtW+PAAaxeTXafIAimMZEIfX3x+TMuX/7644ABb4TChc+fj/L1HctA7USFtWgRTpxAdDQqVZKzhMGDh3C5Ie3bm5mbm3+7ePQo3NxQo4aEZ6UfJvymVy9s2oTZs8dpaNQfM2aa7PESBCEPJhIhl4spU7B06dcfq1Wrpq2dC6yysnJjoHaiYlq/HmFhOHMGVarIX8iLF84DBlyKjNzI5X5/p4iZL/qjpk3B5yMlRbYa/f2LOJzTAsHxHTtOyRgsQRByYiIRAvD3x61bSE4GAB6P9+VLeosWUQUF/zBTO1HR7NmDJUtw9uz31Q7yiYxE794/Xfn4EVevomdPqR7v2VPmRqGWllarVk483oCiohV37sj2LEEQ8mEoEWprY/x4rFz5/cr69VVv38bbt8zUT1QU27bt6dhx+pQpb0/btKAQAAAgAElEQVSdgqWlQkV9+IDERHTs+NPFiAh07gwDA6lKkKN3FMCFCweKi295e7d3c8O7dzI/ThCErBhKhADGjMHJk993GW7eHDVrYsoUxuon1N/Lly8nT9565kw9V9dQW1tFS4uKQocO0Nf/6eLevRLmi/7IzQ2fPuHBA3lqv3gRVlawsyMbsBEE7ZhLhIaGCAj4aV+PSZPw779kcymCMpUrV/3ypdjQMKxPH1fFS/u1X/T5c6SloX17aUvgcNCjBw4fljOAGzegoQEHBzkfJwhCSswlQgCTJmHPHvz339cfJ08Gh4N/yEAhQZF163TbtLny+vXZwYMlLW6QJCcHV66U3u1l71706wctLRnKka93VERLC3fvIiMDrVvLWQJBENJgNBHWqAFfX6xf//2Kjw9CQ5kMgVBbjx9j5Ur88w/HQMoRPLGiotCyJQwNf7q4f78M/aIi3t548QLp6XKGUbUqbtzAtWsICJCzBIIgJGI0EQKYOhXr1iEv7+uPq1fj1Stcv85wFIQaGjMGc+YoOkHmm1/7RW/fRk4O3N1lK0dDAz4+OHJE/khsbHDqFHbulPPQDIIgJGI6EdrYwMsL27d//dHUFE2akCkzhKJ27cLnz/jtN2pKy8/HhQvo1u2ni6JpMnLs2a1I76hI69YIC8PSpdi6VaFyCIIoE9OJEMCMGVi2DMXFX39cvBhXrnxvIxKErD58wIwZ2LRJzj21f3XiBFxdYWLy/YpAgPBwDBwoT2lt2yItDRkZCoU0dCjmzEFQEM6eVagcgiB+xUIibNECQuGfVlbtrl6NA9C1KwwNMXs284EQamLSJAwaROWRfr/2i8bGwsQETZrIU5qmJrp0wdGjikY1bx78/NC1K1JTFS2KIIgfsZAIs7KysrPPv3q1bP78r9NmRo7Ejh3MB0Kog5gYxMTgzz8pK/DLF5w5gx49froo0/LBX8m6AXd5du9Gs2ZwdkZmJgWlEQQhwkIirFy5cocO9lzuuMDAIaIrCxYgLw8REczHQqi2wkKMHo21a6Xd6kUap0/D0RHVq3+/UlSEI0fEHcMrUadOuHkT798rHh2uXoWpKeztyQHXBEEZFhIhgIMH13t5XTM07CD6UVsbrVph3jxWYiFU2Pz5cHT8fnA8JX7tF42ORpMm+OH8CZnp6KBdu6+nciqIy8WdOyguLjAyGufs3CU3N5eCQgmiYmMnEQLw9v5+MBOAtWuRloanT9kKh1A9KSnYvBkrVlBZZnExoqPh6/vTRSmPmxBP8bmj3+jpwdf3j4ICXlKS1R9/zKWmUIKowFhLhF5eiI39/mPDhqhTB5MnsxUOoWIEAowejUWLUKsWlcWeO4eGDX8qMzsbZ8+iVy9FS+7aFVeu4NMnRcsR6dOnHZd7FLi8f/+Q7GxqyiSICou1ROjhgVu3ftpQeNYsnDiBkhK2IiJUyaZN4HIxciTFxZbqF01JSXFw6GZiMtvERKhgyQYGaNUK0dEKFvNVp06dcnLuvXkTp6vrZGZG5pEShEJYS4R6erC1RWLi9yuBgdDU/H5+L0GU5+1bhIRg40Z5lreLwefj2LGfGn+LF29+9uz3Dx9uvXjxQvHye/fG7t0vSij6rqenp2dqqvfkCZo1g6MjTpygpFSCqIhYS4QAvL1/6h0F0K8f1q5lKRpCdYwbh9GjofhBS6VcuoQ6dVCnzvcrAQG9udy5Tk5GZmZmipd/9+7SM2dG29i48/l8xUsT4XJx6RKCg+HjgwULqCqVICoWNhOhl9dP82UA/P033r/HxYssBUSogpMncfcuZs6kvuRf54taWLS0sLh26dI+Ho+nePmPHj0UCrt/+FBUWFioeGk/WrUK69cjJEShNR4EUWGx3CKMj/++1xoAIyM0bYoZM9iLiVBueXkYOxYbNkBHh+KSBQIcOVJ6vmhiIpydKatiy5YlLi45Q4as09PTo6zQ/wsKwoULOHwYjo5koJ0gZMNmIqxcGZaWuHXrp4srVuD6dXz8yFJMhHKbPRutW6NNG+pLvnoV1avDxuani0lJaNGCsiqqVas2atS0nBxPykr8WcuWSEvDixeoVQtv3tBUCUGoITYTIcoaJvT2RtWqmDaNpYAIJXb9OsLD6ZpO9Wu/KKhuEQJwcMDt21QWWIqlJTIyUK0arK2RkEBjRQShTlhOhL8OEwIYOxb79rERDaGsioqKFi9e07t32IoVwqpVqS9fKMThw+jT56eLxcW4e5fKvbwBNGmCR4/o3R1NRwf37qFTJ3h4YM8eGisiCLXBciJs2RJXrkAg+OninDng87FtG0sxEcpnz559c+c+zcg4Y2p6gY7yExOhq4tGjX66eOcOrKyo3MUUgLY26tZFWhqVZZbp0CFMmgR//0gNjYb29h1or48gVBnLibB6dVSvjpSUny5yuejYEYsWsRQToXxq1rQRCOKMjB7WrVuXjvIjI9G3b+mL1A4QfkN37+g3y5ZBU3O+QHD67t1XN26Q+TMEUS6WEyHKGiYEsGYN0tNx7x4bARHK5+RJ98DAU8+eXba2tqaj/EOHyhggVPVECGDq1P6amq0rVfJu0YLn5kZm0BBE2dhPhGUOE9atCxsbTJzIRkCEknnyBBERWLjQWF9fn47yb92CQABHx9LXKZ8pI8JkIly4cFZR0dPs7I0JCcjKgrk5hgwh5zcRRGnsJ8JWrcpoEQL46y9cuPDTZqRExTR1KiZPBh1zZETK7BfNzUV6upxH0ovHZCL8pkULpKZi+3YcPw5DQyrPMSYINcB+IqxdGzo6ePiw9PXevaGvj5AQFkIilEdcHK5fx2+/0VhFmf2i16/D0RGamtRXV7MmOBx2eimHDMGnT5g/H0uXokoVHDzIQgwEoYTYT4QoZ5gQgL8/tmxhPBpCaQiFmDIFS5ZAV5euKu7dQ05OGWOBNPWLitjbs9Ao/GbaNHz8iC5dMGAArK2RnMxaJAShJJQiEZY5TAhgyRJkZSEqivGACOUQHo6CAnr3z4yMRJ8+ZZxiQdNMGRFWekd/pKeH3buRkQFzczRtCj09fy0t67//XsdmTATBHqVIhOW1CPX14e6OP/5gPCBCCRQV4Y8/sHw5uHT+Iy1zQxnQ3CJkPRGK1KiBS5dw+PCbgoK04uJDkydfaNwYgwdjyxa8e8d2cATBIKVIhDY2KCpCmSe+jRyZeOdOlxo1WmRlZTEeF8Gm1athZ0fLtqLfPHqEzEy4u5e+npmJvDxYWdFVr5IkQpEePWpaW2vzeL18fDo1bYqUFEyZAlNTaGjA1BStWmHSJBw9iszM/6ytmzZu7ErVeYoEM1auXFm7tv1MOo5rUSNKkQgBeHqW3SjcvXsR0C8z03HXrl2MB0Ww5tMnLF9O+6YKkZHo1auMFmdCAlq0oPjU3x81aoRnz1BQQFf5snr8+HJx8ZMjR0bt3o3kZGRlobAQUVEYMACFhdizB717w9S0T0aG16NHZr17/zLFllBikyevyM4OXbKE7FopjrIkwvKGCX//PZDHWwq89fEhb78KJCQE/fuX3vOMcuX1iyYl0dgvCkBTE/XrIzWVxioUpKWFrl2xahXi4pCZiZISdO5sCsQA6devz752je34CCnk54vGlbSBMKDJ2LH49IntmJSVsiTC8oYJu3btWlx839Dw2OrVpowHRbDj6VPs24fZs+mt5dkzPHsGb+8yXkpMpHGmjIhS9Y5K48SJ8EOHlhsbnwgJae7nB39/Mo6o1I4dQ5MmSE5Gx45PnZ37tG4dxeOhQQOsXg0+n+3glI+yJEI7O7x/X+7iqvHjsWkThEJmYyJYMnUqpkxBjRr01nLwYFHPnvj15HmhENevk0RYhs6d2w0eXPvVK6SmwsoK9vZYvZocAqx0Xr1C376YPBkbNmDTJty4gcjIfvfuaUyahDNncPAgXFwQH892lEpGWRIhhwMPD1y9Wvar8+aBz8eaNczGRLAhLg5JSfSuoAcwbNik2bNbPX8++deXHj9GpUq0p2FVTIQAAgIEYWHQ1kZICK5cwalTaN687EENgnnFxVi9Gk2bonFj3L2Ljh2xaRMGDkTlysL+/bFjBxwdcfkyxo9Hr17w98f792xHrDSUJRGi/GFCAFwuBg4k51GoP9EK+sWLoadHb0Xnzl0uLl734EEZ3fG0Lpz4xtERt2+rXidH48bCmjVx7hwA1K+PkyexYAGGDEH37nj5ku3gKrYLF+DggHPnkJiIkBBoa6OkBNu2YdQoAAgMxLZt4PPB4cDfH6mpMDZG48akp/QrJUqE5Q0Tiqxdiw8fEBnJYEAE4yIiUFCAgQNpr2jIkJVmZht37fr715doXUr/TZUq0NNTyeQxYsRPZ4V2747UVDRrBkdHhISQHb2ZxufzMzLg749Ro7B8OY4dw7eTyg4fRv36X/fLtbND9eq48P/TPI2MsHo1zp5FZCScnREXx0rsSkSJEmHTpkhPx8ePZb9qYIB27TBjBrMxEQwqKsKcOQgNpXcFvcizZ17z529u2dLz15eYaRFCZXtHBw7E2bM/9arp6SEkBAkJSEyEvT3++efGwoUr3r59y16MFUJ2dnbDht7Gxq62tjH16iElBV26/HTDhg0YM+b7jyNGICzspxscHBAbi5AQDBiAIUP4y5dv37lzn1DluimooESJkMeDiwvEzMzevBlPniApicGYCAatWYMmTdC2Le0V8fk4dw4dO5bxUnEx7t6FkxPtMUBlE6GhIXx8sO+XZWn16iE6Gn/9JZgwYdjcuUYDB45nI7oK5OnTpy9eVM/JCezfP2buXOjo/PRqWhrS0uDr+/2Knx9OncJ//5Uup3t33L2LDx8OTJsWFxx89Jyo47uCUaJECLHDhAAsLNC0KcaOZTAggimfPiE0lKFh4MREmJnBzKyMl+7ehaUlKlViIgwVTYQARozA1q1lv9S7N9fa2pDHi87Lo21jHgIAYGLiwOW26N49ed68Ub++umEDRo786fiUypXRrRv27y+jKENDzJ5tUbny3S9fHpeUlPXGUHfKlQjFDxMC2LAB16+XvRkbodLmz0ffvrC1ZaKukyfRuXPZLzHWLwpVToTe3igpKbdv5t69mJiYJRkZS86cYTasCmbqVM6MGdOjotabmpZeY11QgH37MGJE6UdGjCj3PB8PD4+0tCPLlp2dPt02P5+GcJWbciVCFxekpCAnp9wbWrSApSWCgxmMiaDf06fYs4e53dXFJEJmZsqI2NggIwO5uQxVR62hQ3+aMvMjHo/n5lZ/zx4MG4bXr5kNq8KIi0NcHCZNKvvVvXvh6Yk6dUpfb9kSBQW4caPsp2rUqDF5sknTpl8nmlYoypUItbXRtKmExZ4rVuDUKWRnMxUTQb9p0zB5Mu1L90Tev8fjx3BzK/tVJluEGhpo1AgpKQxVR61hw/DvvxDTdGjVCmPHok8fMo+UekIhfv8dS5eWu8po48afpsl8w+Fg+PByv8GIbNiAlBQJ96gf5UqEALy9JazP7dkTxsblfhUiVMujR486dZpw6VLUeKamVpw8iXbtyj56PjcX6elfp5szQ3V7R01N4e6Of/8Vd8+sWahRA7NmMRVThbF7N7jccs/pTEjA589o167sV4cORUSEuG8wuro4cACzZuHmTQpCVRVKlwi9vCQMEwKYPh179pB1oOpg5MhZp093KCz8Q1OzmJkaxfSLXr8OR8eycyRNVDcR4pcFhb/icLB9Ow4fJst/qVRQ8PWczvJORxGtmihvDZKZGdzccOiQuCpsbLB2Lfr3x+fPikarKpQuEbq748YNfPki7p4pU8DjYeFCpmIiaFOjhhePt9jWtrYmI/lHzMIJMNsvKqLSibBrVzx6hLQ0cfcYGyM8HKNHS7iNkN6SJfD2hodH2a9++oSoKPj7iyshIEByz2e/fujcGf7+qrf5kXyULhEaGMDWVvJiwYAArFrFSEAEbQoKkJg4Pjr6dELCCWZqTEhA7dplL5wAszNlROztcfcuBAJGK6UKj4fBg7Fzp4TbWrTAvHno109cdxwhpVevsH49/vqr3BvCwuDjg2rVxBXi44MHD/D4sYS6VqzAx49YsUKeOFWO0iVCSLGIAkBoKHJzK9yIrppZuhTu7mjfXp+xGsX0i4KNRGhkBGNjpKczWimFRo7Ejh0oltSrHRwMJycwNgysxqZPx9ixsLAo+1WhEJs3lz1N5kc8Hvz8sH27hNs0NbF/P1asqBCbqitjIhS/rF5EWxs9eiAkhIl4CDq8fo1//sHixYxWKiYRZmYiJwfW1ozGAxXvHbWxgbU1Tp2SfOf69YiLk/zhS4iRkIDYWEydWu4NZ89CTw8uLpKLCgzEzp2Sp1mYm2PPHgwYUO4BeWpDSRNhXJzkc842bEBGBsiiXRU1fTqCg8tY6kSf9+/x5AlcXct+NSEBzs7lTkCgj0onQkgxZUZEXx+HD2PGDNy6RX9M6kgoxOTJWLwY+uV3oGzYIO3GWw0awMJCqm8wbdtixAgMGqTmkxOVMREaG6NuXclvmGrV4O6OyWWcKEcou/h4xMRg2jRGKxWzcAJAUhLTM2VEVD0R9uuHy5elajHY2GDNGvTrh6ws+sNSO/v2oaAAfn7l3vDqFS5fluHkFim/wQAICYGmppp3vyljIoR0w4QAtmzBvXtITaU/III6okMHFy0S992WDuIHCBMTmR4gFFH1RKivj169sHu3VDf3748OHSrQXESqFBRg9mysXi3uYJbNmzFokAzvqf79cekSpDkjhMvF7t3YsQNHj0pbuMpR0kQozTAhgIYN0agR2YZbxezdi4ICDBrEaKV8Ps6eRadOZb8qFOL6dXYSobU1PnxQ7QVboraFlLlt1Sq8f0+mfMsmNBTu7vAs49Cwr0pKEBYm29ZoBgbw9cWePVLdXL06Dh5EUJAKT+wST0kToWh/GWmmla9Zg9jYn05HI5SZNN9t6ZCQAHNz1KpV9quPH6NSJYb2eCuFy0Xjxqq60ZqIqyt0dHD1qlQ3a2oiPBzLllWIuYiUeP0aa9dKmFZ2+DBsbNC4sWwli04RkfIbjKsrpk5F//4oLJStFpUg1afR69evY2JiUlJSGDuz0dQUVavi3j3Jd7ZtC1NT/PYb/TERVFi6FJ6e4r7b0kTiwglWBghFVL13FGL34P6VhQW2bcOAAVL1yxGzZmH0aAnTykqdwSsld3doaMhwPP2kSbCwUM9pGZITYUBAgJOTU0hIiI+Pj5ub22emOnGkHCYEMH8+IiPV83uKmnn1Cv/8w9Chg6Uo2wrCH6lHIoyKkmEr/C5d4OFxsm7dFi4u3YolrkOswG7dwtmz4pZM4P9n8PbsKU/5w4eXPrZeDNGeebt2jTIyar5ixSZ56lNWkhNhYGDg69evL168+PDhQw6Hs3r1agbCgtTDhABGjIC+PmbOpDkgQmGi5cBMLpkQyczE06flLpwAezNlRNQgEVapgtatceCADI9oap4sLJydklL0Ru0XqSlgwgQsWgRDQ3H3bNiAwEA598j198ehQ+JOvivFwECooRGflbVt2zax25WqGsmJ0M3NTbQPJI/Hs7W1/fjxI/1RAYC3N2JipL35t9+waROZiqbU4uMlLAemj/iFE8XFuHMHTZsyG9MP7O1x757Kr9OSfjq+yLx541u1iuTxOqSllbNRSoUXEYHcXAkbh+bnY98+BATIWUX16mjVCgcPSns/h8MJDZ1uY7OosHC+Ou2Zx5F+2O/ly5dOTk5HjhzxLGeEx9bWduTIkZaWlj9etLOzs7Kyki84S0vumTOC+vUl38nnw9CQu2iRcMIE9pNhTk5OpUqVWAyAz+cXFhbqlXdYmWJKSkqKi4t1dXVlekoohLs7d9w44aBBLPyCBg7kdOyIYcPKrvrmTQQEcJOTFdrxs7i4mM/n6+joyPe4jQ332DFBgwaKhEC7L1++aGholLc9ukAAa2vu0aMCe3sZyoyJwaBB3MREQXnzmGRC61uP4fd1QsKtAQOsdu6s7O0t7ratWzknTyIyUtzbSnzk0dGcxYs5ly/L9u9/+HCOhga2bqXg7SwUCvPz8/XpWU0lEAiEQqHEPf15ov/7/PnzkCFDfn15xowZHh4eALKysnr16jV69OjysiCAoqKiI0eOGP7cjB82bJipqanM4QMA3N11zp/nm5lJNYTQu7fO4sUagYHsf0spKCjQ0NBgMQBRIqSpcFEilHXa1P79vOJirR498pn/Fsnn49w5/YUL8/Pzy4756lVNJydufr5Cf2OiRCiQd//sxo11kpJKzM0lbafEKvGJEMCAAVphYZxFi2T4m2zRAiNHag0cqHHsWAGPp2iEtL71mHxfL1q0etWqO0JhWqNGZ/Lzy80Q8fHxixZlrlzZNT9fXGDiI/fywpgx+jdvfmnYUIZ/vaGhnJYtdbdsKR40SNEhXqFQWFBQwKFnVyeBQKCtrS3xtq//9PT09IKDg3992draGkBubm7Xrl1dXFwWLFggpiw9Pb21a9c6ODjIFXAZWrfGtWu84GDJfwwAGzfCyAinTxv07k1V/XISCoUGBgYsBsDn8zU1NZWnRVhQgIULsXcvDA1Z+Gu5dg0WFqhfv9xPkzt34O4OAwOFzoFSsEXYrBkePOCx+q9GMh6PJz4Rjh4NFxeEhmpK8cnzXUgIkpIQGmog5lAFKdH61mPyff3pU15BQT1j48c8Hq+8Sh89etSnT0h2dtO0tMxevX4XU5rEyP39ERGhFxoqQ4QGBjhyBF5e2q6u2gp+5AuFQi6XS1+LkC/FqMPXRKilpdW5nEl1+fn5Pj4+DRo0WLNmDU1Juzze3liyRNqbDQxgaBjcp89lT0+by5fJSaBKZPFieHmxsGRC5ORJdOki7obERIwbx1Q05XBwwJYtLMegOEtL2Nnh+HHI9GWUy8XevWjWDC4u8PGhLTiV4uExJybm3z17elepUqW8e/T19UtK8rS1n1Wr1kzB6kaOhLs7Fi2SbcZNgwZYvRr9+uH6dbA6FkQByZNlhg0blpKSUq9evdDQ0KVLlx4+fJiBsEQaNEB+Pl68kPb+7OzjwOlr127QGRQhG4knqNEtOlrcwoncXKSno0kTBgMqixpMHBWRdcqMiIkJ9u/HqFF4/pyGmFTQrl16c+f6NxU7g6tWrVoODqdWrforMHCogtVZWqJhQxw/LvODAwfCy0u2HW2Uk+RE6OXlFRAQkJWV9enTp0+fPuXm5jIQlgiHI8MiCgCBgQM0NHpqai4k00eVx7RpGDeOhSUTIpmZSE8XdzDNjRtwcJBz6jmF6tRBXh7++4/lMBTXqxcSE2X48vqNqysmTUL//igqoiEslZKejuRkyesCCwtx924NPz87SiqV5tj6Mq1bhwcPVL8/Q0gdOzu75ORkCgsUCoWrVwuDgmS4v7BQqKcnnDiR2ihkk52dzWb1QmFJSUleXh5NhRcXF+fn50t587VrQnNzYW4uTbFItmOHsG9fcTcsWyb8/XcKKioqKiooKFCkBC8v4YULFERCn4KCgqKiIom3NW26pFq1ltHRZ2QtXyAQ9uwpnDJFruCEQiHNbz3G3tezZwsnTZJ8W2ys0MVFqgKlifzZs7c8Xk8Pj/5ZWVlSFfqDhw+F1asLb96U9bmvBAJBLm2fEXw+X5p/tEq61+g30u8vI6KlhVWrsGaN+p8kqfyEQvz+u4QT1OgmfkMZsL2nzI/Uo3dUIBA8ebL//fvtc+fKvPMGh4OwMERGgsHhF6VTUoIdOzBihOQ7Y2IgfmWFTM6ePc7neyclWV+8eFHWZ+vXx9q1qn3AlrInwgYNCl+82LBliwyTXwIDUbcufH3pC4qQ7NmzZ6NH7xYIPos5QY1ufD7OnUPHjuLuSUxkc5fRH6lHIuRyuf36deJyB44cOVKOx42NER6OMWPU9pQDiU6cgKUlbG0l3xkbi5YtKau3Y8cO5ubRenrJYhbIidGvH9q3R2AgZfEwTNkT4caNWwsKHk+cuD0pKUn6p44dQ1KSOp+epfw8PHy3bHkrFI5i/sz3b+LjYWFR7okTADIzkZMDa2sGYyqfeiRCAJs3L+nSJb5KFbn2vgScnTFjhtqeciDRli1SpZOSEiQmwsODsnrNzc0TEs5yOCeMjcudpyre33/j6VP88w9lITFJ2ROhjY2Vrm4i8LqGLGfkNGoEX18MG0Y2XWNNbq6WltZzc3M2Z1VL7BdNSICzM1hM1T9q0gQPHkA9NqD29JT2VKYyTZgACwtMm0ZdQCri1SvEx6NPH8l3Xr8Oa2sYGVFZu6kpTExw/76cj2trIzISf/2FGyo4bV/ZE2HXrp137jzQpEmMhYVsGxLu24eiIkyYQFNchDhRUTAxOXfkSL+DB9ncol6aAUIl6RcFoKsLCws8eMB2HFTw9MSVK/I/zuFg61acOIF//6UuJlWwbRsGDoQ0O2FQO0D4jZeXQr+4OnWwaRN69wZTO1JTRtkTIYBOnWqmpBjK2k+ipYU1a7B+PV6/picsohzv32PMGOzcWalTJ2+e4rtmyevtWzx7Ju7ECbB96MSv1KZ3tHlzpKXJcKbBr4yMEBmJ4GA1+WYgDYEA27dLNU0GVA8QfuPhoVBTHkD37ujZE0OHqlhvnAokQn191K8vzwfEiBGwspLzmC5CbmPGYPhwWr6uyuTUKbRrBzGJWCjE9eskEdJCWxtOTkhIUKgQBweEhKBfPxQUUBSWcjt9GjVqwNFR8p18PuLiaNmqScGmvEhoKD5+xKpVVATEFBVIhABcXREfL8+Dx47h5s0KPRubYVu34ulTzJ3LdhxS9Is+eQIDA8gy9Ew7tUmEoOgjNTgY9vYYO/bTfblHrlSHlNNkACQno3ZtVK1KfQwNGiAvD69eKVSIpibCw7FsGfbuTWXsIHcFqXkibNAAffpg+HAVa6erqPR0zJqFnTuhpcVyJHw+zp9Hhw7i7lGehRPfODggOZntICiieCebyMKFH3bvbuXsPHXduq0UFKes3r3DpUvo31+qm2kaIATA4cDdnYJfnLk5unbd4POZfu4AACAASURBVO8/o0GDliqRC9U8EQLYswclJRg7ltKAiF8IBBg+HLNnw46aLZ8UInHhBJRpKf03ZmYQCPDuHdtxUMHDAwkJKFH4XCkOJ8/AQD83t8mFC2rx91KOsDD07Svt1tU0DRCKUPUNpnLldxyO3YcPOs+fs380nkSqkQhtbJCTg7dv5XlWUxNr12LzZnn2PySkt2QJeDz89hvbcQCQol80Ovrs9u1+QuE5piKSlr29mvSOGhmhbl0KGrgWFhbHji0LCWmYmDhp504qIlM+QiHCwqTtFxUKceUKjWe5UNKnDWDBgmn//GMdFLSqf/9aCva1MkA1EiGHgxYt5B97Hz4c9eqRWTM0unULq1YhLAxc5fgHJTERBgXNzMoKWbt2BlMRSYsME5ZVjueffw47f173jz9Udb22eBcuQE8PzZtLdfPdu6hWDTVr0hVMs2Z4/JiCzdL09PSCggL++cdt5Eh4eeHpUyqCo41yfG5JQZHeUQBRUbh9G/v3UxcQ8X+FhRg6FKtWQcalnnSRZuGEi0tHHm9wjx5isyUb1KZFCOo62UQaNEBsLP7+W8WmI0pjyxYEBUl7M30DhCKammjWTKEP21KmTMGMGWjTBo8eUVYm5SpKIrSxQf/+CAqCFIcVE7KZORP164PFPUVLkbhwAkCrVn8NHx6/bt0CpoKSlpq1CKU/Q00adeviwgX88w8WLqSyWHZ9+IAzZzBwoLT30zpAKEJVU/6boCAsXIjWrZGSQmWxFFKlRHjjhkJj77t2QSBAcDB1MRHA5cv491/lOo1MYr8ogIQEuLoq4z9+W1s8eYIvX9iOgwp16kBLi+J2gIUFLl9GRARmKF2vtpx27ECPHjA2lupmoRCXL8PLi96QPDwoToQABg/G8uXo2FFJv+cp42dBmQwNYWGh0BcKHg9btnxd6EZQIisLQ4Zg0yaYmLAdyv9t2rT76NFpjo4SJlbFx0voO2WLtjasreXf71HZUN62AGBqivPncfIkpk+nuGRWSD9NBkBaGgwMYG5OZ0CAuzuuX6f+hOQBA7BxIzp1UnSnBTqoTCIE4OqKuDiFShg4EA0bkhOaKDN2LHx8JDe/GPPy5cvp07cVFtqsXRsq5rYPH/Dff2jYkLG4ZKNOvaPUDhN+U706Ll3CpUsIDlbtJcKXL0MggJubtPfTPUAoYmiIevVw8yb1JXfvjm3b0KMHrl2jvnBFqFgiVPyrxIkTSEnBnj1UBFSxHTqEGzewZAnbcfygatWqGhol+vrb2rYV99ESFwdnZ2WZ4Pqrxo2/nDt3W6jSH/D/R0eLUMTYGGfO4PZtBAVBIKClCgaIpslIf/4JAwOEIvT94rp0QXg4evXCOWVau6SsHwZlUbxFCKBuXQwahDFjyKwZhWRkIDgYO3dKtVM+Y3R1datVu3zmzDk/P3En2SQkKGm/qMiaNV4HDoQGBk5lOxAK2Nnh3TtkZtJSeOXKOHMGT59iyBAKVu4z7/NnHDuGQYNkeISxREhTU16kVSscPAg/P0RF0VWFrFQpEdra4t07/PefouVs3w4OB6NGURFThSQUIjAQY8cq3RZlaWnIzeW4uemLvy0+Hi4uzEQkMz6fX1JSWFzs/OqVXPtHKBkuF25uNH6k6usjKgr//YdBg1TvKMfdu9GlC6pVk/b+x4/B5aJuXRpD+kZ0HhN9vRJeXoiORlCQsmwErUqJkMtFixZITFS0HA0NrF9fFBY2zczM85Eyr21RSlevXg0MPJuZqYzT9o4eha+vhI4mgUC5jiEsRUND4+LFfZUqaS1fvprtWKhBa9sCgJ4ejh1DURF69sxfvXr9a9U5d02maTIAYmLQqhVdwZRiZgYDA3rPwGreHNHRCA7GokW3t2/fLmC1g1uVEiEUXk34TXb2VuBNRsaAsWOV7+NciZ09e7Z168B9+7bUrz9DU5PtaH5x9Ch69JBwT2oqatSgZed+qjRp0sTVdfTLl1XYDoQa9I02faOlhYgInD/v9scfTy0tlbjX+wcJCcjJka2fMzaW0dPNGPjFOTlh69ZHc+b4jB9/vmNH6XYcp0cFTYTt27fX1LwMbE1OHqdyPSosys3NEwr1AAMOR+EeaqplZiItTfInhdIunPiRvT3u3GE7CIo4OyMlBfk0b7yspQUdnSLAks9XplHr8m3ZglGjZJgmAyAmhqEBQhG6m/IideoUcjhaQqH58+ds/uJULxEmJlIwSax+/fqFhelv3ybx+a0dHFR41hnDMjJ6Wln9Nn261fbt69mOpbSoKHTqJPkEqIQE5R0g/MbODnfvsh0ERXR1YW/PxNKxu3fPDByYYWR0NjWV9roUlJuLQ4cwZIgMj7x8iaIi1K9PW0y/YKBFCKBJkyYHDy4dPNggN3dzeDjt1ZVHxRJhlSqoXp2a5cYcDqdGDc3UVLx+DXt7kgsly8jA/Pk4cmTYn3/O0WL9yMFfSNMvCtIiZAMzH6nm5uabNv01b57FyJHK/nbeuxdt28q2cfbFi4w2BwE0boxPn/DmDe0V9erVa+PG2adOaU+ciKNHaa+uTCqWCAG4uVG5IWyNGnjwAC9fgrQLJQoOxm+/oVEjZVzflp+Py5fRsaOE23Jy8OyZUpyYKJ6tLR4/RmEh23FQxNOTiU42keBgaGhg0yaGqpPDiBFTx493qVs3QqanGB4gBMDh0DvjtxR7e5w8idGjcfw4QzX+SPUSoYsLlYkQgKkp0tLw/DkcHVV7lwpa7d+Pp08xbRrbcZTj9Gm4uMDISMJtCQlo2hRKOM2nFG1tWFkhLY3tOCji6Yn4eIZW7nK52LoVISF4+ZKJ6mQlFAoPHTpTVBR+/vwOmR5keIBQhMlvMAAcHXHiBEaORHQ0c5WKqF4ipGq+zI9q1sT9+0hPl2GvowrlwwdMnoxNmySPwLFFnfpFRdTpPCYTE5iZMdfZ26ABxo7F6NEMVScTDofTqNHw6tUDFy+eIv1Tb94gKwuNGtEXV9mY6dP+UdOmOHwYw4fj4kVG61W9RGhvjxcvKDg3shQzMyQn484duLtTXLIamDQJAwYo77cEPh/R0ejeXfKdKjFTRkSd5suA8Y/UWbOQkYEDB5irUUolJXj27PfY2HMdO7aV/qlLl+DtLdsUU0o0b460NOTkMFqpmxsOHcKAAYiNZa5S1UuEPB6cnChYVv8ra2vcuYPkZHh6Ul+46rpwAbGxmD+f7TjKd/UqLCyk2pI/MVFlEqGDg1rNl2FmLv43PB42bcL48Xj/nrlKpXH6NKys0KCBbE8xs9f2r7S14eTEwmERHh7Yuxf9+jFXteolQtDTOypSrx6Sk3HjBjv/7JRQfj5GjcLGjTAwYDuU8knZL/rkCXR0YGZGf0BUUL+Jo9Qe0iuRszMGDFC6Ue1du+DvL/NTjG0x+ivme0dF2rXDrl3o0QPXrzNRHUmEpdnYIDkZSUlo356uKlTI7Nnw9JQ8G5Ndx46p2wAhgNq1UVxM13bVzLOyApeL9HRGK120CLGxOHOG0UrFyMrC2bPo10+2p96/x7t3rE11Zrgp/6MOHbB9O9q3v+js7Hf48Ala61LJROjujoQEGmd4NmiA+HjExqJbN7qqUAlJSdi/H6HijvZj3717KCqCvb3kO1VogFDEzk6tGoXu7ky3LfT0vp5zlJvLaL3liYhA+/aS5zaXEhMDT0/WTg3z8EBCAmuHe3TuDC2tGUlJf44YMYfWilQyEVavjkqVQOt22Q4OiI/H2bPw8aGxFmVWUoKgIKxeLcPu+Kw4ehQ9e0p1p2q1CKF2vaOstC3atIGXF/78k+l6yyRfvyhbA4QiRkaoWxfJyawFMHhwNxOTIQJB2z59aFzdr5KJEDT3joo4OSE+HqdOSfs5q2YWL0bNmujP5ka4UpFygPDLF6SmwsmJ/oCoQyaOUuLvv7FvHwVHmSro8WM8fowOHWR+kMUBQhG2fnEiK1b88f59/Pv3y1u0gKMjNm+mpS+QJEJxnJwQG4tjx05oaDjY2ramvT6l8eAB1q3Dxo1sxyFJRgYePZJqlu+NG2jcGLq69MdEHTVrETo44PVrfPjAdL1VqmDlSgQFsXxg4e7d8POTeTOHjx+/7vXBIg8PNhMhAC6Xq6mJ6dNx7hy2bUPLltSfD0USoeSKDAxCBIJd9+/n9++f+ewZE5Wy6969+8OHZ//5p1QLEth17Bi6dpXqw0Xl+kUB2NkhLU31zpstj4YGnJ3ZmXkxcCAsLTFzZvq7d+9YqB4QCrFnjzz9opcvw90dPB4NMUmN3Rbhj+zsEBeHwYPh5YWQECrfGqqaCJ2c8PAhQ2Pg06cP4PF6GhtXiY+vbmWFKlUwfDgyMpiomnnLlq1zcZlx44aXn18227FIdvSotIO4KjdTBoCuLmrXxsOHbMdBHRanIPbpc3rlymG2tl3S2Ni5LjYWenryNOzYHSAUqVMHWlr0zsmQHpeLUaOQkIC4OLRogevXUVxcLFS4t1RVE6GWFuztceMGE3XNmjW5uDj948fo58+RkQF/f5w7h9q1YW6OiRPx+TMTMTDm1at3eXkOenraX77QfIKcwnJzcfWqtEs74uNVLxFC7XpHWWxb8PmZPJ5NXl7Njx8/Ml/7rl0YNkyeB1kfIBRRnkahiKUlTp3ChAno0CHG2Ni9Xj2PT58+KVKgqiZCAG5uLAyAm5ri77/x8iVSU9GyJXbtgokJ6tTBggVqclZA9+4zatasEx292tTUlO1YJDh1Cu7uMDSUfGdGBr58gbU1/TFRzd5erebLuLri9m3aD+kt05AhfkuXttbSGtO8OdObKBYU4OhR+PnJ/GB2Nh4+RLNmNMQkIxab8uXhcDB8OIYOvZafPzg9vZaT05OhQ7FhA27dkmexhwonQhcXFvb++aZhQ+zZgw8fcOsWXF2xdCn09NCoEfT0HCpXbhwQMIq1yBQTGak/YcIIN6XdV/QHUs4XBRAXp3oDhCJq1iLU00PjxgztFVKKhobGxIl+jo5dT59muurDh+HiItvpgyJXrsDFRSl2ule2FuE3c+aM8vd/Pm+ex/nzzdq2RVoaAgNhYgJPT0yYgB07+I0addu795LEclQ4Ebq54do1toMAHBwQEYHcXJw4gcqVnxUU5AuFx3fsiGE7LnkUFyMyUgWWTADg83HqlLQ7HqjiAKGImiVCsP2R6ueHffuYrlS+5YNQjgFCETs7vHunjPscValSZceOlX/8MdHamuPvj9Wrcf06Xr1CSAiqVcP69ckPH+qmp7eSWI4KJ0Jzc2hqQnmmcXbqhPj4ulxuETAAWMDMpFZqnTqFRo1Qty7bcUghNhaWlqhdW6qbVXHKqEjdusjOZmHJAX3Y7WTr1w8nTyKbwXlgGRlISpJzXw4lGSAEwOUyekivggwN0a4d5szB4cP1udwUgeCSxEdUOBECcHVlf51sKXz+8+zshJ49+3l7K11sEu3dK89IBiuk7xctLsatW2jRguaA6MHhoEkTpKSwHQd1vLxw9SpDh/T+ysQE3t44epS5GvfuRe/e8ixgzctDSgqcnWmISS5KOEwo0bJlhsOGpc2d20rinaqdCNkdJhTj0CH4+KBlSyXtWC9TXh5On0afPmzHIZ2oKGkT4Z07sLSUak6NclKz3tFq1VCjBu7dYy0AhntHd++Ws1/02jU0awYdHaoDkpfSDhOW58YNhIdj8WKpVlaodiJkZeKolP79F76+aNVKiTa/F+/QIXh5oWpVtuOQwp07X5tK0lDdflERNdt6G2x/pPr4ID4ezCysv3kTubnw8JD5wZiYy/36eb55EyAQCGiISx7OzkhJYWfGrxz4fAQFITRU2g801U6EzZrh3j0UFLAdRzkiIjB0KLp0wcmTbIciBbXsF4Uqz5QRUbMWIdjuZNPTQ/fuDB1ev2sXhg6V52T59evDP39e9O5dRobS7Nyhqwt7eyXtgfvVmjUwMMCQIdLer9qJUFcXjRrh1i224yjftm0ICED37jh+nO1QxHr/HvHx6N6d7TikI1MiVPUWob09UlOhNA0DCnh6IjaWzQCY6R0tKUFEBAYPlufZwMARGhoLevWyM1Omg6RVpXf05Uv89Rc2bpThK4hqJ0IwuOmo3DZvxqhR6NEDERFsh1K+8HD4+EBfn+04pPD6NZ49k7a76cMHZGaiUSOaY6JTpUqoVg2PH7MdB3VsbMDn48UL1gJo1w7PntG+Z1h0NOrXl3Mbh+fPm/btezYsLJQjR3OSNp6eqjFf5rffMH48GjaU4RGSCJmwfj0mT4afH/bvZzuUcqhWv2jXrtJuQxwfD2dn1g41pYqDg7r1jjJ/SO+PeDz07YvwcHprkXv5IGTs82CMpyfi41mb8Sul6Gjcv49p02R7SsU/IVQkEQJYtgzTp2PwYOzZw3Yov3jyBOnpaNeO7TikI+sAoUr3i4qo2UZrYHuYEPT3jn76hPPn5ZyDXVCA2Fh06kR1TAozMYGZmVJ/J8vLw9ix2LBB5tm2Kp8Ira1RVITXr9mOQwqLFmHmTAwdqnTn/O3bh/79WT7qRUo5OYiPR/v20t6vonttl0ImjlLO1RV8Pm7epKv88HB07gwjI3mePXMGLVrI+SzdWP/FiTdnDlq3Rps2Mj+o8okQgIuLajQKASxciEWLMHYs/vmH7VB+EB6uMv2i0dHw8kKlSlLdLBTi+nV1SITqN3HUyQnPnoGNcyC+698fe/fSVbiC/aLy7UTDANab8mLcuYP9+7F0qTzPkkTItOnTsWQJxo/HmjVshwIAuHEDhYUqky2kX0cPIDUVVauqxspI8erVQ2YmoxuD0Y3HQ4sWLL9tBw3C/v20jHg9eoRnz+QcaxAIEB0t7Sa6zPP0xOXLbAdRFoEAQUFYuhTVqsnzuDokQlUZJvxm6lSEhmLixGJX11khIQvYDWbfPvj5ybPUiXnFxTh9WobPCFVfOPENlwtbW7XaaA1KMAWxYUPUrIkYGrbH37kTgwfLOdYQF4eaNWFpSXVMFLGyApeL9HS24/jF2rXQ0pK/Fa4OidDFBcnJKCpiOw5ZTJqE5s2nJyRkzZ9/KSwsjK0wBAJERKjGcRMAYmNRv74Mx9mo+lL6H6lf76iHB/ujTX5+1PeOCoXYt0+GpdylyNTnwQp2Z/yW6c0bmRcOlqIOiVBfH9bWuH2b7Thk1KWLCYeTIBR+MDSsw1YMFy+iRg00bsxW/bKRdU652rQIoY7zZays3sbFTd+8mfFTkX4wcCAOH6Z4a6pLl2BkBHt7OR9X5gFCEROTswsXTkpLS2M7kK/evn07enR+cLBCy4XVIRFCBXtHAfz555wLF5bXrn1iyZK2bMUg6hdVFcePy5AIc3KQni7/55GyUb8W4V9/LSsurj916sbX7M35rlULTk4U74CoyDSZR4+QmwsnJyrjoVx4+KSHDzsOGTKZ7UAAIDLyaIMG/aKj3QMC3ipSjvokQlXZBO9HrVq1iokxS07G1q0s1P7lC44cUY1+0S9fvoweHZqTs75BA2m3GktMhJMTNDVpjYs5oqWEQql20lcNrVq56OuH8XiCKlWqsBgGtb2jV6/ePHTozoABcj5++DB69lT2AfsmTepxuX81berJdiAAcP/+45wcLx2dGllZ7xUpR30SodIeQyGelRXGjcO4ccjJYbrqEyfg6Cjt2bbs2rNn39atmdnZcRcuXJDyEXXqFwVgYoLKlZXoGGrF+fv3j44+U6PGZR1Wjxrq2xcXLuDzZwqKOnv2bMeOMwsKxj99ek2+EqKilL1fFMDly4f9/Y83bDiT7UAAwNd3jL5+jS1bguzs7BQpR00SYYMG+PyZoaNVKLdqFUxMWBghV6F+0SZNGgOXDA3vWUu9daM6zZQRUb/eUS8vg5wczsOHbMZgaIg2bXDoEAVFFReXFBZqa2trl5SUyPF4ZiZSUpTlSHrx+vY1ZPJwYzEOH9YbNmz8gAG9FCxHTRIhhwNnZ5XsHRWJjkZMDCIjmasxOxsXLqB3b+ZqVER2tou9/cWXL69ZSj2vPDFRrVqEUMdEyOGge3dGz4svE1W9o9Wqda5WbcqpUyHe3t5yPH78ODp2hLY2BZHQrW1b3L6N//5jOw4gPBwDB1JQjpokQqjmfJlvHB3h54ehQ1FYyFCNBw+iTRsl3cbpV1u2ICjIQPo+tCdPoKUFZTrBhgJ2duq24yiAHj3YT4RduyI5mYJtGpcvx4wZ3h4ebvI9rvwLJ77R1kbbtjhxguUwbt9GTg7c5Pz7/glJhMpixw5oaVHz7UYaKtQv+t9/uHBBtr8ZNRsgFFG/FiGANm1w/z7eKjTjT1E6OhSckvb6Nc6fx7Bhcj5eUIBLl9C5s0IxMEkZvsHs349Bg6iZW6Q+idDFBdevQ67OeaWgoYHISBw9iosXaa8rIwO3bqFLF9orosT27ejZE4aGMjyifgOEABo0wMuXyMtjOw5KaWqifXtER7McxqBBivaOrl6NYcNk+1f6o7Nn0awZjI0VioFJXbvi/Hnk57MWgFCIiAjKWg7qkwiNjFC7Nu7dYzsOBbRujc6d0bs37ceR79pV4OsLXV16a6HK9u0IDJTtEbVsEfJ4aNAAqalsx0E1ZWhbtG6Nd+9w86acS+vz87FjB4KD5Q9AhfpFRUxM0KwZpJ7ETb2rV6GrC8Xmin6nPokQgJubqi6i+ObQIRQXy/y5LxNf35F//NEhM1MpZj9LJNoKUqas9uULUlPRtClNEbHJ3l71dlCSqEsXxMSw3NIVCEoKC1t7e7fZtk2ehuH27WjZElZWcteOEyfQvbucj7OF3W8w+/dj8GDKSlOrROjigitXVHuLfi0t7N6NHTuQnPy/9s40rKlr6+MrkhCFgMgoCFZUFEUsLURBMYl6K1BbRKtoW6stDtWiWNHr9IgDl8ei17FUqZdirxODYh0AJ55HE6wWjAOCcoVAQJSigpUhmJHk/RDNSxWQJGdIjvv3KQfOWWslh7DO3nvt/8LLxdWrt1Sq7ffv5+PlAFN++QW+/Va/S27ehOHDzWa8qxfU69ALAL17A5sNeXlkxtDY2KhSyVtbY86d01tGU6OBn36C77833HthITg7m67QdmdEREB2NjkN61UqOHECy4oKSiXCEycWZ2ZGREYuJjsQo4iIgDFjYPJkvOwPGrTL3/9QWpppdIHqksZGyMmBL7/U7ypKzotqoWS9DJA9tgAAR0fHHTuW9OtX+MEH6/S9NjcXrK1h7FjDvesromsivPceuLiQs2ktLw8GDsTy0YFSibCs7JZKtbmwELe200Rx9iw8ewb//Cf2lm/dgpoa7u+/72Wz/bG3jjWHD8PHH+vdUJCSlTJa3n+fmolw6lTIySG50i0q6ovTp3fv3++h7xamXbsgNtYo12YhKNMhERHkPMGkp2NcYE+pRJiZ+VPfvr+tXGlK3d8NwsYG9u6FnTtBJMLY8ubNsGYNkKpppQcHDhiyXErhEaGTE1hawqNHZMeBNf36Qf/+5C/w+/uDjw/89796XHL3LpSVwfTphjutqICmJvA3g+fSDpgyBU6eJNqpTAa5uRAZiaVNSiXC0aPZX3+9q7ExgOxAMGDePPD1xbhRdVERCIUwfz6WNvGjsBBaWvRTnFIqlWz2lNpadnk5qStOeIJmR3Fl40bYskWP5qY7d0J0NFhaGu7x1CmYMsXUhbY744MPQCaDsjJCnZ45AwEB4OKCpU1KJUIA8xZae42LF0EshsREzAzGx8OqVWZTRZKSAgsX6vcPora2tqxMrlbHZWSQLXqBGxROhKdOkR0EQGAgDB0Khw516+T6ejh9GhYuNMqj+c6LAgCNBuHhRN84zOdFgXqJMDAQrl8nOwiMcHaGhASIi8NGTPzePSgoMPZLSxgSCfz2m95tvgcMGDB0KHfw4PS4uCX4xEU+lBRaA4D33we12iT2AW/aBAkJoFS+/cx9+yAyEozpIvXsGZSUAI9nuAXSIXgo39wMfD5MnYqxWaolQldXYDJBLCY7DoxYvRreew8bCZjNm2HFCrCywsAUAaSlwYQJ4Oqq94Uq1dojR9IHDx6MQ1AmASW3EmoxBQFuABgzBgYNervQjFwO+/dDTIxRvrKz4aOPzGbNvkN4PCgrI04k7/hxmDgRevfG2CzVEiEAjB5NnUEhAOTlwZ078PPPRhkpLQWBQO8NeSTyyy+GrGU2NYFYTM2t9DqGD4eqKpDJyI4DB0xkmRAA4uMhIeEtVazp6eDnB8OGGeXIrOdFtTAYEBICOTkEucNjXhQomQiptEwIAJ6esHw5fP89NBshFfCvf0FsLLBY2IWFJ8XF8PgxfPSR3hfm50NgIHW60neIpSUMHAj/+x/ZceAAlwtisUnUxI4dC+7ukJ7e1Tl79hi1iR4A5HK4fNmchLY7g7AnmMeP4dYtXESSKZgIR4+mVCIEgH//GwDi7eyGc7mGlGmLRHDpklFCiATzn//A/PlgYaH3hVeugEGd4MwMqtbLWFhAaCj5zX20xMVBfHyng8JLl0CpNORZrT15eeDnZ9QSo4kQFgb5+dDSgrujjAyYMgWXcj8KJsKAACgp0aMA2ixQqQ5oNL9fuVJigNz75s2wbBnY2OAQFg5IpZCZaWA7G4HgXUmElKyXAVOaHZ04EVxd4dixjn+7ezcsW2bsngezE9ruDFtbCAwkQiQPp3lRoGQitLaGgQOp9sj8zTfTLS1HM5lLHBwgLU2PCysq4MIFcxoOHj8Oo0dD//56XyiRQGkpsNk4xGRiUHVECAChoXD1qlGrABiyfj0kJHTQCqaiAv74Q2/lv9dQqyEnB+ONwiRCwBNMZSVUV8OECbgYxzIRqtVqNd4NhLoHifUyarW6FQch/ZSU7XK5SCpdungxfPUVBAV1pdavVCql0pcNZRISICYGy0707Y3jQUqKgc03rl0Df3/SCvAUCoWMqAqWkSNx1GTvAplMpsB5poXFgrFj4cIFAy+XSCQajQarYCZNAgcHOH78deN79sCiRcYWYF+/Dg4OQEx1cwv+s5ZTpkBuroEieSqVH7hSeAAAC69JREFU6kU3ZrrS0mDmTKDTcTHe3USYlpZGo9G2b9/exTmVlZUVFRXdNIgrJNbLFBcXT5o0CT/7O3fC7dtQXQ3OzpCV1fE56enpMTExAFBZCTk5sHQplgEcPnx4xYoVWFpsR1kZVFQYuBien6+fDA227Nu3b8OGDcT46tcPNBpsdpfqRVxcXHJyMt5ejBlbfPjhh3V1dRgGs24dxMe/HBT6+vo2NDQ0N0N6OgYF2ETOi3p5eUkkElxd9OsHnp5w9aoh1164cOGrbmwZNqwNb25ubneMdysRNjQ0JCYmBpqPgCP16mXaM3Ik1NXBt9/CzJnA40EXw7MtW2DJEiyHg3iTkgLffGNg2Wd+/juxQKjF15eys6MREXDuXLf2sxNAWBiwWH+T09y/Hz7+GNzdjbV8+rTZb5x4DVxnR4uKQCLBUUO4W4kwOjp6/fr19vb2eEWBNT4+8PgxPH9Odhx4snMnXL8OpaXg6AgnTnRwQk0NnDkDy5YRHpmhKBRw5AhERRlyrUwGRUWU1dp+EwovE7q4gJcX5JtMu8y4ONi48eWgsK0NkpMxmGKprITGRqqtZ+MqwJ2eDl9+iaMi69snXLOzsyUSSWRk5MGDB9968h9vCMg7OzvbkTEk8fIalJn5NCgI/5Lev1NeXi6VSu8QIv5Bp8PFi7Bhg8eMGfZBQS1JSWLtloOHDx82NjYuX/5s6lRVTc3jmhosnT569Oivv/7C4w1evGg3YIBDa2ulAbZv3GB5evYlcWa+rq6uoaGBmPsOAHZ29gIB6x//wPTWvo2GhgZLS0sC3uOoUc6pqQxHx1p9L1QqlaWlpfX19RgG4+EBAF4//vhUpVL9/HO9vb2VpWWFkZ/B4cNOY8Ywi4sJ2jKpVqtLSkqs8JeVUqmGZWVVeXnpt1heVVXV0tLSxd+VRgNHjgz78ceqO3f0XoavqqrqTuUKTbv8m5GRce3atdd+5+7uvnDhwqCgoIsXL3p4eEyePHn8+PErV67szNann36al5dH+3vWdnJyIiURSqXD6fR6BgPLr0R3aGtrk0qlLGL3rkulvs3NwS4uL9dvlEqlUqnUaIJ79hRZWGD8KKBQKFQqFR5fKoWiv0ZjwWRWGXCtSuWkVDr36kWaVKVCoWhra+tFlKK5SuWgUPSzsiJ0VCiVSi0sLCyNabXQPVQqB4XC1crqrr4XNjc329jY0LAeOMjlngBtcvldJnMEgIF/ou2RyQb16CG3tCQoETY1Ndna2mL+sbzJixcjLS0f0un6zcWpVCq5XG5tbd3ZCRqNhUQSZGPzuwEhKZXKjRs3zpo1q+vTXiZCgUBQWlr62u+cnJyEQmFRUdFnn30GAElJSd7e3gsWLMC1GASBQCAQCCJ5OTXK5XK5HZXcSaXSxsbGmzdvAkBTU9OjR48qKysJDRCBQCAQCDyhdX/bzVunRhEIBAKBMDv02FAfHBw8dOhQ/EJBIBAIBIJ49BgRIhAIBAJBPbCRWLv/d5qamjAxa15UV1e3F1d7/vz5n3/+SYzr1tbWhw8f6g6fPXtWW/v/pef19fUNDQ3YupBIJGKxWGki255Jpa6urv0ff/tPiTJoN+SQHUWnPH78uH14Uqn0wYMHxpt9+PCh7hvd1tYmFot14oJKpVIsFpv4KEKpVN6/f5/cINVqtVgsVr3SXtN+jM1viMmuWbOmvULhmx/va7e4m37lcrnuJ/X19U+fPu30Ao3RaHdpDB482OcVp06dMt6s2REbGzt58mTta5lM5uPjc+DAAWJcFxcXM5nMFy9eaA/Dw8Pt7OxUKpX2kMvl7tmzx0gXJ06cGDp0qPb1gwcPhgwZsnr1arVabaRZCjB79mxHR0fdH39UVBTZEWEPj8dLSkoiO4pO+eSTT3744QfdYV5enouLi/FmQ0JCtm3bpn2dn58PAPv379ceZmdne3h4GO8CV0QiEQDIZDISY9AOimpqajQajVwunz59OpfLbWpqeu00e3v7xsZG3WF5eTkASKVS3U9eu8XdYfr06dHR0drX9fX1Li4u58+f7+xkzES3T58+ffcVU6jRXERPtmzZUlVVpZUdWL9+vZub29eGNRPSnxEjRtjY2BQWFgKAWq2+fv26t7d3cXExAMhkssLCQh6Ph5UvsVg8fvz4L774IjExkYCdSWZBVFSU7o8/NTWV7HAQ2MDlcgUCgfY1n8/n8Xi6Q4FAMH78ePJCMz9evHgRHh4ulUrPnTtna2tLgMfk5OSsrKy8vDwA+O6776ZNmxYSEtLZyRRsw0QWTCYzNTU1Njb22LFjBw4cSE1NJSxP0Gg0LpfL5/MB4M6dO15eXiEhIdrDgoICa2vrESNGYOLo3r17HA4nJiZm48aNmBhEIEwWHo935cqVtrY2ABAIBOvWrct/pfzG5/M73G+G6JCmpqZJkyb16dPn5MmThIlOODo67tmzZ+HChSkpKUKhcOvWrV2cjBIhlgQGBs6ZM2fWrFk7duzw8PAg0rXu6VX7FeVwOLpDHo/XowcGN/rJkyccDmfr1q3LzEjAFIEwFDabrdFobt++rVAoSkpKeDyei4uLSCRqaWkpKirCcJaF8oSHh3t7ex89epRhmKC+ocycOdPPz2/RokWpqak2XbYmR4kQSzQazf379xkMhobwBWoej1dQUCCTyQQCAZfLDQoKKiwsVKvV2kNMXDCZzD59+giFQuLfHQJBPHQ6fcyYMXw+XygU+vn5MRgMDofD5/OvXLni6uo6cOBAsgM0G9zc3EpLS/FuBfUmCoWisrKyO/+QUSLEktTU1Orq6gsXLqxcufLRI4KEBLWMGDHC1ta2oKCgoKAgMDCwV69egwYNEgqFhYWFWC1m2NnZ8fn83Nzc5cuXo1yIMB2srKzaN19tbW3FSgtXuy7I5/M5HA4AaCdaBALBBJwapVOUQ4cOubi4hIaGvlkv2iHa29e+CF8ikXQhRtoZ8fHxTk5Ov/7667x587ruTowSIWbU1tauXbv24MGDPB5v7ty5ixYtItI7jUYbN25cUlKSp6enVvKbw+Fs27aNxWL5+Phg5cXd3f3y5cs5OTmxsbFY2UQgjGTIkCG3bt3SHd68edPb2xsTy9plwkuXLmmnVbQjQrRAqC8MBuPYsWPdz4Vubm42Nja6eyqXy+/du6fvPb19+/bevXtTUlI+//zzgICAtWvXdnGynn3vEZ2g0Wjmz58/f/78UaNGAUBCQoKfn9+RI0dmz55NWAw8Hi8mJmbVqlXaQw6Hk5iYOG3aNGxrdrS5UDvK3LVrF4aWEabM2bNndbtR/fz8IiIiyI2nPdHR0SNHjly6dGlwcHB5eXlSUtKZM2cwsRwQEKBWqwsKCthsNgDY29vb29sLhcKMjAxM7BNAQkKChbY9G0BUVFT//v1JCYPBYGRmZkZGRoaFhZ0/f77rFTsajRYXF7dgwYLVq1fb2toePXrUy8tr4sSJ3XenUCjmzp2bmJioncFOTk729fWdMWNGZ08wFps2bdLn7XQcNJ1O53K5BPS7Mlnq6uokEsmqVavodDoAMBiM4ODgmpoaf39/wmJwc3NzdnaOjIx0cnICgL59+9ra2s6aNQuTsh0LCws3Nzft2+ndu3dERER1dbWbm1ufPn2MN27W0On04cOHe3p6kh0IjtDpdBsbmx6vcHV1NSm1RRaLNWfOHJFIdOPGjV69eu3atSsgIAATyz169HB3dw8NDdUZdHd3Z7PZYWFhmNjHFRqNxmKx6HS67sb5+fn17t2b4DBYLNa4ceOYTKaFhcW0adOkUqlKpXrt+7Jt27Zly5b17NlT95OxY8cOGzbsxo0bIpGIy+Xu3r1br0IbkUhkZ2e3ePFi7TDA2tqazWY/efKks+kxJLGGQCAQCDJxcHAQi8XEJ2kdaI0QgUAgEO80KBEiEAgEgkyysrIMKArFEDQ1ikAgEIh3GjQiRCAQCMQ7DUqECAQCgXinQYkQgaA42dnZ9vb2YrGY7EAQCBMFJUIEguIoFIrnz59ruyggEIg3+T/uljbgToVbOgAAAABJRU5ErkJggg==",
"text/html": [
"\n",
"\n"
],
"image/svg+xml": [
"\n",
"\n"
]
},
"metadata": {},
"execution_count": 6
}
],
"cell_type": "code",
"source": [
"plot_bandstructure(scfres, kline_density=5, unit=:eV)"
],
"metadata": {},
"execution_count": 6
},
{
"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 Array{StaticArrays.SArray{Tuple{3},Float64,1,3},1}:\n [0.00045164925040840857, 0.0004516268445082516, 0.0004516246685735173]\n [-0.00045162684450825397, -0.00045164925040845773, -0.00045162466857404755]"
},
"metadata": {},
"execution_count": 7
}
],
"cell_type": "code",
"source": [
"compute_forces_cart(scfres)[1] # Select silicon forces"
],
"metadata": {},
"execution_count": 7
},
{
"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.5.3"
},
"kernelspec": {
"name": "julia-1.5",
"display_name": "Julia 1.5.3",
"language": "julia"
}
},
"nbformat": 4
}