{ "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.900435944863 -0.70 4.6\n", " 2 -7.905007303589 -2.34 -1.52 1.0\n", " 3 -7.905210899359 -3.69 -2.52 3.1\n", " 4 -7.905211520829 -6.21 -3.37 2.5\n", " 5 -7.905211530867 -8.00 -4.81 2.0\n", " 6 -7.905211531389 -9.28 -5.28 3.4\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.1020962 \n AtomicLocal -2.1987824\n AtomicNonlocal 1.7296077 \n Ewald -8.3979253\n PspCorrection -0.2946254\n Hartree 0.5530387 \n Xc -2.3986211\n\n total -7.905211531389" }, "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.00482543 0.0611647 -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.436926 0.285692 0.417554\n 0.354532 0.389829 0.384601 0.449283 0.627532 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///8vSugAAHI5JREFU/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+41JTD9k0leXt6VK1eCgoJavJ+YmNi01wKJUTd78lgA8/65v2YKKxrhkxwWj5JBsvL/7d19TBR3GgfwH+4KomTFXSj7wssSr+WqtYpSrHqla4+1gCjBhWLDSXkRUFFaq43x7HnYalII2tq0akkv9kULqUWlUIEFpYBaFF9KrWhAgYC7BSkvAgss7M7cH5PbM4Cv0Z1h5/v5azCP5psnyzzu/H4z8++5gh2XKdP/TmzbL1JbZk9wxOogS55W47dv367RaBoaGq5duyYUClesWEEIaWhoePHFF7u6upydnePj4/fv36/RaGQy2XfffVdWVkYIuXTpUmJi4oIFC2iazs3NjYmJmTdvHiEkPj6+r6/P09Pz2rVr58+fLyoqekqxwTq2zpkw56hp4yyB6ySy+4p5hXKCNzaLAp8sltnJJ5OcBuoff5nwexd97jZ9WGWzV0S47/8LdU9cfX19cXGxi4tLWFgYs12lr6+vvLz89ddfFwqFhJCenp5jx4719/eHhoZ6eHgQQoaGhk6fPn39+nWhUDh37lw/Pz/mn2ppaamsrGxvb5fJZGq12rIEOFpqaurDXxrt7e0d805/uI8n1bSUs+YptHG1Z8/fyiXnwoS2PQjxSXsMNt+00610fIWpdGF76u+Sv3tM3DDzCVwUsfmmPSVPcRCyAoPwaXtSTfvPkYKkTdscxDKRg6C58phtb3XDJ+0x2HzT2tvbvf62zM7JxdihO3vsK3/f2Q/+Ow9i8017SrAwA+zYuXMn9d6pgfX5PYp5R++9+wnAVqXv3Teo2tC/9pg5/qv3/rWT7Ti8hkEI7DCbzURoTwih7CcPDtrizYMA9zVoHCL2UwghxGHKoE3ePzt+YJcSsOOfG1O27w2mpX916bweocHuJ+CdzSmJx4NXGG8UT2iqTj+wh+04vIZBCOxYkxC7PHiJXq+fPXs27nwCHlIqlXUXT1+5cuXZZ3eJxWK24/AaBiGwRi6Xy+VytlMAsGby5Mnz589nOwXwe41Qq9Xefas+PAzcxPkY0LRHZTabRz+pH+5vcHCQuSEbHhWvB2FKSsrt27fZTjHOrFq1yoiF/UdhMBhiY2PZTjHO/PHHH6mpqWynGGfq6uq2bdvGdopxideDEAAAAIMQAAB4DYMQAAB4zdZ2jd64cUOr1Vrea3F/3d3dUVFRtv1wrydueHg4JCQE7zx7eBRFGY1GtVrNdpDxxGg0dnR0oGmPpK+vr6mpCU0bITw8fN26dfevsbVnjVZXV9+6deshn7bX2Njo7e394Dq4C5r2GNC0R0XTdFNTE5r2SCiKamlpYd4CCxbe3t7Tp0+/f42tDUIAAIBHggtcAADAaxiEAADAaxiEAADAaxiEAADAa4K0tDS2M1hDZ2fn8ePH6+rqvLy8xrxfwmQyabXaM2fOiMXiqVOnWj8hBxkMhvz8/JqaGoVC4ejoOLrgxo0bJSUl169fnzZtGt6LzWhqasrLy2tra/P29r7PTSa1tbW1tbVKpdKK0biroaEhLy+vvb3d29vbzs5uzJpLly4VFRXpdDo3N7dJkyZZOSEH1dfX5+fnd3Z2KpXKMZum0+kKCwtra2vFYjF+PR+A5oGGhgY3N7eoqKiQkBAfH5/Ozs4RBWazecmSJS+99FJCQoJYLD516hQrOTmlq6vLx8cnKCho5cqVUqn05s2bIwo+//xzmUwWGRm5YsUKkUiUn5/PSk5O0Wq1YrE4ISHBz88vODiYoqgxy/R6vaurq0QisXI8bvrpp58kEsnq1at9fX3DwsJGF1AUlZSU5OHhsWrVqmXLlmVmZlo/JNfk5uZKJJLExMRZs2a9+eabowuKioqcnZ0TEhLi4uKcnZ1xTrs/XgzCtWvXJiUlMcdBQUHp6ekjCgoLC5VKZX9/P03T+/btW7RokbUjck9GRoZarWZO5WvXrk1OTh5R0NzcbDQamePMzExfX19rR+Se+fPnHzhwgKbp/v5+T0/PkpKSMcvCw8M3b96MQciYM2fOwYMHaZru7e2VyWQVFRUjCg4dOjR9+vSuri4WwnESRVE+Pj45OTk0TXd3d0skkgsXLoyoCQ0NTUtLY463bt2q0WisnXJc4cUaYX5+fkREBHOs0WgKCgpGFBQUFCxdupS5+hcREXHmzJnOzk5rp+SYgoICjUbDXHKJiIgY3TQPDw/LRWaZTIZXUty+ffvcuXMajYYQ4ujoGBISMrpphJDs7GwHB4fly5dbPSAXtbS01NTUME1zcnIKCgoa3bTs7Ow1a9Z0dHScOnWqo6ODjZjcUl9f39DQEBYWRgiZOnWqWq0e3TSJRGIwGJjj/v5+FxcXa6ccV2ztEWujURTV2tqqUCiYHxUKhU6nG1Gj0+n8/f2ZY1dXV3t7e51Ox/N3Rut0urub1traajabBQLB6MrBwcGMjIzVq1dbNyDn6PV6BwcHyxlHoVDU1NSMqPnzzz8/+OCDn3/+ua6uzuoBuUiv14tEIssKlkKhaGxsHFFz8+ZNg8Fw5MgRd3f38vLynJycwMBAqyflEL1e7+LiYlkoHfOclp6eHh0dHRwcbDab7ezsDh8+bPWY44ntfyOkKIqiKMtiskAgMJlMI2rMZvPd+xrGrOGbu3siEAhomjabzWOWxcTEKJVKvD3uYT5F69ev37Jli5ubm3WjcRdzmrb8OGbTBgcHBQJBVVVVbm7uhx9++Pbbb1s3I+c8TNNOnDih0+neeOONlStXNjY2FhcXWzfjOGP73wiFQqGrq2t7e/vzzz9PCGlra5PL5SNqZDKZ5Q29PT09AwMDo2v45u6etLW1ubi4jN5tS1FUXFzcnTt3fvzxxzG/LPKKVCodGBjo6+tzcnIihLS1tclksrsLmpub8/LyRCLRL7/80traajAYkpOT09LSRpTxilQq7enpMRqNDg4OZKymEULkcnlAQABz6lepVOvXrzeZTEKh7Z+77kUqlXZ2dlqu0LS1tVku3li8//77WVlZS5cuJYRMmzZt8+bN0dHRLGQdJ2z/GyEhZPHixVqtljnWarUqlYo57ujoYL7lqFQqZl8DIaS4uHjGjBn4P7tKpbI0rbi42NK07u7uoaEhQghN0+vWrWtqajp69ChzFuM5uVz+3HPPMU2jKKq0tHTx4sWEEJPJxKxsicXib775Rq1WBwYGzps3b+LEiYGBgVOmTGE5N6uUSqWXl1dJSQkhxGw2nzx50tI0yzr9a6+9Vl9fzxzX19dLpVI+T0FCiI+Pj1gsLi8vJ4QMDw+XlZUxTRseHu7q6mJqBAIB83tKCDEajfh/6gOwvFnHKi5cuCASidLS0t59912JRNLc3EzTNPMpqa6upmnaaDTOmDEjKioqMzPTzc3t8OHDbEdmX3Nzs0Qi2bhxY1pamkgkYhpF07S7u/v3339P0/Qnn3xiZ2cXHR2dlJSUlJSUkpLCal5O+Prrr6VSaWZmZmRk5AsvvDA0NETTdFVVFSFkxK0UFRUV2DXKyMrKUigUu3fvDg8P9/X1NZlMNE2Xl5dPnDiRKdDr9XK5fNOmTR9//LG7u/sXX3zBal5O2Lt3r5eX1549e0JDQxcsWMB8ugoLC0UiEVOwa9cuhUKRmZmZkZEhlUp3797Nal6u48vbJ65evXrkyBF7e/vo6GjmNSU0TX/55Zfh4eHM7obu7u6DBw92dHQsWbIkICCA7byc0NzcfOjQoaGhocjIyJkzZzJ/mJ2d/fLLL3t7e1dXV1++fNlSLBQK4+PjWUrKIWVlZaWlpc8880xsbCzzZIb29vbjx48nJibeXdba2qrVamNiYliKyS2lpaVlZWVubm5xcXHMxpnW1tYTJ05YPlF6vf7bb781Go2BgYELFy5kNSxXFBUVVVRUKBSK2NhY5rrCrVu3Tp48+dZbbzEFZWVllZWVdnZ2KpXqlVdeYTUs1/FlEAIAAIyJF2uEAAAA94JBCAAAvIZBCAAAvIZBCAAAvIZBCAAAvIZBCAAAvIZBCAAAvIZBCAAAvIZBCAAAvIZBCAAAvIZBCGA7kpOTPT09m5qamB+HhoYCAgL8/f0HBgZYzQXAaXjWKIDt6O3t9fPzc3Z2rqystLe337hx4/79+8+ePTt37ly2owFwFwYhgE25ePHiokWLUlNTX3311WXLlu3du3fDhg1shwLgNAxCAFvz6aefvvPOO05OToGBgbm5ucy73QHgXjAIAWzNnTt3PD09e3p6fvvtt1mzZrEdB4DrsFkGwNasWbNmwoQJHh4eKSkpJpOJ7TgAXIdBCGBTsrKycnJy9u3b98MPP1RVVe3YsYPtRABch0ujALbj6tWr/v7+cXFxn332GSHko48+2rZtW1FRkVqtZjsaAHdhEALYCIPB4O/vLxAIzp075+joSAihaXr58uXnz5//9ddfZTIZ2wEBOAqDEAAAeA1rhAAAwGsYhAAAwGsYhAAAwGsYhAAAwGsYhAAAwGsYhAAAwGsYhAAAwGv/BY5B5RxjyIdSAAAAAElFTkSuQmCC", "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: 9%|█▌ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 100%|████████████████| Time: 0:00:01\u001b[K\n" ] }, { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=126}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOydeVxM6x/HPzPtK+2UUChCaSNEsl1bCGVLSmS7slwke9asV5SEyE6ua8nPLlxrCCGRLTfiCqm0z8zz+yNbaZnqLDPTeb/u676mxznP9zNnzpnvPM/zfb5fHiEEHBwcHBwcNRU+2wI4ODg4ODjYhHOEHBwcHBw1Gs4RcnBwcHDUaDhHyMHBwcFRo+EcIQcHBwdHjYZzhBwcHBwcNRrOEXJwcHBw1Gg4R8jBwcHBUaPhHCEHBwcHR42Gc4QcHBwcHDUa2XSEERERDx8+ZNjo7t0YO5Zhm0wjEonYysknFAqLXkRG4t49ViQU4+RJ/O9/xVry8jBjBktqvrFiBVJTKejn+9WuMgIBpk+nQInkU/1rVcSzZwgJoaQn2qHqLTODOF9ZsukIz5w5w6QjFAgEc+fOnThx8bNnLxkzygr5+flsPQM5OTlFL9asgUDAioRiaGhg8uRiSpSUcPIkzpxhTxOQkoKICAr6+X61q0xsLC5coECJ5FP9a1XEoUN49oySnmiHqrfMACKRSCDG94VsOkKGady43bJl7zMzdS5e7PT06VO25cgyaWl4/RrW1mzrABwdUb8+du/+0cLjwd8fy5ezpwnw8sL27ZCERPpnz6JbN7ZFSBWXLsHJiW0RNRXOEVJASkoaIZuBCYT4hoWFsS1Hlrl4ER06QE6ObR0AgMWLsWRJsUHhsGF4/RpXr7Imyc4OGhq4fJk1Ad/hHGGlEAhw7Ro6dGBbR02Fc4SUoAj8CwiByy1atGBbjCxz4QKcndkW8Y327VGvHvbu/dEiJ4dp0xAUxJ4mwNMTkZFsCgCQlYUHD9CuHcsypIg7d9CgAXR12dZRU+EcYXXZuxdAKNAFMJGTe9a6dWu2FckyEuUIASxYgMWLiw0KR43CnTu4c4c1SSNG4PBhZGWxJgBATAzatoWKCpsapItLl9CpE9sixCM6OrpLly6zZ89mWwiVcI6wWmzejBkzoKqaBDwF/pWXH5eQkMC2KJnl7Vu8fw9LS7Z1/ISzM4yMsH//jxYlJUyZglWrWJOkr4+OHXHoEGsCwM2LVh5pWSAMDg7u23f8rVs9ly8/Z2Jiw7YcyuAcYdUJDcXy5bh0CQoKnXi8Y8D/1NT2Ozo6sq1LZrl4EU5O4EvYPTt/PhYtKjYonDABMTFISmJNkpcXy7OjnCOsFEIhrl2DVHxzzJ8fDPwJBAIXk5Pfsy2HMiTsS0V6WLECGzbg8mXIy4OQZpqaKXZ213r12mlkZMS2NJlF0uZFi+jcGXXrIirqR4uaGsaPZ3NQaG39+vp1N2NjhwULVjJv/d9/kZ6Oli2Ztyyt3L0LIyPo67OtQww0NJSBop94yewqoRbOEVaFBQuwYwdiYlCvHiIj0aABb8qUiX//vfTUKYu8PLbFyS4xMZLoCAHMm4eFC/HzBks/P/z9N/79lx09w4dPLCjwe/36cnBwwqlTp5g0HRX1t5vbpIYNtwMiJu1KNdKyQPj2LdTULgJhQFNg+Lx5o9hWRBmcI6wchGDqVBw/jn/+gaEhCMGuXXj1Ct7eMDaGtTUOH2Zboozy+jUvMxPNm7OtozS6dkWdOsUGhdra8PbGunXs6ElJSQU6AAqZmd3j4x8zZjcq6vDYsXtu3hz24EHckiV/MmZX2pGKBcLERLRtCz5fPzDw9e3bcTze3enTF7EtijI4R1gJRCL4+uLqVZw9+zXQ+dIlFBSgXTs0aAAAY8ZgyxZ2Ncome/fuHTRoat26B3k8tqWUwdy5WLwYop9GQdOnY8eOT/v3n0xOTmZYTI8endTU/gCOqqmt79v3N8bsHj58/vPnmUDbvLyFR46cY8yuVCMU4soVdOzIto5yiY1F587o1g2Kipg1C2ZmInV1hIayLYs6OEcoLkIhRo3C06c4fx7a2l8bd+yAoiLGjPn6Z9++ePQIXG4ZatmwYaOHR9CjR50fPNg1cKCEzsZ07w4dHRw8+KPly5dneXldvLyu2dt7/PXXUSbFhIUFhYbajBhxT0Vla8OGzRiz26GDtaLiXuCTklJku3ayE1JIK/HxMDSU6AXC6Gi4uODPP/G//2HbNigqAkDbtjhwgG1lFEJkEXd393379lHYYX4+GTiQ9OxJcnJ+NH75QjQ1iYEBKSj40Th9Opk1i0LLkkVOTk5hYSHDRps06QBcBAjwSVm5EcPWxefkSWJhQYTCr39OnTofOAYQIK1ly66sSHJxIRs2VOXEzMzMKp0lVFcPMjPrNnFiQM7Pj4pMU7Vr9Z01a8iECVRpoZ7ISGJgQK5dIwMHkjlzvjZmZmYeOkTk5VlVJh5CobDg5y/oMuBGhBWTnw83NxQU4PDhYnuEo6KgpwcfHygo/GgcMwaRkSgsZF6mzNKsmTFwFgAQo6urxbKasunRA1paPzbwaWtrKii8AQC80dKqxYqkhQuxfDlycxkyt2kTv3dv/ydPzoSELFPhttOLhyQvEK5YgcBA/PMPUlKQkIC5c3/804ABAHCU0ZkOGuEcYQVkZ6NPH6ip4e+/oaRU7J8iI/HhA7y9izWamaFJk5IFejiqw759W/T1z/F49bW1F5w/v4dtOeUxezYCA7+uFE6ZMtbK6oiubjt5ed+xYxezosfGBjY22LaNCVvZ2Vi7FnPmMGFLZhCJcPWqJDpCoRATJmDPHly+DB0dTJmCiAgoKxc7xsICMpNZmXOE5ZGRge7dYWKC3bshL1/sn5KTcfcu7OzQuHHJs0aPxtatjGmUfVRVVceOveHv/+zjx4dmZmZsyymPXr2gpvY1clhdXf3WrVOvX1/YtSs2KKgZW5MEgYFYtgwM7OoJDYWzM7d9sHLEx0NfHwYGbOsoTn4+hg3D48e4cgVGRvj9d3h4lJI51s0N16+zoY8GOEdYJunp6N4drVohPLyUbCaRkahVq/RKvG5uiI1lbQ+ZTHLzJmxspKMW6Pz5WLDgR/iokpLSkCGoXx/BwezosbGBtTXtg8LsbPz5J2Qr/SS9CAQCD49Jzs526en9UikppkwRnz+je3cIBDhxApqaiI7GnTsIDCzlSD8/ZGWxmUGJQjhHWDr//QcnJzg7IzQUv4bsE4Jt25CTg759SzlXRQWDB7NfAUBmIAS3b0uNI+zdG6qqJddOQkKwciVeslS2OTAQQUHIz6fRxIYN6NwZXOUV8YmI2Hn4sGZGxu3//ps0erQ/23K+8vYtnJ1hZYWDB6GsjIwMTJyIrVtLz5+uqYm6dbF2LeMqaYBzhMUoLCwsLCxMSUGHDhg4sMx6OhcuIC8Po0eXXDX8jq8vtm4FS7XcZY1nz6ChAX19Cag2Kx7z5mH+/GJ7Chs2xOTJmDiRHT22tmjZksZBYXY2goMxbx5d/cskiYnJOTntARDS/vnzZLblAMDjx2jXDv36Yf36r3Ngfn4YMKC8Kom//YYTJxgTSCOcI/zB0qXBhob2deq0trJaNX48Fiwo88jISBQUwMenzAMsLWFggHPclmIquHkT0lXbysUFKio4dqxY44wZePWKtcRDixZh+XK6BoXr16NrVzRtSkvnskr//v3l5AKBXbVqeY0a5c62HNy8CWdnzJ2LhQu/tpw8icuXsWRJeWdNn47Xr/HlC/366IaBnRzMU4V9hOnp6bq6bQAhIFJX7/ju3buyjszKImpqpF27CjrctIkMGlQpCVIAK/sIJ08mq1ZVd7cWwxw5Qlq1IiJRscZ//iGGhuTzZ3Yk9exJNm0S92Dxr3ZWFqlThyQmVlGVDFCFO1MoJD16EB+fhBUr/jx37hwdqirFmTPEwID8738/WjIySP365MyZ0o//+S2rqZHVq2nWVw24fYSVIycnh8/XAfgAT1VV90vZP3KioqCmVvE01/DhOH8e795RrLMGInUjQgB9+0JeHtHRxRo7dMBvv/34xc0wixdj2TIUFFDc7fr16NaNGw5WjvnzkZ+PTZssZs6c0qVLF3bF7NyJESPw99/o1etH4x9/oHdvsWpptWlTVJxcuuEc4VcMDQ1r1aqtqOitozOmRQu+qalpWUdu3oy8vK/7SctBXR2urti1i2KdNY3CQjx4ABtpS9fF42H2bCxcCFJ8ZXP1ahw4wE79eltbNGtGcQzXly/YsIHbO1g5oqOxaxcOHCi5I4thCCEAgoOxcCEuXSq2OyImBufOYcUKsfrx9cWDB/RIZBDOEX4lMREZGbujosZGR486dy6KV0Z255cv8fAhvL1L7i0tlaINhURqgjwkkfv3YWoKdXW2dVSeXr0KXrzw0dGxtbfv9eZNUYoZaGtj2TKMHctOINXChVi6lMpB4bp16N4d5uaUdSjzPHsGHx/s3Qs9PdY03Llzt2HD1gYG9g0bDtmypfDy5WKfYHY2fH0RFgYNDbF6GzwYhEh9yAznCAFAKMTIkVi2jNevn0Pbtm3L8oIAIiLA4/3Isl0+bdtCRQWXL1OmswYijfOiRWzaFJGba5qeHhcXN93XN+B7+8iR0NRkJyWHgwOaNsWOHdT0lpmJ9eu5vYOVIDcX7u5YvBjt27Mpw8trxqtXB9LSbr9+3ez33/eXKCU+cyY6dUKPHpXo0MxM6lPMcI4QAFavhoYGRlVU2IAQbN2KJk0qURXPy4vLMlMtbt2CvT3bIqrE06evCwrsABBin5yc8r2dx0NYGBYtwrdRIqMsXIglS6gZFAYHo3dvbjhYCcaPR/PmpWfhYJIvX7IBQwAiUcOMjI8//9O1azh2DKtWVa5DNzdcuUKhQBbgHCGePMHq1diypZSN8yWIiUFuLqZMqUTnI0YgOhrp6dURWKOR3hHhqFGDdHXn8nhbebzhbm4eP/+TmRnGjsXUqSyoatsW5ubYubO6/WRmYsMGBARUfCRHEcHBuHcP4eFs6wD69fPm8/uqqQXWrbvOw8Pte3tODry8sGEDtCqZ2d7PDxkZ0l1+rqY7QpEIo0dj0SKUHRzzg7AwCAQYNKgS/evooGdP7JHoTNGSS1YWkpMltCp9hdjYWF+7tjckRDB+/Kzz531+3l8PYO5c3L+P48dZEBYYSMGgcN069OkDyc78KkFcv47ly/H331BVZVlJTg5OnPBdvXrtnj2tEhIuGP00MTp3Lhwc0L9/pfvU1oa+Ptavp1In09C/kYMFxN9HuGYN6dSp5H6vUsnKIioqxMen0mLOnyeWlpU+SzJheB9hTAxxdPz6Wrr2Ef6MUEg6diTBwSXbz5whDRqQL19YkNS1K9m6tbwDyr/anz8TPT2SlESxKimlwjvz3TtibExOnGBGTgVMnEg8PUtpv36d1KlD3r8Xq5Nf37KHB2nYsNriaIDbR1gxL19ixQqxJkUB7N8PHg+//15pK87OyM3FrVtVEFjTkd550Z/h8xEZicWLkZBQrL1bN7RrV0HmDpoIDMTixVUfFP75J1xc0KQJpZpkFIEA7u4YMwY9e7ItBTh/HkePYt26ku35+fDxQUhI1WNZp03Dq1dMFDmhiZrrCEUieHtj9uxS6iiVyrp1qFcPrVpV2hCPB29vLmSmKkhvpEwJTEwQGIiRI0tWbP7zT0REID6eaT3t2qFx4yrO2GdkICyM2zsoLtOnQ1NTIi5XRgZ8fLBlSylLgAsWoEULDBxY9c6traGigs2bqyOQTWquI9y4EQUF4o7wnj3Ds2f4448q2vL2xsGDyMqq4uk1FtkYERYxfjz09UvG4xkYIDAQY8eixAoiAyxahKVLIRBU+sS1a9G3r1hr6hz79+PoUURGllLHjXkmTYKLSyn7Iu7exY4dFKzw2dtj9+7qdsIWEvD5sMGrVwgMxNatkJMT6/iiYkxDh1bRXJ06cHJCVFQVT6+ZvHuH7GyYmLCtgyJ4PGzZguBgxMUVax87Fnw+Q0Xkf6ZdO9SvX+lvrqLhIBcsKg6PH2PyZBw6BB0dtqUAR4/i2jUsX16yvaAAI0fizz8pKA7s48PC3AZV1ERHSAh8feHvDwsLsY4XiRAZiT59xE21UCpc2frKcvMm2rQRa/lWWjAywqpV8PQstpTC5yM8HHPm4P17pvUsXowlSyo3KFyzBv37c8PBisnKwoABWLFCIrIDpqVh/Hhs315KhqYlS9CwIYYMocDK8OEQCnH2LAVdMU9NdIRbtuDz50rs4jp/Hrm5mDWrWkZ79kRqqhT/YmKeW7dkZ170O56esLDAokXFGlu2xIgRmD6daTHt28PYuBIrhZ8/Y9MmSRwOCgQCL68JjRs7LliwqOKj6YcQeHujSxd4ebEtBQAwfjy8vEopKxgfj02bKEsKw+ejcWOEhlLTG8PUOEeYmor58xERIe6kKICgIOjqwta2Wnb5fIwcie3bq9VJjeLmTRmJlClBWBgiI0sm3lu0CFeu4Px5psUsXozFi8UdFK5ZA1dXSZys7tFj6M6duc+fL1u8+MKcOYFsy8Hy5XjzBmvWsK0DALBjBx4/xvz5xRpPnDgRFLTOw+PZ6tUokWKtOgwYILUZJRnYycE85ewj7NmTLF1aia4yMoiCAlm1igJV//5L9PRIbi4FXbEFY/sIRSKio0P+++9Hi/TuI/yVI0dIo0YkK6tY44kTpEkTFm4PJyeyY0fJxl+v9sePREeHvHjBkKpKoajYHEgDCHBdVbXfzp2MXsYS1+rcOWJoSF6/Zk5AObx+TQwMSHx8scYZMxZpanoBOxQVbZ48eVKFbst6GNPSCI9HkpOr0CVdcPsISyEyEqmpmDGjEqcUZdkePZoC68bGsLHB339T0JXM8/QpNDWhr8+2Dnro1w9t25ZMV92zJ5o3F7f2DYUsWCDWoHDNGri5SdZw8NMnjB2LWrVQUNAV2A9kA3vl5e18faGqCi0tdOiAwEC8esWcpJQUjBiBvXupHGZVGUIwejSmTIGlZbH2qKjjmZkRgGdh4YyDB6PLOLsq6OpCV7eUfYqSTw1yhG/fwt8f27ZBQaESZ61fj44dUbs2NRq4kBkxkaWNE6USEoKjR3H6dLHG0FCEhuLxY0aVODvDyAj79pV3zMePCA+Hvz9TmirixAm0bQs9Pfz1F3x8kJw8u3HjfYqKzR0cXnz8OCs3F2/fYulSaGggLAwmJlBURKNGGDECR4/SuFMlPx8DB2L6dDg50WWiUoSEICOjlN/9amoGQBxA1NWvNG1K8U+bLl1w7Bi1XTICA4NT5il1atTVlSxYULl+kpKIvDy5fJkqXSQ/nxgYSHFuKsamRv38yOrVxVpkaWq0iLNnibEx+fSpWOPatcTJSaycfxRy/jxp3Jj8/MGWuNoBAWT8eEYllUp6Opk6lWhrEzk54uAg7oOZn0+OHCEeHsTUlCgoED6fGBgQJyeyYgVJTyfHjx9XUDDm8020tMwzMjKqoOr7tRo9mvTvz/RnVxbPnhF9ffLrxOfJk0RL64WZWU8jI9tx4/xFVZJbzsN44wbh8Uh+fhV6pQUxp0ZriiPcs4e0bFnpj2f4cKKrS6UwQsj06WTWLIr7ZAzGHKGDA/nnn2ItsucICSHjxxMvr2ItAgGxsSG7djGtpGNHsnv3jz9/vtofPhAdHZYXfqKiSIsWhMcjBgZk/vxqfc/euEF+/51YWRFVVQIQHq8r8AAgwPL+/UdUocOia7VzJzE3J1XypNQjFBJHRxISUrL9+HFiYECuXatu/+U/jCoqJCysuiaognOEPxxhWhqpW5fcvFm5ToRCoqpKvdN68oTUqUPE+GgkEWYcYUEBUVcvmY1aJh3hly+kSRNy6FCxxlu3iL4+SUtjVMm5c6RJkx+Dwp+vtr8/mTiRUTHfefeOjBxJ1NSIvDxxdia3b1Pfv5ycOZAFEOAQEKikRAwNiYMDGT2abNsmlvvPzMy8e5fo65OEBIrlVZmlS0mXLiXHplFRRF+fxMZS0H/5D6OjI3FwoMAKJXCO8IcjHDSIBARUupODB4mcHPn4kUphRXToQA4fpr5bBmDGEd66RaysSjbKpCMkhFy9SurUKRYfSwiZOJGMGcO0kg4dyJ49X19/v9ofPhBdXZKSQrv1wsLCQYO89PRs3N29hULhsWPEwYHw+URbm0ydSmMUqIfHOB7PnsdbyOcbnzx54eJFMncu6dmTNGtGNDUJj0d4PKKqSkxNSc+eZOpUcuTIj/FoUFAQYAQY83iT9u8X0CWxkty7R/T0yKtXxRr37iWGhiXDR6tM+Q/jli1EWZkaQ9WHc4RfHeGRI8TcvCoPko0Nad2aYmFF7NxJevWipWe6YcYRhoaW4gZk1RESQmbMIP37F2vJyCB1674ZPnzR/PlBH+n4LVYaZ86QJk2IQEDIT1d75kzy++9MWO/b14PH8wYeAN7y8mvl5UnXruTOHSZMHzlyZMaMGQllDOhu3yarVpHBg4mVFdHRIfLyBCAqKqR+fQLUA+4BIsDH0bEDE1orIi+PWFqW3A+zZQsxNCQPH1JmpfyHsbCQ8PkkJoYyc9WBc4T7CCEfPpC6dasS7ZKRQeTkSHQ09doIITk5REen5E82qYAZR+jlRTZvLtkow44wL4+0bPljNEYIyc3N1dW1BfbLy29r2rQ9Y0o6dCB79xLy7WqnpTE0HCSEaGlZAkkAAZ6oqLRnsOplpfnvP7JnD/n9dwIYASKAANs1NduwUlqyBLNmEVfXYi1hYaRBA/L0KZVWKnwYmzQhgwZRabHKcPsIAcDPD8OGwdGx0icuXQpVVfTuTYMmQEUFQ4YgMpKWzmUAmd87UQIlJezYgSlTkJLytSUxMZHHswYGCwTe6ek6qampzCiZOxeLF//YYLByJYYORb16TJg2NbUCNgOvgRAnJyN5eSaMVg19fQwbhg0boKqqCAwHtgMLTE0D6teHpyfOnWNN2PXriIwsljJt9WqsWIGYGHGLzVFF//64cIFRi9WFAZ/MPEUjwuhoYmpaxQrgenpk9GiqZf1EfDwxNv46DSVFMDAizMwkGhrkVyMyPCIsIjCQdOv2NcDh06dP+vp2wDvgpZpaK2YidYtwdCT79pHMzMyi4SAzGVLu3CHy8oI6dXzV1Zt27jwgX3Ki78ulsLDQ2dlZX99k48aNhJC3b0lQEGnUiFhYkKAg8uEDo2Kys4mZGfnrrx8tQUHE3JyWAX2FD+O7d4THI6mp1JuuLDV9ajQi4i9jY3LpUlVOv32b8Pm0f4p2duTUKXpNUA4DjvD8edKhtNUWmXeEhYWkdWsSHv71z1OnzrZo0dnKqkfTprEzZjAn49Qp0qwZ+fw5848/yOTJTFj891+iokJ692bCFh2UuDNFInL5MvH1JVpaxM2NnD3L0M7C8ePJyJE//lywgDRrRt68ocWWOA+jjg5h8r4tixrtCLt2HaKsbOPqWiU3SEj37qRpU2oVlUJ4uKRMo4sPA45w+XLyxx+ltMu8IySEJCYSXd2Syznv35OWLSudC6LKZGdnKypa8ngmPJ7TX39dodvcx49EU5O0aUO3HRop685MTyfh4cTSkpiZkaCgkoHB1PJzcgaRiEybRqytadyBI87DOGgQadKELgHiU6MdoZ1dOPCfkZFtFc4VComCAtm6lXJRJcnKIlpa5O1b2g1RCAOOcMAAsn9/Ke01wRESQtasIe3bl5wzf/+etGhBAgOZEODu7gksAAhwo3ZtC1ptZWcTPT1iYkKEQlrt0EuFd+bt28TXl2hrEzc3cuwY9Qsinz+T+vW/Ti+JRMTPj9ja0jsxK87DePUq4fNLWeNgmBodLJObqwPoFxbyRJVPLBgSAj4fo0bRoasY6upwdcWuXbQbki5qWqRMCaZMgbw8goOLNerp4fx5REVhyRLaBbx8+R9QVEy2RU5ODn2GRCJYWYHPx6NH4Mvm99BXbG0RHo6XL9G1KxYuRMOGmDULL18KPT396tWzc3Ts//bt2+r0//vv6NcPv/0GkQijR+POHcTEQEeHKvlVpF07KChIT0ggAz6Zebp0maCqOtrdfWylzsrOzp41a1bt2lv69MmmSVgJrl0jZmaSkpxQHOgeEb59W2ZOuxoyIiSEvHhBdHXJgwcl2//7j1hYkCVL6LV+7NgxPt8EWM/jOQ0a5E2fIVtboqHBdAIdOqjsnXnrFhk7lqiqRsjJBQCExzvTs2dVUrsVcegQMTMj2dlEICAjRxInJ8LAgyLmW3ZwII6OdGupADFHhOwEKefl5V2+fJkQ0rFjR2Vl5VKPeffu3Zs3b6ysrOR/iqTOzMy8fPmyhoZG+/bt5couraup+dbPz3zZsmWVUqWv3yI72xVIPneupUDwRJ7+CG49vWevXo3S1s62t29y7FhkWZei5hAbW6OHg0WYmGDpUowciVOnPisry2loaBS16+sjJgadO4PHK1nCiUJcXFyuXtXbuHFjly7eI0eOpMlK165ISMCjR9DVpcmC5GJnBzs7yMklb9zoCIAQx1OnFjVtChMTNGyIhg1/vCinDJlIJLp69eqXLwoTJzocPgwFBQwZgi9fcPIkVFSYey/lM3Ikpk1jW4SY0O+SS/Lhw4dmzZo5OTk5Ozs3bdo07ZffhK9evTIyMiryCv/9tMT8+PHjOnXq9O3b19ra2tnZuZwY63IK85bFjRs3+PxOAAEIn9/79OnTlTq9anTsOBC4DRAlpVVr1vySIlfyoHtEOGdOmVEhNWdEWISx8Qw1tfZ6eq3nzStWFfrdO9KsGVm+nF7rtF7t4cOJggK5e5c+C4xStWsVFxeno9Ma2FWrlvuSJeufPydnz5LwcOLvT9zcSPv2xNSUKCsTU1PStSvx9SVBQSQqity+Td68ISKRqEOHftra4xQVvRs1Gp6fT/r3J717M1eOWMy3nJ9P+Hxy9SrdcspDckeEoaGhpqamx48fB9C/f/+QkJCFCxf+fICuru6pU6f09fUNDAx+bl+0aNGwYcPWrFlTUFBgb29/6NChoUOHUqWqfv36hKQCeQAfeG7CSAXStLQPgBmA/Hzz169vM2BRwrl5E1OmsC1CAnj58mVOzsPs7CvZ2SQszGHGjLHfx4UGBjhzBs7O4PEkqECg+MyejX37cPIkWrViWwqr2Ef7pBwAACAASURBVNjYXLq07fjxM7a2Y7p27QrA1LTkMVlZSE7Gy5df/3/jxtfXItHT3FwlgSAMQEZGj99++09f3yAqqnKVVhlAURENG2LDBrRrx7aUCmHAJ5fA2tp6z7d0Uvv377f6Nb8yIYSQDx8+oPiIUFVV9fa3/POBgYFubm5lmajCiDA1lcjJbeXxGvP5xoMHV25xscqEh+/Q0enO56/W0LBOTExkxmh1oHVEKBIRLa0yo8xr1IgwMTFRT8+taH5CVdUpNbVkutGUFNK4MQkOpksATVd77VrC53/N4iYzMH9nJiS8UVd3BESAQF7eetCgDIaDM8V/y1OmED09WrVUgOSOCFNSUoyNjYte169f//Xr1+KclZ6enpOT8/OJp0uU9/6JjIyM69ev83i8oj+VlJRcXFzK6Vwkgr29XIsWPnFxXkUtQqFQHFXVxMdneOvWlgcOPDpyJLpJkzrMGK0OQqGQz+d/v7DU8uQJtLTkdHRKvwxCoVDyrw9VNG7cWEXlJZ/fC8hVVjbp0EE7PFzUqRP5fkDdujh7Fl27ygGiiRNJOV1VDTqu9l9/8aZP5//5p8jdncjSJ8n8nVm3rgoh74FWgIjHK9ixQ5HHY1SC+G956lQEB8ulpgqLz+4xR1Hl4QoPY8ERFhQUfI9DUVBQyM/PF+esosN+PjEvL6+sg9PT069fv/7dxWpoaHTp0qWc4JqBA5U+f8adOzniaaESMzOzuXPNoqJUYmPzW7Wq9GYPhsnPzxeJRDQ9c9evy9vaoqz7oaCgQMxbRQZ4/PjfzEwtkSgIUODzPQMDP3t71+rWTbR8eaG6+tenWl8f//sfr2dPZaFQMHasgFoBlF/tCxf4w4crT55cOHp0oYx9jMzfmc+ePVNVbZedvRnga2m5paS8+j5CYAbx37KuLmrXVlm3TrRwYSHdqkpFJBKV883/HRYcYd26dYumPQGkpaXVrVtXnLP09PTk5OQ+fPigra1ddKKhoWFZBzds2NDV1XXIkCHi9LxlC86cwaVL0NdXFed4OvD2xp49ypI/k87j8RQUFGiKp42Ph4MDVFVL/xSEQmFZ/yTJ5OTklC+7sBAvXiAp6et/T58iKQkfPxKBQAewBJCfr9Gjh7BPH97MmfI2NvKhoejX7+u55ua4eBHOzopKSooTJlApm9qrffcu+veHhwfWrlUAJGwhq9owf2daWlryeA+BI3JyBZqab83MzGiapCmLSr3ljh1x7JjCypXsfO5i/nBnYSNr+/btY2Jiil7HxMS0b9++6HX5A1g5Obm2bduWemJ1ePgQ48dj8eKqVKigEG9v7N8POrcvSwEytpU+MzPT2rqriUmPBg3sEhISihrT03HlCjZvxqxZcHeHnR00NeHsjPXrkZQEU1NMnYp//sGXL00dHQt0dLy0tAYrKpr5+OgIBAgPx759mDkT7u749ksS9esjJgarVxerOSBRpKSgfXv06CE9e6slHqFQSVn5qJtbwsyZydeuHWPYC1aWyZPx9Ckqn9qEWWhfrPyF+Ph4DQ2NFStWrFq1SlNT8969e0XtampqZ86cKXodEBDg5+cHYNKkSbNmzSpqPHbsmI6OTlhYmL+/v56e3vv378syIWawTH4+0dEhzs7VfktU0Ls32bWLbREVQV+wTH4+UVMrr1SI1AXLLFq0Sl4+HCDAA13dAVZWREWFGBqSTp2Iry9ZtYocPUoePyZlbQISiUQ3b968d+9efj6ZPp0YG3+tdJqTQ/z9Sd26xYqvJieThg3Jpk2UiafqaqelSX0q0Qph/s6cNq1Yfm3mqexbVlQkO3fSpKUCJDdYxtLS8p9//tm5cych5OLFi1ZWVkXtS5YsMTMzK3qtpaVVq1atoKCgn090cXGJioo6fPiwurp6bGysnp5eNZU4OYEQnDpVzW6owccH69fDw4NtHSwRH48mTaCmxrYO6khPzxIImgAA9NXUsrZtQ5Mm+LYJomJ4PJ69vX3R61Wr0K0bRoyApycCAxEUhMGD4eODqChs2oR69dCgAc6c+brX3teXlrdTBXJyYGEBXV1cu8a2FBni/n3s2YP799nWURmsrBARgREj2NZRDvS7ZBYQZ0Q4fz6RlydPnjCjqGIKC4mhIUlKYltHudA3IgwJIb6+5R0gdSPCTZuey8m10tScpavb9vDh6Op3+P49cXEhdnZfb5KCAhIURPT1SXj41yx9SUnE2LjYSLHKVP9qC4WkUSNiYMDcLm+2YPLOFApJ27YkIoIxg6VT2bccHExUVWnSUgE1Oul2hVy7hiVLsGEDvg1B2UdeHh4e2LaNbR0scesWvo1/ZIELFzB/vumFCxf++qvznTsH+vfvU/0+9fRw9CjGjEGHDti8GQoK8PfH+fPYtg1OTkhKQpMmOHcOc+Zg9+7qW6su9vZ4/x4PH6LG5w2kko0bIS8Pb2+2dVSSceOQl4dbt9jWUQ4M+GTmKX9EmJ5OVFWJuzuTisTiyRNSpw4R4+cLa9A3ImzWjMTHl3eAFI0I4+NJnTpVLAotDgkJxNKSDBr0tf6cUEjCw4meHgkKIgIBefyYGBjc0dNrradn7eIyUlClqj9Vu9rp6ekGBlZycvXl5R0VFVNfvKhCH9IHY3fm27dEX58kJDBjrTyq8JYbNiQeHnRoqQBuRFgmbdtCVxcHDrCt4xfMzNCkCU6cYFsH42Rl4fVrWFiwrYMKnj9Hz54ICUHHjnSZsLBAbCwMDWFtjStXwOfD1xexsTh3Dvb2yMmBhsYfaWkH0tLuxMQYRkUdpEvHLwwY4P3ff8OFwlcCwR86OoMYSVNYg/Dzw5gx0vqY9OmDs2fZFlE2Nc4R/v47XryQ3NV7Hx9ERLAtgnFu3oS1Neiv9kE779+jVy8EBmLgQHoNKSsjOBjBwRg0CAsXQiiEiQnOnMGUKejZE+/eZQD1AGRnN3ryJI1eKT+RlJQOFG2GdcjI+FDB0RyV4fRp3L2LOXPY1lFVZs7E+/c/tv1IGjXLEf79N8LCEBUFIyO2pZSBuzuuX0dqKts6mEU2dhBmZqJHD3h5YfRohiz264f4eMTGwtERL1+Cx4OnJ+LiIBCIgL7AMh5v9erVfVesAN35t0QijB2Lt28XAJOBSB5v8JAhv9FrsiaRm4sJExAcLEElliqLsTFq1UJICNs6yqAGOcK3bzF0KMaM+ZGYQwJRUcGgQTVu67EMRMrk5cHFBW3bIiCAUbsGBjhxAkOGoE0b7NsHAEZG0NdXABYB5oDf4MEXz5+HnR1iY+nScPYsdHWxezd27XI+cmSeq+vJzZtHRESsp8tezSMwEG3aoFcvtnVUD0dH7N796v3792wLKQ0GliuZ59dgGaGQ1KtHLC3ZUlQJbt0iJiYSWraepmAZIyPy8mUFx0hysIxAQAYOJEOGEKGQNQ23bpEmTciIEeTLFzJ9emCtWp58/mZNTWsDgxdubiQ4mNSpQ3x9SUaGWL2JebXT04mTE+HziatrmZkBZB6678yHD4muLnnzhlYjlaMKb7mwsFBDwxrowuebjBgxkQ5VpcIFyxTD1RXp6bh8mW0dYmBnh1q1cPEi2zqY4u1bFBSgYUO2dVQVQjBuHD5/RmQk+Ow9T3Z2iIuDigrs7DBs2PxZs5q5uf1z+vTG589NbG2xaBEcHJCfDwsL7NxJjcX162FggMePceMG/v4biorUdMvxM4Rg7FgsWYKyMytLB8HBwV++WAHnRKLEPXuOsS2nJDXCEW7diuPHceIENDXZliIeo0bVoJCZ2FjpXiCcOxfx8ThyBEpKLCvR0EB4eFHi3CWBgY8OHnTs12/8+/cv/f3x8iUcHPC//6FNG6xZg86d8eRJ1Q09fw4LC0ybhqlT8e6d1E9rSzJbtqCwEGPGsK2j2uTn5wPqAAB5gCeSsNyjsu8IHz3CuHFYuJDGcHbKGTECJ04gPZ1tHYwg1QuE69fj0CGcOAF1dbalfGPQINSufSwvL1IkGpuWNnXr1mgAGhrw98fjx2jSBP/+i/x8ODpi4UJUtnxQUVCMmRkUFfHmDYrnQOSgmA8fMH8+wsPZnGmgCj8/PxWV4zyeO9DZ2dmBL2FvSbLUUI5AgI4d0b495s1jW0plqF0bPXtizx62dTCC9IaM7tuH1atx+jR0ddmWUpw6dbSBRwCRl78VGmpsaYmAAFy5gtq1ERSEpCR06AChEPv3o3lznDsnbrfnz0NPD7t3IzIS9+6BrVKrNYcpUzByJFq1YlsHFairq2dkPN29eyCfv2PePInbxC3jjtDJCSKRRG/kLAsfH2zZwrYI+iEEcXFSOSKMicGUKYiORoMGbEv5hb17g1u2nFK3rp2Xl/LHj/23b0ft2pg1CwYGcHfH4cOYOhXx8fjtN7x/j4ED4eqK8kP5cnLQqxe6dUObNkhPl+zsybLCxYu4ckXKfsGXj7y8/LBhg83NG65dy7aUX5BlRxgYiJs3ce2aVC7jOzsjJwdxcWzroJknT6CtLXEjqgq5dQtDhuDQIXwrnSJZmJub379/LjU1bvPmVXJyPFtb+PvjyhU8eoQ+fXDuHJo2hasrtLSwdy8GDcLZszA1xcqVpReNCw2Ftjbu3MH16zhxQiqfJqkjPx/jxmHDBgmacqeKYcMkMRJQNh1henrmkiW7Fy0SbdiApk3ZVlMleDx4ecl+yIw0zos+e4b+/REeznIx5yqgrw9PT0RF4e1bBAUhLw/TpyMmBr17o2lTzJ0LU1PExWHChAl16hjb2dm9fAkLC0yejEmT8O4d2rSp2ERMzIUpU+YfOnSE/ncjyyxdipYt4eLCtg4amDYNWVm4d49tHSVgYCcH8zRtuhrwk5OTjJK7VeXtW6KjQ7Kz2dbxE5TvI/z9d7J2rVhHSsg+wjdviIkJ+3VwKOT5c7JuHenalairk1q1CI83AzABNgE9gWArK5KaKm5X0dEntbS6Ayc0NT3WrNlIp2oJgvI788kToqdHXr+mtlcqqeZbbtCAuQTcYu4j5BFC2PbF1KOikpmXpwk0ys9PVJTmqZy+fTFoEDw92dbxjdzcXAUFBXnqsoK2aYO1a9G+fcVHZmVlaYhf1pYe0tPRsSNGjMDMmewKoYVPn3D6NIYPNyJkE+AC5AItunZ9Tgg+f67g3IICZGfj/ftJX74MBxyANBubkXFxNSJ/PLV3JiHo2hX9+2PSJKq6pJ5qvuVp07BrF9IYSYIrEomEQqGCgkL5h0l/nuPSMDS8/eJFMyD3/n1FOzu21VQDHx+sWSNBjpBaCgqQkCA1QXG5uXBxwW+/yaYXBKCtjaFDMXq0ICfnBuACxAEF/v7g8VC7dgXnKihAXR0REc3Wro3Oy7Pn8Y7k5lrk5XHFCCvNjh3IyMCECWzroBN/f6xbh1evJCjQTDYdoapqBHCpTx9PFxd064Y//4SODtuaqkTv3pgwAY8fS+tKZ/ncuwczM6ipsa1DDIRCeHigYUOsXMm2FJp5/jy+bl1bYCcgCggY2bVrJc4NDPRNT5976lQ7a2sbefmVtrbYtw+WlrRplTk+fcLs2YiOhpwc21LoxMAA+vpYvRobNrAt5RuyGSxjYSHYt291dHRQYiK0tGBlhZ07IY1zwPLyGDFCZnNwS36kzIcPH0aM8Gvd2qVLlwOZmdi2TRa2NpdPnTp1CHmTmfmIkDfLli2r1Lny8vIbNwa9eBF76FDYgQMa/v7o0gXBwVL56LHC9OkYPBi2tmzroJ9evXBEkgKqZPyxrl0bwcE4dgwhIXB2RmIi24Iqz5gx2LEDhYVs66AByc8p4+Y2fu/eDrduhV29unP27BvSvNzMAp6euHkT+/fD1RUfP7KtRuK5fBmnT2PhQrZ1MEJAAN68qXjtmTFk3BEWYWODGzcwbBg6dMCsWcjLY1tQZWjUCObm+N//2NZBA5I/Inz6NFkkcgPqCYVu8fG32ZYjfZiY4PJltGqFli1x+jTbaiSYggKMG4eQENSqxbYURmjSBLVq4c8/2dbxjdId4YcPHw4dOrRgwYKJEydOnjx56dKl586dy8nJYVgchfD58PXFgwdITUWLFjh1im1BlUEmy9ZnZOD1azRrxraOcmnatAWPtwG4pa29vXNn6UlWK0nIy2PhQuzaBR8fTJ4sm3Mb1WflSjRoAFdXtnUwSOfO2L+fbRHfKbGd4syZM3379pX7tlarpqam9C2pvpqa2ujRox88eEDLdg9K+bUe4c+cP0/MzUmfPuTff5kUVXVycoiOjkSopXAf4dmzxMmpEsczv49QKCSOjl+6dl02YIDv+fMxDFtnFzqu9vv3pHdvYm9Pnj6lvG82qf61evqU6OiQFy8okcMElNweN28SPp/k5la/p/KodD3C5OTk7t279+zZMy8vLzQ0ND4+Pi8v78uXL3l5eZmZmVeuXJk/f35cXFyrVq3GjBmTnZ3NgtOmiM6dcfcubG3RqhVWrIBQyLagilBRweDBlJWRkxAkf140KAh8vtrp0wGHDoV37uzMthypR08P0dEYPhzt2mH3brbVSBKTJyMgACYmbOtgFnt7KCtj0ya2dQD4eWo0Nja2UaNGL1++PH369NixYy0tLb+PBTU0NNq3bz9z5sw7d+7Exsa+e/fuzZs3LAmmBhUVLFyI2FjExMDODjdusC2oIopmRyWshle1kPBImbt3ERyMXbtkP0yUSXg8TJ6MCxewciU8PfHlC9uC2Objx4979ohSUuDnx7YUNmjXTlJC4n885YMHDw4LCzM2Ni7/BFtb2+joaDMzM5qFMUHjxjh9GrNnw9UVnp64ezfZx+cPL6+pz58/Z1taSWxsULs2LlxgWwd1SPKIMC8Pnp4IDkb9+mxLkUWaN0dsLLS0YGmJ69fZVsMSGRkZLVt2atZs2MiRdlOn3qko84ls4ueHhw8l4vc993MXbm5ITEStWgJ7+0Hbt3fbubNnp07u+ZWtWEo/shQy8+YNBAIJyitRgunT0aoVhgxhW4fsoqKC4GCsXo3+/bFwoUR8FTJMSMjWx4+Hp6WdFgqjwsLmsy2HHVxcICeHvXvZ1iGOI7x69aqHh0e3bzCgiXlq18a0aSm1apkS0oOQ7gUFLZ49e8a2qJJ4eODkSXz4wLYOKpDk4eCZMzh2DOvXs62jBjBgAG7fxoUL6N4dqakQCoWpqalCyV+0p4Ls7DyhsGirRK3cXKna0UUpNjYSsUxYQYo1gUAwZsyYoKAgIyMjZgSxhZGRkZLSU+AeIJef/8BE8laua9VCnz7Yu1cWlhMkdoHwwweMGoVdu6ClxbaUmoGxMWJisHo1bGxeiETufL6hklLqxYsHJfABpJbRo0euXt1XWfmqsnLs0qU1dEQIYMwYiUisWsGIMDc3187Orm/fvrbfYEYW8ygqKp46taNr15Vt2izl8yNevVJlW1Ep+Phg61a2RVCBxI4IJ0yAhwecuRBRBpGTg78/mjdfnZa2+r//jv3776rZs1ezLYp2IiLq9e59+cQJt/v3j/Tt24ttOazh5QWBgP2N3RWMCDU0NLS0tF69etVAYtdzqMPS0vLs2b0Atm6FuztiY6EqYd7QyQl5eZI7nBITQnDnjiS+hYgIPH3KRfazg7p6AVD0vKl+/FjAshqauX8fmzfj3j01IyNpK+5MNXw+LCywfj169GBVRjn/1qhRI21t7e3btzds2FD7G4wpY5HRo2FjI4kzkDwevL2lPmTm8WPo6EhcPZCXLxEQgJ07wSUUZYVFiybVqTNWV3dirVrjb93yk6CcI1QjFMLHBytXQtaXm8TF0xOXL7OsoTxH+Pz580+fPhUlEfj0DcaUscvGjbh+HTt2sK3jF7y9ERUl3RuwJHBeVCSClxfmzkXLlmxLqalYWVklJf1z/LhnSsql8+dbzpsHT09Ic1bHMvnzT2howMuLbR0Sw++/IycHsbFsaqjE9gmhUJibm0ufFIlCTQ1RUZgxA48esS2lOHXqoEMH/PUX2zqqgQRO7S5dCiUlia4JXhPQ0NBo06aNhoaGjQ3u3EFhIVq3RkIC27IoJTkZK1diyxbweGxLkRiUlWFigrVr2dRQzBE6ODiEhoYWvSaEDB8+/MZPOVcOHDigKmmLZnTSvDlWrIC7u8T9LJX2DYWSNiKMi0NICCIiuO8mCUJDA/v2YeZMdO4sEfvMKIEQ+PrC3x+NGrEtRcJwc8O5c2wKKOYI3759m5mZWfSaELJ3797k5GQWREkM3t6wt8eYMWzrKE6vXnjxQiprKwLIz8ejR7C2ZlvHN3JyMHw41q9HRSmVOFjA0xMxMVi6VEamSbdvR1oaJk9mW4fkMWMG0tORlMSaAC6zTAVs3IiEBElJiFeEvDxGjsT27WzrqBL37sHcHCoqbOv4xh9/wMEBgwezrYOjDJo3x40bEAphZ4eHD9lWUw3evUNAALZtg3wFofo1EW1t1K2LNWtYE8A5wgpQUUFUFGbORHw821J+YvRo7NollaXdJGpe9PRpnDnDJZGRdDQ0sGcPZs1Cp05SvI/2998xZowEzYVIGi4uiI5mzTrnCCvGzAzBwXB3R1YW21K+YWJCVFTWNW/eb9asJQUF0rTpSnIiZdLSMGoUtm2DpibbUjjEoCjIPjgYnp6Quipwx4/jwQPMncu2DgkmIADv3uH9e3asc45QLIYOhaMjfH3Z1vGNsLBtb98+ffo0bP164bx5K9mWUwkkZ0Q4fjy8vODkxLYODrFp1gw3b0JJCXZ2ePCAbTVik5mJCROwdSuUldmWIsE0aABtbdZiR0s6wjlz5sjJycnJySkoKAAYNmyY3Dc8PDzYUCgphIbi8WNs28a2DgBATMytvLxRgGFuru/FizfZliMuGRl4+xbNmrGtA9i8GS9fYsECtnVwVBIVFWzZgoAAdO6M4GC21YjH9OlwcUGHDmzrkHi6dcPBg+yYLrZu6+7u/kE2qhvQgLIyoqLg6AgbG7RqxbKYPn06nj27LjNzppxcRO/eUjOouXkTNjaQk2NZxosXmDsXMTFcEhlpxdMTrVvD3R1xcdi4EerqbAsqm0uXcOKEdIf5MEZAAFq1wpcvLHygxRzhqlWrmLYvVTRpgvXr4e6O27dZXljy8homFAoPHVp761abdu0kbHtH2UjCvKhAgOHDsWABWrRgWQlHdWjaFLGxmDULdnaIikLdumlJSUmWlpYaGhpsS/tBfj7GjUNoKGrXZluKNGBpCTU1hIbC359p09waYeUYPBhOThKxWOjjM+LEie0bNoybOVNOWuqaSkKkzJIl0NSUiMovHNWkqLrvnDno2PGSqWmvvn2jzM07SlQl0fnz0aoV+vVjW4f00LEjdu1iwe4PR/js2bO0tDRxzklMTExPT6dNkqSzYQOSkrBlC9s6AACDB3+ds5Vwbt+Oa9y4fXS0zcGD0wgh7MnA5s2IjOSSyMgOI0bAxGTdly8HPn0Kfvs2aOXKzWwr+kp8PHbswLp1bOuQKqZNQ2IiBAKm7f5whPHx8aampjNmzEgoI7sfIeTSpUvDhg1r1apVRkYGUwoljiLHM3cu7t5lWwrA4yEoCAEByM9nW0q5eHhMff58r0h05/TpnNOnT7OiITsbw4cjOBh167Jin4MudHSUgaJvpIzbt5X/+49lPQAEAvj4YNUqGBiwLUWq6NIFCgosJDD5sUY4cOBANTW1mTNnrl69umnTpg4ODmZmZtra2gKB4NOnT/Hx8devX09NTe3Ro0dcXFzDhg2ZVipJNG6MDRsweDD7i4UAOnaEpSVCQvDHHywrKYesrBygPoDs7JavX6eyomHKFHToADc3Voxz0Mj69fN++21EYaG+hkamo+OR5s0xfDhmz2bTCa1Zg1q1ULMD7auIvT22bsXo0cxaJcURiUQxMTEjR46sX7/+92N4PF6LFi2mTJny8OFDQgV79uypX7++urp63759P378+OsBiYmJbdq0UVVVbdGixbVr14oao6KiTH/i/v37ZfXv7u6+b98+SqSWw9ixxN2dbiNi8fgx0dMjHz7QbignJ6ewsLAKJ06ePE9Obri8/BojI+t3795VoYeicmBV4OrVa05Obra2nsbGj6vaR42jylebLYRC4fv374tep6QQPz+io0P8/Ul6Ou2mf71WSUlER4c8f067abag9fbYs4coKFDWm1AoLCgoqPCwko7wZz5//vzkyZPk5OTs7GzKdBGSnJysrq5+9erV3NzcIUOG+Pj4/HqMra1tYGCgQCDYvn173bp1i95JREREt27dnn8jLy+vLBPMOMK8PGJjQ8LC6LYjFuPHk2nTaLdSZUd44QJp0ODijh0709LSqma6as/ex48f9fVtgUfADX19W5FIVDXrNQ2pc4S/kpxMfH2Jvj5ZsIBkZNBoqMS1EolIly5k3ToaLbIO3beHvDw5coSarihwhDSxaNGi/v37F71+8OCBqqpqTk7Ozwfcu3dPTU3te6OpqenRo0cJIREREd9PLB9mHCEh5OlToq9P4uIYMFUB798TXV3y7Bm9VqrsCAcNIqGh1TJdtWfvxo0bOjoTAQIQA4Oe3wcNHOUjA46wiMREMmIEMTAgQUGk+NcMZZS4VuHhpHVrIhDQYktCoPv2sLYmnTtT05WYjpCF7RNJSUktvxUCt7CwyM/Pf/369c8HPH36tFGjRirfKhS0bNky6Vt9jkuXLtWvX9/Ozi4kJISwF3z4ncaNERKCwYPx5k1Waio7S19F6Olh8mQJTWaYmoqYGHbWSywsLPLyrgMXeLyjGhqf9fT0WBDBwR5Nm2LnTpw/j7g4mJsjOJjesLK3bzFvHiIi2M8aIdV4e+P6dUYtslARJD09/bsj5PP5ampqnz59+vmAT58+qf+UWkBTU7PoAAcHh9OnTxsbG9+7d2/UqFFycnLjx48v1cT9+/ejoqKGDh1a9Keurm5SUpI8PeVPevRASMiRRo2CNTV1zc2VoqN3yLH0EPj6wtZW/fz53NathTSZyM3NVVBQqOyVXL9eyc0NPF5+dbKWf/nypQpn7dypoKe3q1u3DdraKv7+Ur75HwAAIABJREFUkVmSkzddsqna1ZZY6tdHRATi4uRWrFBcs4Y/fXrBiBGFVH0f/Hytxo5V8fERNWhQrVtd8qH79vDwwJQpGqdP57RrV92vMpFIpKCgUJQxtDyoGX9WBk9Pz7lz5xa9FggEfD7/efFl5UOHDrVs2fL7ny4uLmvXri3Rybp16zp16lSWCcamRouoW9cayAVIrVrTTp48yZjdX4mIIB070th/FaZGCwqIoSGpfpRVFWZjLl0ideqQp0+ra7oGIjNTo79y9Srp3Jk0bUp27KBmAvP7tTp4kJibk9xcCvqUcBi4PczMiHjrYBUguVOj5ubm8d+K+z18+FBVVdXQ0LDEAc+ePcv5VpH6wYMH5ubmJTqRk5MjEjA1WgQhAOQAFBYqFbJaJNDLC1lZOHaMRQklOXgQFhZo3pxpuy9fYsgQ7NmDxo2ZNs0hybRrh/PnsWEDQkJgZYWQkOctWnSuW9eub1+v6jy8GRmYOpUrMUEZQ4fiwgUG7VHgcytJSkqKurr66dOnMzIyXF1dx40bV9S+bNmyvXv3Fr12cHCYNWtWdnZ2SEiIsbFx0Sjk4MGDiYmJGRkZFy5cMDY2Xld2YBbDI8JVqzbq6XWoVctDQaHrq1f5jNktlZMnibk5EeM3UFWowoiwXTty+DAFpiv1IzQzkzRvLikBvdKIDI8IvyMSkSNHiJqaK3ALIMrKizZu3FKFfoqulbc3mTSJaomSCgO3R0YG4fFIfHx1+5HcEWG9evV27949derUxo0bq6iorFixoqg9LS3te8KaPXv23L5929jYePfu3UePHi1alHrw4EGfPn1MTEz8/PwmT548adIk5sWXyvTp4+/d23/58kx//9MjRyoynx/oZ3r0QP36klLIOz4eKSno04dRo0Ihhg1D584YN45RuxzSBY+Hfv1gaPgeaA4gL89y+vTUrl0xeTK2bMGNG5UoxH3hAmJisHQpjWprGpqaqFePufKEPCIxE4wUMnjwYFdX1yFDhjBsVyRCnz5o0QIrWa2VGx+PHj3w5An1WW8qGywzZgxMTDB7NgWms7KyxCwsMHUqEhJw4gToiY6qEYh/taWd4ODNixYdy8zsrq0dGR29t6CgaVwcHj1CQgLi46GhgebNv87tW1jA1hbfgtkBQCAQTJo0LS7ueUrKss2brVxc2HsbzMLM7eHnhwMHUM2EeSKRSCgUVhgsw31VUAmfjz17YG8Pe3s2U3lZWaF7d6xahcWLWdMA4PNnHDqExERGjUZG4sQJ3LjBeUEOsZg82bddO+tHjx516XKsXr16ABwdf/xraiqK/OKVK9i8GU+eQF//h18MDOzy4kVzYBSP56mktAroztrbkEXmzEFICN68gZER7ba4bwuK0dLCoUPo3p2dCJHvLFsGKyv4+sLYmDUNERHo04fRfI9XrmDmTFy4AC0t5oxySDv29vb2ZZQHMzSEoSG+D/UKCvD4MRIScP8+/voLL168BC4BICRnw4at3btzjpBKDAygq4tVq5io4PFjjfDJkycHDx6k3WANwMoKa9ZgwABkZrKmwcgIvr5YuJA1AYRg82ZMnMicxeRkuLkhMpLN3x8cso2iIiwtMXQoli/HsWNQV1cGrgECHi+6TRtLttXJID174vBhJgz9cITXrl1b+O2Ls06dOteuXWPCvozi4YFOneDpCRZXYAMCcPIka7WiTp6EmhratGHI3Jcv6NsXs2ejVy+GLHJwnD27U1NzlLy8abt2eXMlM6uTlOPvj5QUfP5Mu6EfjlBLS+vTp08iaSl2LvGEhODDB6xZw5oADQ3Mnk1NoEoVCA2Fnx9DtkQiDB8Oe3tITBwxR43AwcEhI+Pxp08JV65I0tZdGcLCApqaCA6m3dCPNUI7O7vPnz937969UaNGmZmZK1euNChteSc8PJx2UTKBggL27UPr1rCyQrdu7GgYOxYbNuDsWaYFvHqFmzfB2ER7QADS05kzx8HBwRidOmHfPixYQK+VH46wXr16Bw8eXLlyZXR0dH5+fkxMTKlR8pwjFB9jY+zZgxEjEBuLevVYEKCggOXLMWMG7twBn8EtoyEhGDUKqqpM2Nq1CwcPIjYWiopMmOPg4GCS6dPh5ISCAnof8GLfjn369Pnnn39SU1P19PROnTr1qTRo1CKLdO6MSZMwaBC9Oe/LYcAAaGhg927mLObmYscOhjazX7+OP/7AsWPgqkpwcMgkjo5QUsKWLfRaKX2YsGXLll/Te/7Ms2fPNm3aRI8kWcPfH8bGmDaNNQGrV2POHHxL3Uo7e/fCwQEmJrQbSk3F4MHYtg0tWtBui4ODgy0cHLBtG70mSneELi4uOjo65Zz25MmTdQxs7pAJeDxs24YLF7B9OzsC2rSBgwMTC85FhIUxsWsiNxf9+mHyZKbzt3FwcDBMr17xd+/OsLLqcuvWLZpMsJBrtAaioYHDhzFrFu7cYUfAihVYs6a6yYrE4do1ZGTQHptDCLy90bw5/viDXkMcHBzskpeXFxDQm5C29+9PbNfONYeeqS3OETKEuTnWr8fAgfj4kQXrpqYYPhxLltBuKDQUEyfSHpgzfz5SUsCFbXFwyDzXr18XiVoCA4ABIpHtxYsX6bDCOULmGDwYrq4YOhRCugrIl8fChTh4EElJNJpIS8PJk/D0pNEEgL/+wu7dOHwYSkr0GuLg4GAda2trIAFIAp4S8qh169Z0WOEcIaOsXImCAiZGZr+ipYWpUxEQQKOJ8HC4uUFbm0YTd+5g4kQcPQp9fRqtcHBwSAi1a9feunWxhkZfZWV3IDIpSZcOK5wjZBR5eRw4gIgInDjBgvXJkxEXhytXaOlcIEB4OMaOpaXzzMzM6Ojoe/feubpi82ZYcmkdOThqDN7eIzMzH+fm3u3WrX2vXigooN4E5wiZxsAAUVEYNQovXjBtWlkZS5Zg+nRaMqAePQoTE9jYUN/zrVu3tLWbeXgcsrbubmt7pF8/6k1wcHBIPv/7H+TlQUeRjyo6Qj6fL351Vo4SODhg9mwMGMDc3r7vDBtGnj/319GxadGicyKlpQKLwmTowNc3QCgMIeQQcPjcOTrndjk4OCQYeXn88w+uXMHq1RT3XAlHmJOTk5WVVfS6Z8+eDx8+pFhLTcLPD61awdeXabsnTvwvJycjPT0uISF0yBDKsmInJuLxY7i6UtVfMQQCOaCw6CWPx6PFBgcHhzRgYYFVqzBrFsVxf6U7wiFDhmzdurVE419//dWoUSPCYmEh2SIsDAkJBa6uwa6uvn//fZQZo8nJr/Py7AAe0PTjR8qqm4SEwNeXlmSAx44hNXUVnz+Nz+/B53cNCWGpmgYHB4dkMHUq7O3h6AgKSyWV4ggLCwsPHz5sZGQE4N9//z158mRRe7du3dLS0pKTkykzXrNRUUHz5guPHPl05Mg4H5+tp0+fZcBo37699PRC+fzNCgo+amq9BAIK+szKwv79GDOGgq5+RiDArFnw88Px45bZ2c+OH5/88eODESM8KDbDwcEhbVy4gLw8DBhAWYelOMIP/2fvzuNizv84gL9muqNEdymqFSFtmxxJKYQIseTIukOO+G3YsORYwjpyJ8cmt4RyRxFyJJH7CEWnu+ie+fz+mEWbqWaa71z1eT72j3zn8/183jM7zbv5fj+f9+ft25KSEhMTEwCXLl3y8/PjHdfT02OxWG/fvmVs8Drv+vWLwFzgl48fp0ZHX5DAiKampjdvRq9fT44cGdKs2UJvbwYWNYaFoVs3GBszEd9Xr17B2RkPH+LWLXTsCFVVVUdHRy0tLSbHoChKPqmqIiYG0dGM1SDlkwiVlZUB8G4Hvnv37tuOE7m5uYSQqmuQUkLp0MFOSSkUyGKxwhs2tJfMoI0bN/b1ndC7t9uhQ3j7FuPGiXqFISSE4Wky0dFo2xZ9++LoUfGuSqQoSk61b4///Q8TJyIjg4He+CRCbW1tIyOjv//+++7du9u2bVNUVDx48CCAtWvXamtrN5bKxnq11ObNS8eNe9227Xg/v847dvSX5GZJANTUEBWF9HSRcuH58yAEnTszE1JZGQIDMXUqjhzB7NmgM2MoiqrMypVo1gwdOzLRF+Fnx44dbDYbQPfu3UNDQwE0bNgQwMqVK/m2lzWDBw/et2+ftKMQzv37pGlTsmCBpMf98oW4uJAxYwiHU33jgoKC0tLS8kc8PcmWLcxEkp5OOnUiHh7k3Ts+j+bl5TEzDCUA+moLrg6+VrLzlD98IGpqZMyYShtwOJySkpJq++G/FnD06NEdO3ZMT093cXFRUlLS0tJKTk52cHDo3bs3E8mX4qNlSyQkoE8fvH6NLVsgsVWa6uo4fhzu7pgwAVu3Cvcl7NUrXLiAsDAGwjh+HOPGYeJEzJ8v9prdFEXVDlpaOHAA/ftjwACIkp1YpDYuh/Dy8vL09BwyZIi0AxHa588YPBhKSti3D+rqkhv3yxe4u6N5c4SEVJULCwsLlZSUvtVSmDsXX75AxI0py8qwZAl27sTevejUqdJm+fn5GhoaIo1ECYy+2oKrg6+VrD3lkSNx8CAyMvhMKeByuRwOR0lJqeoe6N/esqV+fURFQV8fLi7IzZXcuPXqIToaKSmYMUPQU0pKsGMHJk0SadzXr+HigitXcONGVVmQoiiqMmFhMDCAs3PNe6CJUOYoKiIkBL16wcFBvLsmVaCpiTNncPWqoLnw4EG0aYPmzWs+4vnz6NABXbvizBno69e8H4qi6rhr1/DkCabVtFgWTYSyiMVCYCDmzIGTExISJDdugwY4exZXrgi087soxUU5HAQGYswY7N+PwEB6U5CiKJHo62PbNmzcWMPddegnkOwaMwa7dmHgQBw/LrlBebkwPh7+/lU1S05GZmYN707n5KBnT1y+jBs34OhYszApiqL+Y8QIuLujZ8+abGZAE6FMc3NDVBR8fLBli+QG1dJCTAwuXMCsWZW2Wb8evr5QUBC0z0ePHtnaupma2g8bttTODp064exZejmUoigmHTuG+vXRvbvQJ1aTCL98+bJt27Z3797VMC5KZPb2uHwZa9fCz08s+wjypaWF06dx+jQWLuTz6IcPOHoUo0cL0eGgQZNv317z6tW1AwcezJwZRy+HUhTFODYbFy7g+nUsXy7kiVU//O7du/Hjx6enp9c8NEpk5uZISEBSEkaNQmmphAbV0cH584iIwOLFFR/atg19+0JPr/pOysqQkIBFi/DkySegFaDAZjuqq0t8S2KKouqGFi2wejXmzEFKihBn0T/L5UOjRjh3DgUF6NULeXkSGlRXF+fP48ABLFny/SCXiy1bqpkm8/w5tm7F4MHQ1cWIEcjIgIND53r1/FisMG3trT17imGHaYqiKADAtGlwcEDXrhB8dx26y7zcUFXF/v2YNg2Ojjh5EpKp+aqnh9hYuLiAzcacOQBw6hRLSwv2P1QIz8lBfDzOncPp0ygthaMjunXDmjX/7krB5a46fDjyxYvXXl5HeBubUBRFicn58zAwQL9+iI4WqH01iVBRUdHQ0LDaZfmUZCgoYONGBAfDwYEMGxaalHSxS5e2s2dPVRRnQbZvufDRo7MnT04sKBi2ePEUwADAly+4ehXnzuHcOaSloUMHODoiMhJ2dhU7YbPZgwb9Kr4gKYqivlFWxpkz6NChtHPngJUrhzs42FbdvpoPUCMjo8zMTObCoxjg54eUlPAVK64TMj8hYUdx8epFiyqf38kEfX1s3PjQ1XUksAR46+9vX1b26tw5JCbC3h7duiEkBLa2dP4LRVGywt4eamquCQkuZ8/qOThU05h+dMmlT58SCPEFmhcUTN+z5/KnT2Ifcf36hUAfYCwwG9B5+vRZQACysxETg9mzYWdHsyBFUbKlsDADWBQWZlRtS/rpJZfc3DpqaGwBUlVU1ikpOTRpgr59ER4OxjNiTg5CQuDmhjNn3ICrwBcgA8jZsqWpqytUVRkejqIoiikKChzgyaRJ1e+2ShOhXBo//reFC39xdp4zZ472vXv+WVkYPx4xMTA3h4cHdu0SNSO+eoWtW+HhAUtLHD8Ob29kZ49p3lwDaAa0Hzq0h1jvSlIURYlu7941amq9NDXPV9uSbsNUqxQW4tw5HDqEEyfg4IBBg9CvHxo0EPT0tDQcPYpDh/DoEdzdMWgQ3NygolK+//9swyRJsrbzS+1GX23B1cHXSo6esoDbMNG/62sVNTV4eMDD43tGnDGjYkYsLCxMTEw0MTExMzPjnfX8OaKjcegQHj9Gr16YPRs9e4LOFKYoqo4QKBG+fv06JSXFxsbGmLcojJJ5lWVEd/f3y5f3zM/vyGan+Ph4KymNPXQI+fn/5r9evUAveVIUVdfw/9jz8vKysrIKDAwEcOHCBXd398LCQhUVlf379/fv31+iAVKi+ZYR8/MRHY2VK4+kpY0ApgLFK1Y4/f772J07YW9f1a70FEVRtRufyTJlZWXHjh1z+LryIiAgoFmzZpcuXRo+fLifnx+Hw5FshBQzNDQwbBhmzaqnpsYroZ7ftKliUBDataNZkKKoOo1PInz//n1xcbGFhQWA3NzcGzduzJ4929HRccmSJenp6bQAt1z79deBbdve1tNzNjDotmXLX9IOh6IoSvr4XBrlTbApLi4GcOrUKUJI165dATRq1AjAu3fvvk2yoOSOkpJSfPzR/Pz8evXqsekaeIqiKL7fCBs2bGhsbLxt27ZPnz5t27bN1tZWX18fQFpaGgBdXV3RR01ISBg/fvz48eOvXLnCt0F2dvYff/wxfPjwbdu2cbnfl0NGRUWNHDly6tSpDx8+FD2MOktDQ4NmQYqiKB7+n4ZLliwJDg7W0tJKSEiYO3cu7+Dx48d1dXVNTU1FHPLWrVs9e/a0tbW1s7Nzd3dPSkqq0KC0tNTZ2fnDhw/9+/cPDg7+669/r+BFRESMHz++e/fuenp6jo6O2dnZIkZCURRFUfxnjY4aNcrOzi45Ofnnn39u06YN76CRkdG6detYIs+sCA4OnjBhgq+vL4C0tLS1a9eGh4eXbxAVFcVisbZs2cJisUxMTDw8PGbNmqWiovL333//9ddf3t7eAJKTk7dt2zZv3jwRg6EoiqLquEpXjVlbW1tbW5c/wlShlmvXrq1evZr3s7Oz89SpU39s4OTkxMu47dq1+/z587Nnz1q0aHHz5s3du3d/OzEmJoaReCiKoqi6rNJE+Pr1623btt2/f7+wsPD48eMAjh8/rqGh4ezsLOKQ2dnZ2travJ91dXWzsrJ+bND467azbDZbW1s7KytLW1ubw+FUfeI3z58/X7169cGDB3n/VFFR2bx5My2PKTopllj78uWL6FcjKAHRV1twdfC1kqOnzOVylZSUalhiLTEx0c3NjXdl8t073rIz3Lt3LywsTPRZKqqqqiUlJbyfi4qK1NXVKzRQU1MrLS399k9eGzU1NQBVn/iNnp7eTz/91KlTJ94/NTU1NTU1RQybAsBisaSVCDkcThX/xylm0VdbcHXwtZKjp8zlcgWpp83/E23ixIm2trZHjhxJTk7m3ZMD0KdPn4CAgNzcXD09PVEiMzEx+bYYMT09/ceybcbGxvfu3eP9nJ+f//Hjx8aNGzdo0KB+/frp6em8Kax8T/ymfv36HTt2HDx4sChxUj9ifyWtoSU/bt1EX23B1cHXSr6esiBFYPg8mQ8fPty6dWvx4sUNGjQo//23SZMmADIyMkQMa8CAAeHh4YQQQsju3bsHDhzIOx4REZGTkwNg4MCBMTExvEmhe/fu/eWXX3hTVQcOHLhr1y4ARUVFhw4d+nYiRVEURdUYn2+EvKX0P+6y8eHDBwCiXxbz9fWNiIho164di8UqLS2dPHky7/ioUaOOHDnSvXv31q1bjx492t7evnXr1klJSREREbwGf/75p4uLy71797Kzsy0sLDw9PUWMhKIoiqL4ZDU9PT1dXd2TJ0+2adOm/DfC/fv316tXz9LSUsQhtbS0EhMTExMTAdjb2ysoKPCOP3z48Ntq/TVr1kyaNCkjI8PW1lZLS4t30MLC4smTJ4mJiRoaGjY2NvJyt5aiKIqSZXwSIZvNnjFjRmBgIIfDMTY25nK5Dx48OHDgQFBQ0IwZM1TK79NaUwoKCh06dKhw0MTEpPw/LS0tf0y6qqqqnTt3Fj0AiqIoiuLhf51z9uzZb9++DQwMLCsrA9CqVSsWizVy5MhFixZJNjyKoiiKEi/+iZDNZq9atcrPzy8uLi47O1tLS8vJycnKykrCwVEURVGUuFU188XU1HTkyJESC4WiKIqiJI9/IkxLS6ts7YW5ubk446EoiqIoieKfCNu3b89b0vcjQVbpUxRFUZS84J8IQ0NDi4qKvv0zPz//4sWLx44dW7p0qaQCE8mXL18uXLjAVJVwiqpTLl686OfnN2zYsFmzZkk7FoqSBJbg3/CCgoJOnz594cIFccbDjNatA+7fv2BkVJyRcUvasdQqUiy6nZ+f/2ORB4pxCxYsWLRoK9AHiDcyQkbGY2lHJOvq4DtTjp4yl8vlcDjVFt0WIhG+evXK1NQ0NTVV9m8TGhjcy8lpCZgVFz9VVlaWdji1B02EtVhBAZ4/h42NMZe7GvAC3gE2hLyWdlyyrg6+M+XoKQuYCIX4RPv06ROAgoICkeKSCEJYQAnAolmQqoNKSko2bNh2717quHG/Ojh0rPAol4uMDDx/jufP8eLF9x8+foSZGQBlgLfhzEeghbc3Vq+GaGX2KUrWCTRrtLS0NDU1dcGCBVpaWqKXWJMAHZ0NublngCXLliEgQNrRUJRkjRnze2SkTmFhn6iombt2bVJWbsPLdrz/Hj2CsjLMzf/9z8kJo0bB3BxNm4LNxtixO3fsGALsADIsLUdcvAgDA7Rogb//hru7tJ8YRYmHELNGjYyMwsPD5eI7lpXVm7Fjly5cOGTuXHz6hKAgaQdEURIUG3utsDARwLt3vmPHnu/QoY2Z2b85z8wMZmZQVeV/4u7dOHOmS8OG2QsXbg8J8d60ScXJCcnJmDEDHh7Q0oKvLxYuhPzswENRAhFo1qiiomLjxo2tra1VK/sFkjEKCgpGRli5EitWYOVKfPyILVukHRNFid/Ll1i+HG/fNmWxjhPipKV1dM+eia6uAp174gRmzoSbGwwMMGrUYEBl82Y4OcHWFhcuIC8PM2di1SosX44+fbB5M/T1xfxkKEpiSG00ePDgffv2cbnEzY2MGUMUFMjIkdKOqVYoKCgoLS2VytB5eXlSGVdepKaSadOIjg6ZNo3cvZvTv/9YK6sua9duEfD0a9eIri45c4Y0bEjS00leXl5eHtHWJtnZFVuGhxNzc8JiESsrcvw4w89CHtXBd6YcPWUOh1NSUlJts9p8jYPFwtatiIrCjh3Yuxd0+0KqVrp/H7/9BgcHNGyIJ08QHIzWrfWOHNn24EGcn98EQXp48AD9+mHXLty4AU9P8LaB0dDAgAHYsaNiY29vpKYiKQl6eujbF9ra+PNPcLlMPyuKkqDvl0ZPnToVJMDNtIsXL4ozHoY1aYIFC7BlC+Li4OICFxfExUk7JopiSEoK/v4bMTGYMAGPH6NBg5p0kpEBd3csXw4XF4wZg5iY7w9Nnoy+fTFrFr7uGfrdj9dLu3XDxo3c0aMHJibeb9bM9MKFiG87iVKUjPv+jVBJSam+AKQYa834+kJZGUlJuHUL166hY8XJ5BQlfxIS4OGBHj3QqhVSUxEYWMMs+OkT3N0xbRpGjsQ//8DeHq1afX/UxgZGRjh1qtLTNTUREoKCAuzYgcePYWGx7eJFw4KCy3fudHV2/rUmAVGUNAixoF6OeHl5eXp6fiux9uIF2rdHQgIUFGBtDTMzJCdDGovC5R5dUC91ly9j+XLcu4fp0zFhQqXzPwVRWAg3N3TqhKAgcLlo0QI7dsDRESj3au/ahQMHcOKEQB0aGjpnZ68A2gNv1dQcCgqe1Dw4+VEH35ly9JQFXFBfm+8RfmNmhj/+wPjxaNoUDx8iLQ3W1igrk3ZYFCWMy5fRtStGjoSHB54+hZ+fSFmQw8Hw4WjSBMuWAcCRI2jY8N8sWJ6XF5KS8OKFQH16eXVksRYCicBCQgbk5tY8PIqSpKr+tM/NzX3+/Pnnz5/LH+zWrZuYQxKL6dNx5AhCQjBxIp49Q8uWsLDAw4dQV5d2ZBTFT1JSUkDAahUV5ZUr/3j6tPmSJSgowMyZGD6czx07YRECHx+UlODgQbBYALBqFf74g09LFRV4e2Pr1n/zZdXWrg3icPyPHp3avr11SsqaJk1w+jScnUWNlqLEju9c0pycnK5duwreXtbwlk9UOPjoEdHRIamphBDy6RMxNCT6+uTDBymEJ7/o8gnJyMvLMzD4BbgNXFVS+sXGhnPwIOFyGet/9mzSrh35/Pnff164QCwtCYfznwC+/fzsGdHTI4WFQo/i7U3YbBIUJGKwsq5OvTN55Ogpi7R8wsfH5/79+3v27OnTp8+YMWNOnz49bdo0LS2t8PBwcSZl8WreHL//Dh8fEAJNTTx7BiUlmJujko0XKUpqnj17VlpqB9gAHerXNzlxImvQoH+/uolu0yZERiI6GvXq/Xtk5Ur4+1daL8bCAj//jMhIoQcKD8f69Zg7F+7udH0FJdP4vPe5XO7Zs2dXrVo1bNgwHR0dAwODHj16BAcHBwYGBgUFEXmeXOPvj7w87NwJAOrqSE2Fri6aN0dGhrQjo6hyGjSwzMtLAuJZrDMNGmQaGRkx1fP+/QgKQkzM9zrajx4hKQkjRlR11qRJ2Ly5JsP5+uL6dcTHw8wM9JYhJbP4JMI3b94UFhba2dkBUFFRycvL4x0fMWLE/fv3nz17JtEAGaWoiLAwBATg9WsAUFbGw4cwM4OlJZ7UiQlulBx48QI9e9bz8dk9ePD+ESNOXbhwiMXQl8HYWPj5IToaTZp8P7h8OaZNq2bejYcHXr1CSkpNBrWzQ2YmVFXRtCni42vSA0WJG59EqKVtecP+AAAgAElEQVSlxWKxeJsuGRsbP3r0iHe8sLAQcrINUxWsrDB5MiZO/PefbDaSk2FnhzZtkJws1cgoCrh7F87O8PfHhg1WBw5sCgtb26R81hLBzZsYMgSHD8PG5vvBjAwcOwYfn2rOVVDAmDE1L9irqYnHjzFwIFxcsHx5DTuhKPHhkwhVVFRatWqVmJgIoFevXnFxcStWrDh37ty4ceMaNGjQrFkziQfJsDlzkJWF3bu/H4mPh6sr7O1DFRSaKCiYDB48XnrRUXXXhQvo2hWrV1efmYSVmop+/bBlS8UFEmvXYvRoaGtX34OPDw4cwNfLQzURHo516/69ZUhRsoXvFJrdu3evXbuW9/O0adN4V2bq1au3f/9+BufziA/fWaPl3b5NDAxIVtb3I5mZmUBLoAgoYbNbPHnyROxRyiE6a1R8jhwhurokJob5njMziZkZ2bat4vFPn4iODklL43MK31f711/J5s2iBnPjBqlXj5iakjdvRO1KRtT6d+aP5OgpizRrdPjw4X5+fryfg4ODc3Nzb9y4kZGR4eXlJbkULU42Nhg7FhPKVSROT09ns/UAFUAJMH8h4BJiimLCzp2YMgVnz4LxZbp5eejdG+PHY+zYig9t3gx3d5iaCtrVpEnYsEHUeOztkZkJFRWYmtJbhpSs4J8IM/47jVJHR8fe3r5BzaoZyqr58/H8OQ4e/Pef7du3V1d/BUwHAtjs164C7uFGUSJbvhxLluDCBfz8M8M9l5Rg4EB07IiAgIoPFRdj3TrMmCFEb7zfiStXRI1KUxNPnvx7y3DFClF7oyjR8U+EdnZ2bdu23bp1a35+voQDkhhlZWzfjmnTvs/qfvPm3rx5Gn37anC5N7OyaClSSuwIwf/+hz17cOkSfvqJsW4zMzPd3Ib99JPDL7+s0dHB+vV82uzeDRsboVPv+PE1XEfxI94twzlz6C1DSgbwvWAaEhLStm1bAPXq1fvtt99iY2M55ctOyLxq7xF+M3MmGTKk4sEWLYidHfNR1QL0HiGDiouJlxfp0oV8+sRwz46O/Vmsc0CRouKw48fP/NiAyyUtW5LY2Ep7qOzV/vCBNGxIcnKYivTfW4ZNmsjxLcPa986slhw9ZVEryyQmJj58+NDf3//ixYuurq5NmjT5448/5HoRIV+LFyMlpWLVjOPHkZxck1IaFCWgz5/Rty+Ki3HqFDQ1Ge48Le01IV0BFQ7H/e7d+z82iIqCujpcXITuWUsLnp7/lqRghL09Xr+GsjJMTYvNzXuqqJi1aeNaob4xRYlbVbtPtGjRIjAwMDU19eTJk46OjsHBwZaWlhKLTDJUVLBtG6ZNw/v33w9aWODXXzFmDK0LRYnFu3fo3h1GRjh0SKQdJCrTrVtnBYU/gBPa2hs8PNx+bLByJWbNqmHnU6Zg82ZwOCJFWJ6WFp48Qf36U1686FRScv/u3c69e3sz1jtFCaD6bZgUFBRMTEyMjY01NTWJPNdXq0zHjvj114qzBsLDUVKC33+XUkxU7fXyJRwc4OyMHTvEtSnm9OkrNDSaT516/ezZTa3K77QLALh+HVlZGDCghp3b2kJPD2fPihpkBRzODcALUAeG37v3kuHeKapKVSXCt2/frlu3zs7OztraOjQ01MPD49KlSxKLTJKWLkVCAo4d+35EWRlr1mD9elqSm2LS/ftwcsLEiQgKEuMoGzYo+vuPXrduka2t7Y+PLlsGf3+R9nKqcenRKvTq1ZHFmg3EA/Py8iacP89w/xRVFb53DqOjo/v376+srMxms7t27RoeHv7lyxeGb2KKk+CTZb6JiyPGxuTp0ze5ubnfDpqbEwcHpoOTZ3SyjCiuXiWGhkTcRSnevyeNGpHsbP6PPnpEDAxIQUE1nVT9ahcUEB0d8uJFDSOszPTpMy0sHGfM+IO3f5O3N8P9i0kteGcKS46esoCTZfgnQn19fVNT09mzZ6fytu+TNzVIhIQQa+s/VVUddXWdfHxm8Y7cvUvYbHLqFNPxyS2aCGssKoro6pLTp8U+0LJlZMyYSh8dN44sWlR9J9W+2jNmkDlzhIxMGEePElVV0rgxef5cjKMwQt7fmTUgR09ZpER48+ZNLoPbgEpcDRLhmzdvdHQ6AwQgurrdX716xTvu4UG0tcUQonyiiVAopaWlT548+fLlS1gY0dcnV69KYERiakpu3+b/aHY2adRIoIUK1b7aT58SPT1SVCR8iAL78IHY2BBFRbJhgxhHEZ08vjNFJEdPWaTlE3Z2dkzt/CIvOBwOi6X89V/KZWVlvJ8OHsSXL5g3T1pxUfIqNzfX0tLB0XG+gYFjQMDV+Hh06CD2QSMjYW7+n/0lylu7Ft7e0NFhYKCffoK1NY4eZaCrymhp4fZtLFoEPz907Ag53/aGkmmVzlo7ceLE0aNHMzIySktLyx+PiYkRf1RSoK+v36tXq5Mne797x3ZwaNy0aVPecVVVLFmCgAD4+0NLS6ohUnJl7drQtLTpXO4w4KW5+RRLy+MSGHTdukqnOufnIzQU168zNtakSVi/HuIuPxwQgH794OICPT0cO4auXcU7HFU38f9G+L///a9Pnz6nTp2S990HhRIWFnz9+np397W9e/9n47Xff4e+fs2nm1N1U7krLixFRUmsSL11C69fo29f/o9u3YoePWBhwdhw/fohNRX37jHWYWVatkRWFjw94eaGESPEPhxVB/H5RsjhcLZs2eLr67tu3ToFUSZZyyFzc/OxY7FxI8b/d0fCyEh06ID4eDg5SSkySt506DCeze6jpXVCUfH+6tXrJDDi2rWYOpX/uojSUqxbhyNHmBxOURFjxyIkhH8tU2ax2QgPx+DB8PLChQuIj4eZmdgHpeoOPt8I3759W1hYOGbMmLqWBXnc3ZGcjMzM/xy0t4erK4YOlVJMlLx5+RKTJumfPn01Pj7g2bN4JyfH6s8RTW4uTpzA6NH8H927F82b45dfGB50wgTs3QuJVeb38EBmJrS1YWkpiexL1R18EqGOjo6+vn6FnZjqDhUVeHjg8OGKx48cwdu3WLxYGjFRcoVXR3TuXHTrptiyZcv69etLYNBNm+DlhUaN+DxECNaswcyZzA9qaAhnZ+zbx3zPlfk2g2bGDDqDhmIMn0SooKCwZs2aP//88+XLlxKPRyZ4eWH//ooH69dHYCAWL0ZenjRiouQEIRg7Fr/8gsmTJTdoSQm2bsXUqfwfPXkSAPNb/vJMmoRNm8TScxUCApCSghcvoKcHWoOGEh3/WaNHjhzJzs5u3ry5lZWVrq5u+Ydq66zR8rp3x6hRePkSX6eO/isgAMHBGDoUJ05IJzBK9i1ahFevEBcn0UEPHECbNrCy4v/oypWYPRtiWg/VrRsKC3HtmiQWh5TXsiUyMzFyJNzcMHQoV0vLPy4u2du7Z0DAbInGQdUKlS6faNOmjSTjkCmKivD0xMGDfCr0HzgAV1ckJcHOThqRUbLt2DGEhuLGDaioSHTc9euxcCH/hxITkZaGQYPENTSLhQkTsHmzpBMhys2g6d9/OZf7Hpg7d+6iwsLSRYvosl9KOCwijQ0lCCHXrl3LyclxcHDQ09Pj2+b58+e3b99u0aJFy5YteUfevn2blpb2rYGVlZW6ujrfc728vDw9PYcMGVLjCOPi4O+PpCQ+Dzk54dUrvHhR477lWGFhoZKSkqKYNk2oUn5+voaGhuTHFdyjR+jSBdHRsLeX6LhXrmD0aDx6BDa/xVCDBsHJqdKrppUR6tX++BHm5nj8GP+9eCQ56urNCwuvADrAtcaN/V+9uizJ0WX/nck4OXrKXC6Xw+EoKSlV3az6bZgYRwgZPHjw2LFj9+zZ06pVqytXrvzYZufOnR06dIiIiOjevfuSJUt4B6Oiorp16zbhq/JJkXHOzsjJwcOHfB46ehSvX2PtWvENTsmfDx/g4YFlyySdBQGsW4dp0/hnwefPcfFipVNJmaKlhX79EBYm3lGqYGKiBxwCioGIjIyhHTsiOVlqwVByqbLaa1evXvXy8mrdunXLli15R9auXbt9+3bRi7/FxcUZGxvzqtWtXbvW2dm5QoOioiI9Pb3Y2FhCyNOnT9XU1HJycggh27dv79+/vyBD1KzodgV+fmThQv4P+fsTFRUiVxtyMIPWGuWrrIz06kVmzZLC0K9fE21t8ukT/0cnTiTz59ekW2Ff7evXiYUF4XBqMpbosrOzmzZtr6Rk+ssv3WNjizt0IGw2MTcnhw9LYnRZfmeKiRw9ZZFqjUZHR3fu3PnBgwfm5uafPn3iHVRRUQkMDCQiX0o9evRo7969ed+shwwZcvHixfflt4cHEhIS2Gx2ly5dAPz00082NjYnefPegIKCgoSEBN6eGCKGUS0vr0rnha9cCQ0NWuSC+pe/P8rKsHSpFIbesAG//QZNzYrHCSEvX345cACTJkkijHbtoKWFc+ckMdaP9PX1X7y4VlKSlpR01sVF+epVpKfDxgaDB6NBA/z5J7iSKOxDyTH+N3umT58+ZMiQsLCw+Ph4b29v3kFXV9dJkyZlZmYaGxuLMmRGRobN16rA+vr6KioqGRkZjcotgHr9+nXjxo2/Vf02MTF5/fo17+enT58GBAQ8evSoefPmR48ebcR32RTw4cOHc+fOffz4kfdPTU1NL+FLItrbo6hI4fZtrrU1n6S7ezfc3RVu3eJUVuC4VuJwOGw2WyoF2TkcDofDkfy41dq9m3XyJDshgQNAwgEWFmLHDoVLlyq+MJcuXRk6dGpenmaDBtoNGuzjcJQr6aBSNXi1fXxYmzaxunaViZxjYIBDh/D5M/z92StXsv7+G97eZPVqbiWTCkQis+9M8ZGjp8zbRqnaZnwSYW5u7vPnzw8dOlThI4+X/7Kzs6tNhM+ePRs7duyPxzds2GBtbV1SUlJ+toWSklJxcXH5ZqWlpXwbDB06dMyYMQAKCwv79Okzf/78DRs28A0gPz+//LdGFRUVDw+PGhTKGTBAae9eBAaW/vhQ586wtVX18mLduVMkbLfyq7i4mHfzWfJDl5SUVHifyILbt9kzZ6qcOlWkpsaVfHRhYYr29mjcuOILM2HC3OzsM4ABh7M4PHyvt7fQJZFq8GoPGIA5c9SePi0xNZXC/Du+lJQQHIzgYGzdqrhsmdKOHQqdOnFDQkqaNGEyW8vmO1Os5Ogpc7lcQT75+SRCXvL7MYtmZmYCqGyiZnmGhoZL+V0n4m3pYGho+ObNG96RgoKCz58/GxkZlW9mYGDw9u3bb/988+aNo6MjADU1Nd4RNTU1b2/vzZs3VxaAqampiLNGeby98euvWL5cie9XoOhomJhg7171ceNEHEdusFgsac0a5XA4grz3JCk7G0OHYvt22NmpSiWArVuxZg2fX8mSklJAC0BZmX5hYXENXrcavNrq6hg+HHv3qi1aJOxoYjd9OqZPx4EDmDOH3aqVaps22LKFsfUeMvjOFDc5esoC/uHO5x6hrq6umZlZeHg4viZFno0bN+rq6lpaWlbbab169Trxw7sv2KlTp7i4OF6ijY2NtbCwMDQ05EXM5XIB2Nvbv3r1ijcptKCg4OrVq7xEWN6DBw94Z4mVrS1UVJCYyP9RQ0P4+MDPDyUl4g6EkjnFxfD0xIQJ8PCQTgCxseBw4OLC56FZsyaoqPRSV/dv3HjrsGFiW0L4g0mTsG0bSvlcQJEJXl5ITcWtW9DSQqdOMDFBWBi4XO6MGbM6dOi1d+9eaQdISQ/fKTS7du0CMHz48MDAQD09vYiIiEGDBgEIDg4WfRpPYWGhhYXF2LFjQ0NDTU1NQ0JCeMd79+49d+5c3s+TJk1q167dzp07e/To0bt3b97ByZMnL1y4MCQkZNKkSerq6leuXKlsCEZmjfIsWED+97+qGmhpEW9vRoaSA3TW6DfjxpEBAwiXK7UA+vYlW7fyf6isjDRsmHrgQOznz59r1nnNXu2ysrLGjf80MHDw9p5a46El4/Vr4ulJFBSIgsIfwEjgCIvVateuXTXoStbemRIgR09ZwFmjlS6fCA0NLV9crX79+kFBQVyGfu9zcnICAwMnT5587NixbwcjIyMvX77M+7msrCw0NHTChAlr1qwpKCjgHYyNjZ0zZ86ECROWLFnCuwVYGQYT4cOHxMioqnnhR48SNps8fszIaLKOJkKeNWuIjQ2R4kf9ixdER6fSAC5dIra2IvVfs1d7zZrNKiqzgCIlpY0+PtJYTSKk/HzCZrcASgACHK9ff+yUKWTPHpKTI0QnMvXOlAw5esqiJkJCSHFx8dWrVyMjIy9cuJCfn89cbGLHYCIkhNjYkIsXq2rQvHm6ispwHR2bVavWMDWobKKJkBBy7hzR1ydV/iUmdjNmkD/+qPTRP/4gf/4pUv81e7WHDp0CXAUIkGtr21OkCCRFV9caOA5wgImGhitMTYmaGgGIoiLR1iZt2hAvL7JyJblxg//pvXv3NjOz+OeffyQbtZTJzi9jtRhIhPKL2US4dCnx9a2qgbq6JXAIuMNmt+TVAaitaCJ88YIYGpK4OGnGkJ9PdHRIenqlDVq3JteuiTREzV7tI0eiGzbsDcSqqY0JClonUgSScvfu3UaNWikoNLa17cYpd+Xn1i2yfDnx9iY2NqRRI6KkRACiokIaNyYdOhAfHxIeTurVswA8gQ2A8ZYtW6T4LCRMRn4ZBSFgIuRfazQ+Pr6E3wwQTU3Npk2bVlYdVHaIXmu0vBcv0L49MjNR2WRJBQVzLvc5AGD51KkZ69ZJYjtyqajjtUY/f4aDAyZNktAq9cps2IBLl3DgAP9H09Nhb4+sLP5F1wRU41f71KkzS5acU1S0v3BhkFTWm4pPRgbOncPVq7h7F2lpePcORUWNgXSADew0MPgrK+uZtGOUEFn4ZRSQgLVG+X+iDR48OCcnp7JzOnXqFBYWZmFhIVKA8sPMDObmiI2Fmxv/BmpqCl++HAJaAAf69Vsp2egoCSEEY8bAzk7KWZAQbNyIbdsqbRAVhd69RcqCoujVq4e5eQ83N3Ht+iRFxsYYORIjR34/wmIR4BHQEkjicEalpaFJE+nFR4mA/6/L5s2btbW1fX19T506lZiYGBUVNWzYMGNj42PHjm3atOnFixceHh7yUlmAEV5elf4BDuDmzWhLy+AGDcayWMvZ7K4SjIuSnIULkZGBLVukHMbp01BXR6dOlTY4cQK9e0swoB80bw4Op05sz7JgwXigO2DCYt2cPj2gbVtMmIDKv0FQMozvBdOOHTuuWLGiwsFJkyYNGzaMEJKYmAjg0qVLIl+/FRdm7xESQjIzibY2KSqqppmHB2nYUGqlhyWgDt4jzM/PT09Pj4zkGhuTjAyphPAfPXqQKib5f/5MNDTIx4+ijiLiqz1sGNm2TdQY5MWHD3lWViQmhrx5Q2bPJtraZPbsSsug1w617x4hn2+E79+/v3r1at++fSsc79u374kTJwC0bdvW0NDwRV34k+8rQ0O0bo0zZ6ppFhmJsjKx73pDSUxkZLSFhbOt7RQvL7f9+wv/WwFJCp48we3bVe2yGxODDh3QoIEEY+LHxQVxcVKOQWIUFDB/PgICoK2NoCDcuoUPH2BpieXLISdlyCh+l0YJIQCePat44/fZs2fk68waZWVlVVXplJWSlqqvjvIoKmLfPuzeXWkxGkq+/P77X7m5ce/eHQNcU1MjpB0O1q7FxImo4jdP6tdFeVxdERsr7SAkyMsLpaU4fhwATE0REoLz55GUhObNsXWrpEuxUzXAJxFqa2u3b99+8uTJ8fHxvCOEkKioqD///NPd3R3AmzdvXr9+bWZmJtFIpW3QIJw4gS9fqmnWuzc6d5Za2S2KWVwu4U0oY7OVpX5T/ONHHDiAiRMrbUAITp2SiURobg5lZTx5Iu04JIXFQmAg5s79vt9Tq1Y4eBB792L3btjY4NAhqcZHVYf/ZJldu3YpKCg4OztrampaWFjUq1evX79+5ubmwcHBAJ4+ferj42NrayvZUKVMRwcdOuDEiepbHj+OT58wY4b4Y6LErE2b/ykpuerq/mZufnTw4F+lG8z27ejdGwYGlTZISoKGBn76SYIxVa5Ll7r1pbB/f6ir4/Dh/xx0cEB8PNauxZIlcHTE5ctSCo6qDv/lE5aWlnfv3o2IiEhJScnOzjYxMbGzs+vfvz9vAZmDg4ODg4Nk45QJvKujgwdX06x+fWzahHHjMGECWrSQSGSUGBw7hnv3vFJSXEtLs1u2bFmDbbwYxOFg0ybs319Vm+PH0aePpAKqjosLTp6s6vtr7bNwIaZNg6dnxQXH3bohORmHD2PECPz0E1avhrW1lEKkKiPuSTtSwfisUZ5Pn4iWlqBT8uzsiLk54yFIWd2ZNfroEdHXJ4mJkhyzKpGRpFOnatq0bctYyRvRX+30dKKnJ82i5BJT/rVydq5qTm9xMQkJIQYGZNAg8uIFiY+/NHXq3F279nLkbaJ5nZg1SlVGUxPOzjh2TKDGp0/j1SvMmyfmmCgx+PwZAwZg6VK0bSvtUL5atw7TplXVICsLz55Vtb5QwkxMoKGBBw+kHYdkLVmC+fMr3ZdNWRk+PnjyBHZ2+PnnK25u89av7zR5cvy8ecslGyZV0fdEeODAASMjI96e7zY2NkaVkF6oMkGQuaM8OjpYsQJBQUhLE3NMFKMIwejRcHHBmDHSDuWre/fw9Ck8Patqc+IEevZEdZWkJMrFpW7dJgTg6IhmzfDPP1W10dDA7NkYOPB0UdEfQK/8/LUREQJMPaDE6fvFbDMzs0GDBrVo0QJA37598/LypBeV7OrXD5Mm4c0blNuiqlLTpyMkBH364O5d8UdGMWTpUmRkYM8eacdRzpo1mDKlmiR34gQGDpRUQIJxcUFEBKZOlXYckrV0Kfr3x4gRUFOrqpmzs9XBg8c/f+7GYh1t1cpKUtFR/H1PhO3atWvXrh3v58WLF0spHlmnro6ePXHkCHx8BGp/7hyaNsX69XXu40BOnTuHzZtx/TqUlaUdyldv3+Lo0WqWIhQXIy4OoaGSikkwLi6YOhVcrtQKn0pF27b45ReEhlZzKXvEiKF37z6LjHT88qVVixa0QLGU1aV3KEMEvzoKwNgYc+fC3x/v34szJooJaWn47Tfs2wdjY2mHUk5ICAYOhLZ2VW3i4mBtDR0dScUkGEND6OkhJUXacUjcX39h6VJ8/lxVGxaLtXLl/NTUa7dubd+xo9HNm5IKjuKn0kR47Nixzp07N2rUqHHjxrwjK1asWLt2raQCk13u7rh9G5mZgrYPDISBgUwsc6aqUFiIgQMxbx46d5Z2KOWUlSEkBFOmVNNMRgrK/KhO1Vr7xtoazs7YuFGgxkZGWLUKI0eiqEjMYVGV458Iw8LC+vfvr6qq2q9fv28HDQwMli1bJvX6GlKnooK+fREhTL2tmBjcuIHwcLHFRInM1xdWVvD1lXYc/xURgWbN0KZNNc1OnJChFYTl1c1ECGDxYqxZAwEnWnh7w8oK9H6UFPGvNTpnzhw/P7+YmJhRo0Z9O96pU6fc3NyMjAzJRSerhLo6CsDSEj4+8PFBQYHYYqJEEByM27cREiLtOMo5dSrG2tp13LiePXter7rlvXvgcNC6tWTiEo6rKy5dqovFNi0t0aMH1qwRtP3mzfjnH1yv5n81JS58EmFOTk5mZuboH/ZQMDAwAJCbmyuJuGRbt25ITcXLl0KcsnkztLRQ7gs2JSsSEhAUhMhIqKtLO5Sv3r9/P2rUnHv39n35smXVqgllZWVVND5+HD9sFSMrtLVhYoJbt6QdhzQsXIiNGwWdHKCri40bMWoUCgvFHBbFD59EqKysDKDwh/8hL1++BNBA6lu8yABFRXh64uBB4c46eRKxsYKux6ckIysLXl7YuRMyVUM+LS2NkJ8BfaAp0KTqvz5l9gYhT529Otq0KTw9sVLgCaH9+8PGBn/+Kc6YqErwSYSNGjVq2bLlpk2bCCEsFot3kBCyfPnyxo0b/yQjNX2lzcurmsKPP7K1xeDB8PautPAEJWGlpRg8GFOnomdPaYfyX1ZWVioqycB+RcUdDRu+MzQ0rKzl+/e4exddukgwOCHV2UQIYMECbNsmxJ71Gzdi3z583fWHkhz+k2WCgoL27NnTrVu3Q4cOFRYWrl+/3tnZOTw8fNmyZd9SYx3n5ITcXDx8KNxZe/ZAWRlDh4onJkpIU6dCRwczZ0o7jh+oqqqOHBnVps2zefPeXrkSVcUv3cmTcHGpaodCqevSBVeu1NE//oyMMHw4goIEba+tjS1bMHp0NUsvKOZVVoQ0Ojq6RbmtExo3bhweHs5cKVTxElPR7QqmTyeBgUKfdfYsYbPJhQtiCEj8alPR7V27SIsW5NMnZntlBpdLmjUTqOT3kCEkNJT5AJh9tX/5hVy5wmB/sqXq1yo3l+jokPR0ITr87TcyZYqoUYlVHSq63adPn4cPH7569ermzZtPnjxJT0/39vaWSGqWGzW4Ogqge3f07AlPz+97eFKSl5wMf39ERkJTU9qh8HP6NOrXr77kN4eDmBj06iWRmERQl6+O6upi7FgsXSrEKevWISoKZ8+KLSbqB9VUlmncuLGdnV2zZs3oFdEftW+P4mLcuSP0iUeOoKQE48eLISZKAO/eYeBAbNoEK1kt8bh5s0A1+S5dgpmZbNXB4asuJ0IAs2cjMhLPnwvavkEDbN+OCRMEXYZIiY6WWKs5FkvoBYU8ysr45x/88w+SksQQFlUlDgfDhmHECJkrUf1NejoSEuDlVX1LGZ8v+o2TE65fr7uVUxo2xKRJWLRIiFO6dUP37vD3F1tM1H/RRCgS3tVRQoQ+8ddf4eAgH59itcysWSAE8+dLO47KhYTgt98EWtQoU1vSV0FDA61a1enV4v/7H06eFG5u3Zo1iI3FqVNii4kqhyZCkfz8M9TUcONGTc49dQofP2Lq1PyEhISSujmpTuIiI3HkCPbtg4KCtEOpREkJduzAhHkCEtQAACAASURBVAnVt3z+HJ8+wc5O/DExoY5fHdXUxP/+h4ULhTilXj2EhmL8eHz4ILawqK9oIhTVoEE1uToKoH59DBoUvmGDY+fOCzQ0LJ4+fcp0aNR/3L2LCRMQEVHNTg7SdfgwWrdG8+bVt4yKQu/ekJd793U8EQKYOhXx8bh9W4hTXFzg6Ynp08UWE/UVTYSiGjoU+/fXsJpiZOQiIJrLjSkpWTx+PL0hIBb5+fnbt+/cseOAp2fpqlX45RdpB1SlzZsxaZJALeXlBiGPoyNu3cKXL9KOQ3rq1cPs2QgMFO6s5ctx9SoiI8USEvUNTYSiat4c+vq4fLkm5xJCAN6+4yplZXWvMrH4lZSU2Nm5+fp+8PF5UFg44LffpB1QlR48QGoqPDyqb/n5M27cQNeu4o+JIerq+PlnXL0q7TikauJE3LqFa9eEOEVdHWFhmDwZtMazWNFEyICazR0FsGDBRDa7E5vtBay0slrNdFwUUlJSPnywKSn5H4ezsKys+L1s74+8cSN8fKCkVH3LM2fg4AANDfHHxBx6dVRFBXPnYsEC4c7q2BHDh1ez3z0lIpoIGTBsGPbti9+372CekAt/Zs/2f/487sCBQaGhl3bssKxZNqWqYGBgUFj4ACgGPgI5mrK5fh4A8Pkz9u/HuHECNZav66I8NBECGDMGqam4cEG4s/76C/fuCV3lnxKcorQDqA3++mvW58/Zo0a1MDR0vXPnvFAbdDRp0qRJkyYAnjzB8OEwN4e9vdgCrXsOH25cr95oTU0HFRWl9etXKCrK7ht+9264ugq0Op7LxalT8rdNQceOuHcPeXkyWs1HMpSU8OefmDdPuJspKirYtQu9e8PJCQYGYguuDqPfCBlw7Nj5srJdJSVz3rzxPHfuXM06WbECbm5wdhaiVj1Vte3bsWYNbtwYnZmZ9OLFtT59ekg7oqps3SroNJnEROjoyNa+UYJQVUXbtrhyRdpxSJu3N96/F7qC2i+/YMwYgdbVUDVAEyEDVFUVgFwAKiqP9fT0atzPyZMwM4O1dR0t1c+ssDAsWoS4ODRpIu1QBHDlCj5/houLQI1PnJCPdfQ/oldHASgoIDAQc+cKXYhjwQK8fIndu8UTVt1GEyEDdu7829S0t5KSraOjUefOnUXpKjERXC4cHJgKrY6KiMCcOTh9Wm6+Nm3ejMmTBV0UePy4/N0g5KGJkGfQIJSVISpKuLOUlbFrF/73P7x6JZ6w6jCaCBng4uKUlpY4f36yhYXAO49VQl0dd+7g/n2BSk1SfB09imnTcPas7NbUruDtW5w8iREjBGqcmYn0dHToIOaYxKNdOzx5QkulgMXCwoWYNw9lZcLtQWNjgylTMHZsTco6UlWgiZAx/fvjyBEG3qDGxoiJQUSE0NOsKQBnz2LiRERHo1UraYcisG3bMGAAGjUSqHF0NHr1ggxP+qmKsjI6dMClS9KOQwZ061aQltZXW7uDhUXHR48eCX7inDn4+BE7dogvtLqIJkLGtG4NFRWkpDDQlaMjduzAkiWIiGCgt7rj/HmMGIGoKLmpwAmAy8XWrZg4UdD28rhwojx6dZRn06btRUVd8/JuPH++dfz4OYKfqKiIsDAEBHwOCtoRFhZeVGc39WAUTYRM6tsXR48y09XIkfDzw5AhdKsmQSUkYNgwHDyIdu2kHYowTp2Crm71e/DyFBYiPh49ZHr2azVoIuTJzf1QWsqbx2Xy7p1wF4ubNSsDes6d+2HixNcdO8rnvCkZQxMhk/r1w7FjjPW2ejW6dIGTE62uVL3kZAwYgH/+gbOztEMRkuDFRQHExsLWFg0bijMgMbOzw8uXePtW2nFIm4/PMAODhSoqi1RU+s+ZI/A7AADw+PFjwILL/b2oKCArq/7r16/FFGTdQRMhkxwckJWFFy8Y6/DcORgawtYWZWWM9Vn7pKTA3R3btqFXL2mHIqT0dFy/jsGDBW0v79dFASgqwtERFy9KOw5p++mnn+7ejdm+3UZZOdTDQ+B3AADAwMAAeAwUAHlcbpq2LG+nIidoImQSm43evYWeFV21lBQUFkK0RRm12ePHcHdHcLBcLq3bvBkjRwq0By/PyZNy+TQroFdHeXR0dIYP79e1a7MjR4Q7UVtbe+XKGaamXerVc+vVa7Gampp4AqxDaCJkGLNXRwGoqyMxEbduYdgwJrutHZ49Q7duWLFCiC9VsqOkBP/8Ax8fQdvfuQMFBbRoIc6YJIImwvK8vWuyRn7kSK+0tBu3bl07darPx49iCKuOoYmQYW5uSE5m+BaIhQXOnsXBg1ixgslu5V16Otzc8Oef8vonQkQE2rSBpaWg7Y8fR9++4gxIUn7+GTk5yM6WdhyywcMDKSk1XCNvaYk+fbBmDdMx1T1SSISEkJkzZzZq1Khhw4b+/v5cbsUlpe/fv/fz8+vUqZOFhcW7d+++HS8sLBw+fHiDBg309fWDg4MlG7WgVFTQtStOnGC4W2dnrFuHgADGZqXKu4wMuLrC31+Ib1SyRqhpMqgVNwh52Gx07iz0Dgy1lbIyPD2xb18NTw8MxMaNdD6dqKSQCPfs2RMVFfXw4cMnT56cOHFi9w/XBUpLS7W1tX19fZ8/f84pt/X7smXLsrKyMjMz4+PjFy9efP36dckGLijGr47y+Ppi6lQMHox795jvXL7k5qJ7d/j4wNdX2qHU1IMHSEsT4obfmzd48KD23CqmV0fLGzECYWE1PNfUFMOHY/lyRgOqe6SQCHfu3Onr66uvr6+rqztlypSdO3dWaKCvrz9//vyePXtWOL5jx45Zs2bVq1evefPmw4YN+/FEGeHhgdhYFBQw3/PatejYER061OnZ52/fomtXDB+OWbOkHYoINmzA+PFCFIg5eRLdu0NFRZwxSRBNhOV16oSiIty+XcPT581DeDgtQCoSKSTCJ0+eWFtb835u3br106dPBTmrsLAwIyNDwBMJIV++fPnwVWFhoehhC05LC23bIiZGLJ3HxUFPD7a2+OGKcp3w6RN69UKPHpg7V9qhiODzZxw4gLFjhTil1lwX5WndGp8+IT1d2nHIBhYLw4cjPLyGp+vqYvx4LFnCaEx1jFhKFqakpMTHx1c4yGazfX19AXz48KF+/fq8gxoaGuXvAlbh/fv3AAQ88d69eydPnvT39+f9s379+nfv3lVQUBDyedRcz57KERFsV1exVD+6cgVWVvV0dOaUlUVbW5tFRm5XF3wCvmgKCwuVlJQkv73tgwcP9u3b5+TU4++/3eztOQsWFOfnSzgEJoWGKnXpoqihUSjgsygtRUxM/WXLvuTnS6jW8ufPn8U9RKdOamfOlA0ZUirugcSNkddq4EC2m5v6vHmfa/a75evLsrWt5+NTYGkpiT+QJfD2YAqXy1VSUlJSUqq6mVg+0T5+/Pjs2bMKB7/lIR0dnby8PN7Pnz59EnADPx0dHRaLlZeXx9v/veoTra2t58+fP2TIkJpEz4ShQ7F8OdTUlMSRMjQ04OTkGx3NAU4mJGz29BybkHCc+WH4UVRUlHwijI6O7t9/CiEjgoMXtmmTvXHjaBZLWZIBMO6ff7BuHTQ0NARsf/48WrSAhUV9sUZVgeDh1Uz37rh6VXH8eFWxjiIZor9WNjYwM8ONGxo1q5+noYEZM7BqVb29e0UMRPARxfv2YAqXyy0/0aQyYvlEc3JycnJyquzRFi1a3Llzp1u3bgBSUlKaN28uSJ8qKipNmjS5c+eOiYmJUCdKhbExmjRBQgIqfxlEcutWErAJaAL4paS4imUMmfHHH6u53FDADZjy7JkTizVa2hGJ5NIllJUJVwqull0X5XFxQZCou5bVKrwFhTUuJDt9Oiwtcfs2fv6Z0bDqBincIxw3btyGDRtSU1NfvHixbt26cePG8Y4PHTr0zp07vJ9v3brF+/nOnTtJX8tOjxs3bunSpW/evElMTNy/f/+YMWMkH7zgxDR3lMfNzZ7FWgE8BP4uLh5Yuwtzczj1AV7ZujQ1Nbn/AsFbNSHgHrw8x4/XhoIyFbRogbIyJusRyruhQ3HiBGp80VFNDbNnY/58RmOqO4g0LFmyxNTU1MTEZPHixd8OdunS5fr167yf27VrZ/dV27ZteQeLi4snT55sYGDQrFmznTt3VtH/4MGD9+3bJ7bwBZKSQpo2FWP/Q4aM09a2cXPzcnMrYbPJ3LliHOubgoKC0tJSSYxECCEkK4v070+srDLV1MzZbEslJZPY2FiJjS4OubmkYUPy/r0Qpzx6RIyMCJcrtpj4ycvLk8Aow4aR7dslMI54MfhaeXiQXbtqfnpJCTE3J1euMBVOpSTz9mAEh8MpKSmptpl0EqG4yUIiJIQ0a0bu3JHEQDt2ECUlYm5OcnLEO5AkE+HBg0Rfn8yeTYqLCSHk7du3khlXrP76i4wfL9wpf/9NJk4UTzSVk8wnXWgo8faWwDjixeBrdfAg6d5dpB527iSdOzMUTeVqXyKkJdbEyMNDQoVgRo/+9xKTqSkOHJDEiGKVlgY3N/z1F06dQlAQlJUBQFlZvifIAOByERqKCROEO6tW3iDkcXFBbKy0g5AlHh5ISoIouyqNGIF373D2LHMx1Q00EYqRWG8TVmBsjNRU+Ppi2DC4u8vrtk2EYOtW2NvDwQGJibC1lXZAjDpxAgYGsLMT4pRPn5CUBNdaOh3KwgJKSnjyRNpxyAxVVQwYgP37a96DggIWLcKcOSASWmhTS9BEKEaOjsjIkOh0gNWrcekSrlyBvj6SkyU3LiNevED37ti5ExcvIjAQ1a38kT/CFhcFcPo0nJyE2KdJ7nTpQkvM/MeIEfjnH5F6GDAAiooQdmunOo4mQjHibU8YHS3RQR0ckJMDa2u0bYsFCyQ6dI3xvgi2b4/u3XH5MqyspB2QGKSlISlJ6O2iavF1UR5aa62Czp1RUICUlJr3wGJh4ULMmSOvl4WkgiZC8ZLk1dFvVFVx4QJCQrBsGVq2lPXCpKmpcHXFrl24dAmzZ0OC9X8katMmjBwJVYFXf2RlZQUE/HX48EoHhw/ijEvKXF0RF0ev433HYmHo0JrsUFhejx4wMsKePQzFVAfQRChe3bsjKUk6qWjcOKSmorAQjRsjIkIKAVSLy8XWrejYET17Ij4eMlwgQVTFxQgLw/jxgrYvKirq0MFj+XKzwsJGXl61YhPCSpiYoH59PHgg7ThkyciR2LMHApRDqUpQEAIDUVzMUEy1HU2E4qWmhq5dcfKkdEY3McGLFxg5El5eGDBAtup0P3gABwfs3o0rVzB7Nti1+p146BBsbdGsmaDtHz16VFhoS8gwQsZ++NAoMzNTnNFJGb06WoGlJYyMRH1N2rVDq1bYto2hmGq7Wv3xIxukcnW0vJAQxMXh3Dno64t074EpZWVYvhxdumDoUFy4IER6kF/CTpMxNTUFbgO5QDqLlSZgPV45RRPhj0aMqPlmFN8sW4alS/HlCxMB1XY0EYpdnz44d04s2xMKzskJublo1Qq2tkRb211BoYmKStPIyEiJBZCScrdVKxcjI7sRI+Z37Ii4ONy8CT+/Wv5FkCclBa9eCTfnpVGjRgsX/qWk5PXzz+MPH94i+e0+JMnVFRcvytblCqkbOhRRUTUvt8ZjbQ0nJ6xfz1BMtVod+ByStkaN0LYtzp+Xchi8GTQuLivfvzfjctNKSmKGD/eX2OhDhkx98GBrVtbNvXtfd+kSc/o0TE0lNriUbdqECROEngSUm+s2dWpccvIZB4cO4olLVhgaQkdHJq5VyA5dXXTqxMCVpEWLsGoVPtTm6VbMoIlQEqR+dfSbRo1SABsAgEVRkVbLlvD1RWKiWMYiBA8fYvt2jB6Np08/Az8BLBarraVlXdmP9dGjR0FBG/btix8t/IYZkZEYMEAMMckkenX0RyNGiDp3FECzZujfH6tWMRFQrUYToSTwEqGI08AYMXOmH5u9GFjDYg2wsjKzs8OZM+jYEQoKMDHBiBGiJuzSUiQlITgYgwdDTw9du+LMGfz8M/r27aqpOUFBIURXN7R3714MPRuZdvNmUufOIwMC6hcWrj10aKtQ5z57htxcdOwoptBkDk2EP+rXD4mJyM4WtZ+FCxESgpwcJmKqvWgilIQmTWBigoQEaccB2NvbJyefGDv2/tq1rg8eHA4PR2oqysoQGwt3dyQl4ddf/02KAwZg9+7/3LkpKyubMWPG8uXLK/SZk4PoaPzxBxwd0bAhfvsNDx6gTx/cvInMTBw8CD8/REQE7d/vuX4969atk0ZGRhJ9zlLyzz+Rb98uBkaVloaFhOwT6tzDhzFgQJ24gcrj4oJLl2TiL0XZoaqKvn2xT7g3Dh9GRhgxAsuWMRGTHLp8+XJWVlb17SRQ/1vyZGT3ifICA8nvv0s7CMHcukVmzCBWVkRFhbDZRF+f9OpFdu4sZLEMAU+gk6qqWWoqCQsjPj6kZUuiqUm6dSMLFpCYGPLli7iikqOC9zzBwVuUlQMBwmKd6d37N6HOtbcn58+LKS6BSP7Vbt2aJCZKeExmiO+1io0ltrYM9PPmDdHWJs+fM9AVj7z8Mrq6DmKz7dasSa62ZZ35m1Pa+vfH4cPSDkIwtrZYvRoPHqCoCJcuoX9/pKZi7NilhDgCkcDloiLdbt3y4+Jgb49Dh/DxI2JiEBiIbt1qc1VMYU2cOEZdPatRo7YdO24MDRViL/bXr/HiBZycxBeaLKI7UfyoSxd8+IC7d0XtR0cHvr746y8mYpIrFy5c53KvFxS0qbYlTYQSYmMDNpuB97SEOThgyxY8foxZs0qBjwCAMiDnxo3inTsxbhxathRus/W648QJJTOzLW/f3rxy5ZihoaHgJ0ZEoG9f1OoVE3zQ24Q/YrEwZAgzldJ+/x3R0Xj4kIGu5EVREQhRA0oEueROE6Hk9Osnoe0JxWHZsmWKio+BVoClurq1jo6OtCOSaYRg6VIEBtbkr4TISAwcKIaYZJuzMxISUFoq7ThkzG+/YfduBu6eNmiA339HYCADIcmFixehrw9lZX8Wq6WGxtVq29NEKDmys4iiZkpL02JjNyQnnzE2jpbN4qWyIyoKZWXw8BD6xJz/t3fncVHV6x/APzMDA8MqgzkwsYi8BH8IiqgoKgICRnINKdNMDFGoNJUkcU2RyCy6SobhkqW45Yq7loCyiAy4ggJiKSaYLCKCgGwz5/fHuZdLhiwzB84wfN9/cc6c85wHXuLDOef7/T4luHVLZRsQtkEohKUlrl5lOw8l83//B5EIyckMhFq0COnpuH6dgVDKTCbDrFmYMAHu7qitDXz2LMvf36bds0gh7D7jxqGwEA978iS60aNH29paxMZi4UIGBnarsHXrEBYmz+3gsWPw9u5EkwqVQVGUSHRw3rzgffsOsZ2LcmFkQiEATU2sWIHVqxkIpbSys2FsjGPHcO4c4uLA5UJHR0dHR6fdE0kh7D48Hry9cfIk23kozMkJc+ciIIB0z2nd6dN48QJvydU0olfNo29p27ZdKSlnsrKmffLJyW3bYtlOR4m8/z5OnGBmjcbAQOTnIymJgVBKKCQEw4ZhwACUlmLixM6dSwpht+rpT0ebrV2Lp0/J2vati4jA2rXyzAJ89gyZmXjjjS7ISekdORJfW/s5MLaycs2RI+fZTkeJ9OsHR0dm/oBWV0dYGD7/nIFQSqWwEDY2iIlBbCzS0+UZu04KYbeaOBFXr6rC0n9qaoiNxcqVuHuX7VSUzNmzqKmBr6885x4/Dg8PaGsznVNP4OQ0VCD4Bajm8/ePHWvPdjrKhZFmFLSZM1FdjXPnmImmDCIjMWAAADx8CD8/OYOQQtitBAK4uuLMGbbzYMKgQVizBjNnkpF+f/PVVwgLk3NRmKNHe+N4Udrq1YuDgpr69p08dKhs1apP2U5Hufj6QiJh5q08l4sVK2pnz/7yjTc+OH78FAMR2fPsGZycsHIlPv8cublQpFkZKYTdTWWejgJYsACvvYZ/rLnWe/32Gyoq5Cxmz58jNbVz3ZpUCZ/P37Qp4vDhixzOF+rq6myno1wEAkyejIMHmYl29GhIWZnu+fNL58zZnJ7e/tQC5XTkCIyNUVSE/HyEhSkajRTC7jZ5MhISUFfHdh5M4HCwYwc2b0ZmJtupKIcvv8SaNXLeDp4+jXHjoKfHdE49ypgxuHMHT5+ynYfyYfDpqERyjaKCAduKirnnz6cyE7Qb1dVh0iRMn47Zs1FYCEtLBmKSQtjdDA1hb4+EBLbzYIhYjJgY+Puz3HlYGcTHo7wc774r5+m9cx79S/h8jBmjssMaFeHmhtJS3L7NQCgHBzs1tR1AgUCwx82th3W7TEmBkREyMyGRYMsWxsKSQsgCVXo6CuDttzFiBJYtYzsPtkVE4PPP5bwdfPEC8fHyTMBXPe7u7HexVkJcLmbMYKAZBYDdu78LCrrn6LhUTc3P1LQHrGk7d+6namqmamqmVlbhbm4YMwbFxRg5kslLkELIAl9fnDypUk1nNm/GqVM4e5btPNiTmIjiYkybJufpv/6KkSNB1q0D4OGhOs9LmOXv/3JnNPno6enFxKzPyDi8bNn0kBAmMutKjx492rXrtFT6u1T6+++/n/z554dnzzK/Ei8phCwwN4dYjB77lroV+vrYswcffYTycrZTYUlEBNaskf/3k25ASAAYOhSVlXjwgO08lI+NDQwMkMrcS73QUOTl4fRpxgJ2hYKCAsAC0AQ0uVyRWJzfFVchhZAdKvZ0FICzM6ZPR1AQ23mwISkJRUV47z05T29sxLlz8PFhNKcei8MhLZleicEhMwD4fERH49NPlXrsnpXVGIqqA1YBKzU0ct3c3LriKqQQsmPKFBw7xnYSTFu3DvfuMfmL2lOEh2P1avlvBxMSYGMDsZjRnHoy8prwVd5/H3FxePGCsYCenhg2THlnQNXWwtaW279//OLF9Z991lBaeluta/qT9bKmZ0rD3h4yGXJyMHgw26kwR0MD+/fD3R3jx8PcnO1suktaGgoLMXOm/BF68zz6Vnl6YtUqUBRpdfkyY2OMHIlTp+R/G/1PUVGwt8eMGbCyYiwmI2QyDBkCHg+5uZqamv/u0muRO0LW9Oj2hK8yeDBCQjBrlkoNBWpbWBhWrZL/dlAqxalT5Lno35ibQ0+v53Wx7h7MPh0FYGKC0FAsWsRkTEaMHInycuTldUczFlIIWaN6rwlpS5aAx0NUFNt5dIv0dPzxh0K3g8nJMDODhQVzOakEMnb0VXx9cekSSkqYjLl4MR4+VK7GOO7uyM3FzZvo06c7LkcKIWucnfHgAQoL2c6DaVwu9uzBhg3IzmY7la5Hr+XP58sfgcyjbxV5Tfgq2tqYPBmHGG3ayOdj61YsWoSaGibDym3WLKSmIiOj+96wkELIGh4PkyYp119hTDExwddf4/33lXo0muIkEuTn44MP5I9AUThxQs5WFarN3R2XLqGhge08lJKfHzOtelsaPx5OTvj6a4bDymHlSuzfjzNnMGRI912UFEI2qerTUQD+/hg8GGvWsJ1HV1q7FqtWKXQ7mJ4OAwNYWzOXk6owMICVFTIy2M5DKXl44K+/cOcOw2GjorBtG/K7ZJ5eJ3L45hvs3QtPz269LimEbHrjDWRmqkJ7wlb98AP278fFi2zn0TWuXUNeHmbPVigImUffBg8P8nS0dVwupk/Hvn0MhzUywrJlWLiQ4bAdd+QIlizB999jxozuvjQphGwqK/sT8LSwGD5lSkCjyrX169sXO3fC3181K/2aNVixQqHbQQDHj5MXhK/k7k7Gy7wSn79z/XoHE5Phe/ceZjBscDBKShAXx2DIjrp4Ee+9h88+wyefsHB1UgjZNHfussrK8MrKa+fPW27d+jPb6TDP0xOTJ+NTlWuzev06bt1CQICiQXg82NkxlJPKGTcO2dmoqmI7D+VTVla2Y8ePUqnk0aNLISHra5nr/KKmhh9+wOLFqK5mKmSH3LiBiRPh54fIyG69bjNSCNlUVPQYGArgxYth9+49YjudLvHtt8jIYHiQG+vCwrBsGTQ0FApC5tG3TVMTo0YhOZntPJRPWVkZh2MJ8AEBh2NSwegjl3Hj4OyMr75iMGQ7CgsxdizeeAO7dnXfRV9CCiGb5s/3Ewrf43Bi9PXDAgOZWytCmWhpYd8+LFyIwkLU19eznQ4DbtzA9euYM0fROHFx5AVhO8gkilZZW1uLREVaWmE83gpNTdnrr7/ObPwNG7BjB/LymI3auqdPYWeHIUNYXvubFEI2LVoUdPbs51Onavr4HLG1tWU7na4yfDi8vJIGDrQ3NXX19p4l7eGrzoSHY/lyCAQKBcnJQXU1RoxgKCcVRabVt4rH41258mtsrN22baMaGo5nZTEcXyTCqlXdMWqmoQF2dhAKcflyl1+rbaQQsmzUqFEREXMuXuxPUWyn0pVSU1fU1yeUlaWnphqdVvK+L226eROZmQgMVDROXBymTiVrabbDwQElJXikmi8NFKKhoTF16tS5c6dERqp1xYTdBQtQXt61bzRkMtjZQSpFbq6c7awZxPb1CcDaGgIBbt5kO4+u1NDQBPQBUFcnfvr0GdvpyO+LL7B0qaK3gyAvCDuGy4WrK2nJ1JZZs2Bry/yEXR4Pmzdj8eIuHKw0ahT++gvZ2d2xlGi7SCFUCm+9hVOn2E6iK82bN9PQ8C0NjRV8/v4pU95iOx055eRAIsGHHyoap6AAJSVwcmIiJ1VHXhO2i56wm5TEcNixY+HpiS+/ZDgszcsLt28jOxv9+nVJ/M5SzUL49OnTyspKtrPohMmTe0YhvHXrVqFcq6OuWvVpSsq/jxxxFYuTL182kCPC+fPnZTKZHCcyKCwMS5ZAS0vROIcPw9cXPB4TOXWNxMREJZnY6uGB+Hi2k3i12traZLYHtvbti59/xpw5zN+9RUZi585WFg2+cOFCgwLL3338MS5cQHp6d6w1X15eLpFI2j+OUkWmpqbBwcFsZ9EJTU3UyF5+/wAAD+lJREFUa69RhYVs59GeoKCgH374QZEIv/1GWVpSL150+sR+/foVFxcrcmkF5eRQ/fpRz58zEGrUKCo+noE4XWfAgAG///4721n8h4UFlZvLdhKvkJGRMXz4cLazoCiK+ugjKiCA+bCbN1POzpRM9redVlZWeXl5nQ1lZzeByzXjci05nHPnzjGWYduOHTv21ltvtXuYat4RAqB61OATHg9eXiwPIO4gBX+wEyfCwUEp1vbtrPBwhIZCR0fROI8e4d49uLoykFIvQZaY6YgNG5CWhiNHGA47bx7q6/HLL4rG2bFjx+3b2jLZA5ksk8OZ7+XFRHId0MH/r1S2EPY4PeXpqOI2bUJMDMtr+3ZWXh6SkvDxxwyEOnoUkyfL38i3FyKvCTtCWxu7dmHhQhQXMxmWy8XmzQgNhYLvmg4erKAoe4ADCAGK9dccLyGFUFm8+SbS0vD8Odt5dD1jYyxfzubavnL44guEhDBwOwiy0HbneXggKQnK8cpSqTk5Ye5cBgZzvWTkSLz5JsLD5TlXJsPq1dDVRXLyhxzOQQ4njMN519LydS7rEyb+jtOzHiF2kKGhoY6OjpWVFduJdE5e3mKx+Ly+fg7bibxSXl6erq6uiYmJgnEoipeVtXbAgL16eh29MUxJSRk9ejRfwVWuO6+pqammhv/HH5H29qt5vBcKRpNKNW7c+NrBYRmXq9St9i5fvjxs2DCB4tNEGJKdvdrScpe2ttK1sa6qqrpz546joyPbifwHRallZ68xNz/cp88tBsM2NupnZYUNGfIFn/8MQHp6+tChQ7XaGzbW0KB35Up0U5OeWHzO0nJHU1PdgwcPBAKBqakpg7m1raysDMDN9manqWYhPHjwIJfLNTCQZ3Qi0YaysjKBQKDDyJ1RJxUUFFh0wyAzAgDw4MEDc3NzDpnw3x6pVPro0SMzMzO2E+lWPeifR319vZaWlpubW9uHqWYhJAiCIIgOUq4HtQRBEATRzUghJAiCIHo1UggJgiCIXo0UQoIgCKJXU7VCeOXKleQWCgoK2M5IFZSWlib9fU3flJSUYoYm7iYkJNTU1DRv5ufn5/2jJWhlZeW8efP+ee6tW7fu37/fvHnt2rWrV68yklXvJJVKk//uyZMnbCelpB4/ftzuoHyVceXKlUf/7YYlk8kSEhJKSkrozYaGhoSEBGXruS2VShMSElquOJ2Tk5OT8+qZaV26zlv3s7KysrOzm/BfP/74I9sZqYKampqBAwfu3buX3jxw4ED//v2rqqoUj0wvMHH79u3mPcHBwYGBgS8dVlxcLBaL/3n6u+++GxYWRn8dExNjbGx88+ZNxbPqtaqqqgCMGTOm+TcoOTmZ7aSU1LZt28aMGcN2Ft3E398/JCSE/vrGjRscDic8PJzeTEpKEgqFUqmUvexaFxAQ4OfnR39dVFRkaGh46dKlVx2sggs9hYWFvUNavTFKS0srNjb2X//6l4uLi4aGRnBw8J49e3R1ddnO63+++eabLVu2JCUl9bhVFJTQ/v37zc3N2c6CUCKurq6bN2+mv05KSpo0aVJzz42kpCQXFxdlWykGQFRUlJ2dXVxcnK+vb2BgYFBQ0NixY191sAoWQqIrODk5ffDBB/Pnz1dTU3vnnXc8PT3Zzuh/wsLCDh48mJqa2p0rVhBE7zFhwoTAwMBnz5716dMnOTl50aJFs2fPrq+v19DQSE5O9vHxYTvBVujr62/fvn3OnDl37959+PDh8ePH2ziYFEKio9atWzdo0CAej7d79262c/mfLVu26OnppaWl9VOSFp8EoXLMzMxMTU0vXbrk7e2dkZGxd+9eBweHzMxMR0dHiUTy3XffsZ1g67y8vNzc3FatWpWenq6hodHGkUp3P0sorby8vKqqqufPn1cx3gBUAba2tiUlJampqWwnQhCqzNXVNTk5OTs728LCQltb29nZOTk5OSMjQyAQ2Nrasp1d62pqajIyMrS1tYuKito+khRCokPq6+sDAgK+/fbbWbNmfcjc+vYcDkdHR+d5i6YbVVVVnXr76OzsfPTo0dmzZx84cICprAiCeImLiws9kHj8+PEAxo8fT2+6uroq4QtC2vLlywcPHhwXFzdv3jx69e1XUdJvgFA24eHhr732WmBg4Lp16+7evbtv3z6mItvY2EgkEvpriqIkEkln/8D09PSMi4sLCgoitZAgusiECROuX79+4sQJFxcXACNGjMjKyjp//ryrsvaYTktLO3To0Pbt2z08PLy9vUNCQto4mLwjJNp348aNLVu2XLt2jcPhaGlp7dq1a8qUKR4eHiKRSPHgERERM2bMoJfwP3v2rIaGxsyZMzsbhK6FU6dO1dDQ8PX1VTwrgmjX/fv3m6e3qqmpRUdHs5tPlzIzMzMzM0tJSTl27BgAdXX1oUOHJiQkxMTEsJ1aK2pqambPnr1p0yb6/6ioqChbW9vjx49PmTKl1eN5a9eu7dYEu5iRkdGIESP69OnDdiIqJTc318/Pz97ent40NTW1srLicrnGxsaKB7e0tHznnXceP3785MkTV1fX77777p+vtWtqarZt2/bZZ5+9tN/AwGDo0KFisZiO4+bmVlRUZGtry+PxFE+sF+Jyuaampo6Ojt3f97HH0dbWtrS0NP4vsVg8bNgwtpPqWoMGDfLy8hoxYgS9OXDgwFGjRnl5eSlhP6aCggJra+sZM2bQm5qami4uLjU1NdbW1q0eT9owET1ASUmJg4ND89oWBEEQDCLvCAmCIIhejdwREj1AU1NTfn7+4MGD2U6EIAgVRAohQRAE0auRR6MEQRBEr0YKIUEQBNGrkUJIEARB9GqkEBIEQRC9GimEBEF0ws6dO69cucJ2FgTBJFIICYLohODg4LZbuxFEj0MKIUGwgKKo0tLSxsZG+U4vKyt7+vSp4mlIpdLi4uKW3T9aopPseNet0tLSioqKtg+orq7udJYE0cVIISQIBkRHRxsZGdXV1dGbK1asEAqFe/bsoTdv3bolFAovXLgAIC8vb+LEiQKBQCQSaWlpDR8+PC0tjT4sPj5eKBSmpKS0jBwZGSkSiZrL3tatW83MzPr162doaGhra5uUlNRqPlu2bDE0NHypDdtHH31kY2Mjk8kASKXSsLAwkUhkbGysr68/bty43Nzc5iOlUmlERISRkZFIJNLX1zcxMTlw4EBDQ4NQKKyurt64caNQKBQKhc0rFcfGxpqZmYlEIqFQOGTIEPo7pS1YsGDkyJGnTp3q37+/SCRasmSJXD9gguhKFEEQCqNfmyUkJNCbQ4YM4fP5M2fOpDejoqL4fP7z588pikpNTV2yZEliYmJeXl58fPzYsWP19fVLSkooimpsbDQyMgoICGgZ2dra2sfHh/46MjKSy+UuX7786tWrEonEx8dHIBDk5ub+M5/i4mI1NbX169c376murtbR0QkNDaU3g4KCtLS0NmzYkJWVdfHixVGjRhkZGZWXl9Offvjhh1wuNyQkRCKRXL9+ffv27bt27ZJKpfHx8QKBYMaMGfHx8fHx8Xfv3qUo6pdffgEwbdo0iURy4cKF0aNH8/n869ev06ECAgL09fXNzMx++umny5cvSyQSZn7iBMEcUggJggFNTU1CoXDlypUURZWUlHA4nPnz54tEIplMRlGUt7e3s7NzqyeWl5fzeLyffvqJ3gwJCaE7FdOb9M1iXFwcRVGVlZU6OjoLFixoPreurs7c3HzevHmtRvb29raysmre3LlzJ4Ds7GyKonJycjgcTnR0dPOnxcXFAoEgKiqKoqjc3FwOhxMSEtJqWF1dXfrbbGZjY2NjYyOVSpu/I11d3alTp9KbAQEBABITE1uNRhDKgPQjJAgG8Hg8FxeXhISEdevWXbhwQV9ff8mSJTExMbm5udbW1qmpqS37gpaWlh48eLCgoKCmpgaApqbmH3/8QX8UEBCwcePG48eP+/n5AYiNjTU0NJw0aRKAtLS06upqExOThISE5lDm5ua3b99uNSV/f/9p06ZlZmY6OjrSoUaOHGlnZwfg/PnzFEUZGBi0DCUWi+lQdNGaO3duR77xurq6O3furFmzprlNuVAonDhxYssHvLq6uhMmTOhINIJgBSmEBMEMd3f34ODgioqKxMREd3d3CwuLgQMHJiYmVlZWVlVVubu704edOXNm6tSp5ubmzs7OBgYGXC6Xx+M1D0ixtbV1cHCIjY318/Orq6s7dOiQv78/3aCxpKQEwPr165tLDs3CwqLVfHx8fPr27RsbG+vo6Pjnn3+mpKR8//339Ed0qIULF750Cj1q5smTJwBMTEw68l0XFhbKZLKXOlOKxeLy8vLmTUYaOBNE1yGFkCCY4e7uLpVKk5OTExMTQ0ND6T10IdTW1qZvywB89dVXjo6OFy5coLsHy2SyTZs2tYzj7++/ePHihw8fXr58+dmzZ/7+/vR+fX19AMeOHXNzc+tIPnw+f/r06fv379+4cWNsbKyamtr06dNbhrp9+zbd0/gldF/rkpISPT29dq+ira0NoKysrOXOsrIy+hK0lyo3QSgb8g+UIJgxaNAgExOTH3/8saCggL7/c3d3v3jx4m+//TZ+/Pjmnu8FBQX29vZ0FQSQmJjYPNaUNnPmTHV19b1798bGxtra2jb3PXdyclJXVz98+HDHU/L396+oqDh58uTu3bvpG0R6v4uLC4BXhRo/fjyAQ4cOtfqpjo7OixcvmjfFYrGZmdmZM2ea99TW1iYmJo4ePbrjeRIEy9h+SUkQqmPWrFkATE1N6c0nT57QN0Pffvtt8zFvvvmmsbHx1atX6+rq4uPjBwwYoKmp+cknn7SM8/bbb5uYmPB4vA0bNrTcHxoayuVyV69eXVBQUFtbe+fOne+//37nzp1tpGRra9u/f38Ap0+fbrnfx8dHW1s7Ojr60aNH1dXV2dnZERERv/76K/2pr6+vQCCIjo5+/PhxRUVFfHz8iRMn6I88PDwGDhx47ty5q1evFhUVURQVHR0NYOnSpcXFxffv3/f19eVwOBcvXqSPDwgIaDlmhyCUECmEBMGYXbt2AWg5/8HBwQHAjRs3mvfk5+cPGjSI/jNUT09v9+7dZmZmLxVCeukWNTW1x48ft9xPT+9r+cTS3Nz80KFDbaQUGRkJQCQSNTY2ttxfW1s7f/58+u0jbfDgwZcuXaI/rampmTNnjpraf16dCASCrVu30h9lZWWNGzeOfiK6dOlSiqJkMll4eLhAIKAPNjQ03LNnT/OFSCEklB9pzEsQ3a2pqenevXu1tbWDBg1qrh8d19jYmJeXV19fLxaLX3/9dUUyqa2tzc/PpyjKxMSkX79+L31aWVmZn5+vpaXVv39/HR2dtkPV1NTk5uby+XwbGxt1dXVFsiKIbkYKIUEQBNGrkcEyBEEQRK9GCiFBEATRq5FCSBAEQfRqpBASBEEQvRophARBEESvRgohQRAE0av9P6aLe7iFpxhRAAAAAElFTkSuQmCC", "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "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 [-4.387405772893939e-16, -8.999427264780328e-17, 2.551583324994531e-16]\n [-1.4675232268360493e-16, -3.6433503490690063e-16, -5.359761785968941e-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.2" }, "kernelspec": { "name": "julia-1.8", "display_name": "Julia 1.8.2", "language": "julia" } }, "nbformat": 4 }