{ "cells": [ { "cell_type": "markdown", "source": [ "# Tutorial" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "DFTK is a Julia package for playing with plane-wave\n", "density-functional theory algorithms. In its basic formulation it\n", "solves periodic Kohn-Sham equations.\n", "\n", "This document provides an overview of the structure of the code\n", "and how to access basic information about calculations.\n", "Basic familiarity with the concepts of plane-wave density functional theory\n", "is assumed throughout. Feel free to take a look at the\n", "[density-functional theory chapter](https://juliamolsim.github.io/DFTK.jl/dev/#density-functional-theory)\n", "for some introductory material on the topic.\n", "\n", "!!! note \"Convergence parameters in the documentation\"\n", " We use rough parameters in order to be able\n", " to automatically generate this documentation very quickly.\n", " Therefore results are far from converged.\n", " Tighter thresholds and larger grids should be used for more realistic results.\n", "\n", "For our discussion we will use the classic example of\n", "computing the LDA ground state of the\n", "[silicon crystal](https://www.materialsproject.org/materials/mp-149).\n", "Performing such a calculation roughly proceeds in three steps." ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "3×3 Matrix{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}}:\n 0.0 Å 2.7155 Å 2.7155 Å\n 2.7155 Å 0.0 Å 2.7155 Å\n 2.7155 Å 2.7155 Å 0.0 Å" }, "metadata": {}, "execution_count": 1 } ], "cell_type": "code", "source": [ "using DFTK\n", "using Plots\n", "using Unitful\n", "using UnitfulAtomic\n", "\n", "# 1. Define lattice and atomic positions\n", "a = 5.431u\"angstrom\" # Silicon lattice constant\n", "lattice = a / 2 * [[0 1 1.]; # Silicon lattice vectors\n", " [1 0 1.]; # specified column by column\n", " [1 1 0.]]" ], "metadata": {}, "execution_count": 1 }, { "cell_type": "markdown", "source": [ "By default, all numbers passed as arguments are assumed to be in atomic\n", "units. Quantities such as temperature, energy cutoffs, lattice vectors, and\n", "the k-point grid spacing can optionally be annotated with Unitful units,\n", "which are automatically converted to the atomic units used internally. For\n", "more details, see the [Unitful package\n", "documentation](https://painterqubits.github.io/Unitful.jl/stable/) and the\n", "[UnitfulAtomic.jl package](https://github.com/sostock/UnitfulAtomic.jl)." ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n Energy Eₙ-Eₙ₋₁ ρout-ρin α Diag\n", "--- --------------- --------- -------- ---- ----\n", " 1 -7.904938386565 NaN 1.95e-01 0.80 3.9\n", " 2 -7.909781227844 -4.84e-03 2.98e-02 0.80 1.0\n", " 3 -7.910052173409 -2.71e-04 3.04e-03 0.80 3.1\n", " 4 -7.910052839513 -6.66e-07 4.43e-04 0.80 2.4\n", " 5 -7.910052852934 -1.34e-08 3.04e-05 0.80 2.0\n", " 6 -7.910052854037 -1.10e-09 5.59e-06 0.80 3.1\n" ] } ], "cell_type": "code", "source": [ "# Load HGH pseudopotential for Silicon\n", "Si = ElementPsp(:Si, psp=load_psp(\"hgh/lda/Si-q4\"))\n", "\n", "# Specify type and positions of atoms\n", "atoms = [Si => [ones(3)/8, -ones(3)/8]]\n", "\n", "# 2. Select model and basis\n", "model = model_LDA(lattice, atoms)\n", "kgrid = [4, 4, 4] # k-point grid (Regular Monkhorst-Pack grid)\n", "Ecut = 7 # kinetic energy cutoff\n", "# Ecut = 190.5u\"eV\" # Could also use eV or other energy-compatible units\n", "basis = PlaneWaveBasis(model; Ecut, kgrid)\n", "# Note the implicit passing of keyword arguments here:\n", "# this is equivalent to PlaneWaveBasis(model; Ecut=Ecut, kgrid=kgrid)\n", "\n", "# 3. Run the SCF procedure to obtain the ground state\n", "scfres = self_consistent_field(basis, tol=1e-8);" ], "metadata": {}, "execution_count": 2 }, { "cell_type": "markdown", "source": [ "That's it! Now you can get various quantities from the result of the SCF.\n", "For instance, the different components of the energy:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "Energy breakdown:\n Kinetic 3.0795140 \n AtomicLocal -2.1806344\n AtomicNonlocal 1.7339382 \n Ewald -8.3979253\n PspCorrection -0.2946254\n Hartree 0.5417559 \n Xc -2.3920760\n\n total -7.910052854037\n" }, "metadata": {}, "execution_count": 3 } ], "cell_type": "code", "source": [ "scfres.energies" ], "metadata": {}, "execution_count": 3 }, { "cell_type": "markdown", "source": [ "Eigenvalues:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "7×10 Matrix{Float64}:\n -0.170182 -0.1318 -0.0883279 … -0.0562604 -0.11494 -0.0700165\n 0.201343 0.0909039 0.0122927 0.0111148 0.0420616 0.0176498\n 0.249296 0.174773 0.176137 0.132964 0.220115 0.11233\n 0.249296 0.231428 0.20237 0.161043 0.220115 0.190456\n 0.350986 0.360027 0.340134 0.291807 0.320727 0.327376\n 0.36997 0.395898 0.389483 … 0.331814 0.38819 0.460226\n 0.36997 0.401676 0.412474 0.565536 0.38819 0.462717" }, "metadata": {}, "execution_count": 4 } ], "cell_type": "code", "source": [ "hcat(scfres.eigenvalues...)" ], "metadata": {}, "execution_count": 4 }, { "cell_type": "markdown", "source": [ "`eigenvalues` is an array (indexed by kpoints) of arrays (indexed by\n", "eigenvalue number). The \"splatting\" operation `...` calls `hcat`\n", "with all the inner arrays as arguments, which collects them into a\n", "matrix.\n", "\n", "The resulting matrix is 7 (number of computed eigenvalues) by 8\n", "(number of kpoints). There are 7 eigenvalues per kpoint because\n", "there are 4 occupied states in the system (4 valence electrons per\n", "silicon atom, two atoms per unit cell, and paired spins), and the\n", "eigensolver gives itself some breathing room by computing some extra\n", "states (see `n_ep_extra` argument to `self_consistent_field`).\n", "\n", "We can check the occupations ..." ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "7×10 Matrix{Float64}:\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0" }, "metadata": {}, "execution_count": 5 } ], "cell_type": "code", "source": [ "hcat(scfres.occupation...)" ], "metadata": {}, "execution_count": 5 }, { "cell_type": "markdown", "source": [ "... and density, where we use that the density objects in DFTK are\n", "indexed as ρ[iσ, ix, iy, iz], i.e. first in the spin component and then\n", "in the 3-dimensional real-space grid." ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=1}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdeUBUVdsA8OfcmWFHdhg2QRFQUBERkcQNt7RMzTRL3xZzzdQ003yzLLNy6zWzPrNSy0qt3kwtS0VwQ3BfcAMRRdmGfV9n5p7vj/ElQlDAuXO35/eXczsMD7c589x77nnOIZRSQAghhOSK4TsAhBBCiE+YCBFCCMkaJkKEEEKyhokQIYSQrGEiRAghJGuYCBFCCMkaJkKEEEKyhokQIYSQrGEiRAghJGuYCBFCCMma1BJhRkZGXl5eCxvr9XpOg5EkPGltgCetDfCktQGetLaRWiJcs2bNTz/91MLGVVVVnAYjSXjS2gBPWhvgSWsDPGltI7VEiBBCCLUKJkKEEEKyhokQIYSQrGEiRAghJGuYCBFCCMkaJkKEEEKyhokQIYSQrGEiREimsrKy+o4Y+8SEF7D4DMkcJkKE5Ci9nE7/bFeC85ADpY5v/ng8r5rvgBDiDyZChGSksBb+7xob9buu9x6dQ6/hvsn/7VR9s9AjPPAX7cgDuh9ushVavkNEyOSUfAeAEOJcjR5istjvU+mBTLa/O5nXlRntw5gxATD+cKMGs09o+7uTF/wNDfiNGiETwUSIkGSxFBJy6fc32Z9vscEO5AV/ZnN/la2qiZYWChjVnhnVHoprFb/fZb9KZmfG65/wZl7wZwZ7EmLyyBEyJQ4TYXp6emxsrL29/ZNPPmlubn5/g4qKij/++KOmpmb48OHu7u71x+Pj41NSUpycnB5//HELCwsAKCgouHjxYn2DHj16ODs7cxc5QmJ3tZh+f5P97gbraE5e8Geuj1epLVv0gw7m8II/84I/k1FJd92mi07rC2pgrC95KYAJdcKEiKSJq0R4/Pjx0aNHP/vss8nJyZ988snRo0dVqn9ciJaVlUVERPj5+bm5uS1cuPDo0aPBwcEA8Morr5w4cWL06NGXLl166623Tp06ZWdnl5iY+MILL/Tq1cvwsx9++CEmQoTud6eC7kyjW26wtXqY6EeOPqkMsGtj9vK2JvO6knldmavF9Jfb7LhDeksFjO9IXvRnOthiRkSSwlUifP/995cuXbpgwQKdTtejR4/du3ePHz++YYOtW7eq1erff/+dEOLs7Lxy5crvv/++vLx869atycnJAQEBlNLOnTvv37//2WefBYDg4OCYmBiOokVI1Apr4dfb7LZUNrmEjuvAbO6n6Ks22nhmsAMJdlC8GwoJufSX22yfvTr/dmR8B2ZSJ8bZwki/AyFecfI0vLq6Oi4u7umnnwYApVI5atSoP//8s1GbP//8c8yYMYQQABg7duy+ffsAwNzc3N7evqSkBADq6uqqqqrUarWhfUVFxZ49e44fP15TU8NFzAiJTrUOfrnNjjqo67hTeyiLLg5hNJNUm6IUUcbLgvUYAlFqsj5ScXeianEIc66ABvyiHXVQty2VrdIZ+5chZFqc3BFmZ2dTSj08PAwvPT09z50716hNVlZWwwbFxcVVVVVWVla//fbbiy++GBgYmJycPGfOnAEDBhjaEEJ++OGH5OTkqqqqv/76KyAgoMlfrdPpDh48WF5e3vDg1KlTHRwc7m+s1Wq1Wpwt3jp40trAiCdt9fovLl27MXrmW39VuP+RCb2cYVJH+L4fWCtZAKB6vZbjLcoZgMfd4XF3KK2D3zPoz7f0C07qH/eEcO3VQ99tGPfE0EkTxhnlF+EnrQ3wpN1PqVSSh10ZcpIIWZYFgPrfzTCMXt+4d1JKGzYw/BTLsitXrgwNDX3uueeuXr26cePG8ePHd+jQ4Yknnhg1apThp6ZOnbpw4cK9e/c299urqqqKi4sbHtHpdIaQ7o+zyePoAfCktYGxTlpqaurqHfvLuj9zfOWGpcs/XNOLOv5vFprp/5/YKuH5DvB8B8iqgp/Tydsz3qvsPzv+3TfGjxmlVBrhiwU/aW2AJ61tOEmEhlu93NxcLy8vANBoNPU3f/Xc3d3z8vIM/9ZoNLa2tjY2NseOHTt16lReXp5hQDUpKWnTpk0rV640ZEoAIISMGzdu1qxZzf49SuWYMWPmzJnTkjjr6uqanM6KHgBPWhsY66S19+1YW1Vln7Dpq3XLR3U1e/Q3NIqO5vCWAyj/9fj7G/6tVPtbW1sb5W3xk9YGeNLahpNnhNbW1hEREfv37wcASumBAweio6MBQKfTGUZNASA6OtrQAAD2798/ePBgADAzM9NqtXV1dYbjlZWVZmaNe/uZM2d8fX25CBshgfvhrkXfT+LyriSMGjmC71gaWzhnpuZygsPrP8VkUb5jQah1uJo1unTp0hdffFGj0Vy5cqW0tHTixIkAkJqaGhQUVFxcbG9vP23atC+++GL69Onu7u7r168/cOAAAPTq1atbt24jR44cN27clStXjh079sknnwDA7NmzdTpd+/btr1+/vnfv3geMiyIkVSV18MEF/cERSpVKoNUL1uaqVeHs/JP6i2OVSlyVBokHV5/WJ5988sCBA3q9vl+/fidPnrS0tAQADw+P7777zjBy4uLicv78+c6dO6tUqoSEhIiICABQKpVHjhyZPn16YWFh9+7dk5OTO3XqBACzZ8/u1q1bTU1N3759r127NnDgQI7CRkiwPrigH+PLdHMUaBY0GOvLuFvBlhv4mAqJCTEMVErG3Llz/f39W/iMsLy83NbWluuQJAZPWhs8+klLK6N99uouj2vpAjE8ulREH/9LlzxeZfdoDzHxk9YGeNLaBscvEBKBBafYxSEK4WdBAAhxJE+2Zz66yHEZB0LGg4kQIaE7nEOvFNE5QaLprR/2Umy9waaWSmq0CUmYaLoWQvKkp/B6on5tBGOu4DuUFnO1hPndFEvO4pNCJA6YCBEStG9SWEdzGOsrsq76RjfmUiE9hKUUSAxE1rsQkpVyLSw/z66NEM/N4P+YMfBROPPmab0eUyESPPkmwtra2qNHj5aWlvIdCELN+uCC/on2JMxZ0CUTzRnfgWmngu9ScYAUCZ18E+GYya88vzGme//hfAeCUNNuldOtN9jlYeK7Haz3nz6Kd86y5bgKNBI2+SbCqqrqGiu3rDJt9126jy6yt8pxBAcJy6LT7BvdxFEy0ZwwZzLMi6y8hKUUSNDkmwj3bt+yaaT7nfi9W/srcqtp37264P/qVl1is6swIyL+Hcmh5wvo611F30NXhiu+SWHT8UITCZjou1mb2dnZPfPMM56enmHOZH2kIut51aYoRXYVDf1NF/W7bv0VNq+a7xCRXLEUFp7Sr+7NWIh4WPQeN0uYE6RYfAafFCLhkm8ibKTJDbijftd9lcyW4RMOZFpbb7DmChjXQSLdc2F35nQ+Pa7Bm0IkUBLpaUZkroBR7ZltAxU5z6sWhzCHsqj3du2og7ptqWylju/gkAxUaGHZefbTPgpRThVtioUCPurFzEvUs5gKkSBhImyWpRJGtWd+Hqy485xqfAfml9usx4/aCbH63++ydTjMgzjz0UX9ME8S7iKZPAgAMNGPsVbB9zex5yAhwkT4cPZm8II/8/sw5a2Jqifbk8+ush4/al84ov/9LqvDfo2MKqOSbr7BrugltY5JANb3USw9i8MqSIik1t845WQOL/gzMSOUF55WhjmTVZdY35908xL18Rp6My1tx46dFRUVfMeIxG3hKXZusMLDSlK3gwY9nclAd7IaSynQo8nPz//hhx81Go0R3xMTYVt4W5N5XZn4UcpDIxQO5jDlSG3nQWNf/OnaS7MX8B0aErETufRUHl0g/pKJ5nwczvzfdTajEh8VorYbMnbiK39pBj01wYjvKdkuZxqd7cl7PRXJE8ycLRl9iaaKEXPxM+IVS2H+Sf3H4Yylku9QOONlTV4LUrx1Gp8ooLarJJa0KMvSwsKI7yndPmdCDMNcP3n40wOX9qoi9RSkM9sPmdC2VFZJYKKfxK9NF3VnuvxXF6+hUWrsJ6jVyrVQPW3nZ7anX3z8QyO+rcR7nck4ODi8P3GgvYVicwpe7aJWq9LBsvPsp5HSv4iyVMKKXszC07gpBWqL987rn+xoOXNstKWlMYffMBEa0+ePKd45py+s5TsOJDYrL+kHqElvaZVMNGdSJ0ZJYDuWUqBWul5Cf7jJruhl/PWWMBEaU7ADeaYDs/w8zotDrZBZSTdel2DJRHMIwNoIxZIzWEqBWmfBSf3bPRQuxnw4eI9c+p7JrAhT/HSLTSrCgR/UUotOs68FKdrbyOJ20KCPK4lSk7VJeFOIWmr3HfZuBczqwknOwkRoZA7m8G6oYv5JvClELXIyj8Zr6MLusuuJq3szn1/T363AS0b0cLV6WHSaXRepUHHTUWTX/UxgRhemsAb+exuvdtFDUICFp/QfhzPW8pu+7WVNZnZh3jmH3QQ93JokNsSRDPPkatQEE6HxKQh80Vfxxil8BIIe4sebrI7C851k2g2XhCiO5NDT+XhTiB4ks5Kuv6pf05vDbiLTHsi1vm7kMTeyJgkHSFGzqnWw9Cy7trf0SyaaY6WE93syrydiKQV6kDdPs7ODGF9bDjsKJkKurOnNfHENN+ZGzVqdxPZ1IzKvK3/Bn9FR+PkWDpCipp3IpQm59M3u3G5RjYmQK17WZF6w4k1cTQo1JauSfn5N/1G43DsgQ2B9pOLNU2wVPkdA99FTmH1C/0kE5w/R5d4POfVmd+ZSET2YhTeFqLElZ9hZXRgfOZVMNCfSlfRxJeuu4CUjamzTddbJAp7pwHmewkTIIXMFrO7NzE/Ua7GPowbOF9DDOXQRx6M9IrKqN/PpFX12FV4yor8V18LyC/p1fUzRTTARcmuMD9PeBjZex0yI/jbvpH5FL8ZGxXccgtHBlkwLZJaexW6C/rb0nP7Zjkx3R1OMmmAi5Nx/+ig+vKjPr+E7DiQMO9PYSi38S64lE835dw/FwSx6tgBvChEAwNVi+t/b7Ls9TTRqgr2Rc13syeROzNKzWEqBoEYPS86y6yMVDD4c/CcbFbzfk5mHpRQIAABeS9B/EKZwMjfRr8NEaArv9VTsy6BnsHBY9tYmsb1dSD95l0w05+UAplYPv+KSTLK3M40trYNXAk2XnjARmoKtCpaH4dWu3OVWw/qr+o9ks8tEazEE1kYoFp1ma3D0RMaqdfDWGfbTSIUpV5rAPmkiL/kzLMCPuAebjL11Rj+9M+PXDm8HmzXQnYQ6kfVYSiFjH13S91OT/qYdNcFEaCIMgS8eUyw+zZZp+Q4F8eFCIT2QyS4OwZKJh1gTway9rM+p4jsOxIdb5fTL6+xKky80gYnQdMKcyTAv8vFFHPeRo9cT9R+EKdphycTDdLQlLwcwy3B3a1l64yT7RjeFp7WpR00wEZrUynDF5hT2Rik+K5SRtLS0GR9/XajJfCkAu1uLvBOq2JN4/e113+Tm5vIdCzKd2Gx6uZjO78pDN8GeaVJulrAoRPHGKbzalZEhY5/fobGr/mqKfLeZaCVbFdRumvR5msW4F2fwHQsyES0LcxP06/ow5nw8PcBEaGpzg5nUUvgzA28K5cLSzkF5M76Tt5rvQMTEw9EOUk+4unvwHQgykc+vsT62MKo9PylJfhtj882Mgc8fU8w6oR/soeTl2geZ2HOf7r126cKPL4XzHYiYXIw/NPCri8+PDOM7EGQKedXw8UX98VG85SOu0u+NGzf69evn4ODQq1evc+fO3d+gtrZ22rRpzs7OPj4+X375Zf3xb775JiAgwNnZOTw8PDY21nBQr9cvWLDA1dXV09Nz1apVHMVsMkM8SZADWX8V54jLwpFcZlxUN4bB0ZdWMDMze7JPt8MavuNAJrHkjP7lACbQjreHB1x1zueffz46Olqj0UyZMuXpp5/W6RrvNrZ69err16+npqbu3bv37bffPnv2LABcunRp3rx5O3bsKCgomD179tNPP11bWwsAX3311aFDhy5fvnzs2LENGzbs37+fo7BNZl0fZk0SLrcvfTV6OJNPH3PBi55WG+Cqj8UtzGTgXAH9K5P9dw8+x8c4SYQXL15MSUlZsmSJubn5rFmzKKUxMTGN2mzevHnJkiUODg4hISHPPffcli1bAOD27dve3t5hYWEA8Mwzz5SVleXn5xsaL1iwwM3Nzc/Pb9q0aZs3b+YibFPqaEumd2aWnMHvR4k7kUu7ORJbJX6ht1o3B1pUSzMr8dRJGQV4LUG/MlxhZ8ZnGJwkwhs3bgQEBFhYWAAAIaRbt24pKSkNG9TU1Ny5c6d79+6Gl927dzc0GDJkiFKpXLNmzbFjx+bPnz9u3DgvLy8ASElJub+x2L3dQ3E0hx7XYD+XsrhsNtoDZ4u2BQEY4M4czsEOImXf3WAJwL/8eX5wwMnDyaKiIhsbm/qXdnZ2hYWFjRoAgK2tbX2DgoICALC2tp48efKGDRv27dt3586dtWvXAkBdXV1FRUV943bt2hkaN6mmpmbu3Llz585tePDixYt+fn73N66srCSEzy+pd7sxcxOUR4fWiWgvAt5PmrjEZKjeD9HjSWuDysrKKCflgTtkrLrxgxXUHHF90ip05O2zZtujtJUVHG5TZ2Vl9dAn9JwkQkdHx/Ly8vqXpaWlzs7ODRs4OTkRQsrKyuzt7QGgpKTExcUFAHbu3Ll169bk5GQrK6v09PRu3bp17do1MDDQxsamrKys/t0MjZtkYWHx2WefzZkzpyVxUkobJmzTmxIM29J1P2dbTTXhOuuPiPeTJiLlWkgu0w7yMddW6fGktRaldEQHi0+u621sLPiORTTE1T3fP6Uf6Q0D2lvxHQg3Q6MBAQGpqak1NTUAQCm9cuVKQEBAwwbm5ube3t5XrlwxvLxy5Yq/vz8AXL16tWfPnlZWVgDg6+vr4eFx7do1wxve31gaPu2jeOesvqSO7zgQB47ksBGuxAKLZNoqwI4wBFJxJSYpullGv0tlP+gliO7BSSLs0aNHQEDAqlWr6urqvvrqK0rp0KFDAeDQoUNLliwxtHnllVdWrlxZWlp65cqVHTt2TJkyBQB69+596NAhQ87bu3dvRkZGaGioofG6devy8/Nv37799ddfGxpLQw8nMsaXeR9XVpSiuGwa7SGae31hGuhO4vAxoRTNTdT/u4dCbcl3HADAXUH99u3bp06dum7dOj8/v127dimVSgAoLi6+deuWocGiRYvu3r3bsWNHGxubDz74IDw8HACeeuqpRYsWjR07tqioyMfHZ/v27b6+vgAwY8aM1NTU4OBgpVI5d+7ckSNHchQ2Lz7spQj6r3ZKANPNUTSD+6gl4rLp1/0wET6SaA+yL4PO6Mx3HMio9t5h08thdpBQegehVFJXW3PnzvX392/hM8Ly8vL6OTj8+vwa+1s6GztSBAv9COekCVxeNQT+oi34l0pB8KS1heGkZVfRkF263EkqEU0o45EoPml1LHT7VfdZpGK4l1D+pwolIcvcrC5MYQ3sSseyQuk4nMP2d2dwoe1H5GFFHM3J5WJJXa/L3NoktqsDEU4WBEyEAqEg8Pljivkn2SqcKC4Vcdk02l1AXV28BnuQuGxMhBKRVUnXXdGvjRBW6hFWNHIWpSZ9XMnay3hTKBGx2XSwJyZCI4j2ILFZ2C8kYtFp9tUuTAdbYXUNTIQC8kkEs+GqPr0cL35F724FLdfSYAdh9XaRGuTBxOdSLaZC8UvIpcc1dFGIIEomGsJEKCBe1mROsGIxLkAqfrHZNNoDp3cYh5M5dLAlZwvwAlHcWArzEvVrIxhr4U0KxEQoLIu6M2fyKa6vKHZx2XQwLjFqPIM9CO5EIXZfJbMWShjfUYhJR4gxyZmFAtZGMC/857+btnx3/95VSCwO51Bca9uIoj2YuGwcKRGrioqKj9dvfPvHuM8fE+g0akyEguN057gm7vu5247/sH0H37GgtkguoSoGOgpsOoCo9VeTswW0Gq8MxentFauWHs6p+nmpS00O37E0DROh4Li5uVmXZdKc6+29vfiOBbVFLBZOGJuNCro7khO5ODoqSh3ae8OtM+2gWrALggvvqaXsde7cOen4weAdFb2jfPmOBbVFXDYd1wEToZFFe5DDOewQT8FNOEQPNeS5ae7aAcmveAo2EeIdoRC1VzuHBXjhnr1ixFI4pmEH4h2hsUV7MLFYVi9OMVn0ybBOgs2CgIlQsIZ6MjFYRCxCFwqpqwXxsMJEaGSPuZHrxRQ3LBOjmCx2qLAXl8BEKFBDPUkMzhcXoThcUIYbZgxEuJLjGrw6FJk6Fk5o6CBh70cm6ODkLMyZ5FTRnCq+40CtFJvNYuEER3B0VIwScmlne+JoznccD4SJUKAUBAZ5MIewdkpU6lhIzKUD1NitOBGNq2+LUEwWO0xIG000CXuscA31JDGZ2O3F5FQeDbAjDsK++BWvMGeSVUk11XzHgVojJosO9RR6ohF6fHI2zJMczGIxE4pIbDaLDwi5oyDQT80cwWES8SiuhZQSGuEi9E6BiVC4fG2JjYpcKcJUKBpx2TRa2JMCxC7ag8ThSrzicSib7acm5oIv/sROK2g4d1REqnRwoZD2dRP6xa+o4Sa94iKKcVHARChwQz0JVhOKxXENDXMmAtxiRkqCHEilluKenWJxKIsOFfxMGcBEKHDRHkxCLq3V8x0HaoG4bBbHRblGAAZ5MDg6Kgo3y2itHrrYYyJEj8beDLo4kIQ87PYigGttmwYWUYhFTBYd6klE0SUwEQrdME8Sk4mjo0JXUgeppbS3qyh6vbgN9SSHcDa1GBgSId9RtAgmQqEb6sngfBnhO5zN9lUTM+xP3POxIVZKklyCnULQ9BSO5rCDxTBTBjARCl+kK7lZRvNr+I4DPVAsFk6YULQHwbXWBO50PvW2JmpLvuNoGey6QqdkIEpNDmMRsbDF4QNCE8LHhMInonFRwEQoCjg6KnA5VZBXTXs4iabbi91gD+ZIDqvHPiFgMVnsUC/R5BfRBCpnWFYvcLHZ7EAPhsE8aCpuluBhRS4UYqcQqHItXCqkUeJZXAIToQh0sScshdRS7PYCheOipodLzAjZkRy2twuxEs/iEpgIxWEw3hQK2OEc3IzX1KI9SBw+OBcqsaysVk9MscoZjo4KVloZrdNDoB0mQpMa6I6LLgmXuGbKACZCsRjmyRzJYbV4BSw8sdl0MG5Jb3J2ZtDZnpzKx6tDwcmqpAU1Ips7holQHJwtwNeWnMFuLzxx2TQaEyEfBuPoqCAdyKJDPEU2dwwToWjg6KgAUYAjOewgTIR8iPZgcL6MAIluXBQwEYrIUE8Gt2QSmitFtJ0Z8bERWbeXhr5u5GIhrdDyHQdqgBq2YRHbJGpMhKLRT02SimhpHd9xoAZwxwkeWSmhpzOJz8WbQgG5WEgdzIivrcg6BSZC0bBQQIQrOZqDN4UCgg8I+RXtweBjQkER47goYCIUF1xrTVD0FOJz2UG41jZ/BuPq2wITk8ViIkTcwvkygnImn7a3Ji4WfMchYxEuJK2MFuDeLMJQo4dTeXSAu/jSivgilrMeTqS4jqaXYy4UhLhsXFCGZ0oG+rqRYxocHRWE4xoa4kTszPiOo/UwEYoJARjsweBYkEDEZbO4ByHvsIhCOGKyWHGtrFZPlEHLGY6OCkSNHk7n035qvCPkGW7SKxwHM0U5UwYAOFwe/Ny5c/v27XN2dn7++eft7e3vb5Cbm7t9+/aqqqqnn366S5cuAFBaWvrTTz81bBMVFRUUFHTnzp0DBw7UHxwxYoS3tzd3kQvZME+y6LSepQpxLdwgPQm5tKsDaafiOw7ZC3EkhTU0q5J6WmOX4FNBDaRX0F7Oovy/wNUd4b59+4YOHQoAR48ejYyMrK6ubtSgoKCgZ8+eV69erampiYyMPH36NADU1NSc+5+EhIQZM2bk5OQAQFJS0rJly+r/U2lpKUdhC5+nNXG2IBdxJza+xWWz+IBQCBgCA9yZwznYI3h2MIsd5M6oxDnIyNUd4YoVK1auXDl9+nRKaURExE8//fTSSy81bPD111+HhoZ+8803AGBmZrZy5cpdu3a5ublt2rTJ0OCnn346evTooEGDDC/9/Pzq/5PMGUZHe4rzyksy4rLph+EKvqNAAPe2ZKKTO/Edh7yJtILQgJP0XVlZefLkyZEjRwIAIeTxxx+PjY1t1CY2NnbEiBGGf48cOfL+Blu2bHn55ZcZ5l6ERUVFmzZt+vnnnwsLC7mIWUSGehJca41f5Vq4Ukz7uIi120sMPiYUglgxJ0JO7ggN45murq6Gl2q1OiEh4f42bm5uhn+7ubmVlZVVVFTY2NgYjmRmZh4+fPirr74yvDQ3N+/YsePVq1eTk5NfffXVAwcOhIWFNfmrdTrd7t2779y50/Dg/PnznZyc7m9cU1OjUonvIU8fB5iUryiqqOFlA2iRnjTjiskkvZwI0dXU6FrUHk9aG7T8pPmYg55VXM2v8bPlOiih4+uTllwKBBTe5rU1wqvpNDMzq7+hag4nX6WG30rpvWs0vV6vUDQeRGIYhmXv3dYY/tGwzebNm6Ojo318fAwvhw0bNmzYMMO/Fy5c+NZbb8XExDT3262srBwdHRseUalU9wdg+I1NHhc4ewWEOMLJQsVQDx5+u0hPmnEdy4NB7tDy84AnrQ1addIGucPRPEVAE3Py5IWvT9rhXBjq2YoeYUqEPPw+lZNE6O7uTgjJycnx9fUFgJycHHd390ZtPDw8DDeOAJCdne3g4GBpaWl4SSndtm3bxx9/3OSbDx48+JdffmnuVyuVymHDhs2ZM6clcapUKpFepw/zYg9r6EgfHj524j1pRnRYo/sqSqFStXQgCE9aG7TqpA3xYv/MoK8GC/GL2JT4+qTFaXQv+DMqkU6V4egZoaWl5aBBg3bv3g0AOp3ujz/+MDwvrKmpOXv2rF6vB4ARI0bs2bPHcNe4Z8+e+ueFABAbG1tSUjJ69Oj6Izrd3yNQhw4dCgwM5CJsERnqSQ5iNSFPCmvhTgUNw8lKQjLEg8Rlsyz2CT5oWTiuoYNEuLJaPa6eMr4NtNoAACAASURBVC1btmzMmDGpqanJycnW1tZjxowBgNu3b4eHhxcXF9vb27/88stffvnl6NGj3dzcfvvttyNHjtT/7ObNmydPnmxubl5/ZMqUKcXFxT4+PtevX798+XLDmkJ5CnchmZU0pwrcrfgORX5is9gB7oxSxL1egjytiaM5uVJMuzviBYqpJeZR/3bEWcyL7nKVCPv373/u3LnY2Njo6Ognn3zSzMwMANq3b3/w4EHDjBg7O7uzZ8/+8ccfVVVVy5cvbzh2+uqrrza651u9evWJEycKCgqGDBkSHR3drl07jsIWCwWBAe5MXDY7qRN+H5taHO5BKEiGuaOYCE1PpDtONMThvMMOHTpMnTq14RFra2tDlb2BjY3NxIkT7//Bfv36NTqiVqvHjRvHRZDiZagmnIS1UyYXm01nB+H1h+BEe5Btqez8rvi/xtRisujHIq+pxQ+NWBkSIT4TMbGMSlqmpV3xtkN4oj2Y4xqqxQpb0yqpg+vF9DFXcfcITIRi1akdMVfA9RJMhSYVm0WjPXCdVyFyMocOtuRcAfYIk4rLZh9zI+biviHERChmQz1JTCZ2e5PCB4RChkvMmF5MFhXp1ksNif4PkDNca830DufQaA9MhAIV7cHEZWOPMKmYLDrMS/Q9AhOhiA32YI5raK2e7zhkI7mEMgT82om+20tVfzU5k0+rW7buHXp06eW0UkuDHUTfIzARipiDOQTak5N5OBZkIrHZdAjeDgqYrQq6OZIE7BGmciCLDvWUwiNzTITihqOjphSXjeOiQhftQXB01GRisuhQ8Y+LAiZCsRvqycTgWmsmwVI4pmEH4kwZYRvswcThfBmT0FM4ks0O8ZBCEpHC3yBnj7mR5BJaVMt3HDJwsZC6WhBPa0yEgvaYG7lWTEvr+I5DBs4VUHcrIo1VHjERipsZA33V5DCOBXEvFsdFxcCMgd6u5JgGewTnRL0lfSOYCEUPR0dNIy6bxUQoCtE4OmoSMVmsBCoIDSTyZ8jZUE9yABMhx7QsJOTSAWLeaEY+BnsQTIRcq9TB+QLaTy2RS0Ps2KLX1YHoWEgrw57PoVN5NMCOOJo/vCXiXZgzyaikudV8xyFpR3JoL2diI5XdpjERSkG0B8HRUU7hA0IRURDop2aO5OBjQg7FZLFDvaSTPqTzl8iZYScKvqOQsrhsNloS08RlItodR0e5FZMpnZkygIlQGoZ6ModzWD12fG5U6eBCIY2SyuMQORjsiatvcyi7iuZU01An6fQITIRS4GYJXtbkbD72fE4c19BQJ2LN4SbWyMiCHUillqaXY4/gREwWHeLBKKSTBzERSgWOjnIHx0VFhwAM9GAO52CP4ISUKggNsHtLxFBPBhcd5UhcNh2MM2XEBosoOEIB4rLZIZgIkQD1V5PzhbRMy3ccklNSBzdKaW9XSXV7OYj2ILHZLGZCo7tcRG1UpIOtpHoEJkKJsFJCbxdyDMeCjO1wNvuYGzHDjiI2HW2JhYIkl2CPMLKDkhsXBUyEUoKjo1yIy6b4gFCkonF0lAMxmSwmQiRcOF+GC/iAULwwERpdrR5O5tEBaqklDqn9PXIW6kQKamhGJfZ8o8mthtxq2kNC9VKyMtgD62uNLD6XBjsQB8mtNYiJUDoYAoM8mFi8KTSeQ1nsAHeGwTwoTm6W4G5JLhZijzCamCwJjosCJkKJwdFR44rDJUZFbrAnjo4aU0wWlczWSw1J8E+Ss6Ge5BBOGTeeuBx8QChuhiIKvqOQiMJaSCujEVIsJcJEKCk+NsTejCQVYSo0grQyWqunne0l2O3lY6A7k5BLa/V8xyEJhicFKikmDSn+TfKGo6PGEpdNB2PhhMjZm0GgHTmNy/Aag/RWVquH/VxqhnqSmEwcCzICHBeVhsGeODpqHLHZmAiRSER7MIl5tFrHdxwiRwGOZLODMBGKX7QHg/NlHl1KKdWxEGgnzR6BiVBqbFXQzZGcyMWe/0iuFFEbFfGxkWa3l5UoN3KxkFbipeGjicmiwyR6OwiYCCVpqCfBtdYeRUlJybtrv+hafJrvQJARWCnBN33/v9dtrq2t5TsWEZPwA0LARChJQz0ZnC/zKF5b9M6e1IrY1a+WlZXxHQt6VJcuXbr56/r/i0n6/Muv+Y5FrHQsHMthB0uxgtBAsn+YnEW4kNvlNK+a7zhEq4OPN0lLtFWBhYUF37GgR+Xo6GhRlQ9Zl328PfmORaxO5dMOtsRFur0BE6EEKRno787E4Uy5tvrXa2+6T/7w5oVEMzMzvmNBj8rb2/ty4hGzlzaNGj2W71jESqorq9XDRChNWE34KBLzaP/QLpaWlnwHgozD08UhsFMHXHS0zaS6slo9Kf9tcjbMkxzERNhWibk0UorrSMlZpBtJzMMe0RZlWrhcRPu6SblHPCQRVlZW3rhxo7S01DTRIGMJsCNKBnB77rZJzKORku72MhTpiomwjeKy2UhXYqnkOw4uNZsIv/3224CAABsbm8DAQHt7+/bt23/66aeU4idJNAZ74OhoW5RrIa2MdnfERCgpfVxJIhbXtklMFh3qJfGxw6b/vDVr1rz88suOjo5r1qzZtm3bunXrAgMD58+f/8Ybb5g4PtRm+JiwbU7n01BnYibxji87ndqROpZm4rbVrSftCkKDJm536+rqVqxYMWPGjC+//LL+4Ouvv75ixYr33ntvyZIlLi4uLXz30tJSOzu7BzSoqqpSqVQqlaol71ZTU8MwDE7ka6EhnsyMeK2WVUhytXju4ANCqertwpzMo890wP+5rXCngpbWSX+ApInvyPz8/LKystmzZzc6/uqrr+r1+jt37rTkfc+cOePv7+/v7+/l5XXo0KH7G5SVlT3xxBOenp7Ozs7vvfee4eC2bdsc/ykjIwMAamtrJ06cqFarXVxc5s2bhyO0LeFkDn7tyCl8LtJKiXksJkJJwvkybRCTRYd4MJLvD00kQhsbG0JIYWFho+MFBQUKhcLHx+ehb0opfeGFF15//fW8vLwNGzZMmjTp/sWNPvzwQwDIz8+/du3apk2bjh49CgDPPvts2v8sXry4Q4cO3t7eALBhw4aMjAyNRnP79u0///zz119/bdtfKze41lprUYDT+dLcehRFupKTmAhbSQ7jotBkIrSzsxs8ePCcOXPS0tLqD2ZnZ0+fPn327NktGRc9c+aMRqOZPn06AIwdO7Zdu3b79+9v1Gbbtm3z589XKpWenp6TJk3atm0bAJibmzv8z88///zKK68YGn/33XevvfaahYWFo6Pj1KlTDY3RQ+Faa62VWkptVMTDSvo9X4bCXcilQtyktxVYCoez2cEySIRNT4ndunXr4MGDAwICQkJC1Gp1QUHBpUuXzMzMXFxcJkyYYGgzc+bM6OjoJn88LS2tU6dO9U/+Onfu3DCnAkBVVZVGo+nSpUt9gx9++KFhg6SkpGvXrj333HOGl7du3WrY+AGJkFJaWFh469athge9vb1b+BhSYvqpydViWlwLDuZ8hyISiXn4gFCyrJUQYEcuFuIdf0udL6QuFsTbWvqnq+lE6OXlde7cue+///7IkSOZmZksy3br1g0Abt++Xd/mAesRl5aWWllZ1b+0tbUtLi5u2KCkpAQA6tvY2NgYjtT75ptvnn76aQcHBwCoq6urqqp6QOOG6urqPv/88++++67hwT/++KN9+/b3N66oqGjufSQj3En11626UV5GGyCV9kk7lqnqaUfLy428Z4+0TxpHuDhpvRxVRzK0QZaS3ZPJuCftj1vKAa5QXi7uZYutrKwUCsWD2zRbJGljYzNr1qxZs2a14Re7uLg0rMEvKSlxdXVt2MDZ2ZkQUlpaakh1JSUlDUdc6+rqduzYsXPnTsNLMzMzOzu7+rzbqHEj5ubmy5YtmzNnTgtDtbW1bWFLkRrhw8YX0ee7PORz0CoSPmlni3WzuilsbY1/CSzhk8Ydo5+0/p7s73epra2UF88z4kk7lq97ozsn3UFoOJlZHxQUlJKSYrg2YVn2woULXbt2bdjAzMzM39///Pnzhpfnzp1r2GDXrl02NjaDBg2qPxIcHHzu3LkmG6MHG+pJ/rqShTuxtQSW0kseltW3XG5R2enb+f3UsugOnCTCLl26REREvPXWW5mZmcuXL7e3tx84cCAA7NmzZ+rUqYY2M2bMWL58+Y0bNw4ePPjLL7/UHweAzZs3T5kyhWH+jm3mzJlr1qy5fPlyQkLC119/PWPGDC7ClqTE3dvurn+pQ0if8vJyvmMROiyllzwsq2+hGzdudI4YUPfpmKRT8XzHYgpcrR+3ffv2+fPnDxw4MCgoaO/evYQQAFAqlebm96ZtzJs3r7i4eOzYsTY2Nlu2bAkODjYcLysrYxjmxRdfbPhukydPzsnJef75583MzNasWdOvXz+OwpaeK8k3aNCQ8qv7SkpKcHTuwU7m0T4usrj+lTMsq2+JzMzMaqdO1MY1Ne1WVFQU3+FwjkisOH3u3Ln+/v4tfEZYXl4u+dxQWlo6eflXteqgg28+YZQ3lPBJe/KA7pVAZqyv8W8JJXzSuMPRSfv4EltYQ9dGGPOpuXAY66RRSv3nfzvctfrTN6fJYco9DgNJnJ2d3UdvL0z3Hc53IEJHAU5hKb0M4DYULVGjJ5oek9cselUOWRAwEcpBV0dSVENzqviOQ9iwlF4msKy+JU7n066OxErSWy81hIlQ+ghApBs5kYtrrT0IltLLRH1ZPd+BCNpxDe0npy05MRHKQpQbE4+zxh8IN52QD1x9+6Hic9koeRROGGAilIUoNYnXYM9/kMQ82gcToTz0wceED6SncCqPPuYmo+wgoz9VzsJdyI1SWq7lOw6hMpTShzhhIpSFSCyrf6CkIupuRVws+I7DhDARyoIZA6FOeBXcLCyllxUsq3+w4xoqkwVl6mHXl4t+ahKvwfkyTcMHhHJjKKvnOwqBitdQWT0gBEyE8hGlZvAxYXNO4q70MhPphpv0NutErrymjAImQvl4zI2cKaB1eE94HyyllyEsq29OWhklBHxlsONEQ5gI5aKdCjq1I+cLsPM3hqX0MoRl9c2R4QNCwEQoK/3UBKsJ74el9DKEZfXNic+lUTIbFwVMhLIS5YbVhE3AmTLyhGX1TYrHO0Ikbf3dmeMalsW+/09YSi9PWFZ/v/wayK2mwQ6y6w6YCGVEbQn2ZiS5FDv/37CUXrawrP5+xzXsY25EIb/egIlQXvrhWmv/hKX0soVl9feL19AotRw7gxz/ZjnDRUcbwQeEcoZl9Y3Ey6+C0AATobxEuZHjOBzUAJbSyxmW1TdUqYPrJbSXixy7AyZCeQm0J1U6moHDQQCApfSyh2X1DSXm0h5OxELBdxx8wEQoLwSgrxuutXYPltLLHJbVNxSfy8pzXBQwEcpQlBuW1d+DpfQyh2X1DR2X60wZwEQoQzhfph7OlEFYVm+gY+FsvnwLajERyk5PJ5JeTotq+Y5DALCUHmFZvcH5QuprSxzN+Y6DJ5gIZUfJQLgLdn4spUcAWFb/P/Jca7seJkI5isJNerGUHgEAltX/T7xGjmtt18OvATnqp2Zwvgw+IEQGWFZPARLy2L54R4hkpY8ruVhIa+Q9axxL6ZEBltWnlFBrJfG2lm93wEQoR9ZK6GJPzuTLt/NjKT2qh2X1x+U9LgqYCGVL5pv0Yik9qodl9fG5NErG46KAiVC2otxkPV8GS+lRPSyrl+dmvA1hIpSpfmomIZfq5dr3E3OxghD9Tc5l9TlVUKalne1l3R0wEcqUswW4WZIrxTLt/HhHiBqSc1n9MQ3b142ReWfARChfst2kF0vpUSNyLquP18j9ASFgIpSzKLnOl8FSetSInMvqZbsZb0P4ZSBfUW7kuCzvCLGUHt1PnmX1ZVq4WUZDneXeHTARypdfO0Ip3C6XXefHUnp0P3mW1Z/Q0HAcHcFEKHNRatndFGIpPWqSPMvq43NZmRdOGGAilDUZbtKLpfSoSfIsq5fzZrwN4SmQNRlu0ouFE6hJMiyrr9XDhUIsqAXARChz3R1JThXNq+Y7DhPCUnrUHLmV1Z8toIF2xFbFdxwCgIlQ1hQE+riShDwZrbWGd4SoOXIrq5f5ZrwNYSKUuyg1I5/RUSylRw8gt7L6eA0r800n6im5e+t9+/bt27fP3t7+1Vdf9fLyur9BSkrKN998U1NTM378+P79+9cfv3HjxrfffltUVBQQEDBr1ixLS8ubN2/++uuv9Q3Gjx/fsWNH7iKXlSg3sviMXO4IsZQePUB9Wb2XDHbmowCJefTrftgZALi7I9y+ffvUqVNDQ0MrKyv79OlTVlbWqEFmZmZkZKSlpWWXLl1Gjx4dGxtrOH78+PGIiAidTterV6/09PTy8nIAuH79+oYNG4r/R6vVchS2DEW4kqvFtFLHdxwmgaX06MHkU1Z/pYg6WRB3K77jEAau7ghXrVr1ySefPP/88wBw6dKlH3/8cdasWQ0bbNq0adiwYcuXLweA2traNWvWDB48mFI6Y8aM1atXT5s2rdEb+vr6rly5kqNo5cxCAd0dyak8Gu0h/QxxMo99JRAvgVGzDGX1z3TgOw7uHdfgymp/4+RLoaysLCkpKTo62vAyOjo6Pj6+UZv4+PhBgwY1apCRkZGSkhIVFfXJJ5+sX78+Jyenvr1Go3n33XfXr19/+/ZtLmKWM5ls0oul9Oih5FNWj5vxNsTJHaFGowEAZ2dnw0sXF5fDhw/f38bFxaW+QWVlZVlZWXp6uoWFxeTJkydNmpSSktK9e/fz5897e3u3a9du4MCBVlZW58+fX7p06d69e+uTaCN1dXWbN28+cuRIw4Nr1qxRq9X3N66urlYoFI/2t0pBuD3zfymkKrBFtcTiPWmpZWCtVNlDdVWVqX+1eE8aj3g5acHWcKlQVVxeZS7O/10tP2nHc1RLgrSm7wumZ2FhwTAPueXjJBGam5sDgFarVSqVAFBXV2dpaXl/m7q6OsO/6+rqCCFmZmZKpbKqqmrVqlVDhgwBgJycnI0bN3700UcDBgwYMGCAoXHHjh3ffffd48ePN/mrFQpFz549H3/88YYHnZ2dLSws7m+s1WqbPC43A71gSoJeaaZStmCAQLwn7UImjXSlvAQv3pPGI15OmgVAgB2bXGku0pGDFp60OxVUS9lgF1l8Jh+aBYGjRKhWqxUKRWZmpr+/PwBkZmZ6eHg0auPp6ZmZmWn4d0ZGhiFXGSaXGn4KAAICArKyshr9YHh4+NatW5v71QqFIjQ0dMKECS2Jk2GYlpwjyXOyhPY2bFIJ6dWCRejFe9JO5ukj3fgJXrwnjUd8nbRIN3qqgESKc+GxFp60E3lsPzV+Jv/GyYkwNzcfMWLE9u3bAaC6unr37t1jxowBgIqKit9//12n0wHAmDFjfvnlF71eDwA7duwwNGjfvn14eLjhbo9l2fj4+G7dugFA/aRTSul///vfkJAQLsKWMzls0oul9Kgl5FBWH6+hWEHYEFezRpcvX/74449fvHgxLS2tc+fOI0aMAICMjIynnnqquLjY3t5+8uTJW7dujYyMdHV1vXTp0rFjxww/uHbt2vHjx+/bty8tLc3MzOzVV18FgGnTpt28edPHxyc5OVmv1//xxx8chS1bUWqyK52+3pXvODiDpfSohSJdydtSr6yNz6XTO+Pt4N+4SoShoaHJycknTpxwcnKKiIgw3IP7+flduXKlXbt2AGBpaXns2LETJ05UVVX169fPxsbG8IP9+/e/du1aQkKCWq0OCwsz/OB333134cKF3NxcDw+P0NBQlQpXxzOy/mryeqKegkKqiQJL6VELSb6svqgWMipod0dp/nVtw+HKMg4ODk8++WTDI2ZmZsHBwX//bqWyfgpMQ05OTqNGjWp4xMLCIjIykqM4EQB4WRNLJUktpQF20uweWEqPWs5QVv9MB2l+YOI1bB9X0pKZcfKBJwPdI+1qQtyVHrWctHerj8/FPQgbw9OB7olyk+xu9VhKj1pF2mX1uOnE/TARonskvEkv7kqPWkXCu9VX6+ByEe3tgn3hHzARonuCHUhxLc2ukmAuxMIJ1CoS3q3+VD7t5kisOJwcIkqYCNE9BOAxN+aEFB8T4q70qLWkWk2I46JNwkSI/ibV0VG8I0StFekmzUSIm/E2CRMh+luUmwQnjmIpPWoDSe5Wr6dwKp9GuuHXfmN4RtDfermQ1FJaWsd3HEaFpfSoDQxl9VmVksqFlwqppxWRx1LbrYNfD+hvZgz0dJZaBdVJHBdFbSK93erxAWFzMBGif+inJvG5klpoMTEXS+lRW0jvMSFuxtscTIToH6LcGCnNl8FSetRm0iurP5HL4h1hkzARon+IdCNnC2idVO4JsZQetZnEyupvllGGEB8b7AtNwESI/qGdCvzbkXMFErkQxsIJ1GYSK6s/rqH98XawGZgIUWNS2qQXS+nRo5BSWT1uxvsAmAhRY1ES2oYC7wjRo5DSfJn4XJwy2ixMhKixfmomXsOy4u/+WEqPHpFkyurzayC/hgY7YF9oGiZC1JjaEhzMyfUS0fd/LKVHj0gyZfXHctjHXAmDebAZ+CWBmiCNTXqxlB49OmmU1eNmvA+GpwY1IcpNCvNlsJQePTppPCaMxzVlHggTIWpClFr0u9VjKT0yCgmU1VdoIbmUhjljX2gWJkLUhAA7Uq2ndytE3P+xlB4ZhQTK6hPzaKgTsVDwHYeAYSJETSCGtdbE/JgQCyeQUUigrD5egyurPQQmQtQ0sW/Si6X0yFjEXlZ/XEOjcA/CB8Kzg5om9k168Y4QGYuo58toWThbQCNxTZkHwkSImhbqRO6U06JavuNoEyylR0Yk6rL68wW0oy2xN+M7DmHDRIiapmSgtytJEGf/x1J6ZESiLqs/jiurtQB+VaBmRbkxIt2kF0vpkXGJt6w+XoOb8T4cJkLULPFuQ4Gl9Mi4RPqYkAIk5LJ98QHhw2AiRM3q40ouFdFqHd9xtBKW0iOjE2lZfXIJtVURL2vsCw+BiRA1y0oJQfbkjNg26cVSemR0Ii2rP47joi2DiRA9iBhHR7FwAhmdSMvqcTPeFsJEiB4kSk1EN18GS+kRF8RYVo+b8bYQJkL0IFFuTEIu1Yuq++MdIeKC6ObLZFXSci0NtMe+8HCYCNGDOFuAuxW5XCSa/o+l9Igjka5EXBUUhpXVsCe0BCZC9BDi2qT3DJbSI250akdq9WIqq4/PxZkyLYVfGOghxLVJL46LIu6Iq6weN+NtOUyE6CGi1OSYRjTzZbCUHnFHRI8JS+sgrZz2wGcELYOJED1ER1uiIORWuQj6P5bSI06JqKz+RC7t7YLPCFoKzxN6uL5u5LgYRkexlB5xSkRl9bgZb6tgIkQPJ5ZNehPzsIIQcchaCf4iKas/noub8bYCnin0cGLZpDcxF2fKIG6JYnS0Vg8XC/EZQStwmAirq6tv3LhRUVHRXANK6c2bN/Py8u7/T3fu3ElLS2PZf8zRuH37dnZ2tvEDRQ/T3ZFoqmheNd9xPAxOGUVcE8V8mTP5tLMdsVXxHYd4cJUIY2JifHx8JkyY4OPj89NPP93fIDc3Nyws7IknnujWrdv06dMpvffZunnzZs+ePXv37j1ixIhevXoZDpaVlfXv33/w4MG9evUaP368VqvlKGzUJIZApBs5Iey11rCUHpmAKMrqcTPe1uIkEer1+mnTpm3YsOHixYu//fbbzJkzKysrG7VZsWJF586dU1JSUlJSDh069NdffwGATqd76qmnxo4dm5ube+PGjZiYGEPjdevWWVlZ3bx5My0tLTU1dfv27VyEjR4gyo0R8uhodXX1c7MX2f6xlK2r4TsWJGXtavLztsx7/b1V9dfuAhSvYbGUvlU4SYSJiYmVlZXPPPMMAPTv39/T0/PPP/9s1Gb79u2zZs0CAHt7+4kTJ+7YsQMAYmNjKyoq/v3vf1dVVQGAk5NTfePp06czDGNpafniiy8aGiNT6qcW9MTR3bt3H8iCgkrt77//zncsSMrW/d9Xta5dNh84c+HCBb5jaRpLITGP9sWZMq3Bycm6c+dOx44dFQqF4aWfn196enrDBuXl5UVFRZ06dTK87NSp0507dwAgOTnZy8tr+PDhQUFBrq6u3377raFBRkbG/Y2bxLJsRkbGuX+qra018l8oP71dyPUSWiHUMemQnuH0+mGXjPiwsDC+Y0FS9nj0APsz3+pKNH5+fnzH0rQrxdTFgrhZ8h2HqCi5eNOKigoLC4v6l1ZWVuXl5Q0bGF7Wt7GysiorKwOAwsLCkydP/vXXX8OHDz99+vTAgQOjoqJ8fHxqamrub9wknU63Y8eOQ4cONTy4fft2Ly+v+xtXVlYSggMILdXVTnXkblUvayGetJhqz2GfH9sZpWUY5gHzs/iCn7Q2EOZJC+sZmnbq0GMHrRIKdf0UQvykHcpR9nEiFRX4jOAeKysrhnnILR8nidDV1bW4uLj+ZVFRkVqtbtjAxcWFYZji4mIHB4eGDdRqtbe39/DhwwGgd+/eQUFBCQkJnTp1cnBwqH/D+9+tITMzs0WLFs2ZM6clcVJKbWxsWv/3ydRAT/25UtUAV73QTpqewhepuq39Fe3aCfQyGD9pbSDkk7YwhP3sBjOiIyffn4+CUnqmxGy4F7GxsXh4a/Q/nAyNhoSEpKSklJSUAIBOpzt79mxoaGjDBiqVqmvXridPnjS8PHnypKFBz549Kyoq9Pp7KzeUlJTY2toCQGhoaH3jxMTERu+GTCNKzQhzk95fb7OO5tAXd+JGpvKCP3O1GC4IsrL+BK613XqcXNF07Nhx+PDh06dPX7hw4datW/38/CIjIwFgx44dO3fu3LNnDwDMmTPnnXfe8fT0TE9P37Nnj+HJc58+fQIDAxcsWDBlypRdu3bV1tYOGTLE0HjmzJlBQUGlpaVbtmyJi4vjImz0YH3dyHNxVCu8VLg6iV3WE6cGINNRMTAnmFmbbxVQzAAAF9BJREFUxP44SMF3LP9wt5JoWerXDhNh63B1a79t27Zly5a9+eabgYGBe/fuNRx0d3evv5mbOnWqVqt9//337e3t//jjj44dOxqO79mz57333ps3b17nzp3j4+MNd4SjR48uLy9fs2aNmZnZjh07evbsyVHY6AHszcDXllwuYQbY8R1KA4eyaKUOnvDGRIhMamYXpsNObVoZI6isk1DA9FNjX2g1IuRqmDaYO3euv79/C58RlpeXGxItaqHZCXpv87q3wgT0KG7oX7rJnZgX/QXd+fGT1gbCP2lLzugrtLDhMQHdFE45XNPD1WxusKC7gwDh+UKt0Mu2cs9fh0pLS/kO5J5LRTS5BJ7zw48x4sHrXRXb09h8wUzP1Gg0h2IP93EUapGTgOE3CGqFz+ZOPHUyMerx0XwHcs/KS+z8rgxuuoZ44WYJz3RgvrgmiG2ZdDpdaP/hmSf3bXh/Ed+xiA9+haBWUFCWmFsV1QhiwsytchqbxU7rjJ9hxJuF3ZiN19lKHd9xAABAhZZlzK30OmFEIyr4JYJa4dDune8M6QCv/lIlgL62Nomd0YXBJfYRj/ztSH81szmF/0vD3FqFat6edc/02PL5J3zHIj6YCFEr2NvbvzH5qSg/5/9c5rnn51XDzlvs7CABzVNA8rQ4hFmbxPJeVvTWGXZOP5+XnhnVcFUv1EKYCFGrrQxnPr2iz6niM4bPruqf82PUApq+imSqlzPxt4OfbvGZCS8U0thsdmF3vC5sI0yEqNU62JKXA5hl53mbI1Cpg69T2Pld8dOLBGFxiGJ1EstjIdrCU/oPwhT4mKDN8KsEtcXboYq9d9jLRfz0/U3X2WgPppOQCpmRnA3zJCoG/srgpzv8ls4W1MBLAfhl3nZ47lBb2JvB0lDF6yd5uCnUsrD+KvtmN/zoIgF5sxuz6hI/3eGtM+yaCIUCLwsfAX6boDaa2ZnRVMGBTFNfBf94k+1sBz2dsd8jARnfkcmugoRcU3eH/7vGdmoHwzyxOzwSTISojZQMfBzOzD+p15lwlgAFWHuZXRyCkwKQsCgIzO/KrEky6ZSZkjr4+JJ+VW/sDo8KEyFqu6d8GE9r2HrDdJ3/9zuspQKiPfD6FwnOlEDmVD57rcR0N4UfXNCP9WW6OmB3eFSYCNEjWdNb8e45fZmpVjdcncQuDsEPLRIiCwW82kXxialuCm+X022p7LKeeDtoBPidgh5JDycyzItZk2SKaQLxGppTBWN98UOLBGp2ELP7DptRaYqbwkWn2QXdFFhKaxT4nYIe1Ue9mI3X2LsVnHf+VUn6xSEMzo5DguVgDi/6M59d4fym8GQePZVH5+F2S0aC5xE9Kk9rMiuIeecct53/egk9m0//1Qk/sUjQFnRjtt5gS+o4/BUUYOEp/cfhjBVXG6vLDn6tICNY3F1xKIueLeDwpnDlJXZeV4Ul9nwkbF7W5Mn2zMbrHF4X7kxjq3S4Dacx4alERmCjgvd6MgtPcfWkMLOS7rvLzuyCH1ckAm+FMJ9d0Vdzs0NLHQvvnGM/jVQw+IzAePCbBRnHlECmqBb23OHkQviTy+yUQMbejIv3RsjIOtuTcBdm201O+sK6y2x3R9JfjWnQmDARIuNQEPi0j+KNU2ydsbt/US1sS2Xn4LwAJB5v92BWX2L1xn5WkF8Dn1zWrwzHvmBkeEKR0UR7EP92sMnYT0e+uMaO9WW8rfESGIlGhCtxt4Jd6UbuC++d10/qxATYYV8wMkyEyJj+00ex4qK+qNZob1ijh43X9QtwiW0kNotDmI8uGnNvppRS+sstdmkoVtAbH36/IGPqYk/G+DAfXzTarJnNKWwfVybIHi+Bkcg82Z7RsRCXbbRUuPCUfkkPhZO5sd4P/Q0TITKyD8IU36ayN8uM0P/1FNZdYd/sjp9SJD4EYGF3o+3NdCSHXiuGV3HiNDfwtCIjc7WE17sq3j5rhKcjP99ivawh0hVvB5EoPe/H3CiFc49cX8tSWHhKv7o3Y47DotzARIiMb0FX5lQejdc8av9fe5ld3B27PhIrFQPzjLE307ZU1oyBpzvg1zVX8Mwi47NUwge9mIWnH2n2+IFMqmXhcW+8HUQiNr0zczjnkZ4UVOtg2Xl2bQQussshTISIE5M7MXoWfr7V9mvhVZf0i7vj6hlI3KyVML0z85/Lbe8Iq5PYvm7kMTfsChzCRIg4QQDWRijeOsPWtGmuwJl8erMMJnTEzycSvTlBip23WE11W342txq+uKb/qBd2BG7h+UVcGeBOQhzJhqttuRZelcQu7M6o8OOJxM/VEiZ2ZD6/2pZLwn+f0U8JZHxt8XaQW/hNgzj0SQSz9rK+oKZ1P5VWRo/lsFMC8MOJJOLN7sxXyWy5tnU/lVRE92WwS0Jwvhjn8LsGccivHZnYkXn/fOuuhVclsa8GMTYqjoJCyNQ62JJoD+br5NaNjrx5Sv9uqMIO15rnHiZCxK33eip+vs1eK2nprLncavj1Njs7CK+CkaQs6cGsu9KKJen/zKB3K2FaZ/yKNgU8y4hbDuawqLtiyZmWfgGsu6Kf1IlxseA0KIRMLcSRdLGH7S3bm0lPYfFp/doIBT4mNw08zYhzc4KZa8X0UNbDbwrLtLA5hV3QFT+WSIIWhyhWJbEtWYf762TWxQKewCJaU8FvHMQ5MwY+CmfePK1/6FfAl9fZ4V44Rw5J02APYmcGf9x9yE1huRY+uMCuicCnA6aDiRCZwvgOjI0KfnjguFCtHj67yr6BOy4h6XqjG/PRpYckwo8v6od7kTBnvBw0HfzSQSaytrfi7bNspa7ZBj/cZLs7QqgT9n8kWU/7MkW1cCK32bGRzEr6VTL7fhh+M5sUnm5kIhGupK8bWdfMWlMU4D+4xDaSOgWBN7oxq5q/KXzrDPtaMONtjZeDJoWJEJnOynBm/VV9TlUT/2l3OmutggHu2P+RxL3kz5wroFeKm7gpvFBI47LZN7rh5aCpYSJEpuNrS14OYN4910R9/ZokdkkIfhqR9Jkr4LVgZm1TezPNS9Sv6KWwxaUkTA6/epBJvROq+DODXi76x+Xw0RxaVAujffDTiGRhdhCzL4O9W/GPXrArnS2tgxf9sRfwQMnR+1JKN2/evHfvXmdn5wULFnTt2vX+NomJiZ9//nl1dfXEiRMnTJhgOPjee+/l5OQY/h0QEPDGG28AQFJS0hdffFH/g3PmzGnyDZHw2apgSQ/m9ZP62JF/f/ZWJenf7M7glktIJtqp4CV/5tMr7H/63BsF1bLw1hn288cUuOsgL7i6+ti4cePKlStnzpwZGBg4YMCAgoKCRg1u3rw5fPjwqKioF198ce7cub/99pvh+M8//2xvbx8WFhYWFhYQEGA4eOfOndjY2LD/sbOz4yhsZAIzOzO51bA/897l8OUierGQTuqEF8JIRl7vynybyhbW3nv5xTU2oB0M88Q0yA+u7gjXrVu3Zs2akSNHjhw58ujRo999953h3q7exo0bx40bN2vWLADIyclZt27d2LFjDf/pqaee6tu3b6M3VKvV06dP5yhaZEpKBj4OZxac1A95WqlkYFUSO7+rwgLnByA58bQmY32YjdfYpaFMcS2svPSPMRJkYpxchpeUlNy8eTMqKsrwsm/fvmfOnGnU5uzZs/UNoqKizp49W/+fNmzYMHPmzI0bN9bV1dUfvHv37vTp0xcvXnz/WyHRGdWe8bSGLTfYjEq6P4OdjisLI/lZFMJ8dlVfqYMPLuif9mWCHfB2kDecXINoNBpCiKOjo+Glk5OTRqNp1CY3N7dhg+rq6pKSEnt7+zFjxvj6+lJKv/nmm507d8bFxSkUCldX15deesnPzy85OXnQoEHff/99/e1jI3V1dRs2bKgfaDXYuHGjh4fH/Y0rKioe9U+VH2OdtGXBZNy+0l3a1H+FhjO1NeW1D/8R8cJPWhtI/qR5MBBuB//aePioRdiZMebl5S3doeUBJH/S2sDKykqheMiIEyeJ0NramlJaW1trZWUFANXV1TY2NvcHV1Nzb8PW6upqQoih8UcffWQ4+Nxzz/n4+Jw4caJ///4RERERERGG487OzitXrmwuESqVymHDho0ZM6bhQW9vb8Ob38/W1raNf6SMGeWkhTKV5atHHPAI8cs5Ztvv3Ud/Q4HDT1obSP6kFX096WSVk3POso7TEo31npI/aVzgJBGq1WozM7Pbt28HBwcDQHp6ure3d6M23t7e6enphn+np6cbfqRhg3bt2rm7u99/K9m5c+e8vLzmfjXDMP7+/kOGDHn0vwJxSqfTWasYnYOHvqaS71gQ4odSW6VwCrHMSeA7ELnj5NmMSqUaM2bM5s2bAaCkpGTXrl2G6oiSkpJNmzYZnvxNmDBh+/bt1dXVALBlyxZDg7KysuLiYsObxMbG3rp1q1evXgBw9+5dSikA1NbWbt68uU+fPlyEjUzJzs7u6O7t300KXffRe3zHghA/9m7fsuUJ9xP7f3t4U8QlruYprVixYvjw4SdOnMjMzBw5cuTAgQMBICcnZ+bMmc8++6yZmdmECRN++eWX4OBge3t7rVZ76NAhAEhPT+/bt6+fnx+l9O7du1988UXHjh0BYNmyZQcPHvT29r5161ZAQMDPP//MUdjIlIKCgoKCgviOAiHeODg4vPDCv/iOAnGWCP39/VNSUq5evWpvb+/r62s4GBgYWFRUZKgCVKlUu3fvTk1NrampCQ4OZhgGALp3756Tk5OWlqZQKPz8/CwtLQ0/uHXr1vT09MLCQjc3Ny8vL45iRgghJEMcVq6oVKoePXo0PMIwjIODQ8Mj/v7+jX7KxsYmJCTk/nfz9fWtT6gIIYSQsci6fuuvv/6qn7mKWmjPnj0s+5CdRVFDOp3u999/5zsKkamurt6/fz/fUYhMaWlpXFwc31GIkqwT4dy5c0tKSviOQmSmTp2q1Wr5jkJMamtrp02bxncUIlNYWDh//ny+oxCZO3fuLF26lO8oREnWiRAhhBDCRIgQQkjWpLbMq4uLy9dff713796WNNZqtRMnTlSpcB/MVrCwsBg5cqRhli9qCZZlVSrV0KFD+Q5ETOrq6mpqavCktUpVVVVeXh6etEa+/PJLPz+/B7chhkJ1ydBqtUePHuU7CvT/7d1ZTBNfGwbwIoR/taSr2iCWRdwVI6QdIkVwa4iiaDAoCTEKXqgxWMCgBiSKS6IBd4PojZFERQW9Mhi8UUBkCaYUIUhBXHArYhfagrRwvotJiGlB8920Q+f5XZ0zvBdPTqbzkjOTGQAARoiOjv7na+e8rRECAAD8X7DBBQAArIZGCAAArIZGCAAArIZGCAAArMaWRnjjxg2lUhkfH19RUTFhwevXrzdu3EhR1PHjx+kPRUFFRUV8fHxMTExJSYnrXzs7O9VqdWxsbFxc3JkzZ/CyOg6H8/v377y8PIqiNm3a1NTU9JfKgoKCvXv3ui0Ykw0PDx85ckShUGzevLmlpWXCmt7e3oyMDLlcrlKpnj9/7uaEDGS1WnNychQKRVJSklardS0YHh4+ceJEbGxsbGzsqVOncE37B8ICDx8+nDt37qtXr6qqqkQiUUNDg1NBf3+/UCi8detWa2trTExMfn6+R3IySmNjo1AofPr06atXr2QyWXl5uVPB7du3T548WVtb++LFi5UrVx44cMAjORnl6NGjSqWytbW1tLRUJBINDAxMWPbkyZP58+dLJBI3x2OmrKysNWvWaLXaa9euSSQSk8nkVPD9+/c5c+bk5eU1NzdXV1fX1NR4JCej7Nu3LyEhQavVFhcXS6VSm83mVHDs2DGFQqHRaN68eRMZGXnixAlPxJwyWNEI4+Lirl+/To9zc3N3797tVHDhwgWVSkWP6+rqZs6cabfb3ZmQgfbs2ZOTk0OPS0pKVq9e/Zfix48fh4WFuSUXc/3+/VssFtfX19PTdevWXblyxbXMYDBEREQ8evQIjZAQYrPZ+Hx+S0sLPVUqlaWlpU41OTk5O3bscHs05jKZTDNmzOjo6KCnUVFRd+7ccapRqVTjp19RUdHmzZvdGnGqYcXWaFtbm0KhoMcURbnuJGi1Woqi6LFCoRgYGPj69atbIzKPVqsdX7To6OgJt1/GaTSahQsXuiUXc3358sVgMMjlcnpKUVRra6trmVqtzs7Olkql7k3HUB8/frTZbJGRkfR0wkVrbGykKEqtVm/duvXixYsOh8PtMZmF/mLrkiVL6OmE17Tk5OQHDx50dXV1dnY+fPhw+/btbo85lXh/IxwZGTEYDEKhkJ6KRKIfP3441ej1+vECf39/Ho/nWsM29HYxPRaJRCaTabK7gM3NzZcvXz5//rwb0zGRXq/n8Xjjb+wTiUR6vd6ppqqq6vPnz3v27HF3OKbS6/UCgcDHx4eeTrhonz59KioqoigqOzv7/v37hw8fdntMZvnzesWZZNHS0tJ4PN6qVatiYmJmzpyZkpLi3oxTjPc3Qn9//xkzZthsNnpqsVj+PIdoAoHAarXS47GxMZvN5lrDNnw+/89F43K5//33n2vZ27dvk5KS7ty5M+HnlFlFIBAMDQ2Nf6zR9UyzWq2ZmZnnz583Go2Dg4OEEIPBMDo66omwTPHnT48zyc+Tz+enpKSkpaWtWbOmuLi4rKzMvRkZRyAQjP82OZMs2q5duxYsWKDX6/v7+4OCgjIyMtybcYrx/kbI4XBCQkJ0Oh091ul0rl+6Dw0N7e7upsc9PT1+fn5BQUHuTMhAoaGhTos2/m/7uM7OzoSEhEuXLm3bts3tARknKCjIx8ent7eXnnZ3dzudaXq9/tevXwkJCeHh4ampqQaDITw8vKenxwNZGUMmk9nt9r6+PnrqumgcDicsLEwsFtNjiURitVpZvjsaHBxsNBp//vxJTydctJcvX+7cudPX19fX1zc1NbWmpsbdKacWT9+kdIczZ86sXbvWbrcPDg4uXbr07t27hJDR0dGCgoIvX74QQtra2gQCwfv37wkhmZmZuDNPCLl3797ixYvNZrPD4Vi/fv3p06fp48XFxW/fviWEdHV1yWSysrIyj8ZkluTkZLVaTQjp6enh8/n04wx9fX0FBQVOlTU1NXhYhpaYmHjkyBFCyLt37wICArq7uwkhHz58KCwspAsePXoUERFhNpsJIbm5uWvXrvVgWoZYt24d/SBoW1sbj8fr6+sjhOh0urNnz9IFq1atys3NHRsbGxsbU6vVWLS/Y0UjtFgsiYmJs2fPFovF6enpDoeDEGK327lc7vjjaufOnRMKhTKZLDIy8uPHjx7NywgOhyMjI0MsFkul0sTERIvFQh+fP39+ZWUlIeT06dOiPwQGBno0LyP09vauWLEiODhYKBQWFRXRB5uamrhcrlNlfX19eHi42wMykU6nW7ZsWUhIiFAoHH/Qsa6ujs/n0+PR0dGDBw/OmjVr3rx5crlcp9N5LixTdHR0LFq0KDQ0VCQS3bx5kz5YXV0tlUrpsUajWb58eUhIiEwmW7lyZXt7u+fCTgEs+vrEwMCAn5+fQCCYrGBoaMhoNAYGBrozFcOZzWa73S6RSDwdZCr5+vWrWCzmcrmeDjKV/HPRLBbLyMjI+B4pEEK+ffsmkUgmvHlPMxqNPj4+f7noAY1FjRAAAMAVKx6WAQAAmAwaIQAAsBoaIQAAsBoaIQAAsBoaIQAAsBoaIQAAsBoaIQAAsBoaIQAAsBoaIQAAsBoaIQAAsBoaIYD32LdvX3Bw8IcPH+jpyMhIXFwcRVFDQ0MezQXAaHjXKID3GBwclMvlQqGwtrbW398/Ozv7xo0b9fX1UVFRno4GwFxohABepaWlRalUHjp0KD4+fsuWLVeuXMnMzPR0KABGQyME8DZXr17NysoKCAjYsGFDZWWlj4+PpxMBMBoaIYC3MZlMwcHBZrNZq9VGRER4Og4A0+FhGQBvs3///mnTpslksoMHDzocDk/HAWA6NEIAr3Lr1q3y8vKSkpKKioqGhobCwkJPJwJgOmyNAniP9vZ2iqLS09OvX7/O4XDOnTuXn5//7NkzlUrl6WgAzIVGCOAlrFYrRVG+vr6NjY3Tp0/ncDiEkKSkpKamJo1GExgY6OmAAAyFRggAAKyGe4QAAMBqaIQAAMBqaIQAAMBqaIQAAMBqaIQAAMBqaIQAAMBqaIQAAMBq/wMm/gRNESQfqAAAAABJRU5ErkJggg==", "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "metadata": {}, "execution_count": 6 } ], "cell_type": "code", "source": [ "rvecs = collect(r_vectors(basis))[:, 1, 1] # slice along the x axis\n", "x = [r[1] for r in rvecs] # only keep the x coordinate\n", "plot(x, scfres.ρ[1, :, 1, 1], label=\"\", xlabel=\"x\", ylabel=\"ρ\", marker=2)" ], "metadata": {}, "execution_count": 6 }, { "cell_type": "markdown", "source": [ "We can also perform various postprocessing steps:\n", "for instance compute a band structure" ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Computing bands along kpath:\n", " Γ -> X -> W -> K -> Γ -> L -> U -> W -> L -> K and U -> X\n", "\rDiagonalising Hamiltonian kblocks: 6%|▉ | ETA: 0:00:02\u001b[K\rDiagonalising Hamiltonian kblocks: 15%|██▍ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 24%|███▉ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 35%|█████▋ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 44%|███████▏ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 54%|████████▋ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 63%|██████████▏ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 72%|███████████▌ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 83%|█████████████▍ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 94%|███████████████▏| ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 100%|████████████████| Time: 0:00:01\u001b[K\n" ] }, { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=210}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOydeUBM6xvHvzNtKlkKlaIslSVljSyRKyWSLbss0cJ1s95wUddyZd8lu2vPUrZI4V6RfamUS0IKRSnt28z7+6N+baZmO+fMTM3nr86Z9zzPM9OZ88z7vs/CIoRAjhw5cuTIqauwJW2AHDly5MiRI0nkjlCOHDly5NRp5I5Qjhw5cuTUaeSOUI4cOXLk1GnkjlCOHDly5NRp5I5Qjhw5cuTUaSTpCPfv3+/r6ytBA+TIkSNHjhyWpPII//nnn7Fjx6qqqiYkJEjEADly5MiRIweSmhHm5uZ6enr6+PhIRLscOXLkyJFThmQc4bJly6ZNm9a6dWuJaJcjR44cOXLKUGRe5YMHDx4+fLh58+bQ0NCaR6ampnp6ehYUFFQ8OXz48AkTJgio68cP5OaydHUJgPfvWZmZMDcnHA5iYlhmZlJdW47D4SgoKDCvNzUVL16wjY05amqsJk2EuzY2NmnUqCZdu37dsqWZjk69nwfEx3O1tNiNGlFjqiAUFuLNG5apaen/Oj8f796xOnSg/l+fm4tPn1hGRgJJLtmPYLFYfEdmZyM5mdW2rVTfqxXhcrlstnA/rzMykJ7OatWq9D1++QIAurpUWpWUxFJW5jRrxsOwt28zPT1zv3xRu3o1T1dXWyixX76wAI6uLvvJE1b37uTuXVbnzqR+fYqMFoz0dGRlsVq2FPQOEe3BEhPDMjIiysqlh1FRLFNTIuT/uSYKChAXxzU1pWVipqCgwPeeZNoR5ufnu7q6Hj9+XJB/xsWLF79+/Tpr1qyKJ01MTAT8RxKCSZNYvXtjxQoCYPhwdlER4uK4+fmwtmYnJXEZvmWFIjc3V0NDg3m9e/eydu9mOTpyu3RRnD1buOfv0aNHEhI4iYlhISEXtLQ0zM1Jly7o0gWdOxN9faSmpg4ffmXWrL4LF7alyfif+fwZdnbs9++59eoBwLt3cHBgx8dzKf+NERWFhQvZDx5wBRlcVFTE5XJVVFT4jrxzh+Xnh+Bg2XCEhJC8vLz6Qn6vLlxg3bmDY8dK36O/P4vFwqpVVL7lBQvuZ2dfDwz0UlNTe/cOL16wnj/H8+d48YKVkxOfl+fNYrUMDrZwdZ0ulNidO1nq6sUrVii4uLD37OGeO8f680/WlStcNTUKbefD1aus0FAcPy7oxyXCg4XDgY0N+9kzrp4eAKSn45df2F++cJWUhDW2Wt6/z7ay2hwSYturVy/KhAIAuFwuh8Ph/+OMMEtkZKSKikrr1q1bt26to6OjqKjYunXrz58/8xx84MCB6dOni6xr0yZiYUEKCwkhJDWV1KtHFiwofcnamly5IrJgJsjMzJSIXisr0rcvd8uWAjc3oa/Nycnx8ztw7949LpfExZGAALJkCbG1Jc2akaZNibJyJ2CZgoIB9UbXyODB5MSJ8sOuXcmNG9Rryc8n9euTnByBBhcWFubn5wsy0seHLFsmlmFMwuVys7KyhL1q+HBy8mT5YbduJDycSquSk5NZLD3AQ13910aNSIsWxMGBrFxJAgPJhw+EEHL1avDRoyeKioqElWxrS86ezSWE7NtHhg0jXC5xdSWDBpHcXCrtr5nYWNK2rRDjRXiw3LlDunYtPzx9mgwfLqwMPrRvbwW4s9nN09LSqJXM4XAKS3xAjTA9I+zYseOXkrUP4ObNm/Pnz3/y5EnDhg0pVxQdjY0bcf8+Sn62hIVBTQ1DhpS+amuLkBAMHUq5WpknNhZubjA15Z49K/S1ampq7u4uJX+3bYu2beHkVPrS58/o0IFVWFjAYhVTZ6xAeHhg82ZMnFh6OH06Dh+GjQ3FWlRUYGqKx4/Rvz+VYp8/x+TJVAqUNgoK8O+/OHSo9DA1FfHx6NmTShUsFovFAiHF2topDx/i5wV/e/shvK7jT1QUOnbkAHB2ho8PYmPh54epUzFiBC5eRD0emwPU064d0tKQmsrjfVHF1auVHpUhIbC1pVgFi8UGFABl/kPpgelgGQUFhcb/p379+mw2u3HjxsJuKvAlNxdjx2LLFrRqVXrm6lXk5qJfv9JDOztcv06tztoAl4u0NIwZw+3Uibx8CQoza5o3x7t3t728FKOjwygTKhgODkhMxPPnpYcTJ+LaNaSnU6+od29ERFAs89kzdO1KsUyp4tYtmJtDS6v0MCQEAweCwjU3AM2aNbt379yKFTrR0cco9BZpacjLg54eAaCigtmzsWUL2GwcOQItLYwcicqxDXTBYqFrVzx+TKOKK1fKHSEhtDjCx4+v9ujR1dPzqqamJsWiBUOSCfVWVlbh4eF0SJ43Dz17lk8CAFy/jp49UbYpY2aG3FzEx9OhXIa5eRMKCujUCQ0akMaN8f49lcI1NTX/+OOPdu3aUSlUABQU4OKC/fvLzICtLU6epF6RpSXu36dSYGoqMjPLf8zVShiYbQDo1avX4sWL1Sjdu4uMhJkZygKeZs/GxYv4/BkKCjh2DA0bYsIEFBVRqLBaevbEo0d0Cf/4Ed++oUeP0sPoaKiro00birWoqak5OU0qLOxAsVyBkaQjVFNTa9myJeViAwNx6xZ27Cg/ExeHnByMHl1+hsWCra18UliVwMDyaD0zM0RHS9Qa6pg1C2fOIDOz9LBkdZRySmaEFE6jS6aDAsSWyjDBwZVmG2FhtDhCOoiORqdO5YeNG2Py5NLHTokv5HAwcSKK6d8KsLCg0RFevgx7e5St2V2/Djs7WhS1acONi6NFsiDUtlqjnz7BwwPHj6NBg/KToaFgsar+/0q2CeVUJCIC5ualf3fqhKgoiVpDHTo6+OUXHD9eemhjg2/fEBlJsZbmzaGhAQq/zLV+XTQmBhwOOnYsPXz+HI0awcBAojYJTBVHCGDhQhw4gB8/AEBJCWfPIi8PEydyHj9+WkTn3LBnTzx8SOUvsIowM2WH3BFSCJcLZ2d4eqJKCO7581BTg5FRpZM2Nrhzh6F1fFkhPr78Lu/UqfbMCAF4eGD37tKHBZuNqVNx5Aj1WqjdJnz+HF26UCZNCrlyBcOHlx/S95Clg6gomJlVOtOiBWxtcfBg6aGyMs6eRWjoWCurff36jaDPEh0dqKnh3TvqJefl4d698siynBw8eUJxOFgZLVtyv3yR2AO5VjnCdetQXIzff690ksPBgwcYNqzq4MaN0bEj7t5lzDppJycHOTnlcZ5mZrVnRgjA2hosVvm/e/p0HD9O/bfO0PCln9+2srhoMXn2rJY7QsZmG5TD5eLVq/K5bBleXti6FYWFpYeqqtDV/V5Q0OPTpwxa7aFpmzA0FD16oCyo//Zt9OgBmtKvFRXRogXFcQmCU3sc4ZMn2LEDx46hSq70o0dQUMDIkTwukW8TVuT8edSrB+3/19YwMUFiInJzJWoTpbi6ws+v9O9WrdCpEy5fpljFvn2THz1qOHKkq/iiMjPx5QtMTMSXJKWkpyMysnx6kZWF589hZSVRmwTm7Vs0a1Zp/6UEMzO0a4czZ8rPhIYeHzWquHNnGqKzKmBhgYcPqRfL8C+Vtm2p3FkQitrgCFNSUjp3tuvbd/CqVV9+Dr4JDkZBAQYM4HGhPImiIlevVtqeUVSEsTFiYyVnENVMnYqQEKSklB7SETKjrd2IzQ5r1KiFmHJycnJsbacqKIz+9i2ZEsOkkOvXMWAAVFVLD2/ehKUlmKzJIg7R0VXXRctYvBgbN5bv2Onp6R0+7H7/vkFSEo320DEjJKRSKBPojJQpwcgIb9/SKL8GaoMj/Oeff2Ji+hQVWdevf+vnVwMD0aEDeBYV6t4dX7/i40faLZQJnj6FhUWlM7Vsm7BhQ4weXZ67PWYMHjxAYiKVKp4/D50zZ5mh4R4x5dy5c+f586ZZWdbnzgVSYpgUIrvrouAVKVPG4MFQVMSNG+VnNDQwYQL27aPRnu7dER1dviRLCS9eQEUFxsalh+/fIycHpqZUqqiCkZF8RigGPXoMZrGie/Z8ZvfTz5WcHMTFVUqcqAibjUGDKt2ydZnExEqRC6h1jhDAnDnYuxccDgCoqmLMmPJQUkpQUlJatqxjQEB5qoZo9OzZi8WKNTC4YG8vO85BGDgchITA3r78zI0bsuQIo6KqdYQAFi7Exo2Vzvz2G/btQ34+XfaoqaF1a4o39a9erfRAuHYNdnb0JvPIl0bF4u+/G0+bFhARcVarrEDF/7l9GwoKVZ/vFZEnUZSQlISiokoPJtS6eBkA5uZo3hzBwaWHJauj1Mad6+jA2hqnTokl5OnTxsbGwe/f36qtrcru30eLFtDXLz188waFheggsXRqofk5ZLQi48bh7Vs8fVp+xsgInTtDhLKFgkP56ijzU3b50qjo5OZi714sWsT71QsXoKhY02+3IUMQFsZQAQhp5uRJNGhQvmFTgpkZ9cl2EsfDozxkplcvKClRHzns5lauQjT8/eHhQZE1UkmVh+z16+V1gKWfnBwkJ6Nt9T1UFBXh6YlNmyqdnDsX27bRaBW1afXfvuG//8prUhYW4t9/MWgQZfJ5YmiI5GQa5801IPOOcP9+WFmVL2RX4do1DBhQ03S+aVO0bk1LwJVsERbGI0BRVxcsFpJrV7jGuHF49qx8BYaOkBkbG+Tmil7+MTkZt25VKhBY+5D1DcL27atGp1fB1RW3blUq4mhvj+xsPHhAl1XUBo4GB2PQIJQ1ILx3D+3b46cVN4pRUEDLlpLJoJBtR1hUhK1bsXAh71c/f8b37/wfKHZ28tVRREejb18e52tTfZkSVFTg7IwDB0oPnZ0RGMhJScmhUAWLBRcX+PuLePn+/Rg3jkdofq3h6NHrHz+GlUVmFRTg3j1YW0vUJmGoIVKmDHV1uLhUKvTIYsHdHTt30mVVx45ISnoTF/dVfFFJSUm7d//dv39a2RnGfqlIKl5Gth3hyZMwMqq2acu1ayCEf8MdeTYhIfj6lXdIUe2LlwHg7o4jR0pXYBQU0goLexkbD7x+PZRCFS4uOH9elB4XXC4OHkTlXtS1isuXL7u7H8jL23779s2SM+HhMDWFhLoOiIIgjhCApyeOH0dqavmZmTNx/To+faLFqgsXLhQUzLGwsPvw4YOYovr1G/3kSfqhQ85lZ+SOUHohBJs2wcur2gGnTqFFC/7T+d69ER+PrxT8kJJVSrYWLC15vFQrHWHr1ujaFefOAcDHjx8VFJpnZY25d+8pv+uEoEkTDB4sSo+Lq1ehq1ubS4z++NGgsDBNTS29rJG9bK2Lgl+kTBna2hg5Env3lp/R0MC4ceWNUKjl8+evgFF2tub37z/EFJWeXk9Z+b2mZuk/KDkZHz9Wza2iibZtJRMvI8OO8PJlKCril194v0oIHjyAgwN/OYqKGDiwTidRnD2LZs1476TWvsDREspCZrp06bJ69dBmzTJNTNypVeHmhj3C5xP6+8PNjVpDpIicHKxf39/be8+jR4d7/n8lh+40bcp5+VKgGSGAxYuxezfy8srPzJ2LvXtpqag5e7bL/v19e/RYGhhozn909ezahRYtrl26NOratdLUopAQDBrEZ0+UKiSWSihk43tGOXDgwPTp06t7tU8fEhBQ7bVRUURZmTx4IJCiffvI5MnC20cnmZmZjOmysCA2NuWHxcXFOTk5JX/n5hI1NVJURI0iJt9UzRQXE0ND8vRp6eHly6RzZ8LlUqmCyyXGxuTePR4vFRYW5ufn/3w+IYFoaZH/f/ayB5fLzcrKqmGAs3PVL1pSEmnShHA49BpGqLv3EhOJjo4QYu3svs6YceLLly9lZwYNIsePU2ILD75+Jfr6JCys2gE1G/zgAdHWJm/fVjo5YQI5eJAi+6qnxLB370jLllSK5XA4hYWFfIfJ6owwPBzJyRg1qtoB58+DxUL37gJJGzIEISHgcqmyTsZ4/RoDB/J+SVUV+vp484ZZg+hHQQEzZ5YX+xg2DMrKCAqiUgWLBVdX4UJm9u3DlCkyU2ZMWLZvR2Rk1Q8kJAQ2NuXt7qSfGoqr8SQubsLhwykDBowtOzN3Lo0hM02b4tgxzJhRaW9SQL5+hZMT/Pwq9d3lchEWxj/SgipatsTXr5Xm0MwgOzdgZdavh5dXTbP18+fRo4eg03l9fTRrhmfPqLJOligsRGYmxo6tdkBtXR11dcW5c6Xd4wCsWIFVqyhOrp82DZcv4/t3gQYXF+PwYcycSaUB0sP9+1i3DhcuVHXzsrhBKOC6aAlNm6oqKn7OylIpOzNsGL59o7GV7oABmDgRzs7C3cwcDiZPxrRpVfsTPHkCHR20ELd6rqAoKMDQkJaWUjUjk44wOhpPn2LKlGoHFBbi9WuMHy+ETDs7XLsmvmmyx5UrUFJCDQVMamW8DICmTWFjg2PHSg/pmBRqacHeHn//LdDgwEAYG/Po7FMLSEnBuHE4fLjqbcbh4OZNDB4sIbNEQsCQ0TJu3Qo4enRIYeHFshgQNhuzZ2PXLjqsK2X1amRkCKdi6VJwufD2rnqe+V8qEomXkUlHuH495s9HvXrVDggPByG8Wy9VR52ttXbpUnmlK57UVkcIwMMDe/aU/3CmY1Lo5gZ/f4Fk1tYwmeJijBsHFxcetWMePUKLFtDVlYRZoiJgyGgZqqqqEyYMXLJEbdas8tvAxQVXrtBYqkJREcePY80aPH8u0PiLFxEQgNOneSyhMe8IJRIvI3uOMDERISF8HhnHjkFLC82bCyHWygovXwq6ilWbePSIz05qbV0aBWBlBUVF3LlTekjHpLBfPygoIDycz7D4eERHC/fTTVbw8oKyMpYv5/GSzK2LFhUhPh7t2wt94bx5yM7G0aOlh40aYexYuvIoSmjdGrt2Ydw4ZGXxGRkXV7pN0KRJ1ZcyM/HyZXmhNWaQO0KB2LABs2aVN03mSWhotWkV1aGigr59cfOmOKbJJB8+VK21XYXWrZGejgx6O2xLjCp1QemYFM6axT9kxs8PM2ZARYXPMNEorLE9T25urr5+dzU14+s01JUICsKFCzh1ivduvcw5wlevYGhY01pUdSgo4OBBLFlS3g7zt9+wdy/FjZOq4OSEvn0xf35NY3JyMHIk1q7l/Wv4xg307UvXbVkdEulBIWOOMC0Np05h7tyaxqSnIyUF06YJLbwOro5++4b8fIwYUdMYFgsdOuDlS6ZsYhZnZ1y/Hnn8+DUOhwN6JoVTp+LatZoqNhQU4NgxuLhQqbSE9+/fq6gYqqp2HDPm5OLFmDoVw4ahd2907AgDAzRpgvr10aBB9KdP7fPy1jk5BVBbf/z1a7i54dw53kUt0tMRG4vevanUSDfCbhBWxMwMU6eWu6UOHWBiggsXqDKNNzt34t69mgo7eHiga9dqQ7Qk8ktFIj0oZMwRbtuGsWP5bCpcvAgWC1ZWQgsvaVhP7WxAyjlzBurqfKbXAMzMau02YUrK24ICdxeXgF27SheqKJ8UNmqE4cNrCpkJCEDXrjV1MxCZGzf+KSzsyeV6Xrz48uRJPHiAlBSoq6NTJwwfjnnzsGMHzp3r0LBhlKLici2taf37Q1MTbm4UFFrKzsaoUVi3Dt26VWcbBgxgerYhJuI4QgA+Pnj8GJculR7SmkdRgro6AgKwcCHvStbbtiEqqqblCok0iWzZEqmpjGdQUJm7SDVVEuqzs0mzZiQujs9Vv/xCTE1F1NimDYmKEvFaamEm99zBgZibVz1ZMaG+hJ07iYcHBeqkJ6G+jPfv32tqdmexnDZsKM8ZtrAgFy5QqSUigrRpU56wXyWhvk8fEhhIpboSOBxiZFSsqLhIX9/y1atXglySlUWWLyf6+oTFIvr6ZP164VLdyxLquVzi5ETc3GoaPH062b1bCOFiQsm9N2QIuXRJLLG3bxMDA1JyUXExadWKPHokvl182LqV9OhBStLKywyOiCDa2iQ+vtqrYmKIoSHttpVR8ZNs355ER1MjthYm1O/di4ED+f9wfvRI9KCDkklh3eHFC4HWpmpfD4oyDA0NHz06PX78vDdvppedpHxSaGmJBg1w+zaPl2Jj8f59pZ5EVNGjB5KTFT592piYGNGuXTtBLqlfH6tXIzERcXEYMAA+PlBWhqUlbt/G48eP/f39uYJVndi8GQkJ2L69pjFhYTK2QQjhQ0Z/ZsAADByIFSsAQEEB7u7YvZsS02rC0xO6uvDxKT9TXUJLRSTYJFIC8TLUuF16qDgjLCwkLVuW18SqjnfvCJtdtUSQ4Fy+TH75RcRrqYWZyZOiIgkJqXry5xnh9++kQQMKKpBJ4YywhMxMoqdHHj4sP0P5pHD3bjJ2bOnfFWeEc+YQb28qFZXQty9RUyNJSeLKOXaMmJsTNvsL0BEY3anTzG/fqh1cMiO8d49oa5N372oSGxVF2rYV1zahEP/eS0vj8S0QQWxaGtHVLa299/07adyYJCeLaRp/vn0j+vokNJRkZmYWFRErK/Lnn3wuGTyYBAXRblgZFT/JhQvJ+vXUiK1tM8Jjx9ChA/+q/AcPon79SiWChMLaGo8fIztbxMtli5gYcDjVFlerSOPGaNAACQn02yQhNDSwahXmzaMxp3DKFISG4suXSidzc3H6NGbMoExLCY6OePQIT59CT09cUZMn48UL/PvvZxaLCzSNjdVr2hQqKjAwgL09/vqr6lJBTEzauHE4ehStWtUkVuYKbeP/xdVq6PItIJqa2LIF7u4oKkLjxhg9urw7Jn00aYLjxzF1KmJjU3//HRoavBNaysjLw4MHGDCAdsN4wnxOvWw4QkKweXNNHZfKuHgRvXqJrkhdHRYWvJewah9nzkBTE4qKAg2uxdmEJUybBg6nPL6O8vBRDQ2MHl2eSVbCyZPo2xctW1KmBcCMGQgOxt27EGw1VCD69u0aELBqyZLGmZm/5+Tg1CkMGYKvX7F5Mzp3hoICmjZFly6f2Gy9Xr1GcDiL+a55ylziBIQvrlYD48ejVSts2AAAv/6KvXtRXEyN5Bro3x8/fvTq1Wvo1q22R4/yqe/677/o0oV/GB1NML80KhuOMDAQ6ur8f55wOHj9uqbSa4JQd5Io/vkHHToIOrgWB46WwGZj2zYsXozMzNIzlE8K3dywb1+l2u6UV5NZvBh//43gYPToQaVYAGPGjFm37i81NTU1NYwahb178eQJ0tLA5SIiAu7uKCj4F2gOzP/yJaJjRzg5YeVKnD6NFy9KeyCXkZODx48lNtsQGTFDRquwcye2bcOrVzA3h4bG2f79f31L8yToxQvk5HwGNgNv27X7MXAgZs/G7t24davqQgUk/UtF7gh5s2kTli7lPyw8HBxOTS0pBMHODsHBYkmQFWJjhXgY1eJCa2VYWmLQIKxfX3pI+aSwe3doaiI0tPTwxQt8/Uplpc1Nm7BlC/7+m7leASX07InVqxEbO1FPr1hR0fvw4dmnT8PJCWw2AgPh7AxNTbRti+HD4eWFDRviDQwsi4r6pKcnMmql2IgfKVORli2xYgXc3ZGTk5uUtC4iwtbV9adCnxTx9ClGjMDQoRg61IPNnmNq2jwysuGyZWjfHi9fYvVqmJtDUxN9+mDWLGzeDBubBbt2df/wYStN9vBFXx/p6cjJYVAlNTuS9FASLHPrFjExESiMe9IkYmBAgd6WLcnr1xTIEQe640qKiwmLRV6+5PlS1WAZQkhUFGnfXlylUhssU0ZyMmnalLx5U3pIeZ9Cf38ycmRpsMzMmeSvvyiTfPAgYbPJrl2UCRSB6voRFhWRN29IYCBZt4706fM3i+WtrLz43LlzTNom5r3H4RANDZKeTqVYDof07k38/bk9egxRVe3ZvPn+z5/FsZEHDx+SYcOIvj7ZsYPk5RFSvcGpqeTOHeLvT+bNI8rK5kBc27a9KbamRqoY1qEDiYykQKyAwTIy4AhtbQVtC9m8OZk5kwK9M2eS7dspkCMOdPuM69eJoiLvl3g6wsJCoqpKcnPFUir9jpAQsm4dcXQsP6Q2fDQri2hqkg8filJS8ikMFzx/nrDZZNUqaqSJDN/GvISQzMzMceNmT578W66YN5OQiHnvxcXxTqoTU2x0NGnShCQlkby8vG3bSKtWlP0Ej4ggw4aRli3Jtm2VvraCGHz8+Ll+/caEhf1DjSmCUcUwR0dy/jwFYmuJI3R0XKmvTwoK+A/OySFsNnn8mAK9Z88Se3sK5IgD3T7D3Z20aMH7JZ6OkBBiZsY/faVmZMIRFhQQY2Ny7VrpIeWTQnd34u1dvHVr0bhx1Ai8dYsoKBBPT2qkiYMgjlBSiHnvXbhAHByoF0sIWb6cjBlT+vfRo0Rbm9y/L6Kot2/fpqam3rtHhg0jBgZk27bSWWBFpPY7WMWwRYuIry8FYmtJ+sSLF/3mz4eyMv+RJ05ASUnQlvQ1M3gw7t6VQJdkJomIQOfOwl1S6wNHS1BWxqZNWLAARUUADTuFHh7w9/+0Z8938cNkMjMzg4M/29pi4kRs20aFcXKqgdoNwoosX46QkFlNmvQMCLjo7Ix9++DoiLAwoeWcPx/Utatr8+aDJkxIGDECcXHw9BSlPriUwHC8jLQ7wsTEzS1aCNQw9+RJymK6GjSAmRn/1jkyzbt3Qidy1YV4mRIcHGBoWF7yg9rw0YSEy8nJv8TF2b97d0gcOU+fPm3c2HTo0OFGRicF7P0rR2SoDRmtCIeTy2ZHpqX57dlzGsDw4Th/HpMn48wZ4czz9v6SldVBSanJhQvpLi5QUqLFWsZguAeFtDtCRcUOCQkCRRU/fYrRoynTW7uTKLKzkZMj9MdVdxwhgC1bsGZNaetUaieFkZGRLJYZYP7ixQtx5Dx+/JjL7QYMNjK6TI1lcqqnJJueDtTU1H77baym5vJWrUo7U5T0g/v9d2zezP/yxES4uWHQIEycONPPr8f584u7dRNyqUcqYboHBQWrsLRx4MCBPn0G5v28zv0Tt2//x2K9Sk2lTPXjx8TY+HPFysgMQ+tS/uHDRFW12ler2yNMSiJNm4qlV2r3J3gyfz6ZNav07y1bHterN3DixF85QlWh5ixQI5YAACAASURBVAWHw5k8eebIkZMFubFrwMeHsNk7LC1HfKuh6Bmz1NY9wtxcoqZGioooFluRknDlipEyHz4QExPi5VXt/nRqKvHyIlpaxMuLZGQIqkhqv4NVDONyiZoaEd/YWrJHaGxsUI/fOveFCxcGDpxEyIzQ0NNU6X3y5Ojbt5NatbLM4tvgWQa5ehWGhkJfpacHLpeCBj2ygo8Prl7F48cA8PDhkfz8lZcuvUtKShJTLJvNPnRoz6lTB/je2DXw9SvWrMGGDXMjIgKb/NxZXA6lvHyJdu0ErcEkGtraWLQIixaVnzEwQEQE7tzBjBlV687k5GD9erRrh/R0vHwJX1+JlYChDxYLrVvj3TuG1Em7IxSEZ8+eEaILaD1//pwqmTExbwixy8xUy6iNrdmfPYOFhSgX1uI2FD/ToEF5AdLffptsaLiiuLh5fLy+pO0CgMGD0aYNFi6UtB11AwqLq9XAvHl486ZSNY+S8gvJyRg9ujR2r6gI+/bB2BhPn+L+ffj7Q0eHdsMkBZPxMrXBESoorGSxxltbK69cuZIqmWvW/D5yJMvcfFmLFi2okik9JCVh+HBRLqwjgaNlTJ+O4mKcOoXevXu9f3/n0qX9kyezEyVdEeXoUbx8iWsCxZDJoQD6ImUqoqyMnTvh6YmCgvKT6uq4eBHKypna2vY6OlZt2sQGBSE4GAEBtHRyliqYjJeReUdYWIgNG5QXL55069Z5dXV1qsQ2bNhwz57Fr17ZM1AMl2ESElBUBHt7Ua6tlfEynz9/joyM5PlSSQFSL6/ShiQ2Npg3D46OyM1l1MKK5OfDwwOzZ/Np7yCHQujLnaiCjQ06dKiaCaOsDA+PZ4WFhikpznZ214ODYW7OhDESh8l4GZl3hG5uUFbGX39RL1lbG4aGpVtEtYlTp9CwoYgJRrXPESYmJpqZDevbd8X27cd4DrC0xIAB8PUtPVy8GJ07U1wpWyhGjIC6OnbskJgBdZCXL5mYEZawdSs2bcKnT5VO9ulj6eDAsbK67u09jiE7pAD50qigfP2KY8eweTMUFGiRP3hwLUyiuHkTxsYiXmtqiv/+Y6JlDGNEReWlp6vn5uovWfKja1d4eSE0tGothY0b4e9f/p3cswevX2PLFuaNxd27CA3FuXMSUF1n+fwZAHNbca1bw929ass5FRWVs2f9//33nJ74HSZlB7kjFJTx46Gnh5kz6ZJfK7MJo6PRr5+I16qrQ1eX8+YNl/9QWeDxY7i6Gvv4bDhxwjojw33/fmhqYsMG6OrCxgbr1+PpU3C50NHB1Klv+/ad4u29CUC9eggKwtatEtilGzECDg7o359pvXWZFy+Kzcyo68UlAMuW4d49/PsvkzqlkebNkZUFhsL2xU3ToJOSotvVvRoZSdhscusWjQYUFJCGDUlaGo0qqoOmdB8ul7DZJCKipjHV5RESQmJiYlRVu2tpdX///r0I2qUqhyksjOjolNcUrUh2NgkNJV5epFs30rQpcXIi7dqNB4IUFfuVvfF794iODomLE0V1SfcJYa+aOZOoqvKoHik91L48wgcPHjZo0L1hw57J1dRHp+mWDgggpqa8MxfFRKq+gxXhaVinTuTZM7HE1pI8whqYMAFdusDamkYVysro1w83b9KogmEePACAXr1EvPzJk6eFhTZZWb3ErIoicY4fx5QpuHiRd505dXUMGgRfXzx5gqdPMWQI0tLygJXFxWq+vholtdZ694a3Nxwcyhv50kpcHA4dwoEDMlw9Uha5e/dhVtaooiKTV69eManXyQnNm2P/fiZ1SiOMxcvIqiMMDsZ//+HUKdoV2drixg3atTDGuXNo1gwsloiXjxkz2sYG+vpa9qJFnUoHu3Zh6VKEhAiUTNmiBaZPx7VrK9u3b9agQY+QEK3p00vrcbu7o39/ODtXajpPE3Z26N4dEyfSrkhORWbOdG7c+Pu4ce2srKwYVr11K3x8kJrKsFrpgrFtQll1hC4uGDYMRka0Kxo8GNev066FMcLDxQoEV1NT2737r+JiH2VBGoJIH4TAxwd79iAiQrg4wG7dusbGhr58uZoQvHgBe/vSieDOncjIwOrVNNlbytat+PgRly7Rq0XOz6iqNszP37hnz1I2m+lHZYcOGD8eK1YwrFa6YCyVUCYd4Y4dSE0FMxX3jY2hrAxm10Vo5M0bDBokloRWrZCRgfR0igxiEA4H7u4IDsadOxCtTEKLFggNRWoq2Gz07YukJCgp4cwZHDqEs2epNvf/fP8OLy8sWwZtbbpUyKmOV6/QqpXElqNXr8alS7Uwg0tw5Euj1cLhYNkyzJ7NXHk9G5taEjuan4/MTDg5iSWExYKpqezVlykowIQJiI/HzZsQpzankRGCgxEZid69YWWFV6+grY2gIPz6K16+pM7cCgwdCj09/PknLcLl1ExkJEOp9DypWOSvbiJfGq2WuXMBMJrFVWuSKK5cgZKSKOW2q2BuLmOOMDsbDg5QUsK1a9DQEFeamRmCghAYiIkTYW2NO3fQpQu2bMHw4dTv6Jw7h0ePcPUqxWLlCEhUlITLuEyfDg4HJ09K0gYJoquLnBz8+EG7IhlzhFlZOHAAvr50ZdDz5JdfEBFRGxrWX7ok4pJgFczMUE1JMmkkJQVWVjAxwbFjlHUr7dULJ09i/374+MDJCadPY9IkjBqFceMKHj9+xqUoeKa4GDNmYPp0dOhAiTw5QhMZKWFHyGZj9278/jtDwcnSBouFNm0QH0+7IhlzhOPHo0kT/Poro0pLGtbfvcuoUjp4/Bjdu1Mgx9xcNhxhUFBw9+6ju3S5PHgwdu4EteEOv/yCffuwahX274eXF3x8sH49Hj4cZGW1y8mJggpsGRkZDg6pCgrYt098YXJEhLEqozXQrRtsbbF2rYTNkBTMrI7KkiN89w7Xr+PAAQmoHjy4NiRRfPggYq3tKpiZ4dUrcDgUiKKVX39d+fTp1vx8n7JKodTi6AhfX8ydi7NnERSEBQugoZFdUNAjPPxrfr5Yku/cuaOl1fn6dbtp044yHq4op5QvX8DhoHlzSdsBLFv2Y8sWe0PDfrGxsZK2hWmYCRyVpS/ZmDEwNaXmUS4ssr5NyOVyO3a0yc/vmZ9PQe9idXXo6jJXBlBk9PRGqqqOdHYWqeOUYDg7w8sLkybhzBlERUFf38vS8rqFxTZLSyQkiC52+/YYLrc30D81tc4X2pIcUVHo3FnSRgAAkpKeKykZJCQ4X7xYi3K5BIOZwFGZcYT//ovISCYy6HnSvTuSkyF2c3KJ8erVq1evMoH1q1dvpUSg9K+OfvqE+Pg/3rx5um2bN62KZs/GtGkYMwYnT+Lly78iIma+e+fh4QFLS/zzj9DSuFyMGoWgIPd27bT693+3c6ckanvLASAFG4RlWFpa2tlBUfHG0KHjJW0L08iXRisxeTIGDZJY1ACbjYEDERYmGe3iY2RkpKjYisVa5uw8jBKB0t+hd+VKuLlBn5F+8n/8gaFDMWIEOnVqqaGx+8OHrvHx2LMH48cLt5L/6RMMDHDjBm7fZr16tfOffwIbNWpEm9Vy+CDZ3ImKqKioXLjgN2XK2ZAQKVioZRa5Iyzn6FEkJ0s4hlimV0eVlZU1NU/Pnh26di01lSqkfEb433+4fBmLFjGncd06dO0KNnvTggXWz57Nz8/H3LlYuBBbt8LNrbQkW82cPIlWraChgeRkMF7PSw4PJJ47UYVZs+DvX+dyCnV0kJ9PewaFDDhCQuDpienToaUlSTMGD0ZYGBNVJemAy8W3b3B2VqdKoJQ7Qi8vLFmCxo2Z08hiYfduREWNX7Wq0dixLtu348IFnD0LLS3Ex2PgQKSk1HS5mxumTMGcOYiNRf36TBktp3oKChAfj/btJW1HBSwtoaFRF9sztWlD+zahDDjC339HQQF275awGXp60NHBs2cSNkM0bt0Ci4UePSgTaGCA7GwprQh89y6iojBnDtN62WwYGjZQU7v//r1u5854+BDXr2PyZERGghB0746nT3lclZICQ0McO4YbN7CVmg1cORQQG4u2baGiImk7KjNjRl1sScHA6qi0O8LiYqUdO+DtTVkqtDjI7upoQAD09ERvOvEzJYXWoqMpE0ghS5Zg9WrJPMKePw+9e3deZubegwcRHY22bREWht270bo1cnMxYACOH680PjgYhoZQUMDHj/jlFwkYLKc6pCdSpiLOzggOxrdvkraDWeSOEOHhrg0bYskSSdsBABg8WFYd4b171KTSV0Q6C61duICsLIm1K1JRUencuTOLxerWDf7+ePcOgwbhr7/w8CGcnNCsGVxdMXNm6QK7uzscHDBuHOLjq5Y/zc7OtrGZ0KvX8E+fPknkjciRhlT6n2nYEMOHV/05VethIJVQ2h3hhw9/ODlRkPpGCVZWiIxkovAd5bx7B0dHimVKYaE1DgcrVmD9eoqLyIhMo0ZwdcWLFzh5EoQgIwNt2uDIkUv16rWuV6/7gQMxZ89i3z58+IC7d3HyJDZuhKcnRoxA167hN2+2ePTI9uRJee8lySCdM0IAs2Zh3766FTLDQCqhdDwwamJoaOgmSdtQSr16sLTE7duStkNIkpJQUIBRoygWK4XxMocOQVubd9N5yVI2QfTwAJvtCYwEnFisOb/+Cg0N9O+PpUtx5Qq+fkWbNpg6Ff7+vbt0+U9TM+jgwSGJiZK2vk4inTNCAH37QkGhNlR8FJzz5zc+fNjd3X0pfSoU6RNNEbs1NecUFkJKGsGWbBOOGCFpO4Th+HE0aEB9LKKpKf77D8XFUJSOmyg3F3/+icBASdtRPQ0bYvZssNleHh5rAfaECTN9faGjw3P+2vDp00sAtm5F7964dAldujBubh3m0yewWNDRkbQd1TBzJvbvR79+kraDKS5ePM/lngsKGrl37zqaVEj7jHDMmJ7Nm88dOBBfv0raFACAra3sNawPCaElClxdHfr6ePOGesmisXUr+valMjKWJtzd3Z8/vxwRcfrvv1c0b85nFXf+fOzcCVtbXLnClH1ygMhIaSmuxpMpU3D5Mr5/l7QdTLF795oGDRa6uq6kT4W0O0INDYXz5zFoEHr3hjTUm+3QARwOQ02TqSI6GgMH0iJZelZHU1OxbRtWr5a0HYLRsWPHrl27Cjh4xAhcuQI3N8lnENUdpC2VvgpaWhg6tA6FzAwZMmjy5PNNmoykT4W0O0IALBZ8fLBmDaytpaJDqWw1rC8qwvfvmDSJFuHSU2ht7VqMHw8jI0nbQQ8WFrh7F7t2wdNTVks6yBbSU1ytOmbNwt69kjaCQTp1ojdZSwYcYQnjxyMoCK6u2LFDwpbIVjbhlStQVKSrRquUzAg/fMCxY/jjD0nbQSetWiEiApGRGDu2NvSIlnKkNmS0jP79ASAiQtJ2MIWZmdwR/h9LS9y9i337BC3eSBM2Nvj3XxQUSMwAoTh/Hi1b0iVcSjIoli/H3LnSG9pAFY0bIyQE9eph4MA6l1LNJAUFSEhAu3aStoMfdarKTKdOiImhcTlElhwhgFat8OABvnzB0KHIyJCMDY0bo3173L8vGe3C8vAhevakS3jLlsjLk3AcU1QUwsKwYIEkbWAMFRUcOwZbW1ha4vVrSVtTS3n5EkZG0hKmXgPTpyMoCOnpkraDETQ00LQp4uPpki9jjhBA/foIDETXrrCwkNizwNZWZhrWJyRg9Gi6hLNYtK/d8+X337FiBTQ0JGkDk5RsmS9bhv79cejQi1u3bpE6lVxNP9K/LlqClhbs7CTck4dJaH3UyJ4jBKCgAF9f/P47BgyQTHq7rNRai4tDcTGGDqVRhWQLrf3zD+LiMGuWxAyQFDNmwMcnetasuQ4O+/7+u848CxlB+iNlyiipMlNHoDU0TyYdYQkzZyIgABMnws+PadW9eiEhQVpSG2vg2DFoatJbflqC24SEYNEirFsnA6tYdGBlpdCgAScvr+j9ewVJ21KrkPLciYpYWyMvDw8fStoORpDPCKulXz+Eh2Pr1u/6+k79+o1OZaotkIIC+vdHaCgz2kQnLAympvSqkGDg6JkzYLHg5CQZ7RKnQ4cO9+4d2LRpgb//ONlKbJVyoqNlZkbIYsHFpa6EzNS2GWF2dvbx48dnz57t4uJy4MCBIvECQNu2hZdX6JcvZvfvdw9hcL1SJpIoYmMxeDC9KkxN8eYN00G8hJDQ0NteXjG+vlT2lpI5OnTosGBBn9WrWQ4OdSVogm4SE6GsjGbNJG2HwLi4IDAQmZmStoN+jIzw5Quys2kRztsRPn782NfXd+LEiYMHD7a3t582bZqfn99bin523rx589ixYx07duzXr9/27dtnzJghpkAHB2tT0wgW69+2bZlr6VYSLyPNYQq5ucjMpL0hkaoqWrZkOmpp374jw4fvT052adFCaiq8SY6ZM2FvjxEjUFgoaVNkHxnaICyhSRMMHFgnQmYUFNCuHWJiaBFeyREWFRX5+fl17NjRwsJi+fLl9+7d+/btW2Ji4o0bN2bPnm1sbGxtbX3pkrh9YRwcHEJCQubMmTNt2rQjR46cOXOmULxvcLNmzSIjry1bdt3Pj7lUMgMDNGokFVl01XH+PJSVYWhIuyLmV0c5HHZhIVFTA6suzwcrsHEjtLTg7i5pO2QfWQkZrUjdqTJDX1p9uSOMjo42NTVdsWKFtbX1zZs3MzIyEhISnj9/Hh0d/fnz55SUlPPnz2tra48bN87a2jpDjCQ+doUyw1+/fm3UqJESFe3nFy/GjRuMPpGlPIni4kW0bs2EIuYLreXnO/ft6x4RcdiothZVExI2GydOIDYW6+iqzl9XkKFImTJsbJCTg6dPJW0H/dAXL1PeQScxMdHNzc3Dw0NVVfXncc2aNRs5cuTIkSNTUlLWr1+fkZHRqFEjMXVnZWXNnz/f29u7ht/19+7dG1W5k97QoUPHjRvHc/C8eUpLlyoEBOSLaZiAWFkp7typNHs2LQWvcnJyxJzuPHqkNnAgJztb6BI4HA6noKCAK3AVBxMTxb17lbKz+X8O4r8pAHl52LxZ/dy5Hi1acLNp2jGgmaKiIi6XK+bu+M+cOMEaOFBNW7tg7NhiaiULCCEkNzdXIqr5IuC99/y52sKF+dnZgt78lNzS4jNpkvKePazt2/l/2aXE4J8RxLC2bRUCA5UFedSUweVyVQSIm2dJKhs3NzfX3t7exMRk79691b3/gwcPnj59elblNDETE5NOnTrxHF9YiI4d2fv3cwcMoNxeHuTkQE+PnZTEpbzVH4CsrCwN8bLEVVTYly5xbW2FvrDEEaqpqQk4PjERvXqxP33i/+wQ/00B2LSJ9eQJTp+W4u1ZfpQ4QkG+n8ISGwsbG/bZs9zevSmXzR9CSE5OTn06vg9iI8i9l5eHZs3Y379zBV+iouSWFp/kZJiast+/5/K1RUoM/hlBDPv2DR06sL99E6LSGpfLJYTwXXSs1FN1y5YtQ4cONTExEVyNaOTn548cOdLAwMDPz6/mXwEtWrQYO3asgGLr1SspusG+f5+JYEINDfTogfBwNh0Z62w2m11zq7oaefoUXC5sbESRQQgRSruBAYqK8O0bW1ubz0gx3xSAnBxs3YobN8BmS+OvWgEp+RDE/Ch4YmqKI0fg5MQOD0fbtpSL54OAdw4hhPlJiSCGxcTAxAQqKkL8X8S/pSmheXMMGICAADbf4hICGsz8/0gQw7S1oaSEL1/YenpCSOZwOPy1VzzYsWNHu3btunfvvm/fvqysLCFUCUNhYaGTk1ODBg0OHjxI+T00aRIKCnD5MrVSq0VqkyhOnULTpsz1ju/UiaFtwp07YW2NalYE5ACArS1Wr4bUJlS8ffu2VSvLVq0sP378KGlbqiKLG4RlzJqFHTteffnyRXxRzs7zdHQsVq7cJL4oyqEpIqGSHwoPD/f19f3+/bubm1uzZs3Gjh0bFhZG+dppQEDAlStXHj16ZGJi0qZNmzZt2lDyzyuBzcaaNfDyQjEjuyRSGy/zzz+MtthmptBadja2b8eKFbQrknWkOaEiIuJ+YuLgDx/6Dh/+eMsWvHolaYMqIHO5ExVJTz8bGzu/Y8eh79+/F1NUcPA/X78e9/UNnDsXwcGQqm1fJhxhixYtvLy83r59Gx4e7uzsHBwcbGNjY2BgsGTJkg8fPlCl0tHRMT4+/vbt26H/p2nTplQJBzB0KHR1GWrfbGaGzEyIfeNRz+vXsLNjTh0zhda2b8egQWjfnnZFtQCpTaho126ksnLGuHGctWvtExJgbw9dXTg74+xZ/PghYdtkMXeijLS07woKBrm5jcRfzNPT8zY29j5yxLdtW2zdiubNYWOD9eulIjCVrsBRUj0ZGRl79+61tLQEoKCgYGdnFxAQUMN4yjlw4MD06dNFuPDhQ6KvT3JzKbeIB1OmEH9/6sVmZmaKfG16OmGxyJcvIl5eXFyck5Mj1CVPnhBzc/7DxHlTP36Qpk3Jf/+JLECKKCwszM/Pp1tLbi7p2ZMsWPA6PDy8JGSAbrhcblZWVs1jHB3J1q2VzsTHE39/MmwYUVcn3boRb2/y5AnJy8sPCgqKi4ujyja+9x6XSxo3Jt++USyWMYqKig4ePN248T81f0f4Gnz3LjEwIAUF5WdyckhoKPntN2JgQAwNiasrCQggGRkkJibm8uXLRUVFVJgv6Cf57BkxMxNCLIfDKSws5DusJkdYxqtXr3799deS/TwhTBAbkR0hIWTECLJ5M7Xm8Ob4cTJqFPVixfmC+fkRNTXRVYvgCPPziZpapS8PT8R5Uz4+ZNo0ka+WLphxhISQiIg4BQULNTWnffuOMKCOryN89Ijo65O8PN6v5uSQK1fInDmkTRuirr5EUXF+kybmVHkavnI+fCB6etSLZZhVq/h8Tfga3L8/OXy42ldjYsimTWTQIFK//mclpW4qKu4rV64XxVDhDStBwEdNGQI6Qv7xXXfu3NmwYcPRo0e5XG6bNm1omJTSgq8vfH2ZiBdITz8TFNTdzm4K7ZoE5soVMJxorqICQ0P89x9d8n/8wK5dWLqULvm1FS0tboMGCnl5KrdvSyazsAorVuCPP1CvHu9X1dQwdCh27cLbt3ByUmazc7KzCWMxmTK9LlqGpyeuXEFcnIiXh4UhORlTqn+YdeiAhQsRGoroaAU1NU5hYU52NgXlUASHpkdNtTdZUlLS+vXrjY2N+/fvf+rUKTs7u9DQ0DiRP2DGMTGBgwM2b6Zd0fnzF7jcnffvv8rLoyWzXgSeP0e/fkwrpbXQ2pYtcHSEsTFd8msrxsbG//yz98iRKU+eTF+yRMLG3L2L168hYGnhAwdWBgaO1da+9uCBOs12lSLTkTJlNGiAOXOwYYOIl69YAR8fKAjQ18vQsFlkZJCbm0t0tKeIykSFlm3CKjPEvLy8gICAYcOGKSgoAOjWrdu2bdvS0tIEnYhSijhLo4SQpCSipUUSEym0iAdPnjxr02ZM585+1IoVecmFyyVsNrl7V3TVIiyNEkLWrSOLFvEZI9qbSk8nTZuS+HgRLpVSGFsaLSM5mXTtSubMIRwOjVpqXhodMIAcEXKB9to1YmRU7VKqUPC998aMIadOUS+WeTIySNOm5N073q/WYPClS8TUVLg7pKiImJsTSkJHBP8k16whXl6CihVlaXTBggXa2tpjx4598uTJggULXr169eTJE09PT01NTar9LxPo6WH6dKxdS6+Wbt26PH169sMH92/f6FUkIOHhAGBpybRe+maEGzdi5EiG6qbWVrS18e+/iIvDmDHIZ6gEYSVu3MDnz5g0Sbir7OxgaoqNG+mxqTK1Y2kUQMOGmDVL6EkhIfD2xpo1EGopWlER/v6YP5/RiF86MigqvemLFy9aWFgEBAR8/Phxw4YN7dq1o1gb4yxbhvPnady7KqFhQzg44MQJerUIyOnTaN5cuLuZEszM8OIF9WLT0uDvD4kv69UC6tfH5ctQUoK9vQTa13l7Y9UqUSo8bN+O7dvxhuZ2Wzk5+PSJ6Z11+liwAAEBSEgQ4pKzZ8FmY/hwoXX17Ak7O/j4CH2hyNDuCJ89exYaGurk5ERJOwhpoHFjLFgAb2/aFc2ciX37aNciCOHh6NpVAnpLih4lJ1MsduNGjBuHVq0oFls3UVbGyZMwMcHAgfj6lTm9ly4hOxtOTqJc26IFlizB3LlU21SZ6Gi0b89cJSa60dLCzJlCREhwOPDxwdq1Ilam3LABp0/j+XNRrhWBli2Rk4O0NCplVnKEDRs2LPv7w4cPy5Ytc3R0tLe3Lzlz6tSpa9euUamcETw9cf8+Hj6kV4uVFQA8eECvFkGIj4eDg2RUd+pE8epoaioOHpQHi1KJggL8/ODkBEtLUNRpmw8la25r14q+SjFvHr5+xdmzlJpVGZkursaTxYtx8iQ+fRJo8IkT0NSECAX6S9DUxNq1cHODwB1rxILFgqkpxfEyvO/NBw8emJubHzhw4NOnT1H/n4UmJibOnz+fSuWMoKqK5cuZWFubNg0HDtCupWaSk5Gfj9GjJaOd8kJr69djwgTo61MpUw4ALy8sXIj+/WlZza5CQAAUFcX6ccbARlSt2SAso0kTTJ0q0KSwqAirVmHNGrHUTZ8OdXXmVsUo79DL2xF6eHj07NkzPj5+06byuqtDhw59/fp1SkoKlfoZwcUFKSkIC6NXy9SpuHBBArsvFTlxAvXro3FjyWinttBacjIOHYKXF2UC5VRk9mzs2oUhQ0qjq2iCw8Gff4q+5laGhQWGDKFxj6N25E5UYfFiHD0KvoWcDx9G69YQs3Udi4Vdu7BiBX91lEB5BgUPR/j9+/cXL178+eefGhoaFTtxGBgYAPj8+TOV+hlBQQGrVmHRInpn7trasLbGmTM0quDLtWuQYIQTtTPC9evh7Ayh+q3IEYqRI3HyJEaPxpYtL69duyZItxphOX4cmpoYPJgCUevXIyAAz55RIKoKhODly1rY0kRHB5MmYevW1UrMvAAAIABJREFUmsYUFsLXF6tWUaCuY0fMmMHQL1fK293wcISFhYUA1NWr5rGmpaUBUJTNDeXRo6GqSu82A4CZM3HwIL0qaiYqCtbWEtPeoQPi4mJTUjLEF/XlC/7+G4sXiy9JTk1YW+PgwQ+LF7uMHHl2yxY/aoUXFWH1anHX3Moo24ii3F8/eZKipvZOS4tisdLA0qU4dAg1ZHb5+aFTJ/TqRY06b2/cvYtbt6iRVgNmZoiJoXJiw8MRNmvWTFtb+9KlSwAqzghPnDhRv359Btr20gGLBV9fLF1Kb2MaW1skJzPRh4EnxcVIS8PEiZLRDmD3br+ioiWmpv0zMsT1hevWYcYMNG9OiV1yaqJzZ4VGjTjFxXn+/kqvX1Mp+dAhtGkj7ppbRaZNg4YG9u+nTCCAuLg4Gxv7tDSXy5evUilXOtDVxdix2LaN96t5edi0CX/+SZk6NTXs3g0PDxQUUCaTJxoaaNIE795RJpCHI2Sz2QsXLly9erW3t/ebN2+4XO6LFy+WLl3q4+Pj6emprKxMmXJm6d8fRkb0ztjYbEybhkOHaFRRA9evg82W5J7/p08pQKfc3Hq54nUw+/IFp07Jp4MM0aJFi0ePzty+PWfxYtc+fbBkCTVPsfx8rF1LzZpbGSwWdu6keCMqIyMjP78pYJyczGBCCYMsXYq9e5GayuOlHTvQuzfFjUuHDEG7dthEf09firMJedab4XK5S5YsqZhNyGazXV1dBalVQyFillj7maCg/xQVrXv1GsG3WYzIfPxItLQo6AAlQummadNIq1bi6iWillgjhOTk5Hh7H2zSJKK620SQNxUVFTVo0PGFCxnpoSUhmC+xJiCfP5MxY4iREbl1S5TLK5ZY27qVODpSaVsZXl5k8mShr6ru3ktPJ/XrX9u162RxcbEIxkhhibWfmTWLrFxZ+neZwVlZREeHxMZSry4hQZSaiMJ+ksuXE29v/sPE6j7BYrHWrVuXkJBw4sSJTZs2HT58+M2bN/7+/rKeaB8Tc6W4eMqLF82f0tZiskULdO+OwECaxNfE/fuwsJCA3jLU1NR8fGa0b2955YqIEjIyMqytp4aFxaSlUbdeI0dgdHVx9iw2b8a0aXB2xvfvIsrJycGGDXRVG1m5Evfu4eZNaqQdPoyRI+3mzJmgIEipadnkjz+wZ0/VVjybN2PwYFraXLdsiUWLMHs29ZIrQm3gaE05rrq6uhMnTly4cOG0adNkqAFTDUyaNMbc/AKXm92yJY0ew8VFMiEz79/D0VECeqvg4QE/UaMulJSUcnI4qqofmzdvyH+0HHpwcEBUFBo3hpkZ/v5bFAk7d8LKiuI1tzIo3IgiBP7+8PCgwiwpxsAAw4Zh587yMxkZ2LMHK1fSpXHBAiQn4/x5uuSD8lTCsrlhSkqKgCuf379/z87OFmSkmFC+NFrC4sXE1ZVyqeUUFBBtbSJmb21hFwri4wmLRURa0ayKyEujJRQUEB0d3q3k+b6pGzeIoWFqePhDZjqqSwqpXRqtwr17pGNHMnQoSUgQaHzJ0uiPH6RpU1rW3Cri6EjWrBFiPM9778YNYm4ulhkysTRKCHn7lmhpkYyMUoOXLKH3GUgICQ8nzZuTjAxBxwv7SRYXE3V1wnePS+il0Rs3bpiYmOzbt+9H9fUbEhMTvb29W7duLYvZhGUsW4aLF2noaPV/lJUxZQrTITPHj6NRI6ipMaqUJ8rKmDED/v5CX1hYiLlzsWuXVt++FiwxE7DlUEHv3qW9Lbt04Vpa/tGv35i3ApRl27IF9va0rLlVZOdObN+O+HixhPj50b6CJyW0aYMhQ7BnDwB8+4YDB7B8Ob0a+/aFrS3F0VIVUVBAu3aIjaVIXJlLLCws3LFjh5aWlqqqqoODw9q1a8+ePRsWFhYSEnLy5MmlS5f269ePzWabmJhcvHhRKNctMjTNCAkh27eTIUPoEFzKf/8RHR0iTmiRsL+P+vYlffqIrq4iYs4ICSEJCURLi8f0tOY3tXo1GTlSHLUyg6zMCMsIDo5VUhoJXDAy8tq2jdy9y3vtgcvlJiZmM9Y8cv16Ymsr6OCf771Pn0iTJkTMtS1ZmRESQt68Ic2akc+fs+bNI56eTGhMTSU6OuTZM4EGixYeuH8/nzECzgirRo1mZWX5+/v36dOnSuJ8/fr1hw0bFhQUJFpslWjQ5wgLC4mxMQkJoUN2KX37kqAg0S8X9rZo3Jj4+IiuriLiO0JCyPDh5ODBqidreFMJCaRJk2q7idYyZM4R5ufn9+7tqKPTc8mS8DlziIUFUVcnZmbExYX4+ZEnT0p/88XGxtraXnZ3Z8iqoiLStu39Pn3m3Lt3n+/gn++9lSvJb7+Ja4MMOUJCyNixRba2+xs1Sv/0iSGNe/bkamt7//GHb1FRUc0jRfgkN2/m/x8U0RGWkZ2d/fjx46tXr964cSMqKorv26AD+hwhIeTCBWJmRuhz64cPEwcH0S8X6rbIzycsFnnzRnR1FaHEEV67Rrp0qXqyhjfl6EjWrhVTp8wgc47wZwoLyZMnxM+PuLgQMzOirk46d/7AYukDY/v0mciYGTo6PYD7+vo9+I6scu8VFRE9PfLypbgGyJYjNDDoC0xhs40Ye54fPHiYxVqiojLr6tWrNY8U4ZMMDSUDBvAZI1b6BAB1dfXu3bvb29vb2Nh06tRJRiur1cDIkWjSBEeO0CV/3DhERCAxkS75FQkKgpKSdLUVtbVFdjYePRJocEgIYmKwcCHNNsmhDiUldOsGd3ccOIDISHz9ikWLcgAVoHl+Pq/kbXoYOLCnurpXfn5PYSNIL1yAsTE6dqTHLGmlqCgPaAvkFtJaYasC3bp1ado0rLj4WVpaB8qFU1jcmPFG5tLEpk1YsQJZWbQIV1WFk5OI0efCEhgIQ0MmFAkOi4VZswTKoygowG+/YccOqKjQb5YcelBTw6RJHU6eXOPqWnDrFs0lfStw4sTO5OTgwYN3OjuDECEu9POr/VkTP/PkyeWpUxOuXTukxlRYnbm5+cePd2/derBwoeHLlxQLb9oUSkqC9lysmTrtCLt0wcCBQvRxFpaSGtwMNKt89AiWlrRrEZYZM3DxIv9G0uvWwcwMQ4YwYpMcOhk3btzmzRsaNGjApNL69dUPHsTHj/jrL0EvefUKr19jxAg6zZJKdHV1d+7cNpiSbiACo6KiYmWluGkThg+vqfy3aFBVaK1OO0IAvr7YtYuuBcxu3dCoERO12JOSJNaMtwa0tODgwGfxOSEBu3dj40aGTJJTK6lXD4GB2L8fp08LNN7PD7NmQcbLZMkYzs4YPx6jRlFcj5uqtPq67gj19TFrFpX116vAQJWZ6GgUF8PWll4tolFSZaaGOfGvv2LRIqlb15Ujc+jo4OpVzJvHf1s6NxenT2PWLEbMklOBtWuhrw83NyplUlVora47QgDLliE4mJaGnwAmT8b169QvCFTkxAk0aQLpbArSqxcaNUJYGO9XL13C69eYN49Zm+TUUjp2xOHDGD2azwLP8ePo1w/6+kyZJef/sFg4dAivXmHDBspkUtWhV+4IoaGBP/6gq+lPw4ZwcMCJE7QIL+HWLZiZ0ShfTNzdeYfM5OVh3jx5jIwcKhkyBL/9BkdH5ORUO6YuFBeVWlRVERSE3bsRFESNwI4d8fYtBV1m5Y4QANzckJKC69dpET5zJvbto0UygNTU1JiYG4MG0dwHUwwmTcLdu0hIqHp+3TpYWMDOThI2yam9LF4MCwtMmcJ7QT4iApmZ+OUXxs2S8390dREUBFdXPH5MgTQVFRgYQPyG0nJHCACKivjrLyxYgOJi6oVbWQHAgwfUSyaEGBhY5OaGnDghBV0nqkFVFZMm4cCBSifj47F3L43xunLqMjt34scP/PEHj5dKiovKC9lKli5dsHcvRo1CUhIF0igJHJU7wlKGD4eeHg4fpkX4tGlVPQElcLncvDwFQD0x8SP10qljzhwcPFhp+cLTE7//Dj09ydkkp/aipISzZ3H+PPbvr3Q+NRVXr2LqVAmZJacCo0bBwwOOjsjNFVcUJfEy5Y4wLS0tJiZGXHmyzMaNWLkSmZnUS546FRcuUC85L0+Bzb6mqfng1KltFIumFCMjdOhQ3qz4wgW8fw9PT4naJKdWo6mJy5exfDlu3y4/efAgRo2CpqbkzJJTgWXL0LUrnJ3FzbS+cWPutm09/P2PiSOk3BFeunRp7NixJX/37NnzGU1hlFJM584YPBibNlEvWVsb1tY4c4Zise7u0NBo++3bjSFDGM2QFYGybr15eVi0CDt3yrO45NCLiQkCAjBxIko6R3G58PeHu7ukzZJTgT17kJYGHx/RJRBCoqP/LSj429//lDiWlDvC+vXr5/w/1iohISE/P18cuTLKX39hzx58pGGhsaTKDIWkpuLMGfj6gi0Ly9uOjoiPR0wMe80a9O6NgQMlbZCcOkD//lizBg4OyMhAaKhi06bo3l3SNsmpgJISzp3DqVM4flxECSwWa8WKXxs0+H3o0GXiWFJeStvc3DwxMdHNzc3c3DwvL+/ixYsvedWGc3V1FUeflKOnB3d3eHtTv1loawsPD0RGwtycGoEzZkBTE7Ly31BUhKPjf2PHbsrMnB8bW8dKHcuRHC4uiIqCldWxb9+SvbymAU0lbZGcSmhp4dIl9O//6ciRZTNm2E6cOPF/7N15XE3pHwfwz73t+4oIqSRrm0hCWYfIkiEklL2x72PfxjKYQWMXEYbsQslaREX2fZlkzRZp0Xbv+f1x/VK5+z33nHvref8xL84993m+mep7z7N8H1lbmDhxRO3aI9auVSyO0kdRhIWFWVlZSX+/sin1GCZRsrOpGjWoa9fob3ny5Lfe3usePXokzc3iDyV5+ZLS0KD+/ZemyMqi5RimnxkZ1QfCOBx72ltWOxXgGCah+Hx+dnY221GUl5b2gsOxBebb27emvXH1OoaJUtWAq1VrASzkcOpkZGTI8XbBoVp37gh5SZ5jmMaMGfP27ducnJwqVarExsZmCqNY2lUDhoaYMwdTptDfclxccHy8drt2/RVvKigINWqgXz/FW2KOqakRlxuvpaW6Wx6JCsnY2IDLzedwkqysTNmOhRDOxsaMw7lGUQbFxQZyvF1TEyEhCm3XFjK/ZGBgMHfu3MaNG5sJI7jn7du3MTEx8ner2oYP/77Sml62tla6upfz8kwUbOfJE1y8WH5puOr777/LW7Z0efeOjsqABCE1c3Pz9PSrO3YMuHTpGNuxEMIlJ5+MiOg9fHj8gAGG8lXlHjECe/aIqygknvCFFmPGjLEWu8nr+vXrEydOlLNPlaehgcDAhF69PH18+hQVFdHV7OHD4efOjatSJVbB5aMDBsDBQUWrbIuhqanZp08fU1PyqZxgmrW1dc9KeOqSWunVq9fGjRYWFggOlu1oSYGaNeHlJe3xIz9ThxWHbHj8+GhR0eyrV/PTfy4OJi8ul+vp6bJnj/a4cfIvTL16FdevK2V7PkEQBIu4XOzejbQ0LF4sz9tHjsSmTfJ2Lef7Krpp04Y3bbqFy3W6ccOe3pabNsXYsSJrIUo0eDBcXeHlRW9QBEEQ7NPTw7Fj2LkTkbLvj+/cGR8/IjVVnn5JIhSufv36164dSUj4Y8wYDu31dmbOBIeDv/+W+Y3nz+PxY2XVgSMIgmBdlSo4dgxTpsh8pDmXi+HDsXGjPJ2SRCiOqytWrYK/P758obNZLheRkVixAjdvyvbGkBB4eaFJEzqDIQiCUCkNGiAqCoGBMh8rMWwYDh7E588y90gSoQQDB6JLFwQEgMejs9latbBiBQID8e2btG85dgwvXyIigs4wCIIgVJC3N5Ysga8v3r+X4V1VqqBTJ3nOfyWJULJVq8DhYO5cmpsNCoKzM6ZPl/b+UaPQpQtsbWkOgyAIQgUFB6N/f3TrJtsJFSNHYv16mdedkkQomYYG9uzBvn30V81evx7Hjkm1YXHLFrx/r357BwmCIOS2aBHq18fgwcjOlnaHYNu24HCQmChbR3Imwnbt2p09e1a+96ojc3McPIgxY0DvmRympoiMxPDhePdO3G18PmbMQN++kFT/jiAIouLgcLBlCy5dmmZl5duz51Ap3zVsmMxLZoQnwkWLFkVHR5e7uH379kaNvpdL1tPTE7/jvuJxdsamTejdGx8/0tls69YICpKwh/Svv5CdjXXr6OyXIAhC9enoQF8/MS9vaVzcjexsqd4SHIyTJ2WbXBSSCHk83rJly3g8HoDMzMySMyjatWt3//79V69eydB8xeLvj7590a8fiovpbHbxYnz8KHLYs6gICxYgJAT/L29HEARRiezbF9a377++vmudnJCUJPl+U1P06IEdO2ToQkgi/PjxY15enp2dHYCYmBh/f3/B9Ro1anA4nIyMDBmar3CWLoWODmbMoLNNLS3s2IFZs4SvFV64EEVF+OsvOnskCIJQF+7ubvv2hR040Oqvv9CrF9askbwWRlBlRvqiJSLnCIuLiwFkZ2dn//9xVHDAirGxsbRtV0RcLnbtwtGjNO9qb9AACxZgwAAUFpa5np+PVaswbhz09ensjiAIQu306oWkJERFoWdPfPok7s4WLWBiAunXsQhJhFWqVDEzM9u2bVteXt7evXtzcnISExMBREZGGhgY1KpVS+bwKxYzM0RHY8YMXL1KZ7OjR8PaGgsWlLk4eTI4HCxaRGdHBEEQasrGBvHxcHVF06a4eFHcnTJVmRGSCLlc7uTJk9etW2dkZPT58+fff/+9U6dOHh4ekydPHjFihJ6enuzBVzT162PTJvTpI9t8rHgcDsLDERGBCxe+X8nJwdatmD0bOjq09UIQBKHWNDUxfz62bkX//pg/X2Spk4EDceECXr+Wrk2hV2fNmuXq6pqWltavXz8TExMOh3Pjxo3+/fuPHTtW3uArmp49cf06/P1x7hy0telps0oVhIdj0CDcugVNTYweDT09TJ1KT+MEQRAVRocOSE3FoEGws1uqoXHyr79+79nTt/QNhoYICMC2bZyZMyW3xqHkOPqJKeHh4YmJidu2bWM7EOH4fPTsCRsbhIXR2WxoKL5+xaJFOfXqGa5Zg9BQOhuXBo/HKygo0FfCtGR2draRkRHtzaqdoqIiPp+vU+Ge9CmKys3NNTQ0ZDsQIZT0vad239IqG7B8gRUXUyYmLnl5x42NR8yaFePjAzc3aP7/+e72bXTrhsePi3R1tcS3QyrLyE9welZMzD0np+BNm3bS1exff+Ho0YF2dm00NMJHj6arVYIgiIpGU5MzfHivOnUCxo0blZmJ0FCYm6NjR8yfjzNn4ODAy84etWTJLcntMBBrBWZkBEvLRcnJIyZOnNKtWx9raxomUHV1kZNzHnhYUNCIw5G2mAJBEEQltHr1/NWr55f89dMnXLyICxcwZQqePr2Um/syPd1RYiPkiVBRgwd3sbScbmpq5+SkN368/EfPC5w/j+7dAXCBppqaHHpCJAiCqBwsLNCzJ1avxs2bSEqy19R80LDhHYnvIolQUaNHD87ISHrz5sCDBzAzg5sb/PxkPiW5sBD796NFC4wcifbt8exZyty5/d++leusZYIgCAJo3Ljm48dnQ0IcJN5JEiENNDQ0AFStivnzkZaGDh3QowdatcJP5VqFyMrCmjWoWxdr1nwvLjN+POzsqk+ZMsXS0lLpoRMEQVRcNjY2pqamEm+TkAgLCgri4+NzcnJoiqriMzLC+PH47z+MGIFp0+Dujp07he90efYM48fD1hapqThxApcuwc8PHDIaShAEwSwJifDdu3c+Pj5PnjxhJpoKQ1sbgwbh3j3Mm4f161G/PtasQX4+Hjx4wOPxLl1C375o1Qp6erhzBzt3okkTtiMmCIKorMiqUSXicuHnBz8/nDmDP//EzJkjgFyKem9re3rCBOzYAVKlhyAIgnUkETKhQwd06IAGDdIfPhxkYbHq7l0yBEoQBKEqyGIZ5sTFbV28OCM+PpJkQYIgCNUh4YmwVq1aOTk5pNA2LWrVqjVr1mS2oyAIgiDKkJAIORyOgYEBM6EQBEEQBPPI0ChBEARRqZFESBAEQVRqJBESBEEQlRpJhARBEESlJlUiLCwsfPr0aXFxsbKjIQiCIAiGCU+EU6dO/fvvvwV/vnv3bp06dRwcHKpXr37lyhUGYyMIgiAIpROSCPl8/saNG+3s7AR/nTx5spaW1vbt293c3EaMGEFRFLMREgRBEIQSCdlH+Pnz55ycnPr16wPIyso6f/78P//8M2TIEG9vbzs7u9evX9esWVPBXu/evbthw4bc3Nw+ffp07dpVwdYIgiAIQm5CnggFz3xcLhfA2bNni4qKOnXqBMDa2hpARkaGgl2+evWqdevWNWrUaN++fXBw8LFjxxRskCAIgiDkJuSJ0MLCwszMLDo6etKkSREREY6OjnXq1AHw+vVrAGZmZgp2uXnz5o4dO86aNQtAbm7uypUru3fvrmCbBEEQBCEfIU+EHA5n8uTJU6ZMqVatWnR09NixYwXXz507Z2xsbGNjo2CXSUlJPj4+gj97e3snJyeTeUeCIAiCLcJrjc6aNatRo0bXr193c3Pr2bOn4OK3b9/mzJmjqanoyU0ZGRkWFhaCP1taWhYWFn769MnS0lLozYmJif7+/qWvdO3aNSAgQMEYVF9ubi6HpVMqeDxeQUEBn8+nvWUWvyiVUlRUxOfzi4qK2A6EZhRF5eXlsR2FcEr63lO7b2mVDVhJgfH5fB0dHYm3icxqPXv2LEmBAmPGjKEhLkBXV7ewsFDw54KCAgBiTreoXbt2v379Sl9xdHTU19enJRJVxuPx2PoyeTyehoaGMnpn8YtSKYJEKM3Pp3qhKIqiKNX8X6yk7z21+5ZW2YCVFBifz5dmxFFkIszLyzt8+PC9e/d4PN7y5csBXL161cjISLCaVBE1a9Z88eKF4M8vXrwwMzMTc8BFrVq1+vbtq2CP6ojL5QrWKzGPoigl9c7iF6VSBP8IFe+fQnnfOYoj39ICKhuw8gLj8XiSexd6NS0tzcnJKSgoaPPmzbt37xZcPHz4cFBQkOJh9e7de9++fYJnwcjIyF9//VXxNgmCIAhCPsIT4bBhw/T09B4/fnzgwIGSi717905NTf38+bOCXfbp06d69epubm4+Pj6xsbGC5aMEQRAEwQohQ6PZ2dkXLlyIjY2tW7euYMuEgL29PUVRr169UnAHhba2dmxs7I0bN/Ly8po1a1bxZkoIgiAINSI8EfL5fMH2+dK+ffsGgJbS2xwOx83NTfF2CIIgCEJBQoZGq1atamxsfPHiRQCl17OePHlSS0vLwcGBueiA3bsvXbp0ickeCUKp7t+/r6Njp6tbNy4uju1YCIIAhD4RampqBgcH//7772ZmZqampgDy8vIOHTo0derUgQMHGhoaMhlfYeGkwYPHPXt2nclOCYJexcV49gx37+LBAyxZ0p+i+gKW/v5zc3I6sR0aQRAitk8sW7YsLS0tICBAsJ7VxMSkuLi4devWJWczMeiYnt4MxjslCHkcPnxi165jM2aMNjZ2uXcP9+9/T35PnsDaGo0bo0EDeHoGnzu3EjDNyztsZ4fwcLRty3bcBFG5ccRsNoyPj4+Li/v48aOhoaGPj0/Xrl0Z3oASHh4+f37T169dDhxA2fIyFV92draRkRErXQsqyyhjcyuLXxQD+HyYmjplZ2/icufY2Z1p1AgNG0Lw3wYNoKsLADk5cHLC588FurrYv19n+XKcPIk6dbBrFzw92f4CFEZRVG5uLsODRlJS0vee2n1Lq2zASgqMz+fzeDwtLS3xt4mrl+bt7e3t7U1rVDJr0WLPiRMu/fvjxQtUq8ZuLAQhUlwcpk2DpqaLqenEESP8ly8XfltoKMzM0L27VmEhlZCA6Gikp2PwYHh5oVEj7NmDJk2YjZsgCFH7CFWHkdHHnj1hYoKWLdkOhSCEuX8ffn4YOxazZuHTp53v3iUsXz5N6J07dyI1FV+/on9/fq9e/IMHAcDGBhcu4NYt6OjA2Rmenvh/2SWCIBgiPBFaWVlxRGA4PgAzZkBLC2/eIDiY+c4JQqTXrzFyJNq2RatWuHMHffqAw4G2trbQm589w9SpmDEDmppwc6NateK/fYtnz76/2qQJrl3D5cvIyoKtLTp2xOzZK6dMmU7LbiWCIMQTPjQ6e/bs3Nzckr9++fIlISHh3r17EydOZCqwH5yc4OQER0eEhaFbN/TuzXwIBFFGTg5WrsQ//2DYMDx+DBMTCfcXFSEwEAsWIDkZgYEAwOWie3ccOoSpU3/c1qIF7t/HyZMYODDmzJkzgFlOzoSNG/9R4ldCEISoRCj0oImxY8fev39fyfEIN306Ro3C0KEYMIBMFhJsKirC9u2YPx+tWiE1FVKezjljBqytMWwYatbE5cvfL/bujTlzyiRCAV9fbN9e2KtXOkUVPnjgRWf0BEEII8Mc4fjx4/fv31+66BpjfHxgYYFffoGNDZksJFjw9OnTzp0H9e+/qlEj7N+PmBhERUmbBWNjceAANm9GbCzq1oWd3ffrbdviv/+Ezwj26NHjypWIWbPmXro0ZulS2r4KgiCEkiERamhoUBT18eNH5UUjxtSpWL4cV67gzRuEhLASAlF5zZ695tQp/wMHYhcufHn6NJydpX3ju3cYOhQ7dsDCArt3fx8XFdDUhJ8fDh0S/kYPD4/Fi33Wr8fs2YiMVDR+giDEkCoR8ni8J0+ejBs3TkdHp169esqOSagePZCXh5s3ceAAduwQ+euDIGhXWIinT3319P5o2lS3d28r6d/I5yMoCKNGwccHubk4dQp9+pS5oXdvCNaOijJyJGbOxJAhiI6WK3SCIKQgfI7Qysrq3bt35S7q6ur+/fffYk6TVyoOB1OmYPlyxMVh2DCQnYUEM3g8DB6MmjW7XLnSWUtLtlXTf/6JggLMnAkABw+iVStYWpa5oWNHDBqEN28peR3bAAAgAElEQVRQo4bIRhYtwqdP6NUL586hTRs5vgKCICSQatWohoZGjRo1fHx8aoj5eVW+wEDMn4/UVGzahHPn0LLlj9XnBKEMFIVRo/DxI6KjIWsWvHoVq1cjJQUaGgCwe7eQIX0tLXTpgiNHEBoqrqn165GZiQ4dcO0anJxkioIgCMlkWDXKOi0tTJiAP//Evn1ISkLNmhg6FOHhbIdFVFxTpuDePcTFfS+QJr2sLPTrh02bULs2ALx/j2vXcPiwkDt790ZYmIRECGDvXrRvDw8P3L0Le3vZgiEIQjxVryxTzogRiI/HkyewsMDBg4iIIJOFhLLMnIlz53DiBOSonRkaiq5d0aPH97/u3o3u3SG0emvnzrhxA+/fS27z7Fk0aQIXF/w0a0EQhEJ+PBEmJiZu375d4hu2bt2qzHgk0NfHqFFYuRKbNsHX9/tk4cuXqFqVxaCICmjJEkRH48IFmJnJ/N7wcNy+jZSUH1d278ayZcJv1tVFp044ehTDh0tuOSkJDRqgYUOkpcHYWObACIIQ6kcifPfuXUJCAouhSGncONSrh3nzUKPG98lCT08yWUjQad06bN+O+HhYWMj2xo8fPx46lDRnTvsLF36sKnvyBG/fijtrqXdvhIdLlQi5XNy5A3t7NGiAtDSIqOZGEIRsfiRCf39/f3U468jcHAMHYs0aCAr8CyYLAwJeL1mSb08mTwiF7dyJP//EhQviVnKK0qyZ74sXHRs02NOgwZ7SDfbv/33JjFBdu2L4cGRmwtxcchfa2njwALa2aNgQjx+D2YPRCKJiUssfoylTEB6OL18AwMIC48YdiYrq7uDgu3v3brZDI9Tb4cOYMQOxsbC1left374B0KxenV9yhaKwZ0+ZffQ/09dH+/Y4dkzaXgwN8eABPnxAtWpre/Ua+PXrV3liJQji/yQkwk+fPv1XFjNhiVezJrp1w4YN3/+qrZ0KOFJUi1u3brEaF6He4uIwejROnECDBnK20K7d8eBgtyNHfixlvnIFOjpwdZXwRok768uxtMT06ds/fow7csQiJOQ3uYIlCOI74dsn8vLyJk+e/O+//2ZlZZV7ScyJ9kyaMQNt22LCBOjpYd68eQ8ehB465N65MzmoiZBTYiKCgnDokOSkJUpREeLiql6/3sPA4MfF3bsxcKDk9/r5ITQUWVmSD7Io4ePjyOXe5fM/vXw5RPZgCYL4QfgT4bhx43bu3Dlx4kRvb29/f/+wsLDOnTsbGhqGhYUxHJ8o9evDwwMREQCgqal54MBmN7cRs2ZpsRwWoYZSU1NXrDjm71+8bx+8FDjs4dw51Kv3feOgQFER9u9H//6S32tkhDZtcPy4DN21bNnyxYvErVu3paaOXLJE5mgJgighPBHu379/xYoV8+bNs7W1rV+//pgxY2JiYkJDQ/fv389wfGLMnIkVK1BycOmqVUhJwadPrMZEqJsnT5507Dhu+vToHj02+/go1NTBg+UPy4yNhaOjtNONso6OArC2th46tMGWLZgzB3v2SL6fIAihhCTCDx8+fP361cfHB4CWllZJrbUJEyYkJCSkp6czGZ8YzZujVi2UpGZvb1haCjndjSDE0NTUycnJ19PLdHcXtt1dajwejh5Fr15lLpY7bkK8Hj1w9iyys2XuOjgYM2ciKAinT8v8XoIgIDQRGhgYACgoKABQvXr1tLQ0wXUdHR0An1TpmWv6dCxbhpJZy7FjsXcvqwER6ubYsdpubofOnZs+YsQQRdo5fx42Nj/OGgSQnY3Y2PLHTYhhaoqWLRETI0/vixZhwAB07QqWTs4mCPUmJBHq6+vb2dndvn0bQOvWrePi4g4fPpyenj59+nQdHZ26desyHqRIvr7Q0EBs7Pe/zpwJHg+bNrEaE6E+XrzAokXYvt3Gw6O5gk39PC568OD3A6WlJ8foaInISHh4oHlzqaq1EQRRmvA5wtDQ0GfPngHo0KFD+/bt/f3969Sps23btsWLFxurWGUnwdlMAlwuunYFOdGbkNLYsZg0Sf7NEiX4fBw9Wj4RyjQuKtCzJ06dQl6enGHEx8PaGk2aID9fzhYIonISvn1i8uTJJX8+duxYcnJyWlqai4tLw4YNmQpMWgEBmDcPly+jZUsACAtDrVq4fh1ubmxHRqi2vXuRloYDB2hoKiEBVlYoPVby9i1SU9Gtm2ztWFqiaVOcOlV+rlFKggJsNjZwccHDh/K0QBCVk/AnQsEE4fc7uFxPT88BAwaoYBYEoKGBiROxYsX3v1pbo2FDlMrjBCFEVhamTsXGjdCiY8fNwYP49dcyV/79Fz17Qo5DrBUZHQWgrY07d/DmjbjSpgRBlCM8ETZp0mTQoEFnzpxRke3z4nXokB4d7Vmzpqeg8M2yZbh4EaXOFSaI8iZNgr//91EEBfH5OHQI5cr0yjEuKtC7N06cUGhs09IS167h8mUEBMjfCEFUKsITYe/evWNjYzt27FivXr3Fixe/ePGC4bBkcvNmMofjnZHRPikpCUC3bjA2xu+/sx0Woari43HmDBYvpqe1y5dhaYn69X9cefgQGRmQb1ditWpo0gRnzyoUUr16iIvDwYOYMUOhdgiikhCeCJcuXfr27dvTp0+7urouWrTI1ta2VatWmzdvzsnJYTg+aXTr1q1Dh3xT0+we/z8Iddgw7NjBblCEiioowKhRWLsWRkb0NPjzetFduzBggLjjJsRTcHRUwNsbERFYsQIbNvDz5F5+QxCVg8ii2xoaGh06dIiKinr9+vXq1avz8/NHjhxpbW3NZHBS0tfXP3ZsdWHhmuLi70UeFy9GXh4iI9mNi1BFCxbAyenH2fEKoigcOlRmglCa4ybE+/VXHD2KwkJFYxs4EOPHvwgNbWFk5BQWtl7R5gii4pJ8DJOlpWVQUNDw4cNtbGxU9sAXLS00a4bLl7//VVsb7drRNvZFVBh37mDLFvz9N20NJidDTw+ll5ElJkJXFy4u8rdpbY169XDhgsLBAU2anAVq8Pm/bdiQSENzBFFBiUuEPB7v5MmT/fr1q169+ujRo21sbLZv385YZLJq0wYJCT/+GhaGJ0/w6BF7AREqhs/HyJFYtkyeE3dFOXiwzJqUzMzMZcti+/X7pmCztIyOAhg8eLC3N9fC4uyjR2vGjaOhQYKokIQnwgcPHsyfP9/e3r5r166JiYnjx49//PhxfHz8kCFDmA1PBt7eZRJhvXqws8OkSewFRKiYsDBoayMkhM42Dx0qM0Ho7t755MkLZ88OU7DZX3/F4cM/CsrLjcvlXrhw6OPH4/v3W27YIOf6HYKo8IRvqG/btm1WVpafn19YWJivr6+G3PP+DGrRAnfuIDcXJafBLVqEQYNQWAhtbVYjI1TAixdYvBgJCeBwaGvz2jVoaMDJ6ceVb9/4HI6+pqaiGczCIjs/f7iTU97p0xtomZj390dKClq1gqMjbtyAvkIFxgmiohH+RBgWFvb+/fuoqCg/Pz+1yIIAdHXh7IykpB9X+veHri4WLmQvJkJljB2LCRNoqKZW2sGD5Wtq//ZbdIcOjY4eDRfxDmldunSpsLD2w4edDx8+pmBTJVxdkZaGr19RqxZevqSrVYKoCIQnwj59+hjRtbqcQeVGRwEEBWHjRpaiIVTGvn347z/6j+g6cKD8xolHj6r369fb0NBQwZZbtmxpZ/dQX/9wt25dFGyqtKpVkZ6OWrXg4PBjZRlBEMKHRgE8ePAgOjr65cuXhWXXcW9S4cMd2rQpX3F7xQps3ozoaPj5sRQTwbasLEyZgn37aB4hv3kTPB5cXctcTE6mZw+7iYnJxYvH6taFjQ0NrZWmrY2bNxEUhDZtEB6OwYNpbp8g1JHwRLhmzZrJkydzOBwrKytt9Zlha9kS164hPx+6ut+vGBigZUvMnk0SYeU1aRJ69aKnmlppBw6gT58yM46ZmXj3rkyJGUVYWEBfHy9e0J8LAURGokEDhIQgNRVr19LfPkGoFyGJkKKoefPmdevWbdu2bebm5szHJDdDQzRsiKtX0br1j4tr18LNDa9eoWZN9iIj2PDgwYNFiw6cPev/+HEj2hs/eLB89aKrV+HuLn9BmZ+5uODWLaUkQgAzZ6JePfTvj9u36dmzSBDqS8gc4fv377OysmbOnKleWVDg52lCFxfUrEn2UVRG3boN/fdfZw5nqIkJzS0L1ic3a1bmYkoKmit6vm8Zzs64dYvOBsv59VekpODqVb6u7khLS9fk5GQldkYQKkxIIjQ3Nzc2NlbZIjLildtWLzBtGo4eBY/HRkAEe7S17XV0djg729HesmC9aLmdGGqXCAG4uuLvvw8WFHz89GlSQMAKyW8giIpISCLU0tKaPXv24sWLVbPEtnitWyMpCUVFZS6OGQMNDaxcyVJMBBs+f0ZmZmR09F+nTu2hvfGf14sCuHqV/kR48yadDQrVs6e3nt4tDuePly/H1aiBEyeU3iNBqBrhi2Xev3//6NEjBwcHLy8vCwuL0i+p8qpRACYmsLPD9evw8ChzvU8frF6N6dNZCotg3Lx56NsXHTvSP8P26BGystCiRZmLz5+DywW9RekdHJCRga9fYWxMZ7PlVK1aNS/vKYCcHIweje7dUaMGdu9GmzZK7JQgVIrwRJiQkKCvrw/gxo0bzMZDA8E0YblE+NdfqFoV8fHw9mYpLIJBjx5h927cv6+UxqOi0Ls3uGUHU5KTy6dGxWlooFEj3LkDLy+aWxbK0BCRkVi5EsHB8PFB/fo4dIi2RbAEocqEb6hPTk5+JgLD8clB6DShhQXc3DBtGhsBEYybNg0zZqBaNaU0/vMBhFDCuKiAYOEok6pVw8mTuH8f2tpo1Ai+vsjMZDQAgmCe5GOY1E6bNrh0ScjSmFWrcPUqPnxgIyaCQfHxuHMHSjps4ckTvHsnZFdicrJSEiED62WEql8fN2/ixAncuYOqVREUhJs3H+7bt4+FUAhC+UQmwhcvXsyePdvf39/X11dwZc+ePTExMUwFJj9LS9Sogdu3y19v0wZVq5KHwgqOz8eUKVi2DDo6SmlfsEym3GZBHg+3bsHdnf7umFkvI0rnznj5Eps349ixt66uv/bvv7N589FquISOICQQOTTq5OS0adOm9PT02/9PKa9evZo4cSKDscnv592EAmPG8Hftit6/fz/jEREMiYyEpmb5Wtg0EjoueucOatVSypIWJyfcu8fyzp+QEMTFveBwOBTlmJpaxcgIenqws0O3bli6tMwnzuLi4qdPn7IXKSHZI3JGqzDCE+Ho0aPd3d2fPXu2atWqkotdu3Z99OjRu3fvmIpNfkKnCQGkpY0pLt4TELDk2DHaivoTquPbN8ydi5Ur6TxrqbTnz/HiRZm6RQIpKeUXZ9HFyAhWVnjyRCmNS8/Dw2Pr1okjRxZ9/jwlKwvh4WjbFm/eYMUKuLhAQwOWlnByeqGtXdvVta+TUyuWwyVEqFKlSbNmQRwOfSdTVxRCEmFmZuaNGzcWLlxobGzMKfUbxcbGBsCbN2+Yi05e3t64eBEUVf66pqYGkE1RXC63Ak6OEitWwNNTiWss9++Hvz80f1pqnZJSvsoMjZhfLyNUSEjIxo1hxsbGxsYYMADh4bh+HZmZ4PORkoLx4wFcAayA6ffvm6tnNY4K7uNHfPqUDSwHLI4cUYPnGSYJyQeC4yZ+Pkrm06dPADR//jWgeqpXh4mJkNXzGzeuXb68C5e7E+jGRlyEEr17h7Aw/PGHErsQOi4KZT4Rgr31MtJr2hRz5uD27QA7Ow0OZ2fz5osbNsTmzeDz2Y6MAADw+di5E05OcHCYxeHMMjDoPXFitb59ybGUPwhJhFWrVq1atWp0dDSA0k+Ee/bsMTQ0dHR0ZC46BQidJuRwONOm/ebt3Yj2o+kI1s2ahZAQ2Nsrq/1Xr/DsGdq2LX89NxdpaWjSRFn9qn4iLPH06dXFi+Pq13c6dgwREWjRAikpbMdU6aWmwssLW7ciNhYmJsN3707R0Zl/9iwaNoSrK+bPR9lz9iopIYmQy+VOnjx54cKFCxYsePLkCZ/Pv3379qxZs+bNmzdu3Dh1OZVJ1DQhgK1b8egRUlOZDYhQptu3cfw4Zs5UVvsXLya6urbW1w/hcMovXLl6FU5O0NJSVtfsLhyV1cCBxUePom5dJCZizBj06IFBg8ieJXZkZmL8ePTogdGjER+P4mJ8+gRf3+LevbF/P+bPR1ISrl1DkyaIjWU7VtZRwvD5/GnTpmmV+uHmcDhDhw4tLCwUer+SbN26NTg4WL73pqVR1auLfNXZmWrRQs6omPH161e2ui4uLs7NzVVGy8r7ojp1ov75R0ltUxRF9ev3G3DexMQ3PT293EvLl1MTJsjWWmFhYX5+vvT3m5lR797J1gUr+Hx+dnb2r79SmzZ9v/L5MzVuHFWlCrV6NVVczGZsSvreY/HnVAwej9qxg7KyosaNo758+X5xyBDqzz+pr1+/JiVR9vYUn//9+rFjlJ0d1a0blZbGUrgURSntX5LH40mTtoSvGeFwOMuXL09LS9u1a9eff/65devWR48ebd26VUt5n3vpVqcOtLVFLrfbsAHJyXjxgtmYCOU4eRLp6RgxQoldjBoVwuX+4efnWKtWrXIvKammTGlOTkL2xaqsoUMRHv79z6amWLMGCQk4cQLNmiExkdXIKjo+nw8gNRUtWyI8HHFxWLMGgjPIPn/G0aMYMgQAPDxgaPhjwMzPDw8eoEMHeHhg/nzk57MUPavErXyxtrYODAxkLBTaCUZHHRyEvOTpCRsbhIbi+HHGwyJoxeNh+nSsXKnEwUkA7965+fqejowU8lJyMpYtU2LX+P/C0Q4dlNsLXTp1wqhRuHkTLi7fr9Svj7g4REcjMBCtW/MaNtxuYaE5bNggsnibRklJyf7+Y/PydPT0jq5caT5gQJlNROHh6NEDVaogOxsAgoMRHv6j8LK2NsaPx6+/4vff0bgxRo26/eXLiaFD+9na2rLwlbBBeCJMT0/nCdvEa2ZmZmZmpuSQaCNIhEOHCn91xQr066f00v6Esm3eDEtLdFPyKuCYGHTpIuR6Rgby8mBH/4mHZTg74/x55XZBIy4XQ4Zg+3asWVPmup8f2rbFwIGHZs1K1tXNt7au0rVrV5ZirIBOn76YkTFQW/tqVNSDTp3KbCGiKGzejF27flwJCsL8+fj8GaV/nVtbY+dOnDmDLl2GFBfPPHp02J07Z5kKn2XCP5F5eHjYC2Nubm5qahocHCzYSqHi2rQR9+vj119hZgY1KZVDCJedjUWLlH7SJEXh1Cl07izkpeRkeHgoa/9+CfVaLwNg2DD8+6+QQTZDQ8yebWtqerOg4EFRUW02QquwLC2HmJo+/O03hw4dPMu9FBsLQ8MyA/jm5ujcGf/+K6SdDh3g5lZDR2dvfr6wwbQKSngiXLp0qZmZWceOHdeuXbt3796//vqrefPmNWvW3LRpU0hISFRUlL+/P8OByqFePfB4SE8XecOMGdi1C8XFDMZE0GrJEnTujKZNldtLaipMTYU/9il1K32JRo3w9CkKCpTeEV1q1oSbGw4fFvKSu7v706enNm8+M316E7Lvni7fvmHZMstjx9avWjX35wHn9esxdmz5twwdClFny165cuzatRWamhuFzgVUTEKX0Pj5+U2aNKn0FT6f7+/vP3bsWIqiBKW3U1NTaVnVI4Yiq0YFAgKonTvF3WBgQM2apUgPykJWjUr08iVlYUG9eEFjk8ItXEhNniz8pQ4dqJMnZW5Q1lWjFEU1bkxdvy5zRwwTrBoV/DkqimrfXtzNI0dSffsyEZVAxV41OmcOFRgo/KX0dMrSkir5gS4JmM+n7O0pMb/Fb9+mqlal7t2jN1KRVG7VaHZ29vHjx0NCQkpf5HA4wcHBe/fuBdC5c2cLC4sHDx4wk6oVIWY3ocCIEVi7lqloCFpNm4YxY/DTKk76iZogpCikpjLxRAiVKbQmvR49cPcuxJTgXrMG//2HsDAGY6qgXr3Chg1YskT4qxs2YPBg6OuXv87hYMiQH+t7f9akCf76C/7+qAwP7kISYV5eHkVRP88Cfvz4MTc3V/BnExMTjXJH0aikNm0QHy/uhuXLUVCADRuYCoigw9Wr16pVa37gQPuRI7OU3VdmJu7fRythdaQfPYK5OSwtlR0CoFb1ZQS0tTFwICIiRN6go4ODB7FkCdlToagpUzB2LGoLm3ItKEBEhMidRcHB2LsXeXkiWw4MROvWyt2YpCKEJMJq1ao5OjpOnDjxRal9dvfu3Zs3b16bNm0AZGdnv3r16ucNVSqoUSNkZUFMnXAtLfTujYULGYyJUNjx42ffvx+qrV3z+fOf6snSLSYG7doJP91QSYfxCqV2iRDA8OHYtk3cHHzt2ti6FQMGkNIz8rt8GZcvY/Jk4a9GRcHFBfXqCX/V2hqenjh4UFz7YWF4+hTr1ikap4oTvlhm27ZtT58+rVu3rru7e+fOnZ2dnZs0aVJcXLx27VoA165da9u2rbsyziGlG4cDLy9cvCjunnXr8P49yLlMaqRatcEmJklBQdbNlZ+IRI2LgpGt9CVcXHDjhpADVVSZoyNsbSWU7+raFYGB6N+f5TMX1RSfjwkT8OefMDAQfsOGDQgNFddC6eoHQunq4uBBLF5cwR/chSfCli1b3rt3b/r06TVr1vzy5YuDg8Mff/xx9+5dBwcHAG3bto2NjdVR0hHgdJM4TWhmhnbtQMpwq4tv3/Dnn1bR0ds3bFii7PF5Ph+nT+OXX4S/yuQTYZUq0NPDq1cMdUcXib9nASxeDC4XixYxElDFEhEBTU0EBAh/9eZNvH4NX19xLfj54fFjPHwo7h4bG2zZggED8PGj/KGqOJGVZWrWrLmoQnxvtmmDrVsl3LN1K2xtcfUqQwsfCEUsWwYvLyGn4ypDSgqsrETOvjx4AFdXJsIQEOwmVIcZiR8CAjB1Kt68QQ3RZ8Fyudi1C+7u8PAQ+fBN/Cw7G3Pn4tAhkdtY163DqFEQ/1lRUxNBQdixA0uXirutWzdcuoR+/XDqlIQG1VTFL3Hk7Iw3b/D+vbh7bGzg6ooxY5iKiZDXq1dYv17kAjnaxcSI/EB94wbq14eeHkORQA0XjgIwMMCvv0LidrSqVREVheBgPH/ORFQVwx9/oFMnkWMSX77g4EGUXfsv3NChiIhAUZHk7rhc5Z73yaIfT4QxMTHLli0bMWJEYGBgr169MjMzhb4hXvwqTNWjoYGWLZGYiF69xN22cSM8PJCWhkpTXU8tTZmCMWNgY8NQdydPYtUq4S+lpDA3Lirg7CxhXYNqGjoUAwZg2jQJ9XdatMDUqQgIQEKC8KVJRGn//YfwcHGl2LdvR9euqFZNclP16qFePZw4gZ49xd2moYHISDRrhubNhVdZUms/ngi1tLQMDAwE50vo6+sbiMBeqPKTOE0IoFkz2NpKmFgm2CVYIDdlCkPdvX+Pp0/hWb5e1XfM1JQpTe0KrQk0b17mrAMxJk1C7drM/f9Va1OmYPJkVK8u/FWKwqZNGD1a2takmcoFUK0a9u3DkCEV8cFdGZv56aJ4ZRmBK1coV1fJtx05Qmlo/Di+i12kskw5PB7l7k7t3UtjOBJERFB9+oh81cFB/qIbclSWoSiquJgyMKCysuTslAGlK8uUtno1FRQkVQtfv1L160uoBiWfilRZ5tw5ytaW+vZN5A1xcZSzs/CXhAaclydDkaY//6Q8PKiCAqlulp7KVZapeNzd8d9/+PJFwm09esDcHOPGMRITIaNt26Ctjb59metRzMaJzEy8e4f69ZkLBoCGBho2xN27jHZKi6AgREfj82fJdxoZISoKU6bgvtI3iKorHg8TJmDlSujqirxn/Xr89psMberpoW9f7Nwp1c1TpqBmzYq2zF5kIrx582ZISIiHh0fLli0FV8LDw6OiohTvMjs7e+fOnaNGjQoODt68eXORxFlahWlqonlzqfbBzJqFf/+VPG9MMCw7G/PmYfVqpZ/zUILHw5kzIjdOpKTA3R3Mn6anjutlIPasg59VqspectiyBebmEHPqwcuXuHgRAwbI1qxgdJTPl3wnh4Pt2xEXV+ZcJ3Un/Ec5Li7Ow8Pj4sWLpqamJfVlvn37Nn36dErhPb3nzp3bu3evs7Ozj4/PP//8ExwcrGCD0pBmmhDA+PHQ1cWcOcoPiJDFwoXw9WV0Tu7yZdSpI3LRP/MrZQTUsb6MgJizDn4WGIg6dU5Xr96ia9fBfGl+N1caX79i4UIJ545t2oSgIJFb7EVp2hSmprhwQaqbjYywe3fh0KF9atb0unTpsmw9qSThiXDs2LE9evS4f//+77//XnKxQ4cOz58/f/v2rYJd+vn5nTx5cvTo0YMHD46IiIiKiiosLFSwTYkkFh0tMXJkxa8npF6ePUNEBNNl8MSMi4IkQtm1b4/cXFy/Lu39Wlr78vJWJiS8/FiBd3HLbv58dO8u7tyxoiJERGDUKHkal3LJjICBQZq2dvHr1zPDw4WdtqVuhCTCDx8+PH78eNq0aVpaWpxSQ1GC4qKKJ8LSx2VlZGSYmppqa2sr2KZEHh64dw/Z2ZLvXLoUhYUkF6qQiRMxfbrIBXJKcvKkuETIZHG10pyccPeuVONXqkbiWQflLFjwW+PGS/j8Nk+fVlVmXOrk6VNERmL+fHH3HDiABg3g6ChP+wMH4uRJacvHODg4BATUrVJl4/v3Q+XpTMUIqSwjGPz8+XTHjIwMAHrSbSH+8OHDz2MahoaGpTdgfP36deLEiQsWLBDTTmJiYrlDgLt27RogqqaQWC4ueufPF7ZrJ7mmYa9eugsWaAwenCtHL3TJzc3lMDYhVhaPxysoKFDGkJQcX9SFCxp37+ps356Xk0N7OCJlZHBevtRv1ChXaKfp6VwNDT0TE+GvSqOoqIjP58sxO87lwtJS//bt/Lp1VTEZUhSVJ/osg759OS1a6M+bl/vzkUA/q1fP4ZZWTl0AACAASURBVMqVqNhYzX79qIsX8ywsFJ2RUdIPFGM/p0eORC9ZkjV+fKChIUfMN15YmN7YsUU5OSIrnYsJWEMDnTvrbt/OGz1aqu/M1avnLVuGjh31V68uGDZM0YUVSvqX5PP50lQDFZIIq1SpUrt27X379rm5uZWObPPmzebm5vVEVTIvq0OHDj+PaUyaNGny/8uk5+Xl+fn5tWvXbpTYx/jatWv369ev9BVHR0d9aX6SfuLjw0lJ0e3WTfJP1KZNsLDgxMbqi5mRVjYejyffl0lL1xoaGsroXdYvqrgYv//OXb2aMjNj9J/i/HnOL7/AyEh4p3fucFq0gCL/PoJEKF+1XhcXzqNHek5Oqlh+W7ASXdS/TN268PREbKzBwIHSBu/vj+RkhITox8ZSClb2UtIPFDM/pykpKaNH78rNddDX36WvP1LoPe/fvx84cM7t29V79JijqytyjE18wMOHY8IEzcmTtaQMTF8f+/bBy0unWTMtUTtupaSkf0k+ny/NuhYhiZDD4cyePXvkyJFfvnyxs7MrKiqKjY2NioqKiIhYvny5pqbI8qSl3RI7lfHt27cePXrY2dmtW7dO/KeAWrVq9aVpyby3NxYuBJcr+UOHiQnateOHhkbb2dVyc3OjpXdZcbncnx/KmUFRlJJ6l7XZjRtRrRq6d+cAjD4cx8aiVy+R3yqCcVFpvpFEEfwjyPcv7OKC27cREMDOaIF4Er9zhg3DmjUYNEiG4JcuRfv2WLaMo+ASNhX5lpZP9erVi4reGBhk16/fWVR3W7ZEnjvnrKOTcuNGipeXl6imxAfcti3y85GaypV+YVq9eggPR2Ag99o1hc7mVN6/JE+ak01EbTBcs2aNiYlJyW26urpz587l8XiK73AsKCjo1q1bnz59ioqKxN9J14Z6gdxcysiIysuT6uZWrQKAwVxuncePH9MVgEzIhvpPn6hq1ag7d5QRiDiFhZSZGfXuncgbWrWizp1TsAt5NtQLHDlCde2qUO/KI2pDfYmiIqpGDerBA9mazcigatakTp1SKDa13lD/8iVlZvbl6dPXYu5JTk7R1Gxma9smMzNTzG0SA16yhBo5UuYIJ02ifH0pRfKDim6oHzdu3OvXr8+cObN3794TJ068fv16wYIFtGTsqKio48ePJycnOzo62tvb29vbK74ARxr6+mjcGCkpUt3M4+UA2ny+YX5+vpLjIoSbMwcBAWjcmOl+L11CvXqoKmKJBo+HW7fELdtTNjUttCYgOOtAzLH1QlWrhshIDBmifqdQ0WXrVgQGmtjbiz7CA6hTp5mhYdKTJ/FmZmaK9BUcjKgoyDr/vXw5srOxbJkiPbNJ3DingYFB+/btae+yR48ez549K32lSpUqtPcilGAThbe35DtPn943Y8aiDRv2nDzZpEkT5UdGlHX/Pg4cYKe8iPiNE3fuoHZtGBszGFBZNjbIy8PHjwoNQ7EoJAQtW6b4+2s2by7DpIOPD8aMQZ8+iI+H8teYqxY+HxEROHJEwm3x8WjThqv4GUlWVmjRInPOnCuLF/tIX1xaUxNRUXB3R7Nm6NhR0RiYx8IslJGRkV1ZUs47Kk7KbfUADAwMwsKWzZrVZO5cqTZdEPSaOBHz5sHCgoWuT54Ud5Ypk4fxCsXhoEkTdd1NCOD587isrHnt20+9dOmSTG/8/XdYWaHUxubK4uRJVK8OFxcJt0n5EV8at275hYUld+kSJNO7rKwQGYnBg/H6NT1hMKlS1Bot0aoVkpPz8vJEri0uZ8ECmJszWt+SoChqz57sV68wYgQLvb98iQ8fxI18srWDsDT13VYPgKIoLS1uURGXkrFGFYeDbdtw5AgOHFBSaCpqyxapfhZoTIQGBuDzNQsLZd6i07YtQkPRpw+UXyKFZpUrEV6/nlBQ4GNr6yHYEymN/fsRF4eLF5UaF/FdYWFh/fqtBw3q3LHjTqaGCco4fhydO4srIsr6EyHUPBH+8ssvGzfOsrL6o3Xr1rK+18wMhw4hNBQPHyojNFX09i0SEyV/Fv/0Ca9eSX5qlNLly0ednJymTJGnluisWahaFbNm0RMJYypXIrx27SaP1y07u265SUoxWrVC+/bkoZAhWVlZL1/yebzh799Lt6iJbuInCHNykJYG1ueM1bT0domBA1t+/dr8/Xt53uvsjAUL0LcvRG/cr1C2bEG/fpILh8bHw8sLik8QClhaWvr797x2zVCO9wpKch88qGaHSFeuRDhyZHCXLnx7+zYlR2pI4/BhfP6MuXOVFxfx3aNHVXR0Jg4c+GDlypnM915QgIQEcVP9167BxQVa0u42VpbGjfHkifqNPpXgcuHpKdVpMEKNHg03N3ZGzhnG52PbNgyVooQZjeOiAq1ayT8MZmaGffswahQePaIzJKWqXInQyMho/fr579//JlMtHwMDLFmCpUvx6ZPyQiOQlYWgIERG9omMXF5D1LkPyhQfjyZNxK3QYavWdjk6OrC1xYMHbMehgFatIONamTLWrcPNmzJULlVTsbGwsoKrq+Q7L1ygORG2aIFbt+R/7G7WDPPnq9ODe+VKhABq14aWFtLSZHvXpEmoUUPcGWCE4saOha8vunVjLQDx46JQmUQINd9NCMUeOAAYGODQIcycKcNxFupo82YMHy75tsxMPH8OektgCXZdX70qfwu//QZnZ7U557zSJUIALVrgyhWZ33X0KC5dQlycEgIigEOHcOUKli9nMwbxGyegYolQracJmzfHgwcKbUyqVw9hYfD3L4iJSSwoKKAvNFXx9i0SEqRanZCQgJYtQfviMgWf2gFs2IArV7BkycM0WZ88GFcZE6GnJ5KSZH6Xiwu6dUP//lD4ZGKivNevERqK3bthKM/0PD3S0pCdDWdnkTe8fYu8PNjaMhiTaOq+XkZHBy4uSE5WqJG+fZGf79+z5/42bSrgWM3WrejXD0ZGku+kfYJQQPFEaGCASZPOzJkT6ubW/6Zqj2BU0kQoxxMhgP37kZ+PqVPpDqhyoygMH46xY1l+2Dp+HL6+EDN3nJwMDw9xNzDJxUW9h0ZBx+9ZAObmuUVFDo8e5RRLuzdYPQiWyQwbJtXN8fHw8aE/hlatcOUKpClYLYaJSZaOjlVWlsXVq19pikspKmMibNoUDx/KM4urrY3Vq7F6NRipjVpZrF6NrCzMmMFyGOInCK9fvzFyZKe3byfLug1cSapWhba2etfepCURnj377/r1+p6ee3v1QkWqCnzqFKpWlWraLysLT5/SPEEoYGmJGjVw+7ZCjfTu7f/vv/0XLpw4Z04budcJM6AyJkIdHTRsiNRUed47fDjs7NCjB90xVVb372PJEkRE0LYFSj7fviExER06iLxhxYrw9+9n/vffg5cvXzIYlzgVYL3M1auQ/XDiMqpXrz5qVPCxY9X19ODrK3OpaJW1ZYtUy2QAJCTA01NZ9Vdbt1b0wwqHw+nRw2/27A5796J3b5w+TVNkdKuMiRAKjI4CiI5GaioOH6Y1oEqpoAADBmDFCjg4sBzJuXNo2hSljh0rb9So/lzu7BYtqllbWzMYlzjqvl7GxAR16uDGDRqa0tLCv/+iTh34+uKrSo/ASSUjAxcuICBAqpuVNEEoQMtTu4CPDw4cQGCg5OrhrKi8iVCO9TICjo4ICMDgwYqOnhOzZ8PWFkOGsB2HFBsnzM29HBwuxcZu12D30bUUdU+EoOOBo4SGBsLD4eqKdu3Ufr9veDgCAqRaJgPlJ0IpTymQsrWYGISGqmK12EqaCFu0wOXL8r89MhIUhdBQ+gKqfC5exO7d2LSJ7TgAPp8fE/NNfCK8cgWenkwFJB11XzgKWh84AHA4WL0abdvC21uNZ/EpCtu3Szsu+vUrHj2Cu7uygrG1hZYWpK5HKVnTpjhzBhMnYvt22tqkRSVNhDY20NTE8+dyvl1DA5s3Y+tW/PcfnVFVHoIiMlu3ijz/ljGZmZk2Np7p6T6vX58Rc1tSEjw8GAtKKo6OeP1avWfFWrfGxYt07kficLBiBYKC0K6duq4kOnUKZmbSLn65dAnNm0NHR4nxeHnRfORAw4Y4cwbz5mHtWjqbVVAlTYQAPDzknyYE0L8/GjZEz570BVSZhIaie3cJu9eZkZ6enplZnc8PuHRJXBWNpCSVeyLU0ECDBrh7l+04FGBtDUND+itSTp+OkBC0bk3nowxjpF8mAyWPiwrQ+9Qu4OiIixcRFoZVq2huWW6VNxEqMk0ocOIE7t/Hzp00BVRpHDyIGzdYLiJTwtXVtVq1Tl26vJ80aZSoe758watXaNSIybgky8/P//Bh6MCBvV6p6bMPAFqnCUubOhXTp8PbG/fu0d+48mRk4Px59Osn7f20lxj9mZL+B9nYICEB27ezv29KoFInQkWeCAHUro2hQzF6tBqfA8A8QRGZHTugp8d2KACA9++RmRl66NAyMzMzUfckJcHdnf4SVgq6cuXK+/d6z5513b17P9uxyE/BoqNijBqF5cvRqZOiO+GYtG0b+vSBsbFUN+fk4P59pZehaNwYHz5AvjOzxKteHefOITYW06bR37isKm8ibNoU9+8rWhx940ZoaHxo0GDm0aNHaYqrIvv4MXPQIEyciGbN2A7l/44dQ5cuEmZZVHBcFEDTpk3r1EnX09vVs6cKDDHLSxkjbyUCA/H33+jYEcnJyMzMVFY3NBEsk5H+eKnERLi7Q1dXmTEBXC5atFDW/6OqVXHhAi5exOjR+PiRzf9BlTcR6uqiYUNFq9dzONDS6vnff1V69Rqfo9brFpSvZcse9vb+CQldVapG3dGjkssjJCWhRQtGopGFsbHx5cvRmpoXHBwc2Y5Ffg0aIDtbiQtb+vbFunX5np42tradg4J+U1Y3dDh9GgYGaNpU2vsZmCAUUOqHFVNTnDqFAweG29t31tOrU8jS8FrlTYSgY3QUgJWVNpBAUaa3bytz8Zb6u3HjPkWt4PHuqMxOPOTmIiEBnTuLu4eiVOjQiXJMTWFmJvOZYiqFw0HLlvIf0isNZ+eXgAFFzTl7lo7d+0qzZQtGj5bhfgYmCAUEi3uVx9gYGhrXKWpBQYHJ27fvlNiTaJU9ESq4XgbA7dtnIiJ6du161sdHS8Fq+hUYnw97+7/MzGbOnSvLz7qSxcbC0xOmpuLuefgQ5uaoVo2pmGTk7KxOc2BCKfWBA4CDg8PQob+YmkY3aRKhxG4U8+4dzp5F//7S3p+Xh7t3GRqoaNYMDx8qdGaWRIcOhdnYrK1RY9z8+bX4fCV2JEqlToTyHUxYjoaGxuDBg48ft/DzQ+vWip4sU1GtWwdjY7+0tNPz5//Odiw/SDMuqoJb6UtzclL7bfVKWpdY2pYtfz95svnly3qqWd8LwLZt+PVXaZfJAEhMhIsLQyvOaDkzS7yWLVveuRPz+PHQ9HQMHw7mc2GlToR16oDDQXo6Pa0dPIj27dG6tULHOldIL15g4UKEh4OrSt9uPB5iYuDnJ+E2FdxKX5qTk9o/Ebq54dkzfPmi3F50dLBhA8aMUXpHcvDzC547t3nt2oekf4uSjl4ShYEPKwD09XH8OJ49w8iRTB/7qkq/mdig4Lb6cmJi0K4dvLxILixj5EhMnowGDdiOo6z4eNjZoWZNCbep5pLREhVgaFRLC+7udP4YiuLtDT8/Vdm4ViIvLy8h4W5x8brz52XYBsPYShkB5e1yKUdfHydO4NEjTJjARHclKnsipGWasLTYWLRti9atcecOnc2qrx078PYtJk9mO46fSDMump2NtDQ4OTESkFzq1kVGhtofucDMAweA5csRE4Mz4qrpMU1fX9/WNrB27UVLl0r7Q/LtG27eZHQls5cXUlIY2jBtYIDoaCQlYeJEJroTIImQ/o+ip06hTRs0a0ZyIT5+xO+/IzwcWlpsh/KT6GjJFfKSk+HqqorBl9DQQMOG6l1oDQw+cBgbY8MGjBiB3FwmupNGTg7S0yekpBxr3lza4tlXrsDZGQYGSo2rDBMT2Nszd/6liQni4pCYyNwH6MqeCJs2xb17+PaN5mbj4tC6NZo1U7MKT7QLDUVwsAxboxhz48b3FCKeio+LClSA0VFPT1y/ztAR876+8PTEggVM9CWN/fvh7S3bsmSGx0UFGPuwIiDIhfHxYGbbcWVPhHp6aNBA0W31Qp0+jVatvifayun4cdy6hTlz2I5DmKNH0auX5NtUcyt9ORVgvYyhIerXR2oqQ92tXYtdu1RlIn/HDgweLNtb2EqEzAxflzA1RWwsTp3CggXIU7AGmCSVPRFCOaOjAmfOwMsLTZviwQOltK/KsrIQGoqtW5VeAko+R45IniCkKCQnq/SSUYEKkAih/F3bpVlYYNUqDB3KfpXg589x755sx7AUFCA1lYWBCsE8LsOLOS0tceYM1qyZVrVq5169himvI5II6V8vU9rZs2jZEq6ulS4XTpqEnj3RujXbcQjz/DkyMiQ/6j19CgMD1KjBSEwKcHbGnTtM/4aiHcMPHP37w96e/SNQduxA//6yHSiYlITGjaU9v55G1tYwMsLDh0z3W7UqTE0Tc3OXx8Zef/FCWb2QREjPtnoxzp2Dpyfc3Fj4HmLL+fM4exZ//MF2HCIcOYLu3SGx0puKb6UvYWoKU1P1LrQGoE0bXL7M6E7qdeuwbh2bMxcUhchImcdFGaus9jPGFveWs2/f2r59dw8atMbVFcuXg8ejvwuSCGFrCz4fyvusAeD8eTRtCmfnV927j7jI5IwzG/LyMHw4Nm5k4UOrlKTZOAGV30pfWgUYHa1SBVWqMJqWatTAggUYOlQpv1ilcfEidHVlXkrGygShAPPThALNmjXdt++fTZtap6Tg9Gk0a0b/qg6SCAHlPxQCuHgRPF7v6OgWPj7BqrN0WxlmzUKrVhIqWbPo0yfcvIn27SXfqRZLRgWcndW+0BoYX5cIYMQIGBkhLIzRTkvs2IEhQ2R7S2EhUlPh5aWUeCRi/n9QOfb2OH0aEybA1xczZtC5zJgkQkDJ04QCHA7Mzb8BZ/j8KsbGaNAAf/3F/lw9vS5cSKhe3XP9+l+XLCliOxaRjh9Hhw6Sl/Dk5uLxY7i4MBKTwirAEyHYeODgcLBixeepU3+pU6f1/fv3mez62zccPYrAQNnelZwMR0cZSpLSq3595Obi5Ut2ehfgcDBoEG7exJs3aNwYZ8/S0yxJhIAyF46Wlp6evGVLu0+fTqSmomlTzJ//ffPG8uUoLlZ67wzYufNoRsYsHZ2CvDya6rcqgZTjolevwslJtlUMLKoYibB1ayQkMN1pZuYNTc26L14MOHnyNJP9HjwIDw9Ury7buxguMVoOA2dmScnKCjt3Ys0aDB0KZ+eVjo5tT5yIU6RBkggBwN0d9+4pfT+vnp7esGHDzM3NXVywaxe+fsXp06hTBwsWQFcXrq7YvVu91/4FBY3Q0gofONDZ3t6e7ViE+/YN586ha1fJd6rRuCgABwe8fav2hdbq1gWfj+fPGe3Uy8urSxdoaCT06hXAZL9ybB8EqxOEAmxNEwrVtStu3aIeP458/DhiwIC/583DyZP49EmepkgiBAA9PTg6KmVbvXjt2iEmBnl5iImBlRVCQqCjgxYt4OYWaGbmFBgYwnRAiklMdBwx4vD69Ys5HA7bsQgXFwd3d5iZSb5TLbbSlxBUyakApRu8vJj+Paujo3Po0Lpffvn3yhUrxjp9/RrXr0s++aScoiKkpLA2QSjA1sJRUUxMOIMH+9aqFRAaOozDwYYNqFcPNWqgb1+sWYNLl+DnNyQ6WvK8F0mE3zEzOipKx46IicG3b9i4EYWFuHHjAo93Ys8eVaoNLIW9e2We82CYlOOiAJKT1SkRoqKMjrL1wBEYiN27metuxw707SvzaYJXr8LBQcI50srm5oa0NNU6ymrjxqUvXiQtXdp7/nxER+PDB5w+jc6dce8eBg26c/Lkp2vXJE/1k0T4HQPrZSTichESguvXoampBfQFBsTHsxyS9FJTkZ+v0smDx8OJE1J9DE9LA5eLWrWUHxN9KkYiZLK+TGk9eyIlBRkZDHUnx/ZBqMC4KABNTbi74/JllsMQg8tFo0YICcHmzbh40YrDuc/nSx7rI4nwuxYtVOj/blHR8zt39gYHL2vfHidPsh2NdHbvRmAgVHVMFAASE1GzJurUkXznlSto2VLp8dCrAhxVD8DZGW/e4MMHpvvV00O3bti3j4m+kpLA58uzRVUVEiFUb3RUjE2bqvTq9XT2bMlbNUki/M7ODnw+yyuDS7Oxsdm2DSEh6N4d0dFsRyMJn4+oKAwYwHYcYkk/LqpGW+lLuLhUhEJrGho0n5UtPcZGRwXbB2X9yFhURCUlFbM7QSjA+m5CKT16hE2bsHo1R0eKxd8kEf7A1k+gGJs3Y/x49OyJPXvYDkWss2dRowYcHdmOQyxpDiAUUK8lowKmpjAxYXrJpTKwNU3YoQPevMGjR8rtJT8f+/dj4EDZ3vXhwwdb25Z5eZ5Pn6YoJy4ZeHrixg2GzsySG0Vh9GgsWABra6nuJ4nwB1WYJvzZqlWYOhVBQYiMZDsU0QTjoqrszh0UFUl11vy3b7h/H25uyo+JbhVjdJStBw4uF336KP0T57FjcHOTefr54cOHnz45FBf3jY9nf/7GwAANGuDaNbbjEGvLFuTlYcQIae8nifAHdheOirFsGebOxZAhWLeO7VCEyc9HdDQCGN2FJbMjR6Q6gBBAaioaNpR5RZ8qqAAn9ALw8MCdO1Dy8XPCCUZHlTq8LN/2QS8vLzOz+p06vR0xYgj9McmOrTVNUnr3DrNnY9MmcKXObyQR/tCsGe7eVdFH/nnzsHQpxo3DmjVsh/KTI0fQrBmsmNuFJQ+ZJgjVblxUoGIsHNXTg5MTUtgYAnR3h44OkpOV1f67d0hKkvYDWWlv33KLimYeP/6XKbubJ/5PpbbV/2z8eIwYAWdnGd5CEuEPenqoVw83brAdhwjTpmHNGkyahMWL2Q6lLNUfF339Gunp0u5EVq+t9KX9r70zDWjq2hbwysAgZVbCLCqDguDEoHgZRagWZwSKVXFGr7629FVb67Va6+yzXpXigCBcKFKpCChWgQhBBKyCFBGrtco8SKLIHCA578epuYiIGU5ykri/X8k5e6+9yMDKXnsNyuEaBVLDMT7+WIohM/HxsHAhaGiIPDEtDfz9gU6Xgk5i4e4OBQWkde0YmqtXobgYtm0TbRYyhK8ht95RnE2b4ORJ2LEDduwgW5VXPH8O+fnCBqGQxcWLMGeOsP9HFNcQ2thAQwO0t5Oth8SQuOFYvhzOn4de6RSNFy99EETxZ8gGAwMwNJTHSkadnbBxI5w4IfLRBjKEryGf8TL9WbsWzpyB3bshPJxsVQAAICkJZs+W39aDOML/H6mpgd5eGD1aygpJBxoNbG2hvJxsPSTGzQ2KisipRG9hAVZWkClRAefBKS6G1lZwcxN54suXUFQEM2cSr5IkyGcSxbffgru7OK8VMoSvMW2aXNRWH5qVK+Hnn+H4cfif/yFbFUXwi758Cbdvg6+vUIMVMZW+P8rhHdXXB3Nz0v4QKSUUipc+CAC//gqennL3W1MOjwnLyiAhAQ4dEmcuMoSvYWkJPB7U1pKtx7tYvBh++QVOnIDly3tJbHlfVQWPH4OfH1nrC8Xly+DtDR98INRgRUyl749yxMsAqeVLPv4Yfv0V2tqIlNnTAz//DMuWiTNX3vyiOPIWOMrnQ1gYHDwIBgbiTEeGcCBymFY/KAsWwKVLEB/v7ul5wN5+Bik6xMdDUBCoqJCyuLCI9H9EcUNGcZQjgwJI3XDo64ObG6SmEikzIwNsbWHMGJEn9vZCZqZQjcNkjKUlYBg8fUq2Hq84dgzU1MT8qQHIEL6J/B8TCvD17aNQ6jFswdOn5FSDl/92E1wuZGcL+3+Ey4WyMnB8d2FC+QXfESp6oTUAcHLqyso6/+DBA1JWJ9w7Kl76IADk5ICtrZzmJtnblx45cqGnp4dsRaCykrdnD5w8KX6tY2QIByLngaP9odPp0dG7jIwKebz0jg5Zr463m5BzRyKTCRMmCOstuXsXbGxAU1PKOkkTfX3Q1laGQms//LCjtbXEzS2klYx2wwsWwO3bhDWj4HAgLw8WLxZnrnz6RQGgvr6+oGDdiROs7747TK4ms2cvGzdumqNj5Lhx4gtBhnAgzs5QVgZcLtl6CMfKlSsaGqJ1dMzEyNKVEPlvN/Ho0aPt2/e7uv4h5HhF94viKId3VEvrAxWVpr4+jE5G9py6OsydS0wziu7u7p07WbNnd4kR7YJhkJ4O8+YRoAbhqKioqKn18PnNWlrCHb9LBwzDCgt/53KPsdkZkshBhnAgGhpynVY/KJcuAZMJaWmyW1Eh2k34+68sKbFKTFwp5HjFzSDsj3IEjh48uH3XrtUTJmRriJF/TgREeUe9vQNPnEj/7TdxtnXFxaClJae17A0MDMrKrpibf+bjQ2bwOoVCcXbebmNzKipqjyRykCEcBAXyjuK4uEBwMCxdKq1E4DfJzlaAdhPa2qNUVM5bW48UcrzSGEIl2BFSqdTPP3crLzdgs8lRwMeHmGYUzc2dfL4FjdYlxly59YvimJiYhIRMS08n0ynU1wf37gVevhw7efK729APATKEg6BA8TIC4uNBRUV2la/lP30QAMzMftq1a19WVpIwg/GaLFZW0lZK6iiHaxQA1NTAxwcyJPJ4iQ+VCkFBBDSjsLM7t2SJbm5ushhzU1Pl2hACwPz5BIfXikpWFowZA9bWkspBhnAQpk1TsB0hANBocOECpKXB9etSX6urCy5flvd2Ew0NcOMGbNpkSRHuGLOwEFxd5frIU0hsbKC+XhkKrQHA/PkydfgPQPJmFFVVUFjIOHlyuZHocZ+VldDcDC4u4q8uA6ZOhZYWePSINAUSE4k5oEGGcBCsrKCnRwHS6gfg7Q3+/hAQAHy+dBdKS1OAdhNnzkBIiAgh9FdWvAAAGbRJREFUoIqeSi+ARoNx4+SxDqQYzJ0L16+T05IJABwdQV1dIufQ0aOwapWYccgpKTB/vgiNhEiBQoE5cyA9nZzVOzvh8mUIDCRAlHy/zOTh4qJ43lEAuHAB+vrEzFgSHvn3i/L5EBMDa9aIMEU5QkZxJk5UhngZANDVBScnyMoiTQFJmlG0tUF8PGzcKOZ0OT8gFEDirj0tDaZNA0NDAkQhQzg4inhMCAAqKpCYCImJ0nLtYhgWF3cxJ+ey7LM1RCIzEwwMYPJkYcf39sLdu+DsLE2dZIhyxMvgkOsdXbKEFx//86+/imOKo6LA1xdGChuq9RocDpSVwQxyCkaJho8PVFQQlnMpEufOERa4jgzh4Chc4KiAuXPBywvmzZNKeZGUlIthYck9PdFFRdnESyeOqChYu1aE8b//DqNHg7a21BSSLcpkCBctgkuXyOlEAQDZ2XHt7dlBQYfu3Lkj0kQeD378ET77TMx1L12CmTNBXV3M6bJERQX8/ODyZVmv+/w55OURtmlGhnBw7O27iov/c+OGAu4KATIy/u7LRTh6eiP6+mo0NBr09fWJl04QjY2QkwMffyzCFGXyi8Ir16gSFFoDAFNTsLAgrSeMoaGBhsZTLpejo6Mj0sSUFDA1Ff/UWVH8ojik7NrPn4fZswn78YoM4eDs23egp6di3rzP6uvrydZFZNTV4exZOHWK+IOiykqPqVPj7t1LnjJlCsGiiSMmBgIDRWtboxwZhALwQmtVVWTrQRAkekfnz5/722+RFhaXq6tFi9A/ckT8jqFdXZCTI4+Ftt+Gvz/cuEFwv453QlS8KA4yhINjYsJQU3vI53cPE7XVsXwQFAQuLgR/l3g82L8f9u0bY25uTqRcQsEwOHtWNL8oKJ0hBOXyji5YABcvkra6ra3Njh3G334rwhS8TqnYpdGyssDREfT0xJwue7S0wNUVrl2T3Yo1NVBRAR9+SJhAZAgH53//95//+tcuP79cPQX6PL7OtWvAZsOWLYQJ/M9/wNwcPDwIEygNsrPhgw/AyUmEKadP/1Jf/28zM5mXLZcmShM4CgAODqCiQuafExICHI4IGbr/93/w+edAo4m5nGL5RXFkvGtPTITFi0FVlTCByBC+lYAAh5ISRbWCAKCtDRERcPgw/PknAdLw7aBIv4tJISoK1q8XYfy9e/e+/PJMd/fzI0cipaYUCTg4KM+OEADmzSMzdpRGg23bhP3wV1VBdjasWCHmWnw+XLkCc+eKOZ0sFi6EK1dkV+KRWL8oIEM4BGPHwosX8OwZ2XpIwJo14OAAH31EgKj4eDA1BU9PAkRJDzYbmEzRviEMBgPDmjQ0cu3t5btwqogoTaE1HHKTKABgyRJgsyEn590jjx+H1avFD+IoKAAjIxg9WszpZGFoCDY2kJcni7UqKuDFC3BzI1ImMoRvhUIBZ2e4fZtsPSQjKwuqquC77yQSwuPBvn2wYwdBOkmNmBhYuFC0/0GGhobjxhX89FPyokVy2e1GXMaOhbo6JSm0BgD/+AfU15PZD51Gg2++ge3b3zGsrQ3i4mDTJvEXUkS/KI7MfqwkJEBICME1d5AhHIqpU+HWLbKVkAwDAzhwAHbvlqhiXEICmJjI+3YQwyA6WuQwmfZ2+OOPYR9+SER1CnmCRoOxY5Wk0BoAUKng7w+XLpGpwyefAJsNublDjYmOhpkzxUyix0lPhwULxJ9OIgEBkJIi9aQdDIOkJOIbwCFDOBRTpypkfZkBhIeDpSXMni3mdEXZDubkgJqayJlbN2+Ck5NiZC6LCvKOEguNBlu3wr/+9dYBPB5ERMDnn4u/REUFdHfDxIniSyARa2vQ0oKSEumucvMmqKsT/xIhQzgU06bB7dtSr2EtA65e7SsvD6bRxu7de1jUuT/9BCNGgJeXFNQilKgoWLdO5FkslrzvdMVGmTIoAMDPD0pKgKz2hDhLl0JzM7BYg99NTQUjI4lKt6elwYIFCtwCZcECqf9YSUyEpUuJF4sM4VAMHw7Dh5PZZIQoWlsrKJQqPv/0vn2ZIk3k8WDvXvj+eynpRRgcDly9Kk4pcOU2hEqTQQGv2hNeuUKmDvim8G0n7pIk0eMo7gEhjrTbE/b1QUoKhIQQLxkZwnegHN5RBweHCRO0VVU3cbm7TExE2CgkJsKIEeDtLU3liCA2FhYsEDkHubMTysqUpPvSm+CuUeUotIZDuncUAJYtg7q6QTaFd+5AXZ1EZqy+Hh49And3SbQjmalT4cULYvK1BuXaNbCykkpILZmG8NSpU/v37x9iAJvNZpPrCiEpXobH46UT2uOLQqGUlmZyufc4nKk2NjB5Mqxa9dbB1dXVt2/fBgAeD/bsgV27CFPj2rVrndJpLhcTI3KYDAAUFsKkSaChIQWFhuThw4fl5eXSXmX4cNDUhOpqaa/zXxobGwsKCqQnf+5cYDLFaU/IZrPzCArtp9Hg66//+6XIy8vD/0cdPgyffQZ0uviSL10Cf39QUSFCy7dz/fr1ly9fSkk4hQL+/mK2J2SxWBwOZ+gxYqQP1tXVFQrTPwEjiezsbAaDMXLkyCHGhIaG2tjYyEylQSkqwiZPlvWiXV1d6urq0pMfH4+pq2MMBnbnzqB344OCgvBh06cTua6ZmVlVVRWREjEMw7DcXMzWFuPzRZ64fTu2bRvh6ryb3bt3b968WQYLzZ6NpaXJYJ2/uXLliq+vr1SXmDFDnL8oLy9vOnEf5b4+zNoaY7EwDMNcXV1v3LhRW4sNH461tEgkdvZsLDmZEAWHYvz48WVlZdKTf+UK5u4uzkRnZ+fCwsIhBrS3Y7q6WFOTaGJTUlLmzZv3zmHk7Ag7OzvDw8N37txJyuoiMXkyPHoEHUpVfuvvM/8JE8DFBZYtG9x7Rvh2UHpERUFYmDghBkp8QIijZIGjIB/e0QGbQgA4ehRWrAARG1S8Rns7FBQQWTyTLGbMgHv3pNKeMDUVpk8HBoN4yUCWa/Trr79es2bNaEUon6CqCvb2Uo8Jlj2ampCVBefPQ0oKGBrCb78NHJCUBPr64ONDhnKi0NICv/4Ky5aJPJHLhZISZau1PQAlCxwFgEWLID2dtPaEAkJDobr670IqnZ20mBhJu55duQL/+IdoLVPkEzU1+PBDyMggXjLhZdX6Q8GkcJje2dnZ2to6cCUKxdDQEAAKCwu/+OKL/Pz8rKyssLCwqrd3i/nmm2+OHDmi8foZDoPBMDExIVznIaiqCtHQqDEwyJfZinw+/+7du46OjjJYi8cb9vvvezo6Rk2f/jGN1gMAHA6npaWFxztqYJCnp1dK4FqlpaV2dnaqBNbKBXj+fMrz585WVqdEndjZaVZZuczObh+ByghJQ0NDX1+fDJp4dHaaVVYut7PbK+2FcFpbWxsaGsaOlW6xuvLy7ZaWZ4YNaxB+SltbW01NjZ2dHYFqNDV5d3WZdnZ+o6s7v73dz8bmR0mkVVcvVlN7bmgodGFvcbl3756lpaWGNA/G2WzXtjar0aPjRZp1//59CwsLTU3NQe9iGL20dP+ECdtptC6RxLa0tEyaNCk6OnroYVIxhDExMd++UaGWTqdXVlZ2dXW5uLgkJiY6ODhcvXp1aEMIAJGRke2v14kyMjKSsSFEIBAIhIIydepUrXfttaViCIegrKzMxcXF1NQUADo7O9ls9siRI/Pz842NjWWpBgKBQCAQOLI2hDweT+A1ZTKZ4eHhZWVlOjo6VGJLqCIQCAQCIRwSpL2IBY1GE7S61dTUpFKpitv5FoFAIBBKgKx3hP0RuEbJUgCBQCAQCDINIQKBQCAQXC63q6tLV1eXLAXk9GQuOjo6oB+bN28mWyMZcerUqaNHjwqenj9/XgZlB/74448VK1YInkZFRW3YsEHw9OjRo4mJiZLI/+qrr27cuIE/rq2tXbp0ae7QXd2UlyNHjvT/YO/evZtsjYjh4sWLQ5dLJIXu7u6goKBnz54Jrpw6ders2bOSyDx27Jjg69Db2xsSEpL6qs50e3t7UFAQ6VUhB4XP5wcEBNTV1cl+6ZqamuDgYMHTCxcurFmzpqWlpf+Y7Ozs0NDQARMxDAsKCqqpqRFciY2NPXnypDCL8vn8jRs3svrVhN22bdvFixffNl5ODWFpaWlHR8fyV8yZM4dsjWTE3Llzv//++/z8fACoq6vbuHHjRx99JO1FzczMEhMTHz58iD89efLk+fPnG18Vhzh+/Li6ZP36srKy8CSZp0+fenl5mZmZecl/VyfpcOvWLTqdLvhgz5w5k2yNiOHhw4c3b94kW4uB9PX1JScnd/SrC1VcXFxaKlFqLJ/Pj4mJwR+XlJRkZmbGx/+dMFdQUHDjxo3hw4dLIl9K8Pn8lJSUtrY22S/98uXL5ORk/HF8fPyGDRvWrFkj5OYvOTm5f0p6aWlpcXGxMBOpVKqfn19oaCj+J587dy4xMXGIr5usg2WEZ9SoUfMVuiWJWJiYmBw6dGjNmjUlJSWrV69et26di4uLtBfV1NR0cnJisVhjx459+fLls2fPFi9ezGKxgoODa2trnz596uHhIfkqDx488PPzCw8P/+KLLySXprjY2Ni8hx9spcHT03Pbtm1cLldNTS0nJ2fDhg0xMTF8Pp9Kpebm5np5eVEUt52glPnxxx/37t3LZDIdHBxksNz8+fMTExM3b978/fffh4eHJyUlDZFNKKc7wveZlStXWlpaurm51dXVvVmXQEp4eXnhboS8vDx3d3cPDw/8aW5uroODw4gRIySUf/v2bS8vr927d7/nVhCh6EycOHHYsGF4exYWizVr1qwxY8bcv38fAHBDSLJ+8srOnTt/+OGH/Px82VhBnIiIiNTUVD8/v8DAwKHfGmQI5ZFVq1bdvXt38+bNampqslnR09MzJycHAFgslqenp6enJ36Mx2KxCPlunzhxYtKkScuXL5dcFAJBIlQq1c3NjcVi9fX1lZSUODk54b8aOzs7i4uLkSEcFAzD9u7du3PnThnXlzYwMJg7d25ZWdnWrVuHHokModzR0dGxZcuWwMDAvXv3dnWJVlhPbNzd3Tkczp9//pmbm+vh4WFmZtbT09PU1ETUj9z9+/fX1tauW7eOz+dLLg2BEAa8TEf/jxyPx5Pcdenp6clisUpKSsaPH6+uro4bwps3bw4fPlzadVYVFAqFkpSU9M9//jNDxGrcFApFkrfv9u3bKSkpvr6+7ww5RIZQ7tiyZYuDg8P58+fHjBmzY8cO2SyqoaHh6OiYnp5eV1c3btw4AHB3d09KSnry5Ik7ET2zGQzG9evXCwoKwsLCkC1EyAYNDQ1dXd2Ghv+W566vr8frO0qCt7f3zZs3MzMz8bNzNze3mzdv5uTkeHt7SyhZiVm0aFFcXNySJUuys7OFnEKhUIyNjcV++7hc7qpVqw4cOJCQkJCenp6ZmTnEYGQI5YucnJzk5ORTp04BwJkzZ6Kjo2UWjOft7X348GE3Nzf8N5eHh8eBAwcmTJhAVBScoaEhbgvXr1+PsleVDwzDel/RR3qfpFfMmjXr2LFjvb29AHD37l0Wi+Xn5yehzAkTJqirq0dGRuKGUFNT08zMLDY21lPum1v29fUJ3iPZfwcXLVoUExMTEBAgvC2cNWtWREQEl8sFgIqKimvXrn0odM/GHTt2GBoarl69esSIEceOHVu3bt0QQbNyaghVVVWJbdajEHR3d3/55ZfHjx/H+1WZmJgcPHjwq6++ks2/lRkzZnR3dwsijD09Pbu7u4X/2A2BtrY2/m4aGhpmZWXdunVLZkFA8oaqqqqKigrZWhCPiooKk8nUfoW1tTXZGv3Nv//9766uLhMTE0tLy48++igiImLixIkSyqRSqbNmzeLxeK6urvgVX1/f7u5ued4RUigUdXV1Z2dnwXuUlJQkm6X7l9UMCAiIiopauXLlb292QB2MgwcP0ul0U1NTKysrb2/vQ4cOTROug+i9e/fOnTsXExOD/6wPCgqaPn36EKmuqLIMAoFQcnp6elpaWhhS6m6OkJiMjIzTp0+npaUNelcGb5+c7ggRCASCKFRVVZEVVFxk8PYhQ4hAIBAIMhk+fPj48eNJVAC5RhEIBALxXoN2hAgEAoF4r0GGEIFAIBDvNcgQIhAIBOK9BhlCBAKBQLzXIEOIQCAQiPcaZAgRCMRrVFZWnj59ekAPcQRCiUGGEIFAvEZpaWlYWFj/YscIhHKDDCECIS24XG5zc7N4c3t6ehobGwnpw9XV1dXY2NjT0zPoXS6XO8TdQQd3d3e/bQCutvwU3UYghAEZQgRCWNzd3detW4c/xjDMysrK0NCwo6MDv7J9+/axY8fiFSqio6PxZnUMBkNbW3vp0qUCT+PatWvt7e0H9KJycnL65JNP8McvXrxYsWKFrq6usbGxjo5OUFDQ8+fPB9XHy8tr3rx5/a90dHSYmpoKundVVlbOmzdPW1vb2NhYT09v06ZNeCF/nKqqqoULF+ro6BgbG6urq7u6ujY0NCQkJCxbtgwApk2bpq+vr6+vX1FRAQDt7e2rV6/W09PDtQoODmaz2QJRDAZj//794eHhuNq3bt0S7xVGIMgBQyAQwvHpp5+amJjgj0tLSwFARUXl6tWr+JUpU6YEBwfjj/fs2XPixImioqL79+9HRUXp6ekFBgbity5fvgwA169fF4gtLCwEgOTkZAzDuFyuk5OTmZlZQkLC/fv3f/nlFzMzM3d3dz6f/6Y+hw8fplKpVVVVgiuxsbEAUFpaimFYc3OzmZmZg4NDenr6/fv3z5w5o6Ojs2rVKnxkc3Ozubm5kZFRbGxseXk5i8XasmXLkydP6urqvvvuOwCIjo7OysrKyspqa2vDMMzf319NTe3o0aNlZWXR0dFaWlrOzs59fX24NDqdzmAwfH19MzIycnJy6urqCHvREQjpgwwhAiEseHX8iooKDMMOHz5sa2vr7e29efNmDMPYbDaVSj19+vSgEyMiImg0WldXF4Zhvb29xsbGoaGhgrvr16/X19fv7u7GMOzs2bMAUFRUJLiLNxQtKCh4U2xTU5OKisqePXsEV7y9vadMmYI/3rp16wcffNDfJuFqNDU1YRi2bds2KpVaXFz8ptiLFy8K/kycO3fuAED/hWJiYgAgNTUVf0qn0y0sLLhc7lteOQRCrqGTsw9FIBQQb29vOp3OZDJtbW2ZTKaPj4+RkVFKSgoA5OTk8Pl8Hx8fweCioqKsrCz8wKy+vp7H4z19+tTW1pZOpy9dujQyMjIiIkJTU7Onpyc5OXnJkiVqamoAkJmZqaen19bWJmhe2tvbS6FQysvLBd3vBDAYjNmzZ8fGxm7dupVCoVRVVbFYrCNHjuB3r127NmbMmIqKCty3CQCqqqo8Hu/BgwcMBiMrK2vSpElTpkwR5g/Ht7/BwcGCK8HBwatXr87Ly5s/fz5+xd/f/z3sIYpQDpAhRCCERUtLy9HRkclkrl+/Pj8/f+3atcbGxt9++y2bzWYymRYWFmPGjMFHrly5MiEhYcaMGVZWVrhhA4CXL1/id1esWHHo0KELFy6EhoampqZyOJzQ0FD8VlNTU1tbW1BQUP91dXV1X7x4MahKoaGhAQEBRUVFrq6u8fHxdDo9JCQEv/Xs2bNnz54NEKWnp4ef7bHZbHt7eyH/8OrqagAwNjYWXNHQ0NDR0eFwOIIreDdpBEIRQYYQgRABHx+fyMjI/Pz8jo4OT09PvNl3bm4uk8mcOXMmPuavv/6KjY09efJkWFgYfuXnn38+d+6cQIidnZ2zs3NcXFxoaGhcXNz48eOdnJzwWzo6OkZGRjU1NULqM2fOHAMDg7i4OFdX14SEBPwpfktbW9vGxobJZA46UVdXt6mpSchVNDQ0AIDNZo8cORK/wuVyW1tbdXR0BGPwVuAIhCKCokYRCBHw8fFpaWk5dOiQk5OTnp4ejUbz8PCIi4v7888/BX7RyspKAHB0dBTMunLlygA5oaGhubm5RUVFmZmZK1asEFz38PCora3Fw2eEQVVVNSQkJCkpiclkPnz4ULCzBABPT8+ioqK32VRPT8+7d+8+fvz4zVuampoA0D9zA/fK4mE+OBkZGXw+/01vLQKhkJB9SIlAKBJdXV3Dhg0DgK1bt+JXjh07BgAUCqWhoQG/UldXp6qqGhgY2NzczOFw9u3bp62tDQCFhYUCORwOR01NbdSoUXQ6vb6+XnC9tbXV2tra3Nw8JSXl+fPnHA6nsLDw008/ffr06dtUKikpAYBRo0YxGIyenh7B9crKSj09vYkTJ16/fr2tra2xsZHJZIaGhvJ4PAzDqqur9fT07OzscnJy2tvbq6urT5w4UVNTg2FYTU0NnU5ftWpVfn7+nTt3Ojs7+Xw+nk2RlpbW2tqanZ1tampqbW2NB/hgGEan03ft2kXUi4xAyBhkCBEI0cB3fkwmE39aXl4OAPb29v3HnD59GncnAsCkSZPwGMv+hhDDsMWLFwOAv7//APk1NTX+/v5U6t/eGiqV6ubm1tjYOIRKEydOBIAvvvhiwPXff/+9/6ZNVVV11qxZgkyMkpKSSZMmCe6am5sLMjEiIyNHjx6toqICAGVlZRiGNTQ0CHy/ADB16tTHjx8LFkKGEKHQoA71CIRUaG1tffz4sZaWlrW1tRjTORzOX3/9paGhYW5u3v8oTgzq6+tra2u1tLQsLCwE5lnAkydP2Gy2vr6+paXlO8/5amtr6+vrDQwMRo8eLYlKCIRcgQwhAoFAIN5rULAMAoFAIN5rkCFEIBAIxHsNMoQIBAKBeK9BhhCBQCAQ7zXIECIQCATivQYZQgQCgUC81/w/SI2oiq/1y/gAAAAASUVORK5CYII=", "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "metadata": {}, "execution_count": 7 } ], "cell_type": "code", "source": [ "plot_bandstructure(scfres, kline_density=5)" ], "metadata": {}, "execution_count": 7 }, { "cell_type": "markdown", "source": [ "or get the cartesian forces (in Hartree / Bohr)" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "2-element Vector{StaticArrays.SVector{3, Float64}}:\n [0.00044876924247816386, 0.00044877116961621576, 0.0004487694128655534]\n [-0.0004487711696164536, -0.0004487692424776881, -0.0004487694128654898]" }, "metadata": {}, "execution_count": 8 } ], "cell_type": "code", "source": [ "compute_forces_cart(scfres)[1] # Select silicon forces" ], "metadata": {}, "execution_count": 8 }, { "cell_type": "markdown", "source": [ "The `[1]` extracts the forces for the first kind of atoms,\n", "i.e. `Si` (silicon) in the setup of the `atoms` list of step 1 above.\n", "As expected, they are almost zero in this highly symmetric configuration." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Where to go from here\n", "Take a look at the\n", "[example index](https://juliamolsim.github.io/DFTK.jl/dev/#example-index-1)\n", "to continue exploring DFTK." ], "metadata": {} } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.6.2" }, "kernelspec": { "name": "julia-1.6", "display_name": "Julia 1.6.2", "language": "julia" } }, "nbformat": 4 }