{ "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", "[Periodic problems](https://docs.dftk.org/stable/guide/periodic_problems/)\n", "or the\n", "[Introductory resources](https://docs.dftk.org/stable/guide/introductory_resources/)\n", "chapters 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": [], "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 log10(ΔE) log10(Δρ) Diag\n", "--- --------------- --------- --------- ----\n", " 1 -7.900397381616 -0.70 4.5\n", " 2 -7.905003697155 -2.34 -1.52 1.0\n", " 3 -7.905210885294 -3.68 -2.52 3.0\n", " 4 -7.905211520839 -6.20 -3.36 2.5\n", " 5 -7.905211530841 -8.00 -4.81 1.8\n", " 6 -7.905211531388 -9.26 -5.25 3.6\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, Si]\n", "positions = [ones(3)/8, -ones(3)/8]\n", "\n", "# 2. Select model and basis\n", "model = model_LDA(lattice, atoms, positions)\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 (in Ha):\n Kinetic 3.1020961 \n AtomicLocal -2.1987822\n AtomicNonlocal 1.7296077 \n Ewald -8.3979253\n PspCorrection -0.2946254\n Hartree 0.5530387 \n Xc -2.3986211\n\n total -7.905211531388" }, "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×8 Matrix{Float64}:\n -0.176941 -0.14744 -0.0911688 … -0.101219 -0.0239766 -0.0184075\n 0.261074 0.116915 0.00482545 0.0611648 -0.0239766 -0.0184075\n 0.261074 0.23299 0.216734 0.121636 0.155532 0.117748\n 0.261074 0.23299 0.216734 0.212135 0.155532 0.117748\n 0.354532 0.33511 0.317103 0.350436 0.285692 0.417258\n 0.354532 0.389829 0.384601 … 0.436925 0.285692 0.417471\n 0.354532 0.389829 0.384601 0.44922 0.627531 0.443806" }, "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 k-points) 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 irreducible k-points). There are 7 eigenvalues per\n", "k-point because there are 4 occupied states in the system (4 valence\n", "electrons per silicon atom, two atoms per unit cell, and paired\n", "spins), and the eigensolver gives itself some breathing room by\n", "computing some extra states (see the `bands` argument to\n", "`self_consistent_field` as well as the `AdaptiveBands` documentation).\n", "There are only 8 k-points (instead of 4x4x4) because symmetry has been used to reduce the\n", "amount of computations to just the irreducible k-points (see\n", "[Crystal symmetries](https://docs.dftk.org/stable/developer/symmetries/)\n", "for details).\n", "\n", "We can check the occupations ..." ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "7×8 Matrix{Float64}:\n 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\n 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\n 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\n 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+gvaeTAAAgAElEQVR4nOzdeUBUVdsA8OfcmWHftwFURBZ3ZBFQcEUtTUvbtHxd6rWyxTQtKyv3rHxfP83lbc9K21QylzZLBdxAVBCVzT1RBhh2hn1m7vn+mCJEUIS5c7fn9xczHu4cjnPOc896CaUUEEIIIbli+M4AQgghxCcMhAghhGQNAyFCCCFZw0CIEEJI1jAQIoQQkjUMhAghhGQNAyFCCCFZw0CIEEJI1jAQIoQQkjUMhAghhGRNaoHwhx9+SE5Obmdio9HIaWYkCQutA7DQOgALrQOw0DpGaoHw8OHDaWlp7UxcW1vLaWYkCQutA7DQOgALrQOw0DpGaoEQIYQQuisYCBFCCMkaBkKEEEKyhoEQIYSQrGEgRAghJGsYCBFCCMmafAPhN99vf+zpF/+7fpPBYOA7LwjxICsr6/Gn58ye/5pGo+E7LwjxSaaB8Mfde17cFP9b73krEm68ufJdvrODkKVVVlbe8+iM7W6PfG6MiXtgMt/ZQYhPMg2Ee/cfqhw2B7oNqL3ntZ2/J/GdHYQs7dTZ7KrusdBzGA29v4i1Lywp4ztHCPFGpoEwLibS8eTXUFGgPPpZqU9U5G7D5+fZGhwiRTJwtoy+cMw4ObOH8dJx0OTA5ZQGXVnYPoc3Thqv6ijfuUOIBzINhE9M/9fSBwcO/PX5uYF1BdtWrI5S/HGD+n2vf/ao8WwZtgVIghqMEH+Vvec3w337jE5WcGam+uDXm0acWnX/n1+e/e37oxNVFGDwXsM9vxnir7IGlu/sImRBhFJJtfvz5s0LDg6eO3duexLrdDpHR8emlwW1sPUi+2EO62kDs3szM4IYWyVnGRWtFoWG2oPfQrtcRT87z351gQ1xI7N7Mw91Z5Rt3AA3GGFvHvtpLptdDjOCyZy+TDd7YtnM/gO/aR2AhdYxMu0RtsrHDl4PZa4+plwdpfgpjw3Yrl+Eg0VItFgKB/LplIPGQXsMdQZInqjcf59yco82oyAAWCtgcg9m/33KA+MVdQYI2Wl44A/DgXxp3SwjdAvs8rTEEBjThYzporxURT8/zw7eaxhwp/tohASlsA62XGA/ymE9bGB2b2bLCNXdjm30cSEbYhSrIhXfX2YXphrrjDCrJ/N0b8bdmpscI8QrbNrbFOREVkcp8h5Xze7NfJrL+m83LDppvF6DN8dIuI4W0ikHjf1+0F/R0T33Kk49qJzdu+Mj/I4qmN2byXhY+V2c4oqOBu/QTzloPJCPVQBJDfYI78A0WDS5B5NTQT/OYcN+NIz2ZWb3ZkZ3IbxNniB0s8pG2H6F3ZjFshSe7c18MVzloDLn9Qd6kE+GKv4TrdhygX3umNHFCmb3ZqYFMfbYfiBJwB5he5kGi/58XDWmC1mYauwdb/jPGba04a9/NRqNZWW4EwtxrrGxsbKysullWgl99qixx3b9gXy6frAi+1HlS/0Z80bBJi5W8FJ/5sJk5eooxYH8v1ZZn2u2yrqsrAyfkI4soLi42LzLPPGO7u6YBotm92bSSuinuWzwDv0YXya2KuW/b8w12rt3dWCO/Lbbzs6O72wiafp2+w8LlrxN7Fwigv0eX/XlhmzQ6eHpXsyFySoPGwvl4e9JdEVBrWLrRfb+P4yeNvCEX+1nLz6o1SsZXfGPWz4ePGiQhXKDZObGjRsjJjxSY+1qW1easHdHjx49zHJZ3D7RqaXGxfXw5QV2ydTRjf/eAi4+Vgkb3h/u9MJzsztzTYHD9dkdYK5C69InQjM3EazsyPfzh46btPRf9/A+RK9nYc819s3/bLqkIzTueSi7HvrTsxlH9nf+yvhN6wDJF9qTLyzYajWK9h8HFw4/Whwf/+XHZrkshz3CP/744+eff3Zzc5s9e7avr++tCS5fvrx58+a6urrJkyfHxsaa3qytrf3qq69ycnK8vLxmzpzZvXt3ALh06dLOnTubfnHy5MkBAQHc5bz9PG3gtQHMZzb6SzYOAGCwca6ureU7U0iyjJSC0hoA7BydF/aqH9OF/3lqFQOP9mCy/RuXX/YCALBxrG9o5DtTSLJq6uqpszMAgK1zTW2duS7L1Rzh9u3bZ8yY0adPH61WO3jwYJ1O1yKBRqOJjo6mlPr7+48fPz4pKcn0/oQJE/bs2RMTE1NZWRkeHp6fnw8AOTk5mzZtKv+bXq/nKNsds+K1+R4fT7DdMd/28MdPTp/Kd3aQZE1+fBqz4X63HS8EFKXce++9fGfnH888Mb1r8kbnnQuYjQ8sWvgy39lBkrV4wfOqbS857HxFvX32ilfbNfLXLpQbYWFhX3/9tennYcOGffzxxy0SLF26dPLkyaaf16xZc99991FKTetNrl+/3uIie/fuHTJkSHs+d+7cuRs3bmxnJquqqtqZ8o7y8/MPHjradYvuVDFrrmsKkxkLTT7MVWgP7zcs+uXi8ePH9Xq9WS5oRrW1tUeOHJm668abJw1muSB+0zpA8oX23SVjxLfFh48cKSsrM+NlOekR6nS6jIyM0aNHm16OHj368OHDLdIcOXKkKcGYMWOOHDkCAE5OTv7+/idOnACAgoKC/Pz8kJAQU5qioqJVq1Z9+OGHeXl5XOS5k3x9fUcNH/JWlN0bJ3HVHOLEyWJ6opguuTdo0KBBSqXglrnZ2toOHTp0zT2+n+SwN3C7LeKAnoWlaeyake7Dhg51dXU145U5qU4FBQUA4OHhYXrp5eV16NChFmkKCwubJ6iurq6qqnJyctq9e/e4ceNefvnlkpKStWvXhoaGAoCjo+PQoUMppUePHl20aNHPP/88fPjwVj/6zz//PHnyZEZGhumlQqFYtWqVg4NDq4nr6uoUCkWn/9x/TPOD9eeUv1ypj/OW7KHFZi80OTBLob16XPFWCIXGxloBz8G5EpgRwKw4xW6I6uwdIX7TOkDahfbJBSbAgQx2abyrlRhWVlZ3vHHkJBBaWVkBgNFoVKlUANDY2Ght3fJoJpVK1fRo+MbGRkKIlZVVXV3dY489NmfOnOnTp2dlZT311FPh4eHR0dEjR44cOXKkKfGyZcuWLFlya2Q1cXFxsbW1jYqKMr20tbV1dnZuqxRazVgnrRhIl50hY7ur+F/GwA0uCk3yOl9o+25AYT39dy8i/HP+FkdAn53sS/0VfVw6VQnwm9YBEi60aj2syWZ/uoextr67sMUwd64znARCHx8fhmFu3LgRFBQEAPn5+V26dGmRpkuXLjdu3DD9fOPGDXd3dxsbmwMHDpSXly9evBgA/P39J0yY8P3330dHRzf/xUGDBn311VdtfbSLi0twcPBzzz3XnnwqFAqz3z1NCYT1WYYfr5EpAYJvsTqEi0KTvE4WGkthSbphdTRjrRLBl8rdFl4OISsz6I7Rnfqe4DetAyRcaBvOsqN8mYGenPx1nNQra2vrcePGbdu2DQDq6+v37NkzceJEAKipqdm3b5+pIzhp0qQffvjBdA7F9u3bJ02aBACenp6VlZVFRUWm61y4cEGtVgNAdXW16R1K6a5duwYMGMBFts2CAKyOUixOY/WSHRxFlvb9ZVbFwKTuIoiCJvP7M8e19LgWZwqReZTUw6Ys44oIrqoAV1PuK1euHDdu3NmzZy9evBgQEDBhwgQAyMvLu++++8rLy11cXGbMmLF58+bhw4d7enqePHnStJomNDT0X//6V3R09JgxY7KzsysrK5955hkAePrpp69du+bn55ebm1tdXf3rr79ylG2zGOFDejjA5vPsc31E03IhwWpkYVk6+9kwhYgG220UsCScWXTSmDRBcIt6kBityjD+K5AJdOKqEnB4skxpaemRI0fc3NyGDh1qGqVtaGjIyckJCQkxdd71ev2hQ4fq6upGjBjh5OTU9IuZmZkXL1709PQcPHiwaXqvpqYmPT29sLDQ19c3KirKNAfZKgufLNOWM2V03G+GC1NUjtyc+sgjyR9dwYXOFNrGLHZ/PvvTvSKLKEYKA3Ya3o9R3NvRXf/4TesASRbanzoauduQ9ahKbcvVR+ARa1x9af6VaOznSt4Kk1qnUJI1jWsdLrRqPfSM1/86VhnmLqIO4V9+/JNdmc6mP6RkOpR3/KZ1gCQLbUaSMciJLONsXBTw6RPceSeSef+cUWu2M4CQHP3fOeM9XRgxRkEAeNifsVPCjis4W4467lwZ3Z/PvhzCbajCQMiVHo7k8UDmP2dxfz3qoOJ6+CCb5fRGmGuroxVvnmIbMRSijlp00rg4XMH1HJOI65jwLYtQbLnAXtFJavAZWczKdOOMICbAUZTdQZPh3qS3C3yWi5EQdcSRQppdAc/04jxOYSDkkKcNzOnLrEzHVgDdtT919PvL7KJQ0e8J+0+UYtVpo05Y5+QjcVh00vhOJGPNfSXAQMitVwco/shnz5RhpxDdnbdOsfP7K7w4WyZnMSFuZHQX5v1MvB1Ed2fXn2y1Hh63yMkkGAi55aCC1wcolpzCVgDdhbNlNEHDvtRfItVz1UBmY6axCBeOoXYzUlh8iv1vtKJjS47vlkRqmpA935fJKqdJBdgpRO312gnj0gjOFwhYjL8jmRbEvHcGF46h9vryAuthA2O7WmiCHAMh56wYWD6QWXTSiJEQtcehAnqhEp7ifoGAJS0JV3x7ib1chZUA3Vm9EVams6ujLTdBLqnKJljTAplGI/x0DQdI0R1QgEUnje9FMVbSqpoeNvBiX8VyXDiG2mFjFjvIi8R4WW69tLRqm1AxBFZFKl4/yRqwHUC3tfMqW2eAyT0kWDFfCWEOatjTpdgpRLdT0QhrzxlXDrRoFZBgfROm8d2Irx18fQkjIWqTkcLSNPb/BllogYCFOajgzTDF4lM4U4hu570M40PdmU4+zPJuYSC0nPeiFEvT2DoD3/lAQvX5eVZtC2M6ekq18D3bmzlfCQka7BSi1mlq6ebz7OJwSwcmDISWE+1JojzJhznYKUStqDPAO6fZ1VGi30F/GyoGVg5k3sCFY6gNy9LY2b2ZrvaWvhfEQGhR70Ux/zljLGvgOx9IeNZnsbFqMsiCCwR4MTWQMVDY/SfeDqKWzlfSPdfYhQN4uBfEQGhRvZzJpO7M/+FJ3Ohm5Q2wPtP4tmUXCPCCALwbqXgDF46hW7xxkn0tVOFmzcNHS7/iCc2KgcynueyNGhwcQv94J8P4aA8m2Fni3UGTsV1JV3v46iJGQvSPE8X0ZDGd04efkISB0NJ87chTvZhVp7EVQH/Jr6FbL7JLwqU8O9jC6mjF8nS2FheOob+9cdK4ciBjq+Tn0zEQ8uCNMMWua2xOBXYKEQDAkjT22T6Mt/jP126/SA8S40X+l423gwgA4NfrVFMLM4J4i0cYCHngYgWvhCiWpmErgCC3gv52nV0YIqPuoMk7kcz/ncWFYwhYCotPGf8TxSj5C0cYCPkxrx+TqqXHtdgplLvXT7KvhSqcrfjOh8X1dCYP+zP/wZO4Ze+7y6yVAh7ozmcwwkDIDxsFLI1gFp3EVkDWUrX0bBl9gacFArxbFqHYfJ69jgvHZKyRheXp7OooBb/rxGRaA4Xg3z2Z4jr4/Qa2AvK16KRxRYQlHsAtTD52MLs3sxJP4paxj7LZvi5kpA/P66UxEPJGQWBVJPPaCSOLoVCWfspjSxtgGn8LBITg9VDF3jw2GxeOyVK1Hv5z1vh2JP9VgP8cyNlD/oyDCrZfwTti2WEpLDnFvhel4HlIiG/OVvDaAMXiU1gF5GjNWeO9XZhQN/7rAAZCnq2OUrx1im3AuUKZ2XqRdbSCCd34bwJ4N6cvk1ZCk4uwUygvxfXwQTa7NEIQMUgQmZCzYd6kjwt8dh7viGWkkYW3pX6+dvvZKGA5LhyTnxXpxid6MgGOgrgXxEDIv9VRilWnjTo93/lAlvK/LHaAGxmiFkQTIAQzg5myBvj1OnYK5eKqjm67zL7Ox/narcJAyL8QN3JPF2bdOewUyoJOD/93zrhKAAsEhENB4N1IZtFJXDgmF2+dYheEKLwEc5oS1kZBeHsgsynLWFTHdz4Q9/5zxnhfV6afK3YHbzKxO+Okgm8v4+2g9J0po4kadl4/AUUfAWVFzvwdyfRg5t0MnCaRuIJa+DiHXSaMBQJCszpKsTQNF45J32upxmURCkcV3/loBiukUCwOU3x3mb1chWNDUrbytHFWL8bPAbuDrRjqTUJcycc52CmUskMF9FIVzOolrNAjrNzImYcNzO2nWI6nbEjXxUr6w1X29VChLBAQoP9GM6vPGKtw4ZhEUYBFJ43vRTFWAos8AsuOvL3cnzmoYU+XYqdQmt46xb4SonDn4wHcYtHbhYzrxvzfWRwelaYfrrIGFiYHCC7uCC5DcuaggrfCFG+dwlZAgk6V0GQtFdQCAWFaNZD5OIctxIVjkmOksCyNfY/v87VbhdVSWGb3Zi5Wws8XqsrLy/nOCzIPo9FYWFj4xgnDsnDGjqcHcItIF3syI5hZeaqhsLCQ77wgs9FqtR+dretqD2O6CDAOAtZLYVExEJq24aE3trnaWU0aFfvZhjV85wh1Sm5u7j2PTKu2864pL96ctBvAl+8cicD9bMaY6U/t8PL1pFXH/tjr5ubGd45QxxkMhnsefCyrUFeqLVq1fCnAI3znqBXYIxSW2traQ3u/Nyw6WvxS0p7jOZcvX+Y7R6hTXln23o0HN1U8u9sw/s13127gOzvi8MayFeyz35c+u+di2FPvf/gp39lBnbJv37504l/83C/s64c+XLea7+y0DgOhsBgMBqKyAcIAANg41NfX850j1Cn1DY1g7QAA1Mq+pg7/N9uloaEBrOwBgLVyqMNCE7mGhgajlT0AgNLKaBToqngcGhUWJyen+4dF7/38kTJqF+YLffv25TtHqFNWLXpp7PQnG3sMdruesnjPNr6zIw7vvbXwifmPlPlEuuanvrR/D9/ZQZ0yfvz4rv/ZePFrjVvFxZdfeJrv7LSOUMrVYv3S0tLc3NzAwEBvb+9WEzQ2Nqanpzs5ObVo7q9cuVJYWBgQEND8F41GY3p6urW1dUhICCFtTrfOmzcvODh47ty57cmhTqdzdHRs319jUefPnx//c238E2ERHoKbWBZsoQnW/btKwuoyF02MdHBw4DsvolFSUrJwT7Znr4g1Q7HQ7oIwq+fXFxq/3H9666QuXbt25TsvreOqR/jDDz88++yz4eHhGRkZ77zzzrPPPtsiwbVr10aNGuXt7V1UVDRgwIAdO3YolUqj0fj444+fOnWqf//+qampL7744tKlSwGgpKQkLi7O2tq6trbW29v7l19+sbUVzHGtHOjVq9d9pcaEAirAQIjuCkshtcZ5/ViMgnfHw8PjsZEDV2ThpkspOFykeDguqmtX4c7EcZIzvV4/b968rVu3Hjhw4I8//li4cGFlZWWLNG+//fbo0aOPHTt29uzZ7OzsvXv3AkBiYmJSUlJmZuZPP/2UkJCwfPly0y6CtWvX9uzZ89SpU2fOnKmtrd26dSsX2RaUOB+SqBHoeDpqv4xS6mVD1DZ4SMJdi/Zgs8tpRSPf+UCdlqChcT6CvqfnJBAeO3aMZdnx48cDQERERGBg4K+//toiTXx8/KxZswDAzs7usccei4+PBwBCiK2trY2NDQC4ubkpFAqWZQFgx44dTz75JACoVKrp06ebEktbnC9zrIjqMRSKXEIBHeUr6CZAsKwYGORFjhbiPYS45VXTagPtK+zHrXAyNHr9+nV/f/+mmbwePXrk5eU1T1BZWVlVVeXv72966e/vf/DgQQAYNWrUgw8+OHHixPDw8KSkpI0bN7q7u1NK8/Pzmye+fv16Wx9dU1Nz/vz5AwcOmF4qlcoRI0bcZk5RsNysoYcjOVVCY7zEl3nUJFHDPiWw84VFJM6HSSxg7/fD01lF7KCGjvJlBN6KcRIIa2trraysml6a5vZaJACApjQ2NjY1NTUAoNPpcnNzXVxc7OzslErl2bNnKaUGg6GxsfHWxK3Kz8+/cuVKdna26aWtrW2vXr3amj2uqakRcowc5qn8/U9DiJ2B74zcROCFJigGFo4WWn0Y2YCF1gE1NTWDXRULTimr++F5a+0lwG/a/jxVrCdbXc3byZE2NjZK5R0iHSeB0Nvbu7S0tOllaWmpj49P8wReXl4KhaKsrMx0ZkRJSYkpweeff97Y2Lhjxw4AePnll7t27Tp16tThw4e7u7uXlZWZfrcpcat69uw5YcKEdq4apZQKeQnDvX50Q5ZxuYMN3xm5icALTVBStDTAydjd3V6nY7HQ7haldLiX3bVD+nqlg4ewKoFwCbB6HisxrIy2chD2o8c4GbQZOHDgpUuXiouLAaChoeHEiRPR0dHNEygUioiIiCNHjpheHj161JSgurra1dXV9Ka1tbW9vX11dTUAREdHHz16tEViyRvhQ1K1tB6P4BatRA1OEHaKkoEhanK4EKfKxepiJaUUgpyEXgs46RF27dr1kUcemT59+ksvvbR169aBAwdGREQAwJdffvnVV18dOnQIAF5++eWXX37Zycnpzz///OOPP9auXQsAkyZNWr169dq1a8PDw3/88Uej0ThkyBAAWLBgwWOPPdalS5eqqqqvv/46JSWFi2wLjaMK+rqSVC0dIewFV6gtiQXsS/1wfqtT4nyZRA192J/vfKAOEctiMa6m8Tdv3jxs2LAvvvgiICBg165dpjf79Onz0EMPmX5+/PHHN2zYsGPHjtzc3KSkJNNGy9DQ0CNHjly9evWjjz5ydHRMSUlxdnYGgDFjxnzzzTc///xzamrqvn375HPeSpwPSSzA22FRamQhVUuHeougFRCyOB+SoMGFo2KVqKFxYgiEHJ4swwtpnCzT5I98uuq08fD9AjoJT/iFJhCHCujrJ43HJyoBC61DTIXGUlB/qz/3iMpbykdomI2gvmkUwOdbfeokZXdhTxACHrotcEPV5HQprRHWulHULokFrMA3EYsCQ2ComsHDJcQoq5w6qIjwoyBgIBQ4OyWEu5PkIkn12mUiQUPjfLF+mUGcL0kswCogPiJaLIYVVehG+eJZa+JTa4DTpXSIWhytgMCN8iWJOE0oQokFQj9ZrQkGQqGL82ES8HZYbI4V0XB3Yi+guV0R6+dKdHp6rRprgZiwFA4XsGJZ8Y6BUOhi1ASPHhadRA0rlkEh4SMAI32YJLwdFJWMUqq2Jb524qgFGAiFzoqBaE88elhkEgponA9WLrOJw9FRsUkoEMfGCROsqyIQ58vgbkIR0ekhu5wOwtPSzQd3E4pOokZMq6YxEIoALhYQl0MFNNqT2OCRMubT05lQgEtVWAvEwcBCchEdIZ5BEdFkVM6iPMgVHS2p5zsfqH0SC1jcOGF2cT54OygaJ0uovyMR0VHpWF1FAI8eFhcRbZ8SEdxNKCKiqwIYCMXBdPQw37lAd1bWAJeraKSHmFoBURjtSw5qWKwDopBYwIprsZiY8ipnuFhALJIK2KHeRIUVy9z8HIiDkuRUYC0QOjEeN4/1VRzC3Ym2nhbik7oFLxFPVuNMnC/eDopAShHt60pcrPjOx93AGisOePSwWCRoRHOslOjgehlREONx8xgIRQMXCwiftg4K6miYu8haAbEY5cscKmBxnlDgxHjcvMiyK2e4m1D4DmrYEd6MAuMgN3zswNOWnCnDWiBctQbIEOFx8xgIRQOPHha+RFEdKyVGeDsocEeLaJgIj5vHQCgaePSw8Ilu+5ToxPkQPG5QyER63DwGQjHBo4eF7HoNrdLTfq7iawVEJM6XOVJIDRgKhSpRnMfNiy/Hcoa7CYUsQUPjfBgMg5xytwZ/B5JWgrVAiMR73DwGQjHBo4eFLFGDE4SWEOdL8FHVwnSogA7yEuVx8xgIRQY7hYKVVIA7CC0hzofghlphEt3Jak1EmWk5w92EwnSpiupZ6OmMgZBzI3yY41raYOQ7H+gWCaIdFMFAKDKjfUkCHj0sPIkaOlqcTYDoOFtBbxdyohgrgbCUNcAV0R43j4FQZPDoYWHCHYSWNAoPHRUeUR83L85cyxu2AkJDARI14jtfUbzifBjcTSg0oj5uXqz5ljPcTSg0ORXUTkn8HTEQWsgwb5JWQmsNfOcDNZOgoaNEey+IgVB84nyYJDx6WEjEu0ZApOyUEOpGkouwDgiF6bj5UNEeN4+BUHx87MALjx4WkkR89JLFjfLFs9YEROzHzWMgFCU8elg4WAqHC9mRGAgtK86XwSogHGJfLIaBUJTw6GHhOFNGPWxIF3sRtwJiFOtFMstplZ7vfCAAEP9x8xgIRcl09LAeQ6EAiL0JEClrBUR5kqOF2CnknwSOm8dAKEqmo4fT8ehhAUgswI0T/IjzYfCsNSGQwHHzGAjFCo8eFgIjhWNFdIQ4z1cUO6wCAiGB4+axAosVHj0sBKeKqZ898bThOx+yNMiTXK6ipQ1850P2JHDcPAZCscKjh4UgQeSL5URNyUCMFzmMq8Z4JY3j5jEQipXp6OFUPHqYV3iyGr/ifBl8GAu/EiRx3DwGQhHD3YT8amQhVUuH4wQhf7AK8E7sOwhNsA6LGB49zK9ULe3lQlys+M6HjIW7E00tLarjOx9yRQGSJDEogoFQxPDoYX4l4A5CvikIDPNmkvB2kCeSOW6eq0DY0NDw0ksvBQYGRkZG/vTTT62m2bp1a2hoaHBw8JIlS1iWBYAbN27cc7Nff/0VAI4fP978zRMnTnCUbXGxU0KYOx49zJvEAjYOx0X5FueDo6O8kcy9oJKj665atSo9PT0pKencuXOPP/74mTNnevTo0TzBiRMn5s+fv3fvXh8fnwcffNDHx+eFF15wc3N7/fXXTQlu3Ljx1FNPffbZZwBQXFxcUFCwfv160z8FBARwlG3RMZ21NqaLgu+MyE69EdJKaKxaCq2AqMX5kg9ysEfIj0QNfaSHFKoAJ/ezlNLPPvtsxYoV3bp1Gz9+/Pjx47/88ssWaT799NOZM2cOHTo0MDDwzTff/OSTTwDAzs5uzN9u3LgxavmKpYkAACAASURBVNQof39/U3oXF5emf/Lw8OAi22IU58vgQ3p5cbSQhroRRxXf+ZC9EDdS1Uhv1GAtsDTTcfMjvDEQtqGkpKSoqGjgwIGmlxEREdnZ2S3SZGVlRURENCXIyckxjY6aUEq/+uqrWbNmNb2Tm5s7ePDg++6778svv6QUv/R/ifUiWXj0MB8SC1hpDAqJHQEY7o2bKHggpePmORkaLSkpIYQ4OjqaXrq4uGi12hZpSktLnZycTD87Ozvr9frKykpXV1fTO4mJiaWlpQ8++KDpZVBQ0EcffRQUFJSbm/vKK69UVlbOnz+/1Y/OyMj4/PPPly5danppZWV1+vTppg9qoaamhhDR/y9GuKn2X60d62uh0SFpFFrnHbiuWhZqrK5uV7FjoXVA+wst1l3xxzXykDcuG7PoN23fVcUwT1JdXW+Zj+swGxsbpfIOkY6TQOjq6kopra6uNkWgqqoqd3f3W9PodDrTzzqdTqFQNA9XX3zxxfTp021tbU0v+/Tp06dPHwAIDw+vr69fv359W4FwwIABEyZMmD17tulli8u2QCl1cHDo4B8pGGO6ssfLlY/0tNA0oTQKrZOq9ZBdpY/zs7ZtXwXCQuuA9hfafT3o+7lGBwc86c6i37TkMsMTwYw0ip2ToVFPT08nJ6fc3FzTy5ycnFuXtwQGBjZP0KNHD4Xir6a8srJy165dzcdFm3Nxcamra3PfEMMwdnZ2rn+7TRSUjFF49LDFHS6k0Z6knVEQca2XMzFSuKLDWmA5BhaOFtKRUlk1zcmfoVAopk2btmbNGoPBcPHixV27ds2cORMAioqK5syZYwpjTzzxxNatWwsKChoaGtavX29KYPLtt98GBQWFh4c3vXPkyJGSkhIAuHbt2nvvvTdu3Dgusi1S0Xj0sMUlanDjhLCM9CG4asyS0kpodwfiIYXeIAB3+whXrVpVXV3t4eExePDgZcuWhYaGAoBOp9uzZ49erweAsWPHPvnkk71791ar1V5eXgsXLmz63cTExHnz5jW/2v79+wMCAuzt7cPDwyMiIt59912Osi1GSgZi1Xj0sEXhWdtCg7sJLSyhQCI7CE0IpyswGxoarK2tb5OAZVmDwWBl1a5Tqurr621s7nAHMm/evODg4Llz57bngjqdrmlFj6itOcter6EbYywxTSiZQuuwikbw+15fMkNl1e7bSCy0DrirQruqo0N+MuT/SyWdtrlDLPZNu/c3w9x+zAN+EhkX4fbPuH0UBACGYdoZBQHgjlFQtuJ8cVzIcpIK2CFq0v4oiCyghyOxVpDzFVgLLMF03Pwwb+nUAen8JXIW7k4KamkhHj1sEYkaGueLFUdwRuE0oaUcl9xx81ifpcB09PAhnCa0CMmcrygxcb4Et9VbRqLkqgAGQonAxQKWoa2DGzU03F1SrYA0jPIliRqWxUrAPekdNy+pP0bOcDehZSQWsMN9GAXGQeHxtSPuNuRcOdYCbtUaIK2EDpXEEaNNMBBKRH83UtVI86qxFeBWooZK4DGkUoUPrLeA5CIa5k7spXWaBAZCiTAdPZyEnUKOJUpr+5TExPngNCHnEguk8Ej6FjAQSgcuFuCappaW1tP+rlJrBSQjzpc5XMAasRJwKUGKq6al9vfI2SjcTcixgxo6ypdhMA4KlacNdHMg6SVYC7ii00NWOR3sKbU6gIFQOkxHD1+uwlaAK4kaPFlN6OJwNyGXDhfSKCkeN4+BUFJwjoRTiQW4Ukbo4nxJIm6o5YxUj5uX4J8kZ3G4ao4zV3S0wUh7u2AgFLSRPkxKEW3EUMgNqS4Ww0AoKaN9yUENbinmRIKGjpLcGgHpcbGCYGdyQouVwPwqGuFyFY2S3AQhYCCUmO4OxE5JcvHoYQ7gDkKxGIXLp7mRqGFj1UQlxaAhxb9J3nCxAEcOFUpzUEh64nyYRA2OjZpfYgGV5AQhYCCUHtxNyIXcCqog0MMRA6EIDPMmJ0tonYHvfEhOgnRXTWMglJrRvkwSHj1sbgkaOkaiTYD0OKhggBtJwWlCs9LWQb50j5vHQCg1PnbggUcPm1tigWTvhSUpzgc3UZhZYgE7QrrHzWMglCA8Ysa8KMChAnYkrpQRjzhfBquAeUl7sRgGQgnC3YTmdbaMulqTrvaSbQWkZ4ianCujOj3f+ZCQBInuIDTBQChBo3yZI4V49LDZSO953JJno4CBHuRYEdYB89DU0vIG2k+6x81jIJQgd2vo5kDS8OhhM8GT1cQozhc3UZjNgXwa5yPl4+YxEEoTPqHUXIwUjhSyIyS6fUrCRuGGWvOR/GIxrN7ShKvmzCW9hHa1J2pbvvOB7lK0F7lQScsb+M6HJCRJeoIQMBBK1UgfJrmINhj5zof4JUh6sZyEWTEw2IscLsTbwc66XEUbjLSXs5RrAQZCaXK2gl7O5GQxDg11VmIBK+1BIQmL82XwlKXOSyygo6V+3LzE/zw5i/MlCdgKdI6ehZQiOsIbq4ko4Uy5WcjhedRYwyULjx7uvFQtDXYmrtZ85wN1yEAPcr2GFtXxnQ+RS5LBqmkMhJI1zJuk4dHDnSPVx5DKhILAUDWD04SdkVNBVYz0j5vHQChZDioIcSPJePRwJyRqWKk+d0Ym8JSlTkrQ0NEyuBfESi5lo3wJjo52WL0RTpXQId7SbwUkDB/P2UmS30FogoFQyuJ8cNVcxyUX0RA34qTiOx+oE0LdSVkDza/BWtARFOBwASv5CULAQChtsXj0cCckyqMJkDYCMBxvBzvqbBl1syZdZHDcPAZCKbNRQKQHOVqIrUBHJGhonNS3T8lBnA9OE3ZQgmyOm8d6LnFxvgyetdYBNQY4V0ZjvWTRCkjbKF9yEANhh8hhB6EJBkKJs87a9/H8qU88P1+r1fKdF9HYFr9z8PhHHfe80VBdwXdeUGfZVeZpP58z5IHHjx47xndeROPChQuTpj21b8VM18IMvvNiCXcIhDU1NRcuXKisrLRMbpB55eTk/Pe/q3X3r/iGGTpx2lN8Z0ccUlNT56z+NDPu3UK3/o8/PYfv7KDOuveRaXWhDyVHvfHI7AWFhYV8Z0cEWJYd89DUvd1m6scsmDHraZ1Ox3eOONdmIPzqq6969uzp4ODQq1cvFxcXPz+/9evXU4ojDGKSnp6u6z8JvILY8EnX8gv4zo44pKSeLA99DDx7sIOnZeZc4Ds7qFMaGhoqGij0GQ2+fRuD486dO8d3jkSgqKio0dEHgmKhW6ihe9T58+f5zhHnlK2+u2bNmtdee23QoEFr1qxRq9WlpaW//PLLggUL8vLy1q1bZ+Esog6Liopyeu/JkqDhRJMV1L0b39kRh2FDYly/eLWse7jq0tHwAf34zg7qFGtraw87ZfGZn6iT2vrCwdDQV/jOkQh4e3vb1BRB1n6wtlNdO9W791q+c8Q9eouGhgYnJ6dnn322xftvv/22QqHQarW3/opwzJ07d+PGje1MXFVVxWlmhCAhMWnQg0/4Tn6jrKzMLBeUQ6Gt/uZXlxEzFi5eodPpzHJBORSa2Zmr0DQazb9ffFk1YtYfyWlmuaCQmavQrl69anPPnAnTZ2dlZZnlggLXSo+wuLi4qqpqzpyWsyMvvPDCkiVLrl275unp2Z74Gh8fn56e3qdPn2nTpimVrXxQWlranj17HBwcnnjiCbVaDQC5ubmHDx9unmbKlCkuLi4AkJ2dHR8fr1Kppk+f7ufn1/5IL3NxI0fsHzLC9zu9vTPuDG8v2n/sv1ffu2awgu+MIDPw8fH5YtPawt8NNbgZpt0Mbt09/73+56mtDxlKTyvfDAcHB0JIaWlpi/dLSkoUCkX37t3bc92FCxeuWrVKrVZ/8sknTzzxxK0JEhISRo8ebW1tfenSpcjIyPLyctNHpP0tPj7+lVdeMUXQ9PT0mJgYlmWLi4sjIiJu3Lhx13+ojDmqINCRZJTi/G57pWhpjFoWq8blI8aLSSnCKtBeyUU0VlZVoNV+4pgxY/r373/p0qWmd/Lz80eMGDFv3rz2dDNLSkpsbW1Nv15WVmZra3vx4sVbP2LdunWmn++99961a9e2SPDMM8/MmjXL9PPUqVPfeOMN08/Tpk1r+vlWODTaquePGtafM5rlUpIvNJZSz68br1ezZrym5AuNC+YttAP57NCf9Ga8oDCZq9CeM1+LIQqtjxV8+eWXjY2NPXv2jIiIGD9+fHR0dI8ePdLS0jQazZS/JSQktBVcjx8/3rVr18DAQABwdXWNiopKTExsnoBl2UOHDo0dO9b0cuzYsUlJSc0T1NXVxcfHz5o1y/QyMTFx3LhxTYlbXA3dUYyapOBjKNrnUiW1VZKuMjhWSlYGe5GMUtqIZ0u0T4rMeoStDwF37do1LS3t66+/TkpKunHjBsuyISEhAHD16tWmNFVVVW1dtKCgwMvLq+mlWq0uKLhp7X5JSYler29Ko1arNRpN8wQ7duzw8vKKjY0FAIPBoNVqmyYmb03cXF5e3smTJ8+cOWN6aWtru3z5cnt7+1YT19fXq1SymDmLcIZFhYr6ejOcOir5QjuUTwZ5kPr6ejNeU/KFxgXzFpoCIMCBSdU0RHlI+Y7QLIWm08NlnaKXfYNZKwFvVCqVQnGH+f4250IdHByef/75559/vgMfrFQqWfafWy+j0dhisYzpZVMag8HQ4j9v8+bNTz31FCEEABiGUSgUt0ncnL29vZ+f38CBA5s+yNbWtq1SUCgUdywgaejpAiyFwgZFF7vOXkryhXaiBGK8wLx/o+QLjQtmL7RYNZwoJYPVZryk4Jil0NKKINwNbFUS+caa4sjtcbIoyMfHJz8/v+mlRqN54IEHmidwdXW1sbHRaDSmTqFGo/Hx8Wn614sXLx4/fnz79u2mlwzDmHqB/fr1uzVxC+7u7sHBwe2M3yqVSj736YO8jCdKyRTnzi6ck3yhpZYYZvdVqFTmHBeSfKFxweyFNsSb/SmPqqTSvrfKLIV2ooyN9ZZ4QbXAyXrioUOHVlZWnjp1CgDy8vIyMjJM04GFhYWmQUtCyIQJE3bu3AkALMvu3r37/vvvb/r1L774Yvz48c2j3f333//jjz+afv7xxx+bJ0btFOOF04R3VqWHKzo6wE1GsyPyEeOFT2Jpl5QiNkZmx81z0iO0t7dfvnz5pEmTJk2a9Pvvv8+fP98U1Xbt2vXJJ59kZGQAwOLFi0ePHn39+vW8vDyj0Th16lTT7xoMhq1bt3700UfNL/jqq6/GxsZWVFTodLrLly9/8803XGRb2mLVZGEqLhW4g1QtHehBrHC/mRQFOhGW0us1tBuuhGobBUgtpl+NkFcd4Gq/5Pz580eOHHn69OmZM2cOHjzY9OaDDz7Y9HNYWFh2dvb+/fsdHR3HjRtnbW1ter+hoWHLli0jR45sfrWAgIDs7Ozff//dyspq3LhxbS1+QbcR6UGyymmtAezkske2I5KLqNzuhWVlsBeTXEQfC8D/4jZll1N3a6K25TsflsVhoxgWFhYWFtb8HR8fn+YDnmq1evr06S1+y97efsyYMbdezc3NranXiDrAVgn9XElaCR3mja1Am1K07At95HUvLCumfUSPBfCdDwGT3VZ6AMDnEcpKLO4mvC0KcKKYxqixUkhWrBdJxvNlbitFK8dBEazzMhLjRfCUqdvIKqceNsTThu98IM5EepLsclpr4DsfApaixR4hkrQhapKsxfUybUoporHyuxeWFRsF9HcjaSV4O9i68gYoqKX9XGVXCzAQykgXe2LNkMtV2Aq0Ds/algMcHb2NZC2N9iQK+VUCDITyEqMmyThN2IZk7BHKAJ67exsy3EFogoFQXnCasC2lDVBYR/vKb1BIboaoybEiFutAq5KLZLpYTI5/s5zhwtG2pBTJdFBIbnztiK0CJwhaYaSQVkIHecqxDmAglJdwd3KpilaZ4SkUUpOiZWW4WE6eYtU4TdiKM6W0mwNxteY7H3zAQCgvKgbC3cnJYmwFWkouojFeWB1kAc/dbZU8N06YYM2XnRhcNXcLAwtpJTRaloNCMhSrxpnyVshzK70JBkLZiVGTFNxNeLMzZbS7XAeFZCjMnVzW0cpGvvMhMHI+aBcDoezEejHHtRSXzTUn50EhGVIxEIETBDcrqoOqRtrLRaa1AAOh7HjZgps1ya3EVuAfKTK+F5Yn3FDbwrEidrBXOx7lLlEYCOUID9doIRnPlJGZGC+SUoQTBP9IkesOQhP5/uVyFoOLBZopqAVdI+3pjIFQRoaomdRinCD4R7K8ZwcwEMoRbqtvLrmIjVHLd1BInjxswN2a5FRgLQAAaDDCmVIa5SHfSoCBUI76uxJNLS2p5zsfwpCixR2EcoTb6pukl9JeLsRBxXc++IP1X44UBKI8yQlcNQcAcn0kN8Jt9U3wAWQYCGUqxgt3EwIANBjhbBmNlPGgkGzhBEETfAAZBkKZivFicFwIANJKaG95DwrJVn9XUoATBAAg7zNlTDAQylSMmpwspgbZ9wlxK71sMQSiPUmq7CcIrlVTA0t7OMq6FmAglCkXK+hmT86Vy70VwHthOcPdhACQXESHyHgHoYnc/345w1VzAHAcA6GMxahxggAnCAEwEMpZjOwXC/ypo5SCv7wHheRssBc5VSL3CYJk2S8ZBQyEchbrJffzZfBkNZlzsQI/B3K2TL61oM4A5ytphOxXTWMglK9eLqRKTwtq+c4Hf/CsbRQr79O3U4tpiCuxUfCdD75hIJQvAjDIkxyX8W5CmZ+viOCv9TLyDYS4atoEA6GsxagZ2U4T1hjgfAWNcMdWQNZkvq0eB0VMMBDKWoyMn8d0spiGuhNr2Q8KyVxPZ6LTU02tHGsBBTiuZQdjIMRAKHODvEhGKW0w8p0PPuBiOQQABGCwFzkuy07hxUpqryJd7LEWYCCUN3sl9HQmp0vl2AqkaFlcMooAIMaLkec0YTKOi/4NA6HcyXOOhAIc19JBntgKIIiR68JRPFapCQZCuZPnqrkLldQBB4UQAABEe5IzspwgwCWjTTAQyl2MmhyTXyDECULUxF4JvVxIuswmCKr08KeODnDDWgCAgRAFOBIKNK9aXq0Anq+ImpPhKUvHtXSgB1FhBAAADIQIAAbL79mE+Ehu1JwMz91NLmJxgrAJBkJkelq9jFqBKj1cq8ZBIfQPGW6oTSnCQZF/YCBEsnseU0oRjfQgSvzuo7/1cCSEwDXZTBCwFE4U08FeWAf+ggWBINKD5FbSWgPf+bAU3EGIbjXIU0a3g1nl1MuWeNrwnQ/B4DAQGo3GzMxMjUZzmzTXrl3LycmhtOX3T6PRnD17trq62vRSr9eXN6PX67nKtCxZK6C/KzlZLJdWILmIxuC9MLqZrKYJ8bj5FrhqDq5evdq7d+/p06eHh4e/+OKLtyYwGAxTpkwZMmTIo48+OnDgwJKSEtP71dXVkyZN6t+//8yZM/38/AoKCgBg3759Xl5egX9LSkriKNuyJZ9t9X8PCmErgG4SK6dpQjxruwWuAuHixYvHjRuXkZGRnZ29a9euw4cPt0iwc+fOzMzMCxcuZGVlBQUFrV692vT+3LlzGYYpKCjIyMjIz893d3c3vT9o0KCyv91zzz0cZVu25LNeJrOcetsSDxwUQjcb6EFyK2i1PAabcCt9C5wEQr1ev3PnztmzZwOAu7v7ww8/vH379hZptm3bNm3aNDs7OwCYPXv2tm3bAKCqquqbb75Zs2aNXq9vaGiwtbW1srJq+pXi4uLGxkYuMoxi1SS5iJVDJEwuwiYAtcJaAQPcyKkS6VeCknooqqN9XLAW/IOTQFhYWNjQ0BAYGGh6GRAQcO3atRZprl27FhAQ0JRAo9E0NjZevnzZxsbm3XffjYiI8PX1nT17ttH418FHp06dCg8Pd3Z2fvTRR8vLy9v66Pr6+uvXr6f9LSMj49YJSHQrXztirySXKqVfVni+ImqLTCYIUrTsYC+iwErQjJKLi9bU1ACAtbW16aWdnZ1Op7s1jY3NX+NTtra2lNLa2tqysrLq6mo/P78LFy6Ul5dHR0dv2bJl1qxZQ4YM0Wq1Tk5OpaWlDz/88Ouvv/7pp5+2+tFXr17Nzc3dv39/05Xj4+OdnJxaTdy0GAcBQJS7VUJeo7f/HY5cFHuhHSu0nhPYqNNZtL0Te6HxwvKFFubIfPen8sUAEQ87tafQDl1XRjiDTldngfwIgY2NjUqlun0aTgKhWq0GgPLycg8PDwAoLS319va+NU1Tx66srMzGxsbZ2dmUbObMmQDg6uo6adKko0ePzpo1y83NzZTS3d39lVdeWbBgQVsf3adPn4kTJ86dO7edWXV0dLy7v026hvmyp8uVzzre+Um14i20knooa9RHdnFgLH47LN5C45GFC22UP52fZnBwdBR1Z+mOhXaqwvBmmELkf6WZcTI06urqGhAQkJKSYnqZkpIycODAFmkiIiKSk5NNPycnJ0dERBBCAgMDXV1dmwdIZ2fnFr9YUFDg4uLCRbZlLlYGD6NJLmIHeRLLR0EkCr52xNGKXJT0BIGBhdMl+ACyljjpEQLAvHnzXn31VVtb28zMzGPHjm3evBkALl26NHr06HPnzjk5OT3//PPR0dFDhgxRq9UrVqxYt24dANjY2MyZM2f+/Pnvvffe5cuXf/jhB9Ny0w0bNjg7O/v5+eXk5Cxbtuy9997jKNtyFuZGruloZSM4W905sUilaGmMGncQojaZzlrr6SzZOJFRRrs7EgnX8Y7hMBCqVKp169a5uromJCR4eXkBgIODw/jx403Dtb169fr55583bdpUV1e3evXqKVOmmH5x+fLlGzdufPfddz08PPbv3x8WFgYA/v7+8fHxRUVFPj4+W7ZsmTBhAkfZljMlA+EeJLWY3ttFsq1AchFdHI6BELXJtI/oyZ5854Mz+ACyVhGJLaqcN29ecHBwO+cIdTodztw09+ZJo7WCLIu4XagQb6HpWXD7Wn9jqsryt8PiLTQe8VJoaSX0yUPGc49w1UPg2h0LbWqicVxX8kQw3g7eBIsD/SNGTVK0LN+54EpGKQ3AQSF0W6Fu5Fo1LW/gOx+cwQeQtQoDIfpHjBdzXEuluq8eT9NAd6RkYKAHOSnRbfWaWlpjoEHSnQHtMAyE6B8eNqC2JdkV0mwFcCs9ag/TKUt854ITyUU0Vo2LpluBgRDdRMJPKMXD1VB7xHgxKRKtAngv2BYMhOgmUn0YjaaW1hpooBO2AugOYtQktZgaJVgJTA8gwyrQCgyE6CaxXkSSt8PHcFAItY+7NahtSXa51GpBgxHOldFI3ErfGgyE6Cb9XElRHS2u5zsf5oYPYEPtJ8lTlk6V0D4uxF6sG0O4hYEQ3YQhEO1JUiXXCuAjuVH7xUhxXARXTd8GBkLUkvR2EzYYIaucRnpgK4DaRZLPY8JBkdvAQIhaivFiJLZw9GQx7eNC7HBQCLVPXxdSXC+1CYIULRuDPcI2YCBELcWoSVoJNUioT4iDQuiumCYIjktoXOSqjhIg3R2wFrQOAyFqyUkF3R3ImTLpdApx+xS6WxLbTYibaG8PAyFqRaxaUtvqU4pYDITorsRIa+Eo3gveHgZC1ArTw2j4zoV5XNFRhhA/HBRCd2OwF0kroXqpDI5ij/D2MBCiVkhp1VxyER2CTQC6S04q6OEokQmCGgNcqKTh7lgL2oSBELUi2JnU6Gl+jRRagRQtxcVyqANipXLubqqWhrkTawXf+RAwDISoFQRgsBdzXBKdQnwkN+oYyZy7i6um7wgDIWqdNFqBGgNcqqLhuJUe3T3JnLuLi8XuCAMhap00nsd0XEvD3IkVfs3R3QtyJrUGekPkEwQUILWYDvbCOnA7WDqodYM8yblyWm/kOx+dk4LjoqijCECMWvQTBOcrqJOK+NjxnQ9hw0CIWmerhF7OJL1E3K0AHiuFOkMC+4iScbFYO2AgRG0S+8NoKMBxLQ4KoY6TwAQBnrXdHthGoDaJ/WE0uRXU1Zp42/KdDyRa0Z4kU+QTBLiVvj0wEKI2xarJsSIRH62RjPfCqHNsldDbmaSJdoKgohGu19AQV6wFd4CBELWpuwNRMuSqTqytAJ6viDpP1OfuphTRKE+ixGb+TrCE0O2IerEADgqhzhN1FUjRslgF2gMDIbod8U4TVjTCjRraHweFUOfEqkmyaCcIUrQ0BheLtQOWEbod8T6MJqWIRuOgEOo0PweiEucEAUvhZDGN9sR7wTvDdgLdToQ7uVBJq/V85+Pu4aAQMheRbqI4V0597YiHDd/5EAMMhOh2rBUwwI2cFOGqueQiHBRC5iHSc3dx1XT7YUuB7kCMRw8bKZwqoYOwFUDmINIeYUoRninTXhgI0R3EqEmKVmSLBc6VUV874mbNdz6QJAz0IBerxDdBgE9faj8MhOgOhqiZlCIqrvth3DiBzEjFQKgbOVEspkqgrYPSBtrbGWtBu2AgRHegtgUnK3KhUkytAG6lR+YVK7ZpwhQtO9iLMFgJ2gcDIboz0c2R4KAQMq8YL5Iiqt2EuIPwrmBJoTsT16o5bR2UNdBeOCiEzCdWzSRrKSuaSoBLRu8OBkJ0Z7Gi6hEma9kYHBRCZqW2BVfxTBDoWThdSqMxELYbBkJ0Z6HuJK+aljfwnY/2ScEdhIgDIjpl6XQpDXQkTiq+8yEe2F6gO1MQiPQgqSJZNYeP5EZcENG5u7hq+m5hIETtEqsWx2IBPQsZpTQKz1dE5iai5zGl4L3gXeIqEBYVFc2YMSMkJOSRRx65evXqrQlYll21alVYWNjw4cN//vnnpvcbGxvfeeedqKiosLCwV199ten9DRs2DBw4MCYm5vvvv+coz+g2YtSMKNbLpJfQICccFELmN8CN3KihZWKYIEgporE4QXg3lBxdd/r06f7+/rt27dq8efPEiRPPnj1LyE3/MR988MH333+/bdu2q1evTps2LTU1tXfv3gDw9NNP5+XlrVu3zsXFJTMz05R427Zt77///8vSugAAHI1JREFU/g8//FBVVTV58mR/f/+YmBiOco5aNdiL/EtLjYIPhcm4cQJxQ0Eg0pOkaul93QT9BcuvoQ0sDXQSdCaFhpMe4fnz548cObJu3bqgoKBVq1ZptdrDhw+3SPPhhx8uX748JCRk4sSJjz766GeffQYAZ8+e3bVr186dO4cNGxYSEjJ16lRT4g8++OD111+PjIwcNWrUM88889FHH3GRbXQb7tbga0+yyoUeCVNw1TjiTKwYjhs8VkRjcbHYXeKkvDIzM3v16uXo6AgACoUiIiLi3LlzzRM0NjaeP38+KirK9DIqKsqU4MSJE1FRUd99992ECROeeeaZK1euNF0wMjLS9HNkZGSLqyHLEMW2etxKj7gT4yWCCQKcIOwAToZGi4uLXVxcml66uroWFRU1T1BSUkIpdXZ2Nr10cXHRarUAcP369ePHj/fr12/FihW7d+8eOnRobm6ujY1NRUVFU2JXV1dT4ladPXv2u+++W7dunemlUqlMSkpq+t0WampqWgzYotuIcFYc0TAPugq30K7XQIPRyovUV1fznZWb4TetAwRYaAPsyQmtqlJXrRBWvv5RU1NzRGP1brixulroPVeLsbGxUSrvEOk4CYTOzs41NTVNL3U6naura/MEpjBZU1Njer+6utr0g5OTk62t7fvvv69UKiMjI7dt23bw4MGHHnrIzs6u6YJNiVvVt2/fkSNHzpw50/TSzs7O29u7rcSUUgcHh47/nTIT50c3nDfa29sLttDOatmh3kL8P8VvWgcIsNAcALrYG/7U24e6CTQS1uppbhUztJu1HVfLP6SJk6HRHj16XLlyxWAwmF5euHChR48ezRPY2dmp1eoLFy40JfD39weAgIAABweHpujt4uJSXV1tuuDFixdbJG6VUql0d3cP+NttoiC6W31dSWkD1dYLtAkAPGsbcU/gEwSny5l+rgSj4N3iJBAOGjTIw8Njy5YtAPDbb7+Vlpbed999AJCamrpp0yZTmmnTpm3cuJFlWa1W+913302fPh0Axo8fX19ff+DAAQBIS0vLycmJjY0FgOnTp3/wwQd6vb6qquqLL74wJUYWRgAGeZKTpcKNNHi+IuJajFrQ2+pTSxicI+8ATgIhIWTLli0rV6708/P797///fXXX9vY2ABAVlbWjh07TGmWLFlSU1Pj7e3dq1ev6dOnjx49GgCsra2//fbbWbNmBQUFPfDAA59++mlgYCAAvPTSS25ubr6+vt27dx8xYsSUKVO4yDa6owBd1mebv8jIyOA7Iy3V1dV9sfXbs/u29XMQwz4vJFr9rSt+3/7l7t27WVZwk3CJiYk/fPtVkD6P74yID6GcPXKVZdmysjIXF5fbTFSWl5fb2NjY2to2f9NoNFZUVLi7u7dIXFlZqVQq7e3tb/Oh8+bNCw4Onjt3bntyqNPpTEtbUXv8vv/AlPnLqyJnuJ/Z9tnKlx+a+ADfOfoLy7KhQ0Zd6DLaYGQHFB9OP3JAaIss8JvWAQIstKqqqpDYUddD/mVfmRfnXLH3+6/4ztE/lr27ZuNvaRX+w91TPzv+a3xQUBDfORITDseSGYbx8PC4fZpWl70oFIpboyAAtLX4E1nGh1u2VT26HrqFlvYasXHzcuEEwsuXL2utvBvHLQKAgi3nrl27dptZZIQ6LDk5uSL4Xjp6bjVA6vtDDQbDHZcjWszW7Tsr5iQCoyyzc962c/fi1xfynSMxwX2XqL2CundR5aUBAHPtVGD3rnxn5x+enp6s9jLo66Chhi2+0updFEKd16VLF5XmLLAGqCqyogKKggCg9vKEG5kAYHfjVICfgKqnKAjoPxIJ3PJFCzOmP3V2zSeVjn7Lf/2S7+z8w8XF5Z5/v7J7bZyHDbN65VtCG09DkhESEjLv4VEfvz+k1Gg99931fGfnJv9Zu27M1Ofc9KXjRo94/DFcRXF3OJwj5AXOEXJNp9PNPmUX4U5eHSCU4YR6IwTtMPwyViHY3V34TesAIRfazqvse2fYkw8qhfOFe+qw0c+BvBxcI9hCEzKhtGVIRJaGM/93zqjT852Pv32YzcZ6EcFGQSQ9D/dgDCz8dl0ovYhr1fSnPHZeP2zPOwgLDt21Pi5ktC/zUY4glo/XG2FdJrs4HL/JyHIIwNIIZmmaUB7HsiKdndNX4WrNdz5EC5sP1BHLI5h1wugUfpDNDlGTAdgdRJb1kD9jYOFXAXQKL1dhd7CzsOxQR/R0JqN9mQ+zee4U1hng/Uz2rTD8GiNLM3UKlwmgU/hOBju3H3YHOwVbENRByyOY9zN57hR+kIPdQcSbh/wZI4Vf8vgMhdgdNAssPtRBwc5kTBfmA/46hbUGWHvWuBi7g4gnBGBpOLMsnc9O4aoMdl4/hYsVfzmQBGxEUMctC2fW89cp/DCHHeHDhGB3EPHnQX+GAd46hZer6M957FzsDnYaliDqOFOn8H98dAprDLD2rBFnBxG/CMBb/HUK3z7NvoTdQXPAdgR1yrJwZgMfncIPs9mRvtgdRPyb1J1hAH7Os/Tt4OUq+st19kXsDpoDFiLqlGBnco/FO4U1Blh3zvhmKH57Ef8IwOJwZlkaa+FO4crT7Pz+2B00D2xKUGcti7D0TOEH2WwcdgeRYEzqzigZ+Oma5W4HL1XRX6+zc/piA24eWI6os4KcyNguzKYsC7UCNQZ4/5zxTZwdREKyOIxZnm65TuHKdOwOmhO2JsgMlkYwG7KMVRbpFP4vix3ly/R3xe4gEpCJ3RklA3st0im8VEV/z8fFouaERYnMIMiJjOtqiU5hjQHWZxrxZFEkQIvDmBUW6RSuTGdf6qdwUnH/SbKBDQoyjyXhzPpMY0Ujt5+yKYsd5cv0ccHuIBIcU6dwD8edwouV9A/sDpobliYyjyAnMqEb8z8uO4U1BtiQaVyC3UEkVEvCOe8UrjzNvtRf4YjdQbPCNgWZzZJwZkMWh53CTVns6C5Mb+wOIqF6wI+x4rJTeLGS7s9nX8TFouaGBYrMJpDLTqGpO4gniyKBWxKuWHKK5ahXuOI0Ox+7gxzAZgWZE3edwo1Z7BjsDiLBu9+P2Ck56RRerKQH8nHvICewTJE5BTqR+/3Mv3y0Wg8bcbEoEomlEYqlaebvFC5Px+4gV7BlQWa2OIzZaO5O4cYs9p4uTC9n7A4iEZjQjdgpYbdZO4UXK+lBDXYHuYLFisws0Ik84MdsNF+nsFoPm7KMb2F3EInHMnN3CpenswuwO8gZbFyQ+b0VxmzKMpY3mOdqG7LYe7tidxCJyfhuxEEFu/40z+3ghUp6UMO+gN1BzmDJIvMLdCITzdQprNbD/7LwZFEkPsvCFcvSzdMpXJ7OvhyC3UEOYfuCOLE0gvkg2wydwvWZ7FjsDiIRuq8bcVDBj53uFOZU0IMa9vk+2FZzCAsXcaK7A3nAj9mQZezMRar0sDHL+AZ2B5E4LQtXLO90p3DlaXYhdgc5hk0M4srSCObDbLYzncINmez4btgdRGJ1Xzfiag07O9EpzKmgCdgd5B6WL+JKdwcysXvHO4VVetiEs4NI5N4KU6zoRKdwRTr76gCFA3YHOYatDOLQkvCOdwo3ZLIT/Jie2B1EYjaua8c7hdkVNLGAfa43ttKcwyJGHOruQCZ1Z9Zn3nWn0NQdfCMUv59I9BaHKZZ3aE8hdgctBhsaxK3F4cxHOWzZXXYK12ey92N3EEnC2K7EzQZ+uHp3ncLsCnoIu4OWgqWMuNXdgTzoz2y4m05hZSP8L8u4CLuDSCqWhN/1TCF2By0J2xrEucVhd9cpXJ/JPoDdQSQh93Yh7jYQ3+5OYVY5PVKIi0UtBwsacc7PgTzk396ZwspG+CAbu4NIapaEK1a2u1O4Ip1dGKKwU3KcJ/Q3bG6QJbzV7k7h+5nGid2ZYOwOImm5pwtxt4EdV+7cKcwqp0eL2OewO2hBWNbIEvwcyMP+zPt36hRWNsJHOSzuHUSStCRcsfL0nTuFy9PZVwdgd9CiOCzsM2fO/P777y4uLlOnTnV0dLw1QUlJyfbt2+vq6iZNmhQcHGx68/fff6+qqjL97OnpOXLkSADQaDTHjh1r+sWhQ4f6+Phwl3PEhbfCmPBdhpf6KTxs2kzzfqbxAT8mwBG7g0iC7ulCPG1g+xV2amCbt3pZ5TS5iG4ZobBkxhBXt9779u2Li4urqqr67bffYmNj6+vrWyQoLS2NiIhITU3VarVRUVFpaWmm9xcsWPDRRx/Fx8fHx8cnJSWZ3kxLS3vhhRfi/1ZYWMhRthF3/BzI5B63mynE7iCSvCXhiuXprLHtTuGydHbhAAa7gxbGVXmvWrXq3Xfffe6551iWjYqKio+PnzFjRvMEn3/+ef/+/bdu3QoAdnZ2q1evjo+PN/3T22+/PWTIkBYX7NWr144dOzjKLbKMN02dwv4Kz9Y6hevOGSdidxBJ2pguxNsWdrTRKcwqpylFdCt2By2Ok7vv2traY8eOTZgwAQAYhhk/fvz+/ftbpDlw4MD48eNNP0+YMOGPP/5o+qeEhIQtW7acOXOmefqKioovvvhi9+7dlZWVXOQZWYCfA5kS0PqeQuwOIplYGqFY1kancFk6+yp2B/nASZFrNBoAUKvVppfe3t4pKSm3pvH29jb97OPjU1VVVV1d7eDg4O/vf+nSpYsXLy5YsODJJ59ct24dAKhUKm9v79TU1Nzc3Oeee27fvn1hYWGtfnR+fn5WVlZeXp7ppUqleu2112xtbVtN3NDQYGVl1ek/V146WWgL+0L0T+T5nkYP65tagv+eIQ90A1+rxgYzPddeUPCb1gFSLbSh7uBtQ7453/B4j5vez6qA5ELms8FsZ6qAVAutM5RKpUJxh042J4GQEAIAlP7V0rEse2s+CCHNEwAAwzAA8Ouvv5revHz5cv/+/Z944onQ0NBx48aNGzfO9P78+fNff/3133//vdWPViqVdnZ2rq6uTS+VSqXpyrdiGKatf0Jt6WSh+TnA5B7wv1xYGf7PEGhFI3x2AY5NAIaR5rgoftM6QMKFtjgU5hwnU3qAstnf9/YZWNgfHKw69SdLuNA6zBSPbo+TQOjr60sIKSws7N69OwAUFhY2df6ap2la81JYWOji4mJnZ9c8QWBgYGBgYG5ubmhoaPP3x44du2vXrrY+Wq1WBwcHz507tz35VKlUKhUeYXR3Ol9ob0XQsB8NLw9QNc0UbjxrfMgfgl0lOzWC37QOkHCh3esHXc4ZfrzOTAv6K2hlltMTJcbvRilVnWuSJVxonOLk3sHW1nb48OF79+4FAKPR+Msvv5j6cw0NDefOnTP1/8aNG2dKAAB79+41JTD9k0leXt6VK1dM2yqav5+YmNi01wKJUTd78lgA8/65v2YKKxrhkxwWj5JBsrIs4v/bu9eYKM41DuAv7gqiZMVdKHvhssTTcqq1ClKs2tK1h7WAKMGFYsORchFQUVqrjfHYY7DVpBK0tWnVkp7YixZSi0qhAgtKAbUoXkqtaECBgLsFKReBBRZ2Z86HydljAK/RnWHn//s0mEfzz5NlHnfed2YE2y9Tpv+d2LZdpDbPnuCI1UGWPK3Gb9u2TaPRNDQ0XLt2TSgULl++nBDS0NDw4osvdnV1OTs7JyQk7N+/X6PRyGSy77//vqysjBBy6dKlpKSk+fPn0zSdm5sbGxvr5+dHCElISOjr6/P09Lx27dr58+eLioqeUmywji1zJsw5atowS+A6iey+Yl6unOCNzaLAJ4tkdvLJJKeB+uffJvzRRZ+7TR9W2ewVEe77/0LdE1dfX19cXOzi4hIeHs5sV+nr6ysvL3/jjTeEQiEhpKen59ixY/39/WFhYR4eHoSQoaGh06dPX79+XSgU+vn5+fv7M/9US0tLZWVle3u7TCZTq9WWJcDR0tLSHv7SaG9v75h3+sN9PKmmpZ41T6GNqzx7XimXnAsX2vYgxCftMdh800630gkVptIF7Wl/SP7hMXH9zCdwUcTmm/aUPMVByAoMwqftSTXtP0cKkjdudRDLRA6C5spjtr3VDZ+0x2DzTWtvb/d6Zamdk4uxQ3f22NcBvrMf/HcexOab9pRgYQbYsWPHDur9UwPr8nsUc4/ee/cTgK3atXffoGp9/5pj5oSv3//3Drbj8BoGIbDDbDYToT0hhLKfPDhoizcPAtzXoHGI2E8hhBCHKYM2ef/s+IFdSsCOf21I3bY3hJb+3aXzeqQGu5+AdzalJh0PWW68UTyhqXrXgT1sx+E1DEJgx+rEuGUhi/V6/ezZs3HnE/CQUqmsu3j6ypUrzz67UywWsx2H1zAIgTVyuVwul7OdAoA1kydPnjdvHtspgN9rhFqt9u5b9eFh4CbOx4CmPSqz2Tz6Sf1wf4ODg8wN2fCoeD0IU1NTb9++zXaKcWblypVGLOw/CoPBEBcXx3aKcebPP/9MS0tjO8U4U1dXt3XrVrZTjEu8HoQAAAAYhAAAwGsYhAAAwGu2tmv0xo0bWq3W8l6L++vu7o6Ojrbth3s9ccPDw6GhoXjn2cOjKMpoNKrVaraDjCdGo7GjowNNeyR9fX1NTU1o2ggRERFr1669f42tPWu0urr61q1bD/m0vcbGRm9v7wfXwV3QtMeApj0qmqabmprQtEdCUVRLSwvzFliw8Pb2nj59+v1rbG0QAgAAPBJc4AIAAF7DIAQAAF7DIAQAAF7DIAQAAF4TpKens53BGjo7O48fP15XV+fl5TXm/RImk0mr1Z45c0YsFk+dOtX6CTnIYDDk5+fX1NQoFApHR8fRBTdu3CgpKbl+/fq0adPwXmxGU1NTXl5eW1ubt7f3fW4yqa2tra2tVSqVVozGXQ0NDXl5ee3t7d7e3nZ2dmPWXLp0qaioSKfTubm5TZo0ycoJOai+vj4/P7+zs1OpVI7ZNJ1OV1hYWFtbKxaL8ev5ADQPNDQ0uLm5RUdHh4aG+vj4dHZ2jigwm82LFy9+6aWXEhMTxWLxqVOnWMnJKV1dXT4+PsHBwStWrJBKpTdv3hxR8MUXX8hksqioqOXLl4tEovz8fFZycopWqxWLxYmJif7+/iEhIRRFjVmm1+tdXV0lEomV43HTzz//LJFIVq1a5evrGx4ePrqAoqjk5GQPD4+VK1cuXbo0MzPT+iG5Jjc3VyKRJCUlzZo166233hpdUFRU5OzsnJiYGB8f7+zsjHPa/fFiEK5ZsyY5OZk5Dg4O3rVr14iCwsJCpVLZ399P0/S+ffsWLlxo7Yjck5GRoVarmVP5mjVrUlJSRhQ0NzcbjUbmODMz09fX19oRuWfevHkHDhygabq/v9/T07OkpGTMsoiIiE2bNmEQMubMmXPw4EGapnt7e2UyWUVFxYiCQ4cOTZ8+vauri4VwnERRlI+PT05ODk3T3d3dEonkwoULI2rCwsLS09OZ4y1btmg0GmunHFd4sUaYn58fGRnJHGs0moKCghEFBQUFS5YsYa7+RUZGnjlzprOz09opOaagoECj0TCXXCIjI0c3zcPDw3KRWSaT4ZUUt2/fPnfunEajIYQ4OjqGhoaObhohJDs728HBYdmyZVYPyEUtLS01NTVM05ycnIKDg0c3LTs7e/Xq1R0dHadOnero6GAjJrfU19c3NDSEh4cTQqZOnapWq0c3TSKRGAwG5ri/v9/FxcXaKccVW3vE2mgURbW2tioUCuZHhUKh0+lG1Oh0uoCAAObY1dXV3t5ep9Px/J3ROp3u7qa1traazWaBQDC6cnBwMCMjY9WqVdYNyDl6vd7BwcFyxlEoFDU1NSNq/vrrrw8//PCXX36pq6uzekAu0uv1IpHIsoKlUCgaGxtH1Ny8edNgMBw5csTd3b28vDwnJycoKMjqSTlEr9e7uLhYFkrHPKft2rUrJiYmJCTEbDbb2dkdPnzY6jHHE9v/RkhRFEVRlsVkgUBgMplG1JjN5rv3NYxZwzd390QgENA0bTabxyyLjY1VKpV4e9zDfIrWrVu3efNmNzc360bjLuY0bflxzKYNDg4KBIKqqqrc3NyPPvronXfesW5GznmYpp04cUKn07355psrVqxobGwsLi62bsZxxva/EQqFQldX1/b29ueff54Q0tbWJpfLR9TIZDLLG3p7enoGBgZG1/DN3T1pa2tzcXEZvduWoqj4+Pg7d+789NNPY35Z5BWpVDowMNDX1+fk5EQIaWtrk8lkdxc0Nzfn5eWJRKJff/21tbXVYDCkpKSkp6ePKOMVqVTa09NjNBodHBzIWE0jhMjl8sDAQObUr1Kp1q1bZzKZhELbP3fdi1Qq7ezstFyhaWtrs1y8sfjggw+ysrKWLFlCCJk2bdqmTZtiYmJYyDpO2P43QkLIokWLtFotc6zValUqFXPc0dHBfMtRqVTMvgZCSHFx8YwZM/B/dpVKZWlacXGxpWnd3d1DQ0OEEJqm165d29TUdPToUeYsxnNyufy5555jmkZRVGlp6aJFiwghJpOJWdkSi8XffvutWq0OCgqaO3fuxIkTg4KCpkyZwnJuVimVSi8vr5KSEkKI2Ww+efKkpWmWdfrXX3+9vr6eOa6vr5dKpXyegoQQHx8fsVhcXl5OCBkeHi4rK2OaNjw83NXVxdQIBALm95QQYjQa8f/UB2B5s45VXLhwQSQSpaenv/feexKJpLm5maZp5lNSXV1N07TRaJwxY0Z0dHRmZqabm9vhw4fZjsy+5uZmiUSyYcOG9PR0kUjENIqmaXd39x9++IGm6U8//dTOzi4mJiY5OTk5OTk1NZXVvJzwzTffSKXSzMzMqKioF154YWhoiKbpqqoqQsiIWykqKiqwa5SRlZWlUCh2794dERHh6+trMplomi4vL584cSJToNfr5XL5xo0bP/nkE3d39y+//JLVvJywd+9eLy+vPXv2hIWFzZ8/n/l0FRYWikQipmDnzp0KhSIzMzMjI0Mqle7evZvVvFzHl7dPXL169ciRI/b29jExMcxrSmia/uqrryIiIpjdDd3d3QcPHuzo6Fi8eHFgYCDbeTmhubn50KFDQ0NDUVFRM2fOZP4wOzv75Zdf9vb2rq6uvnz5sqVYKBQmJCSwlJRDysrKSktLn3nmmbi4OObJDO3t7cePH09KSrq7rLW1VavVxsbGshSTW0pLS8vKytzc3OLj45mNM62trSdOnLB8ovR6/XfffWc0GoOCghYsWMBqWK4oKiqqqKhQKBRxcXHMdYVbt26dPHny7bffZgrKysoqKyvt7OxUKtWrr77Kaliu48sgBAAAGBMv1ggBAADuBYMQAAB4DYMQAAB4DYMQAAB4DYMQAAB4DYMQAAB4DYMQAAB4DYMQAAB4DYMQAAB4DYMQAAB4DYMQwHakpKR4eno2NTUxPw4NDQUGBgYEBAwMDLCaC4DT8KxRANvR29vr7+/v7OxcWVlpb2+/YcOG/fv3nz171s/Pj+1oANyFQQhgUy5evLhw4cK0tLTXXntt6dKle/fuXb9+PduhADgNgxDA1nz22Wfvvvuuk5NTUFBQbm4u8253ALgXDEIAW3Pnzh1PT8+enp7ff/991qxZbMcB4DpslgGwNatXr54wYYKHh0dqaqrJZGI7DgDXYRAC2JSsrKycnJx9+/b9+OOPVVVV27dvZzsRANfh0iiA7bh69WpAQEB8fPznn39OCPn444+3bt1aVFSkVqvZjgbAXRiEADbCYDAEBAQIBIJz5845OjoSQmiaXrZs2fnz53/77TeZTMZ2QACOwiAEAABewxohAADwGgYhAADwGgYhAADwGgYhAADwGgYhAADwGgYhAADwGgYhAADw2n8Bo1zlHESSWWMAAAAASUVORK5CYII=", "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 -> U and K -> Γ -> L -> W -> X\n", "\rDiagonalising Hamiltonian kblocks: 14%|██▎ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 93%|██████████████▉ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 100%|████████████████| Time: 0:00:01\u001b[K\n" ] }, { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=44}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOydeVxM6x/HPzNtWmm3hUpJUloQskRcW+hSWVMke+laKnvWrD+pZMsu5GbLRShcaxKyRUIushRp32bm+f2RpaJ9zjnTdN6v+7qv6XHO8/2cMzPnO8/zfJ/vl0MIAQsLCwsLS32Fy7QAFhYWFhYWJmEdIQsLCwtLvYZ1hCwsLCws9RrWEbKwsLCw1GtYR8jCwsLCUq9hHSELCwsLS72GdYQsLCwsLPUa1hGysLCwsNRrWEfIwsLCwlKvYR0hCwsLC0u9RjwdYUhIyKNHj2g2euAAJk+m2SbdCAQCpnLy3b9/v6ioCMCePbh/nxEJpTh7Fv/8U6olPx9z5zKk5jtr1iAlRQj98Pn8WvbA42HOHCEoEX1qf6+KSUpCYKBQeqIcYV0yPVTlkSWejvD8+fN0OkIej7dw4cLp05cnJb2izSgjFBQUMPUd+OOPP758+QJgwwbweIxIKIWiIjw8SimRkcHZszh/njlNwJs3CAkRQj+5ubm17CEmBtHRQlAi+tT+XhVz7BiSkoTSE+UI65JpQCAQ8KrwvBBPR0gzrVt3XbXqU2am6uXLvZ4/f860HHEmNRVv38LUlGkdgJUVWrTAgQM/WzgceHlh9WrmNAHOzti9G6KQSP/CBfTrx7SIOsXly+jZk2kR9RXWEQqBN29SCdkOTCPELTg4mGk54szly+jeHRISTOsAACxfjhUrSg0KR4/G27e4fp0xSRYWUFTE1auMCfjBhQvo25dpEXUHHg83bqB7d6Z11FdYRygUpIH/AD5w1cjIiGkx4sylS7C2ZlrEd7p1Q/PmCA392SIhgb/+gp8fc5oAJyfs2cOkAABZWXj4EF27MiyjDnH3Llq2hJoa0zrqK6wjrC2hoQCCgD6AtoREUqdOnZhWJM6IlCMEsGQJli8vNSicMAF37+LuXcYkjRuH48eRlcWYAADR0ejSBbKyTGqoW1y5gl69mBZRNZ48ebJkyZKTJ08yLUSYsI6wVmzfjrlzISeXCDwH/pOUnPL48WOmRYktHz/i0ycYGzOtowTW1mjWDIcP/2yRkcGsWVi3jjFJGhro0QPh4YwJADsvWn2uXKkbC4RXrlxp377/pk1Kdnarxo6dwrQcocE6wpoTFITVq3HlCqSkenE4p4B/5OUPW1lZMa1LbLlxAz17gitin9nFi7FiBUrG0k6bhuhoJCYyJsnZmeHZUdYRVgs+HzduoE48ORYuXCcQbCBkKSGXwsIimZYjNETsoVJ3WLMGAQG4ehWSkiCkrZLSGwuLGwMH7mvWrBnT0sSW69dFa160mN69oa5ealAoL4+pU5kcFJqavr15015Ly3LJkrX0W//vP6Sno317+i3XVe7dQ7Nm0NBgWkcVUFPTBF4AAJKlpaUYViM8WEdYE5Yswd69iI5G8+bYswctW3JmzZp+7NjKc+cM8/OZFie+XLsmio4QwJIl8PUtNSh0d8exY/jvP2b0jBkzvbDQ/e3bq/7+j8+dO0en6bCwY/b2M1u12g0I6LRbp6krC4Tv3+Pp0yAO5zyHYwCMOXRoA9OKhAbrCKsHIfD0xOnT+PdfNG0KQrB/P16/hosLtLRgaorjx5mWKKbw+cjKQrt2TOv4HTY2aNwYYWE/W1RU4OKCTZuY0fPmTQrQHZDKzOwXH/+UNrthYccnTz54+/bohw/jVqz4H2126zp1YoEwIQFdugBosGRJdGxsHIdzr2dPW6ZFCQ1JpgXUJQQCTJ6M+HhcuAAVFQC4cgWFhejaFS1bAsCkSdiyBaNGMStTDAkNDc3IyNLXj+BwXJnW8nsWLsSsWXB0/LmEOWcO2rX70qlTjKVl21atWtEppn//XqGhs3NyesjLbx4yZA9tdo8fj/r6dR7QOT9f78SJcYsXz6bNdN2Fz8e1a9i5kyZzSUlJGRkZ1T3r0SPMno3u3fH4Mf74A7m5ubKycvPnw8WFCo0Voaenp6SkJPRuWUdYVfh8TJyI5GRERUFR8Vvj3r2QlsakSd/+HDIEM2bg+XPo6TElUwwJCNji4bGVEOmEhKPDh98ID9/FtKLf0K8fVFVx9CgcHb+1ZGcn5efbOzsPVlRcGRw8d8SIobSJCQ7269btcFTU/XPndrZq1ZY2u927mx47FlpYqC8js7drVzPa7NZp4uPRtCl9C4RWVlaNGzeWlKzGk5/HQ14emjXDkyeQl8eMGSCEtGzJiYrCrVvUKf0N79698/T0nDdvnvC7JuKIg4PDoUOHhNhhQQEZPpwMGEByc382ZmcTJSXSuDEpKvrZOGcO8fYWomXRIjc3t6jk1dKCnl534DKgCTxr0ECXZutV5+xZYmhI+Pxvf3p6LgZOAQRIbd/ehhFJtrYkIKAmJ2ZmZtboLL6Cwmo9PZtp07xzS35VxJqa3asfbNxIpk0TlpbKUVNTS01Npc+eUJk/f/7KlSurdQqfzy8sLKz0MHaNsHIKCmBvj8JCHD9eao9wWBjU1TFhAkr+upo0CXv2oKiIfpliS9u2WsAFAMBVNTVlhtWUT//+UFb+uYFPRUVJSuodAOCdsnJDRiQtXYrVq5GXR5O5rVu5gwZ5JyZeCApaLctup68adWKBUOxhHWEl5ORg8GDIy+PYMcjIlPqnPXuQllZ2llxfH3p6ZQv0sNSGQ4d2aGhcBNIaNVoXFXWQaTkVMX8+fH0hEADArFmTTUxOqKl1lZR0mzx5OSN6zMxgZoZdtMwl5+Rg40YsWECHLbFBIMC1a6wjZB7WEVZERgb69YO2Ng4cQJlJ9eRk3LsHCwu0bl32LFdX+pa+6wNycnKTJ9+Sl1d9+vSKvr4+03IqYuBAyMt/ixxWUFCIjT339u2l/ftj/PzaMjVJ4OuLVatAw66eoCBYW7PbB6tHfDw0NKCpybSOeg/rCMslPR39+qFDB2zb9ptsJnv2oGHD31fitbdHTAxje8jEktu3y/4QEVkWL8aSJd8GhQBkZGRGjkSLFvD3Z0aPmRlMTSkfFObk4H//w/z51FoRJ3g83tixM62tLdLTh6YIpZgySy1gHeHv+fgRPXvC2hpBQeBwyv4rIdi1C7m5GDLkN+fKysLRkfkKAGIDIbhzB1J1JIvFoEGQk0OZjMSBgVi7Fq8YKtvs6ws/PxQUUGgiIAC9e4OtvFJ1QkL2HT+ulJFx5+PHma6uXkzLqe+wjrAURUVFRUVFb96ge3cMH15uPZ1Ll5CfD1fXsquGP3Bzw86dYKiWu7iRlARFRZFLMVoBixZh8eKfg0IArVrBwwPu7szoMTdH+/YUDgpzcuDvj0WLqOpfLElISM7N7QaAkG4vXiQzLae+U3eeLtSzcqV/06YdGzfuZGKybupULFlS7pF79qCwEBMnlnuAsTE0NXHxIhUy6x23b6Nu1baytYWsLE6dKtU4dy5evmQs8dCyZVi9mqpB4ebNsLGBgQElnYsrw4YNk5DwBfY3bOg8YYID03JEhYyMjNQSFNG1tM46wm98/fp106ZDaWl3v3y5W1R0evToj+UdmZ2NY8fQrh0qjttgQ2aERWwsOnZkWkQ1WbAAvr4g5GeLtDS2bsWMGah+Wg8hYG4OIyNKpuuzs7F5MxssWj0EAqxZY+bsvHvNms/h4W5eXjOZViQq2Nvb6+rqGn/nLl2FPVlH+I3c3FwuVxXgAhw5ObXs7OzyjgwLg7w8pk+vpMMxYxAVhQ8fhKyzHlLnRoQAhgyBpCQiIko1du+OP/7A0qXMSFq+HKtWobBQyN1u3oy+fdnhYPVYvBgFBdi61XDevFl9+vRhWg6TFBQUZJWuIr106dL33+ncuTM9MlhH+I2mTZs2bNhIWtpFVXWSkRFXR0envCO3b0d+Pv78s5IOFRRgZ4f9+4Wss75RVISHD2FW19J1cTiYPx9Ll5YaFAJYvx5HjjBTv97cHG3bCnlQmJ2NgAB2OFg9IiKwfz+OHKkzgdAUkZubO2zYsKZNmxobGzdp0uTLly8Miqnfb0UJEhKQkXEgLCxGQ4NYWlpyfg0VBQC8eoVHj+DqigYNKu/T1RXOzpgz5zdxpyxV5MED6OhAQYFpHdVn4MBCF5epqqr3dXU1T5zYUVyoUkUFq1Zh8mTcugUJCbolLV0KR0c4O0NaWjgdbtqEfv3Qpo1weqsPJCVh4kQcPw51daalAADeZb0Ljg3mCXhUG5KWkHYxddFupP2j5dKlS69evUpJSZGRkcnNzZX5Hnl46tSpt2/fAlBVVV1A148s1hECAJ+P8eOxahVn6FDLio8MCQGH8zPLdsV06QJZWVy9ih49hCCyflIX50WL2bo1JC9Pp7AwJC4u2s3N559/9hW3jx+P/fsRHIwZM+iWZGkJAwPs3VvVD3DFZGZi82ZcvSqEruoJeXlwcMDy5ejWjWkp37n3/p7fNT8+oTzAncvhtlVvW9IRmpubS0pKBgQEtGzZ8s8//5T4/sNQWVm5uFpLw4b0JSZkHSEArF8PRUVMmFDJYYRg507o6VWjKp6zM3buZB1hzYmNhWUlP05ElOfP3xYW9gBASMfk5J8p1jgcBAfDygp2dmjWjG5VS5di5EiMHy+EQaG/PwYNYoeD1WDqVLRr9/ssHEwxWH8wbzHlw8HfQghp1KjR3bt309LShg0b9sMRdu/e3Z32nUbsGiGePcP69dixo/IJzOho5OVh1qxqdD5uHCIikJ5eG4H1mro7IpwwYYSa2kIOZyeHM8befmzJf9LXx+TJ8PRkQFWXLmjTBvv21bafzEwEBMDHRxia6gf+/rh/H9u2Ma1DZNi7d6+urm5oaKifn58U0/ky6rsjFAjg6oply1B+cMxPgoPB42HEiGr0r6qKAQNwUKQzRYsuWVlIThbRqvSVYmZmeuNGaGAgb+pU76ioiSX31wNYuBAPHuD0aQaE+fpixYraho9u2oTBgyvZQcTyg5s3sXo1jh2DnBzTUkQGExOT8+fPR0RE3Lt379ixY0+ePGFQTH2fGt20CZKSmDKl8iOzs3HmDEaPrvZH2dUVnp4MLAiJAXfuwNS0ziRX+xU9PT09PT2BANbWCAwslVlGRgYBAZg0CdbWkJenVVWXLtDTw/79FWWEqJiMDAQG4vp1ocoSXz5+hKMjdu+u0q9tqqE00161GDBgQEBAQHh4+OfPn5s3b25iYgJg+PDhBkzsxanXjvDVK6xZg+vXqxTVefgwOJya+DNra+Tl1cld4YxTd+dFS8LlYs8edOqEPn1KjW779kXXrlixAqtX0y3J1xejR2PcuBquFP7vf7C1hZ6esGWJIzweHBwwaRIGDGBaChAVJVp5H21tbW1tbUu2TGZoBbX+To0KBHBxwfz5v6mj9Fs2bULz5ujQodqGOBy4uLBZZmqC2Px60NaGry/Gjy9bsfl//0NICOLj6dbTtStat67hjH1GBoKD2b2DVWXOHCgpicTtysjAxInlpkeu59RfR7hlCwoLqzrCS0pCUhJmz66hLRcXHD2K0vkTWCpHPEaExUydCg0NrFtXqlFTE76+mDwZZVYQaWDZMqxcCV71AwY3bsSQISIxyyf6HD6MkyexZ49IpIyfORO2tgzsXq0TiMD7wwSvX8PXFzt3VvVjUVyMadSoGppr3Bg9eyIsrIan108+fEBODrS1Kz+yTsDhYMcO+PsjLq5U++TJ4HJpKiJfkq5d0aIFDhyo3lnFw0E2WLQqPH0KDw+Eh0NVlWkpwMmTuHGDgUn4ukJ9dISEwM0NXl4wNKzS8QIB9uzB4MFQVKy5UTYHd3W5fRudO4tVUp5mzbBuHZycStWL53KxbRsWLMCnT3TrWb4cK1ZUb1C4YQOGDWOHg5WTlYU//8SaNSKRHTA1FVOnYvfuOpmhiR7qoyPcsQNfv1ZjF1dUFPLy4O1dK6MDBiAlhYHVoLpLbKz4zIv+wMkJhoZYtqxUY/v2GDcOc+bQLaZbN2hpVWOl8OtXbN0qisNBHo/n7DytdWurJUuWVX409RACFxf06QNnZ6alAACmToWzM7p3Z1qHCFPvHGFKChYvRkhINebK/fygpgZz81rZ5XIxfjx2765VJ/WK27fFJFKmDMHB2LOnbGayZctw7RqiougWs3w5li+v6qBwwwbY2YniZHX//qP27ct78WLV8uWXFizwZVoOVq/Gu3fYsIFpHQCAvXvx9CkWL2Zah2hT77ZPuLrC3R1GRlU9PjMTV69i1SohmJ40Cebm8POrUsLueg4hiIsTT0eopobAQLi64t69n3tS5eTg749p0/DgAa1xfVZWaN4coaFwcqrkyC9fEByM2FhaZFWTq1cTCLkMqBGyetMmPwMD2Nsz9i27eBFBQbh9W2iZzWvDmzeYNw/nz5e6G6dOnVKszTIPcyQkJFhYWFDRc/1yhHv2ICUFc+dW45TiLNuurkKwrqUFMzMcO4bRo4XQm3jz/DmUlKChwbQOavjzTxw/Dm9vbN78s9HWFrt3w88PS5bQKmbJEri5YfToSqoCbdgAe3vRGg5++QIfHxw+jMJCG+Aw4AKESkpauLlh/Hg0bAgjI9jYwNkZLVvSJOm//zBuHA4fZiCL7K8QggkTMGsWTEx+No4ZM2bhwmX5+R0BDvBBV7dAW7vad4fP50uUM6V25QqaNKEq5ZCkpGQnitZLiDji4OBw6NChMo0pKURDg8TFVa+rVq2IjY3QhB09SqythdYb/eTm5hYVFdFgaP9+4uhYqkVDQ+PDhw80mKaHr19Jixbk3LlSje/eEXV1kpBAt5iePcm+faVaMjMzS/6ZlkZUVcmrV3SKqoh//iGWloTLJSoqxNOTJCd/bN26q7R0S0vLQcWfzw8fSFAQGTCAaGoSDodISREdHTJ2LDlxgvD5Qhbz417l55OOHcmGDULuv8Zs3ky6dCE8Xtl2Q0Nb4DYgUFSc/vff4TXouczHoyQjRxIdnRp0SRV8Pr+wsLDSw+qRI7SzI0uWVK+fxEQiKUmuXhWWLlJQQDQ1SWKi0DqkGdocobs7Wb++VIuYOUJCyIULREuLfPlSqnHjRtKzJxEIaFUSFUVatyYl39gyTzofHzJ1Kq2Sfkt6OvH0JCoqREKCWFpW9YtZUEBOnCBjxxIdHSIlRbhcoqlJevYka9aQ9HRy+vRpKSktLldbWblNRkZGDVT9uFeurmTYMLrfu/JISiIaGuTZs7LtZ88SZeWX+voDmjUznzLFS1AjuRU4wlu3CIdDCgpq0CslsI6wlCM8eJC0b1/tt2fMGKKmJkxhhJA5c4i3t5D7pA3aHKGlJfn331It4ucICSFTphBn51ItPB4xMyP799OtpEcPcuDAzz9LPumKh4PJyXRLKklYGDEyIhwO0dQkixfX6jl76xaZMYOYmBA5OQIQDscGeAgQYPWwYeNq0GHxvdq5kxgakqysmgsTIjwesbQkAQFl20+fJpqa5MaN2vZfgSMkhMjKkuDg2poQFqwj/OkIU1NJkybk9u3qdcLnEzk54TutZ89I48akCm+NKEKPIywsJAoKJDu7VKNYOsLsbKKrS06eLNV4+zZp0oR8/kyrkosXiZ7ez0FhySedlxeZPp1WMT/48IGMH0/k5YmkJLG2JnfuCL9/CYk2QBZAgHDAV0aGNG1KLC2JqyvZtatK7j8zM/POHaKhwcCcdnmsWkVsbMqOTcPCiIYGiYkRQv8VO0IrK2JpKQQrQoF1hD8d4YgRxMen2p0cPUokJCh5HnXvTo4fF363NECPI4yNJSYmZRvF0hESQq5fJ40bk48fSzVOn04mTaJbSffu5ODBb69/POnS0oiaGnnzhnLrRUVFI0Y4q6ubOTi48Pn8U6dKrQLm5VFld+zYKRxORw5nKZerdfbspcuXycKFZMAA0rYtUVIiHA7hcIicHNHRIQMGEE9PcuLEz/Ho+fPnZWV1udyW0tKeR44Ie+2xpty/T9TVyevXpRpDQ0nTpiQ+XjgmKnaEO3aQBg2EY6j2sI7wmyM8cYK0aVOTL5KZGenUScjCitm3jwwcSEnPVEOPIwwK+o0bEFdHSAiZO5eMGFGqJSODNG78bsyYZYsX+32ma2x4/jzR0/sWW/HjSTdvHpkxgw7rQ4aM5XBcgIeAi6TkRklJYmND7t6lw/SJEyfmzp37+PHj3/7rnTtk3Tri6EhMTIiqKpGUJACRlSUtWhAutxNwHxAAExcsWEiH1srIzyfGxmTv3lKNO3aQpk3Jo0dCs1KxIywqIlwuuXxZaOZqA+sIDxFC0tJIkyY1iXbJyCASEiQiQvjaCCG5uURVtexPtjoBPY7Q2Zls3162UYwdYV4eadfu52iMEJKXl6emZg4clpTcZWDQjTYl3buT0FBCvj/pUlNpGg4SQpSVjYFEgADPZGW70bISXUM+fiQHD5IZMwiH0xYQAATY3aGDQ5nJfEbw9iZ2dqVagoNJy5bk+XNhWqnYERJC9PTK/rZjiio6QjHPLOPujtGjYWVV7RNXroScHAYNokATICuLkSOxZw8lnYsB4lR0oio0aID9+zFrFt68+daSkJDA4ZgCjjyeS3q6akpKCj1KFi7E8uU/S2GsXYtRo9C8OR2mtbU7ANuBt0Bgz57NKt7UyCwaGhg9GgEB6NChOYczFtjN4Sxt1GhaixZwcsLFi4wJu3kTe/YgOPhny/r1WLMG0dFVLTYnLIYNw6VLtFqsLTT4ZPopHhFGRBAdnbIxF1VEXZ24ugpbVgni44mW1m+2+Ig4NIwIMzOJoiL51YgYjwiL8fUlfft+C3D48uWLhoYF8AF4JS/fgZ5I3WKsrMihQyQzM7N4OPj2LR1G794lkpK8xo3dFBQMevf+s0B0ou8rhM/nL126dMCAP69evUoIef+e+PkRXV1iaEj8/EhaGq1icnKIvj75+++fLX5+pE0bSgb0lY4IP3wgHA5JSRG+6epS36dGQ0L+1tIiV67U5PQ7dwiXS/m7aGFRdj+16EODI4yKIt27/6Zd7B1hURHp1Ils2/btz3PnLhgZ9TYx6W9gEDN3Ln0yzp0jbduSr18zZ88mHh50WPzvPyIrSwYNosMWFZTxCgIBuXqVuLkRZWVib08uXKBpZ+HUqWT8+J9/LllC2rYl795RYqtSR0gIUVUldH5uy6NeO0Ibm5ENGpjZ2dXIDRLSrx8xMBCuot+wbZuoTKNXHRoc4erVZPbs37SLvSMkhCQkEDW1sss5nz6R9u2rnQuixhQWFqqoOCkqmklK9r5x4yXV5j5/JkpKpHNnqu1QSHleIT2dbNtGjI2Jvj7x8ysbGCxcSiZnEAjIX38RU1OSmkqVuao4whEjiJ4eVQKqTr1eI/z61To//+zt23/V4FyBAJcu0VETZ/RoREXhwwfKDdUtYmPFM9d2VTAwgI8PnJ3B5/9sVFdHVBTCw8sWb6KI7dt3Z2frZWXF8Xir5s+fR6mt3FwYGEBVFTduUGqHGRo1gpsb4uMRGoqXL9G2LRwcEBFR6s0VChkZmDgRO3ZAWRmEYNYsXLmCCxegpiZkQ9XC0xMvXlSv2iWDiKcjzMtTBTSKijiCH+v+VSYwEFwuJkygQlcpFBRgZ4f9+yk3VLeob5EyZZg1C5KS8Pcv1VjsC8PCsGIF5QKSkz8UFhoDAIxSUj5SZ0gggIkJuFw8eQKueD6HvmFujm3b8OoVbGywdClatYK3N1694js5uTdvbmFlNez9+/e16X/GDAwdij/+gEAAV1fcvYvoaKiqCkt+DenaFVJSdSckkIbBKf306TNNTs7VwWFytc7Kycnx9vZu1GjH4ME5FAkrw40bRF9fVJITVgWqp0bfvy83p119mBot5sULoq5Onjwp2/7xIzE0JCtWUGv98ePHGhpmEhIBKiqDNm/eQZ0hMzOipER3RAkVVGWesCSxsWTyZCInFyIh4QMQDuf8gAE1Se1WTHg40dcnOTmExyPjx5OePUk15dSEKl6ypSWxsqJaSyVUcWqUmSDl/Pz84jirHj16NCinbtiHDx/evXtnYmIiWSKSOjMz8+rVq4qKit26dSuvDggAJaX37u5tVlWziqCGhlFOjh2QfPFiex7vmST1Edzq6kmvX09QUcnp2FHv1Kk95d2K+kNMTL0eDhajo4MVK+DkhLNnv8rISPwoHaehgeho9O4NDgfz51Nl3dDQMDb2+D///GNhsaQjZZPUffviyRM8fcr8wIV+LCxgYQEJieQtW6wAEGJ17twyAwNoa6NVK7Rq9fNFBWXIBALB9evXs7Olpk+3PH4cUlIYORLZ2Th7FrKy9F1LxYwfj79qsjzFBNS75LKkpaW1bdu2Z8+e1tbWBgYGqb8s6b5+/bpZs2bFXuFjiSXmp0+fNm7ceMiQIaamptbW1hXEWP+2+kTF3Lp1i8vtBRCAcLmDIiMjq3V6zejRYzhwByAyMus2bAikwWItoXpEuGBBuVEh9WdEWIyW1lx5+W7q6p0WLVpXsv3DB9K2LVm9mlrr1R3lVIsxY4iUFLl3jzoLtFKzexUXF6eq2gnY37Chw4oVm1+8IBcukG3biJcXsbcn3boRHR3SoAHR0SE2NsTNjfj5kbAwcucOefeOCASC7t2HqqhMkZZ20dUdU1BAhg0jgwZRmIiuDFW85IICwuWS69epllMRojsiDAoK0tHROX36NIBhw4YFBgYuXbq05AFqamrnzp3T0NDQ1NQs2b5s2bLRo0dv2LChsLCwY8eO4eHho0aNEpaqFi1aEJIC5ANc4IU2LRVIU1PTAH0ABQVt3r69Q4NFEef2bcyaxbQIEeDVq1e5uY9ycq7l5JDgYMu5cyf/GBdqauL8eVhbg8OBlxezMmuCtzcOH8a5c+jQgWkpjGJmZnblyq7Tp8+bm0+ysbEBoKNT9pisLCQn49Wrb/+/devba4HgeV6eDI8XDCAjo/8ff3zU0NAMC4OUFP3XURHS0mjVCgEB6NqVaSmVQoNPLoOpqenB7+mkDh8+bPJrfmVCCCFpaWkoPSKUk5O78z3/vK+vr729fXkmajAiTEkhEhI7OZzWXK6Wo2P1FhdrzLZte1VV+3G561t7V8oAACAASURBVBUVTRNEJ3d9+VA6IhQIiLJyuVHm9WpEmJCQoK5uXzw/ISfXMyWlbLrRN29I69bE358qARSNCDdsIFwuOXyYir4Zg9LR8295/PidgoIVIAB4kpKmI0Zk0JyRruqXPGsWUVenVEsliO6I8M2bN1paWsWvW7Ro8fbt26qclZ6enpubW/LEyMjI8g7OyMi4efMmh8Mp/lNGRsbW1raCzgUCdOwoYWQ0MS7OubiFL/QY598xceKYTp2Mjxx5cuJEhJ5eY3qM1gY+n8/lcn/cWOHy7BmUlSVUVcu9DXw+X/RvkVBo3bq1rOwrLncgkNeggXb37irbtgl69SI/DmjSBBcuwMZGAhBMn04q6KpmUHGr//6bM3cu93//E4wYQcTpbaT/Y9mkiSwhn4AOgIDDKdy7V5rDoVVC1S/Z0xP+/hIpKfzSs3v0UVx5uNLDGHCEhYWFP+JQpKSkCgoKqnJW8WElT8zPzy/v4PT09Js3b/5wsYqKin369KkguGb4cJmvX3H3bm7VtAgTfX39hQv1w8JkY2IKOnSo9mYPmikoKBAIBBR9527elDQ3R3mfB0JIYWFhFT8tdZ2nT//LzFQWCPwAKS7Xydf3q4tLw759BatXFykofPtWa2jgn384AwY04PN5kycLebuW0G/1pUvcMWMaeHgUuboWidl7SP/HMikpSU6ua07OdoCrrGz/5s3rHyMEeqj6JaupoVEj2U2bBEuXFlGt6rcIBIIKnvw/YMARNmnSpHjaE0BqamqTJk2qcpa6urqEhERaWpqKikrxiU2bNi3v4FatWtnZ2Y0cObIqPe/YgfPnceUKNDTkqnI8Fbi44ODBBqI/k87hcKSkpCiKp42Ph6Ul5OR+/y5wOBxZWdny/lVkyc3NrVhzURFevkRi4rf/nj9HYiI+fyY8nipgDKCgQLF/f/7gwZx58yTNzCSDgjB06Ldz27TB5cuwtpaWkZGeNk2Ysvl8vhBv9b17GDYMY8di40YpQMQWsmqNcO9VVTA2NgYeASckJAqVlN7r6+tTNElTHtW65J49ceqU1Nq1zLzvVfzhzsBG1m7dukVHRxe/jo6O7tatW/HrigewEhISXbp0+e2JteHRI0ydiuXLa1KhQoi4uODwYeTmMqmBccRsK31mZqapqY22dv+WLS0eP35c3JiejmvXsH07vL3h4AALCygpwdoamzcjMRE6OvD0xL//IjvbwMqqUFXVWVnZUVpaf+JEVR4P27bh0CHMmwcHB3z/JYkWLRAdjfXrS9UcECnevEG3bujfv+7srRZ5+HwZWdmT9vaP581LvnHjFM1esLq4u+P5c1Q/tQm9UL1W+Svx8fGKiopr1qxZt26dkpLS/fv3i9vl5eXPnz9f/NrHx8fd3R3AzJkzvb29ixtPnTqlqqoaHBzs5eWlrq7+6dOn8kxUMVimoICoqhJr61pfkjAYNIjs38+0iMqgLlimoIDIy1dUKqTOBcssW7ZOUnIbQICHamp/mpgQWVnStCnp1Yu4uZF168jJk+TpU1LeJiCBQHD79u379+8XFJA5c4iWFomOJoSQ3Fzi5UWaNClVfDU5mbRqRbZuFZp4YQWApKbW+VSilUJ/sMxff5XKr00/1b1kaWmybx9FWipBdINljI2N//3333379hFCLl++bGJiUty+YsUKfX394tfKysoNGzb08/MreaKtrW1YWNjx48cVFBRiYmLU1dVrqaRnTxCCc+dq2Y1wmDgRmzdj7FimdTBEfDz09CAvz7QO4ZGensXj6QEANOTls3btgp4evm+CqBwOh/NjP/u6dejbF+PGwckJvr7w84OjIyZORFgYtm5F8+Zo2RLnz3/ba+/mRsnl1IDcXBgaQk1NPFOJMsX9+zh4EA8fMq2jOpiYICQE48YxraMCqHfJDFCVEeHixURSkjx7Ro+iyikqIk2bksREpnVUCHUjwsBA4uZW0QF1bkS4desLCYkOSkreampdjh+PqH2Hnz4RW1tiYfHtQ1JYSPz8iIYG2bbtW5a+xESipVVqpFhjaj/K4fOJri7R1KRvlzdT0Dki5PNJly4kJIQ2g7+nupfs70/k5CjSUgn1uvpEpdy4gRUrEBCA70NQ5pGUxNix2LWLaR0MIWZFJy5dwuLFOpcuXfr779537x4ZNmxw7ftUV8fJk5g0Cd27Y/t2SEnBywtRUdi1Cz17IjERenq4eBELFuDAgdpbqy0dO+LTJzx6hHqfN1CYbNkCSUm4uDCto5pMmYL8fMTGMq2jAmjwyfRT8YgwPZ3IyREHBzoVVYlnz0jjxqQKP18Yg7oRYdu2JD6+ogPq0IgwPp40blzDotBV4fFjYmxMRoz4Vn+OzyfbthF1deLnR3g88vQp0dS8q67eSV3d1NZ2PI/Hq4GJmo1y0tPTNTVNJCRaSEpaSUunvKS8mqFIQNuI8P17oqFBHj+mx1pF1OCSW7UiY8dSoaUS2BFhuXTpAjU1HDnCtI5f0NeHnh7OnGFaB+1kZeHtWxgaMq1DGLx4gQEDEBiIHj2oMmFoiJgYNG0KU1NcuwYuF25uuHkTkZHo1g08HhQVZ6emHklNvRsd3TQs7ChVOn7hzz9dPn4cw+e/5vFmq6qOoCVNYT3C3R2TJtXVr8ngwbhwgWkR5VPvHOGMGXj5UnRX7ydOREgI0yJo5/ZtmJqC+moflPPpEwYOhK8vhg+n1lCDBvD3h78/RozA0qXg86Gri6goTJwIa2t8/JgBNAeQk6P77FkqtVJKkJiYDhRvhrXMyEir5GiW6hAZiXv3sGAB0zpqypw5+PQJnz8zraMc6pcjPHYMwcEIC0OzZkxLKQcHB9y8iZQUpnXQi3jsIMzMRP/+cHaGqytNFocORXw8YmJgZYVXr8DhYNIk3LuHoiIBMARYxeGsX79+yJo1wq+KXgaBABMm4P37JYAHsIfDcRw58g9qTdYn8vIwbRr8/UWoxFJ1adkSDRsiIIBpHeVQjxzh+/cYNQqTJv1MzCGCyMpixIh6t/VYDCJl8vNha4suXeDjQ6tdTU2cOYORI9G5Mw4dAoBmzaChIQUsA9oA7o6Ol6OiYGGBmBiqNJw9C1VVHD2KQ4esT5xYZGd3dvv2cSEhm6myV//w9UXnzhg4kGkdtcPKCgcOvP706RPTQn4HDcuV9PNrsAyfT5o3J8bGTCmqBrGxRFtbRMvWUxQs06wZefWqkmNEOViGxyPDh5ORIwmfz5iG2Fiip0fGjSPZ2WTOHN+GDZ243O1KSqaami/t7Ym/P2ncmLi5kYyMKvVWxWiI9HTSsyfhcomdXbmZAcQeqoNlHj0iamrk3TtKjVSPGlxyUVGRoqIp0IfL1R43bjoVqn4LGyxTCjs7pKfj6lWmdVQBCws0bIjLl5nWQRfv36OwEK1aMa2jphCCqVORlYW9e8Fl7vtkYYG4OHA4sLDAmDGLvb3b2tv/Gxm5JSlJ28QEK1bAxgZFRTA0xL59wrH4v/9BUxPPnyM2FseOQVpaON2ylEQgwOTJWL4c5WdWrhv4+/tnZ5sAFwWChIMHTzEtpyz1whHu3InTp3HmDJSUmJZSNSZMqEchMzExdXuBcNEixMcjPJx5T6CoiL17sWgRunVb4ev75OhRq6FDp6amvlqwAImJaNUKp06ha1ds3IjevfHsWc0NJSZCXx/z5mHePLx7BzMz4V0DS2lCQiAQiFC2oBpTUFAAKAAAJAGOQMRyj4q/I3zyBFOmYOlSCsPZhc64cThzBunpTOughTq9QLh5M/7+G//8AwUFpqV8Z/RoNGp0Kj9/j0AwOTXVc+fOCACNGmH5cjx7Bh0dvH6NggJYWWHpUlS3fFDxAKVtW8jJ4e1bLF9OySWwFPPpExYuxNatTM40CAt3d3dZ2dMcjgPQ29rakitilyRaaoQOj4cePdCtGxYtYlpKdWjUCAMG4OBBpnXQQt0NGT10COvXIzISampMSylN48YqwBOASErGBgVpGRvDxwfXrqFRI/j5ITER3buDz8fhw2jXDhcvVrXbqCioq+PAAezZg/v3wVSp1frD7NkYPx7GxkzrEAYKCgoZGc8PHBjO5e5dtEjkNnGLuSPs2RMCgUhv5CyPiROxYwfTIqiHEMTF1ckRYXQ0Zs1CRARatmRayi+Ehvq3bz+rSRMLZ+cGnz8P270bjRrB2xuamnBwwPHj8PREfDz++AOfPmH4cNjZoeJQvtxcDByIvn3RuTPS00U7e7K4cPkyrl7F4sVM6xAekpKSo0c7tmnTauNGpqX8gjg7Ql9f3L6NGzeYX7ypAdbWyM1FXBzTOijm2TOoqIjciKpSYmMxciTCw/G9dIpo0aZNmwcPLqakxG3fvk5CgmNuDi8vXLuGJ08weDAuXoSBAezsoKyM0FCMGIELF6Cjg7Vrf180LigIKiq4exc3b+LMmTr5bapzFBRgyhQEBIjQlLuwGD1aFCMBxdMRpqdnrlhxYNkyQUAADAyYVlMjOBw4O4t/yExdnBdNSsKwYdi2jeFizjVAQwNOTggLw/v38PNDfj7mzEF0NAYNgoEBFi6Ejg7i4nDw4MFevWxmz5796hUMDeHhgZkz8eEDOneu3ER09KVZsxaHh5+g/mrEmZUr0b49bG2Z1kEBf/2FrCzcv8+0jjLQsJODfgwM1gPuEhKiUXK3prx/T1RVSU4O0zpKIPR9hDNmkI0bq3SkiOwjfPeOaGszXwdHiLx4QTZtIjY2REGBNGxIOJx7gBGwFRjA4fibmJCUlKp2FRFxVlm5H3BGSWnshg1bqFQtQgh9H+GzZ0Rdnbx9K9xehUktL7llS/oScFdxHyGHEMK0LxY+srKZ+flKgG5BQYJ0XZ7KGTIEI0bAyYlpHd/Jy8uTkpKSFF5WUEtLrF9fpaGVpqbmgwcPNBmN0EhPR48eGDcO8+YxqIIqvnxBZCScnbsWFvoAtkAeYGRj84IQfP1aybmFhcjJwadPM7OzxwCWQKqZ2fi4uHqRPz4rK0ux6gWXK4MQ9O4NOzu4uwurS+FTy0v+6y8cOFDJsrSwEAgEfD5fSkqq4sPqfp7j39G06Z2XL9sCeQ8eSFtYMK2mFkyciA0bRMgRCpfCQjx6VGd2oeXlwdYWf/whnl4QgIoKRo3Cxo2Sd+7cAmyBOAkJvpcXOBw0alTJuVJSUFBASEjbjRsj8vM7cjgn8vIM8/PZYoTVZt8+ZGdj+nSmdVCJlxc2bcKbN9DSYlrKd8TTEcrJhQBXBg92srVF37743/+gqsq0phoxaBCmTcPTp3V1pbNi7t+Hvj7k5JjWUQX4fIwdi1atsHYt01Io5tKlMxoaJnl5+7hchIcH2thU41xfX7f09IXnznU1NTWTlFxrbo5Dh8Qk+p8evnyBjw8iIiAhwbQUKtHUhIYG1q4VoRzc4hksY2jIO3RofUSEX0IClJVhYoJ9+1AX54AlJTFunNjm4Bb9SJm0tLRx49w7dbLt0+dIZiZ27RKHrc0Vo6CgkJv7IjPzCZ//Zmg189NLSkpu2eL38mVMeHjwkSOKXl7o0wf+/nXyq8cIc+bA0RHm5kzroJ6BA3FClAKqxPxr3agR/P1x6hQCA2FtjYQEpgVVn0mTsHcvioqY1kEBop9Txt5+amho99jY4OvX982ff6suLzczgJMTbt/G4cOwsxPdQnSiw9WriIzE0qVM66AFHx+8e1f52jNtiLkjLMbMDLduYfRodO8Ob2/k5zMtqDro6qJNG/zzD9M6KED0R4TPnycLBPZAcz7fPj7+DtNy6h7a2rhyBa1bw9RUFHePiQ6FhZg8GYGBaNiQaSm0oKeHhg0hOjvrf+8I09LSwsPDlyxZMn36dA8Pj5UrV168eDE3N5dmcUKEy4WbGx4+REoKjIxw7hzTgqqDWJatz8jA27do25ZpHRViYGDE4QQAsSoqu3v3rjvJakUJaWmsX4/duzF2LDw8xHNuo/asXYtWrWBnx7QOGundG0dEJtVa2WCZCxcuBAYG/vPPP3w+H4C8vDyPxysoKCh+PWrUKA8PDyMjIwaUCoMmTbBvH6KjMW0a9PSwZYsIhS1VwIgR8PQUrSCr2hMbC3NzCG8jhvARCFBQENinz2YlpZ3Tpy81ZqM+akGfPrh3Dy4u6NYNoaFo3ZppQaJEUhI2bUJsLNM66MXbG5aWEJHQ4p8jwuTk5H79+g0YMCA/Pz8oKCg+Pj4/Pz87Ozs/Pz8zM/PatWuLFy+Oi4vr0KHDpEmTcnJyGBRdS3r3xr17MDdHhw5YswZ8PtOCKkNWFo6OQisjJyKI/ryonx+4XPnISJ/w8G29e1szLafOo66OiAiMGYOuXXHgANNqRAkPD/j4QFubaR300rEjGjTA1q1M6wBQ0hHGxMTo6uq+evUqMjJy8uTJxsbGMjIyxf+kqKjYrVu3efPm3b17NyYm5sOHD+/evWNIsHCQlcXSpYiJQXQ0LCxw6xbTgiqjeHZUxGp41QoRj5S5dw/+/ti/X/zDROmEw4GHBy5dwtq1cHJCdjbTgpjm8+fPBw8K3rwR6e3z1NG1q6iExP/8ljs6OgYHB2tVNvtmbm4eERGhr69PsTA6aN0akZGYPx92dnBywr17yRMnznZ29nzx4gXT0spiZoZGjXDpEtM6hIcojwjz8+HkBH9/tGjBtBRxpF07xMRAWRnGxrh5k2k1DJGRkdG+fa+2bUePH2/h6Xm3sswn4om7Ox49Eonf9+zPXdjbIyEBDRvyOnYcsXt33337BvTq5VBQ3Yql1CNOITPv3oHHE8UCRsXMnQsTE4wcybQO8UVWFv7+WLMGw4Zh/fr6uNEwMHDn06djUlMj+fyw4GAxKrZUHWxtISGB0FCmdVTFEV6/fn3s2LF9v0ODJvpp1Ah//fWmYUMdQvoT0q+w0CgpKYlpUWUZOxZnzyItjWkdwkCUh4PnziEiAoGBTOuoB9jb4/ZtnDiB/v3x4QP4fH5KSgpf9BfthUFOTj6fX7xVomFeXp3a0SVUzMxEYpmwkqA9Ho83adIkPz+/Zs2a0SOIKZo1ayYj8xy4D0gUFDzUFr2V64YNMXgwQkPFYTlBZBcI09Lg6or9+yvPrskiFFq2xJUrWL8eHTq8FAgcuNymMjIply8fFcEvoHBxdR2/fv2QBg2uN2gQs3JlPR0RApg0CdOmMS2i0hFhXl6ehYXFkCFDzL9Djyz6kZaWPndur43N2s6dV3K5Ia9fi2IGzIkTsXMn0yKEgciOCCdPxpgxsGZDRGlEQgJeXmjXbn1q6vqPH0/999+6+fPXMy2Kcnbtaj5o0NUzZ+wfPDgxZMhApuUwhrMzeDzmN3ZXMiJUVFRUVlZ+/fp1S5FdzxEexsbGFy6EAti5Ew4OiIkRuXzQPXsiP190h1NVhBDcvSuKl7BrF16+xKFDTOuolygoFALF3ze5z58LGVZDMfHx2LED8fHyjRvXteLOwobLhaEhNm9G//6Myqjg33R1dVVUVHbv3t2qVSuV79CmjEFcXWFmJoozkBwOXFzqfMjM06dQVRW5eiCvXsHbG/v2gU0oygjLls1s3Hiymtr0hg2nxsa6Hz7MtCDK4PPh6go/PzRuzLQU0cDJCVevMqyhIkf44sWLL1++FBcj/vId2pQxy5YtuHkTe/cyreMXXFwQFla3N2CJ4LyoQABnZyxciPbtmZZSXzExMUlM/Pf0aac3b65ERbVftAhOTqjLWR3L5X//g6IinJ2Z1iEyzJiB3FzExDCpoRrbJ/h8fl5eHnVSRAp5eYSFYe5cPHnCtJTSNG6M7t3x999M66gFIji1u3IlZGQwcybTOuo3ioqKnTt3VlRUNDPD3bsoKkKnTnj8mGlZQiU5GWvXYscOcDhMSxEZGjSAtjbDCbhLOUJLS8ugoKDi14SQMWPG3CqRc+XIkSNyorZoRiXt2mHNGjg4iNzP0rq+oVDURoRxcQgMREgI+2wSIRQVcegQ5s1D794isc9MKBACNzd4eUFXl2kpIoa9PS5eZFJAKUf4/v37zMzM4teEkNDQ0OTkZAZEiQwuLujYEZMmMa2jNAMH4uXLOllbEUBBAZ48gakp0zq+k5uLMWOwebNYJTQXG5ycEB2NlSvFZJp0926kpsLDg2kdosfcuUhPR2IiYwLYzDKVsGULHj8WlYR4xUhKYvx47N7NtI4acf8+2rSBrCzTOr4zezYsLeHoyLQOlnJo1w63boHPh4UFHj1iWk0t+PABPj7YtUukK64whYoKmjTBhg2MCWAdYSXIyiIsDPPmIT6eaSklKN70XRdLu4nUvGhkJM6fx+bNTOtgqRBFRRw8CG9v9OpVh/fRzpiBSZNEaC5E1LC1RUQEY9ZZR1g5+vrw94eDA7KymJbyHW1tIiu7qV27od7eKwoL69KmK9GJlElNxYQJ2LULSkpMS2GpAsVB9v7+cHJCnasCd/o0Hj7EwoVM6xBhfHzw4QM+fWLGOusIq8SoUbCygpsb0zq+Exy86/3758+fB2/ezF+0aC3TcqqB6IwIp06FszN69mRaB0uVadsWt29DRgYWFnj4kGk1VSYzE9OmYedOkahAK7K0bAkVFcZiR8s6wgULFkhISEhISEhJSQEYPXq0xHfGjh3LhEJRISgIT59i1y6mdQAAoqNj8/MnAE3z8twuX77NtJyqkpGBlBS0bcu0DmD7drx6hSVLmNbBUk1kZbFjB3x80Ls3/P2ZVlM15syBrS26d2dah8jTty+OHmXGdKl1WwcHhzTxqG5AAQ0aICwMVlYwM0OHDgyLGTy4x4ULmzIz50lIhAwaVGcGNbdvw9wcEhIMy3j5EgsXIjqaTSJTV3FyQqdOcHBAXBy2bIGCAtOCyufKFZw5U7fDfGjDxwcdOiA7m4E3tJQjXLduHd326xR6eti8GQ4OuHOH4YUlZ+fRfD4/PHxjbGznrl1FbHtH+YjCvCiPhzFjsGQJjIwYVsJSGwwMEBMDb29YWCAsDE2apCYmJhobGysqKjIt7ScFBZgyBUFBbDGTKmFsDHl5BAXBy4tu0+waYfVwdETPniKxWDhx4rgzZ3YHBEyZN09CFEo8VwVRiJRZsQJKSiJR+YWllhRX912wAD16XNHRGThkSFibNj1EqpLo4sXo0AFDhzKto+7Qowf272fA7k9HmJSUlJqaWpVzEhIS0tPTKZMk6gQEIDERO3YwrQMA4Oj4bc5WxLlzJ651624REWZHj/5FmKtHHhODbduwZw+bREZ8GDcO2tqbsrOPfPni//6939q125lW9I27d7F3b51ZyBQRPD2RkAAej267Px1hfHy8jo7O3LlzH5eT3Y8QcuXKldGjR3fo0CEjI4MuhSJHseNZuBD37jEtBeBw4OcHHx8UFDAtpULGjvV88SJUILgbGZkbGRnJiIbsbIwbhy1b0KQJI/ZZqEJVtQFQ/ETKuHOnwcePDOsBwOPB1RVr10JDg2kpdQobG0hJMZDA5Oca4fDhw+Xl5efNm7d+/XoDAwNLS0t9fX0VFRUej/fly5f4+PibN2+mpKT0798/Li6uVatWdCsVJVq3RkAAHB2ZXywE0KMHjI0RGIjZsxlWUgFZWblACwA5Oe3fvk1hRIOnJ3r0gJ0dI8ZZKGTz5kV//DGuqEhDUTHTyupEu3YYMwbz50NTkzFJGzZAWRnjxjEmoO7SsSN27oSrK71WSWkEAkF0dPT48eNbtGjx4xgOh2NkZDRr1qxHjx4RYXDw4MEWLVooKCgMGTLk8+fPvx6QkJDQuXNnOTk5IyOjGzduFDeGhYXplODBgwfl9e/g4HDo0CGhSK2AyZOJgwPVRqrE06dEXZ2kpVFuKDc3t6ioqAYnengskpAYIym5oVkz0w8fPtSgBw0NjZqdeP36jZ497c3NnbS0nmZm1qCD+khmXbtTfD7/06dPxa/fvCHu7kRVlXh5kfR0yk3/eq8SE4mqKnnxgnLTTEHpx+PgQSIlJbTe+Hx+YWFhpYeVdYQl+fr167Nnz5KTk3NycoSmi5Dk5GQFBYXr16/n5eWNHDly4sSJvx5jbm7u6+vL4/F2797dpEmT4isJCQnp27fvi+/k5+eXZ4IeR5ifT8zMSHAw1XaqxNSp5K+/KLdSY0d46RJp2fLy3r37UlNTa2a6Zo7w8+fPGhrmwBPgloaGuUAgqJn1+kadc4S/kpxM3NyIhgZZsoRkZFBoqMy9EghInz5k0yYKLTIO1R8PSUly4oRwuhKCI6SIZcuWDRs2rPj1w4cP5eTkcnNzSx5w//59eXn5H406OjonT54khISEhPw4sWLocYSEkOfPiYYGiYujwVQlfPpE1NRIUhK1VmrsCEeMIEFBtTJdM0d469YtVdXpAAGIpuaAH4MGlooRA0dYTEICcXAgTZqQTZtI+b+ca0WZe7VtG+ncmfB4lNgSEaj+eJiakt69hdNVFR0hA9snEhMT238vBG5oaFhQUPD27duSBzx//lxXV1f2e4WC9u3bJ36vz3HlypUWLVpYWFgEBgYS5oIPf9C6NQID4eiId++yUlKYWfoqRl0dHh4imswwJQXR0WAkMZGhoWF+/k3gEodzUlHxq7q6OgMiWJjDwABHjuDCBVy/Dj09+PtTG1b2/j0WLcLOncxnjajTuLjg5k1aLTJQESQ9Pf2HI+RyufLy8l++fCl5wJcvXxRKpBZQUlIqPsDS0jIyMlJLS+v+/fsTJkyQkJCYOnXqb008ePAgLCxs1KhRxX+qqaklJiZKUlP+pH9/BAae0NX1V1JSa9NGJiJirwRDXwI3N5ibK0RF5XXqxKfIRF5enpSUVHXv5ObNMvb24HAKapO1nBCSnZ1d3dLQ+/ZJqavv79s3QEVF1strT5bo5E0XbbKzs5mWIExatEBICOLiJNaskd6wgTtnTuG4cUXCeh6UvFeTJ8tOnCho2bJWH3XRh+qPx9ixmDVLMTIyt2vX2j7KBAKBlJRUccbQihDO+LM6ODk5LVy4sPg1j8fjcrkvSi8rh4eHt2/f/seftra2GzduLNPJpk2bevXqk2WP4wAAIABJREFUVZ4J2qZGi2nSxBTIA0jDhn+dPXuWNru/EhJCevSgsP8aTI0WFpKmTUnto6xqMDV65Qpp3Jg8f15b0/UQsZka/ZXr10nv3sTAgOzdK5wJzB/36uhR0qYNycsTQp8iDg0fD319UrV1sEoQ3anRNm3axH8v7vfo0SM5ObmmTZuWOSApKSn3e0Xqhw8ftmnTpkwnEhISRASmRoshBIAEgKIimSJGiwQ6OyMrC6dOMSihLEePwtAQ7drRbffVK4wciYMH0bo13aZZRJmuXREVhYAABAbCxASBgS+MjHo3aWIxZIhzbb68GRnw9GRLTAiNUaNw6RKN9oTgc6vJmzdvFBQUIiMjMzIy7OzspkyZUty+atWq0NDQ4teWlpbe3t45OTmBgYFaWlrFo5CjR48mJCRkZGRcunRJS0trU/mBWTSPCNet26Ku3r1hw7FSUjavXxfQZve3nD1L2rQhVfgNVBNqMCLs2pUcPy4E09UaEX79SgwMyNatQrBbPxHjEeEPBAJy4gSRl7cDYgHSoMGyLVt21KCf4ns1YQKZOVPYEkUVGj4eGRmEwxHCTJLojgibN29+4MABT0/P1q1by8rKrlmzprg9NTX1R8KagwcP3rlzR0tL68CBAydPnixelHr48OHgwYO1tbXd3d09PDxmzpxJv/jfMmfO1Pv3D1+9Os/LK3L8eGn68wOVpH9/tGghKoW84+Px5g0GD6bVKJ+PkSPRty8mT6bVLkvdgsPB0KFo2vQT0A5Afr7xnDkpNjbw8MCOHbh1qxqFuKOiEBWFlSspVFvfUFJC8+agrQwEh4jMBKMQcXR0tLOzGzlyJM12BQIMHgwjI6xltFZufDz698ezZ8LPelPdYJlJk6CtjfnzhWBaU1PzwYMHmlVIFuLpicePceYMqImOqhdkZWWJVBkH6vD3375s2anMzH4qKnsiIkILCw3i4vDkCR4/Rnw8FBXRrt23uX1DQ5ib43swOwAIBILQ0CNxcc/Cwx2CggxtbZm7DHqh5+Ph7o4jR1DLhHkCgYDP51caLMM+KoQJl4uDB9GxIzp2hL09YzJMTNCvH9atw/LljGkA8PUrwsORkECr0T17cOYMbt1ivSBLlfDwcOva1fTJkyd9+pxq3rw5ACurn/+akoJiv3jtGrZvx7Nn0ND46RdPn/Y+e7YoO9uqQYPxbdocBPQZuwxxZMECBAbi3Ts0a0a5LfZpIWSUlREejn79mIkQ+cGqVTAxgZsbtLQY0xASgsGDac33+O+/8PLCv/9CWZk+oyx1nY4dO3YspzxY06Zo2hQ/hnqFhXj6FI8f48ED/P03zp69xOPFAigoyD116tycOawjFCaamlBTw7p12LSJcls/1wifPXt29OhRyg3WA0xMsGED/vwTmZmMaWjWDG5uWLqUMQGEYPt2TJ9On8VXr+DoiAMH8EuIMQuLcJCWhrExRo3C6tU4dQrm5o2BGwBPSemikRHrBYXPgAE4fpwOQz8d4Y0bN5Z+f3A2btz4xo0bdNgXU8aORa9ecHICgyuwPj44e5axWlFnz0JeHp0702QuOxtDh2L+fPTtS5NFFpawsCArq3VaWl1mzGjbv39/puWIIV5eePMGX79SbuinI1RWVv7y5YugrhQ7F3kCA5GWhg0bGBOgqIj584UTqFIDgoLg7k6TLYEAY8agY0eITBwxS72gRYsWV68ef/w4esUKb6a1iCeGhlBSoqO48c81QgsLi69fv/br109XVzczM3Pt2rW/jdDbtm0b5aLEAikpHDqETp1gYsLYMGXyZAQE4MIFugW8fo3bt0HbRLuPD9LT6TPHwsJCG7164dAhLFlCrZWfjrB58+ZHjx5du3ZtREREQUFBdHT0b6PkWUdYdbS0cPAgxo1DTAyaN2dAgJQUVq/G3Lm4exdcGreMBgZiwgRUMy1oDdm/H0ePIiYG0tJ0mGNhYaGTOXPQsycKC6n9gpd6Og4ePPjff/9NSUlRV1c/d+7cl99BoRZxpHdvzJyJESOozXlfAX/+CUVFHDhAn8W8POzdiylT6LB18yZmz8apU2CrSrCwiCVWVpCRwY4d1Fr5/TBhx44dv6b3LElSUtLWrVupkSRueHlBSwt//cWYgPXrsWABvqdupZzQUFhaQlubckMpKXB0xK5dMDKi3BYLCwtTWFpi1y5qTfzeEdra2qqqqlZw2rNnzzbRsLlDLOBwsGsXLl3C7t3MCOjcGZaWdCw4FxMcTMeuibw8DB0KDw+687exsLDQzMCB8ffuzTUx6RMbG0uRCQZyjdZDFBVx/Di8vXH3LjMC1qzBhg21TVZUFW7cQEYG5bE5hMDFBe3aYfZsag2xsLAwS35+vo/PIEK6PHgwvWtXu1xqprZYR0gTbdpg82YMH47PnxmwrqODMWOwYgXlhoKCMH065YE5ixfj3Tts306tFRYWFsa5efOmQNAe+BP4UyAwv3z5MhVWWEdIH46OsLPDqFHgU1VAviKWLsXRo0hMpNBEairOnoWTE4UmAPz9Nw4cQHg4GybKwiL+mJqaAo+BROA5IU86depEhRXWEdLK2rUoLKRjZPYrysrw9ISPD4Umtm2DvT1UVCg0cfcupk/HyZPQ0KDQCgsLi4jQqFGjnTuXKyoOadDAAdiTmKhGhRXWEdKKpCSOHEFICM6cYcC6hwfi4nDtGiWd83jYto2qEoCZmZkFBQWPHn2ys8P27TA2psQKCwuLCOLiMj4z82le3r2+fbsNHIjCQuGbYB0h3WhqIiwMEybg5Uu6TTdogBUrMGcOJRlQT56EtjbMzITfc2xsrIpK24yMfBsbR3PzE0OHCt8ECwuL6PPPP5CURL9+wu+5ho6Qy+VWvTorSxksLTF/PkaMQF4e3aZHjyYvXnipqpoZGfVOEGqpwOIwGSpwc/Ph8wOBRsDuixepnNtlYWERYSQlcfkyrl3Dxo1C7rkajjA3NzcrK6v49YABAx49eiRkLfUJd3e0bYtp0+i2e+bMP7m5GenpcY8fB40cKbSs2AkJePoUdnbC6q8UPJ4EUAQA4HM4HEpssLCw1AWMjLB6NebNE3Lc3+8d4ciRI3fu3Fmm8e+//9bV1SUMFhYSL3buxP37hXZ2/nZ2bseOnaTHaHLy2/x8C4ADGHz+LLTqJkFBmDSJkjDOU6eQkrKOy/0LSOdyRwQGMlRNg4WFRTSYOxcdO8LKCkIslfQbR1hUVHT8+PFmzZoB+O+//86ePVvc3rdv39TU1OTkZKEZr9/IyqJdu6UnTnw5cWLKxIk7IyMv0GB0yJCB6upBXO52KamJ8vIDeTwh9JmVhUOH4OYmhK5KwuPB2xvu7jh92jgnJ6lRI7mEhMvjxo0VshkWFpa6xqVLyM/Hn38KrcPfOMK0tLTCwkItLS0AV69e9fDwKG7X0NDgcDhpaWlCM17viYm5AiwAzL5+nRkRcZkGiy3+z96dx8W4/XEA/8y0TCvRHkV1RUi6yZJISYgoW5biWiqyXYSb7WYPl8qesiS7hJItlFxZksi+hKLVnrQ3z++PuZZfppppntnqvF/3j+mZ85zznbmZb/M853yPgcGtWzGbNlHHj49s1WqpuzsNixrDw+HggGbN6Ijvm4wMdO+OJ09w5w66dYOCgoK8vHzjxo3pHIMgCOmkoIC4OMTE0FaDlEsilJeXB8C5Hfj+/fvvO07k5+dTFFVzDVKCL127WsrJhQI5DEZEkyZWohm0efPmPj7eAwY4Hj2Kd+8waZKgVxhCQmieJhMTg86dMWQIoqKgpkZnzwRB1A9dumD2bEyejKwsGnrjkgjV1dX19PT++eefe/fuhYWFycrKHjlyBEBQUJC6unpzsWysV09t27Zq0qQ3nTp5zpzZY9cuF1FulgRAURHR0cjMFCgXXrwIikKPHvSEVFEBf39Mn47jxzF/PsjMGIIgqrNuHVq1QrduNHTFfbLMihUrTpw40aFDBx0dnRUrVri5uTVt2nTt2rV//fWXPClsRR8VFZWtW1cnJ58KDPSKi8PixfD3F2kASkqIicGrV/D0rGMu3LIF06fTk7Fev0avXrh9G7dvw9qahg4Jgqjfrl7Fu3eYOFHQfrivBRw/fny3bt0yMzPt7Ozk5OTU1NRSU1Otra0HDBgg6IBENdq2RVISBg7EmzfYvh0iW6WppIRTp+DkBG9v7NjBX0p7/RoJCQgPpyGMU6cwaRImT8aSJUKv2U0QRP2gpobDh+HigiFDIEh2YtTL5RBubm6urq4jR44UdyB8KyzEiBGQk8PBg1BSEt24X7/CyQmtWyMkpKZcWFxcLCcn972WwsKF+PoVAm5MWVGBFSuwezcOHED37tU209bWTktL09bWFmgwgjdfvnxRVVUVdxTSoQG+V5L2kseNw5EjyMriUuiYzWZXVlbKycnV3AP521uyqKggOhra2rCzQ36+6MZVVkZMDNLSMGsWr6eUlWHXLkyZItC4b97Azg5Xr+LmzZqyIEEQRHXCw6GjA1vbuvdAEqHEkZVFSAj694e1tXB3TaqiUSOcO4dr13jNhUeOoEMHtG5d9xEvXkTXrujdG+fOgXzTIwiizq5fx9OnmFHXYlkkEUoiBgP+/liwAD17IilJdOM2bozz53H1Kk87vwtSXLSyEv7+mDABhw7B35/cFCQIQiDa2ggLw5Ytddxdh3wCSa4JE7B3L4YOxalTohuUkwsTE+HrW1Oz1FRkZ9fx7nReHvr1w7//4uZN2NjULUyCIIj/4+EBJyf064eiIr7PJYlQojk6IjoaXl7Yvl10g6qpIS4OCQmYN6/aNps2wccHMjK89vn48WMLC0cDA6vRo1dZWqJ7d5w/Ty6HEgRBp5MnoaKCPn34PrGWRPj169ewsLD379/XMS5CYFZWSEzE+vVYskQo+whypaaGs2dx9iyWLuXy7MePOHEC48fz0eHw4VPv3Al8/fr64cMP586NJ5dDCYKgHZOJhATcuIE1a/g8sean379/7+npmZmZWffQCIH99huSknD+PCZMQHm5iAbV0MDFi4iMxPLlVZ8KC8OgQdDSqr2TigokJWHZMjx9+hloB8gwmTZKSiLfkpggiIahTRts2IAFC5CWxsdZ5M9y6aCpiYQEFBaif38UFIhu0IsXcfgwVqz4cZDNxvbttUyTefECO3ZgxAhoasLDA1lZsLbuoaw8k8EIV1ff0a+fEHaYJgiCAADMmAFra/TuDd531yG7zEsNBQUcOoQZM2Bjg9OnIZqar1pauHQJdnZgMrFgAQCcOcNQU4PVLxXC8/KQmIgLF3D2LMrLYWMDBwcEBv63KwWbvf7YsaiXL9+4uR3nbGxCEAQhJBcvQkcHgwcjJoan9rUkQllZWV1d3VqX5ROiISODLVsQHAxra2r06NCUlMu9enWaP3+6rDALsn3PhW/fvrx/f/aDB9aLF3sBjQF8/Ypr13DhAi5cQEYGunaFjQ2iomBpWbUTJpM5fPgw4QVJEATxnbw8zp1D167lPXr4rVs3xtraoub2tXyA6unpZWdn0xceQYOZM5GWFrF27Q2KWpKUtKu0dMOyZdXP76SDtjYOHnzbseMQYBrwbs6cTgUFzy5cQHIyrKzg4ICQEFhYkPkvBEFICisrKCraJyXZnT+vVWsRf/LRJZU+f06iKB+gdVHRn/v3//v5s9BHPH/+KIPRGZgIzC8pUcnOzvHzQ24u4uIwfz4sLUkWJAhCshQXZwHLwsP1am1JPr2kkqNjN1XV7UA6i7VRTs66RQsMGoSICNCeEfPyEBICR0csW2ZJUUnAVyCLwXgbGKhtbw8FBZqHIwiCoIuMTCXwdMqU2neYI4lQKnl6jl269Hdb2wULFqjfv++bkwNPT8TFwcgIzs7Yu1fQjPj6NXbsgLMzTExw6hTc3ZGd3WXSJAdZWVMWy2bLFn8m+QJIEIRkO3AgUFGxf6NGF2ttSbZhqleKi3HhAo4eRWwsrK0xfDgGD0bjxryenpGBEydw9CgeP4aTE4YPh6MjWKyf+/+/bZhEiWzDJEqSts+OJGuA75UUvWQet2EiyyfqFUVFODvD2flHRpw1q2pGLC4uTk5O1tfXNzQ05Jz14gViYnD0KJ48Qf/+mD8f/fqBzBQmCKKB4CkRvnnzJi0tzdzcvBlnURgh8arLiE5OH9as6fflSzcmM83Ly11ObuLRo/jy5b/8178/xPFljyAIQpy4f+y5ubmZmpr6+/sDSEhIcHJyKi4uZrFYhw4dcnFxEWmAhGC+Z8QvXxATg3XrjmdkeADTgdK1a3vOmTNx925YWdW0Kz1BEET9xmXKQ0VFxcmTJ62/rbzw8/Nr1arVlStXxowZM3PmzMrKStFGSNBDVRWjR2PePGVFRU4J9S8tW8oGBKBzZ5IFCYJo0Lgkwg8fPpSWlhobGwPIz8+/efPm/PnzbWxsVqxYkZmZSQpwS7Vhw4Z26nRHS8tWR8dh+/aV4g6HIAhC/LhcGuVMsCktLQVw5swZiqJ69+4NoGnTpgDev3//fZIFIXXk5OQSE098+fJFWVmZLIEgCIIA12+ETZo0adasWVhY2OfPn8PCwiwsLDhz1jMyMgBoamoKPmpSUpKnp6enp+fVq1e5NsjNzf3rr7/GjBkTFhbGZv9YDhkdHT1u3Ljp06c/evRI8DAaLFVVVZIFCYIgOLh/Gq5YsSI4OFhNTS0pKWnhwoWcg6dOndLU1DQwMBBwyNu3b/fr18/CwsLS0tLJySklJaVKg/Lycltb248fP7q4uAQHB69c+d8VvMjISE9Pzz59+mhpadnY2OTm5goYCUEQBEFwnzX6xx9/WFpapqamduzYsUOHDpyDenp6GzduZAg8syI4ONjb29vHxwdARkZGUFBQRETEzw2io6MZDMb27dsZDIa+vr6zs/O8efNYLNY///yzcuVKd3d3AKmpqWFhYYsWLRIwGIIgCKKBq3bVmJmZmZmZ2c9H6CrUcv369Q0bNnAe29raTp8+/dcGPXv25GTczp07FxYWPn/+vE2bNrdu3dq3b9/3E+Pi4miJhyAIgmjIqk2Eb968CQsLe/DgQXFx8alTpwCcOnVKVVXV1tZWwCFzc3PV1dU5jzU1NXNycn5t0PzbtrNMJlNdXT0nJ0ddXb2ysrLmE787q3E26n7U2MVjOT8yGAwlJSUBwybE6/OEz613tRb8ggTBC4qiyFvNowb4XknXSzbXNk8Yn1BzG+6JMDk52dHRkXNl8v17zrIz3L9/Pzw8XPBZKgoKCmVlZZzHJSUlv6YoRUXF8vLy7z9y2igqKgKo+cTvDAsNZVVlv2dTeXl5MtOVFhUVFUwmUywTbTZv3uw+3l1ZWVn0QzdAZWVl8vLy4o5COjTA90qKXjJFUa2atqq1GfdEOHnyZAsLi+PHj6empnLuyQEYOHCgn59ffn6+lpaWIJHp6+t/X4yYmZn5a9m2Zs2a3b9/n/P4y5cvnz59at68eePGjVVUVDIzMzlTWLme+F3rktauHRti0W1hE2PR7T3uexbvWEyKbouGFFVVFrsG+F5J0UvmFN2utRmXP+0/fvx4+/bt5cuXN27c+Ofvvy1atACQlZUlYGRDhgyJiIigKIqiqH379g0dOpRzPDIyMi8vD8DQoUPj4uI4k0IPHDjw+++/c6aqDh06dO/evQBKSkqOHj36/USCIAiCqDMuf9pzltL/mvA/fvwIQPBvAz4+PpGRkZ07d2YwGOXl5VOnTuUc/+OPP44fP96nT5/27duPHz/eysqqffv2KSkpkZGRnAaLFy+2s7O7f/9+bm6usbGxq6urgJEQBEEQBJespqWlpampefr06Q4dOvz8jfDQoUPKysomJiYCDqmmppacnJycnAzAyspKRkaGc/zRo0ffV+sHBgZOmTIlKyvLwsJCTU2Nc9DY2Pjp06fJycmqqqrm5uZSdLeWIAiCkFhcEiGTyZw1a5a/v39lZWWzZs3YbPbDhw8PHz4cEBAwa9Ys1s/7tNaVjIxM165dqxzU19f/+UcTE5Nfk66CgkKPHj0ED4AgCIIgOLhf55w/f/67d+/8/f0rKioAtGvXjsFgjBs3btmyZaINjyAIgiCEi3siZDKZ69evnzlzZnx8fG5urpqaWs+ePU1NTUUcHEEQBEEIW00zXwwMDMaNGyeyUAiCIAhC9LgnwoyMjOrWXhgZGQkzHoIgCIIQKe6JsEuXLpwlfb+iKEqY8RAEQRCESHFPhKGhoSUlJd9//PLly+XLl0+ePLlq1SpRBSaQr1+/JiQkkMoyBFEHL1682LBhg4uLi4ODg7hjIQhRYPD+DS8gIODs2bMJCQnCjIce7dv7PXiQoKdXmpV1W9yx1CtiLLGmra2dlpZGSqwJW2xs7MCBk4CBQKKra6eoqP3ijkjSSVG9MbpI0UvmlFiTk5OruRkf1ZPHjBlz+fLlFy9eCBaYKBQUOANXs7Pffy/STRBEzYqKcP8+/vhjMxAEhAJJJ04kijsoghAFPv60//z5M4CioiKhBUMbNlsGKAMY0lIinSBoVFZWtnlz2P376ZMmDbO27lblWTYbWVl48QIvXuDlyx8PPn2CoSGKi82ADwCATxTV2t0dGzZAsDL7BCHpeJo1Wl5enp6e/vfff6upqQleYk0EGjfek5V1DlixejX8/MQdDUGI1oQJc6KiNIqLB0ZHz927d6u8fAdOtuP89/gx5OVhZPTffz174o8/YGSEli3BZMLX12f9+p7ATiB7wIClly9DRwdt2uCff+DkJO4XRhDCwcesUT09vYiICKn4jmVq+nbixFVLl45cuBCfPyMgQNwBEYQIXbp0vbg4GcD79z4TJ17s2rWDoeF/Oc/QEIaGUFDgfuK+fTh0qGWTJpn79/87b17nefPke/bEjRuYPRvOzlBXx4wZWLAA4tiPkiCEiKdZo7Kyss2bNzczM1Oo7h+QhJGRkdHTw7p1WLsW69bh0yds3y7umAhC+F69wpo1ePeuJYNxiqJ6qqmd2L9/sr09T+fGxmLuXDg6QkcHNjbmXl7y27ahZ0906YKrV/HxI+bMwcqVWL4cgwdj82ZyvZSoR6j6aMSIEQcPHmSzKUdHasIESkaGGjdO3DHVC0VFReXl5WIZWktLKzc3VyxDS4X0dGrGDEpDg5oxg7p3L8/FZaKpaa+goO08nn79OqWpSZ07RzVpQmVmUgUFBQUFlLo69etbHhFBGRlRDAZlakqdOkXzq5BGBQUF4g5B1KToJVdWVpaVldXarD5f42AwsGMHoqOxaxcOHADZvpColx48wNixsLZGkyZ4+hTBwWjfXuv48bCHD+NnzvTmpYeHDzF4MPbuxc2bcHUFZxsYVVUMGYJdu6o2dndHejpSUqClhUGDoK6OxYvBZtP9qghChH5cGj1z5kwADzfTLl++LMx4aNaiBf7+G9u3Iz4ednaws0N8vLhjIgiapKXhn38QFwdvbzx5gsaN69JJVhacnLBmDezsMGEC4uJ+PDV1KgYNwrx5+LZn6A8WFkhIQEEB5s7F+vVYswYODtiyhT1+/NDk5AetWhkkJER+30mUICTcj2+EcnJyKjwQY6x14+MDeXmkpOD2bVy/jm5VJ5MThPRJSoKzM/r2Rbt2SE+Hv38ds+Dnz3BywowZGDcOe/bAygrt2v141twceno4c6ba0xs1QkgIioqwaxeePIGxcdjly7pFRf/evdvb1nZYXQIiCHH48Y3QwcGhXlZUYjKxeze6dIGTEx4+hJkZzMyQmgpxVEchCEHFx2PFCrx6hfnzcewYBJnEXVQEJyf064fZs8FmY/16LhdCp0zBtm0YOLCWrtzd4e4OXd39ublrAS3A89mz3XWPjCBEqz7fI/zO0BB//QVPT7RsiUePkJEBMzNUVIg7LILgx7//ondvTJoENzc8eQIvL4GyYGUl3N1haPjf4qLjx9GkCWxsqjZzc0NKCl6+5KlPN7duDMZSIBlYSlFD8vPrHh5BiFJNX4vy8/NfvHhRWFj480Ep/db45584fhwhIZg8Gc+fo21bGBvj0SMoKYk7MoLgJiUlxc9vA4slv27dX8+etV6xAkVFmDsXY8ZwuWPHL4qClxfKynDkCBgMAFi/Hn/9xaUliwV3d+zYgdWra+82KCigstL3xInpXbqYpaUFtmiBs2dhaytotAQhdFznkubl5fXu3Zv39pKGs3yiysHHjykNDSo9naIo6vNnSleX0tamPn4UQ3jSiyyfEI2CggIdnd+BO8A1Obnfzc0rjxyh2Gza+p8/n+rcmSos/O/HhATKxISqrPy/AL4/fv6c0tKiiov5HsXdnWIyqYAAAYOVdFK0loAuUvSSBVo+4eXl9eDBg/379w8cOHDChAlnz56dMWOGmppaRESEMJOycLVujTlz4OUFikKjRnj+HHJyMDJCNRsvEoTYPH/+vLzcEjAHuqqo6MfG5gwf/t9XN8Ft3YqoKMTEQFn5vyPr1sHXt9p6McbG6NgRUVF8DxQRgeBgLFyIgQNBtjElJBmX3302m33+/Pn169ePHj1aQ0NDR0enb9++wcHB/v7+AQEBlDT/Rvv6oqAAu3cDgJIS0tOhqYnWrZGVJe7ICOInjRubFBSkAIkMxrnGjbP19PTo6vnQIQQEIC7uR12Yx4+RkgIPj5rO4kyZqYNp03DjBhIS0LIlyC1DQmJxSYRv374tLi62tLQEwGKxCgoKOMc9PDwePHjw/PlzkQZIK1lZhIfDzw9v3gCAvDwePYKhIUxM8PSpuIMjCADAy5fo10/Zy2vfiBGHPDzOJCQcZdD0ZfDSJcyciZgYtGjx4+CaNZgxo9rqoxzOznj9GmlpdRnU0hLZ2VBQQMuWSCTbOhESiUsiVFNTYzAYnE2XmjVr9vjxY87x4uJiSMk2TDUwNcXUqZg8+b8fmUykpsLSEh06IDVVrJERBHDvHmxt4euLzZtNDx/eGh4e1OLnrCWAW7cwciSOHYO5+Y+DWVk4eRJeXrWcKyODCRPqXrC3USM8eYKhQ2FnhzVr6tgSYqCDAAAgAElEQVQJQQgPl0TIYrHatWuXnJwMoH///vHx8WvXrr1w4cKkSZMaN27cqlUrkQdJswULkJODfft+HElMhL09rKxCZWRayMjojxjhKb7oiIYrIQG9e2PDhtozE7/S0zF4MLZvr7pAIigI48dDXb32Hry8cPgwvl0eqouICGzciIULyXZOhOThOoVm3759QUFBnMczZszgXJlRVlY+dOgQjfN5hIfrrNGf3blD6ehQOTk/jmRnZwNtgRKgjMls8/TpU6FHKYXIrFHhOX6c0tSk4uLo7zk7mzI0pMLCqh7//JnS0KAyMricwnVa4LBh1LZtggZz8yalrEwZGFBv3wralYSQoimUdJGilyzQrNExY8bMnDmT8zg4ODg/P//mzZtZWVlubm6iS9HCZG6OiRPh/VNF4szMTCZTC2ABcoDRSx6XEBMEHXbvxrRpOH8etC/TLSjAgAHw9MTEiVWf2rYNTk4wMOC1qylTsHmzoPFYWSE7GywWDAzILUNCUnBPhFn/P41SQ0PDysqqcd2qGUqqJUvw4gWOHPnvxy5duigpvQb+BPyYzDf2PO7hRhACW7MGK1YgIQEdO9Lcc1kZhg5Ft27w86v6VGkpNm7ErFl89Mb5N3H1qqBRNWqEp0//u2W4dq2gvRGE4LgnQktLy06dOu3YsePLly8iDkhk5OWxcydmzPgxq/vt2/uLFqkOGqTKZt/KySGlSAmhoyjMno39+3HlCn77jbZus7OzHR1H//ab9e+/B2poYNMmLm327YO5Od+p19OzjusofsW5ZbhgAbllSIgf90S4bNkyBoPh7e2tq6s7bty4+Ph4dn3ccKxzZ4wdi2/XgKGgoLB8+fKTJxeYmMiRzQsJYSsrw6hRSE3Fv/+CvoWCAODmNvXChYnp6fFPntxydz//60p5isKGDZg7l++ex43D6dO0rQicOhXXriExES1b4t07evokiDqotrJMcnLyo0ePfH19L1++bG9v36JFi7/++kuqFxFytXw50tKqVs04dQqpqXUppUEQPCosxKBBKC3FmTNo1IjmzjMy3lBUb4BVWel0796DXxtER0NJCXZ2fPespgZX1/9KUtDCygqvX0NWFgYGpUZG/Vgsww4d7KvUNyYIYatp94k2bdr4+/unp6efPn3axsYmODjYxMREZJGJBouFsDDMmIEPH34cNDbGsGGYMIHsu00Ixfv36NMHeno4erSWlex14+DQQ0bmLyBWXX2zs7Pjrw3WrcO8eXXsfNo0bNuGykqBIvxZkyZ4/hwqKtNevuxeVvbg3r0eAwa409Y7QfCg9m2YZGRk9PX1mzVr1qhRI0qa66tVp1s3DBtWddZARATKyjBnjphiIuqvly/RrRv69sWuXcLaFHPmzLWqqq2nT79x/vzWdj/vtAsAuHEDOTkYMqSOnVtYQEsL588LGmQVlZU3ATdACRhz//4rmnsniBrVlAjfvXu3ceNGS0tLMzOz0NBQZ2fnK1euiCwyUVq1CklJOHnyxxF5eQQGYtMmUpKboNODB7C1xZQp8PcX4ihbtsj6+o7fuHGZhYXFr8+uXg1fX4H2cqpz6dEa9O/fjcGYDyQCiwoKvC9epLl/gqgJ19WFMTExLi4u8vLyTCazd+/eERERX79+pXmhozDVuqD+V/HxVLNm1LNnb/Pz878fNDKirK3pDk6akQX1grh2jdLVpYRdlOLDB6ppU6q6t+rxY0pHhyoqqqWTmldMFxVRGhrUy5d1jLA6f/4519jYZtasvzj7N7m709y/kEjR6nK6SNFL5nFBPfdLM5MmTWKxWLNmzfLy8jIyMhJxbhaLXr3QtOkSM7N4VVWmq2vXkJA1AE6ehLk5zp5Fv37ijo+QcjExmDgRERHo21e4A4WEwMUF2trcn/3nH/j4QFFRoCEUFeHhgdBQrFwpUD9VBAauDQz87/GwYRg5EgkJSEyEoSGdoxDEr7hfGo2NjX316lVAQEADyYIA3r17l5OTUFJy5e3by8ePp7558wZA+/YYMADu5M49UScVFRXPnj0rKirauxeenoiOFnoWrKjAtm2YMYP7s3l5iIrClCk0DOTjg7AwlJbS0BVXgwcjJwfq6jAxwZYtwhqFIDiqXVBP184v0qKyspLBkP/2k3xFRQXn0ZEj+PoVixaJKy5CWuXn55uYWNvYLNHRsfHzu5aYiK5dhT5oVBSMjP5vf4mfBQXB3R0aGjQM9NtvMDPDiRM0dFUdNTXcuYNlyzBzJrp1g5Rve0NItGpnrcXGxp44cSIrK6u8vPzn43FxccKPSgy0tbX79293+vSA9++Z1tbNW7ZsyTmuoIAVK+DnB19fqKmJNURCqgQFhWZk/MlmjwZeGRlNMzE5JYJBN26sdqrzly8IDcWNG7SNNWUKNm2CsMsP+/nB2Rn29tDWRkwMevUS7nBEw8T9G+Hs2bMHDhx45swZad99kC/h4cE3bmxycgoaMOD/Nl6bMwfa2nWfbk40TBRFfStqz5CVFcWK1Nu38eYNBg3i/uyOHejbF8bGtA03eDDS03H/Pm0dVqd9e+TmwskJvXtzKR1OEILj8o2wsrJy+/btPj4+GzdulBFkkrUUMjIymjgRW7bA8/93JIyKQteuSExEz55iioyQNl27ejKZA9XUYmVlH2zYsFEEIwYFYfp07usiysuxcSOOH6dzOFlZTJyIkBDutUzpxWTi8GHExMDNDefPkxk0BM24fCN89+5dcXHxhAkTGloW5HByQmoqsrP/76CVFeztMWqUmGIipM2rV5gyRfvs2WuJiX7Pnyf27GlT+zmCyc9HbCzGj+f+7IEDaN0av/9O86De3jhwACKrzO/sjOzs/2bQiCD7Eg0Hl0SooaGhra1dZSemhoPFgrMzjh2revz4cbx7h+XLxRETIVU4dUQXLoSDg2zbtm1VVFREMOjWrXBzQ9OmXJ6iKAQG1qXEdq10dWFri4MH6e+5Ot9n0MyaRWbQELThkghlZGQCAwMXL1786tUrkccjEdzccOhQ1YMqKvD3x/LlKCgQR0yElKAoTJyI33/H1KmiG7SsDDt2YPp07s+ePg2A/i1/OaZMwdatQum5Bn5+SEvDy5fQ0gKpQUMIjvus0ePHj+fm5rZu3drU1FRTU/Pnp+rrrNGf9emDP/7Aq1f4NnX0P35+CA7GqFGIjRVPYITkW7YMr18jPl6kgx4+jA4dYGrK/dl16zB/PoS0HsrBAcXFuH5dFItDfta2LbKzMW4cHB0xahRbTc03Pj7V3b2fn998kcZB1AvVLp/o0KGDKOOQKLKycHXFkSNcKvQfPgx7e6SkwNJSHJERku3kSYSG4uZNsFgiHXfTJixdyv2p5GRkZGD4cGENzWDA2xvbtok6EQJgMhERgREj4OKyhs3+ACxcuHBZcXH5smVk2S/BH+6J8MiRI0IdlaKo69ev5+XlWVtba2lpcW3z4sWLO3futGnTpm3btpwj7969y8jI+N7A1NRUSUlJSBG6ucHXl0sitLVF9+4YNgwvXwppZEJaPXgALy/ExtK8y26tkpLw6VO1NWvWrsXs2cLa5oJjwgQYGeHtW/z/xSMRcXYGi7WnuPgqoEFRKrt3+5JESPCr9m2YaEdR1IgRIyZOnLh///527dpdvXr11za7d+/u2rVrZGRknz59VqxYwTkYHR3t4ODg/c3PSZF2trbIy8OjR1yeOnECb94gKEh4gxPS58MHuLhgwwZ06iTqoTduxPTp+HUbegAvXuDy5WqnktJFTQ2DByM8XLij1EBfXws4CpQCkVlZo7p1Q2qq2IIhpFJ11bivXbvm5ubWvn37tm3bco4EBQXt3LlT8HLg8fHxzZo149QvDwoKsrW1rdKgpKRES0vr0qVLFEU9e/ZMUVExLy+PoqidO3e6uLjwMkQddp/41cyZ1NKl3J/y9aVYLEqqNuSgB9l9gquKCqp/f2rePDEM/eYNpa5Off7M/dnJk6klS+rSLb/bC9y4QRkbU5WVdRlLcLm5uS1bdpGTM/j99z6XLpV27UoxmZSREXXsmChGl6KtGOgiRS+Zx90nuH8jjImJ6dGjx8OHD42MjD5//sw5yGKx/P39KYH35j1x4sSAAQNUVVUBjBw58vLlyx9+3h4eSEpKYjKZvXr1AvDbb7+Zm5uf5sx7A4qKipKSktLT0wUPo1ZubtXOC1+3Dqqq8PAQdgiEdPD1RUUFVq0Sw9CbN2PsWDRqVPU4RVGvXn09fJieEtu16twZamq4cEEUY/1KW1v75cvrZWUZKSnn7ezkr11DZibMzTFiBBo3xuLFYIuisA8hxbjfOvjzzz9HjhwZHh6emJjo/m3zBXt7+ylTpmRnZzdr1kyQIbOyssy/VQXW1tZmsVhZWVlNf1oA9ebNm+bNm3+v+q2vr8/ZCwLAs2fP/Pz8Hj9+3Lp16xMnTjTlumwK+Pjx44ULFz59+sT5sVGjRm78l0S0skJJicydO2wzMy5Jd98+ODnJ3L5dWV2B43qpsrKSyWSKqyB7ZWVlZWWlWIauwb59jNOnmUlJlQBEHF1xMXbtkrlypeq7cuXK1VGjphcUNGrcWL1x44OVlfLVdFCtOrzVXl6MrVsZvXtLRM7R0cHRoygshK8vc906xj//wN2d2rCBLYxJBZL5aylUUvSS2Ww2L9+auCTC/Pz8Fy9eHD16tMpHHif/5ebm1poInz9/PpFbTcDNmzebmZmVlZXJ/nTvXk5OrvT/d3MpLy/n2mDUqFETJkwAUFxcPHDgwCVLlmzevJlrAF++fPn5WyOLxXJ2dq5DoZwhQ+QOHIC/f/mvT/XoAQsLBTc3xt27Jfx2K71KS0vZbLZY/g1QFFVWVlYqvI1/6uTOHebcuawzZ0oUFdmiDy08XNbKCs2bV31XvL0X5uaeA3QqK5dHRBxwd+e7JFId3mpXV/j5KaanlzVvLvSrNTySk0NwMIKDsWOH7OrVcrt2yXTvzg4JKWvRgs5sLYG/lsImRS+ZzWbz8snPJRFykt+vWTQ7OxsALxM1dXV1V3G7TsTZ0kFXV/ft27ecI0VFRYWFhXr/P81OR0fn3bt33398+/atjY0NAMVv24kqKiq6u7tv27atugAMDAxcXV1HjhxZa6g1c3fHsGFYs0aO61egmBjo6+PAAaVJkwQcR2owGAw5OTlZoc5BrH5oRUVF4c0TroPcXIwahZ07YWmpIJYAduxAYCCXf5JlZeWAGoCKCu3i4tI6vGmVlZX8nqWkBHd37NunuGwZv6MJ3Z9/4s8/cfgwFixgtmun0KEDtm+nbb1HHd4raSdFL5nHP9y53CPU1NQ0NDSMiIjAt6TIsWXLFk1NTRMTk1o7VVZW7s4N575g9+7d4+PjOYn20qVLxsbGurq6nIjZbDYAKyur169fcyaFFhUVXbt2jZMIf/bw4UPOWUJlYQEWC8nJ3J/V1YWXF2bORFmZsAMhJE5pKVxd4e0NZ2fxBHDpEiorYWfH5al587xZrP5KSr7Nm+8YPVpoSwh/MWUKwsJQzuUCikRwc0N6Om7fhpoauneHvj7Cw8Fms2fNmte1a/8DBw6IO0BCfLhOodm7dy+AMWPG+Pv7a2lpRUZGDh8+HEBwcLDg03iKi4uNjY0nTpwYGhpqYGAQEhLCOT5gwICFCxdyHk+ZMqVz5867d+/u27fvgAEDOAenTp26dOnSkJCQKVOmKCkpXb16tbohaJk1yvH339Ts2TU1UFOj3N1pGUoKkFmj302aRA0ZQrHZYgtg0CBqxw7uT1VUUE2apB8+fKmwsLBunddtWmBFRUXz5ot1dKzd3afXeWjRePOGcnWlZGQoGZm/gHHAcQaj3d69e+vQlRRNoaSLFL1kHmeNVrt8IjQ09OfiaioqKgEBAWya/t3n5eX5+/tPnTr15MmT3w9GRUX9+++/nMcVFRWhoaHe3t6BgYFFRUWcg5cuXVqwYIG3t/eKFSs4twCrQ2MifPSI0tOraV74iRMUk0k9eULLaJKOJEKOoCCqY0dKjB/16emUpma1C3iuXKEsLATqv26fdIGB21iseUCJnNwWLy9xrCbh05cvFJPZBigDKOCUisrEadOo/fupvDw+OpGirEAXKXrJgiZCiqJKS0uvXbsWFRWVkJDw5csX+mITOhoTIUVR5ubU5cs1NWjdOpPFGqOhYb5+fSBdg0omkggpirpwgdLVpV69EmcMs2dT8+dX+6yfH/Xt2kod1e2TbtSoacA1gALyLSz6CRSBqGhqmgGngEpgsq7uWgMDSlGRAihZWUpdnerQgXJzo9ato27e5H76unXrRo0ac7O6p+up+pcIa5r1IC8v31X0BQQlj5sbDh+uaUve168dSktXlpbOmzt3lIWFuR3X+zZEvfDqFTw8cOAAWrQQWwyFhdi7F7dvV9sgJgZhYSIM6JsRI/qePbvi48c5ior73NycxBAB/y5dOmBrO/Lz58kdOrS5dWvL9wI9qamIi8O9e7h3D3FxWLAA5eVgsaCpiebN0aEDevRAYKBzaqo8RdkfPjz0ypVD1tbWYn0pRN1xT4SJiYll3GaANGrUqGXLltVVB62vRo5Ely4IDq62YGNJSTkwDACbPfb48eMkEdZXhYVwdsaiRejVS5xhhIfDzg76+tyfzcxEfj6srEQbEwDAxWUgiyW3YsVpWdm+8+aJbpKOINq3b//+/f1fj1tYwMLi/45kZeHCBVy7hnv3EBuLvXtRUpIL3ACYbLbS2rUbT5wgiVBacf9oHzFiRF5eXnXndO/ePTw83NjYWGhRSRZDQxgZ4dIlODpyb6CoKPP161GgDXB48OB1oo2OEBGKwtixsLaGj4+Yw9i0qaYvfNHRGDCAe+lREejfv6+RUV9HR2Ht+iRGzZph3DiMG/fjiJLSh+Lix0BbIOXFi94ZGeK8TkAIgvs/l23btqmrq/v4+Jw5cyY5OTk6Onr06NHNmjU7efLk1q1bX7586ezsLC2VBWjBuTpanVu3YkxMghs3nshgrGEye4swLkJ0li5FXh42bhRzGGfPQlkZv6wn+iE2FgMGiDCgX7RujcrKBrE9y9GjG1msAQyGgYLCvVGjJnbqBG9vVP8NgpBgXO8cduvWbe3atVUOTpkyZfTo0RRFJScnA7hy5YrgdzKFhN7JMhRFZWdT6upUSUktzZydqSZNxFZ6WAQa4GSZL1++ZGZmRkWxmzWjsrJEP35VfftSNUzyLyykVFWpT58EHUXA2RCjR1NhYYLGIC0+fiwwNaXi4qi3b6n58yl1dWr+/GrLoNcP9W+yDJdvhB8+fLh27dqgQYOqHB80aFBsbCyATp066erqvmwIf/J9o6uL9u1x7lwtzaKiUFEh9F1vCJGJiooxNra1sJjm5uZ46FCxiDca/NXTp7hzp6ZdduPi0LUrGjcWYUzc2NkhPl7MMYiMjAyWLIGfH9TVERCA27fx8SNMTLBmDaSkDBnB7dIoRVEAnj9/XuX48+fPqW911+Tl5RUUxFNWSlxqvjrKISuLgwexb1+1xWgI6TJnzsr8/Pj3708C9unpkeIOB0FBmDwZNfzLE/t1UQ57e1y6JO4gRMjNDeXlOHUKAAwMEBKCs2eRkIB27XDoEIS/Uw4hKC6JUF1dvUuXLlOnTk1MTOQcoSgqOjp68eLFTk5OAN6+ffvmzRtDQ0ORRipuw4cjNhZfv9bSbMAA9OghtrJbBL3YbIozoYzJlBf7TfFPn3D4MCZPrrYBReHMGYlIhEZGkJfH06fijkNUGAz4+2Phwh/7PXXsiDNnEBaGoCBYWtZ+MYkQL+6TZfbu3SsjI2Nra9uoUSNjY2NlZeXBgwcbGRkFBwcDePbsmZeXl0WVycX1nYYGunZFbGztLU+dwufPmDVL+DERQtahw2w5OXtNzbFGRidGjBgm3mB27sSAAdDRqbZBSgpUVfHbbyKMqXq9ejWsL4UuLlBSwrFj/3ewVy9cv461azFvHmxs8O+/YgqOqA335RMmJib37t2LjIxMS0vLzc3V19e3tLR0cXHhbDtgbW3dMJeOcq6OjhhRSzMVFWzdikmT4O2NNm1EEhkhBCdP4v59t7Q0+/Ly3LZt29ZhGy8aVVZi61YcOlRTm1OnMHCgqAKqjZ0dTp+u6ftr/bN0KWbMgKtr1QXHDg5ITcWxY/DwwG+/YcMGmJmJKUSiOsKetCMWtM8a5fj8mVJT43VKnqUlZWREewhi1nBmjT5+TGlrU8nJIhuwFlFRVPfutbTp1ImKj6dnOMGnBWZmUlpa4ixKLjI/v1e2tjXN6S0tpUJCKB0davhw6uVLKjHxyvTpC/fuPVApbRPNG8SsUaI6jRrB1hYnT/LU+OxZvH6NRYuEHBMhBIWFGDIEq1ahUydxh/LNxo2YMaOmBjk5eP4c3buLKqDa6OtDVRUPH4o7DtFasQJLllS7L5u8PLy88PQpLC3RseNVR8dFmzZ1nzo1cdGiNaINk6jqRyI8fPiwnp4eZ893c3NzvWqIL1SJwMvcUQ4NDaxdi4AAZGQIOSaCVpwKMr16YcIEcYfyzf37ePYMrq41tYmNRb9+kJMTVUw8sLNrWLcJAdjYoFUr7NlTUxtVVcyfj6FDz5aU/AX0//IlKDKSh6kHhDD9uJhtaGg4fPjwNm3aABg0aFBBQYH4opJcgwdjyhS8fYuftqiq1p9/IiQEAwfi3j3hR0bQZPVq5OfXcjdOxAIDMW1aLUkuNhZDh4oqIN7Y2SEyEtOnizsO0Vq1Ci4u8PCAomJNzWxtTY8cOVVY6MBgnGjXzlRU0RHc/UiEnTt37ty5M+fx8uXLxRSPpFNSQr9+OH4cXl48tb9wAS1bYtOmBvdxIKUuXMDWrbhxA/Ly4g7lm3fvcOJELUsRSksRH4/QUFHFxBs7O0yfDjZbbIVPxaJTJ/z+O0JDa7mU7eEx6t6951FRNl+/tmvThhQoFrOG9BtKE96vjgJo1gwLF8LXFx8+CDMmgg4ZGRg7FgcPolkzcYfyk5AQDB0KdfWa2sTHw8wMGhqiiok3urrQ0kJamrjjELmVK7FqFQoLa2rDYDDWrVuSnn799u2du3Y1vXVLVMER3FSbCE+ePNmjR4+mTZs2b96cc2Tt2rVBQUGiCkxyOTnhzh1kZ/Pa3t8fOjoSscyZqEFxMYYOxaJF6NFD3KH8pKICISGYNq2WZhJSUOZXDarW2ndmZrC1xZYtPDXW08P69Rg3DiUlQg6LqB73RBgeHu7i4qKgoDB48ODvB3V0dFavXi32+hpix2Jh0CBE8lNvKy4ON28iIkJoMREC8/GBqamYt1j6VWQkWrVChw61NIuNlaAVhD9rmIkQwPLlCAwEjxMt3N1hagpyP0qMuNcaXbBgwcyZM+Pi4v7444/vx7t3756fn5+VlSW66CQVX1dHAZiYwMsLXl4oKhJaTIQAgoNx5w5CQsQdx0/OnIkzM7OfNKlfv343am55/z4qK9G+vWji4o+9Pa5cQQP849nEBH37IjCQ1/bbtmHPHtyo5X81ISxcEmFeXl52dvb4X/ZQ0NHRAZCfny+KuCSbgwPS0/HqFR+nbNsGNTX89AWbkBRJSQgIQFQUlJTEHco3Hz58+OOPBffvH/z6dfv69d4VFRU1ND51Cr9sFSMp1NWhr4/bt8UdhzgsXYotW3idHKCpiS1b8McfKC4WclgEN1wSoby8PIDiX/6HvHr1CkBjsW/xIgFkZeHqiiNH+Dvr9GlcusTrenxCNHJy4OaG3bshUTXkMzIyKKojoA20BFrU/NenxN4g5GiwV0dbtsSQIVjH84RQFxeYm2PxYmHGRFSDSyJs2rRp27Ztt27dSlEUg8HgHKQoas2aNc2bN/9NQmr6ipubG99LzSwsMGIE3N2rLTxBiFh5OUaMwPTp6NdP3KH8P1NTUxYrFTgkK7urSZP3urq61bX88AH37qFXLxEGx6cGmwgBLFmCsDA+9qzfsgUHD+Lbrj+E6HCfLBMQELB//34HB4ejR48WFxdv2rTJ1tY2IiJi9erV31NjA9ezJ/Lz8egRf2ft3w95eYwaJZyYCD5Nnw4NDcydK+44fqGgoDBuXHSHDs8XLXp39Wp0Df/oTp+GnV1NOxSKXa9euHq1gf7xp6eHMWMQEMBre3V1bN+O8eNrWXpB0K+6IqQxMTFtfto6oXnz5hEREfSVQhUuIRXdruLPPyl/f77POn+eYjKphAQhBCR89ano9t69VJs21OfPNHZJGzabatWKp5LfI0dSoaH0B0BvVeXff6euXqWxP8lS83uVn09paFCZmXx0OHYsNW2aoFEJVQMquj1w4MBHjx69fv361q1bT58+zczMdHd3F0lqlhp1uDoKoE8f9OsHV9cfe3gSopeaCl9fREWhUSNxh8LN2bNQUam95HdlJeLi0L+/SGISQEO+OqqpiYkTsWoVH6ds3IjoaJw/L7SYiF/UUlmmefPmlpaWrVq1IldEf9WlC0pLcfcu3yceP46yMnh6CiEmggfv32PoUGzdClNJLfG4bRtPNfmuXIGhoWTVweGqISdCAPPnIyoKL17w2r5xY+zcCW9vXpchEoIjJdbqjsHge0Ehh7w89uzBnj1ISRFCWESNKisxejQ8PCSuRPV3mZlISoKbW+0tJXy+6Hc9e+LGjYZbOaVJE0yZgmXL+DjFwQF9+sDXV2gxEf+PJEKBcK6OUhTfJw4bBmtr6fgUq2fmzQNFYckSccdRvZAQjB3L06JGidqSvgaqqmjXrkGvFp89G6dP8ze3LjAQly7hzBmhxUT8hCRCgXTsCEVF3LxZl3PPnMGnT5g+/UtSUlJZw5xUJ3JRUTh+HAcPQkZG3KFUo6wMu3bB27v2li9e4PNnWFoKPyY6NPCro40aYfZsLF3KxynKyggNhacnPn4UWljENyQRCmr48LpcHQWgooLhwyM2b7bp0eNvVVXjZ8+e0R0a8X/u3YO3NyIja9nJQbyOHUP79mjduvaW0dEYMADScu++gSdCANOnIzERd+7wcYqdHVxd8eefQouJ+I0y8QEAACAASURBVIYkQkGNGoVDh+pYTTEqahkQw2bHlZUt9/QkNwSE4suXLzt37t6167Cra/n69fj9d3EHVKNt2zBlCk8tpeUGIYeNDW7fxtev4o5DfJSVMX8+/P35O2vNGly7hqgooYREfEcSoaBat4a2Nv79ty7nUhQFcPYdZ1VUNLzKxMJXVlZmaeno4/PRy+thcfGQsWPFHVCNHj5EejqcnWtvWViImzfRu7fwY6KJkhI6dsS1a+KOQ6wmT8bt27h+nY9TlJQQHo6pU0FqPAsVSYQ0qNvcUQB//z2ZyezOZLoB60xNN9AdF4G0tLSPH83LymZXVi6tqCj9INn7I2/ZAi8vyMnV3vLcOVhbQ1VV+DHRh1wdZbGwcCH+/pu/s7p1w5gxtex3TwiIJEIajB6NgwcTDx48UsDnwp/5831fvIg/fHh4aOiVXbtM6pZNiRro6OgUFz8ESoFPQF4jyVw/DwAoLMShQ5g0iafG0nVdlIMkQgATJiA9HQkJ/J21fDnS0nD0qFBCIgDIijuA+mDlynmFhbl//NFGV9f+7t2LfG3Q0aJFixYtWgB4+hRjxsDICFZWQgu04Tl2rLmy8vhGjaxZLLlNm9bKykruL/y+fbC352l1PJuNM2ekb5uCbt1w/z4KCiS0mo9oyMlh8WIsWsTfzRRFRezbhwED0KMHdHSEFlwDRr4R0uDkyYsVFXvLyha8fet64cKFunWydi0cHWFry0eteqJmO3ciMBA3b47Pzk55+fL6wIF9xR1RTXbs4HWaTHIyNDQka98oXigooFMnXL0q7jjEzd0dHz7wXUHt998xYQJP62qIOiCJkAYKCjJAPgAW64mWllad+zl9GoaGMDNroKX66RUejmXLEB+PFi3EHQoPrl5FYSHs7HhqHBsrHevof0WujgKQkYG/PxYu5LsQx99/49Ur7NsnnLAaNpIIabB79z8GBgPk5CxsbPR69OghSFfJyWCzYW1NV2gNVGQkFizA2bNS87Vp2zZMncrrosBTp6TvBiEHSYQcw4ejogLR0fydJS+PvXsxezZevxZOWA0YSYQ0sLPrmZGRvGRJqrExzzuPVUNJCXfv4sEDnkpNElydOIEZM3D+vOTW1K7i3TucPg0PD54aZ2cjMxNduwo5JuHo3BlPn5JSKWAwsHQpFi1CRQV/e9CYm2PaNEycWJeyjkQNSCKkjYsLjh+n4Re0WTPExSEyku9p1gSA8+cxeTJiYtCunbhD4VlYGIYMQdOmPDWOiUH//pDgST81kZdH1664ckXccUgAB4eijIxB6updjY27PX78mPcTFyzAp0/YtUt4oTVEJBHSpn17sFhIS6OhKxsb7NqFFSsQGUlDbw3HxYvw8EB0tNRU4ATAZmPHDkyezGt7aVw48TNydZRj69adJSW9Cwpuvnixw9NzAe8nysoiPBx+foUBAbvCwyNKGuymHrQiiZBOgwbhxAl6uho3DjNnYuRIslUTr5KSMHo0jhxB587iDoUfZ89CU7P2PXg5iouRmIi+Ej37tRYkEXLk538sL+fM49J//56/i8WtWlUA/RYu/Dh58ptu3aRz3pSEIYmQToMH4+RJ2nrbsAG9eqFnT1JdqXY3b8LVFQcOwNZW3KHwaetWXldNALh0CRYWaNJEmAEJmaUlXr3Cu3fijkPcvLxG6+gsZbGWsVguCxbw/BsAAHjy5AlgzGbPKSnxy8lRefPmjZCCbDhIIqSTtTVycvDyJW0dXrgAXV1YWKCigrY+65+0NAwejJ07pan2JkdmJm7cwIgRvLaX9uuiAGRlYWODy5fFHYe4/fbbb/fuxe3caS4vH+rszPNvAABAR0cHeAIUAQVsdoa6JG+nIiVIIqQTk4kBA/ieFV2ztDQUF0OwRRn12ZMncHJCcLBULq3btg3jxvG0By/H6dNS+TKrIFdHOTQ0NMaMGdy7d6vjx/k7UV1dfd26WQYGvZSVHfv3X66oqCicABsQkghpRu/VUQBKSkhOxu3bGD2azm7rh+fP4eCAtWv5+FIlOcrKsGcPvLx4bX/3LmRk0KaNMGMSCZIIf+buXpc18uPGuWVk3Lx9+/qZMwM/fRJCWA0MSYQ0c3REairNt0CMjXH+PI4cwdq1dHYr7TIz4eiIxYul9U+EyEh06AATE17bnzqFQYOEGZCodOyIvDzk5oo7Dsng7Iy0tDqukTcxwcCBCAykO6aGRwyJkKKouXPnNm3atEmTJr6+vmx21SWlHz58mDlzZvfu3Y2Njd+/f//9eHFx8ZgxYxo3bqytrR0cHCzaqHnFYqF3b8TG0tytrS02boSfH22zUqVdVhbs7eHry8c3KknD+x68HPXgBiEHk4kePfjegaG+kpeHqysOHqzj6f7+2LKFzKcTlBgS4f79+6Ojox89evT06dPY2Nh9v1wXKC8vV1dX9/HxefHiReVPW7+vXr06JycnOzs7MTFx+fLlN27cEG3gvKL96iiHjw+mT8eIEbh/n/7OpUt+Pvr0gZcXfHzEHUpdPXyIjAw+bvi9fYuHD+vPrWJydfRnHh4ID6/juQYGGDMGa9bQGlDDI4ZEuHv3bh8fH21tbU1NzWnTpu3evbtKA21t7SVLlvTr16/K8V27ds2bN09ZWbl169ajR4/+9UQJ4eyMS5dQVER/z0FB6NYNXbs26Nnn796hd2+MGYN588QdigA2b4anJx8FYk6fRp8+YLGEGZMIkUT4s+7dUVKCO3fqePqiRYiIIAVIBSKGRPj06VMzMzPO4/bt2z979oyXs4qLi7Oysng8kaKor1+/fvymuLhY8LB5p6aGTp0QFyeUzuPjoaUFCwv8ckW5Qfj8Gf37o29fLFwo7lAEUFiIw4cxcSIfp9Sb66Ic7dvj82dkZoo7DsnAYGDMGERE1PF0TU14emLFClpjamCEUrIwLS0tMTGxykEmk+nj4wPg48ePKioqnIOqqqo/3wWswYcPHwDweOL9+/dPnz7t6+vL+VFFReXevXsyMjJ8vo6669dPPjKSaW8vlOpHV6/C1FRZQ2NBRUWMmZlhVNROJd4n4AumuLhYTk5O9NvbPnz4sLCwMDb2UljYCCuryr//Lv3yRcQh0Ck0VK5XL1lV1WIeX0V5OeLiVFav/vrli4hqLRcWFgp7iO7dFc+dqxg5slzYAwkbLe/V0KFMR0elRYsK6/Zvy8eHYWGh7OVVZGIiij+QRfDrQRc2my0nJycnJ1dzM6F8on369On58+dVDn7PQxoaGgUFBZzHnz9/5nEDPw0NDQaDUVBQwNn/veYTzczMlixZMnLkyLpET4dRo7BmDRQV5YSRMlRV0bOnT0xMJXA6KWmbq+vEpKRT9A/DjaysrOgTYUxMjIvLNDabMXHiPx06lGzZMp7BkBdlALTbswcbN0JVVZXH9hcvok0bGBurCDWqKngPr2769MG1a7KengpCHUU0BH+vzM1haIibN1XrVj9PVRWzZmH9euUDBwQMhPcRhfvrQRc2m/3zRJPqCOUTrWfPnj179qzu2TZt2ty9e9fBwQFAWlpa69ateemTxWK1aNHi7t27+vr6fJ0oFs2aoUULJCWh+rdBILdvpwBbgRbAzLQ0e6GMITH++msDmx0KjAX2PH8+lMEYL+6IBHLlCioq+CsFV8+ui3LY2SFA0F3L6hXOgsI6F5KdMQMmJrh7F+bmtIbVMIjhHuGkSZM2b96cnp7+8uXLjRs3Tpo0iXN81KhRd+/e5Ty+ffs25/Hdu3dTvpWdnjRp0qpVq96+fZucnHzo0KEJEyaIPnjeCWnuKIejoxWDsRZ4BPxTWjokNVVYA0mCykoVgFO27o2iotR/geCsmuBxD16OU6fqQ0GZKtq0QUUFnfUIpd2oUYiNRZ0vOiorw88PixbRGlPDQYnDihUrDAwM9PX1ly9f/v1gr169bty4wXncuXNny286derEOVhaWjp16lQdHZ1WrVrt3r27hv5HjBhx8OBBoYXPk7Q0qmVLIfY/cuQkdXVzR0c3R8cyJpNauFCIY31XVFRUXl4uipEoiqKonBzKxYUyNc1WVDQCZOTk9C5duiSy0YUhP59q0oT68IGPUx4/pvT0KDZbaDFxU1BQIIJRRo+mdu4UwTjCReN75exM7d1b99PLyigjI+rqVbrCqZZofj1oUVlZWVZWVmsz8SRCYZOEREhRVKtW1N27ohho1y5KTo4yMqLy8oQ7kCgT4ZEjlLY2NX8+VVpKURSlqamZm5srmqGFZ+VKytOTv1P++YeaPFk40VRPNJ90oaGUu7sIxhEuGt+rI0eoPn0E6mH3bqpHD5qiqV79S4SkxJoQOTuLqBDM+PH/XWIyMMDhw6IYUagyMuDoiJUrceYMAgIgLw8ADL4uJkokNhuhofD25u+senmDkMPODpcuiTsISeLsjJQUCLKrkocH3r/H+fP0xdQwkEQoREK9TVhFs2ZIT4ePD0aPhpOTtG7bRFHYsQNWVrC2RnIyLCzEHRCtYmOhowNLSz5O+fwZKSmwr6fToYyNISeHp0/FHYfEUFDAkCE4dKjuPcjIYNkyLFgASkQLbeoJkgiFyMYGWVkinQ6wYQOuXMHVq9DWhtTNoHn5En36YPduXL4Mf3/UtvJH+vBbXBTA2bPo2ZOPfZqkTq9epMTM//HwwJ49AvUwZAhkZcHv1k4NHEmEQsTZnjAmRqSDWlsjLw9mZujUCX//LdKh64zzRbBLF/Tpg3//hampuAMSgowMpKTwvV1UPb4uykFqrVXRoweKipCWVvceGAwsXYoFC6T1spBYkEQoXKK8OvqdggISEhASgtWr0batpBcmTU+HvT327sWVK5g/HyKs/yNSW7di3Dgo8Lz6Iycnx89v5bFj66ytPwozLjGzt0d8PLmO9wODgVGj6rJD4c/69oWeHvbvpymmBoAkQuHq0wcpKeJJRZMmIT0dxcVo3hyRkWIIoFZsNnbsQLdu6NcPiYmQ4AIJgiotRXg4PD15bV9SUtK1q/OaNYbFxU3d3OrFJoTV0NeHigoePhR3HJJk3Djs3w8eyqHUJCAA/v4oLaUppvqOJELhUlRE7944fVo8o+vr4+VLjBsHNzcMGSJZdbofPoS1Nfbtw9WrmD8fzHr9m3j0KCws0KoVr+0fP35cXGxBUaMpauLHj02zs7OFGZ2Ycb4UEt+ZmEBPT9D3pHNntGuHsDCaYqrv6vXHj2QQy9XRn4WEID4eFy5AW1ugew90qajAmjXo1QujRiEhgY/0IL34nSZjYGAA3AHygUwGI4PHerxSisyX+ZWHR903o/hu9WqsWoWvX+kIqL4jiVDoBg7EhQtC2Z6Qdz17Ij8f7drBwoJSV3eSkWnBYrWMiooSWQBpaffatbPT07P08FjSrRvi43HrFmbOrOdfBDnS0vD6NX9zXpo2bbp06Uo5ObeOHT2PHdsu+u0+RMneHpcvS9blCrEbNQrR0XUvt8ZhZoaePbFpE00x1WsN4HNI3Jo2RadOuHhRzGFwZtDY2a378MGQzc4oK4sbM8ZXZKOPHDn94cMdOTm3Dhx406tX3NmzMDAQ2eBitnUrvL35ngSUn+84fXp8auo5a+uuwolLUujqQkNDIq5VSA5NTXTvTsOVpGXLsH49Ptbn6Vb0IIlQFMR+dfS7pk3TAE51euOSErW2beHjg+RkoYxFUXj0CDt3Yvx4PHtWCPwGMBiMTiYmDWU/1sePHwcEbD54MHE8/xtmREVhyBAhxCSRyCKKX3l4CDp3FECrVnBxwfr1dARUr5FEKAqcRCjgNDBazJ07k8lcDgQyGENMTQ0tLXHuHLp1g4wM9PXh4SFowi4vR0oKgoMxYgS0tNC7N86dQ8eOGDSod6NG3jIyIZqaoQMG9Kfp1Ui0W7dSevQY5+enUlwcdPToDr7Off4c+fno1k1IoUkckgh/NXgwkpORmytoP0uXIiQEeXl0xFR/kUQoCi1aQF8fSUnijgOwsrJKTY2dOPFBUJD9w4fHIiKQno6KCly6BCcnpKRg2LD/kuKQIdi37//u3LDZ7MOHD8fFxVXpMy8PMTH46y/Y2KBJE4wdi4cPMXAgbt1CdjaOHMHMmYiMDDh0yHXTJsbt26f19PRE+prFZM+eqHfvlgN/lJeHh4Qc5OvcY8cwZEiDuIHKYWeHK1ck4i9FyaGggMGDIfhGu3p68PDA6tV0xCSF/v3335ycnNrbiaD+t+hJyO4TP/P3p+bMEXcQvLl9m5o1izI1pVgsismktLWp/v2p8PDyRo1aMRhDGIwe7dvbpadT4eGUlxfVti3VqBHl4ED9/TcVF0d9/SqsqLS0tKRr94ng4O3y8v4AxWCcGzBgLF/nWllRFy8KKS6eiH57gfbtqeRkEY9JD+G9V5cuURYWNPTz9i2lrk69eEFDVxzSsvuEvf1wJtMyMDC11pYN5m9OcXNxwbFj4g6CNxYW2LABDx+ipARXrsDFBenpmDDhZkGBDUUdo6jE+/e/9u1bFB8PKyscPYpPnxAXB39/ODjU56qY/Jo8eYKSUk7Tpp26ddsSGsrHXuxv3uDlS/TsKbzQJBHZieJXvXrh40fcuydoPxoa8PHBypV0xCRVEhJusNk3ioo61NqyPk/Llijm5mAyce8ezMzEHQo/rK1hbQ0AJ07kDxmSQ1EAKpjMt3fvkpxXi9hYOUPD7Skp/G1GDyAyEoMGoV6vmODCzg47dmDePHHHIUkYDIwcif37EcDH31HczZkDExM8elQ/C/lyVVICilIEyiorWbU2Jt8IRWfwYBFtTygMLi4uLVt+YjLNGAyT5s17KpE0WCOKwqpV8PfnOwsCiIrC0KFCiEmy2doiKQnl5eKOQ8KMHYt9+2i4e9q4MebMgb8/DSFJhcuXoa0NeXlfBqOtquq1WtuTRCg6krOIom5evLj26NGxl/9r707DmjjXPoD/k2AgrBKsgcgicgm8CIqoKCoCApbKsUhrtVYsotBWq1Kpu1Wk1NrSo9Ricamt4FZX3LWVRRaRgLiACmKreASPLCKCgIAk836YczjUIksyMCE8v0/MZOaeGy7xZmae57kLJerqMcq5eKnyOHUKTU2YMqXTJ5aW4uZNlW1A2AahEBYWyM5mOw8l83//B5EIKSkMhFq8GBkZuHaNgVDKTCbD7NmYOBHu7qirC3z2LMff36bds0gh7D7jx6OoCA978iQ6ExOTAQOEsbFYtIiBgd0qbMMGhIbKczt4/Di8vTvRpEJlUBQlEh2aPz94//7DbOeiXBiZUAhAQwOrVmHtWgZCKa3cXBgZ4fhxnD+PuDhwudDW1tbW1m73RFIIuw+PB29vnDrFdh4Kc3LCvHkICCDdc1p35gxevMDbcjWN6FXz6FvasSMmNfVsTs70Tz89tWNHLNvpKJEPPsDJk8ys0RgYiIICJCczEEoJhYRg+HAMGoSyMkya1LlzSSHsVj396Wiz9evx9ClZ27514eFYv16eWYBPnyIrC15eXZCT0jt6NL6u7gtgXFXVuqNHL7CdjhLp3x+Ojsz8Ad2nD0JD8cUXDIRSKkVFsLFBdDRiY5GRIc84PlIIu9WkScjOVoWl/9TUEBuL1atx9y7bqSiZc+dQWwtfX3nOPXWq905BcXIaJhD8CtTw+QfGjbNnOx3lwkgzCtqsWaipwfnzzERTBhERGDQIAB4+hJ+fnEFIIexWAgFcXXH2LNt5MMHaGuvWYdYsMtLvL77+GqGhci4Kc+xYbxwvSlu7dklQUFO/flOGDZOtWfMZ2+koF19fSCTMvJXncrFqVd2cOV+9+eaHJ06cZiAie549g5MTVq/GF18gLw+KNCsjhbC7qczTUQALF+KNN/Dtt2znoTR+/x2VlXIWs+fPkZbWuW5NqoTP52/ZEn7kyEUO58s+ffqwnY5yEQgwZQoOHWIm2rFjIeXlOhcuLJ87d2tGRvtTC5TT0aMwMkJxMQoKEBqqaDRSCLvblClISEB9Pdt5MIHDwa5d2LoVWVlsp6IcvvoK69bJeTt45gzGj4euLtM59Shjx+LOHTx9ynYeyofBp6MSyVWKCgZsKyvnXbiQxkzQblRfj8mTMWMG5sxBUREsLBiISQphdzMwgL09EhLYzoMhYjGio+Hvz3LnYWUQH4+KCrz3npyn98559K/g8zF2rMoOa1SEmxvKynDrFgOhHBzs1NR2AYUCwV43tx7W7TI1FYaGyMqCRIJt2xgLSwohC1Tp6SiAd97ByJFYsYLtPNgWHo4vvpDzdvDFC8THyzMBX/W4u7PfxVoJcbmYORO/dq6LSev27Pk+KOieo+NyNTU/E5MesKbtvHmfqamZqKmZWFqGublh7FiUlGDUKCYvQQohC3x9ceqUSjWd2boVp0/j3Dm282BPYiJKSjB9upyn//YbRo1Cv36M5tQzeXiozvMSZvn7v9oZTT66urrR0RszM4+sWDEjJISJzLrSo0ePYmLOSKV/SKV//PHHqV9+eXjuHPMr8ZJCyAIzM4jF6LFvqVuhp4e9e/Hxx6ioYDsVloSHY906+X8/6QaEBIBhw1BVhQcP2M5D+djYQF8facy91Fu2DPn5OHOGsYBdobCwEDAHNAANLlckFhd0xVVIIWSHij0dBeDsjBkzEBTEdh5sSE5GcTHef1/O01++xPnz8PFhNKcei8MhLZlei8EhMwD4fERF4bPPlHrsnqXlWIqqB9YAq9XV89zc3LriKqQQsmPqVBw/znYSTNuwAffuMfmL2lOEhWHtWvlvBxMSYGMDsZjRnHoy8prwdT74AHFxePGCsYCenhg+XHlnQNXVwdaWO3Bg/JIlDZ9/3lhWdkuta/qT9bKmZ0rD3h4yGW7fxpAhbKfCHHV1HDgAd3dMmAAzM7az6S7p6SgqwqxZ8kfozfPoW+XpiTVrQFHyrFqu2oyMMGoUTp+W/23030VGwt4eM2fC0pKxmIyQyTB0KHg85OVpaGj8s0uvRe4IWdOj2xO+zpAhCAnB7NkqNRSobaGhWLNG/ttBqRSnT5Pnon9hZgZdXQY6s6skZp+OAjA2xrJlWLyYyZiMGDUKFRXIz++OZiykELJG9V4T0pYuBY+HyEi28+gWGRn480+FbgdTUmBqCnNz5nJSCWTs6Ov4+uLSJZSWMhlzyRI8fKhcjXHc3ZGXhxs30Ldvd1yOFELWODvjwQMUFbGdB9O4XOzdi02bkJvLdipdj17Ln8+XPwKZR98q8prwdbS0MGUKDjPatJHPx/btWLwYtbVMhpXb7NlIS0NmZve9YSGFkDU8HiZPVq6/wphibIxvvsEHHyj1aDTFSSQoKMCHH8ofgaJw8qScrSpUm7s7Ll1CYyPbeSglPz9mWvW2NGECnJzwzTcMh5XD6tU4cABnz2Lo0O67KCmEbFLVp6MA/P0xZAjWrWM7j660fj3WrFHodjAjA/r6sLJiLidVoa8PS0tkZrKdh1Ly8MC//407dxgOGxmJHTtQ0CXz9DqRw7ffYt8+eHp263VJIWTTm28iK0sV2hO26scfceAALl5kO4+ucfUq8vMxZ45CQcg8+jZ4eJCno63jcjFjBvbvZzisoSFWrMCiRQyH7bijR7F0KX74ATNndvelSSFkU3n5vwBPc/MRU6cGvFS5tn79+mH3bvj7q2alX7cOq1YpdDsI4MQJ8oLwtdzdyXiZ1+Lzd2/c6GBsPGLfviMMhg0ORmkp4uIYDNlRFy/i/ffx+ef49FMWrk4KIZvmzVtRVRVWVXX1wgWL7dt/YTsd5nl6YsoUfKZybVazs3HzJgICFApy7Rp4PNjZMZSTyhk/Hrm5qK5mOw/lU15evmvXT1Kp5NGjSyEhG+uY6/yipoYff8SSJaipYSpkh1y/jkmT4OeHiIhuvW4zUgjZVFz8GBgG4MWL4ffuPWI7nS7x3XfIzGR4kBvrvvwSK1ZAXV2hIGQefds0NDB6NFJS2M5D+ZSXl3M4FgAfEHA4xpWMPnIZPx7Ozvj6awZDtqOoCOPG4c03ERPTfRd9BSmEbFqwwE8ofJ/DidbTCw0MZG6tCGWiqYn9+7FoEYqK0NDQwHY6DLh+HVevYu5cRePExZEXhO0gkyhaZWVlJRIVa2qG8nirNDRkAwYMYDb+pk3YtQv5+cxGbd3Tp7Czw9ChLK/9TQohmxYvDjp37otp0zR8fI7a2tqynU5XGTECXl7Jgwfbm5i4envPlvbwVWfCwrByJQQChYLcvo2aGowcyVBOKopMq28Vj8e7cuW32Fi7HTtGNzaeyMlhOL5IhDVrumPUTGMj7OwgFOLy5S6/VttIIWTZ6NGjw8PnXrw4kKLYTqUrpaWtamhIKC/PSEszPKPkfV/adOMGsrIQGKhonLg4TJtG1tJsh4MDSkvxSDVfGihEXV192rRp8+ZNjYhQ64oJuwsXoqKia99oyGSws4NUirw8OdtZM4jt6xOAlRUEAty4wXYeXamxsQnoC6C+Xvz06TO205Hfl19i+XJFbwdBXhB2DJcLV1fSkqkts2fD1pb5Cbs8HrZuxZIlXThYafRo/PvfyM3tjqVE20UKoVJ4+22cPs12El1p/vxZBgZvq6uv4vMPTJ36NtvpyOn2bUgk+OgjReMUFqK0FE5OTOSk6shrwnbRE3aTkxkOO24cPD3x1VcMh6V5eeHWLeTmon//LonfWapZCJ8+fVpVVcV2Fp0wZUrPKIQ3b94skmt11DVrPktN/efRo65iccrly/pyRGhsbJTJZHKcyKDQUCxdCk1NReMcOQJfX/B4TOTUNRITE5VkYquHB+Lj2U7i9erq6lLYHtjarx9++QVz5zJ/9xYRgd27W1k0OCkpqVGB5e8++QRJScjI6I615isqKiQSSfvHUarIxMQkODiY7Sw6oamJeuMNqqiI7TzaExQU9OOPPyoS4fffKQsL6sWLTp/I5XJv3rypyKUVdPs21b8/9fw5A6FGj6bi4xmI03UGDRr0xx9/sJ3Ff5ibU3l5bCfxGpmZmSNGjGA7C4qiqI8/pgICmA+7dSvl7EzJZH/ZaWlpmZ+f7/FAFQAAD3dJREFU39lQdnYTuVxTLteCwzl//jxjGbbt+PHjb7/9druHqeYdIQCqRw0+4fHg5cXyAOIOUvAHO2kSHByUYm3fzgoLw7Jl0NZWNM6jR7h3D66uDKTUS5AlZjpi0yakp+PoUYbDzp+Phgb8+quicXbt2nXrlpZM9kAmy+JwFnh5MZFcB3Tw/yuVLYQ9Tk95Oqq4LVsQHc3y2r6dlZ+P5GR88gkDoY4dw5Qp8jfy7YXIa8KO0NJCTAwWLUJJCZNhuVxs3Yply6Dgu6ZDhyopyh7gAEKAYv01xytIIVQWb72F9HQ8f852Hl3PyAgrV7K5tq8cvvwSISEM3A6CLLTdeR4eSE6GcryyVGpOTpg3j4HBXK8YNQpvvYWwMHnOlcmwdi10dJCS8hGHc4jDCeVw3rOwGMBlfcLEX3F61iPEDjIwMNDW1ra0tGQ7kc7Jz18iFl/Q07vNdiKvlZ+fr6OjY2xsrGAciuLl5KwfNGifrm5HbwyTkpLGjBmjqfhIlU5qamqqreX/+WeEvf1aHu+FgtGkUvXr179xcFjB5Sp1q73Lly8PHz5coPg0EYbk5q61sIjR0lK6NtbV1dV37txxdHRkO5H/oCi13Nx1ZmZH+va9yWDYly/1cnJChw79ks9/BiAjI2PYsGHt/jI2NupeuRLV1KQrFp+3sNjV1FT/4MEDgUBgYmLCYG5tKy8vB3CjvdlpqlkIDx06xOVy9fXlGZ1ItKG8vFwgEGgzcmfUSYWFhebdMMiMAAA8ePDAzMyMQyb8t0cqlT569MjU1JTtRLpVD/rn0dDQoKmp6ebm1vZhqlkICYIgCKKDlOtBLUEQBEF0M1IICYIgiF6NFEKCIAiiVyOFkCAIgujVVK0QXrlyJaWFwsJCtjNSBWVlZcl/XdM3NTW1hKGJuwkJCbW1tc2bBQUF+X9rCVpVVTV//vy/n3vz5s379+83b169ejU7O5uRrHonqVSa8ldPnjxhOykl9fjx43YH5auMK1euPPpvNyyZTJaQkFBaWkpvNjY2JiQkKFvPbalUmpCQ0HLF6du3b9++/fqZaV26zlv3s7S0tLOzm/hfP/30E9sZqYLa2trBgwfv27eP3jx48ODAgQOrq6sVj0wvMHHr1q3mPcHBwYGBga8cVlJSIhaL/376e++9FxoaSn8dHR1tZGR048YNxbPqtaqrqwGMHTu2+TcoJSWF7aSU1I4dO8aOHct2Ft3E398/JCSE/vr69escDicsLIzeTE5OFgqFUqmUvexaFxAQ4OfnR39dXFxsYGBw6dKl1x2sggs9hYaGvktavTFKU1MzNjb2H//4h4uLi7q6enBw8N69e3V0dNjO63++/fbbbdu2JScn97hVFJTQgQMHzMzM2M6CUCKurq5bt26lv05OTp48eXJzz43k5GQXFxdlWykGQGRkpJ2dXVxcnK+vb2BgYFBQ0Lhx4153sAoWQqIrODk5ffjhhwsWLFBTU3v33Xc9PT3Zzuh/QkNDDx06lJaW1p0rVhBE7zFx4sTAwMBnz5717ds3JSVl8eLFc+bMaWhoUFdXT0lJ8fHxYTvBVujp6e3cuXPu3Ll37959+PDhiRMn2jiYFEKiozZs2GBtbc3j8fbs2cN2Lv+zbds2XV3d9PT0/krS4pMgVI6pqamJicmlS5e8vb0zMzP37dvn4OCQlZXl6OgokUi+//57thNsnZeXl5ub25o1azIyMtTV1ds4UunuZwmllZ+fX11d/fz582rGG4AqwNbWtrS0NC0tje1ECEKVubq6pqSk5Obmmpuba2lpOTs7p6SkZGZmCgQCW1tbtrNrXW1tbWZmppaWVnFxcdtHkkJIdEhDQ0NAQMB33303e/bsj5hb357D4Whraz9v0XSjurq6U28fnZ2djx07NmfOnIMHDzKVFUEQr3BxcaEHEk+YMAHAhAkT6E1XV1clfEFIW7ly5ZAhQ+Li4ubPn0+vvv06SvoNEMomLCzsjTfeCAwM3LBhw927d/fv389UZBsbG4lEQn9NUZREIunsH5ienp5xcXFBQUGkFhJEF5k4ceK1a9dOnjzp4uICYOTIkTk5ORcuXHBV1h7T6enphw8f3rlzp4eHh7e3d0hISBsHk3eERPuuX7++bdu2q1evcjgcTU3NmJiYqVOnenh4iEQixYOHh4fPnDmTXsL/3Llz6urqs2bN6mwQuhZOmzZNXV3d19dX8awIol33799vnt6qpqYWFRXFbj5dytTU1NTUNDU19fjx4wD69OkzbNiwhISE6OhotlNrRW1t7Zw5c7Zs2UL/HxUZGWlra3vixImpU6e2ejxv/fr13ZpgFzM0NBw5cmTfvn3ZTkSl5OXl+fn52dvb05smJiaWlpZcLtfIyEjx4BYWFu++++7jx4+fPHni6ur6/fff//21dm1t7Y4dOz7//PNX9uvr6w8bNkwsFtNx3NzciouLbW1teTye4on1Qlwu18TExNHRkc/ns52LstPS0rKwsDD6L7FYPHz4cLaT6lrW1tZeXl4jR46kNwcPHjx69GgvLy8l7MdUWFhoZWU1c+ZMelNDQ8PFxaW2ttbKyqrV40kbJqIHKC0tdXBwaF7bgiAIgkHkHSFBEATRq5E7QqIHaGpqKigoGDJkCNuJEAShgkghJAiCIHo18miUIAiC6NVIISQIgiB6NVIICYIgiF6NFEKCIAiiVyOFkCCITti9e/eVK1fYzoIgmEQKIUEQnRAcHNx2azeC6HFIISQIFlAUVVZW9vLlS/lOLy8vf/r0qeJpSKXSkpKSlt0/WqKT7HjXrbKyssrKyrYPqKmp6XSWBNHFSCEkCAZERUUZGhrW19fTm6tWrRIKhXv37qU3b968KRQKk5KSAOTn50+aNEkgEIhEIk1NzREjRqSnp9OHxcfHC4XC1NTUlpEjIiJEIlFz2du+fbupqWn//v0NDAxsbW2Tk5NbzWfbtm0GBgavtGH7+OOPbWxsZDIZAKlUGhoaKhKJjIyM9PT0xo8fn5eX13ykVCoNDw83NDQUiUR6enrGxsYHDx5sbGwUCoU1NTWbN28WCoVCobB5peLY2FhTU1ORSCQUCocOHUp/p7SFCxeOGjXq9OnTAwcOFIlES5culesHTBBdiSIIQmH0a7OEhAR6c+jQoXw+f9asWfRmZGQkn89//vw5RVFpaWlLly5NTEzMz8+Pj48fN26cnp5eaWkpRVEvX740NDQMCAhoGdnKysrHx4f+OiIigsvlrly5Mjs7WyKR+Pj4CASCvLy8v+dTUlKipqa2cePG5j01NTXa2trLli2jN4OCgjQ1NTdt2pSTk3Px4sXRo0cbGhpWVFTQn3700UdcLjckJEQikVy7dm3nzp0xMTFSqTQ+Pl4gEMycOTM+Pj4+Pv7u3bsURf36668Apk+fLpFIkpKSxowZw+fzr127RocKCAjQ09MzNTX9+eefL1++LJFImPmJEwRzSCEkCAY0NTUJhcLVq1dTFFVaWsrhcBYsWCASiWQyGUVR3t7ezs7OrZ5YUVHB4/F+/vlnejMkJITuVExv0jeLcXFxFEVVVVVpa2svXLiw+dz6+nozM7P58+e3Gtnb29vS0rJ5c/fu3QByc3Mpirp9+zaHw4mKimr+tKSkRCAQREZGUhSVl5fH4XBCQkJaDaujo0N/m81sbGxsbGykUmnzd6SjozNt2jR6MyAgAEBiYmKr0QhCGZB+hATBAB6P5+LikpCQsGHDhqSkJD09vaVLl0ZHR+fl5VlZWaWlpbXsC1pWVnbo0KHCwsLa2loAGhoaf/75J/1RQEDA5s2bT5w44efnByA2NtbAwGDy5MkA0tPTa2pqjI2NExISmkOZmZndunWr1ZT8/f2nT5+elZXl6OhIhxo1apSdnR2ACxcuUBSlr6/fMpRYLKZD0UVr3rx5HfnG6+vr79y5s27duuY25UKhcNKkSS0f8Oro6EycOLEj0QiCFaQQEgQz3N3dg4ODKysrExMT3d3dzc3NBw8enJiYWFVVVV1d7e7uTh929uzZadOmmZmZOTs76+vrc7lcHo/XPCDF1tbWwcEhNjbWz8+vvr7+8OHD/v7+dIPG0tJSABs3bmwuOTRzc/NW8/Hx8enXr19sbKyjo+O//vWv1NTUH374gf6IDrVo0aJXTqFHzTx58gSAsbFxR77roqIimUz2SmdKsVhcUVHRvMlIA2eC6DqkEBIEM9zd3aVSaUpKSmJi4rJly+g9dCHU0tKib8sAfP31146OjklJSXT3YJlMtmXLlpZx/P39lyxZ8vDhw8uXLz979szf35/er6enB+D48eNubm4dyYfP58+YMePAgQObN2+OjY1VU1ObMWNGy1C3bt2iexq/gu5rXVpaqqur2+5VtLS0AJSXl7fcWV5eTl+C9krlJghlQ/6BEgQzrK2tjY2Nf/rpp8LCQvr+z93d/eLFi7///vuECROae74XFhba29vTVRBAYmJi81hT2qxZs/r06bNv377Y2FhbW9vmvudOTk59+vQ5cuRIx1Py9/evrKw8derUnj176BtEer+LiwuA14WaMGECgMOHD7f6qba29osXL5o3xWKxqanp2bNnm/fU1dUlJiaOGTOm43kSBMvYfklJEKpj9uzZAExMTOjNJ0+e0DdD3333XfMxb731lpGRUXZ2dn19fXx8/KBBgzQ0ND799NOWcd555x1jY2Mej7dp06aW+5ctW8blcteuXVtYWFhXV3fnzp0ffvhh9+7dbaRka2s7cOBAAGfOnGm538fHR0tLKyoq6tGjRzU1Nbm5ueHh4b/99hv9qa+vr0AgiIqKevz4cWVlZXx8/MmTJ+mPPDw8Bg8efP78+ezs7OLiYoqioqKiACxfvrykpOT+/fu+vr4cDufixYv08QEBAS3H7BCEEiKFkCAYExMTA6Dl/AcHBwcA169fb95TUFBgbW1N/xmqq6u7Z88eU1PTVwohvXSLmpra48ePW+6np/e1fGJpZmZ2+PDhNlKKiIgAIBKJXr582XJ/XV3dggUL6LePtCFDhly6dIn+tLa2du7cuWpq/3l1IhAItm/fTn+Uk5Mzfvx4+ono8uXLKYqSyWRhYWECgYA+2MDAYO/evc0XIoWQUH6kMS9BdLempqZ79+7V1dVZW1s314+Oe/nyZX5+fkNDg1gsHjBggCKZ1NXVFRQUUBRlbGzcv3//Vz6tqqoqKCjQ1NQcOHCgtrZ226Fqa2vz8vL4fL6NjU2fPn0UyYoguhkphARBEESvRgbLEARBEL0aKYQEQRBEr0YKIUEQBNGrkUJIEARB9GqkEBIEQRC9GimEBEEQRK/2/69UzBcxbZvXAAAAAElFTkSuQmCC", "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" ], "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" ] }, "metadata": {}, "execution_count": 7 } ], "cell_type": "code", "source": [ "plot_bandstructure(scfres; kline_density=10)" ], "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{StaticArraysCore.SVector{3, Float64}}:\n [-5.696687678871305e-16, -2.6466373207985906e-16, -7.992410752624469e-17]\n [-1.6456025997912e-16, -3.959837374083753e-17, -6.543052491778762e-16]" }, "metadata": {}, "execution_count": 8 } ], "cell_type": "code", "source": [ "compute_forces_cart(scfres)" ], "metadata": {}, "execution_count": 8 }, { "cell_type": "markdown", "source": [ "As expected, they are numerically zero in this highly symmetric configuration." ], "metadata": {} } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.8.3" }, "kernelspec": { "name": "julia-1.8", "display_name": "Julia 1.8.3", "language": "julia" } }, "nbformat": 4 }