{ "cells": [ { "cell_type": "markdown", "source": [ "# Tutorial" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "DFTK is a Julia package for playing with plane-wave\n", "density-functional theory algorithms. In its basic formulation it\n", "solves periodic Kohn-Sham equations.\n", "\n", "This document provides an overview of the structure of the code\n", "and how to access basic information about calculations.\n", "Basic familiarity with the concepts of plane-wave density functional theory\n", "is assumed throughout. Feel free to take a look at the\n", "[density-functional theory](https://juliamolsim.github.io/DFTK.jl/dev/#density-functional-theory)\n", "chapter for some introductory lectures and introductory material on the topic.\n", "\n", "!!! note \"Convergence parameters in the documentation\"\n", " We use rough parameters in order to be able\n", " to automatically generate this documentation very quickly.\n", " Therefore results are far from converged.\n", " Tighter thresholds and larger grids should be used for more realistic results.\n", "\n", "For our discussion we will use the classic example of\n", "computing the LDA ground state of the\n", "[silicon crystal](https://www.materialsproject.org/materials/mp-149).\n", "Performing such a calculation roughly proceeds in three steps." ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "3Γ—3 Array{Unitful.Quantity{Float64,𝐋,Unitful.FreeUnits{(β„«,),𝐋,nothing}},2}:\n 0.0 β„« 2.7155 β„« 2.7155 β„«\n 2.7155 β„« 0.0 β„« 2.7155 β„«\n 2.7155 β„« 2.7155 β„« 0.0 β„«" }, "metadata": {}, "execution_count": 1 } ], "cell_type": "code", "source": [ "using DFTK\n", "using Plots\n", "using Unitful\n", "using UnitfulAtomic\n", "\n", "# 1. Define lattice and atomic positions\n", "a = 5.431u\"angstrom\" # Silicon lattice constant\n", "lattice = a / 2 * [[0 1 1.]; # Silicon lattice vectors\n", " [1 0 1.]; # specified column by column\n", " [1 1 0.]]" ], "metadata": {}, "execution_count": 1 }, { "cell_type": "markdown", "source": [ "By default, all numbers passed as arguments are assumed to be in atomic\n", "units. Quantities such as temperature, energy cutoffs, lattice vectors, and\n", "the k-point grid spacing can optionally be annotated with Unitful units,\n", "which are automatically converted to the atomic units used internally. For\n", "more details, see the [Unitful package\n", "documentation](https://painterqubits.github.io/Unitful.jl/stable/) and the\n", "[UnitfulAtomic.jl package](https://github.com/sostock/UnitfulAtomic.jl)." ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n Energy Eβ‚™-Eₙ₋₁ ρout-ρin Diag\n", "--- --------------- --------- -------- ----\n", " 1 -7.905198671336 NaN 1.95e-01 4.1 \n", " 2 -7.909818255786 -4.62e-03 2.98e-02 1.0 \n", " 3 -7.910018961017 -2.01e-04 3.00e-03 1.2 \n", " 4 -7.910051946700 -3.30e-05 1.34e-03 2.4 \n", " 5 -7.910052589627 -6.43e-07 7.51e-04 1.0 \n", " 6 -7.910052843503 -2.54e-07 2.03e-05 1.0 \n", " 7 -7.910052853933 -1.04e-08 1.53e-05 3.0 \n", " 8 -7.910052854037 -1.04e-10 3.90e-06 1.0 \n" ] } ], "cell_type": "code", "source": [ "# Load HGH pseudopotential for Silicon\n", "Si = ElementPsp(:Si, psp=load_psp(\"hgh/lda/Si-q4\"))\n", "\n", "# Specify type and positions of atoms\n", "atoms = [Si => [ones(3)/8, -ones(3)/8]]\n", "\n", "# 2. Select model and basis\n", "model = model_LDA(lattice, atoms)\n", "kgrid = [4, 4, 4] # k-point grid (Regular Monkhorst-Pack grid)\n", "Ecut = 7 # kinetic energy cutoff\n", "# Ecut = 190.5u\"eV\" # Could also use eV or other energy-compatible units\n", "basis = PlaneWaveBasis(model, Ecut; kgrid=kgrid)\n", "\n", "# 3. Run the SCF procedure to obtain the ground state\n", "scfres = self_consistent_field(basis, tol=1e-8);" ], "metadata": {}, "execution_count": 2 }, { "cell_type": "markdown", "source": [ "That's it! Now you can get various quantities from the result of the SCF.\n", "For instance, the different components of the energy:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "Energy breakdown:\n Kinetic 3.0795105 \n AtomicLocal -2.1806269\n AtomicNonlocal 1.7339358 \n Ewald -8.3979253\n PspCorrection -0.2946254\n Hartree 0.5417535 \n Xc -2.3920751\n\n total -7.910052854037\n" }, "metadata": {}, "execution_count": 3 } ], "cell_type": "code", "source": [ "scfres.energies" ], "metadata": {}, "execution_count": 3 }, { "cell_type": "markdown", "source": [ "Eigenvalues:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "7Γ—10 Array{Float64,2}:\n -0.170181 -0.1318 -0.0883272 … -0.0562598 -0.114939 -0.0700161\n 0.201344 0.0909045 0.0122929 0.0111153 0.0420625 0.0176509\n 0.249296 0.174774 0.176138 0.132965 0.220115 0.11233\n 0.249296 0.231429 0.202371 0.161044 0.220115 0.190456\n 0.350986 0.360028 0.340134 0.291807 0.320727 0.327376\n 0.369971 0.395898 0.389483 … 0.331814 0.388191 0.460273\n 0.369971 0.401677 0.412475 0.56556 0.388191 0.462721" }, "metadata": {}, "execution_count": 4 } ], "cell_type": "code", "source": [ "hcat(scfres.eigenvalues...)" ], "metadata": {}, "execution_count": 4 }, { "cell_type": "markdown", "source": [ "`eigenvalues` is an array (indexed by kpoints) of arrays (indexed by\n", "eigenvalue number). The \"splatting\" operation `...` calls `hcat`\n", "with all the inner arrays as arguments, which collects them into a\n", "matrix.\n", "\n", "The resulting matrix is 7 (number of computed eigenvalues) by 8\n", "(number of kpoints). There are 7 eigenvalues per kpoint because\n", "there are 4 occupied states in the system (4 valence electrons per\n", "silicon atom, two atoms per unit cell, and paired spins), and the\n", "eigensolver gives itself some breathing room by computing some extra\n", "states (see `n_ep_extra` argument to `self_consistent_field`).\n", "\n", "We can check the occupations:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "7Γ—10 Array{Float64,2}:\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0" }, "metadata": {}, "execution_count": 5 } ], "cell_type": "code", "source": [ "hcat(scfres.occupation...)" ], "metadata": {}, "execution_count": 5 }, { "cell_type": "markdown", "source": [ "And density:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=1}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdd2BUZfY38PPcmfRCeiaZBAIBQg9ISwIEJCERVwQRRAFhwbUtuor+EPsqNtBFUXBdFl/rrgWwr0IgICSQ0CEgJQQIJb2H9EzmPu8fw2YjEJiEmbnt+/mLzAzDyeWee+59KuOcEwAAgFYJUgcAAAAgJRRCAADQNBRCAADQNBRCAADQNBRCAADQNBRCAADQNBRCAADQNBRCAADQNBRCAADQNBRCAADQNFsWwuzs7MbGRus/bzabbfivA46nDXHOsfqgDYmiKHUIqoJkty1bFsK77rorOzvb+s/X19fb8F8HHE8bamlpaW5uljoK9cDJaVs4nraFplEAANA0FEIAANA0FEIAANA0FEIAANA0FEIAANA0FEIAANA0FEIAANA0FEIABTt79mx5ebnUUQAoGwohgFJ9/9N/oif9sXdMwpkzZ6SOBUDBUAgBFIYTZZbwxzLNc7/Pqw3sV+0aPH5t6UsHzNnVWBMOoDP0UgcAANY6WsnX5YpfnOZ6Rnf1YBnL5mf+9KUh6LbucSPW5YoTN5rddDS9B7u3pxDpzaQOFkAxUAgB5O5cLf/+LP/slFjaQHdEsI/G6EYbLHVO13/eXMtn+vvqXhxCGcV8Xa4Y91NLuAe7t6cwI1IwuEkYOIAyoBACyFR+HV+fy9fliieq+J3dhXdjdKMM7BoPegKj0QY22qB7O0b3awH/LEd86YCpvy+b3l2Y1VMIcHVc5ADKov/uu+8GDRoUGRlp+bm5uTkzM/PixYuDBw8ODw8novz8/KKiota/MGTIEEFAzyKAvVQ20U/nxXW5YmYxvzVcWBwt3BImOHUk53SMEo0s0ahrNOs254vrzvAlB01xwWx6d2FqhODpZLfQAZRJv3bt2oceeujVV1+9//77i4qKRo8e3bVr17CwsHnz5r366qsPPfTQe++999lnnxmNRstfSEtLc3d3lzZoAPVpaKH/XBA/yxHTCnl8CJveXfh6vOB+Y002rjqa1FWY1JWqm3U/nBPX5YqPZpjjQ9icXsLkboIzbmgBiIhI/+WXX+7cuXPy5MmzZ8/+7LPPwsPDt27dSkTffffdI4888tBDDxHRH//4xzfeeEPqUAFUqMlMm/LFdWf4j+fF4YHs3p7CFzcLXrZ+aOviTHN6CXN6CRVN9J/z4j9PiA/tMP8hXJjeg00ME/SoiKBteiIaNWqUh4dHenq6p6cn+28fhCAIXl5elj9XVFRkZGRYnhQlixRARcycMov556fEdbliPx82vbuwPMYp0P7deH4ulypiXh3/JpcvyxIf3ineGcGmdxeu3QEJoGKXWl7Cw8Pz8vLmz5+/d+/epKQko9H422+/ffzxx5Z309LSTpw4cejQoYkTJ37++edOTle/XzWZTD/99NOhQ4favjhs2LA+ffpc9fNms9lsNtvud9E6HE8bMpvNoija8Hi2tLT8/Z//z8XFeejtc/99hq3L5X4ubHYkHZ4ihLgTESdy6P9eiCs90pce6cvO1vC1ufz+Heb6FprclaYZaret/bBvrx5Tp0y24T+Hk9O2cDytJwgCu9493qVC6OTk1NTUdOrUqR07dsycOdNoNB4/fnzDhg2xsbEvvfTSsmXLiKisrCw2Nnb16tWPPPLIVb/LZDLt3LnTx8en7Yt+fn7du3e/6uebm5ubmpo6/GtBO3A8bchkMomiaMMvXLdu3dPfHTI11nU9GzR/avKWCeYIz0tvSfufFuJMj0XRY1GUVcnWnmW3Pf1ubTP3+mJl756RvXr1stW/0tTU1N4NNHQCkt16rq6u1hbC0tLS0NDQ119/fcqUKS+//DIR3Xrrrd26dVuwYEFwcLDlMwEBAVOnTt2zZ0973+Xu7r506dLo6Ggr4zObzRh3Y0M4njZkKYQuLi62+kJjz34tuSt9nMRfZvfs189mX2tDse4Ua6Th+f0fXPJeQ+1FY1i4DU8nURRxctoQkt22BCIqLCw8derUiBEjWlpa9PpLpVGv14uieNlN8bFjx0JCQiQIE0DhvqWb7l+98ey+bf369ZM6lmuZede0nNS1I9/e8X2Zz/U/DaAK+k8++WTNmjX33ntvSEjInDlzZs6c6efnFxISsnr16uTk5JCQkGnTpo0cOdLb23vbtm0ZGRmrVq2SOmYAhTlayb85Kx6bFuAtx0fByxkMhvfG8ls3tkzrLvg4Sx0NgP0Ju3btmjt37gcffEBEt912W2pqamVl5Z49e+bOnfvDDz8Q0d13311SUnL48OGhQ4dmZ2d369ZN6pgBFOaRDPOSoTp/JVRBiyH+7PZuwpIDGI4BmsA4t9mK9dHR0Z999pn1fYQ1NTWtMzTgxuF42pAN+wi/PiMuzRL3TdHrFDU7oaKJ+q43bblVP8DXBnHX1tZ6enpe/3NgHSS7bWEmLYAdNbTQ03vFFTE6ZVVBIvJzoWejdY9n4qEQ1A+FEMCOlmaZY4PY2BCllUEiIlrQTyhppB/P2XIaCYAMoRAC2MuFOv734+Ibw5WaZXqBVsToFu4SG/FYCKqm1BQFkL+Fu8TH+uu6eSrycdBifCgb5Mfe+Q0PhaBmKIQAdvFrId9byp8YqPgUWx4j/O2w+XytzUbVAciN4rMUQIbMnB7PNK+IudF9lOSghxd7uK/w3D48FIJqSVYIy8rKPvrk01OnTkkVAID9fHBcDHClOyJUcqP57GBdehFPL8JDIaiTZIl6x733P7mrZeSt09APDypT3kSvHDS/G6uTOhCbcdfT68OFhbvMIkohqJFkhTAoMMC96LdGwT30C9Ocbeb/nOfNaHoBVXhxv3lGD8Em89Dl455IwU1PH51EloIKSdaDse6T1WlpaSNHjqwXnH4+L64+YZ6XxieGCdN7sOQwwVklTUqgOUcr+fpc8dg0tW05xIjei9VhAVJQJckKjiAIQ4cOdXNz83ehOb2En5L0x6c5jTawZVliyL9Nc7aZfzovmnD3CUqjuGVFrTfEn03CAqSgRjJ68gpwpQf6CDsm6bOm6ocG/K4itqAighKsPSNWNdOfomSUVrb1xnDdv0+Lv1WiqxBURY4ZG+bBHhsg7JikP9haEb9ARQS5a2ihxcpcVtR6/i70bLRu4S48FIKqyLEQtgr/b0XcN+V3FTE1n2P0GsjN0ixzjGKXFbXegn5CcQMWIAVVkXUhbNXN81JF3DtZPzSAvXTA3O2rlscyzTuKbLeJFMANsCwrulSxy4paTy/QO1iAFNRFYXkb4XWpIqbdpgt1Zw/uMHf78ncVURRxowoSeEL5y4paLwELkIJ07HGRV1ghbNXdiy2OFo5O0/9yi66LM92Xbu75dcvI+18O6j9y/qP/J3V0oC3bCvn+Mv6k8pcVtd7yGOHtI+aCerTIgOOIohiXdLtxwIh/f73ett+s+NQd4MuWDNVlT9d/k6jLTv+5fPanP23aInVQoCFmTn/JNL81QnBT/rKi1uvhxR7sIyzeg4dCcJzKysojF8qKkv667qeNtv1mxRfCVoP92fq/v9kvbYn3XUvMuE8FR/nguBjoSnd2V08qWQkLkIKDMU9/FjMz7uzaZS8ssu03qyp7E8ePO5ryVdeYif88gRtVcISKJnr1oHlFjHqWFbWeZQHSxzKxACk4yOI95vsefHjn959HRUXZ9ptVVQgtVsXpXjpgLm2UOg7QgOf3mWf0EAb6aWKMzJXuiRS8nOhjLEAK9revjP98QfzrTXa56VRhIezvy2ZFCs/uxeBusK+jlfybs+KL9slMRWBEK2J1z+8zVzVLHQqomshpwU7zmyN0dlrnVoWFkIheHqrbmMczS9BkA3ak4mVFrYcFSMEBVp8QnXU0q6e9CpY6C6GXE70xXHgkw4xRM2Anql9W1HpYgBTsqqKJXj5gXhVnx7ULVZvGs3sKXZwIo2bAHrSwrKj1sAAp2NXiPeaZkUK0PXviVVsIiWglRs2AfWhkWVHrYQFSsBO7jpFppeZCiFEzYA/aWVbUeliAFOyhdYxMFzvvBa3yZMaoGbC5J3aJf9HMsqLWwwKkYHP2HiPTSuWFEKNmwLZ+LeR7SrW1rKj1LAuQ5tch2cAGHDBGppX68xmjZsBWzJwezzS/EyO4a2lZUev18GIP9BGe3otcAxtwwBiZVuovhIRRM2AjHxwXA1xpaoQmsqZznsMCpGALjhkj00oTKW0ZNfMMRs3ADdDysqLWwwKkcOMcNkamlSYKIRG9PFSXglEzcAOe32e+S8PLiloPC5DCDXLYGJlWWimEXk60FKNmoLMsy4o6rKFG0bAAKdyIcgeOkWmllUJIRLMwagY6C8uKdohlAdJXDqIzAjrMkWNkWmmoEBJGzUCnYFnRTnhjuO5fp7AAKXTMvjK+4QJ3fNOLtnIbo2ago7CsaOf4u9AzWIAUOuK/Y2QEh42RaaWtQkgYNQMdtOwwlhXtpEewACl0hGWMzEwHjpFppblCiFEzYL0Ldfz9Y1hWtJOwAClYT5IxMq20mOEYNQNWwrKiNyghlA30YyuwAClcjyRjZFppsRASRs2AFTJKGJYVvXFvxwjLj5gL6qWOA2RMqjEyrTSa5Bg1A9fQ0tKy+KXXpj/24tJBDVhW9Ab18GJ3++Yn/Wnx+//8SOpYQI4kHCPTSqOFkIheGqrbmMd3YdQMXCElJeXvuwormoWGXeuljkUNSr5ffi445vmVnxQUFEgdC8jO6hOiq0RjZFpptxB6O9Gy4cICjJqBKwwcOFA4t98nZ9OwYUOljkUNJiXEu29b6aLXBQQESB0LyItljMx7Eo2RaaXpdp9ZPYX/ly2uPi7+uZ92bwjgSsbwrq4v7Nj3B3OPAA+pY1GD2XdPdx6Y+H6ul7Ozpi84cCVpx8i0Em655Zb33ntPFEUiMplM77zzzqRJk26//fYPPvig9cXXX389OTl57ty5OTk50oZrcyvjdC8fxKgZ+J09pTzcU2f0wlXbZsYZnQ6W81qT1HGAnEg+RqaV8NRTT61evfrdd98lopdffvnTTz998sknH3300RUrVqxYsYKIXnjhhR9//PHZZ5+NiIi4+eabGxtVVTT6+7LZPTFqBn4nNZ8nhEodhLq462loAEvDPoXwX3IYI9NKGD9+/FtvvfXuu+9yznfu3PnAAw+MGzduwoQJc+fOTU9Pb2hoWL169apVq8aOHfvyyy8HBwevX6+24QMv3YS1ZuB3thSI40OkDkJ1Eo3ClgJMKIRLJFxH5koCEcXFxZ07d668vHzChAk//vhjeXl5cXHxL7/8kpycfObMmYaGhqFDLw0ZiI2N3b9/v6QB2x7WmoG26lvoQBkfHSx1HKqTGMpS85FjQCT1OjJX0hNRly5ddDpdUVHRwoULU1NTQ0JCOOe33nrrfffdl56e7uvry9ilaP39/U+ePNned9XW1s6cOdPNza3ti48++uidd9551c/X1dW1frO0Jhvon8ec3jtkur+XgttI5XM8FW1TgTDYV8cb62tF0WRCp5Zt1NXV9XGnC7XOp0rrDG4ohzdK6cn+xB799K4U6dxUW2v3f8vd3V0QrvPcqSeipqYms9ns5eU1Z86ciIiIDRs2mM3muXPnLliw4E9/+lNDQ0Prp+vr6729vdv7Ljc3t6effrpXr15tX4yIiPD09Lzq5znn7b3leH8fw8f/0jK7r1ugq9ShdJasjqdy7aw0J4UzT09PURRdXLADoc14enqOCzXvrtbPCpRFa5iiKTrZ95Xx1CLzsWl6Txn0DlroiejMmTPOzs4Gg2HDhg2bN2+2JP/8+fMffPDBJUuW1NTUlJWVWSYAnTlz5qabbmrvu3Q63cCBA6Ojox0WvQ21jpr5cIz0Q5hAQqn5fPVoXKntIiGUbSngs3pKHQdIR1ZjZFoJRLRmzZo77rjDxcUlMjJy+/btlje2bdsWGRlpMBji4+PXrFlDRGfPnk1NTb3rrrukjNeeMGoGShrofC0fFqDgRic5SzSyTegm1LZ/HJfRGJlW+sGDBzc0NGzYsIGI3n///XvuuWfdunUmk6m+vt4yQPSdd965/fbb169ff+7cuaeeeuqylk81sYyaWbDTvHeKXi59uOBYqQXiuBBBL5BJwZ3F8hXVhekZZVfzqC5IMC0qb6IlB80pE/Vy++/Xf/HFF1FRUTqdjohGjx595syZc+fO6XS68PBwvV5PRIMHDz59+vSJEydCQ0P9/f2lDti+sNaMxm3J54lGuSWpqowPZan5KIQatXiPeVZP6deRuZLQr18/SxW0cHJy6tmzZ/fu3S1VsPXFgQMHqr4KWqwapXv5oLmk4fqfBPXZWohCaF+JRkyi0Ki9pXzDBf7iEDkOwsBzz+X6+WCtGY3KruYmkfCwYleJocK2QrEFE+s1RuT0SIbsxsi0QiG8ipdu0m3Kx6gZzUnN50l4HLSzIDfq6sn2lSG5tEWeY2RayTQsabWOmsFaM5qypYAnhKIQ2h1aR7XGMkZGPuvIXAmF8Opm9RR8nGn1cbTgaIWZ0/ZCcXwoMsLuEkKx6Ki2PLVbpmNkWiHt24VRM5qyr5QbPViIu9RxaMDYELa/jNe1SB0HOMTeUr4xT6ZjZFqhELYLo2Y0JbWAT0AHoUN46OmmAJaOLZk0QOZjZFqhEF7LSzfpfjl45vHXV+Xl5UkdC9hXar6YgHZRR0kIFVLz0TqqcocPH777xb8L9eWyHSPTSu7xScvLicR/zn7vlNvEu+6VOhawo/oW2lfG4w14InSQCUa2GeNlVK2lpSVx6qx1pxr1a5+Uf16hEF5HhCGAndrh7RckdSBgR+lFfIg/83SSOg7NGBbAztfyYnTAq5dOp2POrs7n9vSLMEody/Xpr/8RbduZ8sNdnx0eN3yQ1IGAHaXmi4lG3BQ6jl6gsSHC1gLxnkgcdnVijCW9uy2iLueVKYOljuX6cBZeh16vn3nz4E0FUscB9oSRMo6XaGRbCtA6qlqcaGupy7wEZezKh0J4fYlGIb2IN2H0qEqVNVJuDbZecrTEUHQTqtmRCu7hRD28lJFWKITX5+NM/XxZBlZcU6nN+eK4EMEJqeBYfXwY53SyGmmlTpsUtWAhst8qyWEsJQ+jvdUJK6tJJQFrralXSp6YHKaYtEIhtEpymJCSh4xVpy0F2HpJGgmh6CZUp4YW2lPCx4Yopr4oJlBpjQhk52t5Yb3UcYCt5VTzZjP19UEhlMAEo/BroYil7dVnWyG/KYB5K2c+EgqhVXSMbsZKwWqE8aISCnYjozvbjy2ZVGdTvpgUpqTioqRYpZVsZGgdVZ/UfLSLSmkCugnVKCWPJysqrVAIrTUxnKXkiSJyVkXMnLYVYolRKSUYseio2uTV8bJGPtgfhVCNwjyYvyvLqkAlVI/9ZTzUHVsvSWlcCNuLLZnUZWMeTwoTBCXVQRTCjkgOQ+uoqqBdVHIeehriz3ZgSyYV2ZSnpBmEFiiEHZBkFDZhNqGKbClAu6j0sGG9mpg5bS0QJyht5V6FhSutcSFsXxmvNUkdB9hCo5n2lvL4EIXduqoPxsuoyd5SHuahvO4GFMIOcNfTiEC2rRBJqwZpRXywv5KmOqnViECWW8NLG6WOA2whJY8raEGZViiEHZMcJmCtNXXAlvQyoRcoPkTYitZRVUjJE5MVNYPQQnkRSyspjG1CM44qYKSMfCSEonVUDaqa6WglHxWsvLRCIeyYQX6sxsTP1CBpla28ic7U8BGBystYVUo0YksmNdiSL442MBed1HF0HAphxzCiCUYBSat0qflivAFbL8lFPx9m5nTqItJK2VLyuRLbRQmFsBMwm1AFsOOE3IxH66jypeYrcqQMoRB2QpJR+LVANKFrX8nQQSg3iUZsyaRsJ6q4yCmqiyLTCoWwwwJcKdKb7caG9Yp1+iJvwtZLMjPBKGwpwJZMCqbQiRMWKISdkWRkm7BSsGJtzueJRqbUlFUpgxuFurMD2JJJsTbli4pbWa0VCmFnYMN6RdtSwBNClZqxKpZoZKloHVWmZpF2FvObFTsxV6lxSysumJ2sxloYiiRy2lYojkchlJ+EULYFDS3KlF7E+/syPxep4+gsFMLOcMJaGIp1oJwHubIwDxRC2RkXIuwp5Q3YkkmBNilzQZlWCg5dWklGtgmtowqUms8nKLZLX928nCjan+0oRlopT4oCt15qC4Wwk5LDWEo+R8oqTmq+iA5C2UoIZdiwXnGKGuhCHR+u5HWaUAg7qac3c9XRsUqUQiVpNNOeUj42BKe9TCWGChgvozib8sSEUEGn4DqIQngDkoxYYkZh0ov4ID9svSRfMUHszEUMQ1MYRc8gtEAh7LzkMMwmVJgtBSIWlJEzvUCjDexXDENTDk6UWiBOUHhaoRB23vhQIbMYg9yUJDWfJyp2qpNGJIYKWGtNQQ6W8QAX1tUThVCrLIPc0jHITSHKm+jURT4ySNkZq3qJRmz5qSQp+TxJ4e2ihEJ4g7BhvYJsLRDHGBi2XpK5/r7MJBK2/FQKhW5JfxnF/wLSwpZMCpKazxPQLqoE40OwT68y1LXQwTIeb1D+E+GJEyeqqqqkDkOpbvJnpY38Qh2SVgGw9ZJSJBjZFhRCJdiSL44MYu56qeO4YcKMGTN69OjxzjvvENGXX37p91/e3t6MsbNnzy5evNjNza319bq6OqljlhGBUUKogA1F5e9sDa9r4f19UQgVYIKRbS0QRWSV7G3K50nKbxclIiErKyszM/OFF164cOHCPffcU/FfS5YsiYmJiYiIIKLHH3+89XUPDw+pY5YXzCZUhE35fIJRQBlUhFB3FuTGDpYjreQuJY8nq6KVRSCiqKio+Pj4tWvXtn3jk08+mT9/vuXPnPPKykoJolOCW8KFzfnYUFTuthSgXVRJJhgZGlpkLreG15r4AD81pNWlp9rIyMjz58+3vrp3795Tp07NmDHD8uN7773Xu3dvb2/vl1566RrfZTabjx49uvf3SkpK7Ba8LBjcyOjB9mNDURkTOf1agCVGlSQhlKViWr28bczjt4SrpJXlUi+nu7t724r10UcfTZ8+3dvbm4gWLFjwyiuvODs7HzlyJCEhoXfv3jNnzrzqdzU0NLzyyivu7u5tX3z00UenTZt21c/X1dUxVewTPj5I/9OZln5uEk+tV83xtLmDFczP2cmH19XWWvtXTCaTKIomk8mecWlIR8cWjOjCdhU7lVbVuil/IIY9yCHZfznrNLWrubZW7vcr7u7ugnCdjsxLZ1lFRYXBYLD8uaGh4auvvvrxxx8tP3bt2tXyh4EDB957772pqantFUJPT8/PPvssOjrayvg4556enlZ+WM7+0J0vOWh+xdNV2jBUczxtLuO0mBTesYNjKYQuLordaVR+OnT8PYkG+bdk1XmgQfuqJE/2FpF2lpn+3zgXTzcJo7CZS3UyMzNzyJAhlj+vX78+MDBw9OjRV366tLTU8pgIbY0xsMPlvLpZ6jigHVvQLqpAiaFsC1pH5SqzhEd6sSBVVEEi0mdmZm7cuLGmpmbKlCmWlz766KP58+e3Pnc//fTTsbGxfn5+W7ZsWbduXUZGhnTRypSrjuKC2dYC8Y4INYwkVplGM+0q4esS8F+jMAlG4cld5jeGSx0HXE1Knqj0HSfa0j/77LPh4eFbtmxxdnYmovr6+qioqLlz57Z+wtfXd82aNXV1dZGRkRkZGa0PjtBWcpiwKZ/fESF1HHCFncV8gC/r4ix1HNBBsUEs5yIva6QAifsc4Co25fPlI3VSR2Ez+l9//bXtz+7u7v/4xz/avrJ48eLFixc7NirlSQpj7x5FM44cbclX/B4x2uQk0Ohgtq1QnNYdT/PyUt5EOdU8RkXr1+MMs41+PkzklFONSRSyszmfJxhxnitSghHLNsnR5jxxbIigpvXrVfSrSG2CkaUgaWWmsolOVvORgeq5ddWUxFCWir0J5SclX/Fb0l8GhdBmksMYtmSSmy0F4mgDc1FPX4a2DPBj9S08F1syyUxqPk9SV3cDCqHNJIYKaYW8ySx1HNDGlgJsvaRgjGg8FrWXmSMV3FVHkd4ohHA1vi7Uz5dlliBpZSQ1n2OkjKIlhrItaB2Vk02qaxclFELbSjKidVRGztXyGrUsCqxZE4xsC7ZkkpOUPFFl7aKEQmhbyWECtmSSj835PBFbLymc0YMFuLKsCqSVLDS00O4SPk513Q1q+32kNTKInavlhfVSxwFERJSaz7GymgokhrLN6CaUh+1FfEgA83aSOg5bQyG0JR2jm0MFLJAoB5xoa4GIJZtVIMHItuQjp2RhU56YrIot6S+jwl9JWklGtgmtozJwqJz7ubBwDxRCxRsfKmSW8AaJdzkDIqKUPLVNnLBAIbSx5DC2KR9d+9JLzceW9Crh7UT9MR5bBvLqeGkjH+KvwrRCIbSxbp7Mx5lllSNpJbYF7aIqYhk7KnUUWpeSxycYBVUOP0MhtL3kMIaxo9JqFimzmI814PRWiYRQAeNlJJeSp8IZhBa4UthecpiwCX37ktpRxPv5Ml9sL68WscEsu4pXNEkdh4aZOW0tECeodP16df5W0hoXwvaW8lqT1HFoGNpFVcZZoFEG9itaR6Wzr5SHebAQd6njsA8UQttz19PwQLa9CC05kknNxxKjapMQKmCtNQml5PMklbaLEgqhnSSHCVhrTSpVzXSiiseqaNdQIKIJRmzJJKUUlc4gtFDtLyat5DDMJpTM1gJxFLZeUp2BfqymmZ/FlkxSuGiio5V8dLBqby5RCO1ikB+7aMI+atJAu6gqWbZkQuuoJFLzxVHBar65xPXCLhjRBCMGfEsjtQBT6dUpAa2jEknJ40kqHS9qoebfTVpJRswmlMD5Wl7dzAdh6yU1SjKy1HxsySSBzWrcg7AtFEJ7SQ4TthaIJoyYcazN+TwhVJVrXwCFeTA/F3YYWzI5VnY1N4nUx0fNWYVCaC8BrtTDm+0pRdI61Ba0i6paIlpHHS4lj08MV3lOoRDaUbKRbcIkCgfiRL8WiONDVJ60WpYYii2ZHKvb3dEAACAASURBVG1Tnpis9ptLFEI7SgoTUjBexoEOV3BvZxbhpfKk1bKbQ4WMYt5kljoOzWgWaUcxv1ntw7BV/utJa1QwO1nNy7FAoqNsxtZLaufjTH19WQa2ZHKUHUW8nw/zU/uyvSiEduQk0OhgAS05DrMlX0wIRSFUObSOOpK6F5Rppf7fUFrJYWwTWkcdolmkjGI+LgSntMolGAWMl3GYTWqfOGGBq4Z9YW9Ch8ko5n000IYDo4LZ8UpeiR4H+ytuoPO1fHggCiHcmJ7ezFmgY1WohXa3pUCcgA5CDXAWKC6YbStE66jdbcoXE0IFnQayCoXQ7pLwUOgQqfk8QdWrQEGrBCMWHXWElDw1b73UFi4cdpccxrAlk71VNdPRSh6HrZe0ITGUYSFfe+NEW/LFJG20sqAQ2l1CqJBZzBtapI5D1X4tEONUvTo+tBXtz6qxJZOdHSzjvi6sqycKIdiClxNF+7MdxUhaO9pSgK2XNIQR3RwqbC1ETtmRRsaLWuDa4QhJRmxYb1+p+RwjZTQlMZRtQeuoPaXkiereeqktrfye0sIkCvsxmUzTH3j8zNuz/BsKpI4FHCeoeP+6J6e88tYKqQNRp7oW2l/G4zWzbC8KoSMMDWDFDTy/DrXQ9nbv3r0hp9rUe/zH//pC6ljAcVa8/bbp1mff/vBf9fX1UseiQr8W8JFBzEMvdRyOgkLoCAKjBGxYbx+DBg1yvZjvv//TyROTpI4FHOdPs6a7f/2X8N4D3N3dpY5FhTSyslorzVR8qSUZWUo+/2NvqeNQHW9v7+BnfvnXzbpof6004wAR3TN9at2gKWkYL2MfKfl8fYKGCqGGflVp3RImpOaLItLW1qqa6UIdH+iLKqg5o4KxDYVdnK3hNc18oJ+GcgqF0EFC3MngxvaXIW9tLKOYjwxiepzI2tPHh1U18UJ0EdraxjyeHCZoqAyiEDoSxo7aQ0axiAVltIkRxQSxzBJMTLIxTc0gtEAhdJzkMGET9lGztYxiHhuM01ijYoOFTCxVYVMtIm0rFLW2PIW2fltpjTGwrHJe3Sx1HCrSItL+Mj5SA9vEwFWhm9DmdpXwSC8W5CZ1HI6FQug4rjqKDWa/YvsY28mq4N28mC/2INSqkYHscAVvNEsdh4qk5Isa2XGiLf2CBQtGjBhx7733CoLwyy+/HDlypPU9V1fXxx57jIi+/fbbrVu3hoWFPfjgg76+vtJFq3jJYUJKHp/STeo41GJnMXac0DQ3PfXpwvaX8VHBOA1sIyWP/22k5lavF4YMGfLuu+8+9dRTRFRXV1f5X+vXr//222+JaMWKFYsWLRo0aNDhw4fHjRtnNuPuq/OSw9hGjJexncwSHocroLbFBbMMdBPaSGUT5VTzWA3eXHLOT5486e7uXllZydvo27fvp59+ajKZQkNDU1NTOedms7lXr14//PADb8egQYMOHTrU3rtXunjxovUfVo2uX5pyqkV7fLMGj6f9DmZzc3NjY6M9vlmbampq7PTNX50237G5xU5fLlt2SvavTpsnb9LcweScC0TUq1evgICA/fv3t1bHHTt25OfnT5s27fTp06WlpWPHjiUiQRBuvvnm9PR0yYq2KkwwYhKFbeTX8UYz7+mtvbtXaGN0MNtZjH5320jJ4xrZifcyl5ZYCwoKKiwsbH31o48+uvvuu93d3QsLC319ffX6/30sNze3ve+qr69/9tlnfXx82r54zz33jB8//qqfb2ho0Ok01xg9NoB9dVY3L8L2G/Vq7XhuOc9G+uvstOayyWQSRREdAbZSX18vCHYZmufLyJk5HSmuj/Syx9fLlJ2SfXOe0xNRJpUtY+7q6nrdc+9ShTOZTM7OzpY/19bWrl+/fvPmzUTk7Ozc0vK/S7bJZHJxaXeInpOT08iRI8PDw9u+GBER0d5faW5uvsa3qdXEbvTIHpH0LjbfTl1rx3N/FR9lIDv9yoIgiKKoqeNpV9e+dNyguGC+v8q5X4CGHmXskexHq8hJx/sGqO2cZ+z6J4aeiDjnBQUFRqPR8tJXX33VtWvXkSNHEpHRaKyqqqqtrfX09CSivLy87t27t/ddTk5OkydPjo6OtjI+nU6nqScYiwB36uvL95QL42y915fWjmdGccuKWJ1OZ5fLnyiKjDFNHU+7suvJOcog7irlf4zS0H+WPY7n5gJxYjhp85wXiGjr1q16vX7EiBGWlz766KP58+db/tytW7dBgwatW7eOiKqqqlJSUiZPnixVrKoxtPnEq8tXFhcXSx2IgjW00IlqPlRLDwHQnjhMq79hJ06cWPP3VSOcNXpREmbPnj1z5sylS5c6OTkRUXZ29v79+2fPnt36iaVLly5atGjmzJmxsbETJ04cNmyYdNGqxNpnZm0tdZ465wGpA1GwPaV8gC9z1eLNK1xusB87W4M1m25Iwh0zsxtcVz51v9SBSEM/ZcqUJUuW9OjRw/Kzv7//wYMHg4KCWj+RnJyclZW1e/fuRx99NDY2VqI4VSXQt0v5qZ3BvUKlDkTBdhZjDjVcohdoaADbVaK5paJtyMXDS396R7fuGr0o6adNm9b254CAgICAgMs+ZDQap06d6sCoVO7QjtSxa36bd9tgqQNRsIxicX4UFgiES0YFs4xiMTkMTQSdNPv9lPzTxz+8R6MXJVxKJODs7Hxr3OAdWA6jszjR7lIeG4SzFy6JDRbQTXgjMsr1d44ZbM0AS1XCpUQa8QaWVoS87aQTVdzbiYW4Sx0HyEZcENtdwlswsb5TTCLtKdX0aoUohNKICWK/VfJak9RxKFNGsaaTFq7k60JhHuy3Stxcdsa+Mt6rC/NxljoO6aAQSsNFR0P82S405nQKCiFcaRRW3+6stEIeb9B0QqEQSibewNKL0JTTGRnYdAKuEIvZhJ2VXiSOQSEESYwxCOgm7ITyJiqo4wN8NZ23cKW4IDwRdobIKaOEjw7WdC3Q9C8vrVHBbF8ZNtfusMxiPjKI2WdhNVCwKB9Wa+L5daiFHZNVwUPcWJCb1HFICoVQMp5O1NeH7StF3nZMRrGIdlG4EiOKCRLQ795RaYU83tbrHisOCqGUMImiEzJKeJy2m3GgPegm7IT0Yq7xDkJCIZTWGIyX6SCTSAfK+MhArectXBUGjnYUx0gZIkIhlNYYg5BRjFnAHXConPfwYl00POEJrmF4ADtSwRtsv+m1ap2o4p56Fu6BQgjS8XOhbl7sUAXuYa2FGYRwDW566u/L9pchoayVVoQOQiIUQsnFG1haIfLWWhklPBaFENoXF8x2onXUaulF6CAkQiGU3BgDS8d4GatlYPcluKa4YJaJ8TJWSy/S+poyFiiEEhsbIqQViSIy1wrna7lJ5D28kLfQrlHBbGcx8skquTXcJPKe3kgoFEKpGdzI35Udq0LmXt/OYj4KEyfgmkLdmYee5VQjoa4vrYiPDUFCEaEQygG6Ca2UiQ5CsEIcJlFYJ72Ij0FCEREKoRyMMbB05K0V0EEI1kA3oZUwZLQVCqH04g1seyHmEl5HXQudqOI3+SNv4TrigjBw9PqKGqi8kffzQUIRoRDKQXcv5iSwUxeRuteyu4QP9mcuOqnjANmL9mcXanlFk9RxyNv2QjHeIAiog0SEQigTY7Do6PVgKj1YScdoWCDbjdbRa8IMwrZQCGUBswmvK7MEm06AteKCWWYJuhuuBR2EbaEQygIGjl4bJ9pVwmOCcLqCVeKCBHQTXkNFE52r4YP9UAgvwZVFFvr4sNoWfgF7irbjWCX3c2EGbe8dCtaLDWZ7S7GcfbvSi8S4YKbH5f+/cCRkgRGNMQhoHW0POgihQ3ycqasnO4zl7NuRXsTHGHDx/x8cC7kYE4xuwnZllqAQQseMwia97UvDEqO/h0IoF/Eh6CZs185iHheEvIUOiA3C+jJXV2ui41V8GHa3bgOFUC6i/VhhAy9pkDoO+SlrpOIG3s8XeQsdgIXW2rOzmA8LYK6YktsGCqFcCIzigtiOYvTvXy6jWIwJYjrUQeiIXl1Yg5nnYQDaFdKLRLSLXgaFUEYwXuaqMkt4HDadgA5iRLFBAhYdvVIaRspcAYdDRuJDsL7MVaCDEDoH3YRXajLToXLs4nI5FEIZGRbAcqp5VbPUcciJSaSD5XwECiF0HLoJr7SrhPf3ZR56qeOQGRRCGXESaEQgUvd3DpTxnt7M20nqOECBhgWwo5W8vkXqOOQEEyeuCoVQXsYYhPQijJf5n4wS7EEIneSmp4F+bF8Z7iz/J71IRAfhlXBE5AXdhJfJKOaxaBeFzooLxt6E/9Mi0p5SPhpPhFdAIZSXmEB2uAKNOf+DNWXgRsQGsUwUwv/aX8a7ezEfZ6njkB8UQnlx01O0H9uFMd9ERHS2hptF3t0LhRA6abRByCgWkU4W6CBsDwqh7MQbWBq6CYmIKKOEj0Z/BtwAgxt5O7OT1SiFRNiMt324ysgOptW3Qgch3Dh0E1qInDKKMVLm6nBQZGe0ge0t5c14JiTKKMaQUbhRcegmJCKiI5U80I0FY1PPq0EhlB0vJ+rdhe0r1Xrq1rVQzkU+JACFEG4IptVbpBWig7BdKIRyFG/AJAraVcKH+DNnnKFwYwb6sYJ6XtEkdRxSQwfhNeAyI0djDAzT6ndiV3qwBR2jYYEMq2/vKMamE+1CIZSjsSFCRjE3aztzM4tFrLUNNhEXxDK1vcFZdjV30bGunkioq0MhlCM/FzJ6sKxy7VZCTrSnlMdi9yWwhbhgQePdhOggvDbB3d09Jibm5MmTlp+zsrLGjBnj7Ozs7e29dOlSIlq2bFlkG/X19ZIGrBUa7yY8WskDXFmgq9RxgCrEBrN9Zdyk4WdCdBBem1BbW5ucnHzfffcRUUlJSVJS0syZM6urqwsKCm677TYiqqioSEpK2vxfrq64ODnCGAPT8mxCdBCCDXk7UYQXy6rQbkJhTZlrEwRBWLRo0d69e3NyctasWTNs2LCHH37Yzc3N09NzwIABlg/5+Pj0+C9BQGuVI4wLEdKLtLs0VCY24wWbitPwJr0X6niTmffugoRql0BEnp6eXbt2zcnJOXz4cLdu3SZOnBgWFjZ58uTc3FzLhz788MPg4OARI0Z89dVX1/guznldXV3N75lMJkf8HqoT4k7ezux4lUZTNwNrbYNNxQVrd+DotkIeH4IHmGu5tFGxt7d3RUVFUVHRL7/8smHDhptuuumZZ56ZNm3a/v37p0+fPn/+fD8/v23bts2bN8/f33/ChAlX/a6amprExMTLHhlfffXVP/3pT1f9fF1dHWO42LUr1l+/+aypa0+zlZ9XzfEsb2KlDU7h+sbaWsliMJlMoijiNs5W6urqpA1gsBd7ptC5trZB2jBspUPJvvWCfoQvr6219kqiMu7u7tdtyLxUCKuqqgIDAwMCAv7whz+MHj2aiP7617/6+/sXFBQMGzbM8pnp06dv27btm2++aa8Qent7f//999HR0VbGxzn39PS08sMaND5cTM3nj3nqrPy8ao5narkYGyx6e0n5u1gKoYuLi4QxqIy0J+cgTxK5qZJ5hHuo4WaxQ8meWdayMFrnibkT7ROIqLKy8sKFC1FRUf369eP8UuuB5Q+X3XRwztXxzKEI8Qa2XZPjZTKLeVwQWnLAxmKCtDiJorSRShv5AF9ct69FqK6ufuGFF8aNGxcREXHfffelpKSkp6fX1dW99NJLMTExISEhH3744cmTJ0tKSr766qtPP/106tSpUsesFZHeTCDKrdFc6mLIKNiDNhcd3V4ojjYIAvLpmoSoqKgLFy58/PHHRBQREfH5558/9thjUVFRpaWl69atI6Jdu3b94Q9/GDRo0KpVqz799NP22kXBHkZrbzZhk5kOlfMRgUhcsDFtFkLMILSGvqioqO3PkyZNmjRpUttXPvzwQ8eGBP9jmU04t5fUcTjQgXIe5cM8naSOA1RnWAA7Uc1rTaSpsyutiK8ejY6G68ABkjUNri+DPQjBTlx0NMiP7SvTUEJVN9OZi3yIPxLqOlAIZa2fL6tu5gX1Gkpd7EoP9qO1afXpRTwmiDnhMn89OEKyxohGBQuaWmsts0TESBmwk9hglqGlbSjSi8QxBlzkrw/HSO40tejomRrOiHXDhCewj1HBQkYJ187ShWlFPD4E2XR9KIRyp6luwoxiPhoj3MBugt3Iz4VlV2sioepb6GglBmBbBYVQ7gb7s7w6Xt4kdRwOgQ5CsDftdBNmFPMh/szV2pWpNA2FUO50jGKC2I4iTXRsZGAqPdhZrGZmE6YXiWgXtRIKoQKMMWhivEyNiU7X8MEY6g32NCqYZWhjG4q0Io6RMlbCYVIAjXQT7irhQwOYM05JsKcBvqyonpc1Sh2HnTWZaX8ZOhqshauOAgwPZNlVvEbtOwJlYDNesD+B0YhAtkvtD4V7S3lfH+alpTV0bgQKoQI4CzQ0QP0dGxnFYiw6CMH+tDCbMK2Ix2MAttVQCJUhPoSlq3q8jMhpTymPwe5LYH+W2YRSR2Ff6UUi1tq2Hq47yjDGIKi7m/C3Sm5wZ4GuUscBGhATxPaX8Wb13liaOWWW8FEYKWM1HClliA1ih8p5Q4vUcdgNOgjBYbycKNKLZZWr9s7yYDnv5sn8XaSOQzlQCJXBXU8DfNmeUtWmLmYQgiPFBbOd6u10TytEB2HHoBAqRnyImidRZJSgEILjxAWzTPV2E2Iz3o5CIVSMMQZBreNlihuooolHdUHqgoPEBbMdKr2t5EQ7ikWs2dshKISKMTqY7SpRZw//zmIxLogJyFxwlB5ejBM/V6vCWni0kvu5sFB3pFMHoBAqRhdn6unNDqpxf+3MYh4XjFMRHCo2SFDl3Fx0EHYCrj5KotZuQnQQguOptZswvRgdhB2GQqgkY4JVuElvk5kOV/BhAUhdcCi1DhzFE2EnoBAqSXyIkF4kmtWVvPvKeJ8uzBOLIoJj3eTPsqt4rbqW8D11kQuMIrxQCDsGhVBJAl0pxJ39VqmqSogZhCAJFx0N9lfb3Ny0Ij4OexB2HAqhwsQbWFqhqlI3Ex2EIJFRqtukFzMIOweFUGHGGNTWTZhZLKIQgiRig1lGiaomJKGDsHNQCBVmbAjbXiSqphKeush1Agv3QOqCBEYHC7tKuGrSKb+O15h4lA+yqcNQCBUmzIN56NnJapXkbkYxH43HQZBIgCv5u7DjVSrJpu1FfGwI1qXoDBRC5VFTNyE6CEFaauomRAdhp6EQKo+augl3YsgoSCo2mKlmk150EHYaCqHyxBvYNlU8EV400dkaHu2H1AXJqOaJsKyRCur5QGRTp6AQKk+vLszMSQXrBWcW82EBzAnnIEinvy8ra+SljVLHccPSisTRBqZDHewUXIQUaYxBDYuOZpZg4gRIjBGNCGSZxYqfRJFexMcYcD3vJBw4RVJHN+HOYh6LTSdAarHBggq6CdOK0EHYebgMKZIKBo6aOe0r5bFBSF2QmAq6CS+aKKea34SV6zsLhVCRBvixskZeWC91HDfgSAUPdWd+LlLHAZoXE8QOlit7y+sdRXxkEHPG5byzcOQUiRGNNgg7lNyxgbW2QSY89NRL4VtepxeJ6CC8ETh2SqX0bsKMEh6LQgjyoPS9CdFBeINQCJVK6d2EGcV8FAohyIOid6tvaKHDFXxEILKp81AIlWqIPztbyyuapI6jU4ob6GIz790FqQuyEBek4PEymSV8sD9z10sdh5KhECqVXqCYILZTmd2EO4rEuGCsDgxyEeHFGKOzNYqshWlFItpFbxAKoYKNMQgK7SZEByHITWyQUrsJMZX+xuHwKVi8YteXQQchyI1CuwlNIu3FfNwbhkKoYCMC2dFKXmOSOo4OajLTb5V8GCb/gpwotJtwbymP6sK6OEsdh8KhECqYi46GBrBdSruN3VvK+/mgbx/k5aYAlnORX1TabWVaEY8PwT3ljUIhVLZ4A0svUth4GexBCDLkJNAQf7ZHabeV6UUiNuO9ccIHH3xw7ty51p9FUdy4ceOqVau++eab6upqy4tHjhx5//33v/vuu5aWFonihKsbYxAU102IXelBnkYpbZNeM6fMEj4aK9ffMOHo0aNDhgzZu3cvETU2NiYlJS1evDg7O/uLL77YtGkTEa1bty4hIeH06dNvvvnm5MmTpQ4YficumO0v441mqeOwGifKLGxC3z7IUGwQy8hX0szcrHJudGcBrlLHoXz6VatWGQyG11577fvvv1+2bFlLS8u+ffucnJwsb3POX3zxxVWrVt1111319fU9e/ZMT08fM2aMtEFDKw899fdle0u5IppH6urqBo+7tby87nDw62HJSVKHA/A7mV+u3PTPf9/ySe+N6/8ldSxWQQehrQhENGXKlI0bN3LO161bt2DBggMHDvz888/l5eVElJubm5OTM2nSJCJyd3dPSkrasGGDxCHD7yloEsWFCxcKzJ7msQ/9snmb1LEAXG7Txg18xvLM/Yc4V0ZCpRfxMehlsAU9EYWGhjY1NZWWlubm5r7zzjsBAQGurq7z5s37+eefGxsbfX193dzcLJ8ODQ0tKCho77uamprWrFljMBjavpiYmDh06NCrft5kMplMShukJT9xgfTBCXqqv1n+xzMsIlLXP2GCaf9Tf3lG5qGaTCZRFAUBvS+2If+Tk4j+8bdX71/yvnn2y/IfDGEymZpNph3F9M4IMpmU0zUiBb1ez9h1bhf0RGT5kNlsNplMAwcOXL16NRH99a9/fe65555//vm2XyEIgii2O0aRc15bW3vx4sW2LzY1NbX3V0RRvMa3gZVi/GluidDcIsr/eH6cw+JnPfrNzSIRyTxUUVTA8VQQRRzMQQMH7vj6HwN+EDJL+MgAWT8UiqJ4rFL00gvBLvI/rgqgJ6KioiInJ6fg4GCDwRAbG2t5Iy4u7qOPPgoJCamsrGxubnZ2diaiwsLC0NDQ9r7L1dV14cKF0dHRVv7bzc3NLi7YmPVGGVyou1fLsVrnPq4ucj6eZk7vHW/5ZKzOxUUBjTmWez45H09lMZlMijiYLkT/N0h85xj/boJO6liupbm5eVee09gQ7uLiJHUsaiAQ0S+//DJ+/HhBEJKSko4fP25549ixYxEREZGRkeHh4SkpKUTU3NycmpqakJAgZbxwNfEhCugm/OK0GO5BWFkNZO5PUcKeUv5bpdwTKr1IGUPkFEH/wgsvrFq16ocffiCiRYsWxcfHm81mNze3Dz744KuvvhIE4YUXXnjwwQcPHTqUnp7erVu3xMREqWOGy40xsH+f4g9ESB1H+zjRm1ni2zGyvssGICJXHT3aX1iWJX4+Ttana3oRf2UY+rBtQ3B3d9+5c2d8fDwRRUVF7du3z2g0+vv7Z2RkWGrevHnz1q9fzxi75557Nm3adN1eR3C8sQYhvUgUZXwL+/1Z0UmgRCNOHlCAP/cTNuaJpy7KN6PO1jGRqIcXEso29M8880zbn8PDwxcuXHjZh+Li4uLi4hwYFXRMkBsFurLj1SzGW+pQ2vHWYfH5IdiAEJTB24ke7CMsPyJ+MEqmD4U7S4SxaBe1HTxZq0S8gVJOVkgdxdWl5vPKZprSDScbKMbCgbp1Z8SCepk+FG7NKUcHoQ1hCwA14Jz/8pfxZSanlpm3vPLsIqnDudwbWeZnogU8D4KC+LvQ7F7Cit/EN0fI7qFw4vTZm09WVPYLfejLD6WORSVwk64GJpOpvqaqKfaPaXsOSh3L5XaX8DM1dE8kzjRQmP8bKHyULVbIb/HRg0ezzTf/OTf7mNSBqAcuT2rg7Oy8dvWKoabjgbPekDqWy72eJS4aKDjhRAOlCfNgUyKElUflNV/dzMn7vn8k1qZ9/9k/pY5FPdA0qhKJCeP/c9OIYRtd95TyEYFyaYU8VsX3lIhf3Yw5v6BIz0QLsT+2PDFQ8JLNKbzmhGjsE/3tA729vLykjkU9cKOuHh56/tow4bFMs3z69187KD4+QOeG2y1QpkhvdnOo8M8TcnkovGiiVw6Kfxspu25LpUMhVJU5vQSR05enZZG3Z2p4Sp74UF+cY6Bgzw4Wlh8xy2TLzyUHzLd1ZUMD5NLkoxq4SKkKI1oRq3t6j1gng9Xz38wS/9xP6OIsdRwANyDajw3xZ5/lSH9zefoi/zRHfHkoHgdtD4VQbWKD2CgD+9thifO2qIHW54qP9kfSguK9MES3LEtskboULtwlPjVIZ3CTOAxVQiFUobdGCCuPms/XStlX+LfD5jm9hEBXCUMAsI2YIBbmQWtzpayEWwv40Ur+l/64YtsFDqsKhXmwP/cTnt0nWd5WNNHHJ8XHB+DsApV4ZrDutYOSLedr5vT4LvPbMYILWljsA5cqdXomWpdexHdItDfTe0fNUyOErp7o0geVuCWMuenpP+elublcfVwMcKHJWKTQbnBk1clNT68OE/5vjwRTKepa6IPj4qJBOLVAVRZHC68dkqAQVjXTKwfNK2LxMGhHuFqp1uyeAuf071OOTt0PjovjQ4XeXfA4CKpyZ4RQa6JfCx19b7nkgHlKhDDIDwllRyiEqsWI3o3VPbPXoVMpmsy04jdx0UCcV6A2AqP/GyS8ccihMwpPXeSfnxL/ehMeB+0LFyw1iwli8Qb2ZpbjUvfjk+IQf7oJE35BjWb3FE5WU2aJ4x4KH880Px2NKRN2h0KocstGCH8/Lp5zyFQKM6flR8Sno3H3CurkJNCTA4W3HDVJd0sBP1FNj/TDVdrucIhVLsyDLegnPLPXEan75WnR6EGjgvE4CKp1fx9hdwn/rdLud5ZmTgt3md/BlAmHQCFUv8WDdBnFdp9KwYmWZYnP4HEQVM1VRwv6CW9m2f3O8oPjYqArTeqKS7Qj4Cirn5ueXhsmPLbLbNfpwD+eE/UCJYXhcRBU7pH+wi8XxNMX7ZhOlU306kHzihjcVjoICqEmzOwpuOnoX/acSvHmYfH5wQLKIKietxM92Fd4+zc7ZtPLB813dhcGYsqEPvHtSAAAESdJREFUo6AQagIjWj5S99w+sdZkl+/fUsArmuiOCJxOoAlPDNR9fVosqLfLQ+GJKv7FKfFlTJlwIFy5tGJkEBsXwpYdtstUijcOmZ+OFvA8CBrh70Kzegrv2ueh8Ind5mcH6wKwYL0DoRBqyNLhwj+Oi2drbHwbu6eUn6ymeyJxLoGGLBokrMkWy5ts/LWp+Tynmv6MKROOhcOtIUYP9mh/3dO2nkrx+iHxqWjBGacSaEmYB5vSTXj/mC2zqUWkhbvMK2J1yCYHw/HWlkUDhV0lPN12UymOV/HdJeJ9vXEigeY8O1h4/5jZhv3ufz8uBrvRH8LRx+BouH5pi5ue3hguPJZps6kUrx8SHxugc9Pb5tsAFKSnNxtrENZk2+ahsLKJXj9kfgdTJqSAQqg5d0cKnk70WY4Nsje3hm+4ID7UF2cRaNRzQ4S3DpsbbTEE7a8HzNMxZUIiuIRpDiNaEaN7dp/54g036bx5WHy4n+DjbIuwABQo2o8N9mef3/Bt5fEq/vUZ7DIhGRRCLbopgE0wCstubFeKogZad0b8S3+kLmjaM9G6pVliy42Vwid2mZ/DlAnpoBBq1JsjdB9m39BUiuWHzbN7CYFIXdC2MQYW6kHrcjtfCX++wM/W0sPoYpAODr1GBbvRo/10T+3pZPZWNNFHJ8WFA3D+AFx6KOzcTWWLSIv3mN8eqXNCMkkHx167/m+QsLeMp3VqKsXKo+IdEUI3T3TsA9Ct4cxJoJ/PdyaVVh0TQ9xpIqZMSAqFULtcdbR0uPB4x6dS1LXQB8fNTw3CyQNwyVODhFcPdbjTvaKJlmZhyoT0cC3TtBk9BC8n+qSDY97+cVwcFyL07oJ7WIBLpnUXqppoW2HHbipf3G+e0UMY4ItUkhgKodatiNU935GpFE1meuc3EY+DAG0JjJ6KFt7oyEPh8Sq+Lld8EVMmZACXM60b4s+SwzqQwJ/kiNF+dFMA7mEBfufenkJ2Ne0ttfahcOEu8wtDdP4udg0KrIJCCLR0uO7DbPGUFTtumzktPyI+E417WIDLOQm0cICwNMuqjoafzovna+nBPrgCywL+G4CC3WjhAKt2pfjqtBjiRqMNeBwEuIr7+wg7i8XfKq9zT2kSadFu8e0YTJmQC/w/ABHRkwOFrHKemn+tBOZEyw6LzwzG4yDA1bnr6dH+urcOX+eecuVRsac33RKGG0q5QCEEIiIXHb0+XFi0x2xuvxT+dE7UMUpG9gK079H+wi8XrrVmU0UTLTtsXjYCN5QygkIIl0zvLnRxpo9Ptnszu+yw+NxgAWUQ4Bq8nej+KOFvR9rNo+f3mWdGCv0xZUJOUAjhf1bE6F7cb65uvspbWwt4eSNNjcAJA3Adjw/QfXFaLKy/ylvHqvg3Z8Xnh+BxUF5wXYP/GezPJoYJr19tKsUbWebF0QKeBwGuK8iNZkUK7x69Sh4tzDT/FVMm5AeFEH7njeG6j06KOdW/6+HYW8qzq2hWT5wtAFZ5Klr48IRY2fS7F384J16oowcwZUJ+9H379g0PD3/ttdeGDx/e0NBw++23t743ffr0Bx54YM2aNWvXrm198ccff3Rzc5MiVHCEIDd6cqBu8V7x28T/td68fkhcNEhwRv4CWCfcg03qJrx/THx+yKW0aRbpqT3ie7E6PfJIfvT/+c9/UlNTJ06cmJubyxhLTU3dsGGDXq8nom7duhHRqVOnAgIC7rvvPstfcHbGfuQq98RAYcA3LZvz+QQjI6LjVXxXifjvm52kjgtASZ6OFsb81PL4AMHTiYjovaNiVBeGQdfypI+MjIyMjPzkk0/Wrl07Y8YMIho/fvxl1S4iIiIxMVGiCMHRnAV6Y5iwcJf50B16vUBvHBL/MkDnrpc6LABFierCxoYIH2aLjw8QShvpzSzzjknIIpm69JQ+bNiwI0eOWP58++23T5w4cdmyZY2NjZZXvv/++3Hjxs2ZM2ffvn3ShAmOdWd3IcSdPjop5tbwny+I2DsboBOeGywsPyI2i/TCPvO9vbBhi3xdukPx8fE5ffq0k5PT22+/PXz48PLy8iVLluzbt2/dunVjx46NjY0NDAzcsmXLmDFjMjIyhgwZctXvunjxYmJiopPT79rQnnvuuTlz5lz187W1tbb9ZTTOtsfzpf7sjm8L/1V3dO7o8bqmxpqm6/8VNTGZTKIoNjdfbSoJdFxdXR3nndvCXcEinamXvnHy0m37/UYfuMurpv1Z9h2Fi6f13N3ddbrrzFe5VAhra2v9/PxcXFwWLlxoeaV///69evWqrKy89dZbLa+MGjXq1KlTn3zySXuF0MvL65///Ge/fv3avujj4+Ph4dHeP+/l5WXlLwPWsOHxHEQ1tSvuTO8xOqpyv9fYV2z1tUphKYQuLhjnbhuMMU9PT6mjkEDx+/cc6zI0OHdF2EM7bfvNuHja0KVCmJ2dnZCQ0PaNoKAgIqqtrfX19W19MTAwsLKysr3vYowFBQUZjUb7hAqO5unMzEz0cEJ7DkAn+bjonBj3cUHngqwJRHTo0KG0tLS77747Ozs7Pz+fiJqamp577rk+ffqEhYVt377dbDYTUVZW1ueffz5hwgSJQwaH8PLy2pf60/eLpvztlReljgVAqVK++ff6+4ZnbvpR6kDgWvTBwcEtLS0rV640Go3ffffd/Pnziai5uXnYsGHffPMNY+zFF1/cvXu3u7s7Y+yJJ5645557pI4ZHKRr165du3aVOgoABfP09Lx90iSpo4Dr0Ofk5Li7u1smDt5xxx133HFHdXW1u7t765iX7du3Nzc3NzY2ent7SxoqAACA7emvLG9dunS57BVnZ2fMowcAAFWSsgu3trZWgyOq7UQUxbq6OqmjUA9LK4jUUahHfX29ZagB2ERNTY3UIaiKlIVw6NChBQUFEgagJkePHr1s3C/ciJUrV7788stSR6Eed91117Zt26SOQiVMJpNl/UuwFQzqBQAATUMhBAAATbPlIrBms/nkyZPWf95kMh07dqysrMyGMWhWTk5OY2NjVlaW1IGoRGFhYUVFBY6nrdTW1p4+fTogIEDqQNSgpaWFc46T00pRUVGurq7X/gyz4XCV5cuXf/zxx9dd1a3VxYsXvby8GMPCJTZgNpvr6+ux6pKtNDU1cc6vmz9gpbq6OhcXF8s0Lbhx1dXVVw7vh6tau3ZtVFTUtT9jy0IIAACgOOgjBAAATUMhBAAATUMhBAAATUMhBAAATdO99NJLjvmXmpqaduzYcf78+dDQ0MtGlnLODxw4cOTIkcDAQIzTs1JDQ8P27duLiopCQ0MF4Xc3NOfPn9+zZ09FRYXBYLjsLWjP0aNHDxw44OXlddX9Y0VRzM3NZYzh/LSG2WzevXv3iRMnQkJCWpfvb1VfX5+Wlpabm+vj4+Pm5iZJhMpSV1e3ffv20tLS0NDQy4bZNzQ07N69Oycnx9vb293dXaoIFY87RGlpad++fUePHh0bGxsdHV1ZWdn6liiK06ZN692792233RYUFHTgwAHHhKRo586d69q1a2Ji4pAhQ8aMGdPQ0ND61kMPPRQUFJScnNy3b9/+/fsXFhZKGKdSLFq0KCwsbPLkyQEBAT///POVH1i5cqUgCK+99prjY1OcxsbGsWPHRkdHJyUlhYeH5+bmtn1327ZtwcHBsbGxSUlJcXFxEsWoJKdOnTIajUlJSdHR0ePHj7dM7LHIyckJCwsbP378nXfe6efn99NPP0kYp6I5qBC++OKLt99+O+dcFMXk5OQ33nij9a3NmzeHh4dfvHiRc/7qq6/eeuutjglJ0R5++OH77ruPc24ymYYPH/7hhx+2vrVv377m5mbOudlsnjBhwpNPPilZlApx6tQpDw+PvLw8zvkXX3zRp08fURTbfuDs2bMDBw5MTExEIbTGJ598MnjwYMtJ+PDDD99///2tb1VXVwcFBX399deWHy0Tw+Ha5s2bt2DBAs55U1NTdHT0v/71r9a3HnvssRkzZlj+/O67744cOVKaEJXPQe1m33777ezZs4mIMTZz5sxvv/227VuTJ0+2zASfNWvWxo0b6+vrHROVcrUeT71eP2PGjLbHc+jQoZbGKEEQ+vTpU11dLVmUCvHdd9/Fx8cbjUYimjp16rlz544dO9b6Luf8wQcfXL58+VWbTOFK33777YwZMywn4axZs9qenD///LPBYJg6depvv/1WXl5u/eIbWtaa7M7OztOnT297PF1dXV1cXCx/dnZ2Rrt9pzmoEF64cKF1r/OuXbvm5eW1vpWXl9f6Vnh4OBHl5+c7JiqFam5uLikpae94tiooKPj666/vvvtux0anPHl5ea1r+bu4uBgMhgsXLrS+u2bNmpCQkAkTJkgUnfJcluzl5eWtG4SdOnXKxcVl8ODBf/nLX6KiopYtWyZdmMpw8eLF6urq9pJ90aJFFRUVd9111wMPPPDpp5+uXLlSojAVz0ErHjU1NbX2mTs7O7fd6a3tW4Ig6HS6pqYmx0SlUJZGp9aD5uLicuXOeTU1NVOnTp0zZw72Zrqu5ubmtkt/OTs7t56BBQUFb731VkZGhkShKVJzc3PbZCeipqYmDw8PIqqpqTl06NDhw4f79OmTnZ09ePDgO+64o3fv3lKGK2+WU7G9ZD958uSJEydmz57t7e2dkZGxa9eugQMHShOowjmoEIaEhLQurl1WVhYSEtL6lsFgaH2rqqrKZDK1fReu5Onp6enpWVZWZnmALi0tveyI1dXV3XbbbdHR0W+++aZEMSqJwWA4fvx4649tj+fy5cuDgoKWL19OREePHq2oqAgJCZk3b540gSpE24wuLS11dXX19fW1/BgSEtKzZ88+ffoQUVRUVM+ePQ8ePIhCeA3+/v5OTk5lZWWBgYF0RbI/++yzDz/88BNPPEFEMTExEyZMmDdvHhZ07QQHNY3Gxv7/9u4YJJ0ojgP486Sw5fA8CFGQykLKCy5J6ySIhnAItEhCioLGhkAcheCCAodAIVEojGiIXBuCgyAhcgkSHFqkoYKgFqEl6qjXIBzyh/5Fgyb3/ay/N/x43OMH937vPUl7lvPs7Mzv92shv9+vhQqFgsvlslgsjcmqddXPZ6FQqJ/Pl5eXYDDY3d2dzWZxoflPSJJ0fn5eez+9XC6/v7+73e5aKBAIBINBjuM4jmtra+vo6MBO4bf+WdGSJGnf4djY2NPT09vbGyFEVdXHx8fOzs5m5dkSGIYZHR39arG/vr5q+4Imk0lVVYq7o3+nMT05l5eXLMtub28nk0mWZcvlMqV0fHx8a2vr+fnZbrdHo9GDgwOHw1HfAAlfURSF47jd3d2NjQ2z2Xx7e0spHRgY2N/fX1pa4nl+c3MzkUgkEomjo6NmJ/vXfXx8eDye+fn5w8NDURTj8TildG9vTxCE+mHT09PoGv2J+/t7juPW19dzuRzP8ycnJ5TS5eXllZUVSunk5GStvSsSiXi9XjSOfuv4+Jjn+VwuJ8uyxWJ5eHiglPb09OTz+XQ6bbVad3Z28vm8x+NZXFxsdrKtqkEH6m0228TExOnpabVaTaVSoigSQhiGGRoacjqd4XC4VCpdX1+vrq4uLCw0IJ9W53Q6vV6voiiqqmYymd7eXkKI0Wj0+XwsywqCoJ2jN5vNgiA0Ndm/zmAwzM3NVSqVq6urcDgci8UMBgPDMDabbXh4WBvW3t4+ODhY+x0N/8GybCgUuri4uLm5kWU5EAgQQhiG6e/v7+vrm52dvbu7KxaLoiim02mt6RG+4nK5RFFUFIUQkslkurq6CCFGo3FkZGRqasrtdheLxUqlMjMzs7a2hkbc38EzTAAAoGu4fwsAAHQNhRAAAHQNhRAAAHQNhRAAAHQNhRAAAHQNhRAAAHQNhRAAAHTtE57/X6w8vy8RAAAAAElFTkSuQmCC", "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.ρ.real[:, 1, 1], label=\"\", xlabel=\"x\", ylabel=\"ρ\", marker=2)" ], "metadata": {}, "execution_count": 6 }, { "cell_type": "markdown", "source": [ "We can also perform various postprocessing steps:\n", "for instance compute a band structure" ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Computing bands along kpath:\n", " Ξ“ -> X -> W -> K -> Ξ“ -> L -> U -> W -> L -> K and U -> X\n", "\rDiagonalising Hamiltonian kblocks: 4%|β–‹ | ETA: 0:00:03\u001b[K\rDiagonalising Hamiltonian kblocks: 11%|β–ˆβ–Š | ETA: 0:00:02\u001b[K\rDiagonalising Hamiltonian kblocks: 20%|β–ˆβ–ˆβ–ˆβ–Ž | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 30%|β–ˆβ–ˆβ–ˆβ–ˆβ–Š | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 39%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–Ž | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 46%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ– | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 54%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–‹ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 63%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ– | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 76%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ– | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 83%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ– | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 91%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–Œ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 98%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–Š| ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 100%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆ| Time: 0:00:01\u001b[K\n" ] }, { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=210}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOydeUBM6xvHvzOTFkolERJCIWXnSrbIvvws4VqLRMpyhYSSsmUJ15Li2tdkCbkUV9ZcdCltwiU7RaR9mnl/f9RtU812lpmaz1+dM+95nqeac57zvu+zcAghUKJEiRIlSqorXLYNUKJEiRIlSthE6QiVKFGiREm1RukIlShRokRJtUbpCJUoUaJESbVG6QiVKFGiREm1RukIlShRokRJtUZcR/jx48elS5fevHmTVmuUKFGiRIkShlERc5yzs/PNmzfr1q3bq1cvWg1SokSJEiVKmESsGWFQUJCamlrXrl3ptkaJEiVKlChhGNGO8MuXL15eXlu2bGHAGiVKlChRooRhRC+Nzp07183NrX79+iJHnjhx4siRI1xuKec6ffr0oUOHim/Q27ccXV1SqxYAJCZyhEK0aUOEQiQlcVq1UoxqcAKBgMfjsaU9JQWJiZwmTYRaWhxdXYkvj4t76+BQt337T6tXG+jpqf08gM8nz56hTRsOBbZKSHIyp359oq5eePj5M3g86OnRoksoxNOnnNatxf3KEUIIIWW+/BXx+TMIgRi3lHwhFAo5HA6HI/G//u1bTu3apHbtwsNPn8DjoW5dis3LyEBKirBZs/L/BbGxn+fNU2vUKHvnzpq1i0wRm9evOVpa+bq6vKdPOaam5NUrTt26RFNTZqMlJz0d379zGjeW4GEo9RMpMZFjakqK/uFJSZzmzQnlz7bkZE6dOgItLboiN3k8nsgbU4QjPH/+fGpq6tSpU8XR5+3t7eDg0KxZs5In27ZtK/7/4NMnWFlxT54UWlqCEPTpw23ZEnfvCglB797c+HhhvXpiSmKTrKwsLS0ttrSvWcMJCcGQIcTMTGXuXIlfHQID9zx6pBEdHXrs2NkmTTTbtycdOqBDB3ToQOrUAYCwsJdz5pCXL42pN10US5ZwBg2Cg0PhL3XoEOflS+zeTcvrESHo2ZP76pVQW1us8Xw+XyAQ1KhRQ5zBW7Zw6tTB0qWK8WJXBJ/P5/F4KiriBhYUMWMGZ8ECDBtW+PuuXMlp3x7OzhT/+n//zXF1Pff33zY1a9YE8PIlHj/mPHqER4/w6BEnPf1qdvY9dfV3c+YssrKyklT4nDkcB4fc4cNVBw7kXrggPHKEc/8+59IlIfO+8NYtTmAgLl6U4K8n3RPp5UsMGMB9+1ZY4Aizs2FpyU1OFkr+FiGCBQs+t279yNd3IMVyAQBCoVAgEIh+QyWV4uDgoK+vb2xsbGxsXLNmTT09vdmzZ1c02Nzc/PHjx5ULrAShkAwbRjw9Cw+fPCHa2mTt2sLD//2PHDsmtWxGSU9PZ1F7hw5k8GDhtm25Dg7SXP7jxw8/v50RETf4fBITQw4eJPPnk169SO3apEkT0qvXHQ7HGLByc1tBteGiCQsj7dsXH757R+rUIZmZdKnr3ZtcuSLu4Ly8vOzsbDEH9+tH/vxTSqtYJDs7Oy8vT9Krvn8nWlrkx4/iM0ZGJCmJSsMKMDTsCozk8dpbWxNdXdKoUeHz5OxZ8uoV4fP5AQH7g4PPSie8QQMSH/+DELJ2LZk2jQiFZPZsYmVFMjIo/R3E4MMHoqdHhEIJLpHuifT772T69OLDP/8kvXtLIUYEr1+/5nAacjg9xoyxp146IQKBQJwvrQhHmJmZ+fU/bGxsvL29f5T8RpdGRke4bRvp0oUU2eznR3R1SVRU4aG/P5k2TWrZjMKuI9TSIlu2CK5eze7WjUqxQiF59ozMn3+Mw7EE5owaNYlK6WLbYGJC/v67+MzgweTIEbrUubsTLy9xB0vkCOvWJR8+SGkVi0jnCIODyaBBxYdxcaRpUyqtKkJTsyWwl8Npcfky+fSJSslfvhBtbfL9ezoh5OtXoqtLXr8mQiFxcCA2NkTsfztlNG5MXryQYLx0T6SBA8np08WHCxYUT0so5OHDhxyOGeDZoUN/6qWL7QhFTBhr1qyp+x81atTQ0NDQpGctIC4Oq1fj6FEUrS2FhkIgQPv2hYeDBuHKFSh7RlVOdjYyMjB2LGnblsTFQSikTDKHgxYtsHXrr+7u/SdMyD50aDdloiWxYeZM+PsXn7G3x/79dKnr3h2RkdSLTU6GqioMDKiXLJ+EhqJkkMCVKxg0iBZFt24FDRgQGhLiN3AgqN1DiYmBhQUKVgh1dTFtGnbsAIeD3buhr4///Q+5uVSqE0nXrrh/n14VmZm4exf9+hWfuXIFA2lYvOzUqZOn5xw1Nb3Ll49SL118xHetHz58+PbtWyUDpJ4R5uSQdu3IwYPFZ/h8oqFBbG1LDTM1JY8eSSGeaVicER49StTVSX5+fmZmZpMm5Plz6lXw+fysrCzq5YpHairR1SVfvhQe5uYSfX3J3o7FJyWFaGsTgUCsweLPCM+cIUOHymQYW0gxIxQKScOGpb6HAweSs1IuT4qGpltv2zYyZ06x8FevSJ06pOBZmJ9Pxo0jo0YRPp8OzeWzfj1ZuFCC8VL8WUJCSL9+xYdv3hB9fXHvBUkRCIS1ahGanprUzAhLYmBgoC1m5ICELFmCli1RMiInMhLq6hgxotSwgQNx5Qod+qsOoaEwMir82dwcT56wag0N6Olh2DAcOlR4qKqK8eNx5AgtuurWRb16SEigWOyjR+jYkWKZcktUFLS10bx54WFODiIj0bcvqzZJzpMnMDcvPmzSBAMG4I8/AIDHw+HD4PPx66/Iz2fIHgZmhGXm8ZcvY8AAiBcTLTEcDpo2FT5/TotwMWG/1uiVKzh7FgEBpU6GhyM7u9TEHEpHKAYPH6Jz58KfLSyqoCMEMHs2du4sXiQvWB2lcBG4JJaWuHuXYpmPHqFDB4plyi1lnqc3b6JdO9DzOk0jZRwhgEWLsGUL+HwAUFXF6dPIzISDA13fwzJ07oxHjwq108Tly2UXtOlYFy2iefPq7QhTUjBjBg4eREFofhHnzqFxYzRoUOpk376IikJ6OpMGKhivXxd/favkjBCApSW0tHD9euFhx47Q1gZNRXDp2Caszo6Q7ucpHQiFiI+HmVmpk506oWVLBAUVHqqqIjgYL158a9Zshq3t7MzMTFpN0tJC06aIi6NLfnQ0VFRgYlJ4KBDg+nX070+XOgDNm5Nnz2iULxI2HSEhmDED9vZll0p+/MDTp/jf/8qO19BAt26IiGDKPkXj0yfk5havJ5ubIyaGVYNow9ERu0sE69jZ0RUyQ7kj/PwZWVlo0oRKmXLL58949gw9ehSfUURH+PIl6tSBjk7Z866u2LCheGWiZk1MmhTy5o3x+fO6ly9fptsqWldHQ0MxfHjx4f375UxLqMXYuBrPCHfuxIcP8PQse/6vv6ChgXLL0ShXRyvhxAloaqIoqtfUFG/eIDubVZvoYfJkXLuG9+8LD6dMwfnz+P6dekXp6ZFJST0HD54uEAgoEfjPP+jQAZLXZlFILl2CjU1xHPi7d/j0SfH2R39eFy1gyBAIBLh2rfjMwIG9DA3/FArvdunSnW6r6HaEDM/jW7Qg1dER5uXlhYe/8vYulS9RxKVLyMtD9/K+S4MG4c8/GTBQIQkPR4sWxYcqKmjZkvpYD3lAUxPjxmHfvsJDPT307YvgYOoV7dp1TCj0vnv389u3bykR+PffGdV5XdTGhq6AC/p48gQWFuWc53CwcCE2by4+06xZs9evb//yy42HDxvSbRV9jvDrV8TFoWSTIQYcYXWcEebn55ua9hgyxLl7941Fy9AlCQ1Ft25QVS3nIzMzCARgdzVZbomOhqVlqTNVeHXUyQmBgSiap9GUUOjqam9ouEZTs0Xjxo1ll+bnF7Bu3aBTp/qRapAPy+fj2rVSKYOKuC6KimeEACZNQkwMHj8uddLZGTt20G6VuTlevsSPH9RL/vNP9O0Ltf/KDKelIT6+7IOFcho2JOnptPw6YsKCI8zJyUlJyRcIBqmqlvMO8O4dUlNha1vh5TY2ytXR8vn4EaNGlTpTVeNlAFhYoHFjXLpUeDh4MP79F4mJFGvp1Knj8+dX8/O3vnxJwZ0SHn43N/e39PS0rKws2aXJObdvo0WL4sLiBauINjas2iQVMTEVOkI1NTg7Y+vWUifHjMGzZ7S/gNaoAQsL/PMP9ZJDQzFkSPHh1avo1avYL9IEhwNjY7x4Qa+WSmDBEWpoaOrobJ85M2/37rU/fxoWBi63sjdH5TZhucTFQSBAnz6lTlZhRwhg9uzikBkVFUyahIMHqdeipoZJk4qXYWVh2LBVBgZ3AgM9axV0V6nSlFkXffgQjRqhIe1LhhSTk4PXr9GyZYUDnJxw/jzevCk+U6MGHBxK1T+iCTpWRwUChIdj8ODiM4zN41u0AIuroyw4wqAgNG9uFRDgqldeB50zZ6ClVZyB+zM2Nrh5k+maRvLPyZPQ00OZPh9V2xGOG4eoKLx8WXjo4IDDh0FRUEspZs3Cvn0UpG2dOdN061Y/W9uf4qGrIlUgcQJAXBxMTMrfpilAVxdTp2LnzlInnZxw8iTS0ui1rUsXPHhAsczISBgawtCw+Ex4OEP/uJYt2dzzYsERbtoEN7fyPyIEN2+KKEWoo4O2bXH7Nh2mKTAREWjduuxJQ0Pw+fj8mQ2D6EdNDVOmIDCw8NDUFIaGCAujXpGpKUxMcOGCTEJevEBsbDlJQVWSly/x9WupAFEFdYSVbBAWsWAB9uwpFbRcrx4GD6ZlfaIkdMwIy7y+xMeDy61sQkwh1WtGGBYGPr/U1LskcXHg8yvbICxAuTr6MwkJZddFCzAzQ2ws08YwxuzZ2LeveHmAvhrcjo7FHlc6/P1hb0/7XouccPEihg4tDhBNT0dcHCRvAsg+4jjCpk1hY1NYca0IFxfs2EFvoZnmzZGRgQ8fqJTJWIX0n6lejtDXF4sXV5hHFRqK/PzyH+glUTrCMgiF+PIFY8eW81FVLbRWQPPmaNcOZ88WHk6YgPBwpKZSr2jsWDx6JP2NmpeHw4cxYwalNskrqamphw6FDxyYV3QmPBw9eijkS4A4jhDA4sXFFdcK6N4durr0PqY4HHTqlPvXX5R93S9dinrzJrZr1+IzTM7jq9HS6MOHSErC+PEVDjh9GqamENnoqUsXfPhQnE+t5Pp18HjlZztV7W1CALNnFwcmaGjk8XjDmjfvHBpK8RNITQ2TJ0s/3QwORrt2DC0xsQshpEOHAVFRl44fn1t0UkHXRVFxEmEZOnVCixY4darUyTlz6M2jSE1NvXevu6Pj/06fDpFd2qVLf9rarsjOnvPPP4UbjwxXSG/UCCxmUDDqCDdswKJFFe485+cjOhqjR4uWw+WiXz/lpLCYs2cr7G9XhVMJCxgxAv/+W7j8+/79+/z8/PT05SdPUl92wdFR+pCZgADMmkW1QfLK169CLldLKCyeEV69qpCOMDUVublo1EiswYsWwde3VMPUiRMRFYWkJJqsw6dPnwD97OyBd+/Gyy7twwd+To66mpoa/7+vOMMV0jkcNGuGf/9lSF0ZmHOE//6LiIjKVociI8HliuUIoVwdLc3duxW+t7Zti4QEhoris4KKCmbMKNzAa9q06dy51rq6Z7p0mSvqOokxNYWpKc6fl/jCxEQkJZXtKVZV8ffnNGly8dChdidOFEZSJiZCKISpKbt2SUMlGYQ/U1Bx7a+/is+oqcHevlRRXGoxMzPbuXO6paXakycuMlZo4POxf/8IR8f5ly97W/6XPM/8PJ7F1VHmHOGGDXByqmzZMzi4wvW9nxk0COHhtMTKKyLPn2PAgPI/0tKCvj5r71nMMGsWjh9HQcV/H58lx44d3r27OR2+X7qQmYAATJ9eTinBqsf9+/D2xoULhhMnjqlZs2bBySq/LloAh4PffsPGjfz8Em0JnZxw6BDo60UxefL4iIglP35olcnfkJTffoOeHnbt6tO9RGVL5v9xLMbLMOQIP3/GqVNwcalsTGgoevQQtx6xgQEMDfHwISXWKTZZWcjIKD9SpoAqvzraoAGsrHD8eOHhoEHQ0sKZM9QrGjsWjx9Ldq9mZ+PIETg4UG+MvPH1KyZMwK5dZZOAFdoRij8jBNCxY9zVq5aNGv3y4r8SKUZG6NWLrsbRBaio4OhReHsjOlpKCcePIywMhw6VevayUiG96jvCbdvw66/Q169wwI8fSE7GxIkSyBw0CPR3O1EAzp2DmlplOxlVPl4GgJNTqaRmT094eVG/IKyqKnHIzKlT6NIFzZpRbIm8IRRi0iRMmFB2ayM3F3fvKl5L+gIkWhoFEBf3mMPp8+VLtycl7jcXF/z+O2gtLmtsjM2bMXEipKjcFxuLefNw8mTZvUBWKqRX8aXRzEzs3YvffqtsTHg4OJwK8wvLRblNWEBoKIyMKhtQHRyhjQ2ysooLbQwZglq1cO4c9YokDZmpJmEyXl7IzYWPT9nzt26hbVvo6rJhk2wIhUhIKNuPt3LGjh3j7FxLT69hTk5xLp61NXg83LhBvYUlmTIFHTvC1VWyqzIyMG4cNm0qp1M0K/P4Kj4j3L0b/fpVVjUNwLFjqF8f9epJINbKCvHx+PpVRusUnocP0blzZQOqdiphARxO2QKPHh5YuZL6SaFEITMJCXj5slQJ4ypJeDj27cPRo2Ur/EGR10VfvIC+PmrXluASdXX1rVu9LlxY7upao2R9tdmzmehH4e+Pa9cQFCTueEJgbw9ra0ybVvYjtiqkGxri2zd2Mihod4R8Pn7/HYsWiRh244bEJQxUVdGzZ6nGmNWT16/Lb2JcRMuWePuWxh17OWH6dJw7V/xiNGwYNDSkCfIUifghM7t2YebMKh4m8/o1pk7FiRPldzBXXEco6QZhEV27YtQoLF1afGbqVFy/juRkqkwrH01NHD2KuXPx+rVY4zdvxuvXpfopFsFWhfSCHhSsRPbR7giPHkWrViI2Xd+9w7dvmDJFYuHK1dFPn5Cbi5EjKxujogITk6rZobckenoYOhSHDxefWbECq1ZRvz0jZshMdjZOnMD06RRrl5SoqKiScYzUkpuLMWPg7l5++bSPH/HuHTp1okk5vUjtCAGsXYs//yxOpdDULFUUlz66dMH8+Zg8WXQ4fWQkNm3CyZPll/th8fWFrdVReh0hIZWV2C7i7FlwOOW3pK+cgniZatDltEJOnICmJkR29akOq6P4rzFT0fdhxAioqMhaLPtnCkJmRDZmOnECv/yCJk0o1l7EuXPnQkLKLymSno7Hj3HhAoyMbLt0WaipSVdJm7lz0bw55s0r/9PLl2FjU856qUIgiyOsXRv+/nByQk5O4RkXF+zdW3xIH0uXQlUV69dXNubTJ4wbhwMH0LRp+QNYdIRsxcvQ6wgvXICGBqytRQw7eRJt20qzgtSiBdTVERcnnXVVgfBwtGghelh1iJcB0KMHMjL2Dxq08P1/9feWL4eXF/WvSo6O2L9fRMgMrWEy27dvHz167f/+t97I6GzLlmjQADo60NCAigo4HGhro1MnjB2Lt2+zCJmVm1vL2FiweTOonRkePYqbNyub6ISFKeq6KCRMIvyZoUNhZlbskFq0QMeOEmzgSQ2Xi4MHsWMHIiPLH5Cfj/Hj4ehY4VYUuxXSq+aMcMMGLFkiYgwhiIqqLA2ucqr56mh0tFgz6SqfSljA69ev09IOhYWZr1ixqeDMyJEQChEaSrEikSEzMTF4/16yKGiJ+PtvASH1AW0trW+9esHODt7eOHwYkZH49AmEQCBAbi7+/turW7fDXl6/dejA8/CAhgZ696YmgjEmBgsX4syZCsNJhEJcu4b+/SnQxTzZ2Xj7VtbasDt3wt8f8f+VP3NxwbZtspsmmkaN8McfmDy5VGeoItzdUasWli+v8PKrV9mskM5a4CihDnNz88ePHxcd3rtHjI1Jfr6Iq2JiCI9H/v1XSqUhIcTGRspraSI9PZ0xXSoqJDy87Mn8/PzMzMySZ969I/r61Gjk8/lZWVnUyKKa7OzsVq168njdPDxOF50MDiYdOxKhkGJdR4+SAQOKD/Py8rKzs4sOZ88m3t4UaywiIoLweMTK6o8tW7ZIdOH58+SXXwiXS2rWJJMnk48fJdObnZ2dl5dHCElLI82bk8OHKxv84AFp21Yy+bJD1a334AFp354C4Tt3kp49C797QiExMSH37lFhnxg4ORFb28KfiywPCSFNmpCUlMoudHQk27bRbFxphEJhRkZGwc+vX5OGDakULhAICr60lUOjIxw+nOzeLfoqV1eioyO90owMoqVF/vszygWMOcLYWMLhlPOq8bMjJITo60v84CsXeXaEBRw5kmVmRvj8wkOhkFhYkNBQirXk5pL69cmzZ4WHJR3hjx+kTh3y9i3FGguIiyOqqmT0aOkl5OYSX19iaEg4HGJoSHx9SUTEDTc3t7S0tMovLHCEQiEZPZosWCBCy+rVxNVVeiOlg6pb748/yNSpFAgXCEiPHiQgoPBw0yYyebLMxolHdjaxsCAHDxLyn+VJSaR+ffLggYgLmzUjiYn021eCko5QKCQ1a1L5PGfZESYkEAMDIs4Ds1UrMmSITHr79qX+MScLjDlCDw9St24558t1hH36lDN3lAL5d4SEkD59Sr2EBQWRTp2onxS6uhJ398KfSzrCwEDyv/9RrKuAf/8l6uqULYFER5Nhw4iaWhZgBizU1R0aGVnZ+AJHuHYt+eUXkpsrQnjPniQsjBo7xYeqW2/BArJxIzXCExJIvXqFb0VpaURXl5r3UXGIjSX6+iQxkaSnp2dlkQ4dyJ49Ii5JTCSNGzNiXAlKOkJCiJkZiY6mTLiYjpCuPUJfX8ydCw0NEcP4fDx/jqlTZdJVbbcJIyLQurW4g6tJvEwBO3bA07M4p3DsWOTlUf8lmTmz/JAZmsJkUlPRvj3MzBAWRo1ACwtcuIDUVAGX+wNI//GjUY8e4HKhp4cuXeDkhHPnygY6RkRwduxAcHCFzdQK+PED0dEK2ZK+AFlCRsvQqhUcHbFwIQDo6MDWFnv2UCNZJGZmWLkSkyYhLw9OTjAzE13z9soVGje2xYSVwFFaHOG7dwgJEetZcP06CMGwYTKpq7aOMCEBffqIO7haOUIzM4wahdWrCw85HCxfjlWrKNZSEDJTJoXh8WOkpFBflSMrC23aoF493L9PsWRNTc1Hjy5s2tTqx49tAgFiYuDqCj09hIZi/HhoaEBTE61bw8xsnYaG6YABPWfOvC6yRd+1a+jeXfR7sNxCoSMEsGIFnjwpDK1ydkZAAMXhu5Xg7Izo6I5163Y5dMgxIED0+CtXKuxjwxisxMvQ4gj9/DB9OvT0RI/cvx+NG4tOg6ucdu2Qnl7FOw39jFCIr18liLatJqmERaxdi6NHi1NrbG2RkUHZXKqIn6vM7NqFWbMoTp7Lz0fr1lBRwZMntNRBtrCwcHV1VVdXB9C2LZYtw+XLeP0aubn4+BGbN6NdOzx7FgTMBtr5+l4dNAiurtizB3fuoGQtsSIUt6AMgJQU5OdTWVdFTQ3+/pg3Dz9+wMIC9etHT5u29dOnT5QpqJiYGOTnfwQuEnLF0hITJsDHB8HBiI8vZyUjNxd37qBfPwbsqowq4gi/f8ehQxXm2JYhIoKCG4bDwYAB1D/j5Jzr18HhSJDqVNCht/p0cKxTB0uXYsGCwkMuF8uWwcuLYi1jxyI6uvi+zcjA6dOws6NYS5cu+PED8fFQV6dYskjq18esWThxAhERO3m8gBo1rty8OX/ePDRogHv34OqKZs3QoAH694eLC/z9cf06lizx27evv6bmX6KlyyXR0TJlEJZL797o1w8eHiCEPH8+5fjxWra2zhTrKE10NMaMwaBBMDW15HD6t29vsn8/RoxATg6OH8eoUdDWRps2GDsWHh44cQLXr6d06WLL4Uzh8dio9VkCdnLqKduU/C9YZs0aYmcn1vj0dMLlkocPKVB97Bhd4QlSwEywjLMzMTQs/6Nyg2UIIcbG5OlTWfUqRLBMAXw+aduWXLxYeCgQkLZtqYkYKomrK1m6tDBYZtcuMnYsxfJ79CAaGuTVK4rFSkFR+kQZ3rwhYWFk61bi6Eh69hRwOO2AZ507D2XeQkpuPT8/Mncu9cK/fCENGpC7d4mpaQ81tfENGiz8/l0WeRXy6BEZNYo0aED8/EjBY6Bcy3NzSUwMOXmSrFxJxo0jjRrtBTarqrqfO3eOFrMqpkywTHIyadSIMuHsRI3evx/dsCGJjRVrfEAAUVOjJpYvJYXo6BAxfl8mYMYRduhQYbRtRY7wf/8jwcGy6lUgR0gICQ8nzZuTnJzCwyNHiKUlxSoSE4mBAcnMzMvOzu7QgWJHa2tLatQgT55QKVNqKnKEZbC3dzU07HrkiMxfNcmh5NaztyeBgbQIP3aMmJuTzMy82Ni4hQtJx47k0ycZRZYiJobY2pJ69cj69aTkA0Acy5OTk9u06du+/cCUytMMaaCMIxQIqMygYCdqNCSkTufO4jbxOn4c5ubitqSvnLp10bJlhVWFqiTPn0u8rV2t4mUK6N8frVtj+/bCwwkT8OVLcTVkSjA1RatWuHCB++AB9/t30QUFxcfJCWfOICICbdtSJpMB9u3b9ObN35MmjWHbECmRsbhaJfz6K5o2xfbtNczM2mzejHHj0L07NfthT55g3DjY2KBTJ7x6BTc31KwpmQQjI6O4uL8ePbpct25dCgySAS4XzZrhxQtmlVIoixDu4cP6ImuqFSFLZbWfGTSoGsWOZmcXNtWUiGpSaK0MW7fC1xcfPgAAj4fly+HpSbEKR0fs3cvdu5c3axZlwSweHggMREgILC2pEahEHAQCJCSgTRu65O/YAV/fR66ufh8/fnRzw6JF6NMH0dHSiOLz+aNHz2zdevDgwU8LXODLl3BzU+Bg3SKYj5eh0hGmpJgJBAd69BBrcHIyMjKobFIzcCAuX6ZMmpxz9ixUVcvvAFcJ1XBGCKB5c9jZFTu/iRORmoqICCpVDByYcf26zZEjLjY2b2WXNnfuUh2dfmvXPo6nqKAAACAASURBVDl4UESnSSWU8+IFDAygpUWX/MaNCZ8/bcsW7fHjXQA4OcHPD4MG4c4diUVduPAkNDQvMXFGXt7xgllgFXCBBTAfL0Pt0mjHb992ZWdnizPU3x916kBfnzLdv/yC5OTCF/8qT2ioNP19WrbEhw/IyKDBIPnG0xOXLuHBAwDg8eDuDm9vKuVv375ZIDAWCg03blwmo6i8vLydO09+/75MVdVu8mRKrFMiATExVGYQ/gyHwzE01ObxLquqNis4M24cjh3D2LESlIZ//x6zZmHWLLPGjXNMTHZs2TKa+XBiWlHsGWGNGhfNzAw1xHstuXhRmgaElcDjoU+f6tKw/uFDdO4s8VU8HkxNi8vhVx+0tLBqFRYsKOzHNHky3r0THjwYKxQKKZE/atSoGjWu8Hh7fv3VVkZRKiqqgCmH81vv3mIXDVJCHdSm0leg4vqhQ2uSkjZmZRWe6dsXFy5g5kzRTS4zMuDrCwsLqKsjKUnt+fOTT59GWNC0pckeiu0IdXW/BAauFWekn59/fLzTgAEUz37bt0/auNE3KSmJWrFySHKylItm1XN1FMD06eDzcfIkgIJs97GOjpusrWX1WwVYWFhkZj7/9i12+PDhMopyckLNmpdTU+9evnyEEtuUSAQDjlBFReXXX01++QWbNhWf7NwZt25h7Vps3Fj+VXl5CAxEy5aIisLDh9i2Dbq69NrJIoq+NCoWKSkpixb9TsiAjRspXvrZu9cuJsZ4yBA7asXKGx8+IC8PI0dKc221dYRcLrZuxZIlyMwEAEI+5+X1T0igsrqHioqKjBKePcPevQgMRJ06mpSYpERSGHCEBWzejO3bkZxcfKZ5c9y6hSNHMH9+qVbSQiFOnULr1jh1CpcvIyiows7yVQZDQ6SlFd6qzMCCI6xVqxZQB/Bu3Lg+tZJbtjRSUQmuW9eIWrHyRlAQNDWlrEtX3QqtlcTSEj16FL5xX7ly2M0tFTh8/TrbZpVg8GB07IiJE9m2o7qSmYkPH2TtxysmhoZwcYGbW6mTDRogIgJRUbC3x6dPXwQCwdWr6NwZmzdj3z6Eh6NdOyZsYx0uF02bMlo1U9Z3WCm4cKEmh3PD3/+CvT3FIXFhYcenT39hYtKcWrHyRlgYWrSQ8trqmUFRxMaN6NAB9vZo1qzZ+vULBgzAr7/i7l0YG7NtGfD770hOliZ6UAlVxMXB1JTiOrGVsGQJzMwQEVGqdL6uLq5cQdeuB4yM9vN4OU2a3Fi7Vn3UKIZMkh8KVkeZmZ2D+RkhIXB2hq2tyqxZo1Qrb+UiORwOx9a2RXg4FSn6ckxMjPRxRgYG4HDw8SOlBskT8fHx586dy8vLK/fTgtfwpUsLD62tsXQpRo9GUdgCW6SnY/FiuLmhPsWrJEokgLF10QI0NLBxI1xcyjajqFULHTs+5vMncLm88PBv1dALgvF4GaYd4dq1+PEDe/fSJb9PHzx8WMUzBD5+xOjR0l9ehSeFKSkpPXtOGTfuLwcH34oCQpcswb17uHmz8HDBAnTqhJkzGbOxfIYNQ716xX2jlLACw44QwJgxaNSobAMTAH5+yxcs+HTwoJuhoQGjBskNVdkR8vlYswbz5kGTtlCAmjXRtStu3KBLPus8eQKBAL17Sy+hCsfLxMerfPsmAL5duaJerx7Gj8fevaXiEQBoaGD9eixYUNyIY+dOPHuGzZuZt7eQy5dx5w7OnWPNACUFxMTQVVytErZswapVSE0tdVJfX9/Pz2vMGKki4qoELVowGjjKqCOcMwcqKli/nl4tNjZVuSXTqVPQ04Ms8YlV1RHevIlx43R37Qr96y/HDx9c4+MxdiwePICVFZo3x6xZOHUK374BwPjxEAj+at/eLjz8OgB1dZw7h61b8eefLJgtFGLiRIwfj06dWNCupCSxsUzPCAG0aYMJE7ByJdN65ZyWLZlNJaSmxDch5L82TBV9mpZGVFSIvz+FCsvn0SPSqhXtWiqB1u4TPXuSnj1FjKmo+0QBo0e716jRaeXKzdIZIJ/dJy5cIPr65MqV8j998YIEBBBbW6KjQzp1Im5upFYtMyBWW9usaMydO8TAgDx/LqUBBW2YpLhwyhRSqxbh86XUyyRidp9gEalvvfz8fEvLsVxu59OnL1AuXCRfvxIDA1Lxs1NWmOmHIzVluk8UIBAQDQ0KelCw032iEn79Ffr6mD2bdkXt2uH7d7x6RbsiVkhIkGldFMCdO5f4/LOHDgVRZBH7HDoER0dculRhOw5jYzg6IigI799j3ToIhcjOzgPs09NtijYLLS2xfDlGj2Y0eyk2FkeP4o8/ZJriK5Gd1NTUuLgUoXDdkSPnmdeuqwsvL7i4lEofrOYwnEHBkCN89gxhYTTGyJSEw0G/frh6lQldDCMU4ssX2MpWDmXTphXq6nMXL6a02iZ7bNkCT09cvy5WzTkNDdjYYMMGrFw5R09PqKHRc/RonDhR+KmLC7p2xdSpzD2Phg1Dt24YP54hdUoqon79+p07D23UKGD16gWsGDBzJnJzEVR13k4pgMnVUYYc4fjxMDPDkCHMaKuy24QREeByZd3Pnzx57LBh53R1JWxmKJf4+mLPHty6BVNTyS709FyQmvrQ33+0qipcXYvrWu3ahc+f4etLuaXlsHYt3r1TxsjICw0bLvbxOdWGvg5MlVJQ+WjxYkYXJOQcJuNlmHCEN28iOhrHjzOgqpABA3DtWnFYYJXh9GmJWy+VSxXIoBAIMHs2goNx4wYaN5ZSyNSpWLYMamo4cADz50MoRI0aOHUKu3ZJ0ApAOr5+hZcXVq9GvXr0KlIiJtHRLISMlsTSEj17MvQSphAwmUHBhCOcMgXW1uK2racEAwM0boyHD5nTyAyRkdTcq+3aKbYjzMvDxIlISsK1a7J28nJxweTJ4HLx4AHGj0dODgwMcPIkZsyg9210yBAYGpatsKWELfh8PHtGYz9eMdmwAf7+ePmSZTPkhCq1NHr4MN69w9GjdOspy4ABVbBh/fPnFcaDSISFhQI7wsxMjBiB3FxcuoTatSkQ6O2NAQMgFEIoxJAh+P4d3btj/XqMGIH0dArk/8yZM3jwABcv0iJciRQkJKBZM/Yb2zZqhN9+w6JFLJshJ1SdpVFCMG8epkxhYf1nwACEhzOtlFayspCRIWukTAFNmyI9HV++UCCKYdLSYGMDAwMEB4PCZqSbNsHMDOnpMDWFlRXevoWdHXr3piVwJj8fdnaYNo39+YeSIlhJpS+XRYvw5EkVfImXAiMjfPkC8Rq9ywq9jnDVKmRlwd+fViXlY2WFmBh8/86Capo4dw6qqmjYkAJRHA7atlWwtPpv3779+29W797o1Qv791Ocb8DhIDAQOjr4/BnTp6N7d8TEYPt2vH//wtTUztt7C4W6Jk4El8tQBLUSMZEfR6iqio0bsXAh+Hy2TWGbggyKFy8Y0UWf6NxcrF8Pd3cq39zFR10dlpaQqyY7MhIaiiZNKJOmWNuEd+7cad58oIlJj+HDX6xfDw4NZdV5PBw+jLQ0PHuGdeswYAAiI2FktO3ZsxEbN156/fq17CoSExO7dRsfHHzj6FFwWWiApqRCWI+UKcnIkTAyEi5efOdjFa6OLx6MrY7SeDvOnAkNDTZLBw0YUHWSKPLz82/fftCxY/lNFaRAsbYJExKSvn2zrFGjlbX1K/q0qKsjJAQPHyI+HseOwdYWzZoN0tdfJxCo+vg0kP0NvU+fKffv9yJkVrduqaJHK2GQmBj5avXXpInv778HWFgMTKdpm1pBYCxehi5HmJqKY8fg50fLy7uYDBiAy5dZ004tpqY9X7/ee+WKFVUCFcsR1q49qX79pps29be2tqZVkZYWwsJw6RIePsTVqzh8WA+oaWPTPCVFpV8/fP4sveRLl5CSMgQ4p6KSUZuSIB8lFJGSgrw8GBqybUcJ1NSyeLx6WVm8/DItmqoZjGVQ0OUIx49H48awt6dJvFiYmSE/n6ElZrp58+YT0Ck9nbKlEgsLxMcrRqqlQIBVq1QPHJjv7DyDQ/+LlY4OLl1CQAAiI9Gt25GUlFVXrvzr7v7a2hqWloiNlUamuzuGD8fIkatu3Vr5/XsS5Z04lciCXK2LFrBhwwpv7266ukd0dOqwbQubKLYjTEhARAQOHKBDtmT0719F4q8GDtzM4bz08ppPlcBatWBgwGijE6nZtw/6+tTkjYhJw4YIC4OPD7p2ndaqlXfLlk2HDGkMYNUq9OuHs2clEJWTg+7dsWkT9u7FmTOwsrKqWbMmXXYrkQo5dIRqamru7rb6+m2qUpSDFBT0qWcAWhyhrS06dJC1NjQlVJkkiqysUV26rFuxwpVCmQoRL5OTAx8f2lt3/Uzz5vjzT/z+e5PWrXutWGH999/c+/fh4wNvbyxciKVLUVHj35I8eIB69fDyJZKSWF4dUVIJT57InSMswMEBe/awbQSrGBkhNZWJDArqHeHVq0hIYCGDvlxsbBARURUCkZ88AeW7YwpRaO3339GtG375hQXVbduia9f1Z8/qT5/ul5n5/tIlbNyI9evRpQtu3MCECSLKQm7ejO7d0bEj3r5Fs2ZMGa1EcuRwRljApEm4cgUpKWzbwR5cLpo0YaIHBfWOcOpUDB4scRFkmtDTQ4sWuHePbTtkIz8fqan49VeKxcr/jDAtDZs2Yc0a1gwYObKDru5JHo/fp4+urS3U1PDkCdq0wYsXeP0alpbld/sqqFCzZAnWr0dEhLLFklyTn4+kJEYLQIqPtjZGjMDhw2zbwSrMrI5S7AiDg/VSUuTrP1cFVkevXAGPR/1Lq4UFoqMplkkt69ZhzBiYmLBmgIPD5OfPz6Wl3X39WsPGBt7eaN0aOTkICkLt2khLQ+fOiIgodUlSEurXR2QkoqKUtbIUgMREGBlBbvdtZ87Enj3Vuk8hM/EyFDvCLVsaODhAV5daqTJRBbIJT5+mJbbb2BjfviEtjXrJlPD+Pfbtw/LlLJtRp04dHo+npQVHR9y+XZiQY2sLQjBxIng8DB5cvIV58CDMzNCoET58QPv2ZUVFR8dYW/+6cSMblZaUVIDcrosWYGUFLhd377JtB3soniN88aJPVpbb9u0UiqQAS0skJSlkXc0i7t5F167Ui5XzQmuenpg1S76yuwCYmWH9erx5A0dHREVBRQWtWmHFiqeqqi00NHrY2791dcXjx1BXx48fiI/H5cv44w94eWH6dFhb+16/PtvHZ39GRgbbv4eSQuQ2UqaIGTOqdciM4i2NCoV2hISmpcnX3m6NGrCywrVrbNshA69eYeRIWiTLbVp9YiIuXMDixWzbUQHq6rC1RXg4QkNhaQnAFhgN2KuoTI2JQdu20NFBw4awtcWWLbh7FxwOrKwwY8ZIXd2FfH7zwEBNtn8DJYXI+YwQwLRpCAmR35Ubunn+/OLt20MOHDhJqxZq9/F3AdlfvmjL2CKOcgq2CceNY9sOqUhORl4ejY7w0SNaJMvIsmVYvBg6OmzbIYr27bFzJ9TVh/j5HQU4zZtbOzujSRMYGpZr/LgNG8a9fYuhQ5GcDD8/8HjMm6ykFPJWXO1n9PQwcCCOH8ecOWybwga+vt55eRfc3QfZ2Y2nTwuVM8ImTe64uBzu00e1TPgA6yh0rbWjR6GtjVq1aBEunzPC+/fx4AGcndm2Q2w2b17v57dw+XK7uLh9Q4cWzggrwtAQd+7g2TOMGYOsLAatVPITX74gKwuNG7NthyhmzkRAANtGsMT06RNr1BhmbT2KVi1UOkJV1RoODrWPHMH48Th4kELBsmJigho1kJDAth1SERaG1q3pEm5hgbg4sXLDmWTpUvj4sN8lVSJcXFxWrFjGFa+phKYmzp9H/fro21em+qVKZKRgXZTFeshiYm2NrCw8fMi2HWzg6blgwoQH/ft70qqF+jzC/v1x+zbWrsX8+XL0hLWxUdTY0dhY9O9Pl3AtLejrM1TNT0xCQ/HhAyZPZtsOmlFRQUAARo9G9+54+pRta6or8r9BWACHg+nTq2/IjLk57TF9tJRYa9kSd+8iOhrjxsnL4o+CZhPy+fj6lfpU+pLIVVq9UIgVK+DrW12S0N3c4OGB3r1x6xbbplRL5Kcfr0js7XHqFH78YNsONrCwUExHCEBPD1euoFYt9OiBt29pUiIB/frh1i3k5rJth4RcvAgVFRqXRiFn24RHj0JNDcOHs20Hg9jZ4ehRjB2LEyfYNqX6If+RMkUYGKB3bwQFsW0HG5ib0176g8bGvGpqOHgQdnbo3h1RUfTpEQsdHZiZ4c4dls2QlLNnYWRErwr5cYR5efDyAk0N6OWZfv1w7Rrc3ODlxbYp1Yn8fCQmymlxtXKZORN797JtBBs0bAhC6N1Np9ERFjB/PrZsweDBCAmhW5UIFHF1NDIS3brRq0J+Cq35+8PMDH36sG0HG7Rti8hIhITI18561SYpCYaGdMVj08GgQXj/Xl7uVoZp25be93XaHSGAsWNx8SKcnFh+4R0wQPF6E75+jVH0hg2jeXN8+YJv3+jVIpKMDPj6wseHZTNYpGFD3LyJpCT88kuIra1LQkIi2xZVcRQlUqYILhd2dti3j2072IDueBkmHCGArl3x4AEuXMDMmaw1RerWDcnJ+EhZj3faefECfD6GDaNXC5cLMzMpG69TyIYNGDhQYTZsaEJLC2fOCGJiPIKDR8yYsYxtc6o48l9c7WccHHD0qLxEIDJJFXGEABo1wo0b+PwZgwezM//g8dC3ryLVWjtyBDo6UFenXRHr24QpKdi1C570ZgopBhoavE6dmqmre/340ac69xxgAIWbEQJo3BhduuDsWbbtYBy6A0eZc4QANDVx5gw6d4ahoV2DBt1On77ApHYoWjbh1asM7eSz7ghXrcK0acrutYXcuRPy8eO1unXnrVjBtilVGgUKGS1JQWOm6kbbtkhIgEBAl3zRjjAtLe348eMrV67cvn37R5kXFnk8uLt/53CSPn7cFhh4SkZpkjJwIMLCFKa5V1wcbGyYUMRivExWVpa//+ljx5Ld3NgxQD7R1tYIDkZQEAID2TalivL1K9LT0aQJ23ZIzogRePZMUetkSU2tWjAwwIsXdMkX7QgnTJhw8uRJHo93//59U1PT+Ph4GVVqa2vPn/8/HZ01FhYLZBQlKU2bQktLfhsPlSQ3F9++YcIEJnRZWCA2lp1gxalTF8yd+zg/f4SeHm0ve4qJnh4uXICHB/76i21TqiKKUlztZ1RUMGUK9u9n2w7GMTenceFKtCMMDg4+d+6cp6fn4cOHbWxs/vjjD9m1rl695ObNC4cPd2S+UELBpFD+CQlBjRoMNWfX1oaeHv79lwldZcjJURMKv2tpcTiK+EyimVatEBSESZOQlMS2KVUOBaop8zOOjjh4UPHKg8gIrfEyoh2hlpZW0c8cDkdNTY0SxebmsLHBli2UCJMARdkmPHeO0XUbtgqtaWj4OTqOjI6+Kma56upG795YswaDByNFvrp8KjwK7QiNjWFmhgtMh1iwDK2OUIKSjrdu3QoPD9+0aVNFA7Kzs9etW6dfuhvhyJEjraysyh2/YgXH0lLV3j5PX5+5Xbvu3TF5svrXrzk1a9KlIicnp0aNGjIKuXdPrXt3YU6ONLkmAoEgNzdXItfSpo3KP/9gyJB8kSPz8/P5fD4lE7j4eM7Nm6pxcT00NZGTkyO7QBbh8/kCerbyJ07Ekycqo0ZxQ0PzKHoLlZKcnBwej0fTr0kJ4t96jx+rTpuWn5MjwX4AJfc1VdjZ8QICeMOG5YkzWK4s/xlCSMFXq/JhpqacmBjVnBzJJsJCoVCkZIjvCGNjY8eNG3fgwIEmFc9TeDxeo0aNGjZsWPKkjo5ORXYYG2PCBOLnV2PDBuZuLR0ddOhAIiNVBgygy/vyeDxx/vSV8/YtZ/RoIrUcSW1o145z4gRHnEsIIWJ+t0Sydi1v0SKhtnZV6E4rFAoBUPJn+Zn168n48Zg7V/WPP9h0Qrz/YNGGyhHTPIEAT59yzc3F+sJLKpwZxozB4sXc5GQVY2PRzzG5svxnCCHiWGhigo8fOdnZPE1NCYRzOOJtuxAxSExMbNiw4bFjxyofZm5u/vjxY3EEFvH5M9HTI//+K9FFsrJ6NVm4kEb56enpMkqIjyccDsnNlfLy/Pz8zMxMiS5JTCTGxmKN5PP5WVlZ0phVmthY0rAhkdBM+SUvLy87O5s++VlZpGtXsnYtfRpEk52dnZeXx6YFohDz1ktIIC1a0CWcMRYsIB4eYo2UN8vLIBQKMzIyxBnZqRP5+2/JhAsEAnG+tKJXz54/f96/f38vL69faegGpK+POXOYLqwl//Eyx46hTh2oqjKnsWVLfP6M9HTmNK5ciUWLQN8CdRVDQwPnzmH3bhw/zrYpoggKOmdr6xQbG8e2IRWiiKn0PzNrFvbvpyy1Lisry9l5+aJF3ny2Sn+JAX2Bo6Idob29fVZWVlBQkI2NjY2NjRfVBUMXL8bly4hj8K7p2BGfPslFc6iKuHYN5uaMauRy0aYNc4XW4uIQGYlZsxhSVzVo0AAhIZg3D/fusW1KxQiFwtmzPYKDbYcMWRYVJac5uwodKVNEs2a5GRkTmzWzfkJFDMnJk0EBAdi27eOSJVc+fZJdHi3QFy8jeo9w8+bN6SVmCvXq1aPWAi0tuLrCwwNnzlAruEK4XPTrh6tXYWfHkEZJiY8H8wnmBWn1lpZM6PL0xJIlyumgxLRvj0OHMHYs7tyR02RwLperpmairr7MyGjilCn4+hWDBmHwYNjYoE4dto37j5gYTJ/OthEyExcXl59f480b5wMHTm/eLOuLs4lJJ0J21q6NZ8/cWreGsTEGDcKQIejWDfKzvWhujosXaZEs2hF27dqVFs0lcHbG778jMhLdu9OtqpCCJAr5dIRZWUhPZyiVviQMtIEuIDYW9+7hyBEmdFU9Bg+GqyuGDMHdu9DWZtuan3j+HPn5p1++zDUwUAPw6hXCwnD6NJycYGyM/v0xbBgsLcFusoyCFlcrg7m5ubW1yuXLgd27/y67tEuXzO3t7wYEcHk8nkCAx49x9SqWLEFMDHr3xvDhGDoUjRrJrkcmaCwGKdnOY6VIESxTxN69pFcvCm0Rwdu3pG5dIhDQIlzGrenDh4mamkwGSBEsQwiJiCCWlqKHyR4sM2oU2bZNFgHyCN3BMmVwcPimrz9t6FC779+/M6ZUnGCZyZOJj0855/l8cusWcXMjnToRfX1ia0u8vR+1aWNjb79QKBRSZaE4t15aGtHSIlLolM+Qkx07yMiRIsaItDwlpcKgxc+fSVAQmTKF6OqSNm2ImxuZMmWzqWnvkJA/pTW5LOIHyxBC6tUj799LIFzMYBl5cYT5+aRNG3L5MoXmiKBly6f790fSIVnGG8bWlrRqJZMB0jlCMR8QMjrCmBjSsCGhIuxUvmDYER46dJTLXcXjrTxy5DhjSkU6wvh4Ur8+Efn1f/mS7NpFGjVyAf6qXXvIq1evqLJQnFtPzBc+6YQzT04OadSIREVVNkak5YsWEWdnEYr4fHLzJnF3JyoqFsCbtm1tJLS0QiRyhP36SeYmKIsaZQYeD97ecHNjqOLl06dP376d4uTke+jQMSb0ScKDB8wtEZdERwe6unj5kl4tHh5YuhQaGvRqqfL07durZcu/NDRuHDzYMzubbWv+w8sLrq4oUYqqfJo2hZMTTpz4tVEjj7y8elpahoxYV0jVWBctQk0NCxdi/XrpJXz8iH374O4uYpiKCnr2xNq1mDVrRK1aY42MZkqvUgZoipeRF0cIYMwYaGjgFCMdKfLz81VVVfPyauXlyV2s8Nu3sLVlRzXd/ZgePcKDB3BwoFFFNcHQ0DAxMSIt7bqhYaNBg/D9O9sGAbGxuHEDTk7ijreysnz79vavv+5fu5bRYIyYGKZDsulm9mzcvCl9yPe6dbCzk2D/b8cOn+Tke1FRtqy0rKn6jhDA+vXw8GCihb2ZmdnVq1sbNJhibj6FdmWSEB0NgQD9+7Ojne5+TF5eWLZMOR2kDBUV/PEHunSBlRXevWPZGA8PuLtDoqofADZuxNGjjHYBq2IzQgA1a2LBAvj6SnPt+/c4dgxLlkh2lZ4evLzg4sJCegxNMX3y5Qh790bTpti3jwldnTt3cnEZuH+/fP0Fjh9H3bpgqy4grTPCf/7Bw4dVIWxdruBwsGkTJk9Gz5549ow1M6KicP8+HB0lvlBPD6tWwdmZoUeqUIj4eLRty4QuJnFxQXi4NF1KfHzg4ID69SW+0NERAgEOHpT4QhkxM8PTp8gXXRRZMuTLDQBYtw4+PsjKYkKXnR2Cg5GZyYQuMbl+nc1UX1odoZcXli9XTgdpwc0Nnp7o2xePH7NjgIcHVqyQ8p/r4AChEAcOUGxSuTx7hnr1ULs2E7qYRFMTc+ZIPClMTkZQEBYulEYjl4sdO+DmhtRUaS6XGg0NGBpS/84nd46wUydYWmL7diZ0GRigRw8EBzOhS0wSEzFoEGvaTUzw4QMyMqiX/M8/iI7GjBnUS1ZSgJ0dduzA4MG4dYtp1Q8eID5e+rk+l4uAALi7M/FIrRo1Zcpl/nxcuCBZsNuqVXBxQel2QRLQsSPGj4eHh5SXSw0dhdbkzhECWL0amzbh61cmdM2YASo6DVNDRgZ+/MD48awZwOOhdWtaCq15emLpUrDbRajK87//4dgxjBnDXJGmAtzd4eEh0z/X3Bzjx2P5cupsqoCqt0FYhLY2HB2xYYO4458/x4ULWLBAJqWrV+PiRURGyiREUuiIl5FHR2higlGjUHHfQyoZMgQvXiAhgQldIgkOhpoaGjdm0wY64mWiovDkiXJ3kAn69sWlS3B2Zu717vZtJCdj6lRZ5fj44NIl2h+pVXhGCOC33xAcLG4VZS8vLFgAXV2ZNNaujfXr4eJCWe1vcagujhCAlxcCA5moi62igqlTGdqfEMn58zA2ZtkGOoKyYZ4YUAAAIABJREFUPD2xbJlyOsgQnTvj+nX4+GDjRibUrViBlSspCO+qXRu+vpg1i/o4iJJUbUeopwc7O7GmEElJuHoV8+ZRoHTSJOjoYNcuCkSJCR3PKDl1hA0bYvp0rF3LhK6ZM3HoEBM5GyKJioKVFcs2UD4jjIpCbCzs7amUqaRyWrVCZCQOH8b8+Xj37r2QtioV4eH48IGyurgTJ6JePRofqd+/IzWV/XdNWlm0CIcP48MHEcNWrMCiRaJLH4iJvz98fPD+PTXSRGJsjJQUinvGyakjBODujtOn8fw57YqMjWFigtBQ2hVVDiF49w5jxrBsRrt2iImhMpbdwwPLlzPaW1EJgAYNcP06Tp1a1aLFjHbt+tGkxdMTPj5QEV26X1z8/bF6NV2P1JgYtG3Lcr1vuqlfH5MnY8uWysbExuLOHcyZQ5lSExM4OEicjCg1BT3jqO3cJ79fCl1dzJsHT08mdMlDyExB8zZra9bN+Cszs1u3bqMp6c/58CHi4+W0y0eVR08PxsZPcnKmJiamRkZSPykMDUVGBsaOpVJmy5aYOROLF1Mps4DExMTRo3skJfX78uUL9dLliaVLsW8fUlIqHLBiBdzcKG6C5uGByEj89ReVMiuB8sBR+XWEABYuxK1b+Ocf2hXZ2iIyEm/e0K6oEk6cgL4+lS/X0nH48HmBwDMxMectFTu0K1Yop4Nscvz4FlfXJB8f/3HjuLNm4ccPyiQTglWr4OND/QRrxQrcu4dr1ygWGxb215cv47Ky2jx8+JBi0XJGgwawtcW2beV/GhWFhw8xk+pCoRoa2LEDTk7IzaVYcrlQHi8j145QQwNLluRPn348PPwq3YpsbXH4MK1KRBARIReB3cuXO7Vps19L65emTZvKKOrePSQmYto0KsxSIhWNGzfetGnl0qVWBU+N1q1x9iw1ks+ehUCAkSOpkVYSDQ3s3En9I3XgQFsVlahhw9CnTx8q5col7u4ICEBaWjkfrVghfemDyhk8GK1awc+Pesk/U70cIQBC9sXE3BwzxvfRo0e0KpoxA3v3slA6r4inTzF4MGvaizA1NY2JCa5Rw/PxY44sciZOnNunT+fu3fcqp4PygI4OAgJw7Bjc3TF8uKybcEIhvL3h4wOOTN+RChk0CG3aUJxAdfas/rRph06d2q5WDcKXjYwwfHg5ZUnu3kVCAo2JTNu3w8+P9g42+K8GFoWPa3l3hI0a1atV61l29pfatentxt25M7S1ceMGrUoq5Pt3ZGaymUpfEh4PM2Zgzx6ZhFy6dCM390hc3EmKjFJCAb164dEjdOqEDh0QGCj9cyQoCKqq9L63/f47tm7Fv/9SI00oRGAgZs+mRppCsHw5du4suxju4QFPTxq3KoyM8NtvcHamS34RdetCXZ3KQvPy7gjHjPlfVJR/+/Z/3r9Pe9Tz9OmshcwEBUFdHQ0asKP9Z2bOxIkT0gco5+ejVi1XU9Ol27evpNQuJbKioQEvL4SH448/0KcPnj6VWIJAAG9vrF1L13SwACMjLFxI2SP18mXo66NTJ2qkKQTNm8PGplQuyu3beP2agtIHlbNoEV69wvnz9GoB1dmE8u4IAZiYtNy8uf6yZbRvw06ZgkuXyl9Yp5sLF9CiBQt6K8LAANbWOCZt0+ItW9C27bTExHO9e7OdF6mkPCwsEBmJSZPQsye8vJCXh0yxa88fPQo9PSY6hbm6IjkZISEUiPL3r17TwQKWL4efX3Hp4BUr4OVFezieqip278a8ebQ3M6A2cFQBHCGAXr3Qrh3tlbh1dDBwoPRPf1n45x/07MmC3kpwcoK/vzQXfvyIjRsZKpuuRGq4XDg64u+/cfeuUEenf+PGNlu27BV5lUCANWuwejUDBhY+UufPl/WR+vo1IiPlZd+BSVq3Ru/ehXscYWFUlj6onF690KMH7eVQqI2XUQxHCMDXFxs2gO4UoBkzEBhIr4qfIQQfP7LWlb4irK2Rm4u7dyW+8LffMGsWTExosEkJ1TRrhjNnsni8b2lprosX3+7RA/Pn48gRPH1a/g7iwYPcJk3QuzdD5vXqhZ49ZfW7AQGYMoXitDlFYeVKbNyI7GysXInVq8HjMaTXzw9799Jbw7maOkJTU9jaYs0aerVYWyMzEzQHqJYlMhKEyN2MkMPBrFnYvVuyq27dwt27WLqUHpuU0ICmpmZgoPvkybdiYrx8fWFsjLAwjB4NHR1YWWH+fBw6hLg45ObmLViwdOnSEG9vRs3z9RVs2TKnVauBsVJ1ReHzceAA9WlzioKZGUxMHv3yi/unT0+ZrFpVvz7s7aO7dBkwadJcmir8mZnh+XPk5VEkjlCHubn548ePKRRYhs+fSd265Nkz+jQQQoiPD3F2lklCenq6ROPnzSMNGsiksQz5+fmZmZmyy0lLI7q65PPnUif5fH5WVla54/l8YmFBzpyRXbPikZeXl52dzbYVVJKaSv78k/j4kBEjSIMGpGZNZ2ASMHjnzp1MmhEfH6+paQtcnDFjicjBP996J08Sa2tqLJH0vpYT1NSaAQE8XguG9U6d+htwRVNz5DNRj2yhUJiRkSGFitatyZMnIsYIBIK8vDyRohRmRghAXx8LFtDetMzeHsePIyuLXi0luXkTHTowp058dHQwciQOHhR3/O+/o359jBpFp01KmEJPD4MGYcUKhITg/XvMn2/E4cRxuc8aM9snrEWLFj16qGppbfjwYZwU+R7+/nByosEsxUFTsyaXG6GmxmCfJACAk9OEJk18+HztFy+a0qSCytVRKfxwRdA9IySEZGURIyNy+zatSsiQIeTIEekvl/TNsVYtQu1LNlUzQkLI33+T5s2JQFB8pqIZ4cePRF+fJCZSolbxqHozwp8JDw9/8OABK6qzs0m3bmTNGhHDytx6CQmkQQMixnxALBR0Rpibm7tnzx6qHgiScusWMTAgT59WNkbqGaGPD3F3FzGmCs4IAWhoYNUqLFpEbwkYJmtwf/mCrCyKKxdTSNeu0NHBVTEq3Lm6wsEBpqb026SEJaysrNqxVAZQXR3nziEgAEFBElzl7w8HBwp6JSo0qqqq48ePr8lSsJCVFdatw+DBSE2lXjiFGRQK5ggBTJ2KvDzKSiaWy4gRSExkogMUgJMnUbMm6tVjQpd0zJ4tOo/i9m3cvIllyxgxSEm1xMAA587B2RkPHog1Pjsbx45hxgyazVIiCjs7jBqFMWOoC2z5Dwpz6hXPEXK52LQJixdT/2ctQkUFkyZh/3665JckNBQtWzKhSGomTsTt20hOrnCAQAAXF/j5QVOTQbOUVD86dEBAAEaPFqu21rFj6NEDTZrQb5YSUWzYAG1tKjsgFtC0Kb59w7dvFIhSPEcIoG9fmJggIIBGFQ4O2L8f+fk0qijg0SP06kW7FlmoWRMTJmDfvgoH7NiBunXld3VXSVVi9Gg4OmLkSNHhbLt3V8dqMvIJl4ujR3H/PnbupFIshwMzM0iVVlMWhXSEAPz8sHYtvn+nS76pKZo1w+XLdMkvQCjEp08YN45eLbIzZw727EG5nXo/fcKaNco6MkqYY8UKtGoFO7vKAgUePMDXrxgwgEGzlFSKlhbOn8eaNbh0iUqxVAWOKqojbN0aQ4di/XoaVTAQMnPnDjgcdO9OrxbZad0aLVuWX0h3yRLY26N1a8ZtUlJd4XCwdy/evKmsiNfu3Zg1i/qmwUpkoWlTnDoFOzvExVEms7o7QgDe3tizp7K9KxkZPx43b+LDB7rkR0VF2di0B3qkpdFcOI4KnJzKqTJz9y6uXqU9s1OJkjKoq+PsWQQG4mR5bb6+fcPZs7CzY9oqJSLp0QMbN2LECMqCSKkKHFVgR9iwIebMgYcHXfJr1cLo0TS2rV+zZnNu7nSBoNWhQ4fo0kEdo0cjNhZJScVnBAI4O8PPD7Vrs2eWkuqKgQFCQuDiUk4Q6cGDGDJEriOxqzPTpmHMGMqCSC0sEBtLQTadAjtCAEuW4No1REXRJb9gdZSmnEULCwfgnrr67eHDh9OigFJUVWFvXypAyd8fWloKsMGppKrSvj0CAzFqVNkg0j17qns1GTln/Xro6FDzP6pVK4/LPXbyZISMchTbEWpqwtMTixbRJf+XX6Cujtu3aRF+4IB1jx57vn593EKuWhFWzOzZOHSoMFrvyxesXo2dO+ntzqpESeWMGgUnp1JBpNevgxD06MGqWUoqpSCI9OFDCoLstm7dnZ4e6ei4MkG2VheK7QgBODggJQWhoXTJt7OjJWQmOBhv3+Lw4VoaGhrUS6cHIyN064bgYA6AxYsxeTLMzdm2SUm1Z9kytG6NadMKV278/eHiwrZNSkShqYnz57FunayP7gYN6qmpPefz07W0tGSRo/COkMfD+vVYtIiunL+pUxESQk3OZknmzsXAgWjWjGKxdOPkhIAAblQU98oVGndnlSgRHw4He/bg3TusXo1PnzjXrmHSJLZtUiIGTZogOBj29jIFkU6ZMiEiYou6epiamqEsxii8IwQwbBgMDekqBKOnBxub8oPTpGbfPqSmYq/ofuByR+/eWf/806tPnwkuLk+1tf/P3p0HxLi9cQD/zrSn1Z6sqRCVrUhUsmSNcm3RopD12rdruX725bpkjXCJki1Zs2TN2mJJykVRWRMpJS0z8/tj3JBqtneZmc7nr2vmfc954jbPvOc85xy2oyEIAP8Vka5fv7FFC5fq1deQ6i1F0akT1qzht2vX28jI9sKFC9I1YmPT3NW1lowVh8qQCAGsXYtFi/D5My2N9+iRNHu254YNQZS0xudj9my4ucHIiJL2GLV37z88XkMer3dUFDl7l5AjdeogP39DSUlwSspWtmMhJNCgweXi4oLPn1dOmLBM6kb8/bFtm0xVjUqSCFu3RrduWLeOlsaPHl2Vm+u9eHFQbm6u7K1t3IicHNHbWMsnV1dXTc3bqqorx43zYDsWgvhJnz5dudxuJiaMHpdIyMjOzk5X9x2Hs8jAYIrUjXTqBF1dXLkifRhKkggBrFqFzZuRlkb9+ZPDhvUxMJhXWFgfkGk+FgCPhz//hKcnatSgJDSmNWjQIDf339zcpMGDB7MdC0H85NixXZ8+PXz2LJrtQAgJaGtrf/qUnJl5Q1XVffly6dsZPVqm3aeVJxHWr48GDf5o0aLTuHEUnwbk4zPsw4c7I0aET58u61qBFStQUICAAEriIgiCUAZaWjh2DDt3Yv9+KVvw8kJUFN69k/J25UmEALKyLhQU7A0PP095y1wuNyAA0dE4fFj6RoqKsGoVxo2DbIW+BEEQyqZuXZw+jVmzcPOmNLfr6mLgQOlLJpUqEe7cucrWdllx8Wo6ztStVg0hIZgyRfrdR+fPB59P70bhBEEQCsrCAnv3wt0dT59Kc7u/P3bsAJ8vzb1KlQhdXLrdubN/+fJurq6goq6lrPbtMWaMiPNfKlJQgM2bMX06tLWpD4wgCEIJ9OyJ5cvh6orsbInvtbFB9eqIipKmX6VKhEITJsDeXsp0JdLChcjNxVbJK7SnTYOKCv78k/qQCIIglIafH/r1w8CBKCyU+N6xY6UsmVHCRAhg82a8eoU1a6hvWVUV+/dj8WLJjkXOycHu3Zg/H+rq1IdEEAShTFavRs2aGDdO4hs9PHDlStlN2MWhnIlQQwMREdi8GZGR1DfetCmWLYOXlwTHiEycCC0tzJpFfTAEQRBKhsvF/v14/BiSLqjQ0cGQIdKUzChnIgRgZISDBzFqFFJSqG/c3x/162PxYrEuzsrCwYNYvhyqqtRHQhAEoXxKF1SEhEh24/jx2LEDPAnXkyttIgTQqRMWLoS7O/LzqW981y7s3SvWXgajR0NfHxMmUB8DQRCEshIuqJg8+Ur79h4REWfEvMvKCvXq4exZyfpS5kQIYOJE2NjAy4v6wplatbBrF7y9RVQ3vXqFU6fw99/gKvnfNEEQBMUsLKChMTc+/s/RoxeIf5e/v8QlM8r/8bx5M16+xNq11Lfcqxd698a0aZVd4+OD2rXh6Ul97wRBEEpvxIh+hoaeBQXOc+eKO+A5bBhu3UJamgS9KH8i1NRERAQ2baKlcObvv3H7doWHNKWk4PJlbN5MjnEnCIKQxl9/LcjKuv3ixV/376NbN7x+LfoWLS0MH47duyXoRfkTIegsnNHWRkgIJk9GRkY573p6omFDuLtT3ClBEETVweVya9VCZCTc3NC+vViPNOPHIygIxcVidyFLfAqkUycsWEBL4Uy7dpg0Cb6+ZachHz1CTIxMG6ITBEEQQhwOpkxBWBjGjsWUKSKSXIsWMDHB6dPiNl5VEiGASZPoKpxZsABFRWXPlBg5Eubm6NGD4r4IgiCqLAcH3LuHZ8/QuTOeP6/sSolKZqpQIsR/hTOUn9/L5WLPHqxYgYSEb6/ExuLhQ+zaRXFHBEEQVVzNmjh1Ch4e6NABR46AV0EJzZAhuHsXYh7AULUSobBwJiCA+sKZJk2wciU8PPD1KwB4esLaGnZ2FPdCEARBCIdJjx/H5Mk3DAzsmjd3yMnJKXONhgZGjsQ//4hVqVi1EiEAIyOEhdFSOOPnBwsLLFqES5dUnz2TrGaJIAiCkIidHTw8rubn+z55Us/J6cmiRbh0CQUF3y8YOxZ79nDE2QuzyiVCAPb2mDOn2NZ2Zs+eXm/fvqWw5a1b8fffNgMHOtSpE2JtTWHDBEEQRFlz5oweNOjB9OlWAQHtqlXDypWoVQvt22PuXERFoV694oKCYTt33hfZThXd/rJdu1t5eZ8uXOg0b17IP//MoKrZ/Pw0Hu8TEPHmjSswgqpmCYIgiF/Vrl378OFtwv92cMCcOcjPx40buHIFixbh/v3rBQW579+bi2ynKj4RArCysmrW7EXNmntu3Ojeti0OHZJ4k9YySkoQFgZ3dyPgC+BWu7YGRZESBEEQ4qpWDT17YsUK3LyJhw9bqKn9a2DwUORdVTQRGhgYJCREvX9/+8kT640bERQEMzMEBPw0viymz58REABTU2zejLlz1Z8+vb5+/YSXLxNE30kQBEHQpmnTui9eRI8a1ULklVU0Ef6oc2dcuIDgYERFoUkTLF6MT5/EuvHNGyxejCZNEBWFsDBcv47Bg2Fq2sTPz0+VHLlEEATBtrp16+rp6Ym8jCTCbzp3xsmTuHABqalo2hRTpuDNmwovvnsXXl5o1QrZ2bh7FydPomNHBmMlCIIgqEMS4U8sLREcjPh4ALCwgJcXnj79/i6fj5Mn0aMHBg1Cy5ZISUFAABo2ZCtYgiAIggJkBK8cjRsjIAB//IENG2Bnh9ato/79d7GRUbvs7IAaNTBjBtzdoaLCdpQEQRAEFcgTYYXq1MHKlUhJQXr6npcvtzx4ELdlS+7t2xg8mGRBgiAI5UESoQj6+tixw7958ym+vo49e4qedCUIgiAUCxkaFc3JqUty8hW2oyAIgiBoQZ4ICYIgiCqNJEKCIAiiSiOJkCAIgqjSSCIkCIIgqjSSCAmCIIgqjSRCgiAIokojiZAgCIKo0kgiJAiCIKo0kggJgiCIKk2snWWePHly/fr1evXq9ezZk8sluZMgCIJQHqKzWkREhJ2dXWxs7Pz583/77TcGYiIIgiAIxoh+Ipw/f/6mTZs8PDzy8vJMTU1v3rzZqVMnBiIjCIIgCAaIeCJ8/vz5v//+O3DgQAA6OjouLi6nT59mJDCCIAiCYIKIJ8JXr14ZGhpqa2sL/2hsbPz69euKLi4qKgoODq5Xr96PLzo6OlpbW8seqAIpLi4uLi5mMQAej0dfDCUlJcXFxaqq5NySnxQXF/N4PBWlPqmyuLiYz+ezHUVlaP3VY/33WmpyHrlAIKAvQj6fz+FwRF4m4uOMx+P92AqXyy0pKank4levXhUVFf34orW1NY/HExmHMuHxeOz+yLz/KFzjiqsq/LUIfzp5Lpej9Z9Acf995TxygUBAX4R8Pl+cr6ciEmG9evWys7OLi4vV1NQAvHv3rswD34+0tLTmzZtX1Z7/flVcXKypqcliAMKvLzTFUFJSoqKiwu4PKIdUVFR4PJ7S/7WoqKgIPwrkE62/eqz/XktNziMXCAR8Pp+mCPl8vjgpVsSXOxMTk3r16kVFRQEoKSm5ePFi165dqQmQIAiCIOSAiCdCFRWVP/74Y8yYMVOmTImOjq5du3bPnj2ZiYwgCIIgGCB6uN/f33/v3r05OTm9e/e+ePGiPM8QEARBEISkxKr969atW7du3egOhSAIgiCYRx7vCIICfD6/TL00QRCKgspEmJpqYmMzgMIGCUIhbNmyRVOzib5+q1GjRrEdC0EQEqNyWTSPt7y4+LfXr19XssSCIJTAhw949AhJSUhMRFISrl07BPQCmuzf/zwgAHp6bMdHEIQkqN0fZCMw6/bteu7ulLZKECxJSEjw8ZljZWU1adKqpCROUtK3/PfhA0xNYWGBli3h4gI7u20rVvQCquvqnq5RA+7u2LMHWlpsR08QhHioTISmprGpqVuHDUNGBurUobBhgmDH5MlB9+7Nun//r4SE9HbtGllYoFs3WFjA2Pj7NRkZ8PGxqF8/rVYtwd9/c9PTMXMm9PUxahS2boVS77lGEEqCyjlCDqdk4sR3mpogp1MQii4+Hl27IiNjSJ06C1xcasTG1g8KwrRp6NHjpyxYUoJhw9CyJXx9+W5u/KNH4eWFzEysXYvQUOjq4s8/2fsZCIIQD8VVo+7uHzQ08OoV/PyobZggGPLyJfz90bcvBg7Ekydd3r69GRm5r6LtCv/3P+jo4OlTDB0qcHPjHTsG4a7UU6YgNxczZ2L1ahgYIDCQ0R+BIAiJUJwItbT448ejWzfs2YOICGrbJgh65eVh8WK0bg1DQzx5gilTUPkZG9HR2LULo0ahUSOYmQlMTQWGhrhz59u7HA6WLEFuLkaOxO+/o0YNBAcX/PXXXwkJCQz8LARBiI/6dYRTpuDOHQwZgmHDkJlJefMEQb3iYuzYAXNzpKbi4UOsWiW68vPTJ3h5YccOnDkDD49vLw4ahKNHf7pMXR2bNyMrC05O8PaeO2vWszZt+uTl5dHyYxAEIRXqE2GNGhg2DCYmMDaGkxPlzRMExaKi0KYNDh9GZCSCg2FkJNZd48bB3R3du+P0aQwZ8u1Fd3ccPQqBoOzFeno4ehSNGt0FdPl8/cxMuT7VjyCqGlp2lpkxA9u34/x5pKZiwgQ6eiAImbx69SowMCgy8rWTE6ZMwapVuHAB4h8gFhiIp0+xYgWOH4eNDerW/fa6lRU0NHDvXvl33b9/ctYslYYND9nY6OXmUvBTEARBCVoSYZMm6NkTERHYtw/bt+PCBTo6IQjpOTsPnzCBN2DAcG9vPHyIfv0kuDcpCYsWISQEGhoICcGIET+9O3AgwsPLv9HAwGDNmlVPn7bU1kaLFiA7shGEnKBrr9E5c7BhAwYMwODBGDAAOTk09UMQEisqQmamnrp6Ytu2eqNGQaLzVAoL4eGB1avRvDk+fsS1axjw866Cgwbh8OHKWlBXR3IyCgvRpo00wRMEQTm6EqG1NSwsEBKCAwdQsyYcHGjqhyAkw+PBywtduoRfvDgyOrqCZ7eKzZwJMzMItxQ9dAi9epUtq2nfHoWFSEqqrBEdHTx8iLQ08ntBEHKBxtMn5szB6tUQCHDrFpKTMX06fV0RhFgEAowfj6wsHDqkbm/fUU1NTaLbIyNx4gS2b//2x9DQsuOiADgcuLmVrR39lZER7t5FTAxcXSUKgSAI6tGYCJ2dYWiIkydhbIy9exEQgMuX6euNIESbPRsPHyIiApqaEt/77h3GjEFoKKpXB4D0dCQlwcWlnCt/XURRLnNzREcjMhK+vhIHQxAEheg9j3DGDKxcCQDDh8PVFf37gyygItjy5584fx5nzkBHR+J7+XyMHIlx42Bv/+2V0FAMHgx19XIu7tQJ79/j2TPRzdrY4MwZBAdj9myJQyIIgir0JkJ3d3z8iOvXASA8HHp66NqV1g4JonybNuHAAZw7B0NDaW5fswaFhZg37/sroaHf19GXweViwIAKa0fL6NEDBw5g3Tps3ChNYARByI7eRMjlfttuEQCHg5s3cf8+/viD1j4Joqy9e/HXX7hw4fuCP4nExWH9euzb9/0oiYQE5OR8fzr8lZijo0KDB+OvvzBtGkJDpQmPIAgZ0ZsIAXh74+5dJCYCQOPG2LEDq1fj5k26uyWIb8LDMW8ezp9Ho0YS3/v+/fvff/+fq2v4xo0/3S58HKxk3YWTE168QHq6uB1Nm4Z58+DpichIiYMkCEJGtCdCDQ1Mnoy1a7/9cdQouLigZ8/nS5euzs/Pp7t3ooq7cAHjx+P0aTRrJs3tc+eu2ry5Znb23126vC59USDAwYMVjosKqaigXz8cOyZBX8uWYcwYuLri9m1pQiUIQmq0J0IAEybgzBmkpX37Y0hIbn5+n0WL0m1s+jDQO1Fl3byJESMQHi790vX27dsAB+vUKTb8YWoxOhp6erC0FHGvu7u404SlAgPh6orOncOMjbucPn1G8ngJgpAGE4lQTw+jRmHDhm9/5HD4XG4eoP71KxO9E1XT/ftwd8e+fZXN5IlUt+7ILl0iUlJuamlplb4YEiLicVCoRw8kJuLNG8l6PHCgiMf78/XrhcOHz5IwWIIgpMRQKpo+Hfv2ISsLAAwMDKKjDxoaWjdrdpqZ3omq5skT9OuHoKDy1/mJ7+hRDBtW/cdTeYuKEB6OoUNF36uujl69cOKEZD2qq6vr6QH4o7BwCNmMlCCYwVAirFsXAwdi27Zvf+zUqdOGDT5RUdrFxcz0T1QV6enpISGXunXjrV2L/v1laqqoCGfPlt1KNDISFhZo3FisFiSqHS2Vk/Pvy5fn9fT+FP80DIIgZMHc4OSsWdi0CaX1MV5e0NDAihWM9U8ov48fP7ZvP8DLK8Ta+q/hw2Vt7cIFWFigXr2fXix3W7WK9OqFmJhvAyESMTaunpCA9HT06CGiigIxAAAgAElEQVTxvQRBSIq5RNisGeztsWfP91eGD8eWLYz1Tyg/Pp+fm8tRUdGwtCyRvbWjRzFo0E+vfP6Mc+fKvlgJbW10745Tp6Tp3cgI8fG4ehV+ftLcThCE+BgtV5k3D+vWoeS/z6h16/DhA86fZzIEQpnFxtY0Mgo7eLDv0qVzZGyqpASnTsHN7acXw8Ph5IQaNSRoR7rRUaHmzREZib17ycAJQdCL0URoa4sGDb6f1qanB1tbstEMQY38fEyciF27zN3c+qqqqsrY2uXLMDFBw4Y/vShmveiP+vVDdDSkPo++WzcEBWHhQhw4IGULBEGIxPQChjlzsGoVBIJvf1y/HnfvIjOT4SgIJTR/Ppyd4exMTWu/jotmZiIuTrKz7AHo6qJLF5yWoT561Cj88Qc8Pb/t2UsQBOWYToR9+oDL/T4c2rEj6tTBLLJiipBNXBzCwrBqFTWt8fk4cQLu7j+9eOAAXF2hrS1xa4MGSbyyvoylSzFkCJyd8fixTO0QBFEuFpa0z5iBNWu+/3HqVBw6xHwUhPIoKYG/P/7+GzVrUtPgtWswMkLTpj+9GBIiQb3ojwYORFQUZNxPMDQUbdrAxkaaGlSCICrHQiIcNgxJSXemTNmcnZ0NYPZsCAQIDGQ+EEJJrF+PWrUknr2rxK/jos+eIS1NykPEDAxgY0NBUditW6hTB61agSy0JwhqsZAICwvzP3/237SpaNy4PwBwOOjbl7JBLaKqSUvD2rVUrsMRCHD8eNlx0f37MXw4pC7BcXeXvna0FJeLhASUlEi/dSpBEOViIRFqaGgYGnKB640bNxC+snkz0tNx9y7zsRAKb9IkzJpVdhhTFrduQV8fzZv/9OKBA1KOiwq5ueH0aXz9KmNo0NbGgwd4/lzimh2CICrBQiJUVVV9+vRm06Yrhg79tnLCyAgtWmDmTOZjIRTbgQNIT8fUqVS2+eu4aGws+Hy0by99m3XqwNISly7JGBoAGBvjxg2cOwdHx8NBQUEUtEgQVR475z9oamr26NH82rXvr6xciWvX8OULK+EQCunjR8yYgcBAqKlR2WxERNlEGBKCkSPB4cjUrCwr68to0wYeHpuuXdvl779/3bp11DRKEFUYawchOTggOvr7H11doauLBQvYCodQPLNnY8gQ2NlR2WZcHFRUfjprkMfDwYOQfefSQYNw4sT3bZVk1Lp1CYeTJxCoHTlSm5oWCaIKYy0ROjoiOvr7ynoAvr745x+2wiEUzLVrOHsW//sfxc0ePYrffvv+x8+fP0+atM3AINrcXNaW69eHiQmuXpW1HaFp06b988+YOXNGx8d7tm4NPp+aZgmiamItERoZQV8fycnfX1m+HJ8/UzZ8RCixoiKMG4ctW6CvT3HLx479NC46Z86KHTuyMjJmvnv3TvbGKRwdBeDt7b1q1bBHj/D8ORo1wqdPlLVMEFUNm2fEOzjgx2lCTU04OmLRIvYCIhTEihWwsCh7UqDsEhJQWIi2bb+/YmLSgMOJqVatsFq1arK3X7362R07Og0fPkn2pkqZmSEjA1wuGjbEkycUNkwQVQibibBLl5+mCQFs3IjkZKSlsRQQoQiePMG2bQgIoL5l4bjoj0Ux3t4TqlVb9+zZdR0dHdnbDwsL5vECz56Nz8vLk721Unp6eP4cbduiVStylgtBSIPlJ8IyUyYtW6JhQ0yfzlJAhNwTCDB+PBYvhrEx9Y0fOVK2XjQmBh07muvqUpAFASxcOMHAYGq7di6UpNUfcbm4cgXDh6N3b2zdSm3bBKH8ZD2tRhYmJuBykZoKE5PvL86ZgylTwONBRYW9yAh5tXMn8vPh7099y0+eICcHHTr89GJsLGxtKevC0bHz4sWXnj6lrMEy9u5F8+aYPBnJydi0ia5eCEL5sPlECKBLl5+mCQGMHw81Naxdy1JAhBx79w7z52P7dnBp+N/20CEMGlR2sWBMDGxsqOzF2hoPHlDZYBnz5uHwYQQGwsmJxl4IQsmwnwjLTBMC+O03WmaACIX26NGj0aOfjx0La2ta2v91QxmBALGxFCfC1q2RkPDTqiHKubsjJgaxsWjWjOxQQRBiYTkRlikcFfr7b2RmlpMgiSrrzJlIO7upkZHDBwxIoKP958/x5g3s7X96MTUV2towMqKyIwMD6OvjxQsq2/xVmzZITUVODurXzxo6dNyFCxfo7Y8gFBzLibBFC+Tm4uXLn16sUQNt2mD2bJZiIuRPdvbn/PzaWlr6JSWyHetXgcOH4eZWdlo6JobKCcJSdI+OCtWpgxcvkJvrdeiQiYvL6BKqtrQhCGXEciLkcNC5M65fL/v6338jJgYfP7IREyF/srIGt2497MSJeXbU7qj2n1/HRQHqx0WFmEmEADQ10bDhRyBRIKhubq7y628ZQRBCLCdClFcvA8DBATVqYM4cNgIi5MynT1i5krNnT/+uXZ3oaP/lS6SmllNdQm3JaCnGEiGAJ0+uBwU5PHp0rnlzjqMjLCyQmMhQ1wShQNhPhOVOEwKYOBGhoYxHQ8iflSvh6vrTRtjUOnIErq5lD93l8fDgwU+7zFCFyUSoqqo6evRoC4vaZ84gORn6+rCygp0d3rxhKACCUAjsJ0Jra7x+jczMsq8vWIDiYrINd1X34gV276Z+c+0flTsu+vAhGjaEnh713Zma4v175OZS33LlzM1x6xauX0dmJho0gKcnBQcFE4RyYD8RqqjAzg43bpTzuosLli1jIyZCbvzxByZPprh080dv3yIpCd26lX2dpkoZAFwuLCzw8CEtjYvUqRNSUrBrF06fhr4+pk8nJ1cQhBwkQlSwmhDApk14/hyPHjEeECEfYmJw7RpmzKCxi/Bw9OkDDY2yr9NUKSPE5Ohouby98fEjli5FYCD09eHist7Kqlt8fDybMREEe+QiEVY0Tdi4MczMyNajVdfMmViyBFQc/FChcsdFQecTIeQgEQrNno3cXAwYkHX+/N6HD33s7BZOn44LFyg7PZggFIVcJEIbm287Pf5q1qyPFy4sGTFiPONBESw7dgw5OfD2prGLrCzcvQsXl7Kvf/mC1FQay3PkJBECUFXFnj0G6urZHE6Arq51WBj69oW6OnR1YWGBYcMQGIhXr75dvGbNmr1797IaLyGxT58+jRs37lXpvyJRHjY33S6lpgYbG9y8id69y751/Ph4gaBZaGh6t267fX192YiOYEFxMebOxaZN9O69HhGBnj2hpVX29bg4WFpCXZ2ufq2s8OgR+HxaNk2VlKqqakHB87S0tCZNmghfefUKx4/j2jUkJODUqW/b/3I4w4qKngJ5KSkpK1euZDdmQnzVq7cUCHqGhtoKBCQXVkgOfhEBVDxNaG/fnsM5BmRY07TFJCGXAgPRtCl69qS3F1bGRQHo6qJ2bTx7RmMXEuFyuaVZEICxMSZMQFgYkpKQl4fiYpw6BT29dMAUMAkPb1Du4A0hhy5dgkDAAWyAekeOsB2NHJOXRFjRNOHcubOSk4+rqcXExbVjPCiCHZ8+YflyrF5Ney/lDkKA5koZIfkZHRVJVRU9eyIj45K19RtVVceuXSeYmyMggJSbyrU3b+DlBW9vWFlt19I63KrVX0uXwtmZ1B6WT14SoZ0dHjwof7P8Zs1MBg3SXLKE8ZgIltC9gl7o+HF06wZd3XLeIonwV5qamrdvnx84cHabNoiMxMGD6NgRMTFsh0X8oqQEAQGwtIShIeLi8OpV30uXTr186XjxIgYMgJMTpkzB589sRyln5CURamnBygp37pT/7tatePsWkZHMxkSw4eVL/PMPvSvoAQQGBk+a5GhouP/Xt96/R3Y2zMzoDUDhEqHQqFH8XbvQti1u3MCkSRgwAF5eyMpiOyziP9HRaNsWx48jOhoBAYiIgKMjWrbkd++OY8cwZQoePwYACwsEB9N7HJhikZdEiIpHRwEYGsLRkd71ZIScmDMHkybRuIJeaOnSgLy8U2fPrv/1LeEEYZkTeimnoImwWzdBZiYePACHAy8vJCXB0BAWFmSklH0fPsDfH8OHY+ZMXLqEFi0AYMcOjB8PAH5+2LULAGrUQEAAwsOxZQu6diV7z34jR4mwonoZoR078PixQn52EOK7d4/2FfRCPXr4amp2//13v1/fommv7TIaN0ZuLj58oL0janG58PH5vvGhoSECAnD+PA4fhq0tbt9mNbiqis9HcDBatYKmJpKT4eX17fWbN5Gb+23XpJ498eYNEv47zdPGBrduwdcXPXpgyhQWNvyTN3KUCDt3RkwMiorKf9fUFK1aYdIkZmMimDVzJv73P3pX0AtVqzZx6dI7c+ZM+PWtmBjaJwgBcDiwtGRtozVZ+PkhNBSFhd9fad0a0dH4/Xe4ucHLC+npX+Li4vjkCZFmJSUl3t7TrK0Htm2buns3oqIQEPDTnPe2bZgw4dvYBpcLb2/s2fP9XS4XXl7fameEI6XPnqVkZGQw+jPIDTlKhHp6MDNDJds8bd6MGzfw+jWDMREMOnYMWVn0rqAvdfYs+vQp53WBgIlKGSEFHR1t0ACtW+P48Z9eFI6UJifDwAAmJo5du27x8ZnGUoBVRUJCQnj4+4SE4fXrB1+5gpYtf3o3KwunT//02+Tri/37f/oGA6B6dQQE4MgRLF9+2cLC19r6t/v37zMRvZyRo0SISqcJhe/Wr08eCpVTSQnmz8fatfSuoBd6/BjFxbCwKOet1FRoa9M+QymkoIkQP0w4lWFggA0b+AYGRfn5lomJpISGXg0atCgqyjM23rR4cf9f3929G25uqF79+yuNG8PKCidPltNUx45YsuQDh9MgJ6fO69efaAtZfslXIqx8mhDAqlU4ebL8VRaEQtu2DY0a0b6CXigysvzHQdC/lP5HipsI3dxw/z7S0sp5i8vlXr9+cPZsvZcvNz9/znhkVUlAgJaHR8TLl9fbty+7xlogwM6d38pkflTRNxgAQ4YMCg0d1LPn+L17naiPVe7JVyJ0cMCNG+DxKrzAwwO6upg9m8GYCPp9/oyVK7FmDUPdRUaWv44ejKwgLGVp+e3ZVOGoq2PYsJ8mnH7UvHnzVatGz5tnOHRo2YE4gioZGdi+HRWtrj57Fjo6aN++7OtuboiLQ3p6ObdwOJxBg9wiInqnpGDzZoqjlX/ylQhr1kS9et9Lm8o1fTp27SK12kpl2TL060f7Cnqh/HzcuQNn5/LfZaZkVEhLCw0a4MkThrqj1ujR2LWrsu+sU6eiYUOy5Ikus2bh99/RoEH5727bhsmTy3ldUxPDhqGSjdM1NHD0KJYtw82b1MSpKOQrEULUNCGA+fMBYO1aZsIh6PXly5cdO47u3JlG9wr6Upcuwcam/A1leDw8eIC2bRmKBIo8OmppiTp1cOlShRdwOPjnH1y8iP3lbFpAyOTWLdy8WeGXjPR03LyJoUPLf3fMGOzeXdmDRKNGCAqCh0fV2idB7hKhyGlCYX0aSYTKwctr6oQJ94uKXGvXrvjhglKVjIs+fIiGDaGnx0wggCInQlQ64SSkq4tDhzB9OpKSmIqpChAIMHUqVq+Gtnb5F2zfDi+vCt+1skL16rh8ubIu+vfHb7/B27sKDbzJXSJ0dMS1ayL2/tmwATk55JumMsjL0xAIcvT0OBy6t3L5z7lzFSZCJitlhBQ6EQ4fjvPnRTw3WFpixQoMGYL8fKbCUnbBweByMWxY+e8WFWH3bvj7V9aCyG8wAFatwufPtG98Lz/kLhEaG0NX99uGeBXR0kLv3li4kKmYCHoIBPj8+e/Zswc8eBDFZeRovuTkChdOgNlKGSGFToT6+ujXD6GhIi4bPRo2NhgzhpGYlN2XL/jzTwQEVLgF4NGjaNUKzZpV1siIEYiMFPENRlUVhw5hyxZERUkfrQKRu0QIMaYJAezYgfR0EYOohJzbvx+FhWrLl3erWbMmMz1GRqJv3wrfZf6JsH59lJTg3TtGO6WQnx927hR92datSE4W60qicqtWwcGhsv9Lt20rZ9VEGfr66NsXYWEiLqtbF8HB8PauEnuYyGMiFDlNCKBuXdjYYMoURgIiaFBQgAULsGEDo6e0VzJB+OULUlMZqlz9kaWliDJpeebggK9fERcn4jItLRw6hPnzcfcuI2EpqZcvsW0bli2r8IKkJKSkoH85y+vL8vNDUJDoy5ydMXYsPDxQUiJBnIpIHhOhgwOuXhV92fbtuH8fKSn0B0TQYOVKdOmCzp2Z6zE/HzEx6Nq1/Hfj4mBpCXV15uIRUujRUQ4Ho0aJnnACYGaGTZswdCjI6fZSmz0bkyahYcMKL9iyBWPHQk1NdFNOTvjyRazvJQsXQlMTixZJEKciksdEaGoKPh8it6WwtoaZmehxAEIOvXyJrVuxfDmjnV68CFvb8hdOgI1xUSGFToQARo3C4cNibfY0ZAhcXODpSY7Bk8bt27h+HTNnVnhBXh7CwjB6tFitcTjw8RHrGwyXi9BQHDhQdndZJSOPiRDiTRMC2LABFy8iO5v+gAhKzZqFSZPQqBGjnVYyLgo2KmWEFD0R1q0LOzscPSrWxevX48MHbNhAc0xKRyDAzJlYubKyg1n274ezM4yNxW3TxwdhYWJ9g6leHQcOwN9f9MOJ4pLTRCjONCGA3r1RqxaZKVQwt27h+nXMmsV0v+fPi0iErDwRWljg2TPF3opMnHJ8ITU1hIVhzRpcv05zTMolJAQlJfDwqOya0jN4xWRsjI4dceyYWBd37IjZs6HEe+bJaSIU84kQwOLFCAtTyA0bq6bS5cAMHDr4I+HCCeGx3b96/x7Z2TA1ZTQkIQ0NNG2K5GQWuqZKv354+hRPn4p1cYMG2LULHh54/57msJSFsKzsr78qXDIB4MYNFBRUOP9dEfG/wQCYNg0NG1Y2NqvQxEqEr169unr1amJiooCp0f2WLfHpk1hlu+PGVYm5XKUhXA48fDjT/Z45g379KnxXOEHI1Jr+sqysFHt0VFUVI0Z8P7ZepD59MGIEhg2rbKtSotTq1ejcWURZ2bZtGDdO4v+BXV3x+DGePRPrYuGeeRcuKOdOJqIToa+vb5s2bRYvXuzq6mpnZ/fpExOnVXE4sLcXd5nguHHYsoXmgAgqfPmCRYtEfLelicgJQlbGRYUUfZoQwOjR2LNHgiL7ZctQVJRubT1h/frtdMal8F69wpYtIsrKsrJw5gw8PSVuXFUVHh4VniLyK+GeeRMnRjg7j753T6nO7xWdCMeMGfPq1avLly8/efKEw+EEBAQwEBbEniYEsGIFCguxbRvNAREyW7ECjo6wt2e63/x8xMZWNnAUE8NOpYyQEiRCc3M0aYKzZ8W9XkUFdeuue/Soy+LFBzMyMugMTbHNnYuJE0WUle3cCXf3n87gFZ/wG4z4j+aWlgI+f+Hly76+vvOk6U9eiU6EdnZ2ampqAFRVVS0sLD5+/Eh/VIAk04Sqqhg0qMKjuQg5kZGBwMDKlgPTR7hwQken/HcFAtZKRoWUIBFCwgknAEOGOBsa/v3lC5fHq0NbUIotPh5Xr4ooK+PzsWMHxo2TsovmzdGgAc6dE/d6Dodjb2+tqTmNz+8lZZdySVX8SzMyMo4fPx4REVHRBSUlJRcvXkxNTf3xRUtLSxMTEykis7ZGRgb3/Xt+jRqiL960CbVrc0+f5lcy/MUYPp/PZ3Xbdv5/5KrxGTM4v/+O+vUFzP/dnDnD6dULfH75M9wpKdDW5tapI9NfmCx/57VrQ1WVm5HBF7/2nRV8Pp/D4VT0Mw4ejJkzua9e8Y2MxGpt0KD+rq69/vpL3cODc/kyX5xl4OJESN+vHsO/10VFxVOnaixbJtDSquxX5swZTq1anLZtKwut8sh9fTk7d6JXL3HrP86cCc7NLXZy0ti9W+DjQ0HViEAgoPXzSpzLOML6l0+fPnmWN8Y8d+5ce3t7ADk5Od27d3dxcVlW8Vd6U1PTevXq6f18jI2Pj0+vXlJ+d3B31/LzK+7bV6yZh379tN684cTHi7EuhmZ5eXk6FT19MILH4xUWFmpXdBCLbEpKSoqLi7W0tCS6KyZGxdtbMz4+n56gRGjZstqxYwXm5uX/Shw+rHrihOq+fV9l6aK4uJjH42lqakp3+4ABWhMmFLm4yHX1yNevX1VUVNQqTlkTJ2qYmwumTCkSv02BAMOHa5qb85cskeCuitD6q8fk7/Xhw8enTl1fXGyQnn5IU7PC7Y6+fPkyZEjB8OF1R4yorG6+8sjz8znNm2vHxX2pU0eCrPbsGbdnT61jxwqsrWVNYAKBoKCggKbPKz6fr6GhoaGhUfll354ItbW1J0yY8OvbTZs2BZCXl9e3b98OHTosXbq0kra0tbU3bdpkbW0tbcxlOTkhNlalohMmy9i1C2ZmePJEh8mDVcslEAhYT4Rqamrykwj5fMydi7VrUbs2C38tSUngctG2bYV/Gw8fws4OMv6TyZgI27bFkydagwbJEgLtVFVVK0+E48bBxwd//KEuUTHUvn1o1w4ODuoDB8oaIa2/ekz+XkdF3c7NnaWvv7Wo6GvNmuXP/mVmZlpa9szKMvTxmaCjM7iS1iqPXEcHbm6IiKhW0Um/5WrdGps2YdQo7bg4WY/wFAgEXC63Gj0Lqvh8Pk+MKdBviVBdXb13BaOKX758cXV1bdas2caNGxk7NE7IwQHTp4t7cdOmqF17t63tnrlzey1b9gedcRGS2bsXqqoVHplNt8rrRQHExqLSL3hMsLbGqVMsxyA7OzuoquLmTcnqoQwNERYGV1dYW6NJE9qCUygDB846dWrVkiWe9erVq+iajx8/5uXVADq+eCHrji9+fvDzw/TpkpVzDx2KK1cwZgwOHpSxf/aJLpbx8fFJTEw0NTVdu3bt6tWrj4m5FQEVbG3x77/4/Fnc69+//x+PF7ZqVSCdQRGSycvDwoXYsIG1VXqVJ0IeDw8egPVRBOWolwHg6ytZyYyQrS3mzVPmjUskdeJEw6VLt/7+e2U7hzZv3tzSctqIEfVmz54kY3f29uBycfu2xDdu2ICUFGzeLGP/7BOdCLt06eLr65uTk5OdnZ2dnZ2Xl8dAWELq6mjbFjdvinu9pWULoIeamswjLAR1VqxA9+6srdITuXDi4UM0bCjr2I7sWrRAerpYGz/KOW9vHD+O3FyJb5wyBY0aKe3GJRL59AlnzmDkSBGXCQR49qzfqlUTKZkEEfMUkTI0NHD0KJYtw61bsofAJtFVo5MnT2Ygjoo4OCA6Gi4uYl18//7ZFy94pqYq27fD35/myAgxPH/+7bQstkRFoUOHChdOgL1DJ8pQVYW5OR49YnMVByVq1oSFxWU3t4u7do1u3LixRPfu3g1bW+zfLzoHKLfgYPTpI3pd4MOHqFkTFQ+dSmb48CJz881NmmjMmzeeK8kZoY0aYedODB+O+HiIU+Evn+R0r9FSXbqIu5pQqHFjFT8/TJuGIgpq0AhZzZ6N6dPRoAFrAYgzQSgnuUc5Rkf5fH5Cwu+XLrXx9BR7ev8/wo1LZsxQ7J1XZbd7N8aMEX3Z1atwdKSs00uXDhYWpi9bFntO/EWF/+nXD7/9Bm9vBT5gS94TYadOuHfvzbt3Egy1BAZCQwM+PrTFRIiBx+MdPPjszh1Mm8ZmGGfPikiEcvJECGVJhFwu18TEiMPZ1r59eylut7TE8uUYMkQZRomlc+sWCgrQpYvoK69dozIRtmxpoaNzg8NJMJVq7/lVq5CbizVrKIuHYfKeCC9dOvP16xALC8f09HQxb+FwsG8fDh7Eo0e0hkZUpnXr7iNHzmvSZBorCweFHj0Cl4vmzSu84MsXpKbC0pLBmCqmHIkQwL175zt2DOvTR8rK7dGj0a6dWI9ESikoCGPGiK4sEwgQHS1WvhRTu3btrly5XKvWDTMzMyluV1VFaOi3A2IVkbwnwmfPXggENgUFRpmZmeLf1a8fWreGuzt9cREivHjxgcdzLSpKYTEGkeOicXGwtIR6heuVGWVtjYQEBR5c+lHXrjVlOXFw2zY8eiTBcRZKIycHERHw8hJ95ePH0NGheNKhdWudr1+1pN75tX59hITAywtv3lAZFTPkPRGOH+/n7d3MwmKMpCMtZ84gJYXsxM2OR4+gobF36tS3R46w+Q8gMhHKz7gogBo1oKODtDS246CCvb1MR+9qaeHQIcydi3v3qItJEYSGokcP1K4t+sqrV+HgQHHvwjN/ZPmHc3bG6NEYMULxDtiS90SooaGxdKn/ixdukn5TrlMH48Zhxgx8lWnnLEJihYXw8MD69W3+/nuWMXu7Z+bnIy4OTk6VXSM/lTJCSjM6am+PmBiZCtbMzbFxI4YORU4OdWHJvZ07xR0TpnaCsJS9PW7ckKmFP/+Emhr+9z+KAmKKvCdCAPXqoVo1pEg+xrZpE7S0xBpnICi0YAGaNpXmdDRqXbiAjh0rWzgBto8h/JXSJEJ9fZiYyLpsZuhQdO9ehSYL4+KQnQ1nZ7EupikRdu4s0xMhAC4X+/Zhzx4JzuSSBwqQCAF07CjNrgfCqpmjR/HwIQ0xEeW5fh1hYQgKYjsOMcZF379HdjakKpGji9IkQlDxkQpg/XqkpmLzZpSIf+avwhKWyYizhO/ZM3C5kHCVpljatkVqKmQ8fL12bRw8CB8fpKQozL+aMidCAH36oH17uLlRHRBRnpwceHpi+3a5WFd77pxYE4RsbfxWLiVLhDIOsgHQ0EBw8Jdp07rWrWt3+rTE69sUSH4+jhyBt7dYF1+9KmLMX2pqarCxoWCbGDs7uLicbtHCzsqq21dFmJ1S8kQI4ORJpKUpw2548u/339G3L/r0YTuO/xZONGtW2TXyNi4KwMwMb99Ksz+ZHBLuCSV7EayKSoa2drUPHybt23eFgrDk1YEDcHQUd5sYOiplSlHyDQaAmtrl4uKp//6rmpb2moLmaKYYibBtWzx+jPx8ae6tXRuTJmHWLFI1Q6+ICNy4gVWr2I4DAHDmjIh8/Pbt2+oK1V0AACAASURBVPDwvaamb5mKSCwqKrCwQGIi23FQwdgY2tp4+lTWdpo1azZrVvdOnW5fujTpzh0qIpNLwnFRMdE0QSgkY+FoqT//nDJw4DVz837TppkUFFDQIK0UIxFqaKBVK8THS3n7+vWoVg0eHpTGRPwgMxMTJyIkRERxCmNEThA6Ow9NTMxdsWIYUxGJS8lGRyn5SF2wYOqNG9v27DEeMECCLfgVyMOHePMGPXuKdXFGBoqKINWqd7F06oS7dynYorJBgwbHjm2/f39yrVro3RsMHtYgDcVIhADs7KQfHQUQEoKIiCq3LIkZAgF8fTF2LDp0YDsUAEBeHuLjRXxlVlOrxuE8NzCg5SxQWShTIqTq2UKoTx+EhcHdHVFRlLUpJ3bsgJ8fVFTEuvjyZRofBwHo6MDMTPqnjjJUVLB7N5o2Re/ecj3mrzCJsEMHmRKhiws6dICcnwCuoLZtw9u3+ENuzkKOihK9cGLq1HAnp4HXroUzFZS4lCkRUvVEWMrJCYcPw8MDJ09S2Sy7CgoQFoZRo8S9ntYJQqHOnREdTVlrKirYuRNt28LZGR8+UNYstRQmEcr4RAjg+HFkZGDtWooCIgAAKSn43/8QEgI1NbZD+Y/IcVEAd+9q9uvnoKGhwUhEErC2RmIi+Hy246BCq1b48IHiDbe6dMGZMxg7FgweEE6vw4fRoQMaNhT3emoPnSiX7Mvqy+BwsGEDHB3Rowfev6eyZaooTCJs1AgcjkwbUNWujalTsWCBlEU3xK9KSjBiBBYvFlGfyaTExMRjx4506ybipPPbt9GxIzMRSUZPDzVrSrN9hBzicGBnR/2sXvv2OHMGEyZg/36KW2aFRGUyb94gNxctWtAZEODggOvXKf42xuFg3ToMHQpHR7yWvzJShUmEkHl0FMDatTAwwPDhFAVU5S1bBn19jBvHdhz/ef/+vYODd1bWjeDg1ZVcVlCA5GS0bs1YXJJRptFRyp8thNq0waVLmDcPe/ZQ3ziTHj9GSooEK44uX4aDA+2LX+vWhYEBHj+mvuU5c+DtDWdnvHxJfeOyUKREKMtqwlKHDuH0aShxHTZj4uOxYwf27JGjNemqqqrFxTx19fc6OpUd/hQfj1atoKXFWFySUaZE2KULxdOEpVq0QFQUFi1S7CXCO3bA11eCaQVaF078iPL53VJz5mDyZHTpIl/DHlUuETo6olMnDBlCRUBV2JcvGDECGzbAyIjtUH5gaGhoZha5du2ERYtmVHKZ3I6LCp0547t2bfutW/eyHQgFbGyQnExX6XyzZrh2DRs2YMMGWtqnW1ERQkIkKJMBIxOEQlQtqy/XxImYNw+OjkhKoqsLSSlSIrSxQWIiBeviT5zAmzdYXdngGSHCzJno2FHuvk+8e4e0NKNx4zpxKn1KvXNHXlZ6/IrH4z17Fl9YGLhvn9xVtEpBQwOtW9M4ANO4MS5dwtatWLKEri7oc/QoWrdG06biXp+ZicxMtGpFZ0z/oe+JUGjsWKxejZ495WUjaEVKhFpaMDenYC2goSHmzMH8+QeGDvXPJ5UzEvr69evRo59OncL69WyH8ouICPTuLXqgSZ6fCFVUVJYunaqismb58oVsx0INuj9SGzZEdDSOHMHcuXj16hWNPVFNojIZAFevonNnsXblll3z5vj8md6ZvBEjsG4dunfHrVslb9+yvMeTIiVCULGIQqh58xAeb++hQ9peXv4UNFdlJCYmGhhYeHr2799/t6Eh29H84vhxDBgg4prXr1FUhCZNGAlIKpMmjbK1PaSmJtlJ1HKLpnqZH9Wpg6goBAT0srAYamBAc0klRVJS8OgR+veX4BbGJggBcDjo1In2fXyGDsWaNVn29iZmZv0mTJhJb2eVUrBEKHvhqFDjxo253MfAnUePulDQXJVx69YtHq8VMCA1NZLtWMrKy8ONG3BxEXHZzZuws2MkIBlYWSEhge0gKGJvjzt3QPcxSrVrA3gqEKzLzRWxckZO7NgBHx9ItJCVsQlCIWo3BqpI06aPOZzaAsGUs2eTae+sYgqWCKl6IrS3t//336hdu/5OTfV3d6egwSqiXj0/be3qVlbRu3dvZDuWss6ehb099PREXCbPE4SlLC2VJxEaGqJhQ1kP6RXHkSMbjY3X6Ohsz8ykvS8ZlZRg/374+Ulwy8ePyMhgdM0P3WPa//XS2dvbrkGDs7m5u6ja100KCpYImzZFQQEomQgwNTX19e148SJOnmT/OHWFkJ+PyZO5R4/uvX37uJFcVYsCEG9cFPI9QVhKmZ4IwdRHat++fZOTj44d22NGZSXDciEiAubmMDeX4JZr19Cpk7j7kVKifXs8fYqcHNo72r17U3JySFBQ3f79WTt6RcESIYeDDh2oLELr0gVRUQgLg48PZW0qq7lz4eyM7t1lPmKOBjwezp5F374iLisuxv37aNeOkZhkYGWlPButgZFpwlJLluDWLbnejHTbtr1jxw52cpLsU4yBLUbLUFNDu3bUjMCJw80NmzahZ0921lQoWCIERasJf+ToiAsXEBICX18qm1Uyd+4gPFx+d2q9dg0mJqhfX8RlDx6gaVPRw6es09dHjRpITWU7DopQu4lz5bS1ERSEyZPx+TNDPUpEIBDMn78uO3vB8eP/k+hGJitlStG6mvBXgwZhwwb06EHLpjaVI4kQAJyccO4cgoMlq2auOoqK4OeHjRshh5WiQso0Lipkba08o6ONG0NdHc+eMdRd165wdsb8+Qx1JxEOh2Nk1FFb28/HR4KjcHJy8PQp2ralL67yMVMv86MhQ7B8OVxcmP4WqHiJ0NYW9+6huJjiZp2dERGBf/7BxIkUt6wEli1DkyZyfYjViRNiJUKFqJQRUrJpQoY/Utevx7FjTH+Ii4nH23HlStzvv0tQKnP9Ojp2hLo6fUGVz94esbEUHNIrER8fLFoEZ2e8eMFcp4qXCHV10aQJLZ8R/frh2DFs347p06lvXHElJmLbNgQGsh1HxR48AJeLli1FX6lAT4SWlvKy6QYlGB5k09fH+vUYPZqCjaiodfMmBALY2Eh2F/MThEJ6ejA1ZeE8cz8/zJiBHj2oqYsUh+IlQgB2drh1i5aW+/fHkSPYuBHyX3jGDD4f/v5YuRLGxmyHUjExx0U/fMCHD3J0YlTllOyJkJnC0R/99hssLLBqFaOdihQcLNnmokIMryD8EfP/cEKTJ2PiRHTtSvF5lhVRyERIbeFoGQMH4tAhbNiA2bPp6kKBbNwINTXJFjwxT8xEeOsWbG0Z2qFKdmZm3w6fUw6Wlnj7Fgyv8NuyBdu2sVaR/6uiIoSHY8QIye7Ky0NSksQPkVRhsuK3jKlT4e+Pzp2fz5nz12Oa62cU5FPhZx070vVEKOTujr17sW4d5s2jsRf5l5aGZcsQGChHBy396tUrpKejUyfRV965ozDjogBUVNCiBR49YjsOinC5tBzSWzkjIyxdCj8/8HiM9luR48fRujUaNJDsrhs30L49NDXpiUmULl0QHQ0BS2umZszA58+j1qyp3727F60dKWQibNECHz/i3Tsauxg5Ev/8gzVrsHQpjb3IOX9/zJ6N5s3ZjqNSERHo2xeqqqKvvH1bYSplhJRsdJT5EkQAY8ZAV1deziwMDoaX5J/nrCycKGVsDB0d/PsvawHY2DTU1T2ZmWm0fDn1NZKlFDIRcjiwtUVMDL29eHlh1y4sXowVK+jtSD4FB+PtW0ybxnYcoog5LsrnIzYWtrb0B0QdKytlq5dhPhFyONi2DUuXsr8oMzMT169j4ECJb2RxglCIrWlCoVOn9sbGLkpNPRYfj3bt6PrYV8hECHpWE/7KxwdBQVi48JGBQQ8/v6m09yc3srIwdy527ZLg7GxW5OQgJgY9e4q+MikJdeuiRg36Y6KOkj0R2toiMRFfvjDdr5kZZs3CmDGsje8JhYTAzQ06OpLdVVCABw9YHslgcZoQAIfDadasWf363PBwLFwIV1f4+1N/1DNJhCL4+kJbe2xOzszdu69s3ZrJ7u8SMwoLCydNgre3AmxFFhkJBwdUqyb6SgVaOFFKuKZeaf6X09KClRWNZW6VmDEDOTnYu7eYz962ddKNi966BWtraGvTEJDY2H0i/NHgwUhOBgArK5w/T2XLCpwIY2NpP9tFqG/fVhzOZC5Xc+rUWpqacHJCpNydQUSZgweP1a3b5dgxx9mz5WwFVnnEHBeFQi2lL1W9OnR0kJ7OdhzUYesjVVUVM2bE+/nZNWrU8fXr18wHkJiI7Gxp1gKyO0Eo1LIlsrMZWsYgkqEhtm/H1q0YOxZDhiAri5pmFTURGhjA2Jih7VnDwrZ/+HCHx7v99Stn5058/Yr+/aGpie7d2RwxoMm5c9GfPk3S1FT5/FneD7MpLsb586I32hZSxCdCKN3oKMPL6n/04UMs4Pr+vWUSG5s679kDT09plu6wPkEIgMOBnZ18fdb16oWkJJiYoGXLpPr1u/foMfyrbFsnKGoiBP2LKH5kaGgIgMuFpydu30ZhIYKC8OHDt3G53r3B4kla1PL0nFGtWtySJR4NGzZkOxYRrlxB8+aoW1f0lZ8/Iy0Nlpb0x0Q15UuEt2+zs5jB23vEb7990dAwdXTsynDXJSUIDcXIkRLfWFiI+Hi5OEeaxW8wFdHWxqpVcHEJf/VqzKVL2q6uCevX49YtFEp1MLNiJ0JW5hsAqKjA0xP37uHzZyxahBcvYGuLGjUwciQ6d3atU6f13r172YlMZufOGU+evHHKlNFsByKaROOibduKtcRC3ijTCb0AqleHsTE7P5Guru7Bg6saNZp3+zaDZ/oBAM6fR+PG0mxpFBMDCwvo6tIQk4TkZ5qwjAULhlhZ7eratcTHx/rFC8yahRo10LIl/P0RHIxHjzBx4qxLl0Q/pih2ImTsibAi2tqYMwfJyXj7Fp6eOH/+5Y0bCQUFa3x9JTtjRU4IBDh4EB4ebMchBoFA3I22obDjolCuMyiE2P1I9fBAaCjTnQYHw9tbmhvZ2mL0V+3b4/FjeTzZytzc/MGD81FRez08NAICcP063r/H9u0wM8OJE3ByehwYmHD9uqnIdhQ4EbZqhTdv8PEj23EAAGrVwoYNePmyNlAIrAF8PnxgOybJXb8ObW3FGEK8exeamuJ+y1bEShmh5s2Rns7CkgP6sLKsvpSHBw4fZvQ4hZwcnDuH336T5l55qJQR0tBAmzasjcBJREsLnTtj5kwcOYK4OF0O50lBwVORdylwIlRRoXF9pXTU1dULC9NCQ2eami4yM6N37xs6hIRIU+HNiuPH4eYm1pUCgQInQlVVmJuzc2Y3TTp3xrVrrPXesCEsLHDuHHM9HjyIHj2kWcBaUoKYGNjb0xCTVOR2dLQShw8bd+6cNHeuucgrFTgRgsHVhOJTV1fv16/Xo0eoWxdmZsjIYDsgsRUX4+hRDB3KdhziEX+CMCUF2tqoV4/mgGijZPUyJiZQVWVzn5cRIxASwlx30i0fBBAXh6ZNYWBAdUDSYvdRXgrp6Vi1Ctu3a+nr64m8mCRCWqiqIiHh29fPtDS2oxFPZCQsLNC4MdtxiCEtDW/fivuQp7gThEJKdjAhgE6d2PxIHTIEZ88iJ4eJvp4/x7NncHGR+MbU1FRfXz9t7S00BCUle3vExNC42yflJk/GtGnizp4ofCK8cwfsbRZRGVVVPHyI5s3RsiWeP2c7GjGEhipGmQyAiAi4ukJFvOo/xR0XFbKywoMHbAdBKXa37DI0hIMDjh9noq89e+DhIc0+hQsW/J2c3Dch4cgrxo6mFcXAAI0b4/59tuMQz5Ej+PdfzJwp7vWKnQhr1UKNGmzujF45DgcxMWjZEhYW33YGklufP+PsWSmn9Jkn/rgoFP+JUPkSIeuzTSNGMFE7KhBg/34px0X79OnJ4aw0NdWqXbs21XFJj/V/ODHl5mLaNGzbBg0NcW9R7EQI+VhEUQkOB7dvo317tG0r12fLHTsGR0fF2JP60yfcvQtnZ7EuLihAcjJat6Y5JjrVrQs1NcjNgwEFrK3x+jXev2ctgAEDEBuLt2/p7SU6GtraUv6/p6np2rPn7fj4M2rytO29HC6rL9f8+ejdG10l2ThBGRKhnBf1cjiIjoa9Pdq1k9+v9iEhCjMueuoUnJ3F3YY4Ph4tW0JLi+aYaKZk9TIqKujQgc3vr5qa6NsXhw7R20twMHx8pLz3xAkMHMj0wn+RHBwU4IkwNhaHD0t8dp4yJEJ5fiIsFRWF7t1hYyOP1T2ZmYiNRf/+bMchnio1LiqkZIkQbE8Tgv7R0YICRERI+eWSx8PZs+jXj+qYZFa/PjQ08FT0qjzWlJTA3x/r1qFmTcluVPhEaG2N58+Rm8t2HGI4dQr9+sHBATdvsh3Kz8LC0L8/y0e9iKmwEBcvirvRNhS/UkZI+QpHWZ9t6t4daWl48oSu9o8dQ4cOMDKS5t7oaDRujPr1qY6JCqz/w1Vuwwbo60vz/UPhE6GaGtq0QWws23GIJzwcrq5wdMSVK2yH8gMFGhe9dAmWlhJ83SNPhPKpY0c8eMDmjjkqKhgyBGFhdLUv9fJBACdOwNWV0mioI8+rCYULBwMDweFIfK/CJ0LI8WrCch05gmHD0L07xQdLSi0lBS9eoFs3tuMQj0Tjoq9fo6gIJiZ0BsQICws8eybltvrySUsLrVqx/P3VwwP79tHS8uvXiI2VPpmdOiXB/+QMk+cnwkmTJFg4WAZJhCzYtw+jRqFPH8ydG3758mV2g9m/H8OGKcbJDAIBTp+WYC7z5k1leBwEoKEBExM8fsx2HJRi/SO1QwdwuYiLo77l/fvx229SlmglJqK4WH73+23VCllZtBfcSuHwYTx5IsHCwTKUIRHa2eHWLQgEbMchiaAgtGy5dvXqf7p1+/3kyZMsRnLgAEaMYLF/CcTEQF8fZmbiXq8cE4RCyreakPV6GQDDhtFSMiP18kEAx49j4EBKo6EUl4uOHeWuyiE3F9OnIzBQgoWDZShDIjQyQrVqSElhOw4J2dq+ANQEAv3s7Dy2YoiLQ0kJbGzY6l8yEo2LQlkmCIWsrJStXsbW9uvVq2H37rGZ3j09ERZG8UHB8fHIy0OnTlLeLs8ThEIWFs937Aj+9OkT24F8w+fz//gDffrAyUn6RpQhEUJxFlH8aNu2gIULratVW3LkyHC2YggJwYgR0swts0KiRFhcjPv30b49nQExSPnqZdauXVFQcLdbN9/37C2tNzVFvXqgdnZCuHxQut+p16/x9Ck6d6YyHsrt2eN2/vy7QYP82Q4EAKKjb9St23H7dsc5c2TaPVZ5EqGcL6v/laqq6pIlf54963z6NDuFM3w+Dh9WjHpRPp8fGHji3bvL4j+8PngAExPoid53XjEoXyKsXl1fXf1lcXGJuro6i2FQexjFu3cfwsI+jRwp5e0nT6JvX2n2JmWSnp4G8FxLSy5+tc6du5qV5auhYfThg0xLYZQnESpWvUypzp3Rrx+GDmVh6/CLF2FkJGWRFcMOHjw8deqhvLzVt2+L++CvTOOiAOrXR3ExMjPZjoM6CxdO/+uvqc2aRenr67MYxrBhOH4cBQUUNBUTE9OsWe/s7G5fviRK14L8j4sCuH8/ytnZY+DAQLYDAYDhw/3U1O77+1u2a9dOlnaUJBG2bYvHj5Gfz3YcUjlyBCUlGDWK6X4V6LgJfX29kpJ3mprZOjo6Yt6iTJUyQpaWSvVQyOFw/P1tX7yo9eYNm2EYGaFdO5w+TUFTz5+/yMuzVlFpkZ6eLsXteXm4cUOaM5sYpqur6+fX+eRJudgB7uLFOsOGBa5bN5/LlSmXKUki1NBAq1aIj2c7DqmoqSEkBPv3M7qs6utXnDiBIUOY61EWNWr0rldvw927YZZi15Ur2RMhlHF0VE0NLi5gtWgaADw8qBkdtbQcpK3dYd26Hr1795bi9nPnYGenGIP5ffrg6lW5ePA4cADDqSixUJJECMDOTlFHRwG4usLBgdHdBU+eRNu2MDZmrkdZBAVh0qSWJiZNxLz+wwdkZaF5c1qDYpqSPREKDRjA0NGAlRg0CJcv48MHWdvZuFFl5szREyZ4c6QqlTl+XAHGRYX09WFriwsXWA4jLQ2pqejenYKmlCcRduigwIkQwKlTyM3F1KkMdRcaqjDLB/PyEB4OT08Jbrl1C7a2kG2wRO4o3xMhgD59cOMGPn9mMwY9PfTogfBwmRr5+BGHD8Nf2lJK4UbbirLxPeTjG8z+/Rg8mJrNQJTno0KhnwgBVKuGHTuweTONGwGXys7G5ctwc6O9I0qEhsLZWbL9i+/cUbZxUQCWlvj3XxQXsx0HpXR0YGeHc+dYDsPDQ9aV9YGBcHNDnTpS3n79Oho1ktONtss1cCBOnUJJCZsxhIVRMy4KZUqEjRqBw0FaGttxyMDTE61bQ6r5BckcOYIePcBqsZ4EgoIwZoxkt9y+rWyVMgC0tNCgARPfkxgmD88WffsiMREZGVLeXlyMwEBMnix9AJJuFsE6Y2M0asTmFjMJCcjNlX7jgjKUJxFC8UdHAVy4gIwMLF5MYxczZiydPNlWS2srjX1QJyEB795JNg3A5yM2Fra2tMXEHqUcHXV1xenTLD/pqqryNTU9raw6RkZKM/F16BCaNYO1tfQBnDqlMBOEpdj9BnPgADw8KNsMRKkSoeKuJixlaIjVq7FsGV69oquL4OCjhYWno6P30tUBpbZvx5gxUJGkVDspCXXqSHwyp0JQvoMJAdSrBzMzREezGcP79+/z8l59+rRm+3Zpzq3fuFGm2f3ERBQVwcpK+hZYwWIiFAioHBeFkiXCrKyQ7dudgoKo2yiCDdOmwdSUxgFSU9OZdeuOWLt2AV0dUKegAAcPwttbsruUb+FEKaV8IgTbzxYA6tSpM2qUo47O6o4dJ0p6b3Q0/t/euQc1dW1/fIUkhEeCgAKaglNFFHl4+SlQUEliewdEK2WAInb4qUX0WlEcrcVqtWrrRWRaleFap0V0UKnAhauCKEpvCdhWEPGBhVIjIQKClLSCBBLy4vdHasoP5ZHkPEjYn79ykrPXWslJss7ee+3vfvbMoB9sYeG4FtoejnnzQK2Gn/UUDzCIH38ES0ssbx1MKhFmZx+RSgsPHDhCdiCG8t13UF8PaWnYW25sBIEg9pdfrkdFGUGBWl4eBAbC9Om6tTK9pfRaTG8PCg2aREjuBjJHjuz797+Lz5710VXj6dgx2LbNoBJlI1o4MQSy7mDOnwe9dexeiUklws2b3zc3/3tYGOEaLVjj7AxJSZCUBM+eYWz54EFITARbW4zN4sTJkxAfr3MrE+4Rvv46SCQYrHgbb3h6AoNBfo5fuhRYLCgo0KGJSAQVFfpvugQAHR3w6BEEBelvgURISYRKJRQUYCwGYlKJcNeuzXFxt+bO3Ux2IBiQnAxOThgvsRcKoajIoNo2Ivn1V2hshGXLdGv10UeHGxpira0f4RMUyVAo4OlpgtOEALBiBfm1owCwdy8cOKCD8G9aGsTHg7W1/h4vXYLQ0PEutD0cQUEgEkFrK6FOS0th5kyYNQtLmyaVCAHA39/4tqEYjpISqKzEUho/ORkSEsDODjODuPLNN/D++7r9QYhEooyMMqUyMjn5X7jFRTJ/+xuaJsSR5cvBygouXhzTyT09cPYsbNpkkEejWzgxGCoVQkOJFsnDSlZtMKaWCN94w3QSoYcHrFsH69djI43f3AyXLsHWrRiYIgC5HM6d01mI/LXXXmMy6Uxm8sqVOnYkjQeTLBwFgEWLoK0NRCKy4wDYswc+/3xME5YnT0JwMLi46O/LWIS2R4DgOxiZDIqLsRdJNrVEOHfunzqTpsE33wCLBRERGJj65z/hH/8Ae3sMTBHAhQswb57Oox90Ot3NraigoHrZsmB84iIfUy0cNTODZcvIF+AGgLAwoNGgsHCU01QqOH7c0DvLa9cgIMA4hLaHIyQEbt4EwvarLywEX1/9FXyGw9QSIYUC8+cTuo0D3hQVwfXrUFJikJGWFsjPN5ruIOilJgMA/f1QUwOBgTgENG7w9oa6OlCpyI4DB8bJ6CgA7N4NBw6M0im8dAmcnAytTy4sNOJxUQ3W1hAURJxIHh7jomB6iRAA3ngDbt0iOwjs8PeHiAhYudIgWb9Dh2DDBnBwwC4sPGlqgtpaff4gbt2CuXOBxcIhpnEDiwVOTtDYSHYcOBAcDLdvY18prQfh4aBWj7JJoWbVhCGoVHD1qjEJbQ8HYXcwXV14iSSbYCI0pXoZDTk50N+fZWExPz7+Qz2at7ZCbq6hP1oiyciA1auBwdC5YUUFcLk4BDTOMNXRUUtL4PHgyhWy4wCgUGDPHti/f9hO4Z070Nxs6Cp4oxPaHo533oGSEpDLcXdUUICXSLIJJkJNvQy5i3OxhUoFhWK/SvXdqVP63HelpMC6deDoiHlcuKBUQlYWxMXp07aiAjgcrAMaf5hqIoTxNDoaGQkKxbBTEkeOwJYthm4AVFhorOvoh+DoCO7uUF6OuyOcxkXBJBPh1KnAZJra2BGH429m5kehrHJ0hJoaHRq2t0NODuzYgVtkWHP5Mri6goeHzg2VSqiqwkyNfjxjwokwLAxKS0EmIzsOAAoFPvkEDhx4xUttbXD1KqxbZ6iLoiKjnyDUQsAdTHs73LmDl/YklolQrVardZUnwgcSpwnVanVvby/mZsvKclWqRonkc19f8PcfZZdahUIhfbHkIiUF1q7Fsjs42Dge6FcmAwA1NTBjBmllsXK5XEbU/zdZiVAmk8lxHv+aPBm8vaGsTM/mEolkALuxoKgo6O2F69eHGj9+HGJjDZVnqqsjTmi7B/+NjyMi4OJFPcfhlEplX1/fqKfl5kJ4OFha4mJ8rInw22+/pVAoX3zxxQjnNDY2Pno0LhQ9SFxNWFtbGxyMV+2+pSVcuQJ5efCf/8AIXcPz588nJiYCwNOncPYsbN+OZQxnz5798EN9dc7J8QAACzRJREFUpirHQmsrVFVBZKQ+bcmdIPzqq68+/fRTYny5uoJYDM+fE+PtL/bu3XvixAm8vRjSt5g/f357eztWkZiZwe7df+2J5u3tLRaLpVI4eRISEw01TuQ6ejc3N4lEgrMLYDLhzh192l67du1/R761BwB9x0WLi4vHYnxMiVAsFqekpAQYj4Cj6dXLDCYyEsRi0HQNR55LS02FNWuAzSYqMoPJzIRVq8DKSp+2FRXGKtioK2Zm4OFhmsvqASA8HIqKxsscf3Q0PHsG//3vX8+cOQOBgeDqaqhlE1g4MYTwcBxHRxsbQSSCJUvwsj+mRJiQkLBnzx57Y1mMDbBgAfz8M/T3kx0Hbmi7hrm5MG3aq9WKxWI4c8aYZgfVajh9Ws+pF7UafvxxoiRCMOlpQldXsLMbL0uBqVT45BPQdvUHBiAtDYMC7I4OEAhM7euK6zRhdjbExBhanTQCoxsuKiqSSCTR0dFZWaNv5Xrz5s0hzzg6OtqSsdmBs/Ps/PxWL6/RR4ex5eHDh1Kp9D4hQvqzZkFpqdnmzTPnz7detarzo4/aNM+3tLR0dXXt2NEZHEwRi59gq7PT2tr6xx9/4PEGf/iBxWJNpVAEethuaLC0s5ve3v4rdgNjutHe3i4Wi4m57gBgbz+lrMxi4UJC1Y7FYrG5uTkB7zEwcFpGBjAYOl9LhUJRX1/f2dmJYTCenpS2tjmZma1KpfL06XalkmVr+6uBn0FBweSAAOv6+maMYhwFtVr94MEDK/1GWsaMhQW0tXkUFz9ydtZtIrmpqamnp2fk71VW1pz9+1vu39f5/7ypqWkslSsUzfRvTk7OTz/9NOQ1Z2fnDRs2BAYGXr9+3cXFZfny5UuWLNkxfBdjxYoVpaWlFApl8JMODg6kJMK+vnnm5q002h8E+1WpVFKplMlkEun0+fO/K5X29vZ/bq6tUCgUCgVAgLl5C42G8fpkuVyuVCrx+FHJ5Wy12trCQqBHW6XSTi53sbIirZckl8tVKpWlrlP5+qJQOCiVTpaWhG6KKpVKqVSqubk53o4UCkelcrKl5S+6Nnz+/DmLxRryF2Q4MtlsKlUilTYwGHMGBqwtLAythOjvfx1ggMF4jEV0o9Pd3W1jY4P5x/IyfX3/w2AIqdRunVoplcr+/n7rkbbwMOvpWcRi/QCg84i5QqHYt29fTEzMyKf9mQjLy8vr6+uHvObg4FBdXX3v3r3IyEgASE9Pd3d3X79+PX7FIAgEAoFAEMyfQ6NcLpf7qpI7qVTa1dVVU1MDAN3d3a2trY0mtkAPgUAgEBMbytiX3Yw6NIpAIBAIhNGhw4L6xYsXz5kzB79QEAgEAoEgHh16hAgEAoFAmB7YSKw1/H+6u3WrGjINRCLRYHG1Z8+etbW1EeO6t7e3paVFe/j7778/efJEe9jZ2Sk2eAnFEBcSiUQoFCoUCgPNmgDt7e2Dv/yDPyWTQbMgh+wohuXp06eDw5NKpY8fY1CQ2dLSov1Fq1QqoVCoFRdUKBRCoXCc9yIUCkVDQwO5QarVaqFQqHyxh5zmY3z+kirSxx9/PFih8OWPd8glHqPf/kFryTs7O3/77bdhGwwYjGaVxqxZszxfcPHiRcPNGh3bt29fvny55rFMJvP09Dx16hQxrmtraxkMRl9fn+YwLCzM1tZWqVRqDrlcblpamoEuCgoK5syZo3n8+PHj2bNn79y5U61WG2jWBIiNjZ0yZYr2yx8XF0d2RNjD4/HS09PJjmJY3n777UOHDmkPS0tLnZycDDcbEhKSmpqqeVxRUQEAX3/9teawqKjIxcXFcBe4IhAIAEAmk5EYg6ZT1NzcPDAw0N/fHxUVxeVyu7u7h5xmb2/f1dWlPXz48CEASKVS7TNDLvFYiIqKSkhI0Dzu7Ox0cnIqKSkZ7mTMRLcvXbr08wveMTHtoLGRnJzc1NSkkR3Ys2cPm81eu3YtMa69vLxYLFZVVRUAqNXqW7duubu719bWAoBMJquqquLxeFj5EgqFS5Ysee+991JSUghYmWQUxMXFab/8mZmZZIeDwAYul1v+Ym8hPp/P4/G0h+Xl5Uvw0/syRfr6+sLCwqRS6dWrV21sbAjweOLEifz8/NLSUgDYtGlTRERESEjIcCeb4DZMZMFgMDIzM7dv356Xl3fq1KnMzEzC8gSFQuFyuXw+HwDu37/v5uYWEhKiOaysrLS2tvby8sLEUV1dHYfDSUxM3LdvHyYGEYhxC4/Hu3HjhkqlAoDy8vLdu3dr+oUAwOfzX7neDPFKuru7g4OD7ezsLly4QJjoxJQpU9LS0jZs2JCRkVFdXX348OERTkaJEEsCAgJWr14dExPz5Zdfuri4EOlae/eq+YlyOBztIY/HMzPD4EJ3dHRwOJzDhw9v3brVcGsIxDjHz89vYGDg7t27crn8wYMHPB7PyclJIBD09PTcu3cPw1EWkycsLMzd3T07O5tOpxPpd+XKlT4+Phs3bszMzGSxWCOciRIhlgwMDDQ0NNDp9AHCJ6h5PF5lZaVMJisvL+dyuYGBgVVVVWq1WnOIiQsGg2FnZ1ddXU38u0MgiIdGoy1cuJDP51dXV/v4+NDpdA6Hw+fzb9y4MW3atJkzZ5IdoNHAZrPr6+vx3grqZeRyeWNj41j+kFEixJLMzEyRSHTt2rUdO3a0thIqiOzl5WVjY1NZWVlZWRkQEGBpaenq6lpdXV1VVYXVZIatrS2fzy8uLt62bRvKhYjxg5WV1eDNV3t7e7HSwtXMC/L5fA6HAwCagZby8vI333wTE/sThDNnzjg5OS1duvTletFXorl8g4vwJRLJiGKkr+azzz5zcHA4ffr0unXrRt6dGCVCzHjy5MmuXbuysrJ4PN6aNWs2btxIpHcKhRIUFJSenj5jxgyN5DeHw0lNTWUymZ6enlh5cXZ2Lisru3z58nZsd/tFIAxg9uzZdwbtCVtTU+Pu7o6JZc004ffff68ZVtH0CNEEoa7Q6fS8vLyx50I2m81isbTXtL+/v66uTtdrevfu3ePHj2dkZKxatcrX13fXrl0jnIzb/k4TjIGBgfj4+Pj4eH9/fwA4ePCgj4/PuXPnYmNjCYuBx+MlJiYmJSVpDjkcTkpKSkREBLY1O5pcqOllHj16FEPLiPHMlStXtKtRfXx8wsPDyY1nMAkJCfPmzduyZcvixYsfPnyYnp5eWFiIiWVfX1+1Wl1ZWenn5wcA9vb29vb21dXVOTk5mNgngIMHD1KpVM3juLi46dOnkxIGnU7Pzc2Njo4ODQ0tKSkZecaOQqHs3bt3/fr1O3futLGxyc7OdnNze+utt8buTi6Xr1mzJiUlRTOCfeLECW9v73fffXe4Oxjq/v37dXk7rw6aRqNxuVy897saz7S3t0skkqSkJBqNBgB0On3x4sXNzc0LFiwgLAY2m+3o6BgdHe3g4AAAU6dOtbGxiYmJwaRsh0qlstlszduZNGlSeHi4SCRis9l2dnaGGzdqaDSah4fHjBkzyA4ER2g0GovFMnvBtGnTxpXaIpPJXL16tUAguH37tqWl5dGjR319fTGxbGZm5uzsvHTpUq1BZ2dnPz+/0NBQTOzjCoVCYTKZNBpNe+F8fHwmTZpEcBhMJjMoKIjBYFCp1IiICKlUqlQqh/xeUlNTt27damFhoX1m0aJFc+fOvX37tkAg4HK5x44d06nQRiAQ2NrafvDBB5pugLW1tZ+fX0dHx3DDY0hiDYFAIBBkMnnyZKFQSHyS1oLmCBEIBAIxoUGJEIFAIBBkkp+fr0dRKIagoVEEAoFATGhQjxCBQCAQExqUCBEIBAIxoUGJEIFAIBATGpQIEQgEAjGh+T/b/zKzTXbs5wAAAABJRU5ErkJggg==", "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "metadata": {}, "execution_count": 7 } ], "cell_type": "code", "source": [ "plot_bandstructure(scfres, kline_density=5)" ], "metadata": {}, "execution_count": 7 }, { "cell_type": "markdown", "source": [ "or get the cartesian forces (in Hartree / Bohr)" ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "β”Œ Warning: Forces for shifted k-Grids not tested\n", "β”” @ DFTK /home/runner/work/DFTK.jl/DFTK.jl/src/terms/terms.jl:72\n" ] }, { "output_type": "execute_result", "data": { "text/plain": "2-element Array{StaticArrays.SArray{Tuple{3},Float64,1,3},1}:\n [0.00044924555703957123, 0.00044924534195866773, 0.0004492410964905146]\n [-0.00044924534195892594, -0.00044924555703924543, -0.0004492410964902822]" }, "metadata": {}, "execution_count": 8 } ], "cell_type": "code", "source": [ "compute_forces_cart(scfres)[1] # Select silicon forces" ], "metadata": {}, "execution_count": 8 }, { "cell_type": "markdown", "source": [ "The `[1]` extracts the forces for the first kind of atoms,\n", "i.e. `Si` (silicon) in the setup of the `atoms` list of step 1 above.\n", "As expected, they are almost zero in this highly symmetric configuration." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Where to go from here\n", "Take a look at the\n", "[example index](https://juliamolsim.github.io/DFTK.jl/dev/#example-index-1)\n", "to continue exploring DFTK." ], "metadata": {} } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.5.4" }, "kernelspec": { "name": "julia-1.5", "display_name": "Julia 1.5.4", "language": "julia" } }, "nbformat": 4 }